library(dplyr)
library(ggplot2)
library(knitr)
opts_chunk$set(message = FALSE, eval = TRUE)

Simulate Buffon’s Needle

n <- 20
d     <- runif(n = n, min = 0, max = 1/2)
theta <- runif(n = n, min = 0, max = pi)

(df <- data_frame(d, theta))
## # A tibble: 20 x 2
##              d      theta
##          <dbl>      <dbl>
## 1  0.421537737 0.08101089
## 2  0.080696171 1.10375158
## 3  0.008729169 1.19186020
## 4  0.449464490 0.83817382
## 5  0.328243365 0.88036281
## 6  0.499262640 0.53429573
## 7  0.372272875 0.07884363
## 8  0.498089879 0.73014834
## 9  0.178281839 2.28789951
## 10 0.157958551 1.30119270
## 11 0.466539003 1.23536838
## 12 0.402765643 1.05361257
## 13 0.474079798 2.25967411
## 14 0.366149133 2.10208444
## 15 0.393188542 1.05431343
## 16 0.439498770 3.13983823
## 17 0.244888635 1.29131678
## 18 0.197206581 0.55139906
## 19 0.080525428 1.37555721
## 20 0.327939279 0.58738903
sin(theta)/2 > d
##  [1] FALSE  TRUE  TRUE FALSE  TRUE FALSE FALSE FALSE  TRUE  TRUE  TRUE
## [12]  TRUE FALSE  TRUE  TRUE FALSE  TRUE  TRUE  TRUE FALSE
mean(sin(theta)/2 > d)
## [1] 0.6
2/pi
## [1] 0.6366198

Visualize

(p <- ggplot(data = df, aes(x = theta, y = d)) +
  scale_x_continuous(expand = c(0, 0), limits = c(0, pi)) +
  scale_y_continuous(expand = c(0, 0), limits = c(0, 1/2)) +
  stat_function(fun = function(x) {sin(x)/2}) +
  theme_bw())

df <- mutate(df, cross = sin(theta)/2 > d)
p + geom_point(data = df, aes(color = cross))

Estimating Pi

2 * n / sum(df$cross)
## [1] 3.333333