Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

how to to write a function that demonstrates the central limit theorem with graphics

I'm trying to write a function that creates animated graphics(without using the animation package) where users can control the inputs(sample size and distribution ect..) that demonstrates the central limit theorem. This is what in theory I would want but having trouble with writing the function where users can actually control the inputs as I mentioned above.

msample <- NA # set up empty vector
ns <-3 # sample size
for(i in 1:500){
sam <- runif(ns) * 10 # draw sample
 msample[i] <- mean(sam) # save mean of sample
h <- hist(msample, breaks=seq(0,10, len=50), # histogram of all means
xlim=c(0,10), col=grey(.9),
xlab="", main="Central Limit Theorem", border="blue", las=1)
points(sam, rep(max(h$count), length(sam)),
pch=16, col=grey(.2)) # add sampled values
points(msample[i], max(h$count), # add sample mean value
col="red", pch=15)
text(10, max(h$count), paste("sample no", i))
hist(msample[i], breaks=seq(0,10, len=50), # ovelay sample mean 
xlim=c(0,10), col="red", add=T, # in histogram
xlab="", border="white", las=1)
Sys.sleep(.05)
}
like image 688
geneteics_diva Avatar asked Nov 02 '22 23:11

geneteics_diva


1 Answers

It is not clear what do you want as result. but I think, you can put you code in a function and use dot argument ... as a solution to give extra parameters (distribution parameters for example).

   central.simul <- function(N, ns,type = c("runif", "rnorm", "rbinom"),...){
        type <- match.arg(type)
        msample <- rep(NA,N)  ## EDIT here: intialisation
        for(i in 1:N){
          sam <- switch(type,
                        runif = runif(ns)*10,
                        rnorm = rnorm(ns)*10,
                        rbinom = rbinom(ns,...))
          msample[i] <- mean(sam) # save mean of sample
          add.hist <- i > 1
          h <- hist(msample, breaks=seq(0,10, len=50), # histogram of all means
                    xlim=c(0,10), col=grey(.9),
                    xlab="", main="Central Limit Theorem", border="blue", las=1,add=add.hist)
          points(sam, rep(max(h$count), length(sam)),
                 pch=16, col=grey(.2)) # add sampled values
          points(msample[i], max(h$count), # add sample mean value
                 col="red", pch=15)
          text(10, max(h$count), paste0("sample no ", i))
          hist(msample[i], breaks=seq(0,10, len=50), # ovelay sample mean 
               xlim=c(0,10), col="red", add=T, # in histogram
               xlab="", border="white", las=1)
          Sys.sleep(.1)
        }
    }

You can call it using :

central.simul(10,3,'runif')
central.simul(10,3,'rbinom',size=2,prob=0.5)

The code as it don't work for rnorm for example ( you should modify breaks I think), but it should be a good start.

like image 76
agstudy Avatar answered Nov 08 '22 04:11

agstudy