# weak law of large numbers

# examine probability of being within n*eps of the 
# expected value for sequence of Bernoulli RVs

# histograms of sample means as n increases
S_n <- function(n,p) sum(sample(0:1,n,prob=c(1-p,p),replace=TRUE))
n <- 100
p <- 0.1
sim <- replicate(10000,S_n(n,p)/n)
hist(sim,xlim=c(0,.25))

n <- 1000
p <- 0.1
sim <- replicate(10000,S_n(n,p)/n)
hist(sim,xlim=c(0,.25))

n <- 10000
p <- 0.1
sim <- replicate(10000,S_n(n,p)/n)
hist(sim,xlim=c(0,.25))

# use binomial cdf to calculate prob that
# sample mean is within epsilon of mu
wlln <- function(n,p=0.1,eps=0.01) {
  pbinom(n*p+n*eps,n,p)-pbinom(n*p-n*eps,n,p)
}

p <- 0.1
eps <- 0.01
n <- 100
wlln(n,p,eps)

n <- 1000
wlln(n,p,eps)

n <- 10000
wlln(n,p,eps)

n <- seq(from=100,to=8000,by=100)
P <- sapply(n,wlln)
plot(n,P)