I need to numerically evaluate/approximate non-central high-order moments of high-dimensional Gaussian measures/distributions with given mathematical expectations and covariance matrices. The Gaussian measures have dimension $d$ with say $d>1000$ and the moments have order $k$ with say $k>>10000$.

More precisely, I "only" need to evaluate the ratios of successive moments

$\frac{{\mathbb{E}X_1^{{k_1}}...X_i^{{k_i} + 1}...X_d^{{k_d}}}}{{\mathbb{E}X_1^{{k_1}}...X_i^{{k_i}}...X_d^{{k_d}}}}\quad k = \sum\limits_{i = 1}^d {{k_i}} $

The classical Wick-Isserlis theorem/formula and similar symbolic formulae look helpless because they involve an astronomically large number of terms $k!!$.

As suggested, I tried a Laplace's approximation/expansion but that does not work well: unless I'm mistaken, in standard form, it gives a wrong approximation that does not depend on ${{k_i}}$ and in fully exponential form, we fall on a multivariate system of quadratic equations of dimension $d$ for the maximum of the function...

I'm unable to find a single reference on this problem, that's quite unexpected.

Of course, Monte-Carlo methods is always an option but it is the last one: fast and deterministic algorithms are highly preferred.

So, please, what are my best options in order to numerically evaluate those moments in general? Computational complexity theoretic results welcome too.


PS.: I finally found a nice, ad hoc, deterministic method allowing to evaluate the ratios of successive moments instantaneously. Unfortunately I can't disclose it right now. I really don't understand why I can't find a single paper and why I don't get any answer about such a basic and fundamental problem??????

  • $\begingroup$ Have you tried using the moment generating function for a multivariate Gaussian? $\endgroup$ Jan 11, 2018 at 13:19
  • 1
    $\begingroup$ @NawafBou-Rabee The moment generating function of a Gaussian is a Gaussian. Which leads to Isserlis formula. However one can probably obtain a recursive version (integrating by parts). $\endgroup$
    – lcv
    Jan 11, 2018 at 13:58
  • 1
    $\begingroup$ @NawafBou-Rabee One nice reference is arxiv.org/pdf/1310.2559.pdf $\endgroup$ Jan 11, 2018 at 14:20
  • 1
    $\begingroup$ Thanks, but where in that reference do they state or give the Wick-Isserlis theorem or formula? $\endgroup$ Jan 11, 2018 at 14:27
  • 1
    $\begingroup$ @NawafBou-Rabee Another one thphys.uni-heidelberg.de/~amendola/teaching/compstat-hd.pdf, page 28, but you can find plenty of them by yourself. Unfortunately, the Wick-Isserlis formula for moments of order $k$ involve $k!!$ terms. Hence, with $k>>1000$, a different method is required for evaluating my moments. Astonishingly, for the time being, I'm unable to find numerical methods and algorithms for high-order moments of Gaussian r.v.. $\endgroup$ Jan 11, 2018 at 16:03

1 Answer 1


It seems you may want to turn the large $k$ property to your advantage and use a Laplace-type asymptotic method with controlled bounds. For fixed $d$ you may find some useful tools in the book "Analytic Combinatorics in Several Variables" by Pemantle and Wilson.

If also $d$ becomes large, then this starts looking like a problem in statistical field theory. There is a method people used for this kind of thing called the Lipatov argument (see for instance this article by Spencer).

  • $\begingroup$ Thanks for your input. $d$ is fixed, $k$ is slightly variable but no asymptotics. For sure, good old Laplace is always an option, but I guess there exist dedicated methods and algorithms for those fundamental integrals. $\endgroup$ Jan 12, 2018 at 8:11
  • $\begingroup$ Please, can you elaborate a little bit on the Laplace approximation? I tried the standard and fully exponential forms but, unless I'm mistaken, it doesn't work (see the question). I'm unable to find relevant sections in Permantle & Wilson too... $\endgroup$ Feb 13, 2018 at 15:31

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .