Correlation Function of Squared Standard Wiener Process I'm trying to find the mean and correlation functions of $ X_t:=W_t^2 $, where $W_t$ is a standard Wiener process ($\sigma=1,E\{W_t^2\}=t$ for $t\geq 0$).
Unless I'm missing something large here, the mean is very easy to compute.
$m_x(t) = E[X_t] = E[W_t^2] = \sigma^2 t = t$, for $t\geq 0$
Where I run into trouble is the correlation function. I've omitted some intermediate steps for time's sake that I could include if asked.
$$R_x(t_1,t_2)=E[X_{t_1}X_{t_2}]=E[W_{t_1}^2 W_{t_2}^2]=E[W_{t_1}W_{t_1}W_{t_2}W_{t_2}] = E[W_{t_1}^4] + (\sigma^2t_1)(\sigma^2(t_2-t_1))=E[W_{t_1}^4]+t_1(t_2-t_1) $$
I run into trouble calculating $E[W_{t_1}^4]$. I saw a solution somewhere online that led me to believe it should be $3\sigma^4t_1^2$, but I have no idea how to arrive at this answer. 
I've been trying to use the method for solving $E[V_t^2]$, where $V(t)$ is integrated white noise, $V(t)=\int_0^tX_\tau d\tau$ to give an idea...
$$E[V_t^2] = E[(\int_0^tX_\tau d\tau)(\int_0^tX_\theta d\theta)]=\int_{0}^{t} \int_{0}^{t} R_{X}(\tau-\theta) d \tau d \theta=\int_{0}^{t} \int_{0}^{t} \sigma^{2} \delta(\tau-\theta) d \tau d \theta=\sigma^{2} t$$
... but I seriously doubt it's as easy as squaring the result.
A nudge in the right direction for solving $E[W_{t_1}^4]$ would be much appreciated (assuming it's correct in the first place)
 A: Ito's Lemma: If $dX_t=\mu(t,X_t)dt+\sigma(t,X_t)dW_t$ and $f(t,x)$ is twice continuously differentiable then
$$
df(t,X_t)=\left[f_t(t,X_t)+\mu(t,X_t)f_x(t,X_t)+\frac{1}{2}\sigma^2(t,X_t)f_{xx}(t,X_t)\right]dt+\Big[\sigma(t,X_t) f_{x}(t,X_t)\Big]dW_t
$$
So, in your case take $dX_t=0dt+1dW_t=dW_t$ and $f(t,x)=x^4$. Then you have
$\mu(t,x)=0$, $\sigma(t,x)=1$, $f_t=0$, $f_x=4x^3$ and $f_{xx}=12x^2$. Hence
\begin{eqnarray*}
dW_t^4&=&df(t,W_t)\\
&=&\left[0+0f_x(t,W_t)+\frac{1}{2} 1^2f_{xx}(t,W_t)\right]dt+\Big[1f_{x}(t,W_t)\Big]dW_t\\
&=&
\Big[6W_t^2\Big]dt+\Big[4W_t^3\Big]dW_t 
\end{eqnarray*}
So, we have
$$
W_t^4 = 6\int_0^tW_s^2ds + 4\int_0^tW_s^3dW_s
$$
which means that 
\begin{eqnarray*}
E[W_t^4] 
&=& 
6E\left[\int_0^tW_s^2ds\right] + 4E\left[\int_0^tW_s^3dW_s\right]\\
&=& 
6\int_0^tE\left[W_s^2\right]ds + 0\\
&=& 
6\int_0^t s ds\\
&=& 
3t^2
\end{eqnarray*}
Hence
\begin{eqnarray*}
E\Big[W^2_{t+h}W^2_t\Big] 
&=&
E\Big[E\Big[W^2_{t+h}W^2_t\Big|{\cal F}_t\Big]\Big]\\
&=&
E\Big[E\Big[W^2_{t+h}\Big|{\cal F}_t\Big]W^2_t\Big]\\
&=&
E\Big[E\Big[(W_{t+h}-W_t+W_t)^2\Big|{\cal F}_t\Big]W^2_t\Big]\\
&=&
E\Big[E\Big[(W_{t+h}-W_t)^2+2(W_{t+h}-W_t)W_t+W_t^2\Big|{\cal F}_t\Big]W_t^2\Big]\\
&=&
E\Big[\Big(E\Big[(W_{t+h}-W_t)^2\Big|{\cal F}_t\Big]+2E\Big[(W_{t+h}-W_t)W_t\Big|{\cal F}_t\Big]+E\Big[W_t^2\Big|{\cal F}_t\Big]\Big)W_t^2\Big]\\
&=&
E\Big[\Big(h+2E\Big[W_{t+h}-W_t\Big|{\cal F}_t\Big]W_t+W_t^2\Big)W_t^2\Big]\\
&=&
E\Big[(h+0+W_t^2)W_t^2\Big]\\
&=&
hE\Big[W^2_t\Big]+E\Big[W_t^4\Big]\\
&=&
ht+3t^2
\end{eqnarray*}
where $\cal{F}$ is the filtration generated by the Brownian motion.
