\section{Sparse Linear Regression}

We first provide an example problem that shows some advantages of casting a problem in the form of convex optimization.

Suppose we have a response $y \in \R^n$
(\ie a scalar response for each of $n$ training examples)
and predictors $A = (A_1, \dots, A_p) \in \R^{n \times p}$
(\ie $p$-dimensional features for each of $n$ training examples).
We wish to use a linear model $y \approx A x$, where $x \in \R^p$.

If $n > p$, this is the classical linear regression problem
\begin{equation}
    \min_{x \in \R^n} \norm{y - A x}^2.
    \label{eq:reg-std}
\end{equation}
The solution to this problem is well-defined and easy to find if $n > p$.
But in some applications $n \ll p$, and we further believe that many of the extracted features $A_1, \dots, A_p$ could be irrelevant.
We therefore wish to find a model $x$ with many zero coefficients,
as illustrated in Figure~\ref{fig:regression-model}.

\begin{figure}[ht]
    \begin{center}
        \includegraphics[width=.5\textwidth]{regression}
    \end{center}
    \caption{The linear regression model $y \approx A x$, with $n \ll p$ and a sparse $x$.}
    \label{fig:regression-model}
\end{figure}

For example, if $y$ is the size of a tumor, it might be reasonable to suppose that it can be expressed as a linear combination of genetic information in $A$, but we expect that most genes will be unimportant to this combination: therefore most components of $x$ will be zero.

\subsection{Best Subset Selection}

One natural idea for incorporating sparsity into \eqref{eq:reg-std} is to constrain the number of nonzero components of $x$ to be no more than some $k$:
\begin{equation}
    \min_{x \in \R^n} \norm{y - A x}^2
    \text{ subject to } \norm{x}_0 \le k,
    \label{eq:best-subset}
\end{equation}
where $\norm{x}_0$, the 0-norm, denotes the number of nonzero components in $x$.

This ``norm'' is non-convex, however, which makes solving \eqref{eq:best-subset} difficult (in fact, NP-hard). Figure~\ref{fig:l0-sublevel-set} shows a sublevel set in $\R^2$, illustrating the nonconvexity.
The only known algorithm is essentially brute-force evaluation of all possible subsets; the solution cannot be computed in practice for $p \gtrsim 40$.
We also know very little about statistical properties of the solution.

\begin{figure}[ht]
    \begin{center}
        \begin{tikzpicture}
            %\draw[brown] (current bounding box.south west) rectangle (current bounding box.north east);
            % box
            \draw (-1.2,-1.2) -- coordinate (x axis mid) (1.2,-1.2);
            \draw (-1.2,-1.2) -- coordinate (y axis mid) (-1.2,1.2);
            \draw (-1.2,1.2) -- (1.2,1.2) -- (1.2,-1.2);

            \foreach \x in {-1,-.5,0,.5,1}
                \draw [-,yshift=-1.2cm](\x,0pt) -- (\x cm,-3pt)
                    node[anchor=north] {\scriptsize $\x$};

            \foreach \y in {-1,-.5,0,.5,1}
                \draw [-,xshift=-1.2cm](0pt, \y cm) -- (-3pt, \y cm)
                    node[anchor=east] {\scriptsize $\y$};

            \node[below=.5cm] at (x axis mid) {\small $x_1$};
            \node[rotate=90, above=.5cm] at (y axis mid) {\small $x_2$};

            % set
            \draw[<->,very thick] (0,-1.2) -- (0,1.2);
            \draw[<->,very thick] (-1.2,0) -- (1.2,0);
        \end{tikzpicture}
        % TODO: too much spacing here...
    \end{center}
    \caption{Sublevel set for the constraint of \eqref{eq:best-subset}: $\{  x \in \R^2 : \norm{x}_0 \le 1 \}$.}
    \label{fig:l0-sublevel-set}
\end{figure}

We call this the 0-norm because it is a limit case of the general $p$-norm. A vector $x = (x_1, x_2, \cdots, x_n)$ has $p$-norm $\norm{x}_p = \left( \sum_{i=1}^n \abs{x_i}^p \right) ^{1/p}$. Figure~\ref{fig:2D-norm} shows contours of different $p$-norms in $\R^2$. Note that the non-convex contour in Figure~\ref{fig:2D-norm} is just a sketch for $p$-norm with $p < 1$, and that in fact the $p$-norm is not a norm for $p < 1$.

\begin{figure}[ht]
    \centering
    \includegraphics[scale=0.4]{2dnorm}
    \caption{%
        An illustration for different $p$-norms in 2D.%
    }
    \label{fig:2D-norm}
\end{figure}

\subsection{Forward Stepwise Regression}

Another natural approach is given by the following algorithm:

\begin{enumerate}
    \item Find the variable $j$ such that $\abs{A_j^T y}$ is largest. If the variables have been centered to have zero mean and scaled to have unit variance, then $A_j^T y$ is the correlation between predictors $A_j$ and response $y$.

    \item Update $x_j$ by regressing $y$ onto $A_j$, so that we are solving
        $\displaystyle \min_{x_j \in \R} \norm{y - A_j x_j}^2$.

    \item Find a new variable $k \ne j$ such that $\abs{A_k^T r}$ is largest, where $r$ is the residual $y - A_j x_j$. This finds the variable most correlated to the residuals.

    \item Update $x_j, x_k$ by regressing $y$ onto $A_j, A_k$.

    \item Repeat until the desired number of variables have been selected.
\end{enumerate}

It is unknown whether this algorithm minimizes some criterion function.
A few statistical properties of this estimate are known, such as its consistency and some details about the convergence rate, but the proofs are relatively complicated and many properties remain unknown.

\subsection{Lasso}

Instead of the 0-norm of best subset selection, we can also use the (convex) 1-norm $\norm{x}_1 = \sum_{i=1}^p \abs{x_i}$:
\begin{equation}
    \min_{x \in \R^n} \norm{y - A x}^2
    \text{ subject to } \norm{x}_1 \le t.
    \label{eq:lasso}
\end{equation}
A sublevel set in $\R^2$ is shown in Figure~\ref{fig:l1-sublevel-set}.

\begin{figure}[ht]
    \begin{center}
        \begin{tikzpicture}
            % box
            \draw (-1.2,-1.2) -- coordinate (x axis mid) (1.2,-1.2);
            \draw (-1.2,-1.2) -- coordinate (y axis mid) (-1.2,1.2);
            \draw (-1.2,1.2) -- (1.2,1.2) -- (1.2,-1.2);

            \foreach \x in {-1,-.5,0,.5,1}
                \draw [-,yshift=-1.2cm](\x,0pt) -- (\x cm,-3pt)
                    node[anchor=north] {\scriptsize $\x$};

            \foreach \y in {-1,-.5,0,.5,1}
                \draw [-,xshift=-1.2cm](0pt, \y cm) -- (-3pt, \y cm)
                    node[anchor=east] {\scriptsize $\y$};

            \node[below=.5cm] at (x axis mid) {\small $x_1$};
            \node[rotate=90, above=.5cm] at (y axis mid) {\small $x_2$};

            % set
            \draw[very thick, fill=gray, fill opacity=.8]
                (1,0) -- (0,1) -- (-1,0) -- (0,-1) -- cycle;
        \end{tikzpicture}
        % TODO: too much spacing here...
    \end{center}
    \caption{Sublevel set for the constraint of \eqref{eq:lasso}: $\{ x \in \R^2 : \norm{x}_1 \le 1 \}$.}
    \label{fig:l1-sublevel-set}
\end{figure}

It is less obvious that this is the case than for the two previous methods, but the solutions obtained include exact zeros, with lower values of $t$ resulting in more zeros.
A proof of this fact relies on the Karush-Kuhn-Tucker conditions and will be given later in the course, but Figure~\ref{fig:lasso-intuition} gives some intuition as to why this is the case in $\R^2$.

\begin{figure}[ht]
    \begin{center}
        \begin{tikzpicture}
            % axes
            \draw[<->] (-2,0) -- (5,0);
            \draw[<->] (0,-2) -- (0,4);

            % lasso feasible set
            \draw[very thick, fill=gray, fill opacity=.8]
                (1,0) -- (0,1) -- (-1,0) -- (0,-1) -- cycle;

            % true solution
            \draw[fill=black] (3,2) circle (.1);

            \foreach \a in {3,6.6,10.2,13.8,17.4}
                \draw (3,2) ellipse (\a * .2 and \a * .11);
        \end{tikzpicture}
    \end{center}
    \caption{An illustration of why the lasso method gives solutions with zeros in $\R^2$.}
    \label{fig:lasso-intuition}
\end{figure}

The gray diamond represents the set of solutions allowable by the constraint on $\norm{x}_1$, the black circle is the least-squares solution, and the contours are level sets of the objective $\norm{y - A x}^2$.
The first contour to intersect the feasible set is the optimal solution.
In this case, that corresponds to a sparse solution. In higher dimensions, the sparse edges of the feasible set are even ``spikier'' and so it is very likely that the solution is sparse.

This is known as the Least Absolute Shrinkage and Selection Operator (LASSO) method, and is effectively a convex relaxation of best subset selection.
Thanks to the convexity, this problem can be solved readily, and many properties are known about the solutions; the proof of consistency, for example, follows readily from properties of convexity and the nature of the errors.
