[This article was first published on

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

**Stats raving mad » R**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

It’s π-day today so we gonna have a little fun today with Buffon’s needle and of course R. A well known approximation to the value of $latex \pi$ is the experiment tha Buffon performed using a needle of length,$latex l$. What I do in the next is only to copy from the following file the function estPi and to use an ergodic sample plot… Lame,huh?

estPi<- function(n, l=1, t=2) { m <- 0 for (i in 1:n) { x <- runif(1) theta <- runif(1, min=0, max=pi/2) if (x < l/2 * sin(theta)) { m <- m +1 } } return(2*l*n/(t*m)) }

So, an estimate would be…

```
estPi(2000,l=1,t=2)
# 3.267974
```

Ok, not that great but for the whole scene it’s remarkable good! Now, we set some increasing sample sizes to account for the estimation.

```
n=8000
r=15
mat=rep(NA,r)
size=rep(NA,r)
for (i in 1:r) {
size[i]<-n*i
mat[i]<-estPi(n*i,l=1,t=2)
}
matrix<-expand.grid(size)
matrix[,2]<-mat
names(matrix)<-list("n","pi")
matrix
# n pi
#1 8000 3.182180
#2 16000 3.165809
#3 24000 3.135615
#4 32000 3.145581
#5 40000 3.138486
#6 48000 3.144860
#7 56000 3.162412
#8 64000 3.111932
#9 72000 3.097574
#10 80000 3.155072
#11 88000 3.157404
#12 96000 3.144139
#13 104000 3.126597
#14 112000 3.150226
#15 120000 3.136599
```

Which is the best estimate?

```
matrix[which.min(abs(matrix[,2]-pi)),]
# n pi
# 12 96000 3.144139
plot(matrix,type="b");abline(h=pi,col="red",lty=2)
```

source : [Chiara Sabatti , pdf]

Take a look @

+ Wiki

+ An introduction to geometrical probability: distributional aspects with applications (A. M. Mathai)

To

**leave a comment**for the author, please follow the link and comment on their blog:**Stats raving mad » R**.R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.