Monte Carlo integrate Improper integral

Question

I am trying to use Monte Carlo method to integrate the following improper integral

$${1 \over \sqrt{2\pi}}\int\limits_{-\infty}^{\infty}x^4e^{(-x^2/2)} \ dx$$

Using change of variables, $ y = x^2/2 $, I can transform the integral into gamma function and get the solution as $3\sqrt{2}\pi$
(as shown here)

However, the integral still remains improper at $0$ to $\infty$ and I am not sure how to get that into a proper one so that I can use Monte-Carlo method to arrive at an approximate solution.

Any help is much appreciated.

Thank you

Answer 1

What you have written down is incorrect. $E[X^4]=3$. Thus your answer should be close to 3 instead.

I will denote $p(x)$ as standard normal density $1/\sqrt{2\pi}\exp(-x^2/2)$. Then you want $\int_{-\infty}^{\infty} x^4p(x)dx$ which is really $E[X^4]$ by definition.

Recall weak law of large number for $Z_1,\dots, Z_n$ iid following distribution $f(Z)$, then $\bar{Z}=\frac{Z_1+\dots+Z_n}{n}$ has for any $\epsilon>0$, $P(|\bar{Z}-E[Z]|>\epsilon)\to 0$ as $n\to\infty$. Here you want to choose $Z=X^4$ and by normal MGF, you will see it has $E[Z]<\infty$. The only matter is to compute it.

To simulate $Z_1,\dots, Z_n$, it suffices to simulate $X_1,\dots, X_n$ following $p(x)$ density, where $Z_i=X_i^4$.

Draw $X_1,\dots, X_n$ iid from $p(x)$, then compute $X_i^4=Z_i$. At last take the mean and apply weak law of large number.

Here is the code for simulation in R.

me_calc=c()
for(i in 1:1000){
x=rnorm(100000,mean = 0,sd=1)
z=x^4
me_calc=c(me_calc,mean(z)-3)
}
hist(me_calc)

Here is the simulation histogram and you can see it is fairly close to the true answer 3.