Game 10: Simulation:Law of Large Numbers (rev)

Keywords:

Simulation, Distribution (pmf, pdf), Law of Large Number

Questions

1. How does \(\bar{X} = \frac{1}{n} \sum_{i=1}^n X_i\) behave?
2. How does median(\(X\)) behave?
3. What are the assumptions, conditions and some issues we must take into account to answer these questions?
4. What do we mean by $X_1, ,,,, X_n $ some distribution, say \(Bernoulli(p)\) or \(N(\mu, \sigma^2)\)?
5. What is simulation? What is it for?

die <- 1:6
roll <- function(n) {
  mean(sample(die, size = n, replace = TRUE))
}
plot(sapply(1:1000, roll), type = "l", xlab = "# of dice", ylab = "average")
abline(h = 3.5, col = "red")

# Repeat the experiment 4 times
set.seed(12)
par(mfrow=c(2,2)); 
# Repeat it for 4 times
for (i in 1:4) # For loop, but ....
{plot(sapply(1:500, roll), type = "l", xlab = "# of dice", ylab = "average")
abline(h = 3.5, col = "red")}

For Bernoulli(p)

m<-1; p<-0.3;   
mean_rs <- function(n) {
  mean(rbinom(n,m,p))
}
set.seed(1234)
par(mfrow=c(2,2)); 
# Repeat it for 4 times
for (i in 1:4) # For loop, but ....
{plot(sapply(1:500, mean_rs), type = "l", ylim = c(0,1) ,xlab = "# of Bernoulli tosses", ylab = "average") 
abline(h = m*p, col = "red")}

Bin(m,p)

m<-10; p<-0.4;   
mean_rs <- function(n) {
  mean(rbinom(n,m,p))
}
set.seed(1234)
par(mfrow=c(2,2)); 
# Repeat it for 4 times
for (i in 1:4) # For loop, but ....
{plot(sapply(1:500, mean_rs), type = "l", ylim = c(m*p-2,m*p+2) ,xlab = "# of Binomial experiments", ylab = "average") 
abline(h = m*p, col = "red")}

Exercises:

1. Produce the similar simulation for N(0,1), N(0,10); N(10,1), N(10,10)
2. Produce the similar simulation for U(0,1), \(Y=Z^2\) where \(Z \sim N(0,1)\)

nmu<-0; nsd<-10;
mean_rs <- function(n) {
  mean(rnorm(n,nmu,nsd))
}
# c<-1; summary(rnorm(20,nmu,c*nsd))
set.seed(1234)
par(mfrow=c(2,2)); 
# Repeat it for 4 times
for (i in 1:4) # For loop, but ....
{plot(sapply(1:100, mean_rs), type = "l", ylim=c(-10,10), xlab = "# of normal samples", ylab = "average") 
abline(h = nmu, col = "red")}

Estimate \(E(Z^2)\) using LLN

Note that the theoretical value of \(E(Z^2)=E(\chi^2_1)=1\)

nmu<-0; nsd<-1;
mean_rs <- function(n) {
  mean(rnorm(n,nmu,nsd)^2)
}
# c<-1; summary(rnorm(20,nmu,c*nsd))
set.seed(1234)
par(mfrow=c(2,2)); 
# Repeat it for 4 times
for (i in 1:4) # For loop, but ....
{plot(sapply(1:100, mean_rs), type = "l", ylim=c(0,3), 
       xlab = "# of Z^2 samples", ylab = "average") 
abline(h = 1, col = "red")}

References

1. An LLN example in Stackoverflow
2. Nice R Code
3. Demonstration of the Central Limit Theorem