why is it equal to $E[S^2(t_1)]E[S^2(t_2)]+E^2[S(t_1)S(t_2)]+E[N^2(t_1)]E[N^2(t_2)]+E^2[N(t_1)N(t_2)]$ $Y(t)=a[S(t)+N(t)]^2$,$S(t)$ and $N(t)$ are both Gaussian random process and WSS with zero mean,and  $S(t)$ is independent of $N(t)$
\begin{align}
R_Y(t_1,t_2) & =E[Y(t_1)Y^*(t_2)] \\
& =a^2E[(S(t_1)+N(t_1))^2(S(t_2)+N(t_2))^2]\\
&=a^2(E[S^2(t_1)S^2(t_2)]+E[N^2(t_1)N^2(t_2)]+E[S^2(t_1)N^2(t_2)]+E[S^2(t_2)N^2(t_1)])\\
&=a^2(E[S^2(t_1)]E[S^2(t_2)]+E^2[S(t_1)S(t_2)]+E[N^2(t_1)]E[N^2(t_2)]+E^2[N(t_1)N(t_2)])
\end{align}
Can anyone tell me why is \begin{align} E[S^2(t_1)S^2(t_2)]+E[N^2(t_1)N^2(t_2)]+E[S^2(t_1)N^2(t_2)]+E[S^2(t_2)N^2(t_1)] = E[S^2(t_1)]E[S^2(t_2)]+E^2[S(t_1)S(t_2)]+E[N^2(t_1)]E[N^2(t_2)]+E^2[N(t_1)N(t_2)]\end{align},it seems a little wierd

Seems to be false. Take $S(t)=X,N(t)=Y$ for all $t$ where $\{X,Y\}$ is i.i.d with standard normal distribution. Then the identity you have written does not hold. 
