Ideas for calculation of an expectation I have an $\mathbb R^d$-valued random variable $X$ that is distributed with density $f$. I don't know anything about $f$ except that it is $\mathcal C^2$ and its second derivative is in $L^2$.
Now consider the function $$g(x) = \frac{1}{2}\Sigma^{-0.5}(x - \mu)(x - \mu)'\Sigma^{-0.5}\times
\\\sqrt{2\pi}^{-d}\det(\Sigma^{-1})\exp\left(-\frac{1}{2}(x - \mu)'\Sigma^{-1}(x - \mu)\right),$$ so $g$ is basically the derivative of a multivariate normal distribution with respect to $\Sigma$. Since $g$ is measurable, $g(X)$ is a random variable and I am now interested in $$\mathbb E[g(X)] = \int_{\mathbb R^d}g(x)f(x)\:\mathrm dx$$
My attempt so far was to use a the substitution: $x = \Sigma^{0.5}y + \mu$. Then $$\mathbb E[g(X)] = \int_{\mathbb R^{d}}\frac{1}{2} yy'\exp\left(-\frac{1}{2}y'y\right)f(\mu + \Sigma^{0.5}y)\:\mathrm dy $$
But now I got stuck. Since I don't know $f$, all I can do is to use a second order Taylor approximation: $$f(\mu + \Sigma^{0.5}y) = f(\mu) + y'\Sigma^{0.5}\nabla f(\mu) + y'\Sigma^{0.5}\nabla^2f(\mu)\Sigma^{0.5}y + o(\Vert\Sigma^{0.5}y\Vert^2).$$
If I chose this approach, I have to calculate the following terms:
$$A = \int\frac{1}{2}f(\mu)yy'\exp\left(-\frac{1}{2}{y'y}\right)\:\mathrm dy = \frac{1}{2}f(\mu)\sqrt{2\pi}^d I_{d\times d}$$
This is correct, I hope, because it was the simplest part.
The second term gives me more trouble:
$$ B = \frac{1}{2}\left(\int yy'y\exp\left(-\frac{1}{2}y'y\right) \:\mathrm dy\right)\Sigma^{0.5}\nabla f(\mu)$$
I dont know what to do here. Same goes for the third part:
$$C = \frac{1}{2}\int yy'y'\Sigma^{0.5}\nabla^2f(\mu)y\exp\left(-\frac{1}{2}y'y\right)\:\mathrm dy $$
and the little oh part:
$$D = o\left(\int yy'y'\Sigma y\exp\left(-\frac{1}{2}y'y\right)\:\mathrm dy\right)$$
I have problems with integrating things like $yy'y$, $yy'Cy$. On the one hand, its like computing third moments for a normal distribution but because of the vectors it isn't. I hope that someone has an idea how to process the terms in $B$, $C$ and $D$. 
If nothing helps, I would be also be happy with assuming that $\Sigma$ is a diagonal matrix and then derive upper bounds for the diagonal elements of $\mathbb E[g(X)]$.
 A: We will use the following facts:
\begin{align}
\int_{\mathbb{R}} \exp(-\tfrac{1}{2}x^2)\mathrm{d} x &= \sqrt{2\pi}, \\
\int_{\mathbb{R}} x^2 \exp(-\tfrac{1}{2}x^2)\mathrm{d} x &= \sqrt{2\pi}, \\
\int_{\mathbb{R}} x^4\exp(-\tfrac{1}{2}x^2)\mathrm{d} x &= 3\sqrt{2\pi}, \\
\int_{\mathbb{R}} x^{2k+1} \exp(-\tfrac{1}{2}x^2)\mathrm{d} x &= 0
\end{align}
where $k\ge 0$ is an integer.
It is easy to get
$$I_1 = \int_{\mathbb{R}^d} yy^{\mathsf{T}}\exp(-\tfrac{1}{2}y^{\mathsf{T}}y)\mathrm{d} y
= \sqrt{2\pi}^d I_d$$
and
$$I_2 = \int_{\mathbb{R}^d} y(y^{\mathsf{T}}y)\exp(-\tfrac{1}{2}y^{\mathsf{T}}y)\mathrm{d} y = 0.$$
Let us deal with the third part. Denote $Q = \Sigma^{1/2}\nabla^2 f(\mu)\Sigma^{1/2}$.
Let $Q = U\mathrm{diag}(\lambda_1, \lambda_2, \cdots, \lambda_d)U^{\mathsf{T}}$ be
the eigendecomposition of $Q$ where $\lambda_1, \lambda_2, \cdots, \lambda_d$
are the eigenvalues of $Q$, and $U$ is an orthogonal matrix
whose columns are the eigenvectors of $Q$. We have
\begin{align}
I_3 &= \int_{\mathbb{R}^d} yy^{\mathsf{T}} y^{\mathsf{T}}Q y \exp(-\tfrac{1}{2}y^{\mathsf{T}}y)\mathrm{d} y\\
&= \int_{\mathbb{R}^d} yy^{\mathsf{T}} y^{\mathsf{T}}
U\mathrm{diag}(\lambda_1, \lambda_2, \cdots, \lambda_d)U^{\mathsf{T}} y \exp(-\tfrac{1}{2}y^{\mathsf{T}}y)\mathrm{d} y \\
&= U \left(\int_{\mathbb{R}^d} zz^{\mathsf{T}} z^{\mathsf{T}}\mathrm{diag}(\lambda_1, \lambda_2, \cdots, \lambda_d)z \exp(-\tfrac{1}{2}z^{\mathsf{T}}z)\mathrm{d} z\right) U^{\mathsf{T}}\\
&= U \left(\int_{\mathbb{R}^d} zz^{\mathsf{T}} (\sum\nolimits_k \lambda_k z_k^2) \exp(-\tfrac{1}{2}z^{\mathsf{T}}z)\mathrm{d} z\right) U^{\mathsf{T}}\\
&= \sqrt{2\pi}^d U \left(2\mathrm{diag}(\lambda_1, \lambda_2, \cdots, \lambda_d) + (\sum\nolimits_k \lambda_k)I_d \right) U^{\mathsf{T}}\\
&= 2 \sqrt{2\pi}^d Q + \mathrm{Tr}(Q)\sqrt{2\pi}^d I_d
\end{align}
where we have used the substitution $z = U^{\mathsf{T}} y$ (about the change of variables, see [1], and [2] Ch. 2.1). 
Reference
[1] J. Schwartz, "The Formula for Change in Variables in a Multiple Integral",
The American Mathematical Monthly, Vol. 61, No. 2 (Feb., 1954), pp. 81-85
[2] Robb J. Muirhead, "Aspects of multivariate statistical theory", 2005. 
