Can someone help me figure out how to take the expected value of the quantity $[(x-\mu_x)^2 (y-\mu_y)^2]$?
I know that $E[(x-\mu_x)^2]$ is of course $\sigma^2_x$, $E[(y-\mu_y)^2]$ is of course $\sigma^2_y$, and $E[(x-\mu_x)(y-\mu_y)]$ is $\sigma^2_{xy}$, but it's tripping me up to have four terms inside the brackets instead of two, and I'm not sure how to deal with it.
I tried expanding it to a sum of terms and then taking the expected value of each term, but I think I must be doing it wrong. Expanding yields:
$E[(x^2 -2x\mu_x + \mu^2_x)(y^2 -2y\mu_y + \mu^2_y)]$
$E[x^2 y^2] - 2\mu_y E[x^2 y] + \mu^2_y E[x^2] - 2\mu_x E[xy^2] + 4\mu_x \mu_y E[xy] - 2\mu_x \mu^2_y E[x] + \mu^2_x E[y^2] - 2\mu^2_x \mu_y E[y] + \mu^2_x \mu^2_y E[1]$
Then, because x and y are independent, this simplifies to:
$E[x^2]E[y^2] - 2\mu^2_y E[x^2] + \mu^2_y E[x^2] - 2\mu^2_x E[y^2] + 4\mu^2_x \mu^2_y - 2\mu^2_x \mu^2_y + \mu^2_x E[y^2] - 2\mu^2_x \mu^2_y + \mu^2_x \mu^2_y$
And then to:
$E[x^2]E[y^2] - \mu^2_y E[x^2] - \mu^2_x E[y^2] + \mu^2_x \mu^2_y$
$(E[x^2] - \mu^2_x)(E[y^2] - \mu^2_y)$
Which equals $\sigma^2_x \sigma^2_y$.
However, I feel like I must be doing something wrong. $x$ and $y$ are assumed to be jointly Gaussian, so what about $\sigma^2_{xy}$? Surely that must show up somewhere? Can anyone point out my mistake?
Thanks!