% This is a simple LaTex sample document that gives a submission format
%   for IEEE PAMI-TC conference submissions.  Use at your own risk.

% Make two column format for LaTex 2e.\
\documentclass[11pt,twocolumn]{article} %,twocolumn

%\usepackage{times,amsmath,amsfonts}

% Use following instead for LaTex 2.09 (may need some other mods as well).
%\documentstyle[times,twocolumn]{article}
\usepackage[dvips]{graphicx,graphics}
% Set dimensions of columns, gap between columns, and paragraph indent
\setlength{\textheight}{10in} \setlength{\textwidth}{7in}
%\setlength{\columnsep}{0.3125in} \setlength{\topmargin}{0in}
\setlength{\headheight}{0in} \setlength{\headsep}{-1in}
\setlength{\parindent}{1pc}
\setlength{\oddsidemargin}{-.5in}  % Centers text.
\setlength{\evensidemargin}{-.5in}

% Add the period after section numbers.  Adjust spacing.
\newcommand{\Section}[1]{\vspace{-8pt}\section{\hskip -1em.~~#1}\vspace{-3pt}}
\newcommand{\SubSection}[1]{\vspace{-3pt}\subsection{\hskip -1em.~~#1}
        \vspace{-3pt}}
\newcommand{\bqn}{\begin{eqnarray}}
\newcommand{\eqn}{\end{eqnarray}}
\newcommand {\diff}[1] {\frac{\partial}{\partial #1}}
\newcommand{\jacob}[3]{\frac{\partial^2 #3}{\partial #1 \partial #2}}
\newcommand{\der}[2]{\frac{\partial #2}{\partial #1}}
\begin{document}

% Make title bold and 14 pt font (Latex default is non-bold, 16pt)
\title{Stat 471: Lecture 5\\
Multivariate normal distributions}
% For single author (just remove % characters)
\author{Moo K. Chung\\
mchung@stat.wisc.edu}
% For two authors (default example)
\maketitle \thispagestyle{empty}

\begin{enumerate}
\item For given $z_1,z_2 \sim^{i.i.d.} N(0,1)$, we define a random
vector $z=(z_1,z_2)'$. Let $H=(h_{ij})$. Then we define a
multivariate normal distribution as $w=Hz \sim N(0, HH')$, where
$HH'$ is called the {\em covariance matrix}.

\item The covariance matrix $HH'$ is computed in the following
way. $${\bf Var} (Hz) = \mathbb{E} \big[(Hz)(Hz)'\big] =HH'$$
 assuming the mean is 0. The covariance matrix is {\em symmetric
 positive definite}.

\item Any symmetric positive definite matrix $V=(v_{ij})$ can be
factored as $V= HH',$ where $H=(h_{ij})$ is the lower triangular
matrix and $H'$ is the upper triangular matrix.

\begin{verbatim}
>>V=[2 1
   1 2]
>>Htrans=chol(V)
Htrans =
    1.4142    0.7071
         0    1.2247
>> Htrans'*Htrans
ans =
    2.0000    1.0000
    1.0000    2.0000
\end{verbatim}

It can be shown that $h_{11}=\sqrt{v_{11}}$ and
$h_{i1}=v_{i1}/h_{ii}$ for $i=2,\cdots, n$ (HW 2).

\item Generate 1000 multivariate normal random vectors with zero
mean and covariance $V=\left(
\begin{array}{cc}
  2 & 1 \\
  1 & 2 \\
\end{array}
\right).$ {\em Solution.} From HW 1, I assume that everyone will
have his/her own version of standard normal random number
generator called $\tt{snrnd}$. I assume that $\tt{snrnd}$ has one
argument $n$, the number of random numbers it is generating.
\begin{verbatim}
>>z1 = snrnd(1000); z2=snrnd(1000); 
>>w=Htrans'*[z1 z2]';
>> w(:,1:3)
ans =
    0.1772    0.4068   -1.6214
   -2.4402   -0.7073   -0.5949
\end{verbatim}
\item For given bivariate data $w_i$, test if the data follows
multivariate normal. You may use $\chi^2$ goodness-of-fit test but
a simpler approach would be to check if given bivariate data can
be generated via linear transform $w=Hz +\mu, z \sim N(0, I)$.
First estimate the sample covariance matrix:
$$\hat V = \frac{1}{n}\sum_{i=1}^n (w_i - \bar w)(w_i -\bar w)',$$
where $\bar w$ is the sample mean vector. The sample variance
matrix is an estimator of unknown population covariance $V$. It
can be computed using the {\em Kroneker tensor product} $\otimes$.
For $A=(a_{ij})$ and $B$, $A \otimes B = (a_{ij}B)$.
\begin{verbatim}
>>barw=mean(w,2); temp=w-kron(ones(1,1000),barw);
>>hatV=temp*temp'/1000;
hatV =
    2.0982    1.0232
    1.0232    1.9711
\end{verbatim}
Let $\hat V = \hat H \hat H'$ be the Cholesky factors. Then $z =
\hat H^{-1}w$.
\begin{verbatim}
>> hatH=chol(hatV)'
hatH =
    1.4485         0
    0.7064    1.2133
>> inv(hatH)
ans =
    0.6904         0
   -0.4019    0.8242
>> newz=inv(hatH)*w;
>>qqplot(newz(1,:));qqplot(newz(2,:));
\end{verbatim}
\end{enumerate}
\end{document}