Tuesday, September 18, 2007

Example of Normal Inference


To illustrate inference about a normal mean, suppose we are interested in learning about the mean math ACT score for the students who are currently taking business calculus. I take a random sample of 20 students and put the values in the R vector y.

> y=sample(placedata$ACT,size=20)
> y
[1] 20 22 22 19 27 17 21 21 20 19 17 18 20 21 20 17 18 21 17 20


I compute some summary statistics.

> ybar=mean(y)
> S=sum((y-ybar)^2)

> n=length(y)


The definition of the log posterior of (mean, variance) with a noninformative prior is stored in the R function normchi2post. I construct a contour plot of this density by use of the mycontour function.

> mycontour(normchi2post, c(17,23,1.8,20),y)
> title(xlab="MEAN",ylab="VARIANCE")

To simulate draws of (mean, variance), I first simulate 1000 draws of the variance parameter from the scale times inverse chi-square density, and then simulate draws of the mean parameter. I plot the simulated draws on top of the contour density

> sigma2 = S/rchisq(1000, n - 1)
> mu = rnorm(1000, mean = ybar, sd = sqrt(sigma2)/sqrt(n))

> points(mu,sigma2)


Suppose we are interested in inferences about the 90th percentile of the population curve that is given by Q = mu + sigma z, where z is the 90th percentile of the standard normal. One can obtain a simulated sample from this posterior by simply computing Q on the simulated draws of (mu, variance).

> Q=mu+sqrt(sigma2)*qnorm(.90)

I draw the posterior density of Q by use of a density estimate on the simulated sample.

> plot(density(Q),main="POSTERIOR FOR 90TH PERCENTILE",lwd=3)

No comments: