\documentclass{article}

\usepackage{fullpage}
\usepackage{color}
\usepackage{amsmath}
\usepackage{url}
\usepackage{verbatim}
\usepackage{graphicx}
\usepackage{parskip}
\usepackage{amssymb}
\usepackage{nicefrac}
\usepackage{listings} % For displaying code
\usepackage{algorithm2e} % pseudo-code

% Answers
\def\ans#1{\par\gre{Answer: #1}}

% Colors
\definecolor{blu}{rgb}{0,0,1}
\def\blu#1{{\color{blu}#1}}
\definecolor{gre}{rgb}{0,.5,0}
\def\gre#1{{\color{gre}#1}}
\definecolor{red}{rgb}{1,0,0}
\def\red#1{{\color{red}#1}}
\def\norm#1{\|#1\|}

% Math
\def\R{\mathbb{R}}
\def\argmax{\mathop{\rm arg\,max}}
\newcommand{\argmin}[1]{\mathop{\hbox{argmin}}_{#1}}
\newcommand{\mat}[1]{\begin{bmatrix}#1\end{bmatrix}}
\newcommand{\alignStar}[1]{\begin{align*}#1\end{align*}}
\def\half{\frac 1 2}

% LaTeX
\newcommand{\fig}[2]{\includegraphics[width=#1\textwidth]{a1f/#2}}
\newcommand{\centerfig}[2]{\begin{center}\includegraphics[width=#1\textwidth]{a1f/#2}\end{center}}
\def\items#1{\begin{itemize}#1\end{itemize}}
\def\enum#1{\begin{enumerate}#1\end{enumerate}}


\begin{document}

\title{CPSC 540 Assignment 2 (due January 31 at midnight)}
\author{}
\date{}
\maketitle
\vspace{-4em}

The assignment instructions are the same as for the previous assignment, but for this assignment you can work in groups of 1-3. However, please only hand in one assignment for the group.


\blu{\enum{
\item Name(s):
\item Student ID(s):
}}




\section{Calculation Questions}


\subsection{Convexity}

\blu{Show that the following functions are convex, by only using one of the definitions of convexity (i.e., without using the ``operations that preserve convexity" or using convexity results stated in class)}:\footnote{\red{That $C^0$ convex functions are below their chords, that $C^1$ convex functions are above their tangents, or that $C^2$ convex functions have a positive semidefinite Hessian.}}
\enum{
\item Upper bound on second-order Taylor expansion: $f(v) = f(w) + \nabla f(w)^T(v-w) + \frac{L}{2}\norm{v-w}^2$. \\
\item Probit regression: $f(w) = \sum_{i=1}^n -\log p(y_i | w, x_i)$ (where the probability is defined as in the last assignment).
\item Maximum function: $f(w) = \max_{j \in \{1,2,\dots,d\}}w_d$.
}

Hint: Max are is not differentiable in general, so you cannot use the Hessian for the last one. Hint: Max are is not differentiable in general, so you cannot use the Hessian for the last one. \red{For the second one, you can assume that $PDF(z) \geq - z\cdot \text{CDF}(z)$ (where PDF$(z)$ is the PDF of a standard normal and CDF$(z)$ is the CDF).}


\blu{Show that the following functions are convex (you can use results from class and operations that preserve convexity if they help)}:% (with $\lambda \geq 0$):
\enum{
\setcounter{enumi}{3}
\item 1-class SVMs: $f(w,w_0) = \sum_{i=1}^N[\max\{0,w_0 - w^Tx_i\} - w_0] + \frac{\lambda}{2}\norm{w}_2^2$, where $\lambda \geq 0$ and $w_0$ is a variable.
\item Mixed-norm regularization: $f(w) = \norm{w}_{p,q} = \left(\sum_{g\in\mathcal{G}}\norm{w_g}^p_q\right)^{\frac{1}{p}}$ (\red{for $p \geq 1$ and $q \geq 1$}).
\item Minimum entropy over pairs of variables: $f(w) = \max_{\{i,j \; | \; i \neq j\}}\{w_i\log w_i + w_j\log w_j \}$ subject to $w_i >0$ for all $i$.
}

\subsection{Convergence of Gradient Descent}

For these questions it will be helpful to use the ``convexity inequalities'' notes posted on the webpage.

\enum{
\item In class we showed that if $\nabla f$ is $L$-Lipschitz continuous and $f$ is strongly-convex, then with a step-size of $\alpha_k = 1/L$ gradient descent has a convergence rate of 
\[
f(w^k) - f(w^*) = O(\rho^k).
\]
\blu{Show that under these assumptions that a convergence rate of $O(\rho^k)$ in terms of the function values implies that the iterations have a convergence rate of
\[
\norm{w^k - w^*} = O(\rho^{k/2}).
\]}
\item A basic variation on the Armijo line-search is to set the step-size $\alpha_k$ to be $\alpha_k = (0.5)^p$, where $p$ is the smallest constant such that
\[
f(w^k - \alpha_k\nabla f(w^k)) \leq f(w^k) - \frac{\alpha_k}{2}\norm{\red{\nabla} f(w^k)}^2.
\]
\blu{Show that if $\nabla f$ is Lipschitz-continuous and $f$ is bounded below, that setting $\alpha_k$ in this way means gradient descent requires at most $t=O(1/\epsilon)$ iterations to find at least one iteration $k$ with $\norm{\red{\nabla f(w^k)}}^2 \leq \epsilon$}. Hint: you'll first want to figure out a progress bound of the form
\[
f(w^{k+1}) \leq f(w^k) - \gamma\norm{\nabla f(w^k)}^2,
\]
for some constant $\gamma$ that holds for all iterations. It may also help to recognize that if a function is $L$-Lipschitz continuous that it is also $L'$-Lipschitz continuous for any $L' \geq L$.
\item \blu{Show that if $\nabla f$ is Lipschitz-continuous and $f$ is convex (but not necessarily strongly-convex or PL), that if we run gradient descent for $k$ iterations with $\alpha_k = 1/L$ we have
\[
f(w^k) - f(w^*) = O(1/k).
\]
}
Hint: you will have to use one the \red{$C^1$} definition of convexity in our usual progress bound coming from Lipschitz-continuity of the gradient. You can try to turn this into a telescoping sum in terms of $\norm{w^k-w^*}^2$ (\red{which invovles ``completing the square''}). Finally, note that our usual progress bound also gives that $f(w^{k+1}) \leq f(w^k)$ for all $k$.
}

\subsection{Beyond Gradient Descent}

\enum{
\item We can write the proximal-gradient update as
\begin{align*}
w^{k+\half} & = w^k - \alpha_k \nabla f(w^k)\\
w^{k+1} &= \argmin{v\in\R^d}\left\{\frac{1}{2}\norm{v -w^{k+\half}}^2 + \alpha_kr(v)\right\}.
\end{align*}
\blu{Show that this is equivalent to setting
\[
w^{k+1} \in \argmin{v\in\R^d} \left\{ f(w^k) + \nabla f(w^k)^T(v-w^k) + \frac{1}{2\alpha_k}\norm{v-w^k}^2  +r(v)\right\}.
\]
}
\item In class we showed that if $\nabla f$ is coordinate-wise $L$-Lipschitz and satisfies PL then randomized coordinate optimization with a step-size of $1/L$ satisfies
\[
\mathbb{E}[f(w^{k+1}) - f^*] \leq \left(1 - \frac{\mu}{Ld}\right)[f(w^k) - f^*].
\]
Consider a $C^2$ function satisfying PL where coordinate $i$ has its own Lipschitz constant $L_j$,
\[
\nabla_{jj}^2 f(w) \leq L_j.
\]
Consider using a step-size of $1/L_{j_k}$ and sampling $j_k$ proportional to the $L_j$, $p(j_k = j) = \frac{L_j}{\sum_{j'}L_{j'}}$. \blu{Show that this gives a faster convergence rate}. Hint: the previous result corresponds to using $L = \max_j\{L_j\}$.
\item Consider an SVM problem where the $x^i$ are \emph{sparse}. In the bonus slides, we show how to efficiently apply stochastic subgradient methods in this setting by using the representation $w^k = \beta^k v^k$. Now consider a case where we know an L2 ball that contains the optimal solution. In other words, we know a $\tau$ such that $\norm{w^*} \leq \tau$. If $\tau$ is small enough, we can use a \emph{projected} stochastic subgradient method (projecting onto the L2-ball):
\begin{align*}
w^{k+\frac{1}{2}} & =  w^k - \alpha_k g^k - \alpha_k\lambda w^k\\
w^{k+1} & = \begin{cases}
w^{k + \frac{1}{2}} & \text{if $\norm{w^{k+\frac{1}{2}}} \leq \tau$}\\
\frac{\tau}{\norm{w^{k+\frac{1}{2}}}}w^{k+\frac{1}{2}} & \text{if $\norm{w^{k+\frac{1}{2}}} > \tau$}
\end{cases}
\end{align*}
By constraining the $w^k$ to a smaller set, this can sometimes dramatically improve the performance in the early iterations.\footnote{We can obtain a value for $\tau$ for any L2-regularized problem where the loss is bounded below, by using that $f(w^*) + \frac{\lambda}{2}\norm{w^*}^2 \leq f(0)$ and replacing $f(w^*)$ with its lower bound (which can be 0 for SVMs).} However, the projection operator is a full-vector operation so this would substantially slow down the runtime of the method. \blu{Derive a recursion that implements this algorithm without using full-vector operations}.
Hint: you may need to track more information than $\beta^k$ and $v^k$.
}

\section{Computation Questions}


\subsection{Proximal-Gradient}


If you run the demo \emph{example\_group.jl}, it will load a dataset and fit a multi-class logistic regression (softmax) classifier. This dataset is actually \emph{linearly-separable}, so there exists a set of weights $W$ that can perfectly classify the training data (though it may be difficult to find a $W$ that perfectly classifiers the validation data). However, 90\% of the columns of $X$ are irrelevant. Because of this issue, when you run the demo you find that the training error is $0$ while the test error is something like $0.2980$.

\enum{
\item Write a new function, \emph{logRegSoftmaxL2}, that fits a multi-class logistic regression model with L2-regularization (this only involves modifying the objective function). \blu{Hand in the modified loss function and report \red{the validation error achieved with $\lambda = 10$} (which is the  best value among powers to 10). Also report the number of non-zero parameters in the model and the number of original features that the model uses}.
\item While L2-regularization reduces overfitting a bit, it still uses all the variables even though 90\% of them are irrelevant. In situations like this, L1-regularization may be more suitable. Write a new function, \emph{logRegSoftmaxL1}, that fits a multi-class logistic regression model with L1-regularization. You can use the function \emph{findMinL1}, which minimizes the sum of a differentiable function and an L1-regularization term.  \blu{Report the number of non-zero parameters in the model and the number of original features that the model uses}.
\item L1-regularization achieves sparsity in the \emph{model parameters}, but in this dataset it's actually the \emph{original features} that are irrelevant. We can encourage sparsity in the original features by using \emph{group} L1-regularization. Write a new function, \emph{proxGradGroupL1}, to allow (disjoint) \emph{group} L1-regularization. Use this within a new function, \emph{softmaxClassiferGL1}, to fit a group L1-regularized multi-class logistic regression model (where \emph{rows} of $W$ are grouped together and we use the L2-norm of the groups).  \blu{Hand in both modified  functions (\emph{logRegSoftmaxGL1} and \emph{proxGradGroupL1}) and report the validation error achieved with $\lambda=10$. Also report the number of non-zero parameters in the model and the number of original features that the model uses}.
}


\subsection{Coordinate Optimization}

The function \emph{example\_CD.jl} loads a dataset and tries to fit an L2-regularized logistic regression model using coordinate optimization. Unfortunately, if we use $L_f$ as the Lipschitz constant of $\nabla f$, the runtime of this procedure is $O(d^3 + nd^2\frac{L_f}{\mu}\log(1/\epsilon))$. This comes from spending $O(d^3)$ computing $L_f$, having an iteration cost of $O(nd)$, and requiring a $O(d\frac{L_f}{\mu}\log(1/\epsilon))$ iterations. This non-ideal runtime is also reflected in practice: the algorithm's iterations are relatively slow and even after 500 ``passes'' through the data it isn't particularly close to the optimal function value.

\enum{
\item Modify this code so that the runtime of the algorithm is $O(nd\frac{L_c}{\mu}\log(1/\epsilon))$, where $L_c$ is the Lipschitz constant of \emph{all} partial derivatives $\nabla_i f$. You can do this by modifying the iterations so they have a cost $O(n)$ instead of $O(nd)$, and instead of using a step-size of $1/L_f$ they use a step-size of $1/L_c$ (which is given by $\frac{1}{4}\max_j\{\norm{x_j}^2\} + \lambda$, where $x_j$ is column $j$ of the matrix $X$).
 \blu{Hand in your code and report the final function value and total time}.
 \item To further improve the performance, make a new version of the code which samples the variable to update $j_t$ proportional to the individual Lipschitz constants $L_j$ of the coordinates, and use a step-size of $1/L_{j_t}$. You can use the function \emph{sampleDiscrete} (in \emph{misc.jl}) to sample from a discrete distribution given the probability mass function.  \blu{Hand in your code, and report the final function value as well as the number of passes}.
\item \blu{Report the number of passes the algorithm takes as well as the final function value} if you use \emph{uniform sampling} but use a step-size of $1/L_{j_t}$. 
\item Suppose that when we use a step-size of $1/L_{j_t}$, we see that uniform sampling outperforms Lipschitz sampling. \blu{Why could this be consistent with the bounds we've shown?}
}


\subsection{Stochastic Gradient}

If you run the demo \emph{example\_SG.jl}, it will load a dataset and try to fit an L2-regularized logistic regression model using 10 ``passes'' of stochastic gradient using the step-size of $\alpha_t = 1/\lambda t$ that is suggested in many theory papers. Note that in other high-level languages (like R/Matlab/Python) this demo would run really slowly so you would need to write the inner loop in a low-level language like C, but in Julia you can directly write the stochastic gradient code and have it run fast.

Unfortunately, despite Julia making this code run fast compared to other high-level languages, the performance of this stochastic gradient method is atrocious. It often goes to areas of the parameter space where the objective function overflows and the final value is usually in the range of something like $6.5-7.5 \times 10^4$. This is quite far from the solution of $2.7068 \times 10^4$ and is even worse than just choosing $w=0$ which gives $3.5 \times 10^4$. (This is unlike gradient descent and coordinate optimization, which never increase the objective function if your step-size is small enough.)

\enum{
\item Although $\alpha_t = 1/\lambda t$ gives the best possible convergence rate in the worst case, in practice it's typically horrible (as we're not usually opitmizing the hardest possible $\lambda$-strongly convex function). Experiment with different choices of step-size sequence to see if you can get better performance. \blu{Report the step-size sequence that you found gave the best performance, and the objective function value obtained by this strategy for one run}.
\item Besides tuning the step-size, another strategy that often improves the performance is using a (possibly-weighted) average of the iterations $w^t$. Explore whether this strategy can improve performance.  \blu{Report the performance with an averaging strategy, and the objective function value obtained by this strategy for one run}. (Note that the best step-size sequence with averaging might be different than without averaging.)
\item A popular variation on stochastic is AdaGrad, which uses the iteration
\[
w^{k+1} = w^k - \alpha_k D_k \nabla f(w^k),
\]
where the element in position $(j,j)$ of the diagonal matrix $D_k$ is given by $1/\sqrt{\delta + \sum_{k'=0}^k(\nabla_j f_{i_{k'}}(w^{k'}))^2}$. Here, $i_k$ is the example $i$ selected on iteration $k$ and $\nabla_j$ denotes element $j$ of the gradient (and in AdaGrad we typically don't average the steps).  Implement this algorithm and experiment with the tuning parameters $\alpha_t$ and $\delta$. \blu{Hand in your code as well as the best step-size sequence you found and again report the performance for one run}.
\item Impelement the SAG algorithm with a step-size of $1/L$, where $L$ is the maximum Lipschitz constant across the training examples ($L = \frac{1}{4}\max_i\{\norm{x^i}^2\} + \lambda$).  \blu{Hand in your code and again report the performance for one run}.
}

\section{Very-Short Answer Questions}

Give a short and concise 1-sentence answer to the below questions.
\enum{
\item What is a situation where we can violate the golden rule on our test set, and still have it be a reasonable approximation of test error?
\item Do convex functions necessarily have a minima?
\item What is the purpose of the descent lemma?
\item Why did we show that the gradient norm converges to 0 when we run gradient descent? (As opposed to showing that gradient descent finds the global optimum.)
\item What is the relationship between the projected-gradient method and the proximal-gradient method?
\item Why do we prefer the proximal-gradient method over the subgradient method?
\item Why did we say in 340 that regularizing by the L2-norm doesn't give sparsity, and now we're saying we can use the L2-norm to give sparse solutions?
\item What is an example of a pattern that can be enforced with structured sparsity?
\item Under what condition do we prefer to use coordinate optimization over gradient descent?
\item For a machine learning objective that is the sum over $n$ training examples, why should we never use the subgradient method?
\item What is the advantage of SAG over stochastic gradient methods without a memory?
}



 
\end{document}