\documentclass[11pt]{article}
\usepackage{latexsym}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsthm}
\usepackage[pd,col,3d,turt]{mgtex_23}
\usepackage[tight]{subfigure}

\newcommand{\handout}[5]{
  \noindent
  \begin{center}
  \framebox{
    \vbox{
      \hbox to 5.78in { {#1} \hfill #2 }
      \vspace{4mm}
      \hbox to 5.78in { {\Large \hfill #5  \hfill} }
      \vspace{2mm}
      \hbox to 5.78in { {\em #3 \hfill #4} }
    }
  }
  \end{center}
  \vspace*{4mm}
}

\newcommand{\lecture}[5]{\handout{#1}{#2}{#3}{#4}{#5}}

\newtheorem{theorem}{Theorem}
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{observation}[theorem]{Observation}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{definition}[theorem]{Definition}
\newtheorem{claim}[theorem]{Claim}
\newtheorem{fact}[theorem]{Fact}
\newtheorem{assumption}[theorem]{Assumption}

% 1-inch margins, from fullpage.sty by H.Partl, Version 2, Dec. 15, 1988.
\topmargin 0pt
\advance \topmargin by -\headheight
\advance \topmargin by -\headsep
\textheight 8.9in
\oddsidemargin 0pt
\evensidemargin \oddsidemargin
\marginparwidth 0.5in
\textwidth 6.5in

\parindent 0in
\parskip 1.5ex
%\renewcommand{\baselinestretch}{1.25}

\begin{document}

\lecture{Statistical Techniques in Robotics (16-831, F08)}{Lecture \#14 (October 9)}{Lecturer: Drew Bagnell}{Scribe: Laura Lindzey (2009) Hans Pirnay(2008)}{Bayesian Linear Regression Pt. 2, Gaussian Properties}

\section{Parameterizations for Gaussians}
There are two common parameterizations for Gaussians, the moment parameterization and the natural parameterization. It is often most practical to switch back and forth between representations, depending on which calculations are needed. 

\subsection {Moment Parameterization} 
\begin{eqnarray}
\mathcal{N}(\mu,\Sigma) = p(\theta) &=& \frac{1}{z}\exp\left( -\frac{1}{2}\left(\theta -\mu\right)\Sigma^{-1}\left(\theta -\mu\right)\right)
\end{eqnarray}

Given:
\begin{eqnarray}
\mathcal{N}\left(\left[\begin{array}{c}\mu_1\\\mu_2\end{array}\right],\left[\begin{array}{cc}\Sigma_{11}&\Sigma_{12}\\\Sigma_{21}&\Sigma_{22}\end{array}\right]\right)
\end{eqnarray}
\textbf{Marginal:} computing $p(x_2)$
\begin{eqnarray}
\mu_2^{\mathrm{marg}} &=& \mu_2\\
\Sigma_2^{\mathrm{marg}} &=& \Sigma_{22}
\end{eqnarray}
\textbf{Conditional:} computing $p(x_1|x_2)$
\begin{eqnarray}
\mu_{1|2} &=& \mu_1 + \Sigma_{12}\Sigma_{22}^{-1}\left(x_2 -\mu_2\right)\\
\Sigma_{1|2} &=& \Sigma_{11} - \Sigma_{12}\Sigma_{22}^{-1} \Sigma_{21}
\end{eqnarray}


\subsection{Natural Parameterization}
\begin{eqnarray}
\tilde{\mathcal{N}}(J,P) = \tilde{p}(\theta) &=& \frac{1}{z}\exp\left( J^T\theta-\frac{1}{2}\theta^{T}P\theta\right)
\end{eqnarray}


Given:
\begin{eqnarray}
\mathcal{N}\left(\left[\begin{array}{c}J_1\\J_2\end{array}\right],\left[\begin{array}{cc}P_{11}&P_{12}\\P_{21}&P_{22}\end{array}\right]\right)
\end{eqnarray}
\textbf{Marginal:} computing $p(x_2)$
\begin{eqnarray}
J_2^{\mathrm{marg}} &=& J_2 - P_{21}P_{11}^{-1}J_1\\
P_1^{\mathrm{marg}} &=& P_{12} - P_{21}P_{11}^{-1}P_{12}
\end{eqnarray}
\textbf{Conditional:} computing $p(x_1|x_2)$
\begin{eqnarray}
J_{1|2} &=& J_1 - P_{12}x_2\\
P_{1|2} &=& P_{11}
\end{eqnarray}





\section{Graphical Models}
The matrix $P$ of the natural parameterization has a graphical model interpretation. If and only if there is a non-zero entry for $(z_i,z_j)$, then there is a lik between $z_i$ and $z_j$ in the MRF that corresponds to $\tilde{N}(J,P)$. 
\begin{center}
\InitGraph{8}{3}{0}{0}{1cm}
\TextAt(0.7,2.1)[br]{P = }
\Rectangle(2,0.75,3.8,2.5)
\TextAt(2,2.5)[br]{X}
\TextAt(2.5,2.5)[br]{X}
\TextAt(3,2.5)[br]{0}
\TextAt(2,2)[br]{X}
\TextAt(2.5,2)[br]{X}
\TextAt(3,2)[br]{X}
\TextAt(2,1.5)[br]{0}
\TextAt(2.5,1.5)[br]{X}
\TextAt(3,1.5)[br]{X}
\LineAt(5,2,6,0.8)
\LineTo(7,2)
\TextAt(5,2)[tl]{$z_1$}
\Circle(0.1)
\TextAt(6,0.8)[br]{$z_2$}
\Circle(0.1)
\TextAt(7,2)[tr]{$z_3$}
\Circle(0.1)
\CloseGraph
\end{center}
Following the graphical model interpretation, $P$ is in many cases highly structured. Consider for example the graphical model of a markov chain:
\begin{center}
\InitGraph{5.4}{2}{0}{0}{1cm}
\MultiPut(0.2,1)(1,0){5}{
  \Circle(0.2)
}
\MultiPut(0.4,1)(1,0){4}{
  \RelativeLine(0.6,0)
}
\TextAt(-0.2,1)[tr]{$x_1$}
\TextAt(0.8,1)[tr]{$x_2$}
\TextAt(1.8,1)[tr]{$x_3$}
\TextAt(2.8,1)[tr]{$x_4$}
\TextAt(3.8,1)[tr]{$x_5$}
\CloseGraph
\end{center}
The corresponding matrix $P$ will be non-zero only on the diagonal and off-diagonal:
\begin{eqnarray}
P &=& \left(\begin{array}{ccccc}
X&X&0&0&0\\
X&X&X&0&0\\
0&X&X&X&0\\
0&0&X&X&X\\
0&0&0&X&X
\end{array}\right)
\end{eqnarray}

\textbf{Note:} $P$ describes which variables \textbf{directly} affect each other. 

\textbf{Note:} $P^{-1}$ is, in general, not sparse! (this makes intuitve sense since $P^{-1} = \Sigma$ the covariance matrix, and the covariance of two states along the markov chain are not independent.) 



\section{Update Rules for Bayes Linear Regression}
We want to derive how to update our representation given a series of new data $D = \{x_1, y_1, ..., x_n, y_n\}$. Let the prior be $\theta$, where $\theta \sim N(\mu_0, \Sigma_0)$. Then, we want to calculate 
\begin{eqnarray}
p(y_t | D) & = \int p(y_t | \theta, x_t) p(\theta | D_t)
 & = E_{p(\theta | D)}[p(y_t | \theta, x]
\end{eqnarray}

So, we need to calculate $p(\theta | D)$: 

\begin{eqnarray}
p(y|x,\theta) &=& \mathcal{N}(\theta^Tx_t,\sigma_t^2) = \frac{1}{z}\exp\left(\frac{-(\theta^T x-y)(\theta^T x-y)}{2\sigma^2}\right)
\end{eqnarray}
(Don't worry about the weird notation of $\mathcal{N}$ as a function of $\sigma^2$. This is an arbitrary definition)
\subsection{Deriving the update rules}
Apply Bayes' Rule to the probability of a weight vector $\theta$ given a datapoint $D$.
\begin{eqnarray}
p(\theta|D) &=& \frac{p(D|\theta)p(\theta)}{z}
\end{eqnarray}
This results in the multiplication of two exponential functions. Adding the exponent of the prior to that of the likelihood yields
\begin{eqnarray}
-\frac{1}{2\sigma^2}\left(\theta^Tx-y\right)^2 + J^T\theta - \frac{1}{2}\theta^TP\theta
\end{eqnarray}
collecting terms to find updates $J_\theta'$ and $P_\theta'$:
\begin{eqnarray}
&=&-\frac{1}{2\sigma^2}\left(\theta^Txx^T\theta-2\theta^Txy + y^2\right) + J^T\theta - \frac{1}{2}\theta^TP\theta\\
&=&\left(\frac{x^Ty}{\sigma^2}+J^T\right)\theta  -\frac{1}{2}\theta^T\left(\frac{xx^T}{\sigma^2} + P\right)\theta - \frac{y^2}{2\sigma^2}
\end{eqnarray}
Since this all happens in the exponent of an exponential function, the constat $y^2$-term can be shifted into the regularizing $z$. Thus, the update rules for $J_\theta'$ and $P_\theta'$ are
\begin{eqnarray}
J_\theta' &=& \frac{xy}{\sigma^2}+J\\
P_\theta' &=& \frac{xx^T}{\sigma^2} + P
\end{eqnarray}
\begin{enumerate}
\item in a gaussian model, a new datapoint \underline{always} lowers the variance - this downgrading of the variance does not always make sense 
\item if you believe there are outliers, this model won't work for you
\item the variance is not a function of $y$. The precision if only affected by \underline{input not output}. This is a consequence of having the same $\sigma$ (observation error) everywhere in space. 
\end{enumerate}

\subsection{Transfer to moment parameterization}
The update rules for the natural parameterization at timestep $t$ are
\begin{eqnarray}
J  &+=& \frac{y_tx_t}{\sigma^2}\\
P &+=& \frac{1}{\sigma^2}xx^T.
\end{eqnarray}

Given the transfer rules to the moment parameterization
\begin{eqnarray}
\Sigma &=& P^{-1}\\
\mu &=& P^{-1}J
\end{eqnarray}

the moment parameterization after $N$ timesteps is then
\begin{eqnarray}
\Sigma_{\theta | D} &=& \left[ \Sigma_0^{-1} + \sum_{i=1}^N \frac{x_ix_i^T}{\sigma^2}  \right]^{-1}\\
\mu_\theta &=& \left[ \Sigma_0^{-1} + \sum_{i=1}^N \frac{x_ix_i^T}{\sigma^2}  \right] \sum_{t=1}^N \frac{y_tx_t}{\sigma^2}
\end{eqnarray}

\begin{enumerate}
\item If $\mu_0 \neq 0$, the update for $\mu_\theta$ will be more complicated
\item $\sum_{t=1}^N \frac{y_tx_t}{\sigma^2}$ is the gradient of Bayes online linear regression
\item this looks just like Newton's method
\item Computation time: $o(d^2)$ for update, $o(d^3)$ for mean (that can be reduced to $o(d^2)$ with tricks)
\end{enumerate}

\subsection{Making predictions}
Given all data $D$ up to timestep $t$ and $x_{t+1}$, the probability of an observation $\tilde{y}_{t+1}$ is
\begin{eqnarray}
p(\tilde{y}_{t+1}|x_{t+1},D) &=& \int p(\tilde{y}_{t+1}|x_{t+1},D,\theta)\cdot p(\theta|D)d\theta\\
&=& \int p(\tilde{y}_{t+1}|x_{t+1},\theta)\cdot p(\theta,D)d\theta
\end{eqnarray}
To know $p(\tilde{y}_{t+1}|x_{t+1},D)$, we only need $y$ and $\sigma^2$, because these parameters determine the gaussian.
\begin{center}
\InitGraph{4.5}{3}{0}{0.5}{1cm}
\MultiPut(1,2)(1,0){4}{\Circle(0.2)}
\MultiPut(1,1.8)(1,0){4}{\RelativeLine(0,-0.5)}
\MultiPut(1,1.1)(1,0){4}{\Circle(0.2)}
\MultiPut(1,0.9)(1,0){4}{\LineTo(2.5,0.4)}
\CircleAt(2.5,0.2)(0.2)
\TextAt(0,1.7)[tr]{$x_i$}
\TextAt(0,0.8)[tr]{$y_i$}
\TextAt(0,-.1)[tr]{$\theta$}
\CloseGraph
\end{center}



\end{document}
