\documentclass[11pt]{article}
\usepackage{latexsym}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsthm}
\usepackage{epsfig}
\usepackage[tight]{subfigure}
\usepackage{algorithmic}
\usepackage{algorithm}

\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, F11)}{Lecture \#14 (Monday October 31th)}{Lecturer: Drew Bagnell}{Scribe: Carl Doersch \footnotemark}{Support Vector Machines and Bayes Regression}
\footnotetext{Some content adapted from previous scribes: Alan Kraut, Robert Fisher,  Heather Knight, Ben Xinjilefu, Kevin Lipkin, Alvaro Collet-Romea, Gandalf, Laura Lindzey and Hans Pirnay}

\section{Linear SVMs}
We begin by considering binary, linear classification. In this problem, we are trying to map elements from our feature space into the set $\{-1,1\}$. When we have two features, $F_1,F_2\in \mathbb R$, we can think of the problem graphically. In the following figure, we represent the points with a label of 1 as O's, we the points with a label of $-1$ as X's. We wish to find a linear decision boundary (a straight line) that separates the points from each class. The decision boundary is shown as the dotted line.
\begin{figure}[htbp]
    \centering
        \includegraphics[width=9cm]{binaryClassification.pdf}
    \caption{Linear binary classification}
    \label{fig:bin}  % you can refer to this figure ase robots
\end{figure}

To represent this notion mathematically, we denote the $i^{th}$ training point as $f_i$, and the class of this point as $y_i$. The linear classification problem requires us to find a set of weights $w$ such that:
    \begin{equation}
   \forall i\ y_iw^Tf_i\geq 0
    \end{equation}

Note that the decision boundary need not pass through the origin. We often use a dummy variable $F_0=c$ for constant $c\neq 0$, which allows us to use $w_0$ as an offset.

One issue with this formulation of the problem is that it has a trivial solution. Specifically, if we set all of our weights to 0, we will have $y_iw^Tf_i=0$ for all points. To address this issue, we modify our constraints to be:
    \begin{equation}
   \forall i\ y_iw^Tf_i\geq margin
    \end{equation}
This can have the added benefit of improving generalization. Once again, we can represent this graphically. In the following picture, the space between the dotted line and the decision boundary represents the margin.
\begin{figure}[htbp]
    \centering
        \includegraphics[width=9cm]{binaryClassification2.pdf}
    \caption{Linear binary classification with margins}
    \label{fig:bin}  % you can refer to this figure ase robots
\end{figure}

We further see that the magnitude of our weight vector does not matter, only the orientation affects classification. Therefore we can restrict $w$ to $||w||^2=1$. Taken together, this yields the following formulation:

\[ \begin{array}{rl}
    \text{{\bf Maximize }}        & margin \\
    \mbox{\text{{\bf Such that }}} & \forall i\ y_iw^Tf_i\geq margin \\
                &  ||w||^2\leq 1
    \end{array}
\]
We use $||w||^2\leq 1$ instead of $||w||^2= 1$ to make sure that our problem remain convex.  Note that we do not include a bias term; we can incorporate a bias by simply adding a 1 at the end of each input vector.  This means we are regularizing the bias in this formulation of the SVM, though the original SVM papers did not.

This problem is generally solved in an equivalent form that is easier to optimize, where we hold the margin fixed and minimize the magnitude of the weights:

\[ \begin{array}{rl}
    \text{{\bf Minimize }}        & ||w||^2\ \\
    \mbox{\text{{\bf Such that }}} & \forall i\ y_iw^Tf_i\geq 1 \\
    &  margin>0
    \end{array}
\]
  
  This last formulation is an example of a quadratic program, which we are able to solve efficiently. Unfortunately, the hard margin SVM requires that the data be linearly separable, which is almost never the case. To address this issue, we introduce a slack variable, $\xi$. The problem now becomes:
  \[ \begin{array}{rl}
    \text{{\bf Minimize }}        & \lambda ||w||^2 + \displaystyle\sum_i \xi_i \\
    \mbox{\text{{\bf Such that }}} & y_i w^T f_i \geq 1-\xi_i\\
    &   \xi_i \geq 0 \\
     &  margin>0
    \end{array}
\]
There will be one slack variable $\xi_i$ to correspond to each of the $T$ data-points that we are considering.  Many SVM solvers further transform this optimization into a dual form, where the weight vector $w$ is written as a linear combination of the support vectors, and this linear combination is optimized rather than optimizing $w$ directly.  This allows us to formulate nonlinear versions of the SVM, but this is beyond the scope of the class.

\subsection{Gradient Descent}

To perform gradient descent, we first observe that $\xi = \max (0,1-y_i w^T f_i)$, because the slack variable is 0 if this point is labeled correctly. This allows us to generate the loss function
    \begin{equation}
    l_t = \lambda ||w||^2 + \max (0,1-y_t w^T f_t)
    \end{equation}
    
Our update for $w$ is now as follows.  If the output was correct by at least the margin, the only contribution to the loss is the magnitude of $w$:
    \begin{equation}
    w \leftarrow w - 2 \alpha_t \lambda w
    \end{equation}
And if the output for this time step was incorrect or not outside the margin ($\xi>0$),
    \begin{equation}
    w \leftarrow w - 2 \alpha_t \lambda w + \alpha_t y_t f_t
    \end{equation}
    
We note that this loss function is {\bf not} minimizing the number of mistakes made by the algorithm. In fact, solving this problem with a 0-1 loss function (which would minimize the number of mistakes) is known to be NP-hard. We can visualize the difference between these loss functions with the following figure, in which $Z=w^Tf_i$:

\begin{figure}[htbp]
    \centering
        \includegraphics[width=9cm]{lossFunctions.pdf}
    \caption{Loss functions}
    \label{fig:loss}  % you can refer to this figure ase robots
\end{figure}

Our loss function, the hinge loss function, is convex, while the mistake loss (0-1 loss) function is not. However, we do see intuitively that the loss function we are using is the best ``convexification'' of the NP-hard problem.  In practice the convexified loss is more appropriate anyway, since we generally want to minimize ``how wrong'' we are when we do mislabel a point.

\subsection{Selecting $\alpha_t$}
\begin{itemize}
\item First idea is to set $\alpha_t$ proportional to $\frac{1}{\sqrt{t}}$.
\item If we have $T$ elements, each with a maximum value of $F$, the maximum gradient, $G$, is $\sqrt{TF}$. This results in a regret that is $R\leq \sqrt{FGT}$.
\item This is not as good as we could do.
\item Notice that $l_t$ is an extremely good convex function. It is a quadratic plus a convex function. In the same way all convex functions lie above a line (a subgradient) from every point, $l_t$ lies above a quadratic from every point.
\item Specifically, if it is always the case that
\begin{equation}
f(y) \geq f(x) + \frac{H}{2}(y-x)^2 + \nabla f_x^T(y-x)
\end{equation}
then $f(x)$ is said to be H-strongly convex.
\item In this case $l_t$ is $\lambda$-strongly convex.
\item If $\alpha_t=\frac{G}{Ht}$, then regret$\leq\frac{G^2}{H}(1+\log{t})$. $\log{t}$ is \emph{really} good, and this learning rate and algorithm is essentially the current best for this class of problem.

\end{itemize}

\section{Bayes' Linear Regression}
In linear regression, the goal is to predict a continuous outcome variable.  In particular, let:

\begin{itemize}
\item $\theta$ = parameter vector of the learned model
\item $x_{t} \in \mathbb{R}^{n}$ = set of features at every timestep, used for prediction
\item $y_{t} \in \mathbb{R}$ = true outcome
\end{itemize}

Then our model is as follows:

$$y_{t}=\theta x_{t} +\epsilon_{t}$$

where $\epsilon_{t} \sim \mathcal{N}(0,\sigma^{2})$ is noise that is independent of everything else.  Thus the likelihood if $\theta$ is known may be written:

$$p(y_{t}|x_{t},\theta)=\frac{1}{z} \exp(\theta x (2\sigma^{2})^{-1})$$

\begin{figure}
\centering
\def\svgwidth{.4\columnwidth}
\input{gm.pdf_tex}
\caption{Graphic model of BLR}
\label{F:BLR}
\end{figure}

In BLR, we maintain a distribution over the weight vector $\theta$ to represent our beliefs about what $\theta$ is likely to be.  The math is easiest if we restrict this distribution to be a Gaussian: $\theta \sim N(\mu_0, \Sigma_0)$.  This has the following form:
\begin{align}
\label{E:moment}
p(\theta) = \frac{1}{z} \exp \{-\frac{1}{2}(\theta-\mu_0)^{T}\Sigma_0^{-1}(\theta-\mu_0)\}
\end{align}

This is called the moment parameterization of a Guassian.  This is the most intuitive parameterization, but it is actually more computationally convenient to use the natural parameterization.  In the natural parameterization, the multiplication operation required in Bayes rule becomes a simple sum:

\begin{equation}
p(\theta) = \frac{1}{z} \exp \{J_0^{T}\theta - \frac{1}{2}\theta^{T}P_0 \theta\}
\end{equation}

where
\begin{align}
P_0 & = \Sigma_0^{-1} \\
J_0 & = P_0 \mu_0
\end{align}

Note that for these to be probability distributions, $\Sigma_0$ must be positive-definite; otherwise the distribution will have an infinite integral and cannot be normalized.  It is also generally constrained to be symmetric, as correlations between random variables are symmetric.

Recall that if we have a distribution over $\theta$, and we recieve a new $\tilde{x}_{t+1}$, we can compute the distribution over $\tilde{y}_{t+1}$ by integrating out $\theta$
\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}

We want to derive how to update our estimate of $\theta$ given a series of data $D = \{x, y\}$ up to timestep $t$ and $x_{t+1}$.  We can calculate $p(\theta | D)$ recursively using the Bayes' Theorem:

\begin{align}
    p(\theta | D) = \dfrac{1}{z} p(y_t |\theta,D) p(\theta|D)
\end{align}
where we know the likelihood $p(y_t |D, \theta) = p(y_t |x_t, \theta)$ is given by a Gaussian $\mathcal{N}(\theta^Tx_t,\sigma_t^2)$
\begin{eqnarray}
p(y_t|x_t,\theta) =  \dfrac{1}{z} \exp \{ \dfrac{-(\theta^T x_t-y_t)^2}{2\sigma^2} \}
\end{eqnarray}
The updating rule for $p(\theta|D)$ at timestep $t$ is given by 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_t-y_t\right)^2 + J_{t-1}^T\theta - \frac{1}{2}\theta^TP_{t-1}\theta
\end{eqnarray}
collecting terms to find updates $J_t$ and $P_t$:
\begin{eqnarray}
&=&-\frac{1}{2\sigma^2}\left(\theta^Tx_t x_t^T\theta-2\theta^Tx_t y_t + y_t^2\right) + J_{t-1}^T\theta - \frac{1}{2}\theta^TP_{t-1}\theta\\
&=&\left(\frac{x_t^Ty_t}{\sigma^2}+J_{t-1}^T\right)\theta  -\frac{1}{2}\theta^T\left(\frac{x_t x_t^T}{\sigma^2} + P_{t-1}\right)\theta - \frac{y_t^2}{2\sigma^2}
\end{eqnarray}
Since this all happens in the exponent of an exponential function, the constant $y_t^2/\sigma^2$-term can be shifted into the regularizing constant $z$. Thus, the update rules for $J_t$ and $P_t$ are
\begin{eqnarray}
J_t &=& \frac{y_t x_t}{\sigma^2} + J_{t-1}\\
P_t &=& \frac{x_t x_t^T}{\sigma^2} + P_{t-1}
\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_t$. The precision is only affected by \underline{input not output}. This is a consequence of having the same $\sigma$ (observation error) everywhere in space.
\end{enumerate}

\end{document}
