Predictive Hacks

Estimate Pi with Monte Carlo

One method to estimate the value of \(\pi\) is by applying the Monte Carlo method. Let’s consider a circle inscribed in square. Then for a radious \(r\) we have:

Estimate Pi with Monte Carlo 1

Area of the Circle = \(\pi r^2\)

Are of the Square = \(4 r^2\)

Hence, we can say that \(\frac{AreaCircle}{AreaSquare}=Prob=\frac{\pi}{4}\). This implies that \(\pi = 4Prob\)

We know that the Equation of Circle with center (j,k) and radius (r) is:

\( (x_1-j)^2+(x_2-k)^2 = r^2\) where in our case the center is (0,0) and the radius is 1. This means, that in order to calculate the Area of the Circle we need to consider all the points from the uniform distribution which have a distance from the center (0,0) less than or equal to 1, i.e:

\(\sqrt{ x_1^2+x_2^2} \leq 1\)

Notice that since we chose the center to be (0,0) and the radius to be equal to 1, it means that the co-ordinates of the square will be from -1 to 1. For that reason we simulate from the uniform distribution with min and max values the -1 and 1 respectively.


Example Using R

Below, we represent how we can apply the Monte Carlo Method in R to estimate the pi.


# set the seed for reproducility
set.seed(5)

# number of simulations
n=1000000

# generate the x1 and x2 co-ordinates from uniform 
# distribution taking values from -1 to 1
x1<-runif(n, min=-1, max=1)
x2<-runif(n, min=-1, max=1)

# Distrance of the points from the center (0,0)
z<-sqrt(x1^2+x2^2)

# the area within the circle is all the z
# which are smaller than the radius^2
# in our case radius=1
4*sum((z<=1))/length(z)
[1] 3.14204

Let’s have a look also at the simulated data points.

InOut<-as.factor(ifelse(z<=1, "In", "Out"))

plot(x1,x2, col=InOut, main="Estimate PI with Monte Carlo")
Estimate Pi with Monte Carlo 2

As we can see, with 1M simulations we estimated the PI to be equal to 3.14204. If want to get a more precise estimation we can increase the number of simulations.



Example Using Python

import pandas as pd
import numpy as np
import math
import seaborn as sns
%matplotlib inline

# number of simulations
n=1000000

x1=np.random.uniform(low=-1,high=1, size=n)
x2=np.random.uniform(low=-1,high=1, size=n)

z=np.sqrt((x1**2)+(x2**2))

print(4*sum(z<=1)/n)
3.143288

Let’s have a look also at the simulated data points.

inout=np.where(z<=1,'in','out')

sns.scatterplot(x1,x2,hue=inout,legend=False)
Estimate Pi with Monte Carlo 3

Share This Post

Share on facebook
Share on linkedin
Share on twitter
Share on email

Leave a Comment

Subscribe To Our Newsletter

Get updates and learn from the best

More To Explore

letter frequency
Python

Document Letter Frequency in Python

Letter Frequency We will provide you a walk-through example of how you can easily get the letter frequency in documents

[the_ad_group id="232"]
[the_ad id="2133"]