\documentclass[12pt]{amsart}
\usepackage{pdflscape}
\pagestyle{empty}
\thispagestyle{empty}
\theoremstyle{plain}    
\newtheorem*{definition}{Definition}
\newtheorem*{example}{Example}
\newtheorem*{exercise}{Exercise}
\newtheorem*{conjecture}{Conjecture}
\newtheorem*{theorem}{Theorem}
\newtheorem*{lemma}{Lemma}
\newtheorem*{remark}{Remark}
\begin{document} 
\title{Computational Assignment 1}
\section{computational assignment 1}
\subsection{Overview}
\par\vspace{1 cm}
We covered two methods of construction confidence intervals means of normal populations.  First, in situations 
where the population standard deviation is known, or we have a large enough sample that the sample standard deviation 
can be considered reliable (one rule of thumb: $n>30$; opinions vary a bit on this), the quantity
\[
\frac{\overline{x}-\mu}{\sigma}
\]
has a standard or $N(0,1)$ distribution which leads to a $100(1-\alpha)$ percent confidence interval for $\mu$ of:
\[
\overline{x}\pm z_{\alpha/2}\frac{\sigma}{\sqrt{n}}
\]
where $\sigma$ is the (known) population standard deviation, $n$ is the sample size, and $z_{\alpha/2}$ is the number for which 
\[
P(z\leq z_{\alpha/2})=1-\frac{\alpha}{2}
\] 
(this is called the \textit{large sample} or \textit{$\sigma$ known} situation)  
\par\vspace{1 cm}
If we have neither a large sample nor a reliable estimate of $\sigma$, it is better to use the fact that the quantity
\[
\frac{\overline{x}-\mu}{s}
\]
has a $t$ distribution with $n-1$ degrees of freedom, which leads to a $100(1-\alpha)$ percent confidence interval for $\mu$ of:
\[
\overline{x}\pm t_{\alpha/2}\frac{s}{\sqrt{n}}
\]
where $s$ is the sample standard deviation, $n$ is the sample size, and $t_{\alpha/2,n-1}$ is the number for which 
\[
P(t_{n-1}\leq t_{\alpha/2,n-1})=1-\frac{\alpha}{2}
\] 
(this is called the \textit{small sample} or \textit{$\sigma$ unknown} situation)  
\par\vspace{1 cm}
In this exercise, we will use R to generate simulated samples from a $N(\mu,\sigma)$ population with a specified mean $\mu$ and standard deviation $\sigma$.  Using these samples, we will construct both the "$\sigma$-known" and "$\sigma$-unknown" confidence intervals, and observe the proportion of time that the generated intervals contain the true population mean $\mu$.  Finally, we will construct confidence intervals (incorrectly) using the sample standard deviation $s$ in place of $\sigma$ in the "$\sigma$-known" variation, to see how far off the proportion of intervals that contain $\mu$ is from $1-\alpha$.  In the third case, the confidence intervals will be
\[
\overline{x}\pm z_{\alpha/2}\frac{s}{\sqrt{n}}
\]
Here we treat the sample standard deviation $s$ as if it were the true population standard deviation $\sigma$. 
\subsection{The R code}
One of the objectives of this assignment is to gain experience using R.  
\par\vspace{1 cm}
Like many programming languages, R gives you the ability to define your own functions.  This feature is tailor made for 
and exercise like this one and we define a function called "z\_versus\_t".  The function is written with the following 
arguments (values which are supplied at run time):
\par\vspace{1 cm}
\begin{itemize}
\item \texttt{mu}    is the population mean
\item \texttt{sigma} is the population standard deviation
\item \texttt{ssize} is the size of an IID sample from a $N(\mu,\sigma)$ population
\item \texttt{reps}  is the number of samples to simulate
\item \texttt{alpha} is the alpha level
\end{itemize}
\par\vspace{1 cm}
The first line of the code that defines the function is:
\par\vspace{0.5 cm}
\texttt{z\_versus\_t$<$-function(mu=0,sigma=1,ssize=10,alpha=0.05,reps=10000)} 
\par\vspace{1 cm}
If you invoke the function as 
\par\vspace{1 cm}
\texttt{z\_versus\_t()} 
\par\vspace{1 cm}
the parameter values used are the ones supplied in the function definition (i.e., these are the default values of none are supplied)
\par\vspace{1 cm}
You can override some or all of the defaults when you call the function:
\par\vspace{1 cm}
\texttt{z\_versus\_t(ssize=20,mu=4,sigma=5)}
\par\vspace{1 cm}
would generate samples of size 20 from a $N(4,5)$ distribution.  The \texttt{reps} and \texttt{alpha} parameters would take their default values of 10,000 and 0.05.
\par\vspace{1 cm}
Before you can run the \texttt{z\_versus\_t} function, you need to load it into your R session with the \texttt{source()} command. 
\par\vspace{1 cm}
If you are connected to the internet, you can load it directly from a URL with
\par\vspace{1 cm}
\texttt{source("http://www.sandgquinn.org/stonehill/MTH396/Spring2011/z\_versus\_t.r")}
\par\vspace{1 cm}
Alternatively, you can download the file z\_versus\_t.r from the website, copy it to the directory where R is running, 
and use
\par\vspace{1 cm}
\texttt{source("z\_versus\_t.r")}
\par\vspace{1 cm}
To check if the function z\_versus\_t is loaded, display the contents of the R workspace by typing
\par\vspace{1 cm}
\texttt{ls()}
\par\vspace{1 cm}
If the function is in the workspace, it should appear as
\par\vspace{1 cm}
\texttt{"z\_versus\_t"}
\par\vspace{1 cm}
You can remove it from the workspace by typing
\par\vspace{1 cm}
\texttt{rm(z\_versus\_t)}
\par\vspace{1 cm}
\subsection{The Assignment}
\par\vspace{1 cm}\noindent
\textbf{Question 1} Explain the role of the R function call \texttt{qnorm(1-alpha/2,1,0)} in the program 
\par\vspace{1 cm}\noindent
\textbf{Question 2} Explain the role of the R function call \texttt{qt(1-alpha/2,ssize-1)} in the program 
\par\vspace{1 cm}\noindent
\textbf{Question 3} Explain the role of the R function call \texttt{rnorm(ssize,mu,sigma)} in the program 
\par\vspace{1 cm}\noindent
\textbf{Question 4} Use the \texttt{z\_versus\_t} function to generate 10,000 samples of size 10 from a $N(4,5)$ population and compute confidence intervals.
\par\vspace{0.3 cm}
\begin{itemize}
\item How many of the confidence intervals constructed using the "$\sigma$-known" method contained $\mu$?
\item How many of the confidence intervals constructed using the "$\sigma$-unknown" method contained $\mu$?
\item How many of the confidence intervals constructed treating $s$ as $\sigma$" contained $\mu$?
\end{itemize}
\par\vspace{1 cm}\noindent
\textbf{Question 5} Repeat question 4 using sample sizes of 5, 20, 30, 40, and 50.  If you were designing a "rule of thumb" stating that it is OK to use $s$ in place of $\sigma$ for samples of $k$ or larger, what value of $k$ would 
you pick?
\par\vspace{0.3 cm}
\subsection{The R Code for z\_versus\_t}
\begin{landscape}
\begin{verbatim}
##########################################################################
# z_versus_t.r
#
# This simulation program compares confidence intervals computed two ways:
#   1) Using the large-sample or sigma-known approach (normal distributon)
#   2) Using the small-sample or sigma-unknown approach (t distribution)
#   3) Incorrectly using s in place of sigma in the large sample approach
#########################################################################
# Constants:
#  mu    is the population mean
#  sigma is the population standard deviation
#  ssize is the size of an IID sample from a Normal(mu,sigma) population
#  reps  is the number of samples to simulate
#  alpha is the alpha level
########################################################################
z_versus_t<-function(mu=0,sigma=1,ssize=10,alpha=0.05,reps=10000) {

cat("t versus z simulation: mu=",mu," sigma=",sigma," reps=",reps," alpha=",alpha," ssize=",ssize,"\n")
in_ci_z=0   #number of samples for which CI based on normal contains mu
in_ci_t=0   #number of samples for which CI based on t contains mu
in_ci_w=0   #number of samples for which CI using s in place of sigma contains mu
#
# Generate (reps) samples of size (ssize) from a normal(mu,sigma) population
#
for(i in c(1:reps)){
 x<-rnorm(ssize,mu,sigma)                         #generate sample of size ssize
 xbar<-mean(x)                                    #sample mean xbar
 s<-sd(x)                                         #sample standard deviation s
 LLz<-xbar-qnorm(1-alpha/2,0,1)*sigma/sqrt(ssize) #Upper limit of CI based on normal distribution, known sigma
 ULz<-xbar+qnorm(1-alpha/2,0,1)*sigma/sqrt(ssize) #Lower limit of CI based on normal distribution, known sigma
 if((mu > LLz) && (mu < ULz)) in_ci_z<-1+in_ci_z  #add 1 to counter if CI contains mu
 LLt<-xbar-qt(1-alpha/2,ssize-1)*s/sqrt(ssize)    #lower limit for CI based on t distribution with ssize-1 df
 ULt<-xbar+qt(1-alpha/2,ssize-1)*s/sqrt(ssize)    #upper limit for CI based on t distribution with ssize-1 df
 if((mu > LLt) && (mu < ULt)) in_ci_t<-1+in_ci_t  #add 1 to counter if CI contains mu
 LLw<-xbar-qnorm(1-alpha/2,0,1)*s/sqrt(ssize)     #Upper limit of CI based on normal distribution, known sigma
 ULw<-xbar+qnorm(1-alpha/2,0,1)*s/sqrt(ssize)     #Lower limit of CI based on normal distribution, known sigma
 if((mu > LLw) && (mu < ULw)) in_ci_w<-1+in_ci_w  #add 1 to counter if CI contains mu
 }
cat("Proportion of CIs containing mu - sigma known (normal distribution)",in_ci_z/reps,"\n")
cat("Proportion of CIs containing mu - sigma unknown (t distribution)",in_ci_t/reps,"\n")
cat("Proportion of CIs containing mu - treating s as known sigma",in_ci_w/reps,"\n")
}

\end{verbatim}
\end{landscape}
\end{document}



