Let us consider the following integral: $$\int_{\mathbb{R}^d} |D^{\alpha} P(x;\Sigma)| dx,$$ where $P(x;\Sigma)=((2\pi)^d |\Sigma|)^{-1/2}\exp(-\frac{1}{2}x^T\Sigma^{-1}x)$ denotes the density of a multivariate normal vector with $\vec{0}$ mean and covariance matrix $\Sigma$. Here, $\alpha=(\alpha_1,\dots,\alpha_d)$ is a multiindex and $D^\alpha=\frac{\partial^{\alpha_1}}{\partial x_1^{\alpha_1}}\cdots \frac{\partial^{\alpha_d}}{\partial x_d^{\alpha_d}}$ the derivative operator.
Consider the Cholesky decomposition of $\Sigma = LL^T$. Then if I am correct, the change of variables $L^{-1}x=y$ (or $x=Ly$) will transform the density $P$ in a standard multivariate normal density. We know that the Jacobian is $|dx/dy|=|L|=|\Sigma|^{1/2}$. Altogether,
$$\int_{\mathbb{R}^d} |D^{\alpha} P(x;\Sigma)| dx = \int_{\mathbb{R}^d} |D^{\alpha} P(y;I)| dy.$$ The last integrand is nothing else than: $\prod_{j=1}^d|H_{\alpha_j}(y_j)p(y_j;1)|$, where $p$ are the one-dimensional Gaussian densities with mean 0 and variance 1.
Finally, $$\int_{\mathbb{R}^d} |D^{\alpha} P(x;\Sigma)| dx = \prod_{j=1}^d E\left[|H_{\alpha_j}(Z)|\right],$$ where $Z\sim N(0,1)$.
Is this reasoning correct? Is there any "fixable" flaw in the computations? The whole thing looks too easy to me... I don't think it is correct because the LHS depends on $\Sigma$ while the RHS doesn't. Where is the mistake and in such case can one work out the expressions to get something similar?
I appreciate any help in this matter. Thanks a lot guys! :)