%% This document created by Scientific Word (R) Version 2.5
%% Starting shell: article

\documentclass[12pt,thmsa]{article}%
\usepackage{amssymb}
\usepackage{sw20lart}%
\usepackage{amsmath}%
\setcounter{MaxMatrixCols}{30}%
\usepackage{amsfonts}%
\usepackage{graphicx}
%TCIDATA{OutputFilter=latex2.dll}
%TCIDATA{Version=5.50.0.2890}
%TCIDATA{TCIstyle=article/art4.lat,lart,article}
%TCIDATA{Created=Thu Oct 15 23:57:36 1998}
%TCIDATA{LastRevised=Thursday, February 26, 2009 22:23:59}
%TCIDATA{<META NAME="GraphicsSave" CONTENT="32">}
%TCIDATA{<META NAME="SaveForMode" CONTENT="1">}
%TCIDATA{BibliographyScheme=Manual}
%BeginMSIPreambleData
\providecommand{\U}[1]{\protect\rule{.1in}{.1in}}
%EndMSIPreambleData
\setcounter{section}{4}
\setcounter{page}{93}
\setlength{\oddsidemargin}{0mm}
\setlength{\evensidemargin}{0mm}
\setlength{\topmargin}{-2mm}
\setlength{\headheight}{7mm}
\setlength{\headsep}{4.5mm}
\setlength{\textheight}{220mm}
\setlength{\textwidth}{158mm}
\begin{document}
\section{Linear regression}

\noindent\noindent Regression is all about relationships between variables. In
its simplest form, it constitutes a technique for modelling a relationship
between two variables.%
%TCIMACRO{\FRAME{dhFU}{1.9259in}{2.5175in}{0pt}{\Qcb{Sir Francis Galton 1822 -
%1911}}{}{galton.gif}{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 1.9259in;  height 2.5175in;  depth 0pt;  original-width 6.1246in;
%original-height 8.0315in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'Galton.gif';file-properties "XNPEU";}}}%
%BeginExpansion
\begin{center}
\includegraphics[
natheight=8.031500in,
natwidth=6.124600in,
height=2.5175in,
width=1.9259in
]%
{Galton.gif}%
\\
Sir Francis Galton 1822 - 1911
\end{center}
%EndExpansion
It was the pioneering work of Sir Francis Galton in the 1880s that gave rise
to the technique, the original idea being the direct result of an experiment
on sweet peas. He noticed that the seeds of the progeny of parents with seeds
heavier than average were also heavier than average, but the difference was
not as pronounced; the same effect was true for the seeds of the progeny of
parents with light seeds, where again the differences from the average were
not as great. He called this phenomenon \textit{reversion} and wrote that the
mean weight "reverted, or regressed, toward mediocrity".

\subsection{Data and Questions}

\textbf{Data set 5.1}\quad\emph{Divorces in England and Wales}\medskip\newline
The table below gives the annual number of divorces recorded in England and
Wales between 1975 and 1980 (\textit{Marriage and Divorce Statistics}, Office
of Population Censuses and Surveys, HMSO).\medskip

\begin{center}
$%
\begin{tabular}
[c]{lcccccc}%
\multicolumn{7}{l}{\textbf{Table 5.1}\quad Divorces in England and
Wales}\\\hline
Year & 1975 & 1976 & 1977 & 1978 & 1979 & 1980\\
Divorces (thousands) & 120.5 & 126.7 & 129.1 & 143.7 & 138.7 & 148.3\\\hline
\end{tabular}
$\medskip
\end{center}

\noindent They are plotted below.\medskip\
%TCIMACRO{\FRAME{dhFU}{3.3728in}{2.5988in}{0pt}{\Qcb{\QTR{bf}{Figure 5.1}\quad
%Divorces in England and Wales}}{}{a5regress_f1.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.3728in;  height 2.5988in;  depth 0pt;  original-width 4.1805in;
%original-height 3.2154in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'a5regress_f1.eps';file-properties "XNPEU";}}}%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.5988in,
width=3.3728in
]%
{a5regress_f1.eps}%
\\
\textbf{Figure 5.1}\quad Divorces in England and Wales
\end{center}
%EndExpansion
\medskip The plot suggests a roughly linear trend with some random
scatter.\newline$\square\bigskip$\newline\textbf{Data set 5.2}\quad
\emph{Olympic Sprint Times\medskip\newline}Table 5.2 gives the times in
seconds recorded by the winners in the finals of the men's sprint events (100,
200, 400, 800 and 1500 metres) at each of the 21 Olympic Games from 1900 to
2004 along with the heights above sea level of the different venues. Obviously
the years 1916, 1940, 1944 are missing since the Olympic games were not held
during the World Wars.$\medskip$\newline$\hspace{1.5in}%
\begin{tabular}
[c]{llllllll}%
\multicolumn{8}{l}{\textbf{Table 5.2}\quad Winning times for Olympic running
events}\\\hline
\textit{Year} & \textit{Venue} & \textit{100m} & \textit{200m} & \textit{400m}
& \textit{800m} & \textit{1500m} & \textit{Altitude}\\\hline
\multicolumn{1}{c}{1900} & Paris & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.20} & \multicolumn{1}{c}{49.40} &
\multicolumn{1}{r}{121.40} & 246.00 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1904} & St Louis & \multicolumn{1}{c}{11.00} &
\multicolumn{1}{c}{21.60} & \multicolumn{1}{c}{49.20} &
\multicolumn{1}{r}{116.00} & 245.40 & \multicolumn{1}{r}{455}\\
\multicolumn{1}{c}{1908} & London & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.40} & \multicolumn{1}{c}{50.00} &
\multicolumn{1}{r}{112.80} & 243.40 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1912} & Stockholm & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{21.70} & \multicolumn{1}{c}{48.20} &
\multicolumn{1}{r}{111.90} & 236.80 & \multicolumn{1}{r}{46}\\
\multicolumn{1}{c}{1920} & Antwerp & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.00} & \multicolumn{1}{c}{49.60} &
\multicolumn{1}{r}{113.40} & 241.80 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1924} & Paris & \multicolumn{1}{c}{10.60} &
\multicolumn{1}{c}{21.60} & \multicolumn{1}{c}{47.60} &
\multicolumn{1}{r}{112.40} & 233.60 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1928} & Amsterdam & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{21.80} & \multicolumn{1}{c}{47.80} &
\multicolumn{1}{r}{111.80} & 233.20 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1932} & Los Angeles & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{21.20} & \multicolumn{1}{c}{46.20} &
\multicolumn{1}{r}{109.80} & 231.20 & \multicolumn{1}{r}{340}\\
\multicolumn{1}{c}{1936} & Berlin & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{20.70} & \multicolumn{1}{c}{46.50} &
\multicolumn{1}{r}{112.90} & 227.80 & \multicolumn{1}{r}{115}\\
\multicolumn{1}{c}{1948} & London & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{21.10} & \multicolumn{1}{c}{46.20} &
\multicolumn{1}{r}{109.20} & 225.20 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1952} & Helsinki & \multicolumn{1}{c}{10.40} &
\multicolumn{1}{c}{20.70} & \multicolumn{1}{c}{45.90} &
\multicolumn{1}{r}{109.20} & 225.20 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1956} & Melbourne & \multicolumn{1}{c}{10.50} &
\multicolumn{1}{c}{20.60} & \multicolumn{1}{c}{46.70} &
\multicolumn{1}{r}{107.70} & 221.20 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1960} & Rome & \multicolumn{1}{c}{10.20} &
\multicolumn{1}{c}{20.50} & \multicolumn{1}{c}{44.90} &
\multicolumn{1}{r}{106.30} & 215.60 & \multicolumn{1}{r}{66}\\
\multicolumn{1}{c}{1964} & Tokyo & \multicolumn{1}{c}{10.00} &
\multicolumn{1}{c}{20.30} & \multicolumn{1}{c}{45.10} &
\multicolumn{1}{r}{105.10} & 218.10 & \multicolumn{1}{r}{45}\\
\multicolumn{1}{c}{1968} & Mexico City & \multicolumn{1}{c}{9.95} &
\multicolumn{1}{c}{19.83} & \multicolumn{1}{c}{43.80} &
\multicolumn{1}{r}{104.30} & 214.90 & \multicolumn{1}{r}{7349}\\
\multicolumn{1}{c}{1972} & Munich & \multicolumn{1}{c}{10.14} &
\multicolumn{1}{c}{20.00} & \multicolumn{1}{c}{44.66} &
\multicolumn{1}{r}{105.90} & 216.30 & \multicolumn{1}{r}{1699}\\
\multicolumn{1}{c}{1976} & Montreal & \multicolumn{1}{c}{10.06} &
\multicolumn{1}{c}{20.23} & \multicolumn{1}{c}{44.26} &
\multicolumn{1}{r}{103.50} & 219.20 & \multicolumn{1}{r}{104}\\
\multicolumn{1}{c}{1980} & Moscow & \multicolumn{1}{c}{10.25} &
\multicolumn{1}{c}{20.19} & \multicolumn{1}{c}{44.60} &
\multicolumn{1}{r}{105.40} & 218.40 & \multicolumn{1}{r}{497}\\
\multicolumn{1}{c}{1984} & Los Angeles & \multicolumn{1}{c}{9.99} &
\multicolumn{1}{c}{19.80} & \multicolumn{1}{c}{44.27} &
\multicolumn{1}{r}{103.00} & 212.50 & \multicolumn{1}{r}{340}\\
\multicolumn{1}{c}{1988} & Seoul & \multicolumn{1}{c}{9.92} &
\multicolumn{1}{c}{19.75} & \multicolumn{1}{c}{43.87} &
\multicolumn{1}{r}{103.45} & 215.96 & \multicolumn{1}{r}{111}\\
\multicolumn{1}{c}{1992} & Barcelona & \multicolumn{1}{c}{9.96} &
\multicolumn{1}{c}{20.01} & \multicolumn{1}{c}{43.50} &
\multicolumn{1}{r}{103.66} & 220.12 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1996} & Atlanta & \multicolumn{1}{c}{9.84} &
\multicolumn{1}{c}{19.32} & \multicolumn{1}{c}{43.49} &
\multicolumn{1}{r}{102.58} & 215.78 & \multicolumn{1}{r}{1026}\\
\multicolumn{1}{c}{2000} & Sydney & \multicolumn{1}{c}{9.87} &
\multicolumn{1}{c}{20.09} & \multicolumn{1}{c}{43.84} &
\multicolumn{1}{r}{105.08} & 212.07 & \multicolumn{1}{r}{3}\\\hline
2004 & Athens &  & 19.79 & 44.00 & 104.45 & 214.20 & \multicolumn{1}{r}{505}%
\\\hline
\end{tabular}
\medskip$\newline We would like to construct a predictive model for the
winning times at the 2008 games in Beijing (and the 2012 games in
London).\ Clearly there is now a possibility of more than one explanatory
variable.$\newline\square$

\subsection{Fitting the model}

\subsubsection{\noindent Notation}

The equation of a regression line takes the form
\[
y=\alpha+\beta x,
\]
and gives rise to the regression model%
\begin{equation}
Y_{i}=\alpha+\beta x_{i}+\varepsilon_{i},\quad i=1,2,\ldots,n,\label{line}%
\end{equation}
where $Y_{1},Y_{2},\ldots,Y_{n}$ are observable random variables, the values
$x_{1},x_{2},\ldots,x_{n}$ are specified and $\varepsilon_{1},\varepsilon
_{2},\ldots,\varepsilon_{n}$ non-observable random variables. The data are,
therefore, in the form of $n$ ordered pairs $\left(  x_{i},Y_{i}\right)  $,
$i=1,\ldots,n$.$\medskip$\newline$x_{i}$ is called an \emph{explanatory
variable} or, sometimes, a \emph{predictor};\newline$Y_{i}$ is usually
referred to as the \emph{response};\newline$\varepsilon_{i}$ is called\ a
\emph{departure} or \emph{residual}.\medskip\newline Each $\varepsilon_{i}$
has distribution $N(0,\sigma^{2})$ and they are pairwise uncorrelated.
Clearly
\[
E\left[  Y_{i}\right]  =\alpha+\beta x_{i},\quad V\left[  Y_{i}\right]
=\sigma^{2},\quad i=1,2,\ldots,n.
\]
We are going to be concerned with linear model in its more general form
involving several explanatory variables. The most convenient way of doing that
is to write down the model in matrix notation. Equation $\left(
\ref{line}\right)  $ is a set of simultaneous equations which can be written
more concisely as
\begin{equation}
\mathbf{Y}=\mathbf{X\theta}+\mathbf{\varepsilon}\label{lmod}%
\end{equation}
where
\[
\mathbf{Y}=\left(
\begin{array}
[c]{c}%
Y_{1}\\
\vdots\\
Y_{n}%
\end{array}
\right)  ,\quad\mathbf{X}=\left(
\begin{array}
[c]{cc}%
1 & x_{1}\\
\vdots & \vdots\\
1 & x_{n}%
\end{array}
\right)  ,\quad\mathbf{\theta}=\left(
\begin{array}
[c]{l}%
\alpha\\
\beta
\end{array}
\right)  ,\quad\mathbf{\varepsilon}=\left(
\begin{array}
[c]{c}%
\varepsilon_{1}\\
\vdots\\
\varepsilon_{n}%
\end{array}
\right)  .
\]
\textbf{Data set 5.3}\quad\emph{AIDS data for the USA\medskip}\newline These
data are for AIDS incidence in the USA, adjusted for reporting delays. The
data are taken from Rosenberg, P.S. and Gail, M.H. (1991): Backcalculation of
flexible linear models of the Human Immunodeficiency Virus infection curve.
\textit{Applied Statistics}, \textbf{40}, 269-282.\medskip\newline Newly
reported cases are recorded quarterly and the variable \textit{Time} therefore
counts 3-monthly periods, staring with the first quarter in 1982.\medskip

\begin{center}%
\begin{tabular}
[c]{rrrrrrrrrrr}%
\multicolumn{7}{l}{\textbf{Table 5.3}\quad Reported AIDS cases in the USA} &
&  &  & \\\hline
\textit{Quarter} & \textit{Time} & \textit{Cases} &  & \textit{Quarter} &
\textit{Time} & \textit{Cases} &  & \textit{Quarter} & \textit{Time} &
\textit{Cases}\\\cline{1-3}\cline{5-11}%
1982:1 & 1 & 185 &  & 1984:2 & 10 & 1369 &  & 1986:2 & 18 & 4321\\
1982:2 & 2 & 200 &  & 1984:3 & 11 & 1563 &  & 1986:3 & 19 & 4863\\
1982:3 & 3 & 293 &  & 1984:4 & 12 & 1726 &  & 1986:4 & 20 & 5192\\
1982:4 & 4 & 374 &  & 1985:1 & 13 & 2142 &  & 1987:1 & 21 & 6155\\
1983:1 & 5 & 554 &  & 1985:2 & 14 & 2525 &  & 1987:2 & 22 & 6816\\
1983:2 & 6 & 713 &  & 1985:3 & 15 & 2951 &  & 1987:3 & 23 & 7491\\
1983:3 & 7 & 763 &  & 1985:4 & 16 & 3160 &  & 1987:4 & 24 & 7726\\
1983:4 & 8 & 857 &  & 1986:1 & 17 & 3819 &  & 1988:1 & 25 & 8483\\
1984:1 & 9 & 1147 &  &  &  &  &  &  &  & \\\hline
\end{tabular}
\medskip
\end{center}

\noindent The scatterplot shows that the trend is not linear.
%TCIMACRO{\FRAME{dhFU}{3.5034in}{2.706in}{0pt}{\Qcb{\QTR{bf}{Figure 5.2}\quad
%Incidence of AIDS cases in the USA against time}}{}{a5regress_f2.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.5034in;  height 2.706in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'a5regress_f2.EPS';file-properties "XNPEU";}}}%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.706in,
width=3.5034in
]%
{a5regress_f2.eps}%
\\
\textbf{Figure 5.2}\quad Incidence of AIDS cases in the USA against time
\end{center}
%EndExpansion
The plot has all the appearance of showing a functional relationship between
AIDS incidence and time. One might, for example, try to fit an exponential
function or some kind of power law to model the growth curve. For
example\newline%
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}x_{i}^{2}+\varepsilon_{i},\quad
i=1,2,\ldots,n
\]
$\square\bigskip$\newline At this point you should be careful to note that in
general, a model is \textbf{not} called a linear model because of fitting
linear functions of $x$. \noindent It is called a linear model because it is
\emph{linear in the unknown parameters} $\mathbf{\theta}$. \noindent Thus
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}x_{i}^{2}+\varepsilon_{i},\quad
i=1,2,\ldots,n
\]
still classes as a linear model because it is linear in $\alpha,\beta
_{1},\beta_{2}$. Furthermore, it is clear that all such linear models can be
written in the form of Equation $\left(  \ref{lmod}\right)  $ by suitably
defining $\mathbf{X}$ and $\mathbf{\theta}$. In this notation%
\[
E\left(  \mathbf{Y}\right)  =\mathbf{X\theta},\quad V\left(  \mathbf{Y}%
\right)  =\sigma^{2}\mathbf{I},
\]
where $\mathbf{I}$ is the identity.\noindent\bigskip\newline\textbf{Example
5.1\medskip\newline}The model above,
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}x_{i}^{2}+\varepsilon_{i},\quad
i=1,2,\ldots,n
\]
can be written in the form of Equations $\left(  \ref{lmod}\right)  $ by
defining%
\[
\mathbf{X}=\left(
\begin{array}
[c]{ccc}%
1 & x_{1} & x_{1}^{2}\\
\vdots & \vdots & \vdots\\
1 & x_{n} & x_{n}^{2}%
\end{array}
\right)  ,\quad\mathbf{\theta}=\left(
\begin{array}
[c]{l}%
\alpha\\
\beta_{1}\\
\beta_{2}%
\end{array}
\right)  .
\]
$\blacksquare$

\subsubsection{The likelihood function}

\noindent The likelihood function has the simple form
\[
L\left(  \alpha,\beta,\sigma^{2};\mathbf{y}\right)  =\left(  2\pi\sigma
^{2}\right)  ^{-n\left/  2\right.  }\exp\left[  -\frac{1}{2\sigma^{2}}%
\sum_{i=1}^{n}\varepsilon_{i}^{2}\right]  =\left(  2\pi\sigma^{2}\right)
^{-n\left/  2\right.  }\exp\left[  -\frac{1}{2\sigma^{2}}\mathbf{\varepsilon
}^{T}\mathbf{\varepsilon}\right]  .
\]


\noindent Maximisation of $L$ corresponds to minimisation of $\sum_{i=1}%
^{n}\varepsilon_{i}^{2}=\mathbf{\varepsilon}^{T}\mathbf{\varepsilon}$.
\noindent For this reason, the method of optimisation is also known as the
\emph{method of least squares}.\medskip

\noindent$\mathbf{\varepsilon}^{T}\mathbf{\varepsilon=}\left(  \mathbf{y}%
-\mathbf{X\theta}\right)  ^{T}\left(  \mathbf{y}-\mathbf{X\theta}\right)  $
can be minimised whatever the underlying distribution; it is only optimal in
the sense of maximum likelihood when the underlying distribution is
normal.\bigskip\medskip\newline We need to minimise $\left(  \mathbf{y}%
-\mathbf{X\theta}\right)  ^{T}\left(  \mathbf{y}-\mathbf{X\theta}\right)  $
with respect to $\mathbf{\theta}$. Differentiating with respect to
$\mathbf{\theta}$ and equating to $\mathbf{0}$ gives%
\[
2\mathbf{X}^{T}\left(  \mathbf{y}-\mathbf{X\theta}\right)  =\mathbf{0,}
\]
so the maximum likelihood estimator $\widehat{\mathbf{\theta}}$ of
$\mathbf{\theta}$ is given by%
\begin{equation}
\mathbf{X}^{T}\mathbf{X}\widehat{\mathbf{\theta}}=\mathbf{X}^{T}%
\mathbf{Y.\label{estequ}}%
\end{equation}
Provided $\mathbf{X}^{T}\mathbf{X}$ is non-singular, this can be written%
\[
\widehat{\mathbf{\theta}}=\left(  \mathbf{X}^{T}\mathbf{X}\right)
^{-1}\mathbf{X}^{T}\mathbf{Y.}
\]
It must be a minimum in $\left(  \mathbf{y}-\mathbf{X\theta}\right)
^{T}\left(  \mathbf{y}-\mathbf{X\theta}\right)  $, and hence a maximum in the
likelihood, because $\mathbf{\theta}$ is unconstrained and therefore $\left(
\mathbf{y}-\mathbf{X\theta}\right)  ^{T}\left(  \mathbf{y}-\mathbf{X\theta
}\right)  $ could be arbitrarily large.\bigskip

\noindent Equation $\left(  \ref{estequ}\right)  $ is a set of linear
simultaneous equations called the \emph{normal equations} for the linear
model.\medskip\newline Note that, before you can go ahead and estimate these
coefficients, you need to be extra careful about validating the assumptions
which underlie the method of fitting used. These are

\begin{enumerate}
\item[$\bullet$] the variance of the response is constant (i.e. the variance
does not depend upon $x$);

\item[$\bullet$] the distribution of the response is normal;

\item[$\bullet$] the $\varepsilon_{i}$ are independent;

\item[$\bullet$] the model you are trying to fit is a valid one.\smallskip
\end{enumerate}

\noindent The first three of these assumptions are encapsulated by the
statement $\mathbf{\varepsilon}\sim N\left(  0,\sigma^{2}\mathbf{I}\right)  $.
Under such assumptions the method of least squares and the method of maximum
likelihood are equivalent.

\subsubsection{Simple linear regression}

\noindent The matrix $\mathbf{X}^{T}\mathbf{X}$ is non-singular, so that
\[
\mathbf{X}=\left(
\begin{array}
[c]{cc}%
1 & x_{1}\\
\vdots & \vdots\\
1 & x_{n}%
\end{array}
\right)  ,\quad\mathbf{\theta}=\left(
\begin{array}
[c]{l}%
\alpha\\
\beta
\end{array}
\right)
\]%
\[
\mathbf{X}^{T}\mathbf{X=}\left(
\begin{array}
[c]{ll}%
n & \sum x_{i}\\
\sum x_{i} & \sum x_{i}^{2}%
\end{array}
\right)  ,\quad\left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{=}%
\frac{\left(
\begin{array}
[c]{ll}%
\sum x_{i}^{2} & -\sum x_{i}\\
-\sum x_{i} & n
\end{array}
\right)  }{n\sum x_{i}^{2}-\left(  \sum x_{i}\right)  ^{2}}
\]
and
\[
\mathbf{X}^{T}\mathbf{Y=}\left(
\begin{array}
[c]{l}%
\sum Y_{i}\\
\sum x_{i}Y_{i}%
\end{array}
\right)
\]
so that
\[
\widehat{\mathbf{\theta}}=\left(
\begin{array}
[c]{l}%
\widehat{\alpha}\\
\widehat{\beta}%
\end{array}
\right)  =\frac{1}{n\sum x_{i}^{2}-\left(  \sum x_{i}\right)  ^{2}}\left(
\begin{array}
[c]{l}%
\sum x_{i}^{2}\sum Y_{i}-\sum x_{i}\sum x_{i}Y_{i}\\
-\sum x_{i}\sum Y_{i}+n\sum x_{i}Y_{i}%
\end{array}
\right)  .
\]
\noindent It is more usual to write $n\sum x_{i}^{2}-\left(  \sum
x_{i}\right)  ^{2}$ as $n\sum\left(  x_{i}-\overline{x}\right)  ^{2}$ and
$-\sum x_{i}\sum Y_{i}+n\sum x_{i}Y_{i}$ as $n\sum\left(  x_{i}-\overline
{x}\right)  Y_{i}$ . The expressions for $\widehat{\alpha}$ and $\widehat
{\beta}$ then simplify to
\[%
\begin{array}
[c]{lcl}%
\widehat{\beta} & = & \dfrac{\sum\left(  x_{i}-\overline{x}\right)  Y_{i}%
}{\sum\left(  x_{i}-\overline{x}\right)  ^{2}},\\
&  & \\
\widehat{\alpha} & = & \overline{Y}-\widehat{\beta}\overline{x}.
\end{array}
\]


\noindent\textbf{Example 5.3 }(Example 5.1 revisited)\quad
\emph{Divorces\medskip}\medskip\newline Using the data given in Table 5.1,%
\[%
\begin{tabular}
[c]{lcccccc}\hline
$x_{i}$ & $1975$ & $1976$ & $1977$ & $1978$ & $1979$ & $1980$\\
$y_{i}$ & \multicolumn{1}{r}{$120.50$} & \multicolumn{1}{r}{$126.70$} &
\multicolumn{1}{r}{$129.10$} & \multicolumn{1}{r}{$143.70$} &
\multicolumn{1}{r}{$138.70$} & \multicolumn{1}{r}{$148.30$}\\
$\left(  x_{i}-\overline{x}\right)  $ & \multicolumn{1}{r}{$-2.50$} &
\multicolumn{1}{r}{$-1.50$} & \multicolumn{1}{r}{$-0.50$} &
\multicolumn{1}{r}{$0.50$} & \multicolumn{1}{r}{$1.50$} &
\multicolumn{1}{r}{$2.50$}\\
$\left(  x_{i}-\overline{x}\right)  ^{2}$ & \multicolumn{1}{r}{$6.25$} &
\multicolumn{1}{r}{$2.25$} & \multicolumn{1}{r}{$0.25$} &
\multicolumn{1}{r}{$0.25$} & \multicolumn{1}{r}{$2.25$} &
\multicolumn{1}{r}{$6.25$}\\
$\left(  x_{i}-\overline{x}\right)  y_{i}$ & \multicolumn{1}{r}{$-301.25$} &
\multicolumn{1}{r}{$-190.05$} & \multicolumn{1}{r}{$-64.55$} &
\multicolumn{1}{r}{$71.85$} & \multicolumn{1}{r}{$208.05$} &
\multicolumn{1}{r}{$370.75$}\\\hline
\end{tabular}
\]
The totals for the last two rows are 17.5 and 94.8 respectively. We therefore
obtain
\[%
\begin{array}
[c]{lcl}%
\widehat{\beta} & = & 5.417,\\
\widehat{\alpha} & = & -10577.9
\end{array}
\]%
%TCIMACRO{\FRAME{dhFU}{3.3615in}{2.5901in}{0pt}{\Qcb{\QTR{bf}{Figure 5.3}\quad
%Divorces in England and Wales with fitted line}}{}{a5regress_f5.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.3615in;  height 2.5901in;  depth 0pt;  original-width 4.1658in;
%original-height 3.2033in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'a5regress_f5.eps';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.5901in,
width=3.3615in
]%
{a5regress_f5.eps}%
\\
\textbf{Figure 5.3}\quad Divorces in England and Wales with fitted line
\end{center}
%EndExpansion


\noindent$\square$~

\subsubsection{\noindent Linear transformations - some general results}

\noindent We have seen that the normal equations for the linear model,$\left(
\ref{estequ}\right)  $
\[
\mathbf{X}^{T}\mathbf{X}\widehat{\mathbf{\theta}}=\mathbf{X}^{T}\mathbf{Y}
\]
may be written in the form
\[
\widehat{\mathbf{\theta}}=\left(  \mathbf{X}^{T}\mathbf{X}\right)
^{-1}\mathbf{X}^{T}\mathbf{Y}
\]
when $\mathbf{X}^{T}\mathbf{X}$ is non-singular. $\left(  \mathbf{X}%
^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}$ is a linear transformation of the
random vector $\mathbf{Y}$ and it is worth obtaining some general results for
such transformations.\bigskip\newline The results we shall need are:\medskip

\noindent If $\mathbf{W}=\mathbf{LY}$, then

\begin{enumerate}
\item[(i)]
\[
E(\mathbf{W})=\mathbf{L}E\mathbf{(Y})
\]
This a direct result of expectation being a linear operator.

\item[(ii)]
\[
V(\mathbf{W})=\mathbf{L}V(\mathbf{Y})\mathbf{L}^{T}
\]
To see where this comes from, we first need to be clear in our understanding
of $V(\mathbf{W)}$, the covariance matrix of a random vector.\medskip
\newline\noindent The covariance matrix, $V(\mathbf{W)}$, is
\[
\left(
\begin{array}
[c]{cccc}%
V\left(  W_{1}\right)  & C\left(  W_{1},W_{2}\right)  & \cdots & C\left(
W_{1},W_{n}\right) \\
C\left(  W_{2},W_{1}\right)  & V\left(  W_{2}\right)  & \cdots & \vdots\\
\vdots & \vdots & \ddots & \vdots\\
C\left(  W_{n},W_{1}\right)  & \cdots & \cdots & V\left(  W_{n}\right)
\end{array}
\right)
\]
or, in other words,
\[
V(\mathbf{W)}_{ij}=C\left(  W_{i},W_{j}\right)  =E\left[  \left(
W_{i}-E\left(  W_{i}\right)  \right)  \left(  W_{j}-E\left(  W_{j}\right)
\right)  \right]  .
\]
The matrix may be written in the form
\[
V(\mathbf{W})=E\left[  \left(  \mathbf{W}-E(\mathbf{W})\right)  \left(
\mathbf{W}-E(\mathbf{W})\right)  ^{T}\right]  .
\]
Thus
\[%
\begin{array}
[c]{rcl}%
V(\mathbf{W}) & = & E\left[  \left(  \mathbf{W}-E(\mathbf{W})\right)  \left(
\mathbf{W}-E(\mathbf{W})\right)  ^{T}\right] \\
& = & E\left[  \mathbf{L}\left(  \mathbf{Y}-E(\mathbf{Y})\right)  \left(
\mathbf{Y}-E(\mathbf{Y})\right)  ^{T}\mathbf{L}^{T}\right] \\
& = & \mathbf{L}E\left[  \left(  \mathbf{Y}-E(\mathbf{Y})\right)  \left(
\mathbf{Y}-E(\mathbf{Y})\right)  ^{T}\right]  \mathbf{L}^{T}\\
& = & \mathbf{L}V(\mathbf{Y})\mathbf{L}^{T}.
\end{array}
\]
As we shall see in the next section, these results make it easy to deduce some
of the basic properties of least-squares estimators.
\end{enumerate}

\subsubsection{\noindent Properties of the estimators}

\noindent When $\mathbf{X}^{T}\mathbf{X}$ is non-singular, the estimators have
useful properties.\medskip

\noindent\textbf{Property (i)}$\quad\widehat{\mathbf{\theta}}$ is unbiased for
$\mathbf{\theta.}$\newline$\square\medskip$\newline\noindent
\textit{Proof\medskip\newline}%
\[%
\begin{array}
[c]{rcl}%
E\left[  \widehat{\mathbf{\theta}}\right]  & = & E\left[  \left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}\mathbf{Y}\right] \\
& = & \left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}E\left[
\mathbf{Y}\right] \\
& = & \left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}%
\mathbf{X\theta}\\
& = & \mathbf{\theta}%
\end{array}
\]
$\blacksquare\medskip$\newline\noindent Note that this proof does not require
the departures $\varepsilon_{1},\varepsilon_{2},\ldots,\varepsilon_{n}$ either
to have equal variances or to be uncorrelated; it simply uses the fact that
expectation is a linear operator.\medskip

\noindent\textbf{Property (ii)}\quad The covariance matrix of $\widehat
{\mathbf{\theta}}$ is $\mathbf{\Sigma}_{\theta}=\sigma^{2}\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}$.\newline$\square\medskip$%
\newline\noindent\textit{Proof\medskip\newline}For simplicity, let
$\mathbf{B}=\left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}$;
then $\widehat{\mathbf{\theta}}=\mathbf{BY}$. Thus
\[%
\begin{array}
[c]{rcl}%
\mathbf{\Sigma}_{\theta} & = & V\left(  \mathbf{BY}\right) \\
& = & \mathbf{B}V(\mathbf{Y)B}^{T}\\
& = & \mathbf{B\Sigma}_{Y}\mathbf{B}^{T}\\
& = & \sigma^{2}\mathbf{BIB}^{T}\\
& = & \sigma^{2}\left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}%
^{T}\mathbf{X}\left[  \left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\right]
^{T}\\
& = & \sigma^{2}\left[  \left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\right]
^{T}=\sigma^{2}\left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1},
\end{array}
\]
since $\mathbf{X}^{T}\mathbf{X}$ $\left(  \text{and therefore }\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\right)  $is symmetric.\newline%
$\blacksquare\medskip$\newline\noindent This is an important result because it
shows that the diagonal entries of $\left(  \mathbf{X}^{T}\mathbf{X}\right)
^{-1}$ provide information about the precision of the estimates.\medskip

\noindent\textbf{Property (iii)}$\quad\widehat{\mathbf{\theta}}=\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}\mathbf{Y}$ is a linear
transformation of the normally distributed random vector $\mathbf{Y}$, so
that
\[
\widehat{\mathbf{\theta}}\sim N\left(  \mathbf{\theta,}\sigma^{2}\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\right)  .
\]
$\square$

\subsection{An alternative derivation of the estimating equations}

\noindent We have already seen that the model can be written in the form
\[
\mathbf{Y}=\mathbf{X\theta}+\mathbf{\varepsilon}
\]
where $\mathbf{X}$ can be a $n\times p$ matrix of known constants and
$\mathbf{\theta}$ can be a vector of $p$ unknown parameters. We have also seen
that we need to minimise
\[
\sum_{i=1}^{n}\varepsilon_{i}^{2}=\mathbf{\varepsilon}^{T}\mathbf{\varepsilon
=}\left(  \mathbf{y}-\mathbf{X\theta}\right)  ^{T}\left(  \mathbf{y}%
-\mathbf{X\theta}\right)  =\left\Vert \mathbf{y}-\mathbf{X\theta}\right\Vert
^{2}.
\]
It follows, therefore, that the least squares estimator for $\mathbf{\theta} $
is a statistic $\widehat{\mathbf{\theta}}=\mathbf{h}(\mathbf{Y})$ such that
$\left\Vert \mathbf{y}-\mathbf{X\theta}\right\Vert $ achieves its minimum when
$\mathbf{\theta=h}(\mathbf{y})$.\bigskip

\noindent Let $V_{\mathbf{X}}$ be the vector subspace spanned by columns of
$\mathbf{X}$, and let $\mathrm{rank}\left(  \mathbf{X}\right)  =r\leq p$ so
that $\dim\left(  V_{\mathbf{X}}\right)  =r$. Then
\[
E(\mathbf{Y})=\mathbf{X\theta\in}V_{\mathbf{X}}%
\]%
%TCIMACRO{\FRAME{dhFU}{3.9012in}{2.0418in}{0pt}{\Qcb{\QTR{bf}{Figure 5.4}}}%
%{}{a5regress_f6.eps}{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.9012in;  height 2.0418in;  depth 0pt;  original-width 5.0306in;
%original-height 2.6161in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'a5regress_f6.eps';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.0418in,
width=3.9012in
]%
{a5regress_f6.eps}%
\\
\textbf{Figure 5.4}%
\end{center}
%EndExpansion
\medskip\noindent Clearly $\left\Vert \mathbf{y}-\mathbf{X\theta}\right\Vert $
achieves its minimum over $\mathbf{\theta\in}\mathbb{R}^{p}$ when
$\mathbf{X\theta}$ is the orthogonal projection of $\mathbf{y}$ onto
$V_{\mathbf{X}}$. Thus
\[
\mathbf{X}^{T}\left(  \mathbf{y}-\mathbf{X}\widehat{\mathbf{\theta}}\right)
=\mathbf{0}%
\]
or
\[
\mathbf{X}^{T}\mathbf{X}\widehat{\mathbf{\theta}}=\mathbf{X}^{T}\mathbf{Y}%
\]
is the equation satisfied by the least squares estimator $\widehat
{\mathbf{\theta}}$.\bigskip

\subsection{Estimating the variance}

\noindent We shall need to estimate the covariance matrix. Given an estimate
$S^{2}$ of $\sigma^{2}$ this is
\[
\widehat{\mathbf{\Sigma}}_{\theta}=S^{2}\left(  \mathbf{X}^{T}\mathbf{X}%
\right)  ^{-1}.
\]
Now, since
\[
\left(  \mathbf{y}-\mathbf{X\theta}\right)  ^{T}\left(  \mathbf{y}%
-\mathbf{X\theta}\right)  =\mathbf{\varepsilon}^{T}\mathbf{\varepsilon}%
=\sum_{i=1}^{n}\varepsilon_{i}^{2}
\]
and $\sum_{i=1}^{n}\varepsilon_{i}^{2}$ has expected value $n\sigma^{2}$, we
might expect to construct an estimator of $\sigma^{2}$ based on
\[
Q^{2}=\left(  \mathbf{Y}-\mathbf{X}\widehat{\mathbf{\theta}}\right)
^{T}\left(  \mathbf{Y}-\mathbf{X}\widehat{\mathbf{\theta}}\right)
=\mathbf{Y}^{T}\left(  \mathbf{Y}-\mathbf{X}\widehat{\mathbf{\theta}}\right)
.
\]
It is time for a bit of jargon; $\mathbf{Y}-\mathbf{X}\widehat{\mathbf{\theta
}}$ is called the \emph{vector of residuals} and $Q^{2}$ is called the
\emph{residual sum of squares} (or \emph{RSS} for short).\medskip

\noindent Remember that the $r$-dimensional subspace of $\mathbb{R}^{n}$
spanned by columns of $\mathbf{X}$ was denoted by $V_{\mathbf{X}}$. Therefore
construct an orthonormal basis $\left\{  \mathbf{e}_{1},\mathbf{e}_{2,}%
\ldots,\mathbf{e}_{r}\right\}  $ for $V_{\mathbf{X}}$ and extend it to
$\left\{  \mathbf{e}_{1},\ldots,\mathbf{e}_{r},\mathbf{e}_{r+1},\ldots
\mathbf{e}_{n}\right\}  $ for $\mathbb{R}^{n}$. Since $\mathbf{y\in}%
\mathbb{R}^{n}$, it follows that there exist random variables $Z_{1}%
,Z_{2},\ldots,Z_{n}$ such that
\[
\mathbf{Y=}\sum_{i=1}^{n}Z_{i}\mathbf{e}_{i}.
\]
Since $\mathbf{X}\widehat{\mathbf{\theta}}$ is the orthogonal projection of
$\mathbf{Y}$ onto $V_{\mathbf{X}}$,
\[
\mathbf{X}\widehat{\mathbf{\theta}}\mathbf{=}\sum_{i=1}^{r}Z_{i}\mathbf{e}%
_{i}\quad\text{and\quad}\mathbf{Y}-\mathbf{X}\widehat{\mathbf{\theta}%
}\mathbf{=}\sum_{i=r+1}^{n}Z_{i}\mathbf{e}_{i}.
\]
Hence
\[
Q^{2}\mathbf{=}\left(  \sum_{j=r+1}^{n}Z_{j}\mathbf{e}_{j}\right)  ^{T}\left(
\sum_{i=r+1}^{n}Z_{i}\mathbf{e}_{i}\right)  =\sum_{i=r+1}^{n}Z_{i}^{2},
\]
since $\mathbf{e}_{j}^{T}\mathbf{e}_{i}=\delta_{ij}$.\medskip

\noindent Now
\[%
\begin{array}
[c]{ccccl}%
Z_{j}=\mathbf{e}_{j}^{T}\mathbf{Y} & \Rightarrow & E(Z_{j}) & = &
\mathbf{e}_{j}^{T}E(\mathbf{Y})\\
&  &  &  & \\
&  &  & = & \mathbf{e}_{j}^{T}\mathbf{X\theta}=0\quad\text{for}\quad
j=r+1,\ldots,n,
\end{array}
\]
so that%
\[
E(Z_{j}^{2})=V(Z_{j})=\sigma^{2}.
\]
$Z_{j}$ is a linear function of the normally distributed $Y_{i},$
$i=1,\ldots,n,$ so
\[
\frac{Z_{j}}{\sigma}\sim N\left(  0,1\right)  ,\quad\text{and therefore\quad
}\left(  \frac{Z_{j}}{\sigma}\right)  ^{2}\sim\chi^{2}\left(  1\right)  .
\]
Hence%
\[
E(Q^{2})=\sum_{i=r+1}^{n}E(Z_{i}^{2})=\sum_{i=r+1}^{n}V(Z_{i})=(n-r)\sigma
^{2}.
\]
and%
\[
S^{2}=\frac{Q^{2}}{(n-r)}%
\]
is unbiased for $\sigma^{2}$. Furthermore, the $Z_{j}^{\prime}s$ are
independent, so:
\[
\frac{Q^{2}}{\sigma^{2}}\sim\chi^{2}(n-r)\quad\Rightarrow\quad\frac{\left(
n-r\right)  S^{2}}{\sigma^{2}}\sim\chi^{2}(n-r).
\]
Finally, provided $\mathrm{rank}\left(  \mathbf{X}\right)  =p,$ then
\medskip\ \newline%
\[
\widehat{\mathbf{\theta}}=\left(  \mathbf{X}^{T}\mathbf{X}\right)
^{-1}\mathbf{X}^{T}(\mathbf{X}\widehat{\mathbf{\theta}})=\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{X}^{T}\left(  \sum_{i=1}%
^{r}Z_{i}\mathbf{e}_{i}\right)
\]
is independent of $S^{2}$.

\subsection{Testing the coefficients}

\noindent You have already seen that
\[
\widehat{\mathbf{\theta}}\sim N\left(  \mathbf{\theta,}\sigma^{2}\left(
\mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\right)
\]
so that
\[
\widehat{\theta}_{i}\sim N(\theta_{i},\sigma^{2}\left(  \mathbf{X}%
^{T}\mathbf{X}\right)  _{ii}^{-1})
\]
and
\[
\frac{\widehat{\theta}_{i}-\theta_{i}}{\sigma\sqrt{\left(  \mathbf{X}%
^{T}\mathbf{X}\right)  _{ii}^{-1}}}\sim N(0,1).
\]
Using the result that
\[
\frac{(n-r)S^{2}}{\sigma^{2}}\sim\chi^{2}(n-r)
\]
and remembering the definition of a $t$-distribution we conclude that
\[
\frac{\widehat{\theta}_{i}-\theta_{i}}{S\sqrt{\left(  \mathbf{X}^{T}%
\mathbf{X}\right)  _{ii}^{-1}}}\sim t(n-r).
\]
This enables us to carry out hypothesis tests or calculate confidence
intervals for coefficients.\bigskip

\noindent\textbf{Example 5.4 }(Example 5.3 revisited)\quad
\emph{Divorces\medskip}\medskip\newline Our estimate of the rate of increase
of divorces was $\widehat{\beta}=5.417$ and we would like to answer the
question \textquotedblleft Is the divorce rate changing?\textquotedblright\ In
other words, we would like to test the null hypothesis $H_{0}:\beta
=0$.\medskip\ \newline Using lower case letters for realisations, the $RSS$ is
$q^{2}=65.175$ and $n=6,$ $r=2$ so that
\[
s^{2}=\frac{65.175}{4}=16.294.
\]
We have already obtained
\[
\left(  \mathbf{X}^{T}\mathbf{X}\right)  ^{-1}\mathbf{=}\frac{\left(
\begin{array}
[c]{ll}%
\frac{1}{n}\sum x_{i}^{2} & -\frac{1}{n}\sum x_{i}\\
-\frac{1}{n}\sum x_{i} & 1
\end{array}
\right)  }{\sum\left(  x_{i}-\overline{x}\right)  ^{2}}
\]
and
\[
\sum\left(  x_{i}-\overline{x}\right)  ^{2}=17.5.
\]
~

\noindent Under $H_{0}$%
\[
\frac{\widehat{\beta}}{S(17.5)^{-1\left/  2\right.  }}=T\sim t(4)
\]
giving
\[
t=\frac{5.417}{\sqrt{16.294\left/  17.5\right.  }}=5.614.
\]
The \emph{p}-value for the two-sided test is
\[
2P(T\geq5.614)=0.005.
\]
We therefore have strong evidence to reject the null hypothesis that the
divorce rate is not changing; that is, there is strong evidence of an
increasing divorce rate.\newline$\blacksquare$\ \bigskip\newline%
\textbf{Example 5.5}\quad\emph{AIDS data for the USA\medskip}\newline These
data are for AIDS incidence in the USA, adjusted for reporting delays. The
data are taken from Rosenberg, P.S. and Gail, M.H. (1991): Backcalculation of
flexible linear models of the Human Immunodeficiency Virus infection curve.
\textit{Applied Statistics}, \textbf{40}, 269-282.\medskip\newline Newly
reported cases are recorded quarterly and the variable \textit{Time} therefore
counts 3-monthly periods, staring with the first quarter in 1982.\medskip

\begin{center}%
\begin{tabular}
[c]{rrrrrrrrrrr}%
\multicolumn{7}{l}{\textbf{Table 5.3}\quad Reported AIDS cases in the USA} &
&  &  & \\\hline
\textit{Quarter} & \textit{Time} & \textit{Cases} &  & \textit{Quarter} &
\textit{Time} & \textit{Cases} &  & \textit{Quarter} & \textit{Time} &
\textit{Cases}\\\cline{1-3}\cline{5-11}%
1982:1 & 1 & 185 &  & 1984:2 & 10 & 1369 &  & 1986:2 & 18 & 4321\\
1982:2 & 2 & 200 &  & 1984:3 & 11 & 1563 &  & 1986:3 & 19 & 4863\\
1982:3 & 3 & 293 &  & 1984:4 & 12 & 1726 &  & 1986:4 & 20 & 5192\\
1982:4 & 4 & 374 &  & 1985:1 & 13 & 2142 &  & 1987:1 & 21 & 6155\\
1983:1 & 5 & 554 &  & 1985:2 & 14 & 2525 &  & 1987:2 & 22 & 6816\\
1983:2 & 6 & 713 &  & 1985:3 & 15 & 2951 &  & 1987:3 & 23 & 7491\\
1983:3 & 7 & 763 &  & 1985:4 & 16 & 3160 &  & 1987:4 & 24 & 7726\\
1983:4 & 8 & 857 &  & 1986:1 & 17 & 3819 &  & 1988:1 & 25 & 8483\\
1984:1 & 9 & 1147 &  &  &  &  &  &  &  & \\\hline
\end{tabular}
\medskip
\end{center}

\noindent The scatterplot shows that the trend is not linear.
%TCIMACRO{\FRAME{dhFU}{3.8908in}{3.0044in}{0pt}{\Qcb{Incidence of AIDS cases in
%the USA against time}}{}{aids.eps}{\special{ language "Scientific Word";
%type "GRAPHIC";  maintain-aspect-ratio TRUE;  display "USEDEF";
%valid_file "F";  width 3.8908in;  height 3.0044in;  depth 0pt;
%original-width 7.7236in;  original-height 5.9525in;  cropleft "0";
%croptop "1";  cropright "1";  cropbottom "0";
%filename 'aids.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=3.0044in,
width=3.8908in
]%
{aids.eps}%
\\
Incidence of AIDS cases in the USA against time
\end{center}
%EndExpansion
\medskip The plot has all the appearance of showing a functional relationship
between AIDS incidence and time. One might, for example, try to fit an
exponential function or some kind of power law to model the growth
curve.\medskip\newline The plot suggests that the incidence of AIDS against
time in the USA is not linearly related to time. Try fitting a model%
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}x_{i}^{2}+\varepsilon_{i},\quad
i=1,2,\ldots,n
\]
and test the coefficients $\beta_{1}$, $\beta_{2}$.\newline$\square\medskip
$\newline This time we have 3 coefficients, and the design matrix is%
\[
\mathbf{X}=\left(
\begin{array}
[c]{ccc}%
1 & x_{1} & x_{1}^{2}\\
\vdots & \vdots & \vdots\\
1 & x_{n} & x_{n}^{2}%
\end{array}
\right)  ,\quad\text{with}\quad\mathbf{\theta}=\left(
\begin{array}
[c]{l}%
\alpha\\
\beta_{1}\\
\beta_{2}%
\end{array}
\right)  .
\]
We can use a standard computer package to carry out the regression.\medskip
\newline

$%
\begin{tabular}
[c]{lrrccc}%
\multicolumn{6}{l}{Linear Model}\\
\multicolumn{6}{l}{Response: Cases}\\
Variable & \multicolumn{2}{c}{Coefficient} & s.e. & $t$-value & $p$%
-value\\\hline
Intercept & $\alpha$ & $343.5913$ & \multicolumn{1}{r}{$87.7446$} &
\multicolumn{1}{r}{$3.9158$} & $0.0007$\\
Time & $\beta_{1}$ & $-60.1380$ & \multicolumn{1}{r}{$15.5514$} &
\multicolumn{1}{r}{$-3.8671$} & $0.0008$\\
Time$^{2}$ & $\beta_{2}$ & $15.6277$ & $0.5806$ & \multicolumn{1}{r}{$26.9158$%
} & $0.0000$\\\hline
\multicolumn{6}{l}{$r^{2}=0.9976\quad d.f.=22\quad s=134.7\quad RSS=399155$}%
\end{tabular}
\medskip$

\noindent This is a typical computer package output, and we need to say a few
words about what some of these figures mean.\medskip\ \newline The values for
the coefficients $\alpha$, $\beta_{1}$ and $\beta_{2}$ are simple enough to
interpret - the fitted model is
\[
Cases=343.5913-60.1380\times Time+15.6277\times Time^{2}.
\]


\noindent The column labelled \emph{s.e.} gives the standard errors of the
coefficient estimates\emph{\ i.e.} values calculated from the expression you
saw earlier, namely $s\sqrt{\left(  \mathbf{X}^{T}\mathbf{X}\right)
_{ii}^{-1}}$.\medskip\newline The next column gives the calculated
$t$-statistic and the last column gives \emph{p}-values for tests of the
coefficients being zero. \ There is extremely strong evidence for rejecting
the null hypothesis that either $\beta_{1}$ or $\beta_{2}$ is zero. We can
conclude that $Time$ and $Time^{2}$ are highly significant in modelling the
number of cases.\newline$\blacksquare\medskip$

\noindent Some of the numbers given in the last example have not, as yet, been
explained. \emph{RSS} is the residual sum of squares (as you might expect),
but what is $r^{2}$?\medskip

\noindent$r^{2}$ is a measure of the squared correlation. It can be
interpreted as a measure of the strength of the linear relationship. \noindent
You should think of $r^{2}$ as a way of quantifying the reduction in
variability of the response variable brought about by taking account of the
linear relationship. The value of 0.9976 above means that 99.76\% of the
variability is explained by the model, whilst the remaining 0.24\% is the
unexplained or \emph{residual} variability. Clearly, in Example 5.5 the fit
looks extremely good, but can this model be relied upon? \noindent The crucial
question now is \textquotedblleft How can we assess the quality of our
model?\textquotedblright\newline

\subsection{\noindent Assessing the model}

The first step in looking at the adequacy of a model is to check the
assumptions on which it is based, so let us recall what they were.

\begin{enumerate}
\item[$\bullet$] the variance of the response is constant (i.e. the variance
does not depend upon $x$);

\item[$\bullet$] the trend is described by the equation;

\item[$\bullet$] the distribution of the response is normal.\smallskip\ 
\end{enumerate}

\noindent Well, we are fairly happy about the first of these because the
scatterplot we looked at first shows comparatively little scatter anyway. But
the proper way to perform a preliminary check is to examine the
\emph{residuals}, and the easiest way to check the variance assumption is to
look at a plot of residuals versus fitted values and see whether their spread
has any pattern to it.%
%TCIMACRO{\FRAME{dhFU}{3.5034in}{2.706in}{0pt}{\Qcb{Plot of residuals against
%fitted values for AIDS data}}{}{aidsres.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.5034in;  height 2.706in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/aidsres.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.706in,
width=3.5034in
]%
{../Lecture_b8/aidsres.eps}%
\\
Plot of residuals against fitted values for AIDS data
\end{center}
%EndExpansion
\medskip

\begin{itemize}
\item The points to the right seem a little more spread out but there is no
indication that spread is a function of fitted value.

\item There does seem to be some curvature $-$ re-think the assumption about
the model?

\item We also need to check up on the third assumption. If the residuals show
no dependence upon the explanatory variables and are plausibly normal, we
shall be happy.
\end{itemize}

\noindent We have checked the first two assumptions by examining a plot of
residuals against fitted values $-$ we check the third with a normal
probability plot.\medskip%
%TCIMACRO{\FRAME{dhFU}{3.5042in}{2.7069in}{0pt}{\Qcb{Normal probability plot of
%AIDS residuals}}{}{aidsqq.eps}{\special{ language "Scientific Word";
%type "GRAPHIC";  maintain-aspect-ratio TRUE;  display "USEDEF";
%valid_file "F";  width 3.5042in;  height 2.7069in;  depth 0pt;
%original-width 7.7236in;  original-height 5.9525in;  cropleft "0";
%croptop "1";  cropright "1";  cropbottom "0";
%filename '../Lecture_b8/aidsqq.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.7069in,
width=3.5042in
]%
{../Lecture_b8/aidsqq.eps}%
\\
Normal probability plot of AIDS residuals
\end{center}
%EndExpansion
\medskip This doesn't look too bad. We can therefore conclude that we have a
reasonable model which could perhaps be improved, but which fits the data
pretty well for the most part.\ \bigskip\newline\textbf{Example 5.6}%
\quad\textit{Testicular cancer}\medskip\newline The table below comprises data
from Lee, Hitosugi and Peterson (1973): Rise in mortality from tumors of the
testis in Japan, 1947-70. \textit{J. Nat. Cancer Inst.}, \textbf{51}, 1485-90.
It gives the populations and numbers of deaths from testicular cancer in
5-year age groups and 5-year periods in Japan. The ages refer to the lowest
age in each group and the populations are expressed in millions of persons.

\begin{center}
$%
\begin{tabular}
[c]{ccccccccc}%
\multicolumn{9}{l}{\textbf{Table 5.4}\quad Deaths in Japan from testicular
cancer}\\\hline
& \multicolumn{2}{c}{1951-55} & \multicolumn{2}{c}{1956-60} &
\multicolumn{2}{c}{1961-65} & \multicolumn{2}{c}{1966-70}\\
Age & Popn. & Dths & Popn. & Dths & Popn. & Dths & Popn. & Dths\\\hline
{20 } & \multicolumn{1}{r}{20.4} & \multicolumn{1}{r}{27} &
\multicolumn{1}{r}{21.3} & \multicolumn{1}{r}{39} & \multicolumn{1}{r}{22.2} &
\multicolumn{1}{r}{56} & \multicolumn{1}{r}{24.0} & \multicolumn{1}{r}{\ 83}\\
{25 } & \multicolumn{1}{r}{17.2} & \multicolumn{1}{r}{40} &
\multicolumn{1}{r}{20.0} & \multicolumn{1}{r}{58} & \multicolumn{1}{r}{20.6} &
\multicolumn{1}{r}{97} & \multicolumn{1}{r}{21.8} & \multicolumn{1}{r}{125}\\
{30 } & \multicolumn{1}{r}{12.6} & \multicolumn{1}{r}{18} &
\multicolumn{1}{r}{17.1} & \multicolumn{1}{r}{54} & \multicolumn{1}{r}{19.9} &
\multicolumn{1}{r}{77} & \multicolumn{1}{r}{20.8} & \multicolumn{1}{r}{129}\\
{35 } & \multicolumn{1}{r}{11.7} & \multicolumn{1}{r}{13} &
\multicolumn{1}{r}{12.5} & \multicolumn{1}{r}{36} & \multicolumn{1}{r}{17.0} &
\multicolumn{1}{r}{70} & \multicolumn{1}{r}{19.9} & \multicolumn{1}{r}{101}\\
{40 } & \multicolumn{1}{r}{11.5} & \multicolumn{1}{r}{26} &
\multicolumn{1}{r}{11.5} & \multicolumn{1}{r}{32} & \multicolumn{1}{r}{12.2} &
\multicolumn{1}{r}{29} & \multicolumn{1}{r}{16.8} & \multicolumn{1}{r}{\ 67}\\
{45 } & \multicolumn{1}{r}{10.3} & \multicolumn{1}{r}{16} &
\multicolumn{1}{r}{11.2} & \multicolumn{1}{r}{26} & \multicolumn{1}{r}{11.1} &
\multicolumn{1}{r}{34} & \multicolumn{1}{r}{12.0} & \multicolumn{1}{r}{\ 37}\\
{50 } & \multicolumn{1}{r}{\ 9.3} & \multicolumn{1}{r}{16} &
\multicolumn{1}{r}{\ 9.8} & \multicolumn{1}{r}{27} & \multicolumn{1}{r}{10.7}
& \multicolumn{1}{r}{27} & \multicolumn{1}{r}{10.7} & \multicolumn{1}{r}{\ 29}%
\\
{55 } & \multicolumn{1}{r}{\ 7.6} & \multicolumn{1}{r}{17} &
\multicolumn{1}{r}{\ 8.7} & \multicolumn{1}{r}{19} & \multicolumn{1}{r}{\ 9.2}
& \multicolumn{1}{r}{32} & \multicolumn{1}{r}{10.1} & \multicolumn{1}{r}{\ 39}%
\\
{60 } & \multicolumn{1}{r}{\ 5.9} & \multicolumn{1}{r}{13} &
\multicolumn{1}{r}{\ 6.8} & \multicolumn{1}{r}{21} & \multicolumn{1}{r}{\ 7.9}
& \multicolumn{1}{r}{21} & \multicolumn{1}{r}{\ 8.4} &
\multicolumn{1}{r}{\ 31}\\\hline
\end{tabular}
\medskip$
\end{center}

\noindent The scatterplot shows that the mean death rate from cancer of the
testis in Japan has been rising steadily since 1951.
%TCIMACRO{\FRAME{dhFU}{3.4117in}{2.8245in}{0pt}{\Qcb{Death rate against year
%for cancer of the testis in Japan}}{}{testis.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.4117in;  height 2.8245in;  depth 0pt;  original-width 4.5109in;
%original-height 3.7291in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'testis.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.8245in,
width=3.4117in
]%
{testis.eps}%
\\
Death rate against year for cancer of the testis in Japan
\end{center}
%EndExpansion
\medskip\noindent Note that there is no variability in the explanatory
variable, but there is marked variability in the response. This is clearly
shown in the boxplots.%
%TCIMACRO{\FRAME{dhFU}{3.8908in}{3.0044in}{0pt}{\Qcb{Boxplots of death rates}%
%}{}{testis_box.eps}{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.8908in;  height 3.0044in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";  filename 'testis_box.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=3.0044in,
width=3.8908in
]%
{testis_box.eps}%
\\
Boxplots of death rates
\end{center}
%EndExpansion
\medskip It is clear that the variability in the data increases as the year
variable increases. What should one do about this?\medskip\newline The answer
is to look for a transformation which will stabilise the variance. Here we
need a transformation which compresses large values of the response more than
it compresses smaller values; something like a square root or a cube root or
possibly even a log transformation. Taking the log of the death rate results
in the next scatterplot.%
%TCIMACRO{\FRAME{dhFU}{3.5034in}{2.706in}{0pt}{\Qcb{Testicular cancer: plot of
%log(Death rate) against year}}{}{logtestis.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.5034in;  height 2.706in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/logtestis.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.706in,
width=3.5034in
]%
{../Lecture_b8/logtestis.eps}%
\\
Testicular cancer: plot of log(Death rate) against year
\end{center}
%EndExpansion
\medskip This looks reasonable and we could now go ahead and fit a model of
the form%
\[
\log\left(  Y_{i}\right)  =\alpha+\beta x_{i}+\varepsilon_{i},\quad
i=1,2,\ldots,n.
\]
With the log transformed data, the fitted model turns out to be\medskip\ as
given below.\newline$\qquad%
\begin{tabular}
[c]{lrrccc}%
\multicolumn{6}{l}{Linear Model}\\
\multicolumn{6}{l}{Response: log(Death rate)}\\
Variable & \multicolumn{2}{c}{Coefficient} & s.e. & $t$-value & $p$%
-value\\\hline
Intercept & $\alpha$ & $-105.9198$ & \multicolumn{1}{r}{$14.4887$} &
\multicolumn{1}{r}{$-7.3105$} & $0.0000$\\
Year & $\beta$ & $0.0545$ & \multicolumn{1}{r}{$0.0074$} &
\multicolumn{1}{r}{$7.3808$} & $0.0000$\\\hline
\multicolumn{6}{l}{$r^{2}=0.6157\quad d.f.=34\quad s=0.2479\quad RSS=2.0891$}%
\end{tabular}
\medskip\medskip$\newline Of course stable variance alone is not enough
because the residuals also need to be normally distributed. This can be
checked with a normal q-q plot of the residuals.%
%TCIMACRO{\FRAME{dhFU}{3.5042in}{2.7069in}{0pt}{\Qcb{Normal probability plot of
%testicular cancer residuals}}{}{testisres.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.5042in;  height 2.7069in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/testisres.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.7069in,
width=3.5042in
]%
{../Lecture_b8/testisres.eps}%
\\
Normal probability plot of testicular cancer residuals
\end{center}
%EndExpansion
\medskip The plot seems to show a very rough straight line, but it is not
entirely convincing.\newline$\square\medskip$\newline It is worth noting that,
in doing Example 5.6, we fitted a model with more than two coefficients. How
is this different from fitting more than one explanatory variable? How can the
technique of fitting the model
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}x_{i}^{2}+\varepsilon_{i},\quad
i=1,2,\ldots,n,
\]
be any different from fitting, say%
\[
Y_{i}=\alpha+\beta_{1}x_{i}+\beta_{2}z_{i}+\varepsilon_{i},\quad
i=1,2,\ldots,n,
\]
where $z_{i}$ is some other variable which has been recorded (\emph{e.g.}
annual expenditure on AIDS clinics)? The answer is, of course, that
mathematically they are exactly the same.

\noindent

\subsection{Multiple regression}

\noindent With more than one explanatory variable, the model takes the form
\[
Y_{i}=\alpha+\beta_{1}x_{i1}+\beta_{2}x_{i2}+\cdots+\beta_{p}x_{ip}%
+\varepsilon_{i}
\]
or, as we have seen,
\[
\mathbf{Y}=\mathbf{X\theta+\varepsilon}
\]
where $\mathbf{Y}=(Y_{1},\ldots,Y_{n})^{T}$, $\mathbf{\theta}=(\alpha
,\beta_{1},\ldots,\beta_{p})^{T}$, $\mathbf{\varepsilon}=(\varepsilon
_{1},\ldots,\varepsilon_{n})^{T}$ and $\mathbf{X}$ is the $n\times\left(
1+p\right)  $ matrix
\[
\mathbf{X=}\left(
\begin{array}
[c]{cccc}%
1 & x_{11} & \cdots & x_{1p}\\
1 & x_{21} & \cdots & x_{2p}\\
\vdots & \vdots & \ddots & \vdots\\
1 & x_{n1} & \cdots & x_{np}%
\end{array}
\right)  .
\]
Let us see not only how to apply the general theory but also how to set about
the whole modelling process by looking at a data set.\bigskip

\noindent\textbf{Example 5.7}\quad\emph{Olympic Sprint Times\medskip\newline%
}Table 5.2 gives the times in seconds recorded by the winners in the finals of
the men's sprint events (100, 200, 400, 800 and 1500 metres) at each of the 21
Olympic Games from 1900 to 2004 along with the heights above sea level of the
different venues. Obviously the years 1916, 1940, 1944 are missing since the
Olympic games were not held during the World Wars.$\medskip$\newline%
$\hspace{1.5in}%
\begin{tabular}
[c]{llllllll}%
\multicolumn{8}{l}{\textbf{Table 5.2}\quad Winning times for Olympic running
events}\\\hline
\textit{Year} & \textit{Venue} & \textit{100m} & \textit{200m} & \textit{400m}
& \textit{800m} & \textit{1500m} & \textit{Altitude}\\\hline
\multicolumn{1}{c}{1900} & Paris & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.20} & \multicolumn{1}{c}{49.40} &
\multicolumn{1}{r}{121.40} & 246.00 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1904} & St Louis & \multicolumn{1}{c}{11.00} &
\multicolumn{1}{c}{21.60} & \multicolumn{1}{c}{49.20} &
\multicolumn{1}{r}{116.00} & 245.40 & \multicolumn{1}{r}{455}\\
\multicolumn{1}{c}{1908} & London & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.40} & \multicolumn{1}{c}{50.00} &
\multicolumn{1}{r}{112.80} & 243.40 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1912} & Stockholm & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{21.70} & \multicolumn{1}{c}{48.20} &
\multicolumn{1}{r}{111.90} & 236.80 & \multicolumn{1}{r}{46}\\
\multicolumn{1}{c}{1920} & Antwerp & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{22.00} & \multicolumn{1}{c}{49.60} &
\multicolumn{1}{r}{113.40} & 241.80 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1924} & Paris & \multicolumn{1}{c}{10.60} &
\multicolumn{1}{c}{21.60} & \multicolumn{1}{c}{47.60} &
\multicolumn{1}{r}{112.40} & 233.60 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1928} & Amsterdam & \multicolumn{1}{c}{10.80} &
\multicolumn{1}{c}{21.80} & \multicolumn{1}{c}{47.80} &
\multicolumn{1}{r}{111.80} & 233.20 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1932} & Los Angeles & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{21.20} & \multicolumn{1}{c}{46.20} &
\multicolumn{1}{r}{109.80} & 231.20 & \multicolumn{1}{r}{340}\\
\multicolumn{1}{c}{1936} & Berlin & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{20.70} & \multicolumn{1}{c}{46.50} &
\multicolumn{1}{r}{112.90} & 227.80 & \multicolumn{1}{r}{115}\\
\multicolumn{1}{c}{1948} & London & \multicolumn{1}{c}{10.30} &
\multicolumn{1}{c}{21.10} & \multicolumn{1}{c}{46.20} &
\multicolumn{1}{r}{109.20} & 225.20 & \multicolumn{1}{r}{8}\\
\multicolumn{1}{c}{1952} & Helsinki & \multicolumn{1}{c}{10.40} &
\multicolumn{1}{c}{20.70} & \multicolumn{1}{c}{45.90} &
\multicolumn{1}{r}{109.20} & 225.20 & \multicolumn{1}{r}{25}\\
\multicolumn{1}{c}{1956} & Melbourne & \multicolumn{1}{c}{10.50} &
\multicolumn{1}{c}{20.60} & \multicolumn{1}{c}{46.70} &
\multicolumn{1}{r}{107.70} & 221.20 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1960} & Rome & \multicolumn{1}{c}{10.20} &
\multicolumn{1}{c}{20.50} & \multicolumn{1}{c}{44.90} &
\multicolumn{1}{r}{106.30} & 215.60 & \multicolumn{1}{r}{66}\\
\multicolumn{1}{c}{1964} & Tokyo & \multicolumn{1}{c}{10.00} &
\multicolumn{1}{c}{20.30} & \multicolumn{1}{c}{45.10} &
\multicolumn{1}{r}{105.10} & 218.10 & \multicolumn{1}{r}{45}\\
\multicolumn{1}{c}{1968} & Mexico City & \multicolumn{1}{c}{9.95} &
\multicolumn{1}{c}{19.83} & \multicolumn{1}{c}{43.80} &
\multicolumn{1}{r}{104.30} & 214.90 & \multicolumn{1}{r}{7349}\\
\multicolumn{1}{c}{1972} & Munich & \multicolumn{1}{c}{10.14} &
\multicolumn{1}{c}{20.00} & \multicolumn{1}{c}{44.66} &
\multicolumn{1}{r}{105.90} & 216.30 & \multicolumn{1}{r}{1699}\\
\multicolumn{1}{c}{1976} & Montreal & \multicolumn{1}{c}{10.06} &
\multicolumn{1}{c}{20.23} & \multicolumn{1}{c}{44.26} &
\multicolumn{1}{r}{103.50} & 219.20 & \multicolumn{1}{r}{104}\\
\multicolumn{1}{c}{1980} & Moscow & \multicolumn{1}{c}{10.25} &
\multicolumn{1}{c}{20.19} & \multicolumn{1}{c}{44.60} &
\multicolumn{1}{r}{105.40} & 218.40 & \multicolumn{1}{r}{497}\\
\multicolumn{1}{c}{1984} & Los Angeles & \multicolumn{1}{c}{9.99} &
\multicolumn{1}{c}{19.80} & \multicolumn{1}{c}{44.27} &
\multicolumn{1}{r}{103.00} & 212.50 & \multicolumn{1}{r}{340}\\
\multicolumn{1}{c}{1988} & Seoul & \multicolumn{1}{c}{9.92} &
\multicolumn{1}{c}{19.75} & \multicolumn{1}{c}{43.87} &
\multicolumn{1}{r}{103.45} & 215.96 & \multicolumn{1}{r}{111}\\
\multicolumn{1}{c}{1992} & Barcelona & \multicolumn{1}{c}{9.96} &
\multicolumn{1}{c}{20.01} & \multicolumn{1}{c}{43.50} &
\multicolumn{1}{r}{103.66} & 220.12 & \multicolumn{1}{r}{3}\\
\multicolumn{1}{c}{1996} & Atlanta & \multicolumn{1}{c}{9.84} &
\multicolumn{1}{c}{19.32} & \multicolumn{1}{c}{43.49} &
\multicolumn{1}{r}{102.58} & 215.78 & \multicolumn{1}{r}{1026}\\
\multicolumn{1}{c}{2000} & Sydney & \multicolumn{1}{c}{9.87} &
\multicolumn{1}{c}{20.09} & \multicolumn{1}{c}{43.84} &
\multicolumn{1}{r}{105.08} & 212.07 & \multicolumn{1}{r}{3}\\\hline
2004 & Athens &  & 19.79 & 44.00 & 104.45 & 214.20 & \multicolumn{1}{r}{505}%
\\\hline
\end{tabular}
\medskip$\newline We would like to construct a predictive model for the
winning times at the 2008 games in Beijing (and the 2012 games in
London).\ Clearly there is now a possibility of more than one explanatory
variable.$\medskip$\newline Since the 100m event is largely about getting a
fast start rather than keeping up a sustained speed, we shall leave that event
out of the analysis. Let us start by looking at a plot of average speed in
metres per second against distance of event.%
%TCIMACRO{\FRAME{dhFU}{3.4117in}{2.4509in}{0pt}{\Qcb{Plot of \QTR{it}{Speed}
%against \QTR{it}{Distance}}}{}{sprints1.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.4117in;  height 2.4509in;  depth 0pt;  original-width 4.5109in;
%original-height 3.2309in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/sprints1.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.4509in,
width=3.4117in
]%
{../Lecture_b8/sprints1.eps}%
\\
Plot of \textit{Speed} against \textit{Distance}%
\end{center}
%EndExpansion
\medskip\noindent Clearly the relationship is not linear and we need a
transformation.\medskip\ 

\noindent You can see that the standard deviations of the speeds do not differ
much between events. Therefore the response should \emph{not} be transformed.
Try taking the log of the explanatory variable.%

%TCIMACRO{\FRAME{dhFU}{3.4342in}{2.4509in}{0pt}{\Qcb{Plot of \QTR{it}{Speed}
%against log(\QTR{it}{Distance})}}{}{sprints2.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.4342in;  height 2.4509in;  depth 0pt;  original-width 4.542in;
%original-height 3.2309in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/sprints2.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.4509in,
width=3.4342in
]%
{../Lecture_b8/sprints2.eps}%
\\
Plot of \textit{Speed} against log(\textit{Distance})
\end{center}
%EndExpansion
\medskip\noindent That's better, but it still may not be all that good.
However, we can try a preliminary fit and see where that gets us.\medskip\ 

\noindent The computer produces the following $\medskip$\newline$\qquad%
\begin{tabular}
[c]{lrrccc}%
\multicolumn{6}{l}{Linear regression analysis}\\
\multicolumn{6}{l}{Response: Speed}\\
Variable & \multicolumn{2}{c}{Coefficient} & s.e. & $t$-value & $p$%
-value\\\hline
Constant & $\alpha$ & $17.7837$ & \multicolumn{1}{r}{$0.3233$} &
\multicolumn{1}{r}{$55.01$} & $<2\times10^{-16}$\\
Log(distance) & $\beta$ & $-1.5296$ & \multicolumn{1}{r}{$0.0508$} &
\multicolumn{1}{r}{$-30.12$} & $<2\times10^{-16}$\\\hline
\multicolumn{6}{l}{$r^{2}=0.9061\quad d.f.=94\quad s=0.3749\quad RSS=12.6108$}%
\end{tabular}
$\ \medskip

\noindent On\noindent\ the face of it this does not look too bad with highly
significant coefficients and $r^{2}=0.9061$. However, we need to check the residuals.%

%TCIMACRO{\FRAME{dhFU}{3.1168in}{2.4076in}{0pt}{\Qcb{Plot of residuals against
%fitted values for log(distance) model}}{}{logdist.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.1168in;  height 2.4076in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/logdist.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.4076in,
width=3.1168in
]%
{../Lecture_b8/logdist.eps}%
\\
Plot of residuals against fitted values for log(distance) model
\end{center}
%EndExpansion
\medskip\noindent It is clear that our model will not do. A clear pattern is
present and the log transformation has not worked.\medskip

\noindent The fit can be substantially improved by, instead of taking the log
of distance, adding a quadratic term to it. In order to keep the numbers
sensible, divide \emph{distance} by 100 and re-fit as a quadratic. $\medskip
$\newline$\qquad%
\begin{tabular}
[c]{lcrccc}%
\multicolumn{6}{l}{Linear regression analysis}\\
\multicolumn{6}{l}{Response: Speed}\\
Variable & \multicolumn{2}{c}{Coefficient} & s.e. & $t$-value & $p$%
-value\\\hline
Constant & $\alpha$ & $10.7571$ & \multicolumn{1}{r}{$0.1229$} & $87.53$ &
$<2\times10^{-16}$\\
Distance & $\beta_{1}$ & $-0.5918$ & \multicolumn{1}{r}{$0.0378$} & $-15.68$ &
$<2\times10^{-16}$\\
Distance$^{2}$ & $\beta_{2}$ & $0.0213$ & \multicolumn{1}{r}{$0.0021$} &
$10.03$ & $<2\times10^{-16}$\\\hline
\multicolumn{6}{l}{$r^{2}=0.9134\quad d.f.=93\quad s=0.375\quad RSS=11.6601$}%
\end{tabular}
\medskip$

\noindent The value of $r^{2}$ has improved to $0.9134$, all of the
coefficients have values significantly different from zero and a look at a
plot of residuals against fitted values shows that the pattern in the
residuals has been largely removed.\ %

%TCIMACRO{\FRAME{dhFU}{3.1168in}{2.4076in}{0pt}{\Qcb{Plot of residuals against
%fitted values with a quadratic term included}}{}{sqdistres.eps}%
%{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 3.1168in;  height 2.4076in;  depth 0pt;  original-width 7.7236in;
%original-height 5.9525in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/sqdistres.EPS';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.4076in,
width=3.1168in
]%
{../Lecture_b8/sqdistres.eps}%
\\
Plot of residuals against fitted values with a quadratic term included
\end{center}
%EndExpansion
$\medskip$\newline\noindent We should now go on and construct a normal
probability plot of the residuals, but there is little point because it is
hard to believe that this is going to be our final model. There are obvious
things we haven't taken into account.

\begin{itemize}
\item \noindent Every Olympics, runners go faster, though the improvement has
started to slow down recently. Obviously any proper treatment must include
\emph{year} as a variable (or rather some appropriate functions of \emph{year}).

\item We also need to include height above sea-level of each Olympic venue.
This made a big difference at the Mexico Olympics.
\end{itemize}

\noindent It turns out that we need to fit \emph{speed} against
\emph{distance}, \emph{distance}$^{2}$, \emph{year, year}$^{2},$ and
log(\emph{altitude}).\medskip\newline$\qquad%
\begin{tabular}
[c]{lcrccc}%
\multicolumn{6}{l}{Linear regression analysis}\\
\multicolumn{6}{l}{Response: Speed}\\
Variable & \multicolumn{2}{c}{Coefficient} & s.e. & $t$-value & $p$%
-value\\\hline
Constant & $\alpha$ & $10.001$ & \multicolumn{1}{r}{$0.0471$} &
\multicolumn{1}{r}{$212.394$} & $<2\times10^{-16}$\\
Distance & $\beta_{1}$ & $-0.5918$ & \multicolumn{1}{r}{$0.0109$} &
\multicolumn{1}{r}{$-54.570$} & $<2\times10^{-16}$\\
Distance$^{2}$ & $\beta_{2}$ & $0.0213$ & \multicolumn{1}{r}{$0.0006$} &
\multicolumn{1}{r}{$34.912$} & $<2\times10^{-16}$\\
Year & $\beta_{3}$ & $0.0161$ & \multicolumn{1}{r}{$\ \ \ 0.0013$} &
\multicolumn{1}{r}{$12.375$} & $<2\times10^{-16}$\\
Year$^{2}$ & $\beta_{4}$ & $-0.000059$ & \multicolumn{1}{r}{$0.000012$} &
\multicolumn{1}{r}{$-4.869$} & $4.8\times10^{-6}$\\
log(Altitude) & $\beta_{5}$ & $0.0254$ & \multicolumn{1}{r}{$0.0050$} &
\multicolumn{1}{r}{$5.010$} & $2.7\times10^{-6}$\\\hline
\multicolumn{6}{l}{$r^{2}=0.9931\quad d.f.=90\quad s=0.1040\quad RSS=0.9736$}%
\end{tabular}
$%

%TCIMACRO{\FRAME{dhFU}{5.8211in}{2.1741in}{0pt}{\Qcb{Residuals against fitted
%values and normal probability plot for the Olympic sprints model}}%
%{}{sprintsres.eps}{\special{ language "Scientific Word";  type "GRAPHIC";
%maintain-aspect-ratio TRUE;  display "USEDEF";  valid_file "F";
%width 5.8211in;  height 2.1741in;  depth 0pt;  original-width 7.7236in;
%original-height 2.8617in;  cropleft "0";  croptop "1";  cropright "1";
%cropbottom "0";
%filename '../Lecture_b8/sprintsres.eps';file-properties "XNPEU";}} }%
%BeginExpansion
\begin{center}
\includegraphics[
height=2.1741in,
width=5.8211in
]%
{../Lecture_b8/sprintsres.eps}%
\\
Residuals against fitted values and normal probability plot for the Olympic
sprints model
\end{center}
%EndExpansion
\medskip The plots are satisfactory and the model looks pretty good. We now
have something we can use to make predictions about the times in the 2008
Olympic Games in Beijing. However, there is more to using such models for
prediction than at first meets the eye, as you will see next.\newline$\square$

\subsection{Making predictions}

\noindent Once the coefficients have been estimated, the fitted equation can
be used to obtain values of $y$ for any given values $\mathbf{a=(}%
1,a_{1},\ldots,a_{p})$ of the explanatory variables $x=\mathbf{(}%
1,x_{1},\ldots,x_{p})$. There is a population of such $y$-values and, for a
random observation from this population, $Y\sim N(\mathbf{a\theta},\sigma
^{2})$.\medskip\newline First consider estimating the mean of the population,
$E(Y\mid\mathbf{a})$. Clearly
\[
\widetilde{Y}=\mathbf{a}\widehat{\mathbf{\theta}},
\]
where $\widehat{\mathbf{\theta}}$ is a random vector of estimators, is an
unbiased estimator of $\mathbf{a\theta}$.\medskip\ 

\noindent Furthermore, $\widetilde{Y}$ must be normally distributed and
\[
E(\widetilde{Y})=\mathbf{a\theta},\quad V(\widetilde{Y})=\mathbf{a}%
V\mathbf{\left(  \widehat{%
%TCIMACRO{\QTR{mathbf}{\theta}}%
%BeginExpansion
\theta
%EndExpansion
}\right)  a}^{T}=\sigma^{2}\mathbf{a\left(
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
^{T}\mathbf{X}\right)  ^{-1}a}^{T}.
\]
\textbf{Example 5.8}\quad\emph{Simple linear regression}\medskip\ 

\noindent Suppose we wish to predict $E\left(  Y\right)  $ at the point
$x_{0}$. Then $\mathbf{a=(}1,x_{0})$ so
\[%
\begin{array}
[c]{rcl}%
E(\widetilde{Y}) & = & \alpha+\beta x_{0}\\
V(\widetilde{Y}) & = & V(\widehat{\alpha})+x_{0}^{2}V(\widehat{\beta}%
)+2x_{0}C(\widehat{\alpha},\widehat{\beta})
\end{array}
\]
\noindent As a point of technique, with the simple linear model, it is usually
better to fit
\[
Y_{i}=\alpha+\beta(x_{i}-\overline{x})+\varepsilon_{i}
\]
Then
\begin{align*}
\mathbf{X}^{T}\mathbf{X}  & =\left(
\begin{array}
[c]{cc}%
n & 0\\
0 & \sum(x_{i}-\overline{x})^{2}%
\end{array}
\right)  ,\\
\mathbf{\left(
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
^{T}%
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
\right)  ^{-1}}  & =\left(
\begin{array}
[c]{cc}%
1\left/  n\right.  & 0\\
0 & 1\left/  \sum(x_{i}-\overline{x})^{2}\right.
\end{array}
\right)
\end{align*}
so that $C(\widehat{\alpha},\widehat{\beta})=0$ and $V(\widehat{\alpha})$,
$V(\widehat{\beta})$ are straightforward.\newline$\square$\medskip\ 

\noindent\textbf{Fact 5.9}

\noindent Suppose the two random variables $W$ and $Z$ are jointly normally
distributed. (In other words, $W$ and $Z$ have a bivariate normal distribution.)

Then if $W$ and $Z$ are uncorrelated, i.e. they satisfy $C(W,V)=0,$ they are
independent. \bigskip

\noindent Applying Fact 5.9 to the previous example, $\widehat{\alpha}$ and
$\widehat{\beta}$ are both linear combinations of the normal random variables
$Y_{i}$, so have a joint normal distribution. Further, the new
parameterization gives $C(\widehat{\alpha},\widehat{\beta})=0$, so the
estimators $\widehat{\alpha}$ and $\widehat{\beta}$ are \textit{independent.}
\bigskip

\noindent In general it is easy to construct a confidence interval for
$E(Y\mid\mathbf{a})=\mathbf{a\theta}$. Since \textrm{rank}$(\mathbf{X})=p+1$,
\[
\frac{\mathbf{a}\widehat{\mathbf{\theta}}-\mathbf{a\theta}}{\widehat{\sigma
}\sqrt{\mathbf{a\left(
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
^{T}\mathbf{X}\right)  ^{-1}a}^{T}}}\sim t(n-p-1)
\]
and an inequality can be constructed and re-arranged for $\mathbf{a\theta}$ in
the usual way.\bigskip

\noindent\textbf{Example 5.10}\quad\emph{Simple linear regression}\medskip\ 

\noindent For a given value of $x_{0}$,
\[
E(Y\mid x_{0})=\alpha+\beta(x_{0}-\overline{x}),
\]
so
\[
\frac{\widehat{\alpha}+\widehat{\beta}(x_{0}-\overline{x})-\left(
\alpha+\beta(x_{0}-\overline{x})\right)  }{\widehat{\sigma}\sqrt{(1\left/
n\right.  +(x_{0}-\overline{x})^{2}\left/  \sum(x_{i}-\overline{x}%
)^{2}\right.  }}\sim t(n-2).
\]
If $m$ is such that $P\left(  \left\vert T_{n-2}\right\vert \geq m\right)
=\gamma,$%
\[
\widehat{\alpha}+\widehat{\beta}(x_{0}-\overline{x})\pm m\widehat{\sigma}%
\sqrt{(1\left/  n\right.  +(x_{0}-\overline{x})^{2}\left/  \sum(x_{i}%
-\overline{x})^{2}\right.  }
\]
is a $100(1-\gamma)$\% confidence interval for $E(Y\mid x_{0})$.\newline%
$\square\medskip$

\noindent But suppose we wanted to make a prediction of an as-yet unobserved
random variable $Y$ (\emph{i.e}. not the expectation, but the variable
itself). We are after a \emph{predictive interval} on a random variable $Y$
yet to be observed.\medskip\newline Now $Y\sim N(\mathbf{a\theta},\sigma^{2})$
for a vector $\mathbf{a}$ of given values, and $Y-\mathbf{a}\widehat
{\mathbf{\theta}}$ has a normal distribution with mean
\[
E\left(  Y-\mathbf{a}\widehat{\mathbf{\theta}}\right)  =\mathbf{a\theta
}-\mathbf{a\theta}=0
\]
and, since $Y$ and $\widehat{\mathbf{\theta}}$ are independent, with variance
\[%
\begin{array}
[c]{rcl}%
V\left(  Y-\mathbf{a}\widehat{\mathbf{\theta}}\right)   & = &
V(Y)+V(\mathbf{a}\widehat{\mathbf{\theta}})\\
& = & \sigma^{2}+\sigma^{2}\mathbf{a\left(
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
^{T}\mathbf{X}\right)  ^{-1}a}^{T}.
\end{array}
\]
Note the extra term added to the variance.\medskip\newline Using the
independence of $Y,\widehat{\mathbf{\theta}}$ and $\widehat{\sigma}$,
\begin{equation}
\frac{Y-\mathbf{a}\widehat{\mathbf{\theta}}}{\widehat{\sigma}\sqrt
{1+\mathbf{a\left(
%TCIMACRO{\QTR{mathbf}{X}}%
%BeginExpansion
X%
%EndExpansion
^{T}\mathbf{X}\right)  ^{-1}a}^{T}}}\sim t(n-p-1).\label{predint}%
\end{equation}
\noindent\textbf{Example 5.11}\quad\emph{Simple linear regression\medskip
}\newline\noindent For a given value of $x_{0}$,
\[
\frac{Y-(\widehat{\alpha}+\widehat{\beta}(x_{0}-\overline{x}))}{\widehat
{\sigma}\sqrt{(1+1\left/  n\right.  +(x_{0}-\overline{x})^{2}\left/
\sum(x_{i}-\overline{x})^{2}\right.  }}\sim t(n-2).
\]
If $m$ is such that $P\left(  \left\vert T_{n-2}\right\vert \geq m\right)
=\gamma,$%
\[
\widehat{\alpha}+\widehat{\beta}(x_{0}-\overline{x})\pm m\widehat{\sigma}%
\sqrt{(1+1\left/  n\right.  +(x_{0}-\overline{x})^{2}\left/  \sum
(x_{i}-\overline{x})^{2}\right.  }%
\]
is a $100(1-\gamma)$\% predictive interval for $Y$.\newline$\square\bigskip
$\pagebreak\newline\textbf{Example 5.12}\quad\emph{Predictions for the 2008
Olympic Games.\medskip}

\noindent The model we have arrived at makes prediction of the times for the
Olympic sprint events in 2008 comparatively straightforward. Beijing is 120
feet above sea level so the predicted winning times for the events
are\medskip\ as given below.\newline%
\[%
\begin{tabular}
[c]{lllll}
& \textit{200m} & \textit{400m} & \textit{800m} & \textit{1500m}\\
\textit{Speed} & \multicolumn{1}{r}{10.09} & \multicolumn{1}{r}{9.16} &
\multicolumn{1}{r}{7.81} & \multicolumn{1}{r}{7.11}\\
\textit{Time} & \multicolumn{1}{r}{19.83} & \multicolumn{1}{r}{43.68} &
\multicolumn{1}{r}{102.38} & \multicolumn{1}{r}{211.11}%
\end{tabular}
\]
\medskip\newline Confidence intervals for the predicted average speeds can be
obtained by using Equation \ref{predint}; these can then be converted back to
times for the events. In practice a computer can be used to do the matrix
multiplications, although most statistical packages will do the calculations
directly as one of the options provided with their regression
software.\medskip\newline%
\begin{tabular}
[c]{lcccccccc}
& \multicolumn{2}{c}{\textit{200m}} & \multicolumn{2}{c}{\textit{400m}} &
\multicolumn{2}{c}{\textit{800m}} & \multicolumn{2}{c}{\textit{1500m}}\\
\textit{Speed} & \multicolumn{2}{c}{10.09} & \multicolumn{2}{c}{9.16} &
\multicolumn{2}{c}{7.81} & \multicolumn{2}{c}{7.11}\\
\textit{95\% c.i.} & 9.87 & 10.30 & 8.94 & 9.38 & 7.59 & 8.03 & 6.89 & 7.32\\
\textit{Time} & \multicolumn{2}{c}{19.83} & \multicolumn{2}{c}{43.68} &
\multicolumn{2}{c}{102.38} & \multicolumn{2}{c}{211.11}\\
\textit{95\% c.i.} & 19.41 & 20.27 & 42.67 & 44.74 & 99.58 & 105.33 & 204.78 &
217.86
\end{tabular}
\medskip\newline$\square$

\noindent\textbf{Results from Beijing 2008 and predictions for London
2012\medskip}\newline The actual times were all inside the appropriate
confidence intervals, apart from the 200m:\medskip%
\[%
\begin{tabular}
[c]{lllll}
& \textit{200m} & \textit{400m} & \textit{800m} & \textit{1500m}\\
\textit{Time} & \multicolumn{1}{r}{19.30} & \multicolumn{1}{r}{43.75} &
\multicolumn{1}{r}{104.65} & \multicolumn{1}{r}{212.94}%
\end{tabular}
\]
You might like to try modifying the model by including the Beijing results and
trying to make predictions for London, which is 50 feet above sea level.You
should get the following hostages to fortune:%

\begin{tabular}
[c]{lcccccccc}
& \multicolumn{2}{c}{\textit{200m}} & \multicolumn{2}{c}{\textit{400m}} &
\multicolumn{2}{c}{\textit{800m}} & \multicolumn{2}{c}{\textit{1500m}}\\
\textit{Time} & \multicolumn{2}{c}{19.83} & \multicolumn{2}{c}{43.71} &
\multicolumn{2}{c}{102.58} & \multicolumn{2}{c}{211.40}\\
\textit{95\% c.i.} & 19.39 & 20.28 & 42.66 & 44.80 & 99.70 & 105.63 & 204.88 &
218.35
\end{tabular}
\medskip\newline

Out of interest, the world record times in these events are 19.30, 43.18,
101.11 and 206.00 seconds respectively.


\end{document}