Sunday, August 1, 2010

Estimating Pi with Monte Carlo Integration

There are several tutorials on how to perform this calculation. I like to add the technique I use with MATLAB recently. Basically, Imagine you throw darts at random (really really random) into a board with a circle in it, ie. like this.

The concept can be implement in one MATLAB command:

>> 4*mean(rand(1,5000).^2+rand(1,5000).^2 <= 1)
ans =
   3.141592653589793

or, a slightly different approach, using Monte Carlo Integration with random numbers drawn from uniform probability distribution $x_i \sim Unif(0,1)$.

>> 4*mean(sqrt(1-rand(1,5000).^2))
ans =
   3.141432627301033

Notice that the second approach requires only one series of random numbers, instead of two as the first approach but it has to perform a more complex square root operation.