\documentclass[11pt]{article}
\def\Z{\mathbb{Z}}
\def\Q{\mathbb{Q}}
\def\R{\mathbb{R}}
\def\C{\mathbb{C}}
\usepackage{latexsym,amsmath,amssymb}
\usepackage[all]{xy}
\setlength{\textwidth}{6.5in}
\setlength{\oddsidemargin}{0in}
\setlength{\textheight}{8.0in}
\begin{document}
\centerline{Mathematics 376 -- Mathematical Statistics}
\centerline{Computer Lab Day 1 in {\tt R}}
\centerline{February 6, 2012}
\vskip 10pt
\noindent
{\it Goals}
\vskip 10pt
The goals of today's lab are:
\begin{itemize}
\item To gain some more intuitive understanding of the $\chi^2,t,F$
distributions from \S 7.2 -- their densities, distributions, etc. and
computing with them in {\tt R}.
\item To solidify understanding of how those distributions relate
to sampling from normal populations.  
\end{itemize}
\vskip 10pt
\noindent
{\it Basic R Mechanics -- from Fall Lab Day 1}
\vskip 10pt
Refer to the assignment sheet for the first lab day
if you need to refresh your memory about running 
R in Haberlin 136, about how to save and print your
results, etc.  The Lab 1 assignment is posted on 
the course homepage, so you can open it and look
at the instructions there as needed.  
\vskip 10pt
\noindent
{\it Background on Probability Distributions in R -- from Fall Lab Day 2}
\vskip 10pt
The R package contains built-in functions for computing
probability functions for many of the discrete random 
variables we have discussed, plus functions for the probability
densities and cumulative distributions of many of the 
standard types of continuous random variables, including the 
$\chi^2, t, F$ types that we have just introduced.  Here's how it works.  (Also look at table
on 332 of Dalgaard \emph{Introductory Statistics with R}.) 
Each type of random 
variable is covered by a \emph{family} of 4 functions distinguished by a 
\emph{prefix letter} {\tt d, p, q}, or {\tt r}:
\vskip 10pt
\begin{itemize}
\item {\tt d} -- the probability function in the discrete case, or the 
density function continuous case,
\item {\tt p} -- the cumulative distribution
\item {\tt q} -- the quantile function (essentially the \emph{inverse
function} of the cumulative distribution, but of course the cumulative
distribution can fail to be injective, so this needs to be taken with
 a grain of salt)
\item {\tt r} -- random number generator 
\end{itemize}
Following the prefix letter comes the rest of the function 
name and the inputs needed to compute the corresponding values.
For instance the family of functions for $\chi^2$ random variables 
looks like this, where $df$ represents the number of degrees of
freedom:
\begin{itemize}
\item {\tt dchisq(y,d)} -- the probability density function for a
$\chi^2(d)$ random variable ($d = $ degrees of freedom):
$$f(y) = \begin{cases}
               \frac{y^{(d/2) - 1} e^{-y/2}}{\Gamma(d/2) 2^{d/2}} & \text{if } y > 0\\
                0 & \text{otherwise}
         \end{cases}$$
\item {\tt pchisq(y,df)} -- the cumulative distribution $P(Y \le y)$
\item {\tt qchisq(q,df)} -- computes the $y$ such that $P(Y \le y) = q$
\item {\tt rchisq(N,df)} -- generates a list of $N$ random numbers following
a $\chi^2(n)$ distribution.  
\end{itemize}
The $t$-distribution family of functions are called {\tt dt, pt, qt, rt}.  The 
$F$-distribution family of functions have the names {\tt df, pf, qf, rf}.
\vskip 10pt
\noindent
{\it Today's Lab}
\vskip 10pt 
\noindent
I.  Basic computations.
\begin{enumerate}
\item If $Y \sim \chi^2(12)$, find $P(Y > 10)$.  
\item Find a real value $a$ such that $P(Y > a) = .01$ if $Y \sim \chi^2(10)$.  
Check your answer using the tables in our text.
\item If $Y \sim t(23)$, find $P(-2.2 < Y < 2.2)$.  Compare with the corresponding
probability for the the standard normal $Z$ -- $P(-2.2 < Z < 2.2)$.  Explain
how this matches with our intuitive picture of the $t$-densities and their
dependence on $n$.
\item If $Y \sim F(10,14)$, find $P(5 < Y < 16)$.
\item If $Y \sim F(9,6)$, determine a real number $a$ such that $P(Y \le a) = .05$
\end{enumerate}
\vskip 10pt
\noindent
II. Overall characteristics of the $\chi^2$ family of densities.
\begin{enumerate}
\item Produce a plot showing the $\chi^2$ densities for $n = 4,6,8,10,12$
degrees of freedom together on the interval from $y = 0$ to $y = 30$.  
(Notes: The {\tt R} {\tt curve} command will be useful for this.  Recall
that you can superimpose plots using the {\tt add = TRUE} option.  The
plotting interval can also be controlled using the {\tt from = } and 
{\tt to = } options.)
\item Describe in a short paragraph how the value of $n = $ number of
degrees of freedom affects the shape of the $\chi^2(n)$ density graph.
(Look closely at the shape at $y = 0$ and at the location and height
of the maximum as $n$ changes.) 
\end{enumerate}
\vskip 10pt
\noindent
III.  Overall characteristics of the $F$ family of densities.
\begin{enumerate}
\item Produce a plot showing the $F$ densities for numerator
degrees of freedom $n = 2,3,4,5,6,7,8,9$ and denominator
degrees of freedom $d = 6$ together on the interval from $y = 0$ to $y = 6$.  
\item Produce a plot showing the $F$ densities for numerator
degrees of freedom $n = 6$ and denominator
degrees of freedom $d = 2,3,4,5,6,7,8,9$ together on the interval from $y = 0$ to $y = 6$.  
\end{enumerate}
\vskip 10pt
\noindent
IV.  A special property of the $F$-distributions.  The motivation for
this question is the fact (which you may have noted) that our book's
$F$-tables only include the percentage points $F_\alpha$ for ``small'' values
of $\alpha$ -- namely $\alpha = .005, .01, .025, .05, .1$.  What about
the complementary ``large'' values $\alpha = .995, .99, .975, .95, .9$?  
Why aren't they there?  This problem will show you the answer to those
questions and how to generate the $\alpha = .995, .99, .975, .95, .9$ percentage
points of the $F$-distributions if those are needed.
\begin{enumerate}
\item For instance, our book's tables show that $F_{.1}(10,8) = .254$ (see page 855).
Check this using {\tt R}.  (Note:  The book's table shows the value $F_{.1}$
such that $P(Y \ge F_{.1}) = .1$ (i.e. the .1 is the area in the \emph{upper tail}.
This means that $P(Y < F_{.1}) = .9$.)
\item  What happens if you enter \\
\begin{center}
{\tt qf(.9,10,8)*qf(.1,8,10)}?\\
\end{center} 
 What happens
if you enter $qf(\alpha,n,d)*qf(1 - \alpha,d,n)$ for any particular values of 
$\alpha$, $n$, and $d$?  Try as many examples
as you need to see the answer to the next part.
\item  If we know the percentage point $F_{\alpha}(n,d)$ (say for one of the
``small values'' given in the book's table), how can we find
the complementary percentage point $F_{1 - \alpha}(n,d)$?  State a 
general formula for this, and use it to compute $F_{.95}(15,4)$ using the 
information in the book's tables.  Check with {\tt R}.
\item (for after the lab)  Prove your formula in the previous part
using the definition of $F$-distributed random variables.
\end{enumerate}
\vskip 10pt
\noindent
V. The $T$-statistic.  Recall that if $Y_1, \ldots, Y_n$
are i.i.d. samples from a $N(\mu,\sigma^2)$ population, then 
the closest thing to a standard normal that we can compute
without knowledge of $\sigma$ is the $T$-statistic:
$$T = \frac{\sqrt{n}(\overline{Y} - \mu)}{S}.$$  
Here $S$ is the sample standard deviation:
$$S = \sqrt{\frac{1}{n-1} \sum_{i=1}^n (Y_i - \overline{Y})^2}.$$
In this problem, you will study the empirical distribution of the
$T$-statistics in a case where $n = 8$ is relatively small,
along the lines of something we did on Lab 3 from the fall.
\begin{enumerate}
\item Our main goal is to understand what happens if we generate \emph{lots}
of such collections of $n = 8$ samples, find their sample means $\overline{Y}$ 
and sample variances $S^2$, 
then consider \emph{the distribution of the $T$-statistics of the samples}
$$T = \frac{\sqrt{8}(\overline{Y} - \mu)}{S}.$$  
Do these have an approximately normal distribution, or is it different?
\item Here is one way to do this in {\tt R} using
a simple \emph{for loop}.  First we create space to store the values of 
$T$ computed from 1000 different random samples:
\vskip 10pt
\begin{center}
{\tt ts <- array(1:1000)}
\end{center}
Next, we generate 1000 different random samples of size $n = 8$ from a $N(10,25)$
distribution, find the $T$-statistic values of each of them, 
and store them in the array created above:

\begin{verbatim}
                for (i in 1:1000) 
                   { 
                   sample <- rnorm(8,10,25)
                   ts[i] <- sqrt(8)*(mean(sample) - 10)/sd(sample)
                   } 
\end{verbatim}

(Note: This could all be entered as one input line provided you leave
spaces in the appropriate places.  However, I think it is more readable 
(and it is easier to identify typos if you happen to make one) 
if you press ENTER at the end of each of the lines as
you are typing this.  If you do it that way, you should note that {\tt R} generates a
new input prompt {\tt +} each time, indicating that you 
are still in the body of the {\tt for} loop.  The final {\tt \}}
will close off the loop and take you back to the $>$ prompt.)
\item What does the distribution of the sample means look like?
We can see that by generating a density histogram for the data in the {\tt ts}
array, together with a standard normal density:
\vskip 10pt
\begin{center}
{\tt curve(dnorm(x,0,1),from=-5,to=5)}\\
{\tt hist(ts,breaks=20,freq=FALSE)}
\end{center}
What do you see from this display?  Is the standard normal density a
good match for the distribution of the $T$-statistic?  Is there another 
density that should fit better?  Add that to your plot.  
(Another Question:  Why did I put the first two plotting commands 
in this order and not the other way around?)
\end{enumerate}
\vskip 10pt
\noindent
{\it Assignment}
\vskip 10pt
Lab reports containing input, output, and answers to the question
posed above are due no later than Monday, February 13.  
If you do not finish during the hour today, you can return
to HA 136 any time it is not in use by another class.
 \end{document}
