A demonstration of central limit theorem

simCLT.R

library(ggplot2)
library(dplyr)
library(tidyr)
# integrate(f = function(x)x*(0.6*dnorm(x, -10/6)+0.4*dnorm(x, 10/4)), lower = -Inf, upper = Inf)

scaled_mean <- function(x){sqrt(length(x))* mean(x) / ifelse(length(x)==1L,1,sd(x))}
CLTsim <- function(n, iter){
  mu <- c(-10/6, 10/4)
  mixture <- numeric(iter)
  uniform <-  numeric(iter)
  gamma1 <- numeric(iter)
  gamma2 <- numeric(iter)
  normal <- numeric(iter)
  
  for ( i in seq_len(iter) ){
    z <- rbinom(n, 1, 0.4)
    mixture[i] <- scaled_mean( rnorm(n, mu[z+1]) )
    uniform[i] <-  scaled_mean( runif(n, -2, 2) )
    gamma1[i] <- scaled_mean( rgamma(n, 2, 1) - 2 )
    gamma2[i] <- scaled_mean( rgamma(n, 0.75, 1) - 0.75 )
    normal[i] <- scaled_mean( rnorm(n) )
  }
  data.frame(n, normal, uniform, gamma1, gamma2, mixture)
}
set.seed(1234)
resdf <- dplyr::bind_rows(lapply(c(1, 10, 30), CLTsim, iter = 10000))
resdf <- pivot_longer(resdf, normal:mixture) 
resdf <- mutate(resdf, name = factor(name, levels = c("normal","uniform","gamma1","gamma2","mixture")))

p1 <- ggplot(resdf, aes(x=value)) +
  geom_histogram(bins=101, aes(y=after_stat(density)))+
  facet_grid(n~name)+
  theme_bw(16)+labs(title = "sampling distribution of means (scaled)")
print(p1)

p2 <- ggplot(resdf, aes(x=value)) +
  stat_ecdf()+
  stat_function(fun = pnorm, colour="royalblue", linetype=2)+
  facet_grid(n~name)+
  theme_bw(16)+labs(title = "sampling distribution of means (scaled)")
print(p2)

ggsave(filename = "CLThist.png", plot = p1, width=10, height = 8)
ggsave(filename = "CLTecdf.png", plot = p2, width=10, height = 8)