\magnification=\magstep1
\centerline{Mathematics 376 -- Mathematical Statistics}
\centerline{Lab Project 3 and Problem Set 4 -- More on Confidence Intervals}
\centerline{February 27 and 29, 2012 -- {\it due:} Friday, March 2}
\bigskip
\noindent
{\it Background}
\bigskip
\noindent
We have now discussed constructing:
\bigskip
\item{$\bullet$} Large-sample confidence
intervals for means, differences of subpopulation means,
proportions and differences of subpopulation
proportions (using pivotal quantities that have normal distributions)
\item{$\bullet$}  Small-sample confidence intervals 
for means and differences of means under the assumption
that the population has a normal distribution 
(using pivotal quantities that have $t$ distributions)
\item{$\bullet$}  Confidence intervals for variances
or standard deviations (using pivotal quantities that
have $\chi^2$ distributions).
\item{$\bullet$}  Confidence intervals for a {\it ratio} of subpopulation
variances (using pivotal quantities that have $F$ distributions).
\bigskip
\noindent
In this week's lab we will begin by studying how statisticians might 
test to see whether the assumptions underlying a small-sample case 
are reasonable. Then we will apply some of this to an applied 
problem.  Finally, we will consider the problem of finding {\it shortest
confidence intervals} for standard deviations, given the confidence level
$\alpha$.  
\bigskip
\noindent
{\it Worked Example and Background}
\bigskip
\noindent
In deriving the small-sample formulas for confidence intervals for 
population means and differences of population means, recall that 
we had to make some rather restrictive additional hypotheses in  
order to apply the facts we know about  $t$-distributions: 
\medskip
\item{$\bullet$} we need to assume that the samples are drawn from  
{\it normal} populations 
\item{$\bullet$} for the difference of means case, we need to assume  
the two population {\it variances are equal}. 
\medskip
We have seen how to derive some evidence about whether the  
equality of population variances is reasonable, via our confidence 
intervals for  $\sigma_1^2/\sigma_2^2$ (derived using the $F$-distribution).   
But this still leaves some questions:   
\medskip
\item{$\bullet$}  Can we test whether it is reasonable to assume the samples 
are coming from a normal population? 
\item{$\bullet$}  How much of a difference does it make if they are not? 
\medskip
Thus, it becomes important to know that there are ways--some more 
intuitive, some more precise--to get some indication whether it is 
reasonable to make the assumption that a collection of sampled data
has been drawn from a normal population.

The first way we will discuss is a graphical method called the {\it normal
quantile-quantile plot}.  The idea here is to compare the distribution 
of the sampled data values with the distribution we would expect 
from a normal population in a particular way.  Doing this by hand, 
we would order the sample and compare the quantiles (i.e. order statistics)
for the data to the expected values of the quantiles for
a sample of the same size drawn from a normal population.  A scatter
plot of data points constructed from the two sets of quantiles is
generated to give a visual comparison.    To make life
easier, {\tt R} has a built-in command called {\tt qqnorm} that does
this computation and produces the plot.  To see how this works, 
we will work with the following (simulated) net income data (in units of 
$\$1000$ from a random
selection of $n = 9$ tax returns:
\medskip
\centerline{\tt incomes <- c(20.1,35.5,44.7,50.3,60.8,90.7,140.2,240.9,357.3)}
\centerline{\tt qqnorm(incomes)}
\medskip
What does this plot mean?  Well, for example, the sample maximum (the $Y_{(9)}$) in our
incomes is the $357.3$ (that is one of the {\it ``Sample Quantiles''}).  On the other hand, if we took a random 
sample $Z_1, \ldots, Z_9$ from a standard normal distribution, the 
expected value of the maximum there is $E(Z_{(9)}) \doteq 1.485$ (that is one of the 
{\it ``Theoretical Quantiles''}). (This value can 
be computed by numerically integrating the density function for $Z_{(9)}$ using
the formulas from \S 6.8.)  The point farthest to the right in the QQ Plot
is in fact $(1.485,357.3)$(!)

What are we looking for??   Well, if the data {\it are} coming from a normal
population, then we would expect the scatter plot to be  nearly a {\it straight line}.
Some ``random up and down variation'' from that line would also be expected
because of chance errors coming from the sampling process.

In this example, though, notice that there is an apparent systematic curvature
(concavity) in the scatter plot.   If you see this kind of pattern, it is a 
strong indication that assuming normality of the underlying population is 
probably not justified!!

For now, we will leave this discussion at this somewhat intuitive level.  We will
return to these questions later on, in our discussion of regression (finding a
line that ``best fits a collection of data'')  and the correlation coefficient that
we discussed briefly last semester in connection with covariances.   

For your general ``statistical literacy,'' it is also important to be aware that
there are a number of other much more sophisticated statistical normality
tests that are commonly used to get information in this situation.   {\tt R} has
a very good one (called the Shapiro-Wilk W Test).  Try entering
\bigskip
\centerline{\tt shapiro.test(incomes)}
\bigskip
We will discuss the ideas involved in this sort of test (especially the exact meaning 
of the {\tt p-value}) after Spring Break.  For now, what this is saying is that (with 
a $p$-value less than $.05$ or so), there 
is {\it actually pretty strong evidence to reject the assumption that 
these samples were coming from a normal population}!  On the other hand, {\it larger} $p$-values 
(say $.1$ or larger) could be the result of random variations in the sampling, so we don't 
necessarily want to reject the normality assumption in those cases.  
\bigskip
\noindent
{\it Lab Problems}
\bigskip 
\noindent
A)  We can test for normality now, but does it really matter if we
compute confidence intervals using the small sample formulas
when the normality assumption is violated?  In this problem
you will be working with a new family of distributions in {\tt R}, the lognormal
distributions.  These are defined in problem 4.128 in our book -- the idea is 
that if $Y$ is lognormal with parameters $\mu,\sigma$ then $X = \ln(Y)$ is normal
with the same $\mu$ and $\sigma$.  To get
a feeling for what they look like, plot
\medskip
\centerline{curve(dlnorm(x,0,1),from=0,to=10)}
\medskip
to see the density for the lognormal with $\mu = 0$ and $\sigma = 1$.  This is 
different from a normal density, but not that far off in a way. 
There is a family of functions in {\tt R} for the 
lognormals just like the other cases we have seen:
\medskip
\item{$\bullet$} {\tt dlnorm} gives the lognormal density
\item{$\bullet$} {\tt plnorm} gives the lognormal cumulative distribution
\item{$\bullet$} {\tt qlnorm} gives the lognormal quantiles
\item{$\bullet$} {\tt rlnorm} generates random samples from a lognormal

Following what 
we did in Lab 1, generate 1000 random samples of size 8 from a lognormal distribution with
$\mu = 0$ and $\sigma = 1$.  The expected value of this lognormal is $e^0 = 1$.  
Compute the $T$-statistic 
$$T = \sqrt{8} {\overline{Y} - 1\over S}$$
for each sample, and plot a histogram of the $T$ values.  Overlay a $t(7)$ density
curve.  Discuss your results.  (Comment:  The following statement comes from page 525
of Wackerly, Mendenhall, and Scheaffer -- ``... investigations [of
the empirical distributions of the $T$-statistic] have shown that moderate
departures from normality in the distribution of the population have
little effect on the distribution of the test statistic.''  In order not to be
misled by this, one would have to have a more precise notion of
``moderate departures from normality!'')
 

\bigskip
\noindent
B)  A lab technician tests blood samples from 20 men and computes
the total serum cholesterol level (LDL + HDL) for each.  The data obtained 
were as follows:
$$\matrix{199 & 272 & 261 & 248 & 235 & 192 & 203 & 278 & 268\cr
           230 & 242 & 305 & 286 & 310 & 345 & 289 & 326 & \cr
           335 & 297 & 328 & \cr}$$
(Ouch!!  Not a healthy-eating, regular-exercise group on the whole --
current guidelines say this level should be no higher than 200 in
adult males to minimize risk of heart disease and strokes.)
\bigskip
\item{1)}  Use {\tt R} to compute a $90\%$ confidence interval for
the population mean serum cholesterol level based on this data.
(Start by entering the data as a list:
\medskip
\centerline{\tt choldata <- c(199,272,261,248, etc. )}
\medskip
so that you can work on this with the {\tt R} functions for
mean, standard deviation, etc. 
You may want to refer back to the assignment sheet from Lab 2, 
although that used the large sample formulas.)
\item{2)}  
\bigskip
\noindent
In deriving our formulas for the $(1-\alpha)\times 100\%$
confidence interval for the variance, we used an interval where the 
``tails'' of the $\chi^2$-distributions had equal probabilities 
(areas) $\alpha/2$.  This gives the general form
$$\left[{(n-1)S^2\over \chi^2_{\alpha/2}(n-1)}, {(n-1)S^2\over \chi^2_{1-\alpha/2}(n-1)}\right]$$
for the $(1-\alpha)\times 100\%$ confidence interval,
where, as in the notation of the text book's $\chi^2$ tables,
if $Y$ has a $\chi^2$ distribution,
$$P(Y \ge \chi^2_\beta) = \beta.$$
Of course, this method also gives confidence intervals for 
standard deviations by taking square roots:
$$\left[\sqrt{{(n-1)S^2\over \chi^2_{\alpha/2}(n-1)}}, \sqrt{{(n-1)S^2\over \chi^2_{1-\alpha/2}(n-1)}}\right]$$

\item{4)}  Determine a $95\%$ confidence interval for 
the population standard deviation $\sigma$ for the cholesterol data.
\bigskip
\noindent
This method is something of a {\it compromise} -- it gives a 
definite procedure to get the confidence intervals, but it does not 
always give the {\it shortest possible confidence intervals}.   
\bigskip

\item{5)}  {\it Why} might it be of interest to get the shortest
possible confidence interval $[\widehat{\theta_L},\widehat{\theta_U}]$ where 
$P(\widehat{\theta_L} \le \sigma < \widehat{\theta_U}) = 1 - \alpha$?

\bigskip
\noindent
We will next see how to select $a,b$ (different from 
$b = \chi^2_{\alpha/2}$ and $a = \chi^2_{1-\alpha/2}$ as above)
in order to {\it minimize} the length of the interval
$$\left[\sqrt{{(n-1)S^2\over b}}, \sqrt{{(n-1)S^2\over a}}\right]$$
subject to the constraint that 
$$\int_a^b f(y)\ dy = 1 - \alpha,$$
where $f(y)$ is the $\chi^2(n-1)$ pdf.  
It turns out that finding the appropriate $a,b$ is equivalent to 
solving the two equations:
$$\eqalign{\int_a^b f(y)\ dy &= 1 - \alpha\cr
      a^{n/2}e^{-a/2} - b^{n/2}e^{-b/2} &=0 \cr}\leqno(2)$$
The following Maple procedure does this computation:
\bigskip
\leftline{\tt ShortestSigmaCI:=proc(alpha,n)}
\leftline{\tt local a,b,eq1,eq2,g,x,s;}
\leftline{\tt g:=x\^{}((n-1)/2-1)*exp(-x/2)/(2\^{}((n-1)/2)*GAMMA((n-1)/2));}
\leftline{\tt eq1:=int(g,x=a..b)=1-alpha;}
\leftline{\tt eq2:=evalf(a\^{}(n/2)*exp(-a/2)=b\^{}(n/2)*exp(-b/2));}
\leftline{\tt s:=fsolve($\{$eq1,eq2$\}$,$\{$a,b$\}$,$\{$a=0..n-1,b=n-1..infinity$\}$);}
\leftline{\tt return s;}
\leftline{\tt end:}

\bigskip
\item{6)}  Enter this procedure into your worksheet and 
use it to determine the $a,b$ for the shortest confidence
$95\%$ intervals for the standard deviation with $n = 3,4,5,6,7,8,9,10$.
Compare with the ``compromise'' intervals for this $\alpha$ and $n$.
How much are we gaining?  (Note:  $n$ is the number of samples here -- 
the corresponding $\chi^2$ distribution has $n - 1$ degrees of freedom.)
\item{7)}  (For after the lab)  Show that the solution of 
(2) above does give the minimum length interval. (Hint:
This is a constrained optimization problem.  What method 
from MATH 241 applies to these?)
\bigskip
\noindent
B)  (This one will require some thought before you 
start plugging numbers into Maple!)  A group of pediatric nurses
were interested in the effect of prenatal care on the birthweight
of babies.  Mothers were divided into two groups, and their
babies' weights were compared.  The mothers in group 1
had 5 or fewer prenatal care sessions, and their
babies had birthweights:
$$49,108,110,82,93,114,134,114,96,52,101,114,120,116$$
(all in ounces).
The mothers in group 2
had 6 or more prenatal care sessions, and their
babies had birthweights:
$$133,108,93,119,119,98,106,87,153,116,129,97,110,131$$
(also in ounces).  One question they were trying to address
with this study was:  
\medskip
\noindent
{\it Was there more variation in the birthweights of the babies whose 
mothers had fewer prenatal care visits  than in 
the birthweights of the babies whose mothers had more visits?}
\bigskip
\item{1)}  Explain why $\theta = \sigma_1^2/\sigma_2^2$ is an
appropriate target parameter to consider, and 
why the statistic $\widehat{\theta} = S_1^2/S_2^2$ is a reasonable
estimator ($S_i^2 = $ sample variance for group $i = 1,2$).  
\item{2)} What is the 
value of this statistic here, and what does it suggest?
\item{3)}  But now, the next question is:  How reliable is
that conclusion?  How likely is it that we could 
get a different conclusion for different random samples from the 
same distributions?  To answer this we can try to construct
a confidence interval for the ratio $\theta = \sigma_1^2/\sigma_2^2$.
To do this, explain why $(S_1^2/\sigma_1^2)/(S_2^2/\sigma_2^2)$
is an appropriate ``pivotal quantity'' to use to derive
confidence intervals for $\sigma_1^2/\sigma_2^2$.  Recall this means that:
\itemitem{$\bullet$}  $\sigma_1^2/\sigma_2^2$ must be the 
only unknown part and 
\itemitem{$\bullet$} the distribution of the pivotal quantity 
must be {\it known}.  
\item{3)}  Explain how to get a $(1-\alpha)\times 100\%$ confidence
interval for $\sigma_1^2/\sigma_2^2$, and use your method
to find a $95\%$ confidence interval for $\sigma_1^2/\sigma_2^2$
using the data above.  What does this suggest about your 
conclusion from part 2?
\bigskip
\noindent
{\it Problems from the text:}  
\bigskip
\item{$\bullet$} Chapter 8/88,89,96,100;
\bye

\bigskip
\noindent
{\it Additional Problem}
\bigskip
(A follow-up problem to part 3 of Lab Problem B above.)  
To determine the endpoints of the 
confidence interval for the ratio of variances $\sigma_1^2/\sigma_2^2$,
you needed to work with a formula like this 
$$P\left({S_1^2/S_2^2\over f_{\alpha/2}} \le \sigma_1^2/\sigma_2^2
\le {S_1^2/S_2^2\over f_{1-\alpha/2}}\right) = 1 - \alpha.$$
In the problem from the Lab Project, you were looking for a 
$95\%$ confidence interval, so you needed $f_{.025}, f_{.975}$.  The $F$-table
in the back of the book gives the values $f_{\alpha}$ for
$\alpha = .1, .05, .025, .01, .005$, but not the complementary 
values $\alpha = .9, .95, .975, .99, .995$.  As you might have
wondered, {\it why aren't those other entries there??}
Fortunately, with the {\tt FCDF} function from our Maple 
package, you could compute $f_{.975}$, etc., by the same method
we used for $\chi^2$ and $t$ distributions in the lab, 
so this wasn't a problem.  In this question you will show that 
the other values could be derived from the table entries given with 
a little cleverness, too!
\bigskip
\item{1)} Let the numbers $f_\alpha(n,m)$ be defined by 
$$P(Y \ge f_\alpha(n,m)) = \alpha$$
if $Y$ has an $F$ distribution with $n$ numerator degrees of freedom
and $m$ denominator degrees of freedom.
Show that for all $0 < \alpha < 1$, these numbers satisfy:
$$f_{1-\alpha}(n,m) = {1\over f_{\alpha}(m,n)}$$
(note the reversal in the numbers of degrees of freedom on the 
right.)  {\bf Hint:}  What is the standard way to get a random variable
that has an $F$-distribution?
\item{2)}  Using the table in the text and part 1, determine the 
number $f_{.99}(7,10)$.
\bye
