% \documentclass[serif]{beamer} % Serif for Computer Modern math font.
\documentclass[serif, handout]{beamer} % Handout 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{More Linear Algebra\footnote{See Chapter 2 of \emph{Linear models in statistics} for more detail. This slide show is an open-source document. See last slide for copyright information.}}
\subtitle{STA 302: Fall 2015}
\date{} % To suppress date

\begin{document}
% The text uses parentheses instead of bracket notation, but I'm sticking with brackets; it's clearer. Warn the class.

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

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

\section{Things you already know}

\begin{frame}
\frametitle{You already know about}

  \begin{itemize}
    \item Matrices $\mathbf{A} = [a_{ij}]$ \pause
    \item Matrix addition and subtraction $\mathbf{A}+\mathbf{B} = [a_{ij}+b_{ij}]$ \pause
    \item Scalar multiplication $a \mathbf{B} = [a\,b_{ij}]$ \pause
    \item Matrix multiplication 
          $\mathbf{AB}  = \left[ \displaystyle{\sum_k}a_{ik}b_{kj} \right]$ \pause
    \item Inverse $\mathbf{A}^{-1}\mathbf{A} = \mathbf{AA}^{-1} = \mathbf{I}$ \pause
    \item Transpose $\mathbf{A}^\prime = [a_{ji}]$ \pause
    \item Symmetric matrices $\mathbf{A} = \mathbf{A}^\prime$ \pause
    \item Determinants
    \item Linear independence
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Three mistakes that will get you a zero}
\framesubtitle{Numbers are $1 \times 1$ matrices, but larger matrices are not just numbers.} \pause
You will get a zero if you   \pause

\vspace{2mm}  

\begin{itemize}
    \item Write $\mathbf{AB} = \mathbf{BA}$. It's not true in general. \pause
    \item Write $\mathbf{A}^{-1}$ when $\mathbf{A}$ is not a square matrix. The inverse is not even defined. \pause
    \item Represent the inverse of a matrix (even if it exists) by writing it in the denominator,  \pause like 
    $\mathbf{a}^\prime \mathbf{B}^{-1}\mathbf{a} = \frac{\mathbf{a}^\prime \mathbf{a}}{\mathbf{B}}$.  \pause Matrices are not just numbers.  \pause
\end{itemize}

\vspace{5mm}
If you commit one of these crimes, the mark for the question (or part of a question, like 3c) is zero. The rest of your answer will be ignored.

\end{frame}

\begin{frame}
\frametitle{Half marks off, at least}
%\framesubtitle{At least} 

You will lose \emph{at least} half marks for writing a product like $\mathbf{AB}$ when the number of colmns in $\mathbf{A}$ does not equal the number of rows in $\mathbf{B}$.

\end{frame}

\begin{frame}
\frametitle{Linear combination of vectors} \pause
%\framesubtitle{} 
Let $\mathbf{x}_1, \ldots, \mathbf{x}_p$ be $n \times 1$ vectors and $a_1, \ldots, a_p$ be scalars. \pause A \emph{linear combination} is \pause
\renewcommand{\arraystretch}{2.0}
\begin{displaymath}
\begin{array}{ccccccccc} % 9 columns
\mathbf{c} &=& ~~a_1 \mathbf{x}_1 &+& ~~a_2 \mathbf{x}_2 &+& \cdots &+& ~~a_p \mathbf{x}_p \\ \pause
 &=& a_1 \left( \begin{array}{c}
          x_{11} \\ x_{21} \\ \vdots \\ x_{n1}
          \end{array}\right) &+&
       a_2 \left( \begin{array}{c}
          x_{12} \\ x_{22} \\ \vdots \\ x_{n2}
          \end{array}\right) &+& \cdots &+&
       a_p \left( \begin{array}{c}
          x_{1p} \\ x_{2p} \\ \vdots \\ x_{np}
          \end{array}\right)
\end{array}
\end{displaymath}
\renewcommand{\arraystretch}{1.0}
\end{frame}

\begin{frame}
\frametitle{Linear independence} \pause
%\framesubtitle{} 
A set of vectors $\mathbf{x}_1, \ldots, \mathbf{x}_p$ is said to be \emph{linearly dependent}  \pause if there is a set of scalars $a_1, \ldots, a_p$, not all zero, with  \pause
\begin{displaymath}
    a_1 \left( \begin{array}{c}
          x_{11} \\ x_{21} \\ \vdots \\ x_{n1}
          \end{array}\right) + 
       a_2 \left( \begin{array}{c}
          x_{12} \\ x_{22} \\ \vdots \\ x_{n2}
          \end{array}\right) + \cdots +
       a_p \left( \begin{array}{c}
          x_{1p} \\ x_{2p} \\ \vdots \\ x_{np}
          \end{array}\right) = 
           \left( \begin{array}{c}
          0 \\ 0 \\ \vdots \\ 0
          \end{array}\right)
\end{displaymath} \pause

\vspace{2mm}
If no such constants $a_1, \ldots, a_p$ exist, the vectors are linearly independent.  \pause That is, 

\vspace{2mm}
If $a_1 \mathbf{x}_1 + a_2 \mathbf{x}_2 + \cdots + a_p \mathbf{x}_p = \mathbf{0}$ implies $a_1 = a_2 \cdots = a_p = 0$, then the vectors are said to be \emph{linearly independent}.
\end{frame}

\begin{frame}
\frametitle{Bind the vectors $\mathbf{x}_1, \ldots, \mathbf{x}_p$ into a matrix} \pause
% \framesubtitle{Yields a more convenient definition of linear independence} 
{\footnotesize
%\renewcommand{\arraystretch}{2.0}
\begin{displaymath}
\begin{array}{ccccccccc} % 9 columns
 & & a_1 \mathbf{x}_1~~ &+& a_2 \mathbf{x}_2~~ &+& \cdots &+& a_p \mathbf{x}_p~~ \\
&&&&&& \\ \pause
 &=& \left( \begin{array}{c}
          x_{11} \\ x_{21} \\ \vdots \\ x_{n1}
          \end{array}\right)a_1  &+&
       \left( \begin{array}{c}
          x_{12} \\ x_{22} \\ \vdots \\ x_{n2}
          \end{array}\right)a_2  &+& \cdots &+&
       \left( \begin{array}{c}
          x_{1p} \\ x_{2p} \\ \vdots \\ x_{np}
          \end{array}\right)a_p \\
\end{array}
\end{displaymath} \pause
%\renewcommand{\arraystretch}{1.0}
\begin{eqnarray*}
     & = & \left( \begin{array}{cccc}
    x_{11} & x_{12} & \cdots & x_{1p} \\
    x_{21} & x_{22} & \cdots & x_{2p} \\
    \vdots & \vdots & \vdots & \vdots \\
    x_{n1} & x_{n2} & \cdots & n_{np} \\
\end{array}\right)
\left( \begin{array}{c}
          a_{1} \\ a_{2} \\ \vdots \\ a_{p}
          \end{array}\right) \\
     && \\ \pause
     & = & \mathbf{Xa} \hspace{82mm} % Such an ugly way to get nice alignment!
\end{eqnarray*} 
} % End size
\end{frame}

\begin{frame}
\frametitle{A more convenient definition of linear independence}
\framesubtitle{$a_1 \mathbf{x}_1 + a_2 \mathbf{x}_2 + \cdots + a_p \mathbf{x}_p = \mathbf{Xa}$}  \pause

Let $\mathbf{X}$ be an $n \times p$ matrix of constants. The columns of 
$\mathbf{X}$ are said to be \emph{linearly dependent} if there exists $\mathbf{a} \neq \mathbf{0}$ with $\mathbf{Xa} = \mathbf{0}$.  \pause
We will say that the columns of $\mathbf{X}$ are linearly \emph{independent} if $\mathbf{Xa} = \mathbf{0}$ implies $\mathbf{a} = \mathbf{0}$. \pause

\vspace{10mm}

For example, show that $\mathbf{B}^{-1}$ exists implies that the columns of $\mathbf{B}$ are linearly independent. \pause
\begin{displaymath}
    \mathbf{Ba} = \mathbf{0}  \pause \Rightarrow \mathbf{B}^{-1}\mathbf{Ba} = \mathbf{B}^{-1}\mathbf{0}  \pause \Rightarrow  \mathbf{a} = \mathbf{0}.
\end{displaymath}
\end{frame}

\begin{frame}
\frametitle{How to show $\mathbf{A}^{-1\prime}= \mathbf{A}^{\prime-1}$}  \pause
Suppose $\mathbf{B} = \mathbf{A}^{-1}$, meaning $\mathbf{AB} = \mathbf{BA} = \mathbf{I}$.  \pause Must show two things: $\mathbf{B}^\prime\mathbf{A}^\prime = \mathbf{I}$  \pause and $\mathbf{A}^\prime\mathbf{B}^\prime = \mathbf{I}$.  \pause

\begin{eqnarray*}
    \mathbf{AB} = \mathbf{I}  \pause & \Rightarrow &  \mathbf{B}^\prime\mathbf{A}^\prime = \mathbf{I}^\prime \pause = \mathbf{I} \\ \pause
    \mathbf{BA} = \mathbf{I} \pause & \Rightarrow &  \mathbf{A}^\prime\mathbf{B}^\prime  = \mathbf{I}^\prime \pause = \mathbf{I}
\end{eqnarray*} 
$\blacksquare$
\end{frame}

% \section{Trace} deleted: See STA431s13 for omitted material

\begin{frame}
\frametitle{Extras}
\framesubtitle{You may not know about these, and we may use them occasionally} 
  \begin{itemize}
    \item Trace
    \item Rank
    \item Partitioned matrices
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Trace of a square matrix}
%\framesubtitle{} 
  \begin{itemize}
    \item Sum of diagonal elements \pause
    \item Obvious: $tr(\mathbf{A}+\mathbf{B}) = tr(\mathbf{A})+tr(\mathbf{B})$ \pause
    \item Not obvious: $tr(\mathbf{AB}) = tr(\mathbf{BA})$
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Rank} \pause
%\framesubtitle{} 
  \begin{itemize}
    \item Row rank is the number of linearly independent rows.  \pause
    \item Column rank is the number of linearly independent columns. \pause
    \item Rank of a matrix is the minimum of row rank and column rank. \pause
    \item $rank(\mathbf{AB}) = \min\left(rank(\mathbf{A}),rank(\mathbf{B})\right)$.
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Partitioned matrix} \pause
%\framesubtitle{} 

  \begin{itemize}
    \item A matrix of matrices \pause
    {\LARGE
    \begin{displaymath} \left[ \begin{array}{c|c}
        \mathbf{A} & \mathbf{B} \\ \hline
        \mathbf{C} & \mathbf{D}    
    \end{array} \right] 
    \end{displaymath}  \pause
    } % End size
    \item Row by column (matrix) multiplication works, provided the matrices are the right sizes.
  \end{itemize}
\end{frame}


\section{Spectral decomposition}

\begin{frame}
\frametitle{Eigenvalues and eigenvectors} \pause
%\framesubtitle{} 
Let $\mathbf{A} = [a_{i,j}]$ be an $n \times n$ matrix, so that the following applies to square matrices. \pause $\mathbf{A}$ is said to have an \emph{eigenvalue} $\lambda$ and  (non-zero) \emph{eigenvector} $\mathbf{x}$ corresponding to $\lambda$ if \pause
\begin{displaymath} 
    \mathbf{Ax} =  \lambda\mathbf{x}.
\end{displaymath} \pause
  \begin{itemize}
    \item Eigenvalues are the $\lambda$ values that solve the determinantal equation $|\mathbf{A}-\lambda\mathbf{I}| = 0$. \pause
    \item The determinant is the product of the eigenvalues: $|\mathbf{A}| = \prod_{i=1}^n \lambda_i$
  \end{itemize}
\end{frame}


\begin{frame}
\frametitle{Spectral decomposition of symmetric matrices} \pause
%\framesubtitle{} 
The \emph{Spectral decomposition theorem} says that every square and symmetric matrix $\mathbf{A} = [a_{i,j}]$ may be written
\begin{equation*} 
   \mathbf{A} = \mathbf{CDC}^\prime,
\end{equation*} \pause
where the columns of $\mathbf{C}$ (which may also be denoted $\mathbf{x}_1, \ldots, \mathbf{x}_n$) are the eigenvectors of $\mathbf{A}$, \pause
and the diagonal matrix $\mathbf{D}$ contains the corresponding eigenvalues.   \pause
\begin{displaymath} 
   \mathbf{D} = 
    \left( \begin{array}{c c c c }
        \lambda_1 & 0          & \cdots  & 0            \\
            0     & \lambda_2  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & \lambda_n    \\ 
\end{array} \right) 
\end{displaymath} \pause
The eigenvectors may be chosen to be orthonormal, so that $\mathbf{C}$ is an orthogonal matrix. \pause That is, 
$\mathbf{CC}^\prime = \mathbf{C}^\prime\mathbf{C} = \mathbf{I}$.
\end{frame}

\section{Positive definite matrices}

\begin{frame}
\frametitle{Positive definite matrices} \pause
%\framesubtitle{} 
The $n \times n$ matrix $\mathbf{A}$ is said to be \emph{positive definite} if
\begin{displaymath} 
   \mathbf{y}^\prime \mathbf{A} \mathbf{y} > 0
\end{displaymath}
for \emph{all} $n \times 1$ vectors $\mathbf{y} \neq \mathbf{0}$. \pause 
It is called \emph{non-negative definite} (or sometimes positive semi-definite) if \pause 
$\mathbf{y}^\prime \mathbf{A} \mathbf{y} \geq 0$.
\end{frame}

\begin{frame}
\frametitle{Example: Show $\mathbf{X}^\prime\mathbf{X}$ non-negative definite} \pause
Let $\mathbf{X}$ be an $n \times p$ matrix of real constants and let $\mathbf{y}$ be $p \times 1$. \pause Then $\mathbf{Z} = \mathbf{Xy}$ is $n \times 1$, and  \pause
\begin{eqnarray*}
     &   & \mathbf{y}^\prime \, (\mathbf{X}^\prime\mathbf{X}) \, \mathbf{y} \\ \pause
     & = & (\mathbf{Xy})^\prime (\mathbf{Xy}) \\ \pause
     & = & \mathbf{Z}^\prime \mathbf{Z} \\ \pause
     & = & \sum_{i=1}^n Z_i^2 \pause \geq 0 ~~ \blacksquare
\end{eqnarray*} 
\end{frame}


\begin{frame}
\frametitle{Some properties of symmetric positive definite matrices}
\framesubtitle{Variance-covariance matrices are often assumed positive definite.}  \pause
For a symmetric matrix, \pause
\vspace{5mm}
  \begin{itemize}
    \item[] Positive definite
    \item[]~~~~~~~~$\Downarrow$
    \item[] All eigenvalues positive  \pause
    \item[]~~~~~~~~$\Downarrow$
    \item[] Inverse exists \pause $\Leftrightarrow$ Columns (rows) linearly independent.
  \end{itemize} \pause

\vspace{5mm} 
If a real symmetric matrix is also non-negative definite, as a variance-covariance matrix \emph{must} be, Inverse exists $\Rightarrow$ Positive definite
\end{frame}

\begin{frame}
\frametitle{Showing Positive definite $\Rightarrow$ Eigenvalues positive} \pause
% \framesubtitle{For example}
Let the $p \times p$ matrix $\mathbf{A}$ be positive definite, so that $\mathbf{y}^\prime \mathbf{A} \mathbf{y} > 0$ for all $\mathbf{y} \neq \mathbf{0}$. \pause
  \begin{itemize}
    \item[] $\lambda$ an eigenvalue means $\mathbf{Ax} =  \lambda\mathbf{x}$. \pause
    \item[] $\Rightarrow \mathbf{x}^\prime \mathbf{A} \mathbf{x} = \mathbf{x}^\prime  \lambda\mathbf{x} \pause > 0$, \pause since the eigenvector $\mathbf{x} \neq \mathbf{0}$. \pause
    \item[] $\Rightarrow \lambda \mathbf{x}^\prime \mathbf{x} > 0$. \pause 
    \item[] $\Rightarrow  \frac{\lambda\mathbf{x}^\prime\mathbf{x}} {\mathbf{x}^\prime\mathbf{x}} > \frac{0}{\mathbf{x}^\prime\mathbf{x}} \pause = 0$. \pause
    \item[] $\Rightarrow \lambda > 0$ ~ $\blacksquare$
  \end{itemize}
\end{frame}

\begin{frame}
\frametitle{Inverse of a diagonal matrix}
\framesubtitle{To set things up} \pause
Suppose $\mathbf{D} = [d_{i,j}]$ is a diagonal matrix with non-zero diagonal elements. \pause It is easy to verify that
\begin{displaymath}
       \left( \begin{array}{c c c c }
        d_{1,1} & 0          & \cdots  & 0            \\
            0     & d_{2,2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & d_{n,n}    \\ 
\end{array} \right)
       \left( \begin{array}{c c c c }
        1/d_{1,1} & 0          & \cdots  & 0            \\
            0     & 1/d_{2,2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & 1/d_{n,n}    \\ 
\end{array} \right) = \mathbf{I} \pause
\end{displaymath} 
And
\begin{displaymath}
       \left( \begin{array}{c c c c }
        1/d_{1,1} & 0          & \cdots  & 0            \\
            0     & 1/d_{2,2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & 1/d_{n,n}    \\ 
\end{array} \right)
       \left( \begin{array}{c c c c }
        d_{1,1} & 0          & \cdots  & 0            \\
            0     & d_{2,2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & d_{n,n}    \\ 
\end{array} \right) = \mathbf{I}
\end{displaymath} 
\end{frame}

\begin{frame}
\frametitle{Showing Eigenvalues positive $\Rightarrow$ Inverse exists}
\framesubtitle{For a symmetric, positive definite matrix} \pause

Let $\mathbf{A}$ be symmetric and positive definite. \pause Then $\mathbf{A} = \mathbf{CDC}^\prime$, \pause and its eigenvalues are positive. \pause

\vspace{2mm}

Let $\mathbf{B} = \mathbf{C}\mathbf{D}^{-1}\mathbf{C}^\prime$. Show $\mathbf{B} = \mathbf{A}^{-1}$. \pause

% \vspace{2mm}

\begin{eqnarray*}
    \mathbf{AB} & = & \pause \mathbf{CDC}^\prime \, \mathbf{C}\mathbf{D}^{-1}\mathbf{C}^\prime \pause =  \mathbf{I} \pause \\
    \mathbf{BA} & = & \pause \mathbf{C}\mathbf{D}^{-1}\mathbf{C}^\prime \, \mathbf{CDC}^\prime \pause =  \mathbf{I} \pause
\end{eqnarray*}
So
{\LARGE
\begin{displaymath}
    \mathbf{A}^{-1} = \mathbf{C}\mathbf{D}^{-1}\mathbf{C}^\prime
\end{displaymath} 
} % End size
\end{frame}

\section{Square root matrices}

\begin{frame}
\frametitle{Square root matrices}
\framesubtitle{For symmetric, non-negative definite matrices} \pause % Real if non-neg def 
{\small
To set things up, define
\begin{displaymath}
    \mathbf{D}^{1/2} = 
        \left( \begin{array}{c c c c }
        \sqrt{\lambda_1} & 0          & \cdots  & 0            \\
            0     & \sqrt{\lambda_2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & \sqrt{\lambda_n}    \\ 
\end{array} \right)
\end{displaymath}  \pause
So that
\begin{eqnarray*}
    \mathbf{D}^{1/2} \mathbf{D}^{1/2}  \pause
    & = & 
        \left( \begin{array}{c c c c }
        \sqrt{\lambda_1} & 0          & \cdots  & 0            \\
            0     & \sqrt{\lambda_2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & \sqrt{\lambda_n}    \\ 
\end{array} \right)
        \left( \begin{array}{c c c c }
        \sqrt{\lambda_1} & 0          & \cdots  & 0            \\
            0     & \sqrt{\lambda_2}  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & \sqrt{\lambda_n}    \\ 
\end{array} \right) \\ \pause
    & = &
    \left( \begin{array}{c c c c }
        \lambda_1 & 0          & \cdots  & 0            \\
            0     & \lambda_2  & \cdots  & 0            \\
            \vdots   & \vdots  & \ddots  & \vdots       \\
            0        & 0       & \cdots  & \lambda_n    \\ 
\end{array} \right) \pause = \mathbf{D}
\end{eqnarray*}
} % End size
\end{frame}

\begin{frame}
\frametitle{For a non-negative definite, symmetric matrix $\mathbf{A}$} \pause
Define 
{\LARGE
\begin{displaymath}
    \mathbf{A}^{1/2} = \mathbf{C} \mathbf{D}^{1/2} \mathbf{C}^\prime
\end{displaymath} \pause
} % End size
So that 
\begin{eqnarray*}
    \mathbf{A}^{1/2}\mathbf{A}^{1/2}  \pause
    & = & 
    \mathbf{C} \mathbf{D}^{1/2} \mathbf{C}^\prime
    \mathbf{C} \mathbf{D}^{1/2} \mathbf{C}^\prime \\ \pause
    & = & 
    \mathbf{C} \mathbf{D}^{1/2} \, \mathbf{I} \,
    \mathbf{D}^{1/2} \mathbf{C}^\prime \\ \pause
    & = & 
    \mathbf{C} \mathbf{D}^{1/2}
    \mathbf{D}^{1/2} \mathbf{C}^\prime \\ \pause
    & = & 
    \mathbf{C} \mathbf{D} \mathbf{C}^\prime \\ \pause
    & = &  \mathbf{A} 
\end{eqnarray*} 
\end{frame}

\begin{frame}
\frametitle{The square root of the inverse is the inverse of the square root} \pause

Let $\mathbf{A}$ be symmetric and positive definite, with $\mathbf{A} = \mathbf{CDC}^\prime$. \pause

\vspace{2mm}

Let $\mathbf{B} = \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime $. \pause What is $\mathbf{D}^{-1/2}$? \pause

\vspace{2mm}

Show $\mathbf{B} = \left( \mathbf{A}^{-1} \right)^{1/2}$. \pause
\begin{eqnarray*}
    \mathbf{BB}  & = & \pause \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime 
                       \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime \\ \pause
                 & = &  \mathbf{C}\mathbf{D}^{-1} \mathbf{C}^\prime \pause  = \mathbf{A}^{-1} 
\end{eqnarray*} \pause

\vspace{2mm}

Show $\mathbf{B} = \left( \mathbf{A}^{1/2} \right)^{-1}$ \pause
\begin{eqnarray*}
    \mathbf{A}^{1/2}\mathbf{B} \pause  & = & \mathbf{C}\mathbf{D}^{1/2} \mathbf{C}^\prime 
                       \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime \pause = \mathbf{I} \\ \pause
    \mathbf{B} \mathbf{A}^{1/2} \pause & = &  \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime
                        \mathbf{C}\mathbf{D}^{1/2} \mathbf{C}^\prime \pause  = \mathbf{I}
\end{eqnarray*} \pause
Just write ~~~~~~~ \LARGE{$\mathbf{A}^{-1/2} = \mathbf{C}\mathbf{D}^{-1/2} \mathbf{C}^\prime$}
\end{frame}

\section{R}

\begin{frame}[fragile]
\frametitle{Matrix calculation with R} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> is.matrix(3)  # Is the number 3 a 1x1 matrix?
\end{verbatim} \pause
} % End color
\begin{verbatim}
[1] FALSE
\end{verbatim} \pause
{\color{blue}
\begin{verbatim}
> treecorr = cor(trees); treecorr
\end{verbatim} \pause
} % End color
\begin{verbatim}
           Girth    Height    Volume
Girth  1.0000000 0.5192801 0.9671194
Height 0.5192801 1.0000000 0.5982497
Volume 0.9671194 0.5982497 1.0000000
\end{verbatim} \pause
{\color{blue}
\begin{verbatim}
> is.matrix(treecorr)
\end{verbatim}  \pause
} % End color
\begin{verbatim}
[1] TRUE
\end{verbatim}
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Creating matrices}
\framesubtitle{Bind rows into a matrix}  \pause
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> # Bind rows of a matrix together
> A = rbind(  c(3, 2, 6,8),
+             c(2,10,-7,4),
+             c(6, 6, 9,1)   ); A
\end{verbatim} \pause
} % End color
\begin{verbatim}
     [,1] [,2] [,3] [,4]
[1,]    3    2    6    8
[2,]    2   10   -7    4
[3,]    6    6    9    1
\end{verbatim} \pause
{\color{blue}
\begin{verbatim}
> # Transpose
> t(A)
\end{verbatim} \pause
} % End color
\begin{verbatim}
     [,1] [,2] [,3]
[1,]    3    2    6
[2,]    2   10    6
[3,]    6   -7    9
[4,]    8    4    1
\end{verbatim} 
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Matrix multiplication}
\framesubtitle{Remember, $\mathbf{A}$ is $3\times 4$}  \pause
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> # U = AA' (3x3), V = A'A (4x4)
> U = A %*% t(A)
> V = t(A) %*% A; V
\end{verbatim} \pause
} % End color
\begin{verbatim}
     [,1] [,2] [,3] [,4]
[1,]   49   62   58   38
[2,]   62  140   -4   62
[3,]   58   -4  166   29
[4,]   38   62   29   81
\end{verbatim}
} % End size
\end{frame}


\begin{frame}[fragile]
\frametitle{Determinants} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> # U = AA' (3x3), V = A'A (4x4)
> # So rank(V) cannot exceed 3 and det(V)=0
> det(U); det(V)
\end{verbatim} \pause
} % End color
\begin{verbatim}
[1] 1490273
[1] -3.622862e-09
\end{verbatim} \pause
} % End size

Inverse of $\mathbf{U}$ exists, but inverse of $\mathbf{V}$ does not.
\end{frame}

\begin{frame}[fragile]
\frametitle{Inverses} \pause
%\framesubtitle{} 
  \begin{itemize}
    \item The \texttt{solve} function is for solving systems of linear equations like
$\mathbf{Mx}=\mathbf{b}$. \pause
    \item Just typing \texttt{solve(M)} gives $\mathbf{M}^{-1}$. \pause
  \end{itemize}

{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> # Recall U = AA' (3x3), V = A'A (4x4)
> solve(U)
\end{verbatim} \pause
} % End color
\begin{verbatim}
              [,1]          [,2]          [,3]
[1,]  0.0173505123 -8.508508e-04 -1.029342e-02
[2,] -0.0008508508  5.997559e-03  2.013054e-06
[3,] -0.0102934160  2.013054e-06  1.264265e-02
\end{verbatim} \pause
{\color{blue}
\begin{verbatim}
> solve(V)
\end{verbatim} \pause
} % End color
{\color{red}
\begin{verbatim}
Error in solve.default(V) : 
  system is computationally singular: reciprocal condition 
  number = 6.64193e-18
\end{verbatim}
} % End color
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Eigenvalues and eigenvectors} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> # Recall U = AA' (3x3), V = A'A (4x4)
> eigen(U)
\end{verbatim} \pause
} % End color
\begin{verbatim}
$values
[1] 234.01162 162.89294  39.09544

$vectors
           [,1]       [,2]        [,3]
[1,] -0.6025375  0.1592598  0.78203893
[2,] -0.2964610 -0.9544379 -0.03404605
[3,] -0.7409854  0.2523581 -0.62229894
\end{verbatim}
} % End size
\end{frame}


\begin{frame}[fragile]
\frametitle{$\mathbf{V}$ should have at least one zero eigenvalue}
\framesubtitle{Because $\mathbf{A}$ is $3 \times 4$, $\mathbf{V} = \mathbf{A}^\prime \mathbf{A}$, and the rank of a product is the minimum rank of the matrices.} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> eigen(V)
\end{verbatim} \pause
} % End color
\begin{verbatim}
$values
[1]  2.340116e+02  1.628929e+02  3.909544e+01 -1.012719e-14

$vectors
           [,1]         [,2]       [,3]         [,4]
[1,] -0.4475551  0.006507269 -0.2328249  0.863391352
[2,] -0.5632053 -0.604226296 -0.4014589 -0.395652773
[3,] -0.5366171  0.776297432 -0.1071763 -0.312917928
[4,] -0.4410627 -0.179528649  0.8792818  0.009829883
\end{verbatim}
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Spectral decomposition $\mathbf{V}= \mathbf{CDC}^\prime$} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> eigenV = eigen(V)
> C = eigenV$vectors; D = diag(eigenV$values); D
\end{verbatim} \pause
} % End color
\begin{verbatim}
         [,1]     [,2]     [,3]          [,4]
[1,] 234.0116   0.0000  0.00000  0.000000e+00
[2,]   0.0000 162.8929  0.00000  0.000000e+00
[3,]   0.0000   0.0000 39.09544  0.000000e+00
[4,]   0.0000   0.0000  0.00000 -1.012719e-14
\end{verbatim} \pause
{\color{blue}
\begin{verbatim}
> # C is an orthoganal matrix
> C %*% t(C)
\end{verbatim} \pause
} % End color
\begin{verbatim}
              [,1]         [,2]         [,3]          [,4]
[1,]  1.000000e+00 5.551115e-17 0.000000e+00 -3.989864e-17
[2,]  5.551115e-17 1.000000e+00 2.636780e-16  3.556183e-17
[3,]  0.000000e+00 2.636780e-16 1.000000e+00  2.558717e-16
[4,] -3.989864e-17 3.556183e-17 2.558717e-16  1.000000e+00
\end{verbatim}
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Verify $\mathbf{V}= \mathbf{CDC}^\prime$} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> V;   C %*% D %*% t(C)
\end{verbatim} \pause
} % End color
\begin{verbatim}
     [,1] [,2] [,3] [,4]
[1,]   49   62   58   38
[2,]   62  140   -4   62
[3,]   58   -4  166   29
[4,]   38   62   29   81

     [,1] [,2] [,3] [,4]
[1,]   49   62   58   38
[2,]   62  140   -4   62
[3,]   58   -4  166   29
[4,]   38   62   29   81
\end{verbatim}
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{Square root matrix $\mathbf{V}^{1/2}= \mathbf{CD}^{1/2}\mathbf{C}^\prime$} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> sqrtV = C %*% sqrt(D) %*% t(C)
\end{verbatim} \pause
} % End color
{\color{red}
\begin{verbatim}
Warning message:
In sqrt(D) : NaNs produced
\end{verbatim} \pause
} % End color
{\color{blue}
\begin{verbatim}
> # Multiply to get V
> sqrtV %*% sqrtV;  V
\end{verbatim} \pause
} % End color
\begin{verbatim}
     [,1] [,2] [,3] [,4]
[1,]  NaN  NaN  NaN  NaN
[2,]  NaN  NaN  NaN  NaN
[3,]  NaN  NaN  NaN  NaN
[4,]  NaN  NaN  NaN  NaN
     [,1] [,2] [,3] [,4]
[1,]   49   62   58   38
[2,]   62  140   -4   62
[3,]   58   -4  166   29
[4,]   38   62   29   81
\end{verbatim}
} % End size
\end{frame}

\begin{frame}[fragile]
\frametitle{What happened?} \pause
%\framesubtitle{} 
{\footnotesize % or scriptsize
{\color{blue}
\begin{verbatim}
> D; sqrt(D)
\end{verbatim} \pause
} % End color
\begin{verbatim}
         [,1]     [,2]     [,3]          [,4]
[1,] 234.0116   0.0000  0.00000  0.000000e+00
[2,]   0.0000 162.8929  0.00000  0.000000e+00
[3,]   0.0000   0.0000 39.09544  0.000000e+00
[4,]   0.0000   0.0000  0.00000 -1.012719e-14

         [,1]     [,2]     [,3] [,4]
[1,] 15.29744  0.00000 0.000000    0
[2,]  0.00000 12.76295 0.000000    0
[3,]  0.00000  0.00000 6.252635    0
[4,]  0.00000  0.00000 0.000000  NaN
\end{verbatim}
{\color{red}
\begin{verbatim}
Warning message:
In sqrt(D) : NaNs produced
\end{verbatim}
} % End color
} % End size
\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/302f15} {\small\texttt{http://www.utstat.toronto.edu/$^\sim$brunner/oldclass/302f15}}

\end{frame}


\end{document}

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