What's the value of $\text{cov}(x, x^TAx)$, when $x$ follows a normal distribution When $x\sim N_k(\mu,\Sigma)$ is a multivariate normal distribution, $A$ is a symmetric matrix, how can I show that
$$\text{cov}(x, x^TAx) = 2\Sigma A\mu$$
 A: Writing $$z=x-\mu,$$ we see that $z \sim \mathcal{N}(0,\Sigma).$  Using the bilinearity of the covariance operator repeatedly, make the substitution $x=z+\mu$ and (mindlessly) compute
$$\eqalign{
\operatorname{Cov}(x, x^\prime A x) &= \operatorname{Cov}(z+\mu,\ (z+\mu)^\prime A (z+\mu))\\
&= \operatorname{Cov}(z,\ z^\prime A z + \mu^\prime A z + z^\prime A \mu + \mu^\prime A \mu) \\
&= \operatorname{Cov}(z, z^\prime A z) + \operatorname{Cov}(z, \mu^\prime A z) + \operatorname{Cov}(z, z^\prime A \mu) + 0.
}$$
The first of those three terms must evaluate to $0$ because it is the expectation of a homogeneous odd-order polynomial in the components of $z.$ The symmetry of the distribution ($-z$ also has a $\mathcal{N}(0,\Sigma)$ distribution) shows this expectation equals its negative, whence--since it is finite--it can only be zero. 
For the second term use the identity (often taken as the definition)
$$\operatorname{Cov}(u, v) = E(uv^\prime) - E(u)E(v)^\prime$$
for vector-valued random variables $u$ and $v$, whence (since $E(z)=0$) we obtain
$$\operatorname{Cov}(z,\ \mu^\prime A z) = E(z\ (\mu^\prime A z)^\prime) = E(z\ z^\prime A^\prime \mu) = E(zz^\prime)A^\prime\mu = \Sigma A^\prime \mu.$$
For the third, note that $z^\prime A\mu$ is a $1\times 1$ matrix and therefore equals its own transpose, whence
$$\operatorname{Cov}(z, z^\prime A\mu) = E(z\ z^\prime A\mu) = E(zz^\prime)A\mu = \Sigma A \mu.$$
Finally, the symmetry of $A$ means $A=A^\prime,$ entailing
$$\operatorname{Cov}(x, x^\prime A x) = 0 + \Sigma A^\prime \mu + \Sigma A \mu + 0 =  \Sigma(A^\prime+A)\mu = 2\Sigma A \mu.$$
