An antiquated generator of the normal distribution is:
Generate
Set
The argument being that the CLT normality is sufficiently accurate with 12 terms
(a) Show that and
For any random variable it follows from the definition of expected value and variance of the uniform distribution that
and
. Since
are i.i.d. then we can do:
(b) Using histograms, compare this CLT-normal generator with the Box-Muller algorithm.
(c)Compare both generators with rnorm
Here, I generated samples from three different normal generators and compared them by estimated PDFs, histograms and Q-Q plots. The samples generated from CLT-normal generator lacked fitness to the tail of
while the Box-Muller transformation worked quite well compared with the samples generated from the R function
rnorm
The pictures are kinda crappy as you can see, so here’s the R code for them in case you’d like to reproduce them:
#ex 2.3
#Box-Muller's algorithm
Box.Muller<-function(n){
U1<-runif(n)
U2<-runif(n)
sqrt(-2*log(U1))*cos(2*pi*U2)
}
#CLT normal generator
CLT.Norm<-function(n){
matrix(runif(n*12, -.5, .5), ncol=12) %*% rep(1,12)
}
set.seed(100)
n<-1e5
rv.bm<-Box.Muller(n)
rv.clt<-CLT.Norm(n)
rv.rnorm<-rnorm(n)
#Histograms and Q-Q plots
par(mfrow=c(2,3))
hist(rv.rnorm, 50, main="Histogram of rnorm")
hist(rv.bm, 50, main="Histogram of Box-Muller generator")
hist(rv.clt, 50, main="Histogram of CLT generator")
qqplot(rv.bm, rv.clt, pch=".", main="Q-Q Plot: Box-Muller vs CLT",
xlab="Box-Muller", ylab="CLT")
abline(0,1,col="red")
qqplot(rv.rnorm, rv.bm, pch=".", main="Q-Q Plot: rnorm vs Box-Muller",
xlab="rnorm", ylab="Box-Muller")
abline(0,1,col="red")
qqplot(rv.rnorm, rv.clt, pch=".", main="Q-Q Plot: rnorm vs CLT",
xlab="rnorm", ylab="Box-Muller")
abline(0,1,col="red")
#overlay estimated densities
par(mfrow=c(1,1),lwd=2)
plot(density(rv.rnorm), main="Estimated PDFs", col = "red")
lines(density(rv.bm), col="green")
lines(density(rv.clt), col="blue")
legend(-4,.4,legend=c("rnorm", "Box-Muller", "CLT"), lty=1, lwd=2, col=c("red","green","blue"))