% \documentclass[serif]{beamer} % Serif for Computer Modern math font.
\documentclass[serif, handout]{beamer} % Handout mode to ignore pause statements
\hypersetup{colorlinks,linkcolor=,urlcolor=red}


\usefonttheme{serif} % Looks like Computer Modern for non-math text -- nice!
\setbeamertemplate{navigation symbols}{} % Suppress navigation symbols
% \usetheme{Berlin} % Displays sections on top
\usetheme{Frankfurt}  % Displays section titles on top: Fairly thin but still swallows some material at bottom of crowded slides
%\usetheme{Berkeley}

\usepackage[english]{babel}
\usepackage{amsmath} % for binom
% \usepackage{graphicx} % To include pdf files!
% \definecolor{links}{HTML}{2A1B81}
% \definecolor{links}{red}
\setbeamertemplate{footline}[frame number] 

\mode<presentation>

\title{Double Measurement Regression\footnote{See last slide for copyright information.}}
\subtitle{STA431 Winter/Spring 2015}
\date{} % To suppress date

\begin{document}

\begin{frame}
  \titlepage
\end{frame}

\begin{frame}
\frametitle{Overview}
\tableofcontents
\end{frame}

\section{A Small Example}

\begin{frame}
\frametitle{Seeking identifiability}
%\framesubtitle{} 
We have seen that in simple regression, parameters of a model with measurement error are not identifiable. \pause
\begin{eqnarray*}
    Y_i &=& \alpha_1 + \beta_1 X_i + \epsilon_i  \\
    W_i &=& \nu +  X_i + e_i,  
\end{eqnarray*} \pause

\vspace{5mm}

  \begin{itemize}
    \item For example, $X$ might be number of acres planted and $Y$ might be crop yield. \pause
    \item Plan the statistical analysis in advance. \pause
    \item Take 2 independent measurements of the explanatory variable.
    \item Say, farmer's report and aerial photograph.
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Double measurement}
\framesubtitle{Of the explanatory variable} 
\begin{center}
\includegraphics[width=3in]{DoublePath1}
\end{center}
\end{frame}


\begin{frame}
\frametitle{Model}
%\framesubtitle{Could have written this down based on the path diagram } 
Independently for $i=1, \ldots, n$, let
\begin{eqnarray*}
    W_{i,1} &=& \nu_1 +  X_i + e_{i,1}  \\ 
    W_{i,2} &=& \nu_2 +  X_i + e_{i,2}   \\
    Y_{i~~} &=& \beta_0 + \beta_1 X_i + \epsilon_i , 
\end{eqnarray*}
where  \pause
\begin{itemize}
     \item $X_i$ is normally distributed with mean $\mu_x$ and variance $\phi>0$ \pause
     \item $\epsilon_i$ is normally distributed with mean zero and variance $\psi>0$ \pause
     \item $e_{i,1}$ is normally distributed with mean zero and variance $\omega_1>0$ \pause
     \item $e_{i,2}$ is normally distributed with mean zero and variance $\omega_2>0$ \pause
     \item $X_i, e_{i,1}, e_{i,2}$ and $\epsilon_i$ are all independent.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Does this model pass the test of the Parameter Count Rule?}
%\framesubtitle{} 
\begin{eqnarray*}
    W_{i,1} &=& \nu_1 +  X_i + e_{i,1}  \\ 
    W_{i,2} &=& \nu_2 +  X_i + e_{i,2}   \\
    Y_{i~~} &=& \beta_0 + \beta_1 X_i + \epsilon_i , 
\end{eqnarray*} \pause

\begin{itemize}
    \item[] $\boldsymbol{\theta} = (\nu_1, \nu_2, \beta_0, \mu_x, \beta_1, \phi, \psi, \omega_1, \omega_2)$: 9 parameters.  \pause
    \item Three expected values, three variances and three covariances: 9 moments. \pause
    \item Yes. There are nine moment structure equations in nine unknown parameters.  Identifiability is possible, but not guaranteed.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{What is the distribution of the sample data?}
\framesubtitle{Calculate the moments as a function of the model parameters}  \pause
The model implies that the triples $\mathbf{D}_i = (W_{i,1},W_{i,2},Y_i)^\top$ are independent multivarate normal with  \pause
\begin{displaymath}
    E(\mathbf{D}_i) = E\left( \begin{array}{c}
         W_{i,1} \\ W_{i,1} \\ Y_i 
             \end{array} \right) =
        \left( \begin{array}{c}
         \mu_1 \\ \mu_2 \\ \mu_3
             \end{array} \right) =  
      \left( \begin{array}{c}
         \mu_x+\nu_1 \\ \mu_x+\nu_2 \\\beta_0 + \beta_1\mu_x 
             \end{array} \right),
\end{displaymath} \pause
and variance covariance matrix $V(\mathbf{D}_i) = \boldsymbol{\Sigma} = $  \pause
\begin{displaymath}
    \left( \begin{array}{c c  c }
    \sigma_{11} & \sigma_{12} & \sigma_{13}  \\ 
     &  \sigma_{22}  & \sigma_{23}  \\
     &      &  \sigma_{33}  \\
        \end{array} \right) =   \pause
        \left( \begin{array}{c c c}
        \phi+\omega_1 & \phi            & \beta_1 \phi \\ 
                      & \phi+\omega_2   & \beta_1 \phi \\
                      &                 & \beta_1^2 \phi + \psi
        \end{array} \right).
\end{displaymath}
\end{frame}

\begin{frame}
\frametitle{Are the parameters in the covariance matrix identifiable?}
\framesubtitle{Six equations in five unknowns} 
\begin{displaymath}
    \left( \begin{array}{c c  c }
    \sigma_{11} & \sigma_{12} & \sigma_{13}  \\ 
     &  \sigma_{22}  & {\color{blue} \sigma_{23} } \\
     &      &  \sigma_{33}  \\
        \end{array} \right) =
        \left( \begin{array}{c c c}
        \phi+\omega_1 & \phi            & \beta_1 \phi \\ 
                      & \phi+\omega_2   & {\color{blue}  \beta_1 \phi } \\
                      &                 & \beta_1^2 \phi + \psi
        \end{array} \right).
\end{displaymath} \pause
\begin{eqnarray*}
    \phi     & = & \sigma_{12} \\ \pause
    \omega_1 & = & \sigma_{11} - \sigma_{12} \\ \pause
    \omega_2 & = & \sigma_{22} - \sigma_{12} \\ \pause
    \beta_1  & = & \frac{\sigma_{13}}{\sigma_{12}} \\ \pause
    \psi     & = & \sigma_{33} - \beta_1^2 \phi = \sigma_{33} - \frac{\sigma_{13}^2}{\sigma_{12}}
\end{eqnarray*} 
Yes.
\end{frame}

\begin{frame}
\frametitle{What about the expected values?} \pause
%\framesubtitle{} 
Model equations again:
{\scriptsize
\begin{eqnarray*}
    W_{i,1} &=& \nu_1 +  X_i + e_{i,1}  \\ 
    W_{i,2} &=& \nu_2 +  X_i + e_{i,2}   \\
    Y_{i~~} &=& \beta_0 + \beta_1 X_i + \epsilon_i , 
\end{eqnarray*} \pause
} % End size
Expected values:
\begin{eqnarray*} 
    \mu_1 & = & \nu_1+\mu_x \\ 
    \mu_2 & = & \nu_2+\mu_x  \\ 
    \mu_3 & = & \beta_0 + \beta_1\mu_x  \\ \pause
\end{eqnarray*}
{\small
Four parameters appear only in the expected values: $\nu_1,\nu_2,\mu_x,\beta_0$. \pause
  \begin{itemize}
    \item Three equations in four unknowns, even assuming $\beta_1$ known. \pause
    \item Parameter count rule applies. \pause
    \item But we don't need it because these are linear equations. \pause
    \item Re-parameterize.
  \end{itemize}
} % End size
\end{frame}

\begin{frame}

\frametitle{Re-parameterize}
\framesubtitle{$\mu_1 = \nu_1+\mu_x$ ~~~ $\mu_2 = \nu_2+\mu_x $ ~~~ $\mu_3 = \beta_0 + \beta_1\mu_x$} \pause
{\footnotesize
\begin{itemize}
    \item Absorb $\nu_1, \nu_2, \mu_x, \beta_0$ into $\boldsymbol{\mu}$. \pause
    \item Parameter was $\boldsymbol{\theta} = (\nu_1, \nu_2, \beta_0, \mu_x, ~ \beta_1, \phi, \psi, \omega_1, \omega_2)$
    \item Now it's $\boldsymbol{\theta} = (\mu_1, \mu_2, \mu_3, ~ \beta_1, \phi, \psi, \omega_1, \omega_2)$. \pause
    \item Dimension of the parameter space is now one less. \pause
    \item We haven't lost much.
    \item Especially because the model was already re-parameterized. \pause
    \item Of course there is measurement error in $Y$. Recall \pause
  \end{itemize}
\begin{eqnarray*}
Y &=& \alpha + \beta X  + \epsilon \\
V &=& \nu_0 +  Y +  e \\ \pause
  &=& \nu_0 + (\alpha + \beta X  + \epsilon) +  e \\ \pause
  &=& (\nu_0 + \alpha) + \beta X  
                                + (\epsilon +  e) \\ \pause
  &=& \beta_0 +  \beta X  + \epsilon^\prime  
\end{eqnarray*} 
} % End size
\end{frame}

\begin{frame}
\frametitle{Re-parameterization}
%\framesubtitle{} 
\begin{itemize}
    \item Re-parameterization makes maximum likelihood possible. \pause
    \item Otherwise the maximum is not unique and it's a mess. \pause
    \item Estimate $\boldsymbol{\mu}$ with $\overline{\mathbf{D}}$ and it simply disappears  \pause from 
{\footnotesize
\begin{displaymath}
    L(\boldsymbol{\mu,\Sigma}) = 
    |\boldsymbol{\Sigma}|^{-n/2} (2\pi)^{-np/2} 
    \exp -\frac{n}{2}\left\{ tr(\boldsymbol{\widehat{\Sigma}\Sigma}^{-1}) +
    (\overline{\mathbf{D}}-\boldsymbol{\mu})^\top \boldsymbol{\Sigma}^{-1} 
    (\overline{\mathbf{D}}-\boldsymbol{\mu}) \right\}
\end{displaymath} \pause
} % End size
    \item This step is so common it becomes silent. \pause
    \item Model equations are often written in centered form.  \pause
    \item It's more compact, and calculation of the covariance matrix is easier.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Back to the covariance structure equations}
%\framesubtitle{Six equations in five unknowns} 
\begin{displaymath}
    \left( \begin{array}{c c  c }
    \sigma_{11} & \sigma_{12} & \sigma_{13}  \\ 
     &  \sigma_{22}  & \sigma_{23}  \\
     &      &  \sigma_{33}  \\
        \end{array} \right) =
        \left( \begin{array}{c c c}
        \phi+\omega_1 & \phi            & \beta_1 \phi \\ 
                      & \phi+\omega_2   & \beta_1 \phi  \\
                      &                 & \beta_1^2 \phi + \psi
        \end{array} \right).
\end{displaymath} \pause

\begin{itemize}
     \item Notice that the model dictates $\sigma_{1,3}=\sigma_{2,3}$.  \pause
     \item There are two ways to solve for $\beta_1$: \\ \pause
$\beta_1=\frac{\sigma_{13}}{\sigma_{12}}$ and $\beta_1=\frac{\sigma_{23}}{\sigma_{12}}$. \pause
     \item Does this mean the solution for $\beta_1$ is not ``unique?" \pause
     \item No; everything is okay. Because $\sigma_{1,3}=\sigma_{2,3}$, the two solutions are actually the same. \pause
     \item If a parameter can be recovered from the moments in any way at all, it is identifiable. 
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Testing goodness of fit.} \pause
%\framesubtitle{Six equations in five unknowns} 
\begin{displaymath}
    \left( \begin{array}{c c  c }
    \sigma_{11} & \sigma_{12} & \sigma_{13}  \\ 
     &  \sigma_{22}  & \sigma_{23}  \\
     &      &  \sigma_{33}  \\
        \end{array} \right) =
        \left( \begin{array}{c c c}
        \phi+\omega_1 & \phi            & \beta_1 \phi \\ 
                      & \phi+\omega_2   & \beta_1 \phi  \\
                      &                 & \beta_1^2 \phi + \psi
        \end{array} \right).
\end{displaymath} \pause

\begin{itemize}
     \item $\sigma_{1,3}=\sigma_{2,3}$ is a \emph{model-induced constraint} upon $\boldsymbol{\Sigma}$. \pause
     \item It's a testable null hypothesis.  \pause
     \item If rejected, the model is called into question. \pause
     \item Likelihood ratio test comparing this model to a completely unrestricted multivariate normal model: \pause
{\footnotesize
\begin{displaymath}
G^2 = -2\ln\frac{ L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right) } {L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}})}
\end{displaymath} \pause
} % End size
     \item It's $n$ times the SAS ``objective function" at the MLE. \pause
     \item A likelihood ratio test for goodness of fit.   \pause
     \item Valuable even if the data are not normal.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{The Reproduced Covariance Matrix}
%\framesubtitle{} 
\begin{itemize}
    \item $\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})$ is called the \emph{reproduced covariance matrix}. \pause
    \item It is the covariance matrix of the observable data, written as a function of the model parameters and evaluated at the MLE. \pause
\begin{displaymath}
    \boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}}) = 
        \left( \begin{array}{c c c}
        \widehat{\phi}+\widehat{\omega}_1 & \widehat{\phi} & \widehat{\beta}_1 \widehat{\phi} \\ 
                      & \widehat{\phi}+\widehat{\omega}_2  & \widehat{\beta}_1 \widehat{\phi} \\
                      &                 & \widehat{\beta}_1^2 \widehat{\phi} + \widehat{\psi}
        \end{array} \right)
\end{displaymath}  \pause
    \item The reproduced covariance matrix obeys all model-induced constraints, while $\widehat{\boldsymbol{\Sigma}}$ does not. \pause
    \item But if the model is right they should be close. \pause
    \item This is a way to think about the likelihood ratio test for goodness of fit. 
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{General pattern for testing goodness of fit}
\framesubtitle{Usually works}  \pause
\begin{itemize}
    \item Suppose there are $k$ moment structure equations in $p$ parameters, and all the parameters are identifiable. \pause
    \item If $p<k$, call the parameter vector \emph{over-identifiable}.  \pause
    \item Only needed $p$ equations to solve for $\boldsymbol{\theta}$. \pause
    \item Substituting the solutions (in terms of $\sigma_{ij}$) back into the unused equations would yield $k-p$ equality constraints on $\boldsymbol{\Sigma}$. \pause
    \item Test those constraints with $G^2 = -2\ln\frac{ L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right) } {L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}})}$. \pause
    \item $df=k-p$ \pause
    \item Don't need to actually derive the constraints -- just count them.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{With the same number of equations and parameters} \pause
%\framesubtitle{} 
\begin{itemize}
    \item If the parameter is identifiable, call it \emph{just identifiable}. \pause
    \item Parameters are 1-1 with those of an unrestricted multivariate normal. \pause
    \item Call the model ``saturated." \pause
    \item There are no equality constraints on $\boldsymbol{\Sigma}$. \pause
    \item No likelihood ratio test ($G^2 = -2\ln\frac{ L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right) } {L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}})} = 0$). \pause
    \item This is what happens in regression with all observed variables. 
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{How to proceed} \pause
%\framesubtitle{} 
\begin{itemize}
    \item Verify identifiability. \pause
    \item If the model is over-identified, test goodness of fit. \pause
    \item If it passes (non-significant), \pause proceed. \pause
    \item Now think of your model as a ``full," or unrestricted model. \pause
    \item Compared to some (even more) reduced model that is restricted by a null hypothesis like $\beta_1=0$.  \pause
    \item Fit the reduced model. \pause
    \item Subtract goodness of fit ($G^2$ or ``chi-square") statistics to test $H_0$.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Subtract goodness of fit  statistics} \pause
%\framesubtitle{} 
$G^2$ tests the full model against the saturated model, and $G^2_0$ tests the reduced model against the saturated model. \pause
\begin{eqnarray*}
    G^2_0-G^2 & = & -2\ln\frac{ L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}}_0)\right) } {L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}})} - -2\ln\frac{ L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right) } {L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}})} \\ \pause
    & = & -2 \left(  \ln L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}}_0)\right)
{ \color{blue} - \ln L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}}) }
- \ln L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right) \right. \\ 
&& ~~~~~ \left.
{ \color{blue} + \ln L(\overline{\mathbf{D}},\widehat{\boldsymbol{\Sigma}}) }
\right) \\ \pause
 & = & -2 \ln \frac{L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}}_0)\right)}
                   {L\left(\overline{\mathbf{D}},\boldsymbol{\Sigma}(\widehat{\boldsymbol{\theta}})\right)}
\end{eqnarray*}  \pause
If the software gives you $\frac{n-1}{n}G^2$, use that.
\end{frame}

\begin{frame}
\frametitle{Further comments}
%\framesubtitle{} 
\begin{itemize}
    \item Models with non-identifiable parameters can imply testable equality constraints, but testing them is not automatic.  \pause
    \item Models can imply \emph{inequality} constraints on $\boldsymbol{\Sigma}$, too. \pause 
    \item Recall the solutions  \pause
\end{itemize}
\begin{eqnarray*}
    \phi     & = & \sigma_{12} \\
    \omega_1 & = & \sigma_{11} - \sigma_{12} \\
    \omega_2 & = & \sigma_{22} - \sigma_{12} \\
    \beta_1  & = & \frac{\sigma_{13}}{\sigma_{12}} \\
    \psi     & = & \sigma_{33} - \beta_1^2 \phi = \sigma_{33} - \frac{\sigma_{13}^2}{\sigma_{12}}
\end{eqnarray*}  \pause
We get four inequality constraints.
\end{frame}

\begin{frame}
\frametitle{Four inequality constraints on $\boldsymbol{\Sigma}$} \pause
%\framesubtitle{} 

\begin{displaymath}
    \left( \begin{array}{c c  c }
    \sigma_{11} & \sigma_{12} & \sigma_{13}  \\ 
     &  \sigma_{22}  & \sigma_{23}  \\
     &      &  \sigma_{33}  \\
        \end{array} \right) =
        \left( \begin{array}{c c c}
        \phi+\omega_1 & \phi            & \beta_1 \phi \\ 
                      & \phi+\omega_2   & \beta_1 \phi  \\
                      &                 & \beta_1^2 \phi + \psi
        \end{array} \right).
\end{displaymath} \pause

\begin{eqnarray*}
    \phi     & = & \sigma_{12} > 0\\ \pause
    \omega_1 & = & \sigma_{11} - \sigma_{12} > 0\\ \pause
    \omega_2 & = & \sigma_{22} - \sigma_{12} > 0\\ \pause
    \psi     & = &  \sigma_{33} - \frac{\sigma_{13}^2}{\sigma_{12}} > 0
\end{eqnarray*} 

\end{frame}


\begin{frame}
\frametitle{Inequality constraints} \pause
%\framesubtitle{} 

\begin{itemize}
    \item Inequality constraints arise because variances are positive. \pause
    \item Or more generally, covariance matrices are positive definite. \pause
    \item Could inequality constraints be violated in numerical maximum likelihood? \pause
    \item Definitely. \pause
    \item But only a little by sampling error if the model is correct. \pause
    \item So maybe it's not so dumb to test hypotheses like $H_0: \omega_1=0$. \pause
    \item Since the model says $\omega_1 = \sigma_{11} - \sigma_{12}$.
\end{itemize}
\end{frame}

\begin{frame}[fragile]
\frametitle{Little SAS Example}
%\framesubtitle{} 
{\scriptsize
\begin{verbatim}
/**************************** Babydouble.sas ****************************/
title 'Simple double measurement with proc calis';
title2 'Jerry Brunner:  Student Number 999999999';

data baby;
     infile '/folders/myfolders/431s15/Babydouble.data.txt'
             firstobs=2; 
     input id W1 W2 Y;

proc calis pcorr vardef=n; 
     /* See reproduced covariance matrix, 
        Use true MLE and get exact G^2 */
	 title3 'Fit the centered model';
     var W1 W2 Y;  /* Declare observed variables */
     lineqs        /* Model equations, separated by commas. */
        Y  = beta1*F + epsilon, /* Latent variables begin with the letter F */
        W1 = F + e1,
        W2 = F + e2;
     variance      /* Declare variance parameters. */
        F = phi, epsilon = psi, e1=omega1, e2=omega2;
\end{verbatim} 
} % End size
\end{frame}

\begin{frame}
\frametitle{Results}

Click \href{http://www.utstat.toronto.edu/~brunner/431s15/lectures/BabydoubleResults.html}{\textbf{Here}}
for the output. This link will probably be broken once the term is over. See the course website for another route to the output file:

\vspace{10mm}

\href{http://www.utstat.toronto.edu/~brunner/oldclass/431s15} {\small\texttt{http://www.utstat.toronto.edu/$^\sim$brunner/oldclass/431s15}}

\end{frame}



\section{The general model}

\begin{frame}
\frametitle{An extension of the double measurement design} \pause
%\framesubtitle{} 
Double measurement can help solve a big problem: 
Correlated measurement error. \pause
\begin{center}
\includegraphics[width=3in]{DoublePath2}
\end{center} \pause
\begin{itemize}
    \item The main idea is that $\mathbf{X}$ and $\mathbf{Y}$ are each measured twice, perhaps at different times using different methods.   \pause
    \item Measurement errors may be correlated within sets  but not between sets. 
\end{itemize} 
\end{frame}


\begin{frame}
\frametitle{Double Measurement Regression: A Two-Stage Model}\pause
% \framesubtitle{}  
{\LARGE
\begin{eqnarray*}
    \mathbf{Y}_i &=& \boldsymbol{\beta}_0 + \boldsymbol{\beta}_1 \mathbf{X}_i +   \boldsymbol{\epsilon}_i \\ \pause
    \mathbf{F}_i &=& \left( \begin{array}{c}
         \mathbf{X}_i  \\ \mathbf{Y}_i 
             \end{array} \right)  \\ \pause
    \mathbf{D}_{i,1} &=& \boldsymbol{\nu}_1 + \mathbf{F}_i + \mathbf{e}_{i,1} \\ \pause
    \mathbf{D}_{i,2} &=& \boldsymbol{\nu}_2 + \mathbf{F}_i + \mathbf{e}_{i,2}  
    \end{eqnarray*}  \pause
} % End size

Observable variables are $\mathbf{D}_{i,1}$ and $\mathbf{D}_{i,2}$: both are $(p+q) \times 1$. \pause


\vspace{2mm}

$E(\mathbf{X}_i) = \boldsymbol{\mu}_x$, $V(\mathbf{X}_i)=\boldsymbol{\Phi}_x$, $V(\boldsymbol{\epsilon}_i)=\boldsymbol{\Psi}$, $V(\mathbf{e}_{i,1})=\boldsymbol{\Omega}_1$, $V(\mathbf{e}_{i,2})=\boldsymbol{\Omega}_2$. Also, $\mathbf{X}_i$, $\boldsymbol{\epsilon}_i$, $\mathbf{e}_{i,1}$ and $\mathbf{e}_{i,2}$ are independent.

\end{frame}

\begin{frame}
\frametitle{Measurement errors may be correlated}
\framesubtitle{Look at the measurement model}  \pause
\begin{eqnarray*}
\mathbf{F}_i &=& \left( \begin{array}{c}
         \mathbf{X}_i  \\ \mathbf{Y}_i 
             \end{array} \right)  \\
\mathbf{D}_{i,1} &=& \boldsymbol{\nu}_1 + \mathbf{F}_i + \mathbf{e}_{i,1} \\
\mathbf{D}_{i,2} &=& \boldsymbol{\nu}_2 + \mathbf{F}_i + \mathbf{e}_{i,2} 
\end{eqnarray*}  \pause

\renewcommand{\arraystretch}{1.2}
\begin{eqnarray*}
V(\mathbf{e}_{i,1}) &=& \boldsymbol{\Omega}_1 = 
\left( \begin{array}{c|c}
\boldsymbol{\Omega}_{11} & \boldsymbol{\Omega}_{12} \\ \hline 
\boldsymbol{\Omega}_{12}^\top & \boldsymbol{\Omega}_{22} 
             \end{array} \right) \\ \pause
V(\mathbf{e}_{i,2}) &=& \boldsymbol{\Omega}_2 = 
\left( \begin{array}{c|c}
\boldsymbol{\Omega}_{33} & \boldsymbol{\Omega}_{34} \\ \hline 
\boldsymbol{\Omega}_{34}^\top & \boldsymbol{\Omega}_{44} 
             \end{array} \right)
\end{eqnarray*} 
\renewcommand{\arraystretch}{1.0}

\end{frame}


\begin{frame}
\frametitle{Expected values of the observable variables}
\framesubtitle{$\mathbf{D}_{i,1} = \boldsymbol{\nu}_1 + \mathbf{F}_i + \mathbf{e}_{i,1}$ and $\mathbf{D}_{i,2} = \boldsymbol{\nu}_2 + \mathbf{F}_i + \mathbf{e}_{i,2}$}  \pause
\begin{columns}   % Use Beamer's columns to use more of the margins!
\column{1.2\textwidth}
\begin{eqnarray*}
E(\mathbf{D}_{i,1})  &=& \left( \begin{array}{c}
     \boldsymbol{\mu}_{1,1}   \\ \boldsymbol{\mu}_{1,2} 
             \end{array} \right)    \pause
                    = \left( \begin{array}{c}
\boldsymbol{\nu}_{1,1} +  E(\mathbf{X}_i)  \\ \boldsymbol{\nu}_{1,2} + E(\mathbf{Y}_i)
             \end{array} \right)  \pause
                    = \left( \begin{array}{c}
\boldsymbol{\nu}_{1,1} +  \boldsymbol{\mu}_x  \\ \boldsymbol{\nu}_{1,2} + 
\boldsymbol{\beta}_0 + \boldsymbol{\beta}_1\boldsymbol{\mu}_x
             \end{array} \right) \\ \\   \pause
E(\mathbf{D}_{i,2})  &=& \left( \begin{array}{c}
     \boldsymbol{\mu}_{2,1}   \\ \boldsymbol{\mu}_{2,2} 
             \end{array} \right)  
                    = \left( \begin{array}{c}
\boldsymbol{\nu}_{2,1} +  E(\mathbf{X}_i)  \\ \boldsymbol{\nu}_{2,2} + E(\mathbf{Y}_i)
             \end{array} \right)
                    = \left( \begin{array}{c}
\boldsymbol{\nu}_{2,1} +  \boldsymbol{\mu}_x  \\ \boldsymbol{\nu}_{2,2} + 
\boldsymbol{\beta}_0 + \boldsymbol{\beta}_1\boldsymbol{\mu}_x
             \end{array} \right) \end{eqnarray*}   \pause
\end{columns}

\vspace{5mm}

\begin{itemize}
    \item $\boldsymbol{\nu}_1$, $\boldsymbol{\nu}_2$, $\boldsymbol{\beta}_0$ and $\boldsymbol{\mu}_x$  parameters appear only in expected value, not covariance matrix.  \pause
    \item $\mathbf{X}_i$ is $p \times 1$ and $\mathbf{Y}_i$ is $q \times 1$.   \pause
    \item Even with knowledge of $\boldsymbol{\beta}_0$, $2(p+q)$ equations in $3(p+q)$ unknown parameters.  \pause
    \item Identifying the expected values and intercepts is hopeless.  \pause
    \item Re-parameterize,  \pause swallowing them into $\boldsymbol{\mu} =   \pause
    E\left( \begin{array}{c}
         \mathbf{D}_{i,1}  \\ \mathbf{D}_{i,2} 
             \end{array} \right)$.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Stage One: The latent variable model}
% \framesubtitle{} 
\begin{eqnarray*}
\mathbf{Y}_i &=& \boldsymbol{\beta}_0 + \boldsymbol{\beta}_1 \mathbf{X}_i +   \boldsymbol{\epsilon}_i \\
\mathbf{F}_i &=& \left( \begin{array}{c}
         \mathbf{X}_i  \\ \mathbf{Y}_i 
             \end{array} \right)  
\end{eqnarray*}   \pause

$V(\mathbf{X}_i)=\boldsymbol{\Phi}_x$, $V(\boldsymbol{\epsilon}_i)=\boldsymbol{\Psi}$,   \pause
$\mathbf{X}_i$ and $\boldsymbol{\epsilon}_i$ are independent.   \pause

Proving identifiability, \ldots   \pause
\begin{displaymath}
    V(\mathbf{F}_i)=\boldsymbol{\Phi} =   \pause
\left( \begin{array}{c c}
         \boldsymbol{\Phi}_{11}  & \boldsymbol{\Phi}_{12} \\
         \boldsymbol{\Phi}_{12}^\top  & \boldsymbol{\Phi}_{22} \\
             \end{array} \right)  \pause
= 
\left( \begin{array}{c c}
\boldsymbol{\Phi}_x  & \boldsymbol{\Phi}_x\boldsymbol{\beta}_1^\top \\
\boldsymbol{\beta}_1\boldsymbol{\Phi}_x  & 
\boldsymbol{\beta}_1\boldsymbol{\Phi}_x\boldsymbol{\beta}_1^\top
+ \boldsymbol{\Psi}
             \end{array} \right)
\end{displaymath}  \pause

$\boldsymbol{\Phi}_x$, $\boldsymbol{\beta}_1$ and $\boldsymbol{\Psi}$ can be recovered from $\boldsymbol{\Phi}$.

\end{frame}

\begin{frame}
\frametitle{Stage Two: The measurement model}  \pause
%\framesubtitle{} 
\begin{eqnarray*}
\mathbf{D}_{i,1} &=& \boldsymbol{\nu}_1 + \mathbf{F}_i + \mathbf{e}_{i,1} \\
\mathbf{D}_{i,2} &=& \boldsymbol{\nu}_2 + \mathbf{F}_i + \mathbf{e}_{i,2} 
\end{eqnarray*}  \pause

$V(\mathbf{e}_{i,1})=\boldsymbol{\Omega}_1$, $V(\mathbf{e}_{i,2})=\boldsymbol{\Omega}_2$. Also, $\mathbf{F}_i$,  $\mathbf{e}_{i,1}$ and $\mathbf{e}_{i,2}$ are independent.  \pause

\begin{displaymath}
    \boldsymbol{\Sigma} = V\left( \begin{array}{c}
\mathbf{D}_{i,1}  \\ \mathbf{D}_{i,2}
                                \end{array} \right)   \pause
=
\left( \begin{array}{c c}
\boldsymbol{\Phi}+\boldsymbol{\Omega}_1  & \boldsymbol{\Phi} \\
\boldsymbol{\Phi}  & \boldsymbol{\Phi}+\boldsymbol{\Omega}_2
             \end{array} \right)
\end{displaymath}  \pause

$\boldsymbol{\Phi}$, $\boldsymbol{\Omega}_1$ and $\boldsymbol{\Omega}_2$ can easily be recovered from $\boldsymbol{\Sigma}$. 
\end{frame}

\begin{frame}
\frametitle{All the parameters in the covariance matrix are identifiable} \pause
%\framesubtitle{} 
\begin{itemize}
    \item $\boldsymbol{\Phi}_x$, $\boldsymbol{\beta}_1$ and $\boldsymbol{\Psi}$ can be recovered from $\boldsymbol{\Phi} = V(\mathbf{F}_i)$. \pause
    \item $\boldsymbol{\Phi}$, $\boldsymbol{\Omega}_1$ and $\boldsymbol{\Omega}_2$ can be recovered from $\boldsymbol{\Sigma} = V\left( \begin{array}{c}
\mathbf{D}_{i,1}  \\ \mathbf{D}_{i,2}
                                \end{array} \right)$. \pause
    \item[] 
    \item Correlated measurement error within sets is allowed. \pause
    \item  This is a big plus, because it's reality. \pause
    \item Correlated measurement error between sets must be ruled out by careful data collection. \pause
    \item No need to do the calculations ever again.
\end{itemize}
\end{frame}

\section{The BMI study}

\begin{frame}
\frametitle{The BMI Health Study} \pause
%\framesubtitle{} 
\begin{itemize}
    \item Body Mass Index: Weight in Kilograms divided by Height in Meters Squared. \pause
    \item Under 18 means underweight, Over 25 means overweight, Over 30 means obese. \pause
    \item High BMI is associated with poor health, like high blood pressure and high cholesterol. \pause
    \item People with high BMI tend to be older and fatter. \pause
    \item \emph{But}, what if you have a high BMI but are in good physical shape (low percent body fat)?
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{The Question}
%\framesubtitle{} 
\begin{itemize}
    \item If you control for age and percent body fat, is BMI still associated with indicators for poor health? \pause
    \item But percent body fat (and to a lesser extent, age) are measured with error.  Standard ways of controlling for them with ordinary regression are highly suspect. \pause
    \item Use the double measurement design.
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{True variables (all latent)}
%\framesubtitle{} 
\begin{itemize}
    \item $X_1$ = Age
    \item $X_2$ = BMI
    \item $X_3$ = Percent body fat
    \item $Y_1$ = Cholesterol
    \item $Y_2$ = Diastolic blood pressure
\end{itemize}
\end{frame}

\begin{frame}
\frametitle{Measure twice with different personnel at different locations and by different methods} \pause
%\framesubtitle{} 
\begin{columns}   % Use Beamer's columns to use more of the margins!
\column{1.1\textwidth}
{\footnotesize
\renewcommand{\arraystretch}{1.2}
\begin{tabular}{lll}  \hline
              & Measurement Set One &  Measurement Set Two  \\ \hline
Age           & Self report         &  Passport or birth certificate  \\
BMI           & Dr. Office measurements & Lab technician, no shoes, gown \\
\% Body Fat   & Tape and calipers, Dr. Office & Submerge in water tank \\
Cholesterol   & Lab 1 & Lab 2 \\
Diastolic BP  & Blood pressure cuff, Dr. office & Digital readout, mostly automatic \\ \hline
\end{tabular} \pause
\renewcommand{\arraystretch}{1.0}

\vspace{5mm}

\begin{itemize}
    \item Set two is of generally higher quality. \pause
    \item Correlation of measurement errors is unlikely between sets.
\end{itemize}
} % End size
\end{columns}
\end{frame}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\begin{frame}
\frametitle{Copyright Information}

This slide show was prepared by  \href{http://www.utstat.toronto.edu/~brunner}{Jerry Brunner},
Department of Statistical Sciences, University of Toronto. It is licensed under a 
\href{http://creativecommons.org/licenses/by-sa/3.0/deed.en_US}
     {Creative Commons Attribution - ShareAlike 3.0 Unported License}. Use any part of it as you like and share the result freely. The \LaTeX~source code is available from the course website:
\href{http://www.utstat.toronto.edu/~brunner/oldclass/431s13} {\small\texttt{http://www.utstat.toronto.edu/$^\sim$brunner/oldclass/431s31}}

\end{frame}


\end{document}


# Simulate some data for a simple example
set.seed(9999)
n = 150; beta1 = 1; phi=1; psi = 9; omega1=1; omega2=1
X = rnorm(n,10,sqrt(phi)); epsilon = rnorm(n,0,sqrt(psi))
e1 = rnorm(n,0,sqrt(omega1)); e2 = rnorm(n,0,sqrt(omega2))
Y = 3 + beta1*X + epsilon
W1 = X + e1; W2 = X + e2
datta = round(cbind(W1,W2,Y),2)
cor(datta)







summary(lm(Y~X))
summary(lm(Y~W1+W2))
rnorm(n,0,sqrt())

########### Nicely ambiguous traditional output ###########

> cor(datta)
          W1        W2         Y
W1 1.0000000 0.5748331 0.1714324
W2 0.5748331 1.0000000 0.1791539
Y  0.1714324 0.1791539 1.0000000

> summary(lm(Y~X))

Call:
lm(formula = Y ~ X)

Residuals:
    Min      1Q  Median      3Q     Max 
-6.8778 -2.0571 -0.0718  2.1200  7.4284 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   3.8122     2.5348   1.504 0.134723    
X             0.9288     0.2521   3.684 0.000322 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 3.096 on 148 degrees of freedom
Multiple R-squared: 0.08398,	Adjusted R-squared: 0.07779 
F-statistic: 13.57 on 1 and 148 DF,  p-value: 0.0003218 

> summary(lm(Y~W1+W2))

Call:
lm(formula = Y ~ W1 + W2)

Residuals:
    Min      1Q  Median      3Q     Max 
-7.6711 -2.3889 -0.1303  2.3442  7.6879 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   7.9647     2.1148   3.766 0.000239 ***
W1            0.2369     0.2282   1.038 0.300870    
W2            0.2799     0.2300   1.217 0.225615    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 3.182 on 147 degrees of freedom
Multiple R-squared: 0.03918,	Adjusted R-squared: 0.02611 
F-statistic: 2.997 on 2 and 147 DF,  p-value: 0.053 

> 
> W = W1+W2
> summary(lm(Y~W))

Call:
lm(formula = Y ~ W)

Residuals:
    Min      1Q  Median      3Q     Max 
-7.6787 -2.3707 -0.1431  2.3242  7.6763 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   7.9721     2.1066   3.784 0.000223 ***
W             0.2583     0.1052   2.454 0.015280 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 3.171 on 148 degrees of freedom
Multiple R-squared: 0.0391,	Adjusted R-squared: 0.03261 
F-statistic: 6.023 on 1 and 148 DF,  p-value: 0.01528