# Example: LLN
# This example will demonstrate the law of large numbers theorem. 
# Step 1. Please generate random samples from the standard normal distribution 
#         for 1000 times with different sample size T. Let¡¦s say four cases here, 
#         T=50, T=100, T=300, and T=1000.
# 
# Step 2. For each different sample size T, please calculate the sample average 
#         each time
# Step 3. Plot histogram of sample averages. 
#
# You may observe the average of the results obtained from a large number of trials 
# should be close to the expected value, and will tend to become closer as more trials 
# are performed.


# t: sample size
# n: repetition time

# clean up workspace
rm(list=ls(all=TRUE))

fun_LLN <-function(t,n){
	x <-numeric(n)
	for(i in 1:n){
		x[i] <- mean(rnorm(t))  
	}
	x
}

par(mfrow=c(2,2))
hist(fun_LLN(50,1000), xlim=range(-1,1),freq=F, main='T=50',xlab='Sample Mean')
hist(fun_LLN(100,1000), xlim=range(-1,1),freq=F, main='T=100',xlab='Sample Mean')
hist(fun_LLN(300,1000), xlim=range(-1,1),freq=F, main='T=300',xlab='Sample Mean')
hist(fun_LLN(1000,1000), xlim=range(-1,1),freq=F, main='T=1000',xlab='Sample Mean')




# Example: CLT 
# Step1: Generate random samples with sample sizes T=50, 100, 300, and 1000
#        from the following distribution: Chi-squared £q2(1) distribution
# Step2: compute the normalized sample average for each sample:
#        sqrt(T)(x_bar-£g)/£m
#        where x_bar, £g, and £m are the sample average, mean, and standard deviation,
#        respectively. 
# Step3: Repeat this procedure 1000 times and plot the resulting histograms.

# t: sample size
# n: repetition time

# chiq-square(df=1), mu=1, var=2
# clean up workspace
rm(list=ls(all=TRUE))

fun_CLT <-function(t,n){
  x <- numeric(n)
  s.x <-numeric(n)
  k <- 1
  for (i in 1:n){
    x[i] <-mean(rchisq(t,df=1)) 
    s.x[i] <-sqrt(t)*(x[i]-k)/(sqrt(2*k))
  }
  s.x
}

par(mfrow=c(2,2))
xtemp <-fun_CLT(50,1000)
hist(xtemp,xlim=range (-4,4),ylim=range(0,0.45),freq=FALSE,main='T=50',xlab='Sample Mean')
lines(density(xtemp))
curve(dnorm(x, mean=0, sd=1), add=TRUE,col="red")

xtemp <-fun_CLT(100,1000)
hist(xtemp,xlim=range (-4,4),ylim=range(0,0.45),freq=FALSE,main='T=100',xlab='Sample Mean')
lines(density(xtemp))
curve(dnorm(x, mean=0, sd=1), add=TRUE,col="red")

xtemp <-fun_CLT(300,1000)
hist(xtemp,xlim=range (-4,4),ylim=range(0,0.45),freq=FALSE,main='T=300',xlab='Sample Mean')
lines(density(xtemp))
curve(dnorm(x, mean=0, sd=1), add=TRUE,col="red")

xtemp <-fun_CLT(1000,1000)
hist(xtemp,xlim=range (-4,4),ylim=range(0,0.45),freq=FALSE,main='T=1000',xlab='Sample Mean')
lines(density(xtemp))
curve(dnorm(x, mean=0, sd=1), add=TRUE,col="red")