sizesample = 100
nsample = 10000
xmin = 0
xmax = 6
A = matrix(runif(sizesample*nsample,xmin,xmax),sizesample,nsample)
# generates a (sizesample by nsample) matrix with samples from a uniform random variable between xmin and xmax

mu = (xmax+xmin)/2
sigma2 = 1/12 * (xmax-xmin)**2

y = (colSums(A) - sizesample*mu)/(sqrt(sizesample*sigma2))
myhist = hist(y,breaks=nsample/100,freq=FALSE,xlab="y-variable")
# plots the histogram of the vector y. "breaks" indicates the suggested number of bins, "freq=FALSE" normalizes the histogram so that its area is equal to 1 
curve(dnorm(x,0,1), col="red", add=TRUE)
# "dnorm(x,mu,sigma)" is the density of a normal random variable with mean mu e standard deviation sigma. "curve" plots the density

#in what follows we compare the empirical cumulative distribution of y and the Cumulative Distribution Function of the corresponding normal random variable
#miny=min(y)
#maxy=max(y)
#F=ecdf(y)
#plot(F)
#plot(function(x) pnorm(x),col="red", add=TRUE,miny,maxy)