Calculate π Using a Monte Carlo Method

As this Wikipedia page suggests, we can calculate π using a Monte Carlo method:

• Draw a square on the ground, then inscribe a circle within it.
• Uniformly scatter some objects of uniform size (grains of rice or sand) over the square.
• Count the number of objects inside the circle and the total number of objects.
• The ratio of the two counts is an estimate of the ratio of the two areas, which is π/4. Multiply the result by 4 to estimate π.

Let's write an R function to implement this algorithm. Set the square to be of length 1, in region with 0 < x < 1 and 0 < y < 1. The inscribed circle is the region with 0 < x < 1, 0 < y <1 and x² + y² < 1. We use runif() to simulate the uniform scattering of points. The following function computePiMC() performs a simulation with N scattered points.

computePiMC <- function(N) {
    x <- runif(N)
    y <- runif(N)
    MISSING CODE
}

a. What could be the missing code that replaces MISSING CODE above to complete this program?
4*sum(x^2 + y^2 < 1)/N
4*sum(x+y > 1)
4*sum(x+y < 1)
4*sum(x+y > 1)/N
4*sum(x+y < 1)/N
4*sum(x^2 + y^2 > 1)
4*sum(x^2 + y^2 > 1)/N
4*sum(x^2 + y^2 < 1)

b. The three lines inside the computePiMC() function can be replaced by a single line without changing the result. Which of the following could be the line of code?
4*sum(runif(N) + runif(N) < 1)/N
4*sum(runif(N)*runif(N) + runif(N)^2 < 1)
4*sum(runif(N)^2 + runif(N)*runif(N) > 1)/N
4*sum(runif(N)^2 + runif(N)^2 < 1)/N
4*sum(runif(N)*runif(N) + runif(N)*runif(N) < 1)
4*sum(runif(N)^2 + runif(N)^2 < 1)

c. Now you are going to use the computePiMC(N) to perform a simulation with N=10⁷. Run the following code:

set.seed(600544)
piMC <- computePiMC(1e7)

What is your estimate for π? Give your answer to at least 4 decimal places.

$\hat{\pi}$ =

The estimated value is only accurate within chance variation. It is therefore more useful to provide a confidence interval. Let X be a random variable that counts the proportion of points falling into the circle. Set up a box model consisting of two tickets "0" and "1". The fraction of the "1" ticket in the box is π/4. The value of X is equivalent to drawing N tickets from the box with replacement and taking the average the tickets. The expected value of X is E(X) = π/4 and the standard error is SE(X) = √p (1-p)/N, where p = π/4. (For those of you who have learned binomial distribution in the formal way, NX follows a binomial distribution with n=N and p=π/4.) Since we don't know π (the number we're trying to estimate), we use $\hat{\pi}$ as a substitute. As explained in this Week 7 problem, the 100%*level confidence interval for E(X) can be computed by the formula

$\frac{\hat{\pi}}{4} \pm Z_{\alpha/2} SE$,

where Z_α/2 is the Z score such that the area under the standard normal curve with Z > Z_α/2 = α/2 = (1 - level)/2. It follows that the confidence interval for π is given by

$\hat{\pi} \pm 4 Z_{\alpha/2} SE$.

d. Use the simulation result in part (c) to calculate the 99% confidence interval for π. Give your answers to 4 decimal places.

99% CI for π = ( , ).

e. The value of π is 3.14159265358979.... Is π inside the confidence interval calculated above?
yes
no