Introduction In this homework we will explore some algorithms for gradient based optimization.  These algorithms have been crucial to the development of machine learning in the last few decades.  The most famous example is the backpropagation algorithm used in deep learning, which is in fact just an application of a simple algorithm known as (stochastic) gradient descent.  We will first implement gradient descent from scratch on a deterministic problem (no data), and then extend our implementation to solve a real world regression problem.

Points Allocation There are a total of 28 marks.

• Question 1 a): 2 marks

• Question 1 b): 1 mark

• Question 1 c): 1 mark

• Question 1 d): 2 marks

• Question 1 e): 2 marks

• Question 1 f): 4 marks

• Question 1 g): 3 marks

• Question 1 h): 1 mark

• Question 1 i): 3 marks

• Question 1 j): 4 marks

• Question 2 a): 2 marks

• Question 2 b): 1 mark

• Question 2 c): 2 marks

What to Submit

• A single PDF file which contains solutions to each question. For each question, provide your solution in the form of text and requested plots. For some questions you will be requested to provide screen shots of code used to generate your answer — onlyinclude these when they are explicitly asked for.

•  .py file(s) containing all code you used for the project, which should be provided in a separate .zip file. This code must match the code provided in the report.

• You maybe deducted points for not following these instructions.

• You may be deducted points for poorly presented/formatted work.  Please be neat and make your solutions clear. Start each question on a new page if necessary.

• You cannot submit a Jupyter notebook; this will receive a mark of zero. This does not stop you from developing your code in a notebook and then copying it into a .py file though, or using a tool such as

nbconvert or similar.

• We will setup a Moodle forum for questions about this homework. Please read the existing questions before posting new questions. Please do some basic research online before posting questions. Please only post clarification questions.  Any questions deemed to be fishing for answers will be ignored and/or deleted.

• Please check Moodle announcements for updates to this spec.  It is your responsibility to check for announcements about the spec.

• Please complete your homework on your own, do not discuss your solution with other people in the course.  General discussion of the problems is fine, but you must write out your own solution and acknowledge if you discussed any of the problems in your submission (including their name(s) and zID).

• As usual, we monitor all online forums such as Chegg, StackExchange, etc. Posting homework ques- tions on these site is equivalent to plagiarism and will result in a case of academic misconduct.

• You may not use SymPy or any other symbolic programming toolkits to answer the derivation ques- tions.  This will result in an automatic grade of zero for the relevant question.  You must do the

derivations manually.

When and Where to Submit

• Due date: Week 4, Monday June 19th, 2023 by 5pm.  Please note that the forum will not be actively monitored on weekends.

• Late submissions will incur a penalty of 5% per day from the maximum achievable grade.  For ex- ample, if you achieve a grade of 80/100 but you submitted 3 days late, then your final grade will be 80 _ 3 x 5 = 65. Submissions that are more than 5 days late will receive a mark of zero.

• Submission must be done throughMoodle, no exceptions.

Question 1. Gradient Based Optimization

The general framework for a gradient method for finding a minimizer of a function f  : Rn   →  R is defined by

x(k+1)  = x(k) _ αkVf (xk ),       k = 0, 1, 2, . . . ,                                           (1)

where αk   > 0 is known as the step size,  or learning rate.  Consider the following simple example of minimizing the function g(x) = 2′x3 + 1. We first note that g/ (x) = 3x2 (x3  + 1)- 1/2. We then need to choose a starting value of x, say x(0)  = 1. Let’s also take the step size to be constant, αk  = α = 0.1. Then we have the following iterations:

x(1)  = x(0) _ 0.1 x 3(x(0) )2 ((x(0) )3 + 1)- 1/2  = 0.7878679656440357 x(2)  = x(1) _ 0.1 x 3(x(1) )2 ((x(1) )3 + 1)- 1/2  = 0.6352617090300827 x(3)  = 0.5272505146487477

.

.

.

and this continues until we terminate the algorithm (as a quick exercise for your own benefit, code this up and compare it to the true minimum of the function which is x*   =  _1).  This idea works for functions that have vector valued inputs, which is often the case in machine learning.  For example, when we minimize a loss function we do so with respect to a weight vector, β. When we take the step- size to be constant at each iteration, this algorithm is known as gradient descent. For the entirety of this question, do not use any existing implementations of gradient methods, doing so will result in an automatic mark of zero for the entire question.

(a) Consider the following optimisation problem:

xe(m)Rn(in) f (x),

where

f (x) =

2 |Ax _ b|2 + 2 |x|2 ,

and where A e Rm×n, b e Rm  are defined as

A =  ┌ _11

' 0

1┐

1 ' ,

┌ 3 ┐

'_2' ,

and γ is a positive constant. Run gradient descent on fusing a step size of α = 0.1 and γ = 0.2 and starting point of x(0)   =  (1, 1, 1, 1). You will need to terminate the algorithm when the following condition is met:  |Vf (x(k))|2   <  0.001.  In your answer, clearly write down the version of the gradient steps (1) for this problem.  Also, printout the first 5 and last 5 values of x(k), clearly indicating the value of k, in the form:

k = 0,        x(k)  = [1, 1, 1, 1]

k = 1,        x(k)  = . . .

k = 2,        x(k)  = . . .

.

.

.

What to submit: Your working, and a printout of the value of the ridge solution based on (X train, Y train) . Include a screen shot of any code used for this section and a copy of your python code in solutions.py.

We will now solve the ridge problem but using numerical techniques.  As noted in the lectures, there are a few variants of gradient descent that we will briefly outline here. Recall that in gradient descent our update rule is

β (k+1)  = β (k) _ αkVL(β(k)),       k = 0, 1, 2, . . . ,

where L(β) is the loss function that we are trying to minimize. In machine learning, it is often the case that the loss function takes the form

L(β) =

n

↓ Li (β),

i=1

i.e. the loss is an average of n functions that we havelablled Li. It then follows that the gradient is also an average of the form

VL(β) =   VLi (β).

We can now define some popular variants of gradient descent .

(i) Gradient Descent (GD) (also referred to as batch gradient descent): here we use the full gradi- ent, as in we take the average over all n terms, so our update rule is:

β (k+1)  = β (k) _   VLi (β(k)),        k = 0, 1, 2, . . . .

(ii) Stochastic Gradient Descent (SGD): instead of considering all n terms, at the k-th step we choose an index ik randomly from {1, . . . , n}, and update

β (k+1)  = β (k) _ αkVLik (β(k)),       k = 0, 1, 2, . . . .

Here, we are approximating the full gradient VL(β) using VLik (β).

(iii) Mini-Batch Gradient Descent: GD (using all terms) and SGD (using a single term) represents the two possible extremes. In mini-batch GD we choose batches of size 1 < B  < n randomly at each step, call their indices {ik1 , ik2 , . . . , ikB }, and then we update

β (k+1)  = β (k) _   VLij (β(k)),        k = 0, 1, 2, . . . ,

j=1

so we are still approximating the full gradient but using more than a single element as is done in SGD.

(e) The ridge regression loss is

L(β) =

 |y _ Xβ|2(2) + φ|β|2(2).

Show that we can write

L(β) =   Li (β),

and identify the functions L1 (β), . . . , Ln (β). Further, compute the gradients VL1 (β), . . . , VLn (β)

What to submit: your working.

(f) In this question, you will implement (batch) GD from scratch to solve the ridge regression problem. Use an initial estimate β (0)  = 1p  (the vector of ones), and φ = 0.5 and run the algorithm for 1000 epochs (an epoch is one pass over the entire data, so a single GD step). Repeat this for the following step sizes:

α e {0.000001, 0.000005, 0.00001, 0.00005, 0.0001, 0.0005, 0.001, 0.005, 0.01}

To monitor the performance of the algorithm, we will plot the value

∆ (k)  = L(β(k)) _ L(β(ˆ)),

where β(ˆ) is the true (closed form) ridge solution derived earlier.  Present your results in a 3 x 3

grid plot, with each subplot showing the progression of ∆(k)  when running GD with a specific step-size. State which step-size you think is best and let β(K) denote the estimator achieved when running GD with that choice of step size. Report the following:

(i) The train MSE:  |ytrain _ Xtrain β (K) |2(2)

(ii) The test MSE:  |ytest _ Xtest β (K) |2(2)

What to submit: a single plot, the train and test MSE requested. Include a screen shot of any code used for this section and a copy of your python code in solutions.py.

(g) We will now implement SGD from scratch to solve the ridge regression problem.  Use an initial estimate β(0)  = 1p  (the vector of ones) and φ = 0.5 and run the algorithm for 5 epochs (this means a total of 5nupdates of β, where n is the size of the training set). Repeat this for the following step sizes:

α e {0.000001, 0.000005, 0.00001, 0.00005, 0.0001, 0.0005, 0.001, 0.006, 0.02}

Present an analogous 3 x 3 grid plot as in the previous question.  Instead of choosing an index randomly at each step of SGD, we will cycle through the observations in the order they are stored in X  train to ensure consistent results.  Report the best step-size choice and the corresponding train and test MSEs. In some cases you might observe that the value of ∆(k) jumps up and down, and this is not something you would have seen using batch GD. Why do you think this might be happening?

What to submit: a single plot, the train and test MSE requested and some commentary.  Include a screen shot of any code used for this section and a copy of your python code in solutions.py.

(h) Based on your GD and SGD results, which algorithm do you prefer? When is it a better idea to use GD? When is it a better idea to use SGD? What to submit: some commentary

(i) Note that in GD, SGD and mini-batch GD, we always update the entire p-dimensional vector β at each iteration. An alternative popular approach is to update each of the p parameters individually.

To make this idea more clear, we write the ridge loss L(β) as L(β1 , β2 . . . , βp ).  We initialize β(0), and then solve fork = 1, 2, 3, . . . ,

βk)  = arg β(m)nL(β1 , βk- 1) , βk- 1) , . . . , βk- 1) )

βk)  = arg β(m)nL(βk), β2 , βk- 1) , . . . , βk- 1) )

.

.

.

βk)  = arg β(m)nL(βk), βk), βk), . . . , βp ).

Note that each of the minimizations is over a single (1-dimensional) coordinate of β, and also that as as soon as we update βk), we use the new value when solving the update for β1  and so on.

The idea is then to cycle through these coordinate level updates until convergence. In the next two parts we will implement this algorithm from scratch for the Ridge regression problem:

L(β) =  |y _ Xβ|2(2) + φ|β|2(2)

Note that we can write then x p matrix X = [X1 , . . . , Xp], where Xj  is the j-th column of X. Find the solution of the optimization

β(ˆ)1  = arg β(m)nL(β1 , β2 , . . . , βp ).

Based on this, derive similar expressions for β(ˆ)j  for j = 2, 3, . . . , p.

Hint: Note the expansion: Xβ = Xj βj  + X-j β -j , where X-j  denotes the matrix X but with the j-th column removed, and similarly β -j  is the vector β with the j-th coordinate removed. What to submit: your working out.

(j) Implement the algorithm outlined in the previous question on the training dataset. In your imple- mentation, be sure to update the βj ’s in order and use an initial estimate of β(0)  = 1p  (th vector of ones), and φ = 0.5. Terminate the algorithm after 10 cycles (one cycle here isp updates, one for each βj ), so you will have a total of 10p updates. Report the train and test MSE of your resulting model. Here we would like to compare the three algorithms: new algorithm to batch GD and SGD from your previous answers with optimally chosen step sizes. Create a plot of k vs. ∆ (k)  as before, but this time plot the progression of the three algorithms. Be sure to use the same colors as indicated here in your plot, and add a legend that labels each series clearly.  For your batch GD and SGD include the step-size in the legend.  Your x-axis only needs to range from k  =  1, . . . 10p.  Further, report both train and test MSE for your new algorithm. Note: Some of you may be concerned that we are comparing one step of GD to one step of SGD and the new aglorithm, we will ignore this technicality for the time being. What to submit: a single plot, the train and test MSE requested.

Question 2

Given λ > 0 and v e R, consider the following optimization problem:

βe(mi) , |β| +  (β _ v)2 } .