\documentclass{article}[12pt]

% useful packages
\usepackage{titlesec}
\usepackage{fullpage}
\usepackage{amsmath,amssymb,amsthm,amsfonts}
\usepackage{graphicx}
\usepackage{enumerate}
\usepackage{algorithm,algorithmic}
\usepackage{xcolor}
\usepackage{bbm}
\usepackage{url}
\usepackage{hyperref}
\usepackage{caption,subcaption}

% theorem type environments
\newtheorem{thm}{Theorem}
\newtheorem{prop}{Proposition}
\newtheorem{lemma}{Lemma}
\newtheorem{cor}{Corollary}
\newtheorem{defn}{Definition}
\newtheorem{assump}{Assumption}
\newtheorem{example}{Example}
\newtheorem{conjecture}{Conjecture}

% frequently used symbols
\newcommand{\bE}{\mathbb{E}}
\newcommand{\bP}{\mathbb{P}}
\newcommand{\bQ}{\mathbb{Q}}
\newcommand{\bR}{\mathbb{R}}
\newcommand{\bS}{\mathbb{S}}
\newcommand{\bN}{\mathbb{N}}
\newcommand{\bZ}{\mathbb{Z}}
\newcommand{\sC}{{\mathcal C}} 
\newcommand{\sD}{{\mathcal D}} 
\newcommand{\sE}{{\mathcal E}} 
\newcommand{\sF}{{\mathcal F}} 
\newcommand{\sL}{{\mathcal L}} 
\newcommand{\sH}{{\mathcal H}} 
\newcommand{\sN}{{\mathcal N}} 
\newcommand{\sO}{{\mathcal O}} 
\newcommand{\sP}{{\mathcal P}} 
\newcommand{\sR}{{\mathcal R}} 
\newcommand{\sS}{{\mathcal S}}
\newcommand{\sU}{{\mathcal U}} 
\newcommand{\sX}{{\mathcal X}} 
\newcommand{\sY}{{\mathcal Y}} 
\newcommand{\sZ}{{\mathcal Z}}

% operators
\newcommand{\sign}{\mathop{\mathrm{sign}}}
\newcommand{\supp}{\mathop{\mathrm{supp}}} % support
\newcommand{\argmin}{\operatornamewithlimits{arg\ min}}
\newcommand{\argmax}{\operatornamewithlimits{arg\ max}}
\newcommand{\dist}{\operatorname{dist}}
\newcommand{\tr}{\text{tr}}
\newcommand{\vecop}{\text{vec}}
\newcommand{\st}{\operatorname{s.t.}}
\newcommand{\cut}{\setminus}
\newcommand{\ind}[1]{\mathbbm{1}\left\{#1\right\}} 
\newcommand{\given}{\ | \ }

% grouping operators
\newcommand{\brac}[1]{\left[#1\right]}
\newcommand{\set}[1]{\left\{#1\right\}}
\newcommand{\abs}[1]{\left\lvert #1 \right\rvert}
\newcommand{\paren}[1]{\left(#1\right)}
\newcommand{\norm}[1]{\left\|#1\right\|}
\newcommand{\ip}[2]{\left\langle #1,#2 \right\rangle}

% code commands
\newcommand{\matlab}{\textsc{Matlab }}
\newcommand{\python}{\textsc{Python }}

% header command
\newcommand{\homework}[3]{
    \pagestyle{myheadings}
    \thispagestyle{plain}
    \newpage
    \setcounter{page}{1}
    \setlength{\headsep}{10mm}
    \noindent
    \begin{center}
    \framebox{
        \vbox{\vspace{2mm}
            \hbox to 6.28in { {\bf EE 516: Mathematical Foundations of Machine Learning
            \hfill Winter \the\year} }
        \vspace{4mm}
        \hbox to 6.28in { {\Large \hfill Homework #1 \hfill} }
        \vspace{2mm}
        \hbox to 6.28in { \Large \hfill Due: #2, 11:59PM PT \hfill }
        \vspace{2mm}
        \hbox to 6.28in { {\it Student Name: #3} \hfill {\it Instructor Name: John Lipor}}
        \vspace{2mm}}
   }
   \end{center}
   \markboth{Homework #1}{Homework #1}
   \vspace*{4mm}
}

% For problem titles
\titleformat{\section}{\normalfont\bf}{Problem \thesection}{1em}{}

\begin{document}

% PUT YOUR NAME BELOW WHERE I HAVE \X
\homework{5}{February 17, \the\year}{\X}
% PUT YOUR NAME BELOW WHERE I HAVE \X

\section{\normalfont{(5 pts, 3 pts)}}

Let $A \in \bR^{m \times n}$ with $m \geq n$ and $\text{rank}(A) = n$, i.e., $A$ has full column rank.

\begin{enumerate}[(a)]
    \item Show using an SVD of $A$ that $A^{+} = (A^{T}A)^{-1}A^{T}$. When $A$ is known to have full column rank, one could use this direct formula to compute the pseudoinverse instead of employing an SVD.
    \item Deduce that the orthogonal projection onto the space spanned by the columns of $A$ (i.e., the \emph{range} of $A$) is given by $P = A(A^{T}A)^{-1}A^{T}$ using the SVD (i.e., do not just show that $P = AA^{+}$).
\end{enumerate}

\section{\normalfont{(6 pts, 2 pts, 5 pts)}}

Consider the plane $ax + by + cz = 0$.

\begin{enumerate}[(a)]
    \item Describe how you would use an SVD to find basis vectors for the plane.
    \item How many basis vectors are required to express a point on the plane?
    \item Find the point on the plane that is closest to $(\alpha,\beta,\gamma)$.
\end{enumerate}
\section{\normalfont{(3 pts)}}

In the usual least-squares problem, we found the minimum-norm $x$ that minimized $\norm{r}_{2}$, where $r = Ax - b$ is the residual. We saw that the optimal $x$ can be expressed in terms of $b$ and an SVD of $A$. Let $A$ be an $m \times n$ matrix so that $x$ is an $n \times 1$ vector and $b$ is an $m \times 1$ vector. Thus $r$ is an $m \times 1$ vector. In an upcoming application known as iteratively reweighted least-squares (IRLS), we will instead minimize $\norm{Wr}_{2}$, where $W$ is a
diagonal \emph{weight} matrix; this is referred to as the \emph{weighted least-squares} problem. 
Determine the optimal $x$ for the weighted least-squares problem
\begin{equation*}
    \hat{x} = \argmin_{x} \norm{W(Ax - b)}_{2}^{2}.
\end{equation*}
Your answer should be in terms of the pseudoinverse or SVD of a new matrix that you define---the dependence on $W$, $A$, and $b$ need not be directly apparent, and the solution should be \emph{very} simple.

\section{\normalfont{(5 pts each)}}

\begin{enumerate}[(a)]
    \item Determine the solution to the problem below, either in terms of the pseudoinverse of an appropriately defined matrix or of a set of equations resembling the normal equations.
        \begin{equation*}
            \hat{x} = \argmin_{x} \norm{Ax - b}_{2}^{2} + \norm{Cx - d}_{2}^{2} + \lambda \norm{x}_{2}^{2}.
        \end{equation*}
    \item Assume $A \in \bR^{m \times n}$ with $m > n$ has rank $n$. Find the solution of
        \begin{equation*}
            \hat{x} = \argmin_{x} \norm{Ax - b}_{2}^{2} + \norm{Ax + b}_{2}^{2}.
        \end{equation*}
    \item Does the answer to part (b) change if $A \in \bR^{m \times n}$ with $m > n$ but with rank less than $n$? Justify your answer using the SVD/pseudoinverse viewpoint of the least-squares solution.
\end{enumerate}

\section{\normalfont{(8 pts)}}

Let $X = U\Sigma V^{T}$ be a rank-$r$ matrix. Suppose we are given the estimates $\hat{U}$ and $\hat{V}$ (of $U$ and $V$, respectively) formed from the SVD of the noisy matrix $\hat{X} = X + N$, where $N$ is a noise matrix of the same size as $X$. Consider the diagonal matrix $\hat{D}$ that arises as the solution to the optimization problem
\begin{equation*}
    \hat{D} = \argmin_{D = \text{diag}(d_{1},\dots,d_{r})} \norm{X - \hat{U}(:, 1:r)D\hat{V}(:,1:r)^{T}}_{F},
\end{equation*}
where the notation $A(:,1:r)$ denotes the first $r$ columns of the matrix $A$. Express $\hat{D}$ in closed form in terms of $U, V, \Sigma, r, \hat{U}$ and $\hat{V}$. \textbf{Hint:} Follow a similar approach to that used on page 6.5 of the notes.

\section{\normalfont{(4 pts each)}}

In this problem, you will explore the differences between low-rank approximation (PCA) and least-squares approximation.

\begin{enumerate}[(a)]
    \item First consider one-dimensional data lying with the affine relationship $y = ax + b$. Using the script \texttt{prob6a}, fit the given data in \texttt{xx,yy} using a least-squares estimate for $a, b$ and plot the resulting line. Next, fit the same data using a rank-1 approximation, as described on Pg. 6.12-6.13 of the notes.
        \\ \textbf{Turn in} the resulting plot.
        \\ \textbf{Hint:} To obtain a correct low-rank approximation, you must first \emph{center} the data by removing its mean. This can be done efficiently using the \texttt{mean} and \texttt{repmat} functions. Don't forget to add the mean back when you go to plot.
    \item Next, we'll fit a two-dimensional plane using each method. Complete the least-squares and low-rank approximations in \texttt{prob6b}. Although it is difficult to see, the resulting planes are different. The takeaway here is that we aren't limited to fitting one-dimensional data!
        \\ \textbf{Turn in} the resulting plot.
        \\ \textbf{Note:} There is no offset for this part, so no need to center.
\end{enumerate}

\section{\normalfont{(5 pts)}}

In many applications, we want to fit a system of equations and know that the variables being fit are non-negative. For example, if we are trying to fit $F = ma$, we know that $m$ and $a$ are positive. In this case, we wish to solve the \emph{non-negative least-squares} problem
\begin{eqnarray*}
    x_{\rm NNLS} =& \argmin_{x} \norm{Ax - b}_{2}^{2} \\
    \text{subject to} & x \geq 0,
\end{eqnarray*}
where we have explicitly imposed the constraint that the solution must be non-negative. If we solve the least-squares problem without the non-negativity constraint and the solution happens to be non-negative, then that solution is also the optimal non-negative solution. However, typically the unconstrained solution will have negative values. Setting the negative values to zero after the fact usually does not solve the stated problem.

We can solve optimization problems of the form
\begin{eqnarray*}
    x_{\rm opt} =& \argmin_{x} \norm{Ax - b}_{2}^{2} \\
    \text{subject to} & x \in \sC,
\end{eqnarray*}
where $\sC$ denotes a convex set, via a projected gradient method that combines the usual gradient descent update with an additional projection onto the convex set. For the non-negative least-squares problem, the update iteration is given by
\begin{equation*}
    x_{k+1} = P_{\sC}\left( x_{k} - \mu A^{H}(Ax_{k} - b) \right),
\end{equation*}
and it will converge to the correct solution for $0 < \mu < 2/\sigma_{1}(A)^{2}$ with
\begin{equation*}
    P_{\sC}(z) = \max(0,z).
\end{equation*}

\textbf{Your task} is to implement the non-negative least-squares gradient descent algorithm by completing the \texttt{nnls} file and testing it on the included \texttt{prob7} script. Test your algorithm by plotting $\norm{x_{k} - \hat{x}}_{2}$ vs $k$ as in the previous assignment \textbf{for both} your \texttt{lgsd} function from HW4 and your \texttt{nnls} function from this homework. Note: You can obtain $\hat{x}$ using the \matlab command \texttt{lsqnonneg} or your own result after 500
iterations in \python (I found the \texttt{scipy.optimize.nnls} did not work very well).

\textbf{Turn in} your \texttt{nnls} code and the resulting plot.

\section{\normalfont{(5 pts)}}

\vspace{1em}
\noindent Spend at most one hour reading the paper, ``K-SVD: An Algorithm for Designing Overcomplete Dictionaries for Sparse Representation'' by Aharon, Elad, and Bruckstein again. Answer the below questions:
\begin{itemize}
    \item Review the five bullets from Homework 1, Problem 10. Has your ability to understand any of these improved? If so, which items do you now understand better?
    \item Which topics of the course (e.g., matrix factorizations, least squares) have contributed most to your answer above?
    \item Which aspects of the course (lectures, exercises, homework, demos, exam) have contributed most to your answer above?
\end{itemize}

\section*{Course Project}

The course project description has been listed on the course website. Please turn in the names of your group members and which topic your group has selected.

\end{document}
