# Monte Carlo inference about the mean
n = 1000
p = 2000
X = matrix(0, nrow=p, ncol=n)

for (i in 1:p) 
  X[i,] = pdf3(n,  mu=c(0, 6.5,14.5), 
               sigma=c(0.5,1.25,0.75), 
               p=c(0.2, 0.3, 0.5))

m = rowMeans(X)
m0 = mean(m)
s = sd(m)

alpha = 0.05
ci.lo = m0 - qnorm(1-alpha/2)*sd(X[1,])/sqrt(n)
ci.hi = m0 - qnorm(alpha/2)*sd(X[1,])/sqrt(n)