Determine all $P(X)\in K[X]$ such that $P\big(X^2+1\big)=\big(P(X)\big)^2+1$, for fields $K$ of any characteristic. This question is inspired by this thread.  However, in this question, I take an arbitrary field instead of $\mathbb{R}$ and drop the assumption that $P(0)$ must be $0$.

Let $K$ be a field.  Determine all $P(X)\in K[X]$ such that $$P\big(X^2+1\big)=\big(P(X)\big)^2+1\,.$$
  In other words, if $Q(X):=X^2+1$, we are to determine all $P(X)\in K[X]$ that commute with $Q(X)$ with respect to polynomial composition:
  $$P\big(Q(X)\big)=Q\big(P(X)\big)\,.$$

First suppose that the characteristic of $K$ is not equal to $2$.  I have observed that some solutions are given by 
(i) $P(X)=\frac{1+\sqrt{-3}}{2}$ and $P(X)=\frac{1-\sqrt{-3}}{2}$ if $\sqrt{-3}\in K$, and
(ii) $P(X)=Q_n(X)$ for $n=0,1,2,3,\ldots$,
where $Q_0(X):=X$ and $Q_n(X):=Q\big(Q_{n-1}(X)\big)$ for all $n=1,2,3,\ldots$.   Are these all of the solutions?
If the characteristic of $K$ is equal to $2$, then the given functional equation is equivalent to
$$P\left((X+1)^2\right)=\left(P(X)+1\right)^2\,.$$
If $K=\mathbb{F}_2$, then it follows that
$$P(X+1)=P(X)+1\,,$$
whence $P(X)$ must take the form
$$P(X)=X^{2^{r_1}}+X^{2^{r_2}}+\ldots+X^{2^{r_k}}\text{ or }P(X)=1+X^{2^{r_1}}+X^{2^{r_2}}+\ldots+X^{2^{r_k}}$$
where $k$ is an odd positive integer and $r_1,r_2,\ldots,r_k$ are pairwise distinct nonnegative integers.
Can anybody prove or disprove whether my list is complete when the characteristic of $K$ is not equal to $2$?  Bill Dubuque's answer in this thread shows that, if $K$ is a subfield of $\mathbb{C}$, then the list is indeed complete.  This should imply that, if $K$ is of characteristic $0$ and the cardinality of $K$ is that of the continuum, then the list is complete (as $K$ can be embedded into $\mathbb{C}$).
Can someone try to make a characterization of the solutions when $K$ is an arbitrary field of characteristic $2$?  This case seems to have weird polynomial solutions.  For example, when $K=\mathbb{F}_4=\mathbb{F}_2[\alpha]$ with $\alpha^2+\alpha+1=0$, we have the following solutions $P(X)=\alpha$, $P(X)=X$, $P(X)=X+1$, $P(X)=X^2$, $P(X)=X^2+1$, $P(X)=X^2+X+\alpha$, and $P(X)=X^3+\alpha\,X^2+\alpha\,X$.  
UPDATES: 
(1) By balancing the coefficients, it can be shown that, in characteristic not equal to $2$, there is at most one solution of degree $d>0$, and this solution lies within $\mathbb{F}_p[X]$, where $\mathbb{F}_p$ is the prime field of $K$ (with $\mathbb{F}_0:=\mathbb{Q}$).  Hence, it suffices to show that the list is complete when $K$ is a prime field.  Ergo, for characteristic $0$, we are done.
(2)  In characteristic $2$, via balancing the coefficients, we see that all equations involved are quadratic in the coefficients, whence the solutions must lie within $\mathbb{F}_4[X]$.  Consequently, we may take $K$ to be $\mathbb{F}_2$ or $\mathbb{F}_4$.
(3) In odd characteristic $p$, the list is indeed incomplete.  Note that $P(X)=X^{p^k}$ is a solution for all $k=0,1,2,\ldots$.  Let $R_p(X):=X^p$.  Is it true that all solutions must take the form $$P(X)=\left(Q^{\circ k}\circ R_p^{\circ l}\right)(X)$$
for some $k,l\in\mathbb{Z}_{\geq 0}$?  Here, $\left(F^{\circ 0}\right)(X):=X$ and $$\left(F^{\circ l}\right)(X):=\left(\underbrace{F\circ F\circ \ldots \circ F}_{F\text{ occurs }l\text{ times}}\right)(X)$$ for each $F(X)\in K[X]$ and $l\in\mathbb{Z}_{> 0}$.
(4) It turns out that the characterization in (3) is not complete either, at least not in characteristic $3$.  The polynomials 
$$X^{5}+X^3-X\,,\;\;X^7-X^5-X^3-X\,,\text{ and }X^{11}+X^9-X^7+X^5+X^3+X$$ 
are solutions in characteristic $3$ which do not arise from the formula in (3).    The case of characteristic $3$ may be the only special case where (3) does not give a full characterization.
(5) Combined with Hurkyl's answer, we know all solutions in characteristics $0$, $2$, and $3$.  So far, I have not been able to find any solutions outside (3) in characteristics greater than $3$.
 A: Odd Characteristic
(this section incomplete)
Using the proof of dinoboy, $P$ is either an odd or even function1. If $P$ is even and nonconstant, it must be of the form $P = S \circ Q$, where $S$ and $Q$ commute. Thus, the problem reduces to finding the set of all odd possibilities for $P$. (dinoboy's proof also shows that for all fields of characteristic zero, $P(x) = x$ is the only odd solution)
Note the possible finiteness of $K$ doesn't matter here, since we can pass to any infinite extension field $L/K$
Characteristic 3
Let $T_n(X)$ be the Chebyshev polynomial of the first kind. These polynomials satisfy $T_m \circ T_n = T_{mn}$.
By noting that $Q(X) = -T_2(X)$, this corresponds to one of the special cases of Ritt's theorem (which applies to characteristic zero). Consequently, we get a family of odd solutions $P = T_m$ for odd $m$.
Since there are Chebyshev polynomials in each odd degree, the note in the OP proves this is the complete solution.
The complete solution including even degrees can be more simply stated as $P = (-1)^{n+1} T_n$ where $n$ can be any integer.
Characteristic 2
Let $\alpha = \alpha^2 + 1$ and assume $\alpha \in K$. let $K_s$ be the field isomorphic to $K$, and we fix a bijection $K \to K_s : c \mapsto c + \alpha$.
If $F$ is a polynomial defined over $K$, then $F_s$ is a polynomial defined over $K_s$ given by $F_s(X) = F(X + \alpha) + \alpha$. In particular, $Q_s(X) = X^2$.
The polynomials $P_s$ over $K_s$ that commute with $Q_s$ are precisely the polynomials $P_s$ whose coefficients are in $\mathbf{F}_2$. Consequently, the polynomials over $K$ that commute with $Q$ are precisely the polynomials of the form $P(X) = F(X + \alpha) + \alpha$, where $F \in \mathbf{F}_2[X]$.

