(This article was first published on

**Statistical Research » R**, and kindly contributed to R-bloggers)Recently I’ve been working on some jackknife and bootstrapping problems. While working on those projects I figured it would be a fun distraction to take the process and estimate pi. I’m sure this problem has been tackled countless times but I have never bothered to try it using a Monte Carlo approach. Here is the code that can be used to estimate pi. The idea is to generate a 2 by 2 box and then draw a circle inside of it. Then start randomly selecting points inside the box. This example uses 10,000 iterations using a sample size of 500. An interesting follow-up exercise would be to change the sample size, mean, and standard deviation to and see what it does to the jackknife variance.

library(bootstrap); nsims <- 10000; size <- 500; alpha <- .05; theta <- function(x){sd(x)} pi.estimate <- function(nsims,size,alpha){ out <- as.data.frame(NA); call <- match.call(); pi.sim <- matrix(NA, nrow=nsims); for(i in 1:nsims){ x <- runif(size,-1,1); y <- runif(size,-1,1); radius <- -1; x.y <- cbind(x,y); d.origin <- sqrt(x^2+y^2); x.y.d <- cbind(x.y,d.origin); est.pi <- 4*length(x.y.d[,3][x.y.d[,3]<1])/ length(x.y.d[,3]) pi.sim[i] <- est.pi; } pi.mean <- mean(pi.sim); jk <- jackknife(c(pi.sim), theta); jk.me <- jk$jack.se*qnorm(1-alpha/2); pi.mean <- pi.mean; pi.se <- jk$jack.se; jk.me <- jk.me; pi.me.u <- pi.mean+jk.me; pi.me.l <- pi.mean-jk.me; return(list(pi.mean=pi.mean, pi.se=pi.se,pi.me.u=pi.me.u,pi.me.l=pi.me.l, call=call)); } pi.estimate(nsims,size,alpha);

To

**leave a comment**for the author, please follow the link and comment on his blog:**Statistical Research » R**.R-bloggers.com offers

**daily e-mail updates**about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...