# Are all Monte Carlo algorithms to approximate $\pi$ equivalent?

There are several ways one can approximate $\pi$ using a Monte-Carlo type algorithm. For example, one can draw random points in the unit square, and approximate $\pi$ via the ratio of points that fall in the unit disk. Another way would be to use an algorithm based on Buffon's needle problem: from the proportion of needles that fall across the parallel lines, one can approximate $\pi$. One can presumably devise many other random algorithms to estimate $\pi$.

My question is: in terms of speed of convergence, are all these algorithms equivalent? Or are there some that will yield correct decimals of $\pi$ faster? What is the Monte Carlo method to estimate $\pi$ that has the fastest convergence?

• Different Monte Carlo procedures converge at different rates, so Yes there can be more and less efficient simulation algorithms to approximate $\pi.$ Apr 30, 2018 at 2:37
• @BruceET - since you are the first to answer, I added a more difficult question for you ;-) Apr 30, 2018 at 2:42
• @Frank while it's unlikely this is the best monte carlo integration technique, a good starting point for research is the "Metropolis Hastings algorithm". From the cursory looking I've seen, many modern techniques are directly related to this.
– Mark
Apr 30, 2018 at 2:46
• There is a very simple algorithm that can calculate an approximation of $\pi$ that is accurate to within 47 decimal places using only an average of 2 random coin flips. Apr 30, 2018 at 3:01
• @Michael Flip a coin until heads comes up and then look up and write down $\pi$ to $47$ places? Apr 30, 2018 at 3:05

Trying to keep this simple, answering your original very good question.(+1) Using Riemann approximation by rectangles to get the area under the curve $y = \sqrt{1-x^2},$ for $0 < x < 1,$ and multiplying by 4, you can get pretty close to $\pi$ with only a dozen rectangles. Closer if you use trapezoids. This is not a Monte Carlo method, but the first method below is somewhat analogous (using random points $x,$ instead of a deterministic grid of center points as bases of rectangles).

set.seed(430)  # retain for exactly the same simulation
m = 10^6;  x = runif(m);  h = sqrt(1-x^2)
mean(h)*4
[1] 3.142033           # estimate of pi
2*4*sd(h)/sqrt(m)
[1] 0.001784625        # rough 95% margin of simulation error


So a million random $x$-values gives $3.1420 \pm 0.0018.$

An acceptance-rejection method similar to the one you suggested is somewhat less efficient.

set.seed(430)
m = 10^6;  x = runif(m);  y=runif(m)
acc = y <+ sqrt(1-x^2)
mean(acc)*4
[1] 3.143212           # aprx pi
2*4*sd(acc)/sqrt(m)
[1] 0.003282115        # marg of sim err


Here we get $3.1433 \pm 0.0033,$ which indicates somewhat slower convergence.

Notes: (a) My purpose here is not to give best simulation methods for approximating $\pi,$ but to answer your question by showing two commonly used methods with different rates of convergence. (b) I do not challenge the importance of the Metropolis-Histings algorithm in simulation, but in my experience it is more useful in dimensions higher than one or two. (c) I am not from the 'show-me' state of Missouri, but I have a sister who lives there; if there is a way to get $\pi$ to 47 decimal places by simulating a couple of coin flips, that's news to me and I'd want to see an authoritative reference for that. (d) I have illustrated other Monte Carlo methods on another page; some of them can be adapted to your question.

All random algorithms that approximate $\pi$ by trying something that succeeds with probability $f(\pi)$ lots of times, getting an estimate of $f(\pi)$, and using that to estimate $\pi$ will be approximately equivalent. (Whether the thing we try is "drop a needle that interects a line" or "pick a point that lands inside a circle" or anything else.)

The reason is that their convergence rate will have nothing to do with the method at all. If the probability we're trying to estimate is $p$, then after $n$ trials the number of successes will follow a binomial distribution, and the typical number of successes will be $np \pm C \sqrt{np(1-p)}$. This gives an estimate of $p$ to within $O(\frac1{\sqrt n})$ error, which means millions of trials before we can estimate $p$ to three decimal places.

I only say "approximately" equivalent because the probability we're estimating can depend on $\pi$ in lots of ways. If the probability of success is $\frac12 + \frac{\pi}{1000}$, then it will take many trials before our estimate of this probability becomes significantly better than "it's about $\frac12$". On the other hand, if the probability of success is $10000\pi - 31415$, then estimating that the probability is approximately $0.9$ is enough to give use $\pi$ to five digits after the decimal: $\pi \approx 3.14159$.

• So, according to your last paragraph, can we find a shape in the "pick a point that lands inside a circle" that would have a better relationship between $p$ and $\pi$ than the circle? Apr 30, 2018 at 22:16
• One possibility, using this integral, is to sample a point from $[0,1] \times [0,\frac1{320}]$ and check if it's below the curve $\frac{x^4(1-x)^4}{1+x^2}$: the probability is $320(\frac{22}{7}-\pi)$. But philosophically, all that these methods are doing is offloading a few already-known digits of $\pi$ into the setup of the method itself, not accelerating the speed at which new digits are found. Apr 30, 2018 at 22:29
• Nice - I hear you about the offloading going on here. Just out of curiosity, is there a systematic way to build a series of polynomials that take on more and more decimals of $\pi$ (or any number), or are these polynomials discovered through chance? May 3, 2018 at 15:39
• All I know is what you can find in the references in the Wikipedia article I linked to. So it seems there are some systematic approaches. May 3, 2018 at 19:22