mean.dr <- function (n, m){
mtx <- matrix (nrow = n, ncol = m)
mtx[1:n,] <- replicate (m, rchisq (n = n, df = 3))
M <- colMeans (mtx)
meanM <- mean(M)
sdM <- sd(M)
# save them as "result"
result_toplot <- as.numeric (M)
# plot these means as a histogram
hist (result_toplot, probability = TRUE,
col = "cornflowerblue", border = "darkgreen",
main = paste("histogram of means ( n =", n, ", df = " ,3, ")"))
x <- seq (min(result_toplot), max(result_toplot),
length.out = 1000)
mean_result <- meanM
sd_result <- sqrt (6 / n)
lines (x, dnorm (x, mean = mean_result, sd = sd_result),
col = "red", lwd = 3)
}
Try it:
mean.dr (n = 500, m = 1000)
par (mfrow = c (2, 2))
mean.dr (n = 100, m = 1000)
mean.dr (n = 200, m = 1000)
mean.dr (n = 100, m = 2000)
mean.dr (n = 200, m = 2000)
par (mfrow = c (1, 1))
par (mfrow = c (2, 2))
ns <- c (100, 200)
ps <- c (1000, 2000)
for (n in ns) {
for (p in ps) {
mean.dr(n, p)
}
}
par (mfrow = c (1, 1))