# Monte Carlo Integration for non-square integrable functions

I hope this is the right place to ask, if not feel free to move it to a more appropriate forum.

I've been wondering for quite a while now how to treat non-square integrable functions with Monte Carlo Integration. I know that MC still gives a proper estimate but the error is unrealiable (divergent?) for those kind of functions.

Let's restrict us to one dimension. Monte Carlo integration means that we approximate the integral

$$I = \int_0^1 \mathrm{d}x \, f(x)$$

using the estimate

$$E = \frac{1}{N} \sum_{i=1}^N f(x_i)$$

with $x_i \in [0,1]$ uniformly distributed random points. The law of large numbers makes sure that $E \approx I$. The sample variance

$$S^2 = \frac{1}{N-1} \sum_{i=1}^N (f (x_i) - E)^2$$

approximates the variance $\sigma^2$ of the distribution induced by $f$. However, if $f$ is not square-integrable, i.e. the integral of the squared function diverges, this implies

$$\sigma^2 = \int_0^1 \mathrm{d} x \, \left( f(x) - I \right)^2 = \int_0^1 \mathrm{d} x \, f^2(x) - I^2 \longrightarrow \infty$$

meaning that also the variance diverges.

A simple example is the function

$$f(x) = \frac{1}{\sqrt{x}}$$

for which $I = \int_0^1 \mathrm{d}x \, \frac{1}{\sqrt{x}} = 2$ and $\sigma^2 = \int_0^1 \mathrm{d}x \, \left( \frac{1}{x} - 2 \right) = \left[ \ln x - 2x \right]_0^1 \rightarrow \infty$.

If $\sigma^2$ is finite one can approximate the error of the mean $E$ by $\frac{S}{\sqrt{N}} \approx \frac{\sigma}{\sqrt{N}}$, but what if $f(x)$ is not square-integrable?

• I don't get it: you start out by noting that none of the $E_i$ has a variance and then ask whether the variance of their average would be a reasonable estimator of--that nonexistent variance! Or do I misread this question: perhaps by "statistically independent estimations" you have some different (perhaps robust) estimator of the integral in mind?
– whuber
Jun 21, 2013 at 18:14
• I didn't say $E$ doesn't have a variance, only that I cannot define a variance for it by $S^2$. So the question is whether I can define an error at all and if $\bar{S}^2$ is a reasonable candidate. By statistically independent I mean that the $E_i$ are obtained using different random numbers, e.g. by using differently seeded random number generators (I hope thats the right term then). Jun 22, 2013 at 6:38
• Please explain what you mean by not being able to "define a variance for it by $S^2$." I cannot make sense of this using the standard definitions of variance and $S^2$.
– whuber
Jun 24, 2013 at 14:45
• Well, the function is not square-integrable so, if I am not mistaken, $S^2$ should diverge. If this is the case the definition for $S^2$ makes no sense in the first place, right? By means of the central limit theorem, however, $E$ will still converge to the true value of the integral, but without an error this value alone makes no sense (how 'good' is this result?). Jun 26, 2013 at 13:55
• Sorry, I meant to say "law of large numbers" of course, not CLT. Jun 26, 2013 at 14:04

## 1 Answer

You could just use other scale/dispersion measures such as the interquantile range, which are not affected by tail asymptotics and thus square integrability. With the added benefit that often they are in general more robust anyway.

Obviously one shall apply them to a resampling/bootstrap followed by the mean estimator, not directly just to the raw output from MC sampling of the function before averaging. You can also check in general L-estimators and adapt one of them to merge these two steps into one for performance, but mentally the two distributions shall not be confused, even though the estimator PDF will naturally inherit some characteristics (included maybe lack of square integrability).

• +1, I should add that law of large numbers does not require second moments, so this is a perfectly good advice. Jul 1, 2013 at 13:48
• Thanks for your answer! I have to admit that I read those terms for the first time, but from looking them up at WP I think your answer points me in the right direction. Could you or someone else suggest some articles or books that explains the subjects in more detail? Jul 1, 2013 at 14:30
• I notice now that maybe my answer was a bit unclear. Since you are simulating you don't really need resampling/bootstrapping, in theory you could just add further new samples instead and obtain an empirical distribution for the mean estimator. Only if resources are a concern then you can precalculate partial averages and resample 'em, but the statistics will not be trivial if well done. I'm no boostrap expert so I'll leave advice on that to others, just wanted to point it out if you need to go beyond the straightforward formulation. Concentrate on dispersion measures first, optimize later. Jul 1, 2013 at 15:11
• The mean estimator proposed does not have a finite variance. It does not matter if one adds further samples, the empirical distribution of the estimator will ALSO have non-finite variance. You can confirm this with a few simulations. Jul 2, 2013 at 14:16
• Sure, in fact that's what was being discussed and the reason why one shall use another dispersion measure. Jul 2, 2013 at 15:02