5.10 Gradient descent, momentum and Newton's method

To minimise a function, repeatedly step against its gradient. On a quadratic the error after each step is multiplied by the matrix I minus eta times the Hessian, so the speed is set by the Hessian's eigenvalues and the condition number of 3.5, and the path zig-zags across narrow valleys. Momentum improves the rate from (k-1)/(k+1) to (sqrt k - 1)/(sqrt k + 1), stochastic gradients trade exactness for cost and leave a noise floor proportional to the step, Newton's method divides by the Hessian and converges quadratically near a minimum, and in high dimensions most critical points are saddles.

You are standing somewhere on a hillside in thick fog and want to get to the bottom. You cannot see the valley, but you can feel which way is steepest underfoot. So you step downhill, feel again, and step again. That is gradient descent, and by 5.2 the gradient is exactly the instrument that feels the steepest direction. Almost every large model in machine learning is trained this way. Whether the walk reaches the bottom, how fast, and what goes wrong (overshooting, zig-zagging, getting stuck at a saddle) is decided by the second derivative: the Hessian of 5.6, its eigenvalues from 3.1, and the condition number of 3.5. This chapter proves what can be proved and shows the rest by computation.

Before you start

5.2: the gradient is the direction of steepest ascent. 5.6: the Hessian, positive definiteness and eigenvalues. 3.5: the condition number. 4.2: contractions and fixed points. 5.8: convex functions and the gradient inequality.

Quick checkAt a point

a

with

\nabla f(a) \ne O

, in which direction does

f

decrease fastest, per unit of distance?

Quick checkA symmetric matrix

A

has eigenvalues

1

and

100

. What are its operator norm and condition number?

Quick checkA map

T

\mathbb R^n

satisfies

\lVert T(x) - T(y)\rVert \le 0.9\lVert x - y\rVert

for all

x, y

. What do you know about the iteration

x_{k+1} = T(x_k)

The picture: a ball rolling down a valley

The function is $\frac12(x^2 + 25y^2)$ , a long narrow valley (the ellipses). Drag the start, then try these things in order. (1) Move the step-size slider: at $t = \eta L = 1$ the steep coordinate is cancelled in a single step and the path runs straight along the valley floor; for $t$ between $1$ and $2$ the path zig-zags across the valley, more and more as $t$ approaches $2$ ; above $2$ it diverges. (2) Choose the best step, $t \approx 1.92$ here, and see that even then the path still zig-zags and crawls along the valley floor: the error in the long direction shrinks only by $(\kappa-1)/(\kappa+1) = 24/26$ per step. (3) Turn up the momentum $\beta$ : the orange path damps the zig-zag and runs down the valley much faster. (4) Choose $\kappa = 100$ , then the Rosenbrock banana, where plain descent crawls round a curved valley. (5) Choose the saddle and watch the path approach the saddle point and escape slowly along the downhill direction.

The idea, in words

At a point $x$ the gradient is the best linear summary of $f$ (5.3): $f(x + s) \approx f(x) + \nabla f(x)\cdot s$ . To lower $f$ we want $\nabla f(x)\cdot s$ as negative as possible for a step of given length, and by Cauchy–Schwarz that is $s = -\eta\nabla f(x)$ . The linear model, however, is only trusted for small steps. How small? The quadratic term of the Taylor formula (5.6) says: as small as the curvature allows. Too long a step along a direction of high curvature overshoots the valley floor and lands higher up the other wall; too short a step along a direction of low curvature barely moves. The Hessian’s largest eigenvalue $L$ limits the step; its smallest, $\mu$ , limits the speed; and the ratio $\kappa = L/\mu$ is the whole story on a quadratic.

Where this came from· Cauchy's note of 1847

The method of steepest descent is usually traced to a three-page note by Augustin-Louis Cauchy, “Méthode générale pour la résolution des systèmes d’équations simultanées”, read to the Paris Academy of Sciences in October 1847. His motive was astronomical: to find the orbit of a heavenly body one has to solve simultaneous algebraic equations for the orbital elements, and eliminating the unknowns one at a time is hopeless by hand. He proposed to turn the system into the minimisation of a non-negative function of the unknowns (the sum of the squares of the equations is the natural choice) and to change the unknowns step by step in the direction in which that function decreases fastest, then repeat. It is the iteration $x_{k+1} = x_k - \eta_k\nabla f(x_k)$ of this chapter. The same iteration, with the gradient computed by backpropagation (3.8), trains neural networks today.

Definitions

Definition A

Gradient descent

To minimise a differentiable scalar field $f$ on $\mathbb R^n$ , choose a start $x_0$ and a number $\eta \gt 0$ and iterate $x_{k+1} = x_k - \eta\,\nabla f(x_k)$ . This is gradient descent. Its fixed points are exactly the critical points of $f$ ( $\nabla f = O$ ), and the next sections ask when the iterates converge to one, and how fast.

In one lineThe iteration $x_{k+1} = x_k - \eta\,\nabla f(x_k)$ , which moves from the current point against the gradient by a step proportional to its length.

Why it existsIt is the simplest rule that uses first-order information only and is guaranteed, for a small enough step, to lower the function at every step.

Definition A

Learning rate

The positive number $\eta$ multiplying the gradient is the learning rate (also step size). It may be fixed or change with the step number $k$ .

In one lineThe step-size multiplier $\eta$ in gradient descent; in machine learning it is the most important hyperparameter.

Why it existsIt decides between crawling and exploding: the theory below says the useful range is below $2/L$ , where $L$ is the largest Hessian eigenvalue.

Definition A

Loss function

For a model with parameters $w \in \mathbb R^n$ and data $(x_i, y_i)_{i = 1}^N$ , a loss function is a field $\ell(w) = \frac1N\sum_{i = 1}^N\ell_i(w)$ , the average of the loss $\ell_i$ on each example. For least squares, $\ell_i(w) = \frac12(w\cdot x_i - y_i)^2$ .

In one lineThe scalar field $\ell(w)$ of the parameters $w$ of a model that measures how badly the model fits the data; training means minimising it.

Why it existsIt turns 'fit the data' into the optimisation problem this chapter solves; for least squares it is 2.10's sum of squared errors.

One step lowers the function: the descent lemma

Lemma B standard (Nesterov, Introductory Lectures on Convex Optimization, 1.2.3)

Descent lemma

Let $f$ be continuously differentiable on $\mathbb R^n$ with $\lVert\nabla f(x) - \nabla f(y)\rVert \le L\lVert x - y\rVert$ for all $x, y$ (we say $f$ is $L$ -smooth; for a $C^2$ field this holds when all Hessian eigenvalues lie in $[-L, L]$ , by the mean value inequality of 5.7 applied to $\nabla f$ ). Then $f(y) \le f(x) + \nabla f(x)\cdot(y - x) + \dfrac L2\lVert y - x\rVert^2$ . If $y = x - \eta\nabla f(x)$ with $0 \lt \eta \le \frac1L$ then $f(y) \le f(x) - \dfrac\eta2\lVert\nabla f(x)\rVert^2$ .

In one lineIf $\nabla f$ is $L$ -Lipschitz, then $f(y) \le f(x) + \nabla f(x)\cdot(y - x) + \frac L2\lVert y - x\rVert^2$ ; consequently the step $y = x - \eta\nabla f(x)$ with $\eta \le \frac1L$ gives $f(y) \le f(x) - \frac\eta2\lVert\nabla f(x)\rVert^2$ .

Why it is trueIt is the guarantee that gradient descent makes progress, and it says how much: the decrease is at least half the step times the squared gradient.

ProofNesterov 1.2.3

The idea: Write the change in f along the segment as an integral of the directional derivative, subtract the linear part, and bound the difference of gradients by the Lipschitz constant times the distance.

f(y) - f(x) = \displaystyle\int_0^1\nabla f\big(x + t(y - x)\big)\cdot(y - x)\,dt

Put

g(t) = f(x + t(y - x))

. By the path theorem of 5.3,

g'(t) = \nabla f(x + t(y - x))\cdot(y - x)

, continuous; fundamental theorem of calculus (1.5).

f(y) - f(x) - \nabla f(x)\cdot(y - x) = \displaystyle\int_0^1\big[\nabla f(x + t(y - x)) - \nabla f(x)\big]\cdot(y - x)\,dt

Subtract

\nabla f(x)\cdot(y - x) = \int_0^1\nabla f(x)\cdot(y - x)\,dt

The integrand is at most

Lt\lVert y - x\rVert^2

, and

\int_0^1Lt\,dt = \frac L2

. This is the first inequality.

Cauchy–Schwarz on each integrand and the Lipschitz bound:

\lVert\nabla f(x + t(y - x)) - \nabla f(x)\rVert \le Lt\lVert y - x\rVert

f(y) \le f(x) - \dfrac\eta2\lVert\nabla f(x)\rVert^2

Put

y - x = -\eta\nabla f(x)

: the linear term is

-\eta\lVert\nabla f\rVert^2

and the quadratic term is

\frac L2\eta^2\lVert\nabla f\rVert^2

. So

f(y) \le f(x) - \eta\big(1 - \frac{L\eta}2\big)\lVert\nabla f(x)\rVert^2

, and

1 - \frac{L\eta}2 \ge \frac12

when

\eta \le \frac1L

\blacksquare

Gradient descent on a quadratic: the eigenvalues decide

A quadratic is the model of every smooth function near a minimum (the second-order Taylor formula of 5.6), so this case matters.

Theorem B standard (Boyd and Vandenberghe 9.3; Nocedal and Wright 3.3)

Gradient descent on a quadratic

Let $A$ be symmetric positive definite with eigenvalues $0 \lt \mu = \lambda_1 \le \dots \le \lambda_n = L$ and $f(x) = \frac12x^{\mathsf T}Ax - b^{\mathsf T}x$ , with unique minimiser $x^* = A^{-1}b$ ( $\nabla f(x) = Ax - b$ ). For gradient descent with step $\eta$ , $x_k - x^* = (I - \eta A)^k(x_0 - x^*)$ . (a) The iteration converges to $x^*$ for every $x_0$ if and only if $0 \lt \eta \lt 2/L$ . (b) Then $\lVert x_{k} - x^*\rVert \le \rho^k\lVert x_0 - x^*\rVert$ with $\rho(\eta) = \max\big(|1 - \eta\mu|, |1 - \eta L|\big)$ , and this rate is attained. (c) $\rho$ is smallest at $\eta^* = \dfrac2{L + \mu}$ , where $\rho^* = \dfrac{L - \mu}{L + \mu} = \dfrac{\kappa - 1}{\kappa + 1}$ .

In one lineFor $f(x) = \frac12x^{\mathsf T}Ax - b^{\mathsf T}x$ with $A$ symmetric positive definite with eigenvalues in $[\mu, L]$ , gradient descent converges to $x^* = A^{-1}b$ for every start if and only if $0 \lt \eta \lt 2/L$ ; the error is multiplied by $\max(|1 - \eta\mu|, |1 - \eta L|)$ per step, which is smallest, equal to $\frac{\kappa - 1}{\kappa + 1}$ with $\kappa = L/\mu$ , at $\eta = \frac2{L + \mu}$ .

Why it is trueIt turns the speed of the algorithm into a number computed from the Hessian's extreme eigenvalues, and shows that it is the condition number of 3.5, not the size of the problem, that decides.

ProofBoyd and Vandenberghe 9.3

The idea: Subtract the fixed point from the update to get a linear recursion for the error. Diagonalise A with the spectral theorem: each coordinate of the error in the eigenbasis is multiplied by 1 minus eta times its eigenvalue at every step, so the slowest coordinate sets the rate. Choose eta to balance the two extreme coordinates.

e_{k+1} = (I - \eta A)e_k

for

e_k = x_k - x^*

, hence

e_k = (I - \eta A)^ke_0

\nabla f(x^*) = Ax^* - b = O

, so

x_{k+1} - x^* = x_k - x^* - \eta(Ax_k - b) = x_k - x^* - \eta A(x_k - x^*)

In the eigenbasis, coordinate

i

of the error is multiplied by

(1 - \eta\lambda_i)

at each step:

c_i^{(k)} = (1 - \eta\lambda_i)^kc_i^{(0)}

Spectral theorem (3.1):

A = Q\Lambda Q^{\mathsf T}

with

Q

orthogonal. Then

I - \eta A = Q(I - \eta\Lambda)Q^{\mathsf T}

, and

Q

does not change lengths.

Every coordinate tends to

0

for every start if and only if

|1 - \eta\lambda_i| \lt 1

for all

i

, that is

0 \lt \eta \lt 2/L

. This is (a).

A geometric sequence

r^k

tends to

0

if and only if

|r| \lt 1

. The condition

|1 - \eta\lambda_i| \lt 1

0 \lt \eta\lambda_i \lt 2

. It must hold for the largest eigenvalue

L

; then it holds for all smaller ones.

\lVert e_k\rVert \le \rho^k\lVert e_0\rVert

with

\rho = \max(|1 - \eta\mu|, |1 - \eta L|)

, attained for a suitable start. This is (b).

\lVert e_k\rVert^2 = \sum_i(1 - \eta\lambda_i)^{2k}\big(c_i^{(0)}\big)^2 \le \max_i|1 - \eta\lambda_i|^{2k}\lVert e_0\rVert^2

, with equality when

e_0

is an eigenvector for the worst eigenvalue. The function

\lambda \mapsto |1 - \eta\lambda|

is convex on

[\mu, L]

, so its maximum is at an end.

\eta^* = \dfrac2{L + \mu}

and then

\rho^* = 1 - \dfrac{2\mu}{L + \mu} = \dfrac{L - \mu}{L + \mu}

\blacksquare

\eta

grows from

0

|1 - \eta\mu| = 1 - \eta\mu

decreases and

|1 - \eta L| = \eta L - 1

(once

\eta \gt 1/L

) increases; the maximum of the two is smallest where they are equal:

1 - \eta\mu = \eta L - 1

Why· Why the path zig-zags (the picture behind the proof)

At the best step $\eta^*$ the error coordinate along the eigenvector for $\mu$ is multiplied by $1 - \eta^*\mu = +\frac{\kappa - 1}{\kappa + 1}$ : it shrinks slowly and keeps its sign. The coordinate along the eigenvector for $L$ is multiplied by $1 - \eta^*L = -\frac{\kappa - 1}{\kappa + 1}$ : it shrinks equally slowly but flips sign at every step. So the path moves steadily along the valley while bouncing from one wall to the other, and the same number $\frac{\kappa - 1}{\kappa + 1}$ , close to $1$ when $\kappa$ is large, limits both. Any larger step would make the bouncing grow; any smaller step would slow the long direction. This is the blue path in the picture above.

Proof workshopGive the reason for each line: gradient descent on a quadraticBoyd and Vandenberghe 9.3

Each line of the proof is shown in order. Choose the reason that justifies it, then press Check.

$e_{k+1} = (I - \eta A)e_k$ where $e_k = x_k - x^*$ .
In an orthonormal eigenbasis, coordinate $i$ is multiplied by $1 - \eta\lambda_i$ at each step.
Convergence for all starts if and only if $|1 - \eta\lambda_i| \lt 1$ for every $i$ , that is $0 \lt \eta \lt 2/L$ .
The error shrinks like $\rho^k$ with $\rho = \max(|1 - \eta\mu|, |1 - \eta L|)$ .
The best step solves $1 - \eta\mu = \eta L - 1$ , giving $\eta^* = \frac{2}{L + \mu}$ and $\rho^* = \frac{\kappa - 1}{\kappa + 1}$ .

Worked example new (hand computation of the picture above)

Run gradient descent on

f(x, y) = \frac12(x^2 + 10y^2)

from

(4, 1)

with

\eta = 0.15

, and find the stable range of

\eta

and the best step.

Solution.

Gradient and eigenvalues. $\nabla f = (x, 10y)$ and the Hessian is $\mathrm{diag}(1, 10)$ , so $\mu = 1$ , $L = 10$ , $\kappa = 10$ . The matrix is diagonal, so its eigenvalues are the diagonal entries.
One step. $x_{k+1} = (1 - \eta)x_k = 0.85x_k$ and $y_{k+1} = (1 - 10\eta)y_k = -0.5y_k$ . From $(4, 1)$ : $(3.4, -0.5)$ , then $(2.89, 0.25)$ , then $(2.4565, -0.125)$ . Each coordinate has its own multiplier; the $y$ -multiplier is negative, so $y$ flips sign every step: the zig-zag.
Stable range. Both multipliers must have absolute value below $1$ : $|1 - \eta| \lt 1$ and $|1 - 10\eta| \lt 1$ , so $0 \lt \eta \lt 0.2 = 2/L$ . Theorem (a).
Best step. $\eta^* = \dfrac2{10 + 1} = \dfrac2{11} \approx 0.1818$ , with rate $\rho^* = \dfrac{10 - 1}{10 + 1} = \dfrac9{11} \approx 0.818$ . Both coordinates then shrink by $0.818$ per step in absolute value: $x$ by $+\frac9{11}$ and $y$ by $-\frac9{11}$ .
Compare with $\eta = 0.15$ . The multipliers are $0.85$ and $-0.5$ , so $\rho = 0.85$ , slower than $0.818$ because the $x$ direction now limits. For $\eta \lt \eta^*$ the long direction (small eigenvalue) is the slow one.

Answer. Stable for $0 \lt \eta \lt 0.2$ ; the best step is $\frac2{11}$ with rate $\frac9{11}$ . After $k$ steps at $\eta = 0.15$ the error is about $0.85^k$ of the start.

Your answerFor gradient descent on

f = \frac12(x^2 + 10y^2)

, what is the largest learning rate

\eta

for which the iteration converges from every start?

Your answerA quadratic has Hessian eigenvalues

1

and

100

. What is the best convergence factor per step

(\kappa - 1)/(\kappa + 1)

of gradient descent with a fixed step? (Give a decimal.)

Your answerWith the factor

\rho = 99/101

per step, about how many steps shrink the error by a factor of

10^6

? (Give a decimal.)

Worked example new (links the normal equations of 2.10 to this chapter)

Least squares by gradient descent. Fit

w

to minimise

\frac12\lVert Xw - y\rVert^2

for

X = \begin{bmatrix}1 & 0\\ 0 & 1\\ 1 & 1\end{bmatrix}

and

y = (1, 2, 2)

. Find the best step and the first iterate from

w_0 = O

Solution.

The Hessian is $X^{\mathsf T}X$ . $\nabla\frac12\lVert Xw - y\rVert^2 = X^{\mathsf T}(Xw - y)$ (3.8), so the Hessian is $A = X^{\mathsf T}X = \begin{bmatrix}2 & 1\\ 1 & 2\end{bmatrix}$ . The loss is a quadratic with this matrix.
Eigenvalues. $\mu = 1$ (eigenvector $(1, -1)$ ) and $L = 3$ (eigenvector $(1, 1)$ ), so $\kappa = 3$ . Trace $4$ , determinant $3$ .
The minimiser. The normal equations $Aw = X^{\mathsf T}y = (3, 4)$ give $w^* = \big(\frac23, \frac53\big)$ . 2.10; solve the $2 \times 2$ system.
Best step and rate. $\eta^* = \dfrac2{3 + 1} = \dfrac12$ , $\rho^* = \dfrac{3 - 1}{3 + 1} = \dfrac12$ . The theorem.
First iterate. $w_1 = w_0 - \eta^*X^{\mathsf T}(Xw_0 - y) = O + \frac12(3, 4) = (1.5, 2)$ . The error is $w_1 - w^* = \big(\frac56, \frac13\big)$ , which is $\big(I - \frac12A\big)\big(-\frac23, -\frac53\big) = \big(\frac56, \frac13\big)$ . ✓ Check with the error recursion.

Answer. $\eta^* = \frac12$ , rate $\frac12$ per step (the error halves each time), $w_1 = (1.5, 2)$ , heading for $\big(\frac23, \frac53\big)$ . In a realistic problem with a million columns $X^{\mathsf T}X$ is never formed or inverted; the product $X^{\mathsf T}(Xw - y)$ costs two matrix–vector products.

Same idea elsewhere· The conditioning of 3.5, squared

The Hessian of least squares is $X^{\mathsf T}X$ , whose condition number is $\kappa(X)^2$ (its eigenvalues are the squares of the singular values of $X$ , 3.4). A data matrix with $\kappa(X) = 10$ gives $\kappa = 100$ and a rate $\frac{99}{101}$ ; with $\kappa(X) = 1000$ the rate is $1 - 2\times10^{-6}$ , hopeless. This is why features are standardised (columns scaled to the same size) before training: it reduces $\kappa(X)$ and so speeds up gradient descent by orders of magnitude, and it is also why the normal equations lose precision (3.5, 2.10).

Corollary B standard; ties 4.2 and 5.7 to this chapter

Gradient descent is a contraction

Suppose $f \in C^2(\mathbb R^n)$ and every eigenvalue of the Hessian lies in $[\mu, L]$ at every point, $\mu \gt 0$ . Then for $0 \lt \eta \lt 2/L$ the map $T(x) = x - \eta\nabla f(x)$ satisfies $\lVert T(x) - T(y)\rVert \le \rho\lVert x - y\rVert$ , $\rho = \max(|1 - \eta\mu|, |1 - \eta L|) \lt 1$ . Hence $T$ has exactly one fixed point $x^*$ (the point with $\nabla f(x^*) = O$ , the minimiser) and the iterates satisfy $\lVert x_k - x^*\rVert \le \dfrac{\rho^k}{1 - \rho}\lVert x_1 - x_0\rVert$ .

In one lineIf $f$ is $C^2$ with $\mu I \preceq H_f(x) \preceq LI$ everywhere on $\mathbb R^n$ and $0 \lt \eta \lt 2/L$ , then $T(x) = x - \eta\nabla f(x)$ is a contraction with constant $\rho = \max(|1 - \eta\mu|, |1 - \eta L|) \lt 1$ , so gradient descent converges to the unique minimiser $x^*$ and $\lVert x_k - x^*\rVert \le \dfrac{\rho^k}{1 - \rho}\lVert x_1 - x_0\rVert$ .

Why it is trueIt proves convergence for every smooth strongly convex function, and not only quadratics, by the fixed-point theorem of 4.2, and gives a computable error bound.

Proof4.2 and 5.7

The idea: The Jacobian of T is I minus eta times the Hessian, a symmetric matrix with eigenvalues between 1 minus eta L and 1 minus eta mu; its operator norm is at most rho. The mean value inequality for fields turns that into the Lipschitz bound, and the Banach fixed point theorem finishes.

DT(x) = I - \eta H_f(x)

D(\nabla f) = H_f

, the Hessian (5.6).

\lVert DT(x)\rVert \le \rho

for every

x

For a symmetric matrix the operator norm is the largest absolute eigenvalue (3.5). The eigenvalues of

I - \eta H

are

1 - \eta\lambda

with

\lambda \in [\mu, L]

, and

|1 - \eta\lambda| \le \rho

for all such

\lambda

\lVert T(x) - T(y)\rVert \le \rho\lVert x - y\rVert

, and

\rho \lt 1

because

0 \lt \eta \lt 2/L

Mean value inequality for fields (5.7) on the convex set

\mathbb R^n

T

has a unique fixed point

x^*

, the iterates converge to it, and the bound holds.

Banach fixed point theorem (4.2) in the complete space

\mathbb R^n

, with its a priori error bound

d(x_k, x^*) \le \frac{\rho^k}{1 - \rho}d(x_1, x_0)

. A fixed point of

T

satisfies

\eta\nabla f(x^*) = O

\blacksquare

Exact line search, and why the zig-zag is built in

If at each step we choose $\eta_k$ to minimise $f$ along the line $x_k - \eta g_k$ (an exact line search, where $g_k = \nabla f(x_k)$ ), then $\frac{d}{d\eta}f(x_k - \eta g_k) = -\nabla f(x_{k+1})\cdot g_k = 0$ : each new gradient is perpendicular to the previous one, $g_{k+1}\perp g_k$ . In two variables the path must turn through a right angle at every step. On a long narrow valley, this forces the zig-zag, however clever the step. For a quadratic with the worst start the error in the norm $\sqrt{e^{\mathsf T}Ae}$ is multiplied by at most $\frac{\kappa - 1}{\kappa + 1}$ per step, and by exactly that for the worst start (the Kantorovich inequality). It is the same number as the rate of the best fixed step, which the lab below measures.

Quick checkGradient descent on

f = \frac12(x^2 + 100y^2)

uses

\eta = 0.025

. What happens?

Quick checkWhich describes the zig-zag of gradient descent in a narrow valley?

Convex functions: progress without a second derivative

Theorem B standard (Nesterov, Introductory Lectures 2.1.5; Boyd and Vandenberghe 9.3)

Gradient descent on a convex smooth function

Let $f$ be convex and differentiable on $\mathbb R^n$ with $L$ -Lipschitz gradient and a minimiser $x^*$ . Gradient descent $x_{i+1} = x_i - \frac1L\nabla f(x_i)$ satisfies $f(x_k) - f(x^*) \le \dfrac{L\lVert x_0 - x^*\rVert^2}{2k}$ for $k \ge 1$ .

In one lineIf $f$ is convex and $L$ -smooth with a minimiser $x^*$ , then gradient descent with $\eta = 1/L$ satisfies $f(x_k) - f(x^*) \le \dfrac{L\lVert x_0 - x^*\rVert^2}{2k}$ .

Why it is trueIt guarantees convergence for any convex smooth problem with no assumption on the smallest eigenvalue: the price is a slower rate (inversely proportional to the number of steps) than the geometric rate on a strongly convex problem.

ProofNesterov 2.1.5

The idea: Track the squared distance to the minimiser. The gradient inequality (convexity) says the gradient points away from the minimiser by at least the function gap; the descent lemma says the function gap controls the squared gradient. Together they make the distance decrease and sum up to a telescoping bound.

\lVert x_{i+1} - x^*\rVert^2 = \lVert x_i - x^*\rVert^2 - 2\eta\,g_i\cdot(x_i - x^*) + \eta^2\lVert g_i\rVert^2

Expand

\lVert x_i - \eta g_i - x^*\rVert^2

with

g_i = \nabla f(x_i)

and

\eta = 1/L

\lVert x_{i+1} - x^*\rVert^2 \le \lVert x_i - x^*\rVert^2 - 2\eta\big(f(x_i) - f(x^*)\big) + \eta^2\lVert g_i\rVert^2

Gradient inequality for convex

f

(5.8) with

y = x^*

f(x^*) \ge f(x_i) + g_i\cdot(x^* - x_i)

, that is

g_i\cdot(x_i - x^*) \ge f(x_i) - f(x^*)

\lVert x_{i+1} - x^*\rVert^2 \le \lVert x_i - x^*\rVert^2 - 2\eta\big(f(x_{i+1}) - f(x^*)\big)

Descent lemma with

\eta = 1/L

f(x_{i+1}) \le f(x_i) - \frac\eta2\lVert g_i\rVert^2

, so

\eta^2\lVert g_i\rVert^2 \le 2\eta\big(f(x_i) - f(x_{i+1})\big)

2\eta\displaystyle\sum_{i = 1}^k\big(f(x_i) - f(x^*)\big) \le \lVert x_0 - x^*\rVert^2

Add the inequalities for

i = 0, \dots, k - 1

: the squared distances telescope, and

\lVert x_k - x^*\rVert^2 \ge 0

k\big(f(x_k) - f(x^*)\big) \le \dfrac{\lVert x_0 - x^*\rVert^2}{2\eta}

By the descent lemma

f(x_i)

is non-increasing, so every term is at least

f(x_k) - f(x^*)

. Divide by

2\eta k = 2k/L

\blacksquare

Proof workshopRebuild the proof that gradient descent on a convex smooth function has error at most L times the squared distance over 2kNesterov 2.1.5

Tap the lines in the order a proof would write them. Each line can only use what is already above it.

Python labGradient descent on an ill-conditioned quadratic: the zig-zag, the rate, and momentum

import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt

kappa = 100.0
A = np.diag([1.0, kappa]); x0 = np.array([100.0, 1.0])          # f(x) = x^T A x / 2, minimum at the origin

def run(step, beta=0.0, n=800):
    x, xp, err = x0.copy(), x0.copy(), []
    for _ in range(n):
        err.append(np.linalg.norm(x))
        x, xp = x - step * (A @ x) + beta * (x - xp), x          # heavy-ball step (beta = 0 is plain descent)
    return np.array(err)

L, mu = kappa, 1.0
plain = run(2 / (L + mu))                                           # best fixed step
mom = run(4 / (np.sqrt(L) + np.sqrt(mu)) ** 2, beta=((np.sqrt(kappa) - 1) / (np.sqrt(kappa) + 1)) ** 2)
rate = lambda e, a, b: (e[b] / e[a]) ** (1 / (b - a))
print("plain: measured rate", round(rate(plain, 300, 400), 4), "  theory (k-1)/(k+1) =", round(99 / 101, 4))
print("momentum: measured rate", round(rate(mom, 100, 150), 4), "  theory (sqrt(k)-1)/(sqrt(k)+1) =", round(9 / 11, 4))
for name, e in (("plain", plain), ("momentum", mom)):
    print(name, "steps to shrink the error by 1e-6:", int(np.argmax(e < 1e-6 * e[0])))
print("step 2.02/L, error after 400 steps:", round(run(2.02 / L, n=400)[-1], 1), "(diverges)")
plt.semilogy(plain, label="plain descent"); plt.semilogy(mom, label="momentum")
plt.xlabel("step"); plt.ylabel("distance to minimum"); plt.legend()

Try this. Run it. Plain descent with the best step has measured rate 0.9802 = 99/101 and needs 691 steps to shrink the error a million times; momentum has measured rate 0.8246, close to the theoretical 9/11 = 0.8182 (the extra is a polynomial factor at the edge of the stable region), and needs 83 steps. Step 2.02/L is slightly above the limit 2/L and the error grows. Change kappa to 10 and then to 1000: the plain step count scales like kappa (69, 691 and 6908 steps: about 6.9 kappa, which is (kappa/2) ln 10^6), and the momentum count like the square root of kappa (26, 83 and 264 steps). Print the first few x values of plain descent to see the y-coordinate flipping sign, and change x0 to (1, 1) to see how the start affects the transient but not the rate.

Momentum

Definition A

Momentum

Momentum (Polyak’s heavy ball) replaces gradient descent by $x_{k+1} = x_k - \eta\nabla f(x_k) + \beta(x_k - x_{k-1})$ with $0 \le \beta \lt 1$ . Equivalently the ‘velocity’ $v_k = x_k - x_{k-1}$ satisfies $v_{k+1} = \beta v_k - \eta\nabla f(x_k)$ and $x_{k+1} = x_k + v_{k+1}$ : a ball with friction rolling in the landscape.

In one lineThe modification $x_{k+1} = x_k - \eta\nabla f(x_k) + \beta(x_k - x_{k-1})$ that adds a fraction $\beta$ of the previous step to the gradient step (the heavy-ball method of Polyak).

Why it existsIt lets the iterate build up speed along a valley floor, where the gradients keep pointing the same way, and cancels the sideways bouncing, where the gradients alternate; on quadratics it improves the rate from $\kappa$ to $\sqrt\kappa$ steps.

Theorem B Polyak 1964; Goh, Why Momentum Really Works (2017)

Momentum on a quadratic

Let $f(x) = \frac12x^{\mathsf T}Ax - b^{\mathsf T}x$ with the eigenvalues of $A$ in $[\mu, L]$ . With $\eta = \dfrac4{(\sqrt L + \sqrt\mu)^2}$ and $\beta = \Big(\dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1}\Big)^2$ the heavy-ball errors satisfy $\lVert e_k\rVert \le C(1 + k)\Big(\dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1}\Big)^k$ .

In one lineFor a quadratic with eigenvalues in $[\mu, L]$ , the heavy-ball method with $\eta = \dfrac4{(\sqrt L + \sqrt\mu)^2}$ and $\beta = \Big(\dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1}\Big)^2$ makes the error decay at the rate $\dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1}$ (up to a factor polynomial in $k$ ).

Why it is trueIt shows that the right amount of memory turns the condition number $\kappa$ into $\sqrt\kappa$ in the iteration count: for $\kappa = 10^4$ , about $100$ times fewer steps.

ProofPolyak 1964; Goh 2017

The idea: Diagonalise as before. In each eigen-direction the error obeys a two-term linear recurrence whose characteristic roots are complex, of the same absolute value square root of beta, whenever one plus beta minus eta times the eigenvalue is at most two root beta in absolute value. The chosen eta and beta make this hold exactly across the whole range of eigenvalues.

For eigenvalue

\lambda

c_{k+1} = (1 + \beta - \eta\lambda)\,c_k - \beta\,c_{k-1}

Subtract the minimiser from the update and use

\nabla f(x) = A(x - x^*)

e_{k+1} = e_k - \eta Ae_k + \beta(e_k - e_{k-1})

. In an orthonormal eigenbasis of

A

(spectral theorem) coordinate

i

evolves independently.

The roots are

z = \dfrac{\tau \pm \sqrt{\tau^2 - 4\beta}}2

with

\tau = 1 + \beta - \eta\lambda

A two-term recurrence

c_{k+1} = \tau c_k - \beta c_{k-1}

has solutions built from the roots of

z^2 - \tau z + \beta = 0

Both roots have absolute value exactly

\sqrt\beta

, so

|c_k| \le C(1 + k)(\sqrt\beta)^k

(the factor

1 + k

appears only when the roots coincide).

|\tau| \le 2\sqrt\beta

the square root is imaginary (or zero) and

|z|^2 = \frac{\tau^2 + (4\beta - \tau^2)}{4} = \beta

So the condition

|\tau| \le 2\sqrt\beta

holds for every eigenvalue in

[\mu, L]

, with

\sqrt\beta = \dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1}

\blacksquare

Write

s = \sqrt\kappa

. Then

\beta = \frac{(s - 1)^2}{(s + 1)^2}

\eta L = \frac{4s^2}{(s + 1)^2}

and

\eta\mu = \frac4{(s + 1)^2}

. So

\tau = 1 + \beta - \eta\lambda

ranges over

[1 + \beta - \eta L,\ 1 + \beta - \eta\mu] = \Big[-\dfrac{2(s - 1)}{s + 1},\ \dfrac{2(s - 1)}{s + 1}\Big] = [-2\sqrt\beta, 2\sqrt\beta]

\lambda

runs over

[\mu, L]

.The recipe for

\eta

and

\beta

was designed to make the two ends of the eigenvalue range land exactly on the edge of the oscillating region.

Worked example

Compare plain descent and momentum for

\kappa = 100

Solution.

Plain descent. $\rho^* = \dfrac{99}{101} = 0.9802$ ; steps for a factor $10^{-6}$ : $\dfrac{\ln10^6}{\ln(101/99)} \approx 691$ . Theorem on quadratics.
Momentum parameters. $\sqrt\kappa = 10$ , so $\beta = \Big(\dfrac9{11}\Big)^2 = \dfrac{81}{121} = 0.669$ and, for $\mu = 1$ , $L = 100$ , $\eta = \dfrac4{(10 + 1)^2} = \dfrac4{121} = 0.0331$ . Formulas of the theorem.
Rate. $\dfrac{\sqrt\kappa - 1}{\sqrt\kappa + 1} = \dfrac9{11} = 0.8182$ ; the ideal count is $\dfrac{\ln10^6}{\ln(11/9)} \approx 69$ steps. Same formula as plain descent with $\sqrt\kappa$ in place of $\kappa$ .
In practice. The lab above measures a rate of $0.8246$ and $83$ steps (the factor $1 + k$ and the transient add some). The bound has a polynomial factor in front.

Answer. About $690$ steps against about $70$ to $85$ : a factor of $8$ to $10$ , which is $\sqrt\kappa$ . For $\kappa = 10^4$ the gap is a factor of $100$ .

Your answerFor a quadratic with

\kappa = 100

, what is the per-step convergence rate

(\sqrt\kappa - 1)/(\sqrt\kappa + 1)

of optimally tuned momentum?

Quick checkIn the heavy-ball method, what happens as the momentum

\beta

approaches

1

with the step size fixed?

Beyond the course: for understanding· Nesterov's acceleration, and conjugate gradients

Yurii Nesterov (1983) modified the method by evaluating the gradient at a point ahead of the current one, $x_k + \beta(x_k - x_{k-1})$ , instead of at $x_k$ . For smooth convex functions (no strong convexity) this improves the bound $O(1/k)$ above to $O(1/k^2)$ , which is optimal among first-order methods; for strongly convex functions it gives a rate of order $1 - 1/\sqrt\kappa$ per step without needing a quadratic. On quadratics the method of conjugate gradients (MATH 340) reaches the same rate by choosing the step and the momentum at each iteration from the data, with no tuning, and finds the exact minimiser in at most $n$ steps in exact arithmetic.

Stochastic gradient descent

A loss in machine learning is an average over data, $\ell(w) = \frac1N\sum_{i = 1}^N\ell_i(w)$ , with $N$ in the millions. Its gradient is the average of $N$ per-example gradients, and computing it costs $N$ backward passes. The idea of stochastic gradient descent is to use the gradient of one randomly chosen term.

Definition A

Stochastic gradient descent

Let $\ell = \frac1N\sum\ell_i$ . Stochastic gradient descent (SGD) is $w_{k+1} = w_k - \eta_k\nabla\ell_{i_k}(w_k)$ where $i_k$ is chosen uniformly at random from $\{1, \dots, N\}$ at each step. Because each index has probability $\frac1N$ , the average of the step over the random choice is $-\eta_k\nabla\ell(w_k)$ : the stochastic gradient is an unbiased estimate of the true gradient. A minibatch of $B$ examples uses the average of $B$ such gradients.

In one lineThe iteration $w_{k+1} = w_k - \eta_k\nabla\ell_{i_k}(w_k)$ in which $i_k$ is an example chosen at random; each step costs one example instead of $N$ .

Why it existsIts step is, on average, the true gradient step, but costs $N$ times less, which is what makes training on huge datasets possible.

Where this came from· Robbins and Monro, 1951

Herbert Robbins and Sutton Monro published “A stochastic approximation method” in the Annals of Mathematical Statistics in 1951. Their problem was not machine learning. To find the level $x$ at which an unknown response function $M(x)$ equals a target value $\alpha$ one can only run experiments, and every experiment returns $M(x)$ plus random noise. They proposed to adjust the level after each experiment by a step proportional to the observed discrepancy, with step sizes $a_k$ shrinking so that the sum of the $a_k$ is infinite (the search can travel anywhere) and the sum of the $a_k^2$ is finite (the noise is eventually averaged out), for example $a_k = 1/k$ , and proved convergence in probability. Stochastic gradient descent is the Robbins–Monro scheme with the response function replaced by a gradient. The two conditions on the step sizes are still the textbook ones.

Theorem B standard; the one-dimensional model of Robbins and Monro (1951)

Noise floor of stochastic gradient descent

Let $a_1, \dots, a_N$ be numbers with mean $\bar a$ and variance $\sigma^2 = \frac1N\sum(a_i - \bar a)^2$ , $\ell_i(w) = \frac12(w - a_i)^2$ , so $\ell(w) = \frac1N\sum\ell_i(w)$ is minimised at $\bar a$ . Run $w_{k+1} = w_k - \eta(w_k - a_{i_k})$ with $i_k$ drawn uniformly and independently. Then with $e_k = w_k - \bar a$ : (a) $\mathrm{E}[e_{k+1}] = (1 - \eta)\mathrm{E}[e_k]$ ; (b) $\mathrm{E}[e_{k+1}^2] = (1 - \eta)^2\mathrm{E}[e_k^2] + \eta^2\sigma^2$ , so for $0 \lt \eta \lt 2$ , $\mathrm{E}[e_k^2] \to \dfrac{\eta\sigma^2}{2 - \eta}$ . (c) With $\eta_k = \frac1{k + 1}$ and $w_0$ arbitrary, $w_k$ is the average of the first $k$ sampled values $a_{i_0}, \dots, a_{i_{k-1}}$ , and $\mathrm{E}[e_k^2] = \sigma^2/k$ .

In one lineFor the loss $\frac1N\sum\frac12(w - a_i)^2$ , whose minimiser is the mean $\bar a$ , SGD with a fixed step $\eta \in (0, 2)$ has $\mathrm{E}[w_k] - \bar a \to 0$ but $\mathrm{E}[(w_k - \bar a)^2] \to \dfrac{\eta\sigma^2}{2 - \eta}$ , where $\sigma^2$ is the variance of the $a_i$ ; with steps $\eta_k = \frac1{k + 1}$ the iterate is the running average and the error tends to $0$ .

Why it is trueIt shows the central trade-off of stochastic methods in one formula: the noise floor is proportional to the learning rate, so a constant step cannot converge exactly, and decreasing steps are needed.

Proofstandard

The idea: Write the error update as a contraction of the old error plus a fresh zero-mean noise term. The noise is independent of the old error, so the variances add.

e_{k+1} = (1 - \eta)e_k + \eta\xi_k

with

\xi_k = a_{i_k} - \bar a

, which has average

0

and average square

\sigma^2

and is independent of

e_k

w_{k+1} - \bar a = (1 - \eta)(w_k - \bar a) + \eta(a_{i_k} - \bar a)

\mathrm{E}[e_{k+1}] = (1 - \eta)\mathrm{E}[e_k]

Take averages:

\mathrm{E}[\xi_k] = 0

. This is (a).

\mathrm{E}[e_{k+1}^2] = (1 - \eta)^2\mathrm{E}[e_k^2] + \eta^2\sigma^2

. This is (b)‘s recursion.

Square the update:

e_{k+1}^2 = (1 - \eta)^2e_k^2 + 2\eta(1 - \eta)e_k\xi_k + \eta^2\xi_k^2

. The cross term has average

\mathrm{E}[e_k]\mathrm{E}[\xi_k] = 0

by independence.

\mathrm{E}[e_k^2] \to \dfrac{\eta^2\sigma^2}{\eta(2 - \eta)} = \dfrac{\eta\sigma^2}{2 - \eta}

The recursion

v_{k+1} = \gamma v_k + c

with

\gamma = (1 - \eta)^2 \lt 1

converges to the fixed point

v = \frac c{1 - \gamma}

. Here

1 - \gamma = 1 - (1 - \eta)^2 = \eta(2 - \eta)

This gives (c).

With

\eta_k = \frac1{k + 1}

w_{k+1} = \frac{k}{k + 1}w_k + \frac1{k + 1}a_{i_k}

, which is the running average (induction). The average of

k

independent samples has variance

\sigma^2/k

\blacksquare

Worked example

A loss has per-example optima

a_i

with variance

\sigma^2 = 4

. What noise floor does SGD with

\eta = 0.1

have, and what does a minibatch of

B = 10

do?

Solution.

Single example. $\mathrm{E}[e^2] \to \dfrac{\eta\sigma^2}{2 - \eta} = \dfrac{0.1\cdot4}{1.9} = 0.2105$ , so the typical error is $\sqrt{0.2105} = 0.459$ . The theorem. The iterate fluctuates around $\bar a$ with this standard deviation even though the step is small.
Minibatch of 10. The average of $10$ independent terms has variance $\sigma^2/10 = 0.4$ , so the floor becomes $\dfrac{0.1\cdot0.4}{1.9} = 0.0211$ , ten times lower. Replace $\sigma^2$ by $\sigma^2/B$ in the recursion.
Halve the step. $\eta = 0.05$ halves the floor (to $\approx 0.103$ ) but also halves the speed of forgetting the start, which is $(1 - \eta)^k$ . The trade-off of the theorem.

Answer. $0.2105$ ; with minibatch $10$ , $0.0211$ . The learning-rate schedules used in practice (large first, then decayed, e.g. $\eta_k \propto 1/k$ ) are the Robbins–Monro conditions in disguise: large steps to forget the start quickly, small steps at the end to average the noise out.

Your answerSGD on

\frac1N\sum\frac12(w - a_i)^2

with

\sigma^2 = 4

and

\eta = 0.1

: what is the long-run

\mathrm{E}[(w_k - \bar a)^2]

? (Give a decimal.)

Quick checkWhy does SGD with a constant learning rate not converge exactly to the minimiser?

Newton’s method in many variables

Gradient descent uses a linear model of $f$ and has to guess how far to trust it. If we also use the second derivative, we can ask for the minimum of the quadratic model.

Definition A

Newton's method

For a root-finding problem $F(x) = O$ with $F: \mathbb R^n \to \mathbb R^n$ and Jacobian matrix $J_F$ , Newton’s method is $x_{k+1} = x_k - J_F(x_k)^{-1}F(x_k)$ , that is, solve the linear system $J_F(x_k)\,s = -F(x_k)$ and set $x_{k+1} = x_k + s$ : replace $F$ by its linearization (5.3) and solve exactly. To minimise $f$ take $F = \nabla f$ , so $J_F = H_f$ , the Hessian (5.6). In one variable this is the Newton iteration of 4.2.

In one lineThe iteration $x_{k+1} = x_k - H_f(x_k)^{-1}\nabla f(x_k)$ for minimising $f$ , which is Newton's method for solving the system $\nabla f = O$ , and jumps to the minimum of the second-order Taylor model at each step.

Why it existsIt corrects the step in each direction by the curvature in that direction, so it does not care how badly conditioned the problem is, and near a minimum it converges quadratically.

Why· Two ways to see the same step

(1) Linearise the gradient. $\nabla f(x + s) \approx \nabla f(x) + H_f(x)s$ ; setting it to $O$ gives $s = -H_f(x)^{-1}\nabla f(x)$ . (2) Minimise the quadratic model. The second-order Taylor formula (5.6) says $f(x + s) \approx f(x) + \nabla f(x)\cdot s + \frac12s^{\mathsf T}H_f(x)s$ ; if $H_f(x)$ is positive definite, this is a quadratic in $s$ with a unique minimiser, found where its gradient $\nabla f + H_fs$ is $O$ : the same $s$ . On a quadratic function the model is the function, so Newton’s method finds the minimum in one step, from any start, whatever the condition number: $x - A^{-1}(Ax - b) = A^{-1}b$ .

Theorem B standard (Nocedal and Wright 3.5; Kantorovich); 376 HW7 #3 for one variable

Quadratic convergence of Newton's method

Let $F: \mathbb R^n \to \mathbb R^n$ be $C^1$ near a root $x^*$ , with $J_F(x^*)$ invertible, $\beta_0 = \lVert J_F(x^*)^{-1}\rVert$ , and $\lVert J_F(x) - J_F(y)\rVert \le M\lVert x - y\rVert$ for $x, y$ in a ball around $x^*$ containing the iterates. If $\lVert x_0 - x^*\rVert \le r := \dfrac1{2\beta_0M}$ , then every iterate is defined, and $\lVert x_{k+1} - x^*\rVert \le \beta_0M\lVert x_k - x^*\rVert^2 \le \frac12\lVert x_k - x^*\rVert$ .

In one lineIf $F(x^*) = O$ , $J_F(x^*)$ is invertible with $\lVert J_F(x^*)^{-1}\rVert \le \beta_0$ and $J_F$ is $M$ -Lipschitz near $x^*$ , then for $\lVert x_0 - x^*\rVert \le \dfrac1{2\beta_0M}$ the Newton iterates satisfy $\lVert x_{k+1} - x^*\rVert \le \beta_0M\lVert x_k - x^*\rVert^2$ : the number of correct digits doubles at each step.

Why it is trueIt explains the speed seen in practice (five or six steps to machine precision) and says exactly what 'close enough to start' means, in terms of the curvature of $F$ and the invertibility of its Jacobian.

ProofNocedal and Wright 3.5

The idea: Subtract the root from the Newton update and rewrite the right side as the inverse Jacobian applied to the remainder of the linear approximation. The remainder is at most M over two times the squared error, and the Jacobian is invertible with norm at most twice its norm at the root, close to the root.

The matrix

J(x_k)

is invertible:

\lVert(I + E)v\rVert \ge \lVert v\rVert - \frac12\lVert v\rVert = \frac12\lVert v\rVert

, so

I + E

is one-to-one, hence invertible, and

\lVert(I + E)^{-1}\rVert \le 2

. Therefore

\lVert J(x_k)^{-1}\rVert \le 2\beta_0

Put

e = x_k - x^*

with

\lVert e\rVert \le r

\lVert J(x_k) - J(x^*)\rVert \le M\lVert e\rVert \le \frac1{2\beta_0}

. Then

J(x_k) = J(x^*)\big(I + E\big)

with

E = J(x^*)^{-1}(J(x_k) - J(x^*))

and

\lVert E\rVert \le \beta_0\cdot\frac1{2\beta_0} = \frac12

x_{k+1} - x^* = J(x_k)^{-1}R

with

R = F(x^*) - F(x_k) - J(x_k)(x^* - x_k)

x_{k+1} - x^* = x_k - x^* - J(x_k)^{-1}F(x_k) = J(x_k)^{-1}\big[J(x_k)(x_k - x^*) - F(x_k)\big]

and

F(x^*) = O

, so the bracket equals

+\big[F(x^*) - F(x_k) - J(x_k)(x^* - x_k)\big]

\lVert R\rVert \le \dfrac M2\lVert e\rVert^2

Fundamental theorem of calculus along the segment (5.3):

F(x^*) - F(x_k) = \int_0^1J(x_k + t(x^* - x_k))(x^* - x_k)\,dt

. So

R = \int_0^1\big[J(x_k + t(x^* - x_k)) - J(x_k)\big](x^* - x_k)\,dt

, and the Lipschitz bound gives

\lVert R\rVert \le \int_0^1Mt\lVert e\rVert^2\,dt = \frac M2\lVert e\rVert^2

\lVert x_{k+1} - x^*\rVert \le \beta_0M\lVert e\rVert^2 \le \frac12\lVert e\rVert

Multiply the bounds:

2\beta_0\cdot\frac M2\lVert e\rVert^2 = \beta_0M\lVert e\rVert^2

. And

\beta_0M\lVert e\rVert \le \beta_0Mr = \frac12

. Induction keeps

x_{k+1}

inside the ball.

\blacksquare

Worked example MATH 376 HW7 #3 and Worksheet 19 #2 (one variable)

Heron's method for

\sqrt2

is Newton's method for

x^2 - 2 = 0

. Show that the error is squared at each step, and run four steps from

x_0 = 1.5

Solution.

The Newton map. For $f(x) = x^2 - 2$ : $N(x) = x - \dfrac{x^2 - 2}{2x} = \dfrac12\Big(x + \dfrac2x\Big)$ . One-variable Newton: subtract $f/f'$ . The fixed point is $\sqrt2$ .
The error identity. $N(x) - \sqrt2 = \dfrac{x^2 - 2\sqrt2\,x + 2}{2x} = \dfrac{(x - \sqrt2)^2}{2x}$ . Put everything over the common denominator $2x$ ; the numerator is a perfect square.
Read it. With $e_k = x_k - \sqrt2$ : $e_{k+1} = \dfrac{e_k^2}{2x_k} \approx 0.354\,e_k^2$ once $x_k \approx \sqrt2$ . The error is squared, times the constant $\frac{f''}{2f'} = \frac1{2\sqrt2}$ at the root. This is the theorem above with $\beta_0 = \frac1{2\sqrt2}$ and $M = 2$ .
Run it. $x_0 = 1.5$ has $e_0 = 0.0858$ . Then $x_1 = 1.416667$ ( $e_1 = 2.5\times10^{-3}$ ), $x_2 = 1.41421569$ ( $e_2 = 2.1\times10^{-6}$ ), $x_3 = 1.41421356237469$ ( $e_3 = 1.6\times10^{-12}$ ), and $x_4$ is $\sqrt2$ to machine precision. Check one step: $0.0858^2/3 = 0.00245$ . The number of correct digits goes $1, 2, 5, 11, 16$ : roughly doubling.

Answer. $e_{k+1} = e_k^2/(2x_k)$ ; four steps give $\sqrt2$ to full double precision.

Worked example new (the standard test problem)

Newton's method on the Rosenbrock function

f(x, y) = (1 - x)^2 + 100(y - x^2)^2

from

(-1.2, 1)

Solution.

Gradient and Hessian. $\nabla f = \big(-2(1 - x) - 400x(y - x^2),\ 200(y - x^2)\big)$ and $H_f = \begin{bmatrix}2 - 400(y - x^2) + 800x^2 & -400x\\ -400x & 200\end{bmatrix}$ . Direct differentiation; at the minimum $(1, 1)$ , $H = \begin{bmatrix}802 & -400\\ -400 & 200\end{bmatrix}$ with eigenvalues $0.4$ and $1001.6$ , so $\kappa \approx 2500$ .
First step. At $(-1.2, 1)$ : $f = 24.2$ , $\nabla f = (-215.6, -88)$ , $H = \begin{bmatrix}1330 & 480\\ 480 & 200\end{bmatrix}$ . Solving $Hs = -\nabla f$ gives $s = (0.0247, 0.3807)$ and $x_1 = (-1.1753, 1.3807)$ , $f = 4.73$ . Solve a $2 \times 2$ linear system.
Second step. From $x_1$ Newton jumps to $(0.7631, -3.1750)$ with $f = 1412$ : the function went up. Far from the minimum the quadratic model is a poor guide, and Newton’s method only aims for a critical point, not for descent.
Then convergence. $x_3 = (0.7634, 0.5828)$ , $f = 0.056$ ; $x_4 = (0.99999, 0.944)$ ; $x_5 = (0.999996, 0.999991)$ with $f = 1.9\times10^{-11}$ ; $x_6 = (1, 1)$ to machine precision. Once close, the error is squared at each step: $0.056 \to 0.31 \to 10^{-11} \to 10^{-20}$ in $f$ .
Compare. Plain gradient descent with a step near the stability limit ( $\eta = 10^{-3}$ ; the largest Hessian eigenvalue is $1506$ at the start, so $2/L = 0.0013$ ) has $f = 0.0038$ after $5000$ steps and $1.9\times10^{-8}$ after $20000$ . The valley is curved and narrow ( $\kappa \approx 2500$ at the minimum), and the step is limited by the steepest wall.

Answer. Newton converges in $6$ steps against tens of thousands for gradient descent, but its path is erratic (step 2 lands at $f = 1412$ ). Practical implementations add a line search or a trust region (they cut the step when the function rises) and fall back on a descent direction when the Hessian is not positive definite.

Watch out· Newton's method is attracted to every critical point, saddles and maxima included

For $f(x) = -x^2$ (a maximum at $0$ ), Newton gives $x - \frac{-2x}{-2} = x - x = 0$ in one step: it lands exactly on the maximum. For $f = x^2 - y^2$ it jumps to the saddle at the origin. It is solving $\nabla f = O$ , which has no preference between minima, maxima and saddles. Gradient descent, in contrast, is repelled from maxima and almost never settles at saddles.

Same idea elsewhere· Cost: why neural networks do not use Newton's method

Newton needs the Hessian, $n^2$ numbers, and a solve with it, about $n^3$ operations. A network with $n = 10^6$ weights has a Hessian with $10^{12}$ entries, which in 4-byte floats is $4$ terabytes. Quasi-Newton methods (BFGS and its low-memory form L-BFGS) build an approximation of $H^{-1}$ from the last few gradients and are used for smaller problems; for deep networks, ‘diagonal’ rescalings of the gradient (Adam, Kingma and Ba, 2014) approximate the idea cheaply, dividing each coordinate’s step by a running estimate of its gradient size.

Quick checkNewton's method is applied to

f(x) = -x^2 + 3

starting at

x_0 = 5

. Where does it go?

Training a small neural network

A neural network is a composition of differentiable maps, so by the chain rule of 5.4 its loss is a differentiable field of the weights, and its gradient is computed by backpropagation (3.8): the computational graph is traversed backwards, multiplying Jacobians. Everything in this chapter applies. As a concrete case, a tiny network with $2$ inputs, $3$ tanh hidden units and $1$ linear output ( $13$ parameters) learns XOR: the four inputs $(0, 0), (0, 1), (1, 0), (1, 1)$ with targets $-1, +1, +1, -1$ . No linear model can do this (the points are not linearly separable), but a hidden layer can.

Python labTrain a 2-3-1 network on XOR with backpropagation and gradient descent

import numpy as np

X = np.array([[0., 0.], [0., 1.], [1., 0.], [1., 1.]]); y = np.array([-1., 1., 1., -1.])   # XOR

def loss_grads(P):
    W1, b1, w2, b2 = P
    h = np.tanh(X @ W1 + b1); out = h @ w2 + b2; d_out = (out - y) / len(y)
    d_pre = np.outer(d_out, w2) * (1 - h ** 2)        # backpropagation (3.8): chain rule from the loss backwards
    return 0.5 * np.mean((out - y) ** 2), out, [X.T @ d_pre, d_pre.sum(0), h.T @ d_out, d_out.sum(keepdims=True)]

def init(seed, H=3):
    rng = np.random.default_rng(seed)
    return [rng.normal(0, 1, (2, H)), np.zeros(H), rng.normal(0, 1, H), np.zeros(1)]

P = init(1); G = loss_grads(P)[2]; worst = 0.0
for p, g in zip(P, G):                                 # gradient check by central differences, one parameter at a time
    for i in np.ndindex(p.shape):
        old = p[i]; p[i] = old + 1e-6; up = loss_grads(P)[0]; p[i] = old - 1e-6; dn = loss_grads(P)[0]; p[i] = old
        worst = max(worst, abs((up - dn) / 2e-6 - g[i]))
print("backprop vs central differences, largest error over the 13 parameters:", f"{worst:.1e}")

def train(seed, eta=0.1, steps=3000, H=3, show=False):
    P = init(seed, H)
    for k in range(steps + 1):
        loss, out, G = loss_grads(P)
        if show and k in (0, 10, 100, 200, 300, 400, 500, 1000, 3000):
            print(f"{k:5d}  loss {loss:.3e}  |grad| {np.sqrt(sum(np.sum(g**2) for g in G)):.3e}")
        for p, g in zip(P, G): p -= eta * g
    return loss, out

train(1, show=True)
print("final loss for seeds 0..9:", [f"{train(s)[0]:.0e}" for s in range(10)])
loss, out = train(2); print("seed 2 stalls at loss", round(loss, 3), "with outputs", out.round(2))

Try this. Run it. The gradient check at the top compares the backpropagation gradient with central differences for all 13 parameters and agrees to about 1e-10: backpropagation is the chain rule applied to the whole network. The loss then falls slowly for about 200 steps (a plateau where the gradient is small), drops fast, and ends below 1e-20: the network represents XOR exactly. The last line shows a failure: for seed 2 the same code is still at loss 0.25 after 3000 steps, with two of the four outputs at 0. Training a non-convex loss depends on the starting point. Try eta = 1.0: eight of the ten seeds still converge (and seed 2 now escapes its stall), but seeds 0 and 7 give nan (the step is above the stability limit 2/L near their start). Try H = 2 and H = 8 hidden units and count how many of the ten seeds succeed.

Saddle points in high dimension

Critical points of a smooth function (5.6) are classified by the signs of the Hessian’s eigenvalues. A strict local minimum needs all $n$ of them positive. For a loss with a million parameters this is an extremely demanding requirement, and a critical point with both positive and negative eigenvalues is a saddle.

Proposition B standard (the exact solution of the linear case)

Gradient descent near a saddle

Let $f(x, y) = \frac12(ax^2 - \gamma y^2)$ with $a, \gamma \gt 0$ and $\eta \lt 2/a$ . Gradient descent gives $x_{k+1} = (1 - \eta a)x_k$ and $y_{k+1} = (1 + \eta\gamma)y_k$ . So $x_k \to 0$ and $y_k = (1 + \eta\gamma)^ky_0$ , which grows without bound unless $y_0 = 0$ . The number of steps to grow from $|y_0|$ to $1$ is $\dfrac{\ln(1/|y_0|)}{\ln(1 + \eta\gamma)} \approx \dfrac{\ln(1/|y_0|)}{\eta\gamma}$ .

In one lineFor $f = \frac12(ax^2 - \gamma y^2)$ with $a, \gamma \gt 0$ , gradient descent with step $\eta$ has $y_k = (1 + \eta\gamma)^ky_0$ : it leaves the saddle exponentially fast but with a time of order $\dfrac1{\eta\gamma}\ln\dfrac1{|y_0|}$ , which is long when the negative curvature $\gamma$ is small.

Why it is trueIt shows why training stalls near saddles: where the downhill curvature is small the gradient is tiny and escape takes many steps, though it always happens unless the start is exactly on the stable direction.

Proofstandard

The idea: The gradient is (a x, minus gamma y), so each coordinate is multiplied by a constant at every step. The second constant exceeds one.

x_{k+1} = (1 - \eta a)x_k

and

y_{k+1} = (1 + \eta\gamma)y_k

\nabla f = (ax, -\gamma y)

, so the update

x_{k+1} = x_k - \eta ax_k

and

y_{k+1} = y_k + \eta\gamma y_k

x_k \to 0

and

|y_k| = (1 + \eta\gamma)^k|y_0| \to \infty

y_0 \ne 0

Geometric sequences:

|1 - \eta a| \lt 1

for

0 \lt \eta \lt 2/a

, while

1 + \eta\gamma \gt 1

k = \dfrac{\ln(1/|y_0|)}{\ln(1 + \eta\gamma)}

Solve

(1 + \eta\gamma)^k|y_0| = 1

for

k

, and use

\ln(1 + u) \approx u

for small

u

\blacksquare

Worked example

A saddle with weak negative curvature:

f = \frac12(x^2 - 0.01y^2)

\eta = 0.1

, start at

y_0 = 10^{-6}

. How long does descent take to leave?

Solution.

Setup. $a = 1$ , $\gamma = 0.01$ , $\eta\gamma = 0.001$ , so $y_k = (1.001)^ky_0$ . The proposition.
Steps to reach $|y| = 1$ . $k = \dfrac{\ln10^6}{\ln1.001} = \dfrac{13.8155}{0.0009995} \approx 13{,}822$ . Solve $(1.001)^k\cdot10^{-6} = 1$ .
Meaning. The iterate sits within $10^{-6}$ of the saddle for most of those steps, where the gradient is about $10^{-8}$ : the loss curve looks flat, as if training had converged. The gradient norm in $y$ is $\gamma|y_k|$ , tiny until the very end.

Answer. About $13{,}800$ steps; with a tenfold stronger curvature ( $\gamma = 0.1$ ) it would take about $1380$ .

Your answerIn the saddle example (

\eta\gamma = 0.001

y_0 = 10^{-6}

), about how many steps does gradient descent take for

|y|

to grow from

10^{-6}

1

? (Give a decimal.)

Beyond the course: for understanding· Why saddles dominate in high dimension (a heuristic), and what is known

Suppose that at a critical point each of the $n$ Hessian eigenvalues were independently positive or negative with equal chance. Then the chance that all are positive is $2^{-n}$ : for $n = 1000$ , about $10^{-301}$ . So most critical points would be saddles, and the few minima would be extraordinarily rare. Real Hessians are not like this (their eigenvalues are strongly correlated), but the experience with large networks is consistent with it: the main obstacle in practice is long plateaus near saddles, and the local minima that gradient methods find tend to have loss close to the best known. A theorem in the same direction (proved with the stable manifold theorem, so beyond this book): gradient descent with a small fixed step from a random start almost surely avoids saddles that have at least one strictly negative eigenvalue. The minibatch noise of SGD helps further: the random kicks push the iterate off the unstable ridge sooner than the deterministic computation above.

Quick checkAt a critical point of a loss with

10^6

parameters, the Hessian has 999,999 positive eigenvalues and one eigenvalue equal to

-0.001

. What is it, and what does gradient descent do near it?

Where it is used

Used inside maths· Inside maths: fixed points, ODEs and eigenvalue algorithms

Gradient descent is the fixed-point iteration of 4.2 for $T(x) = x - \eta\nabla f(x)$ , and the contraction corollary above is a direct application of the Banach theorem; the same fixed-point idea proves existence for ODEs by Picard’s iteration (7.1). Its continuous-time limit, $x'(t) = -\nabla f(x(t))$ (the gradient flow), is an ordinary differential equation, along which $f$ decreases: $\frac d{dt}f(x(t)) = -\lVert\nabla f\rVert^2$ by the path theorem. Newton’s method is the multivariable form of the iteration whose local contraction was proved in 4.2, and its convergence proof above uses exactly the mean value and invertibility arguments of 5.7. With exact line search on a quadratic, steepest descent is the template that conjugate gradients improves on.

Used in the real world· AI: training a model, in numbers

A language model with $n = 10^9$ parameters is trained by minibatch SGD or Adam on trillions of tokens. Each step computes a gradient of $10^9$ numbers by backpropagation, at a cost comparable to a few forward passes (3.8). The learning rate is typically chosen by trial, then warmed up from near $0$ over the first thousand steps (to avoid the instability when the gradient is large and the step is above $2/L$ of the starting region) and decayed afterwards (the Robbins–Monro conditions: large steps to travel, small steps to average out the noise). A linear regression with $10^6$ features and $10^7$ rows has $X^{\mathsf T}X$ with $10^{12}$ entries; gradient descent never forms it, and costs two products with $X$ per step ( $2\times10^{13}$ operations) instead of the $10^{18}$ of solving the normal equations. If the features are not standardised, $\kappa(X)$ can be $10^3$ , making $\kappa(X^{\mathsf T}X) = 10^6$ and the rate $1 - 2\times10^{-6}$ : half a million steps per $e$ -fold of error. Standardising can bring $\kappa(X)$ down to $10$ , so $\kappa = 100$ and about $350$ steps give a thousandfold reduction of the error.

Traps

Watch out· A small learning rate is only small relative to the curvature

The stability limit is $2/L$ , where $L$ is the largest Hessian eigenvalue. A learning rate of $0.1$ is safe for $L = 1$ and explosive for $L = 100$ . If the loss goes to $\infty$ or $\mathtt{nan}$ in the first steps, cut the rate by a factor of $10$ .

Watch out· Zig-zag is not a bug

Plain gradient descent on a narrow valley must zig-zag (consecutive exact line-search gradients are perpendicular, and at the best fixed step the steep coordinate flips sign every step). Do not conclude the gradient is wrong. Cures: momentum, preconditioning (rescaling the variables to reduce $\kappa$ ), or second-order information.

Watch out· A flat loss curve is not convergence

Near a saddle (or in a very flat valley) the gradient is tiny and the loss barely moves for thousands of steps, then drops suddenly. Check the gradient norm and the Hessian spectrum before stopping.

Watch out· SGD with a constant rate has a noise floor

The loss curve of SGD fluctuates around a level proportional to the learning rate. It does not mean the model has stopped improving: lowering the rate lowers the floor (the schedule of Robbins and Monro).

Watch out· The loss can go up

Gradient descent with a small enough step lowers $f$ at every step (descent lemma). Newton’s method does not (the Rosenbrock example went from $4.7$ to $1412$ ), SGD does not, and momentum does not. The ‘decrease every step’ check belongs only to the first.

Watch out· Convex and non-convex are different problems

For a convex function, any local minimum is global and gradient descent finds it. A neural network’s loss is not convex, and the final loss depends on the starting weights (the lab: seed 2 is stuck at $0.25$ while nine other seeds reach $10^{-29}$ ).

Problems

check

For $f(x, y) = \frac12(2x^2 + 8y^2)$ , what are $L$ , $\mu$ , $\kappa$ , the stable range of $\eta$ , and the best step?

Solution

The Hessian is $\mathrm{diag}(2, 8)$ , so $\mu = 2$ , $L = 8$ , $\kappa = 4$ . Gradient descent converges iff $0 \lt \eta \lt 2/L = 0.25$ . The best step is $\eta^* = \dfrac2{L + \mu} = \dfrac2{10} = 0.2$ , with rate $\dfrac{\kappa - 1}{\kappa + 1} = \dfrac35$ (check: $|1 - 0.2\cdot2| = 0.6$ and $|1 - 0.2\cdot8| = 0.6$ ).

check

Do one step of gradient descent with $\eta = 0.1$ on $f(x, y) = x^2 + 5y^2$ from $(3, 1)$ .

Solution

$\nabla f = (2x, 10y) = (6, 10)$ , so the new point is $(3 - 0.6,\ 1 - 1.0) = (2.4, 0)$ . The $y$ -coordinate was killed in one step because $1 - 0.1\cdot10 = 0$ (the Hessian is $\mathrm{diag}(2, 10)$ and $\eta = 1/L$ ). $f$ fell from $3^2 + 5 = 14$ to $5.76$ .

check

Newton’s method for $f(x) = x^4$ is $x \mapsto x - \dfrac{4x^3}{12x^2} = \frac23x$ . Is the convergence quadratic? Why or why not?

Solution

The iterates $x_k = (\frac23)^kx_0$ converge only linearly (the error is multiplied by $\frac23$ each time). The quadratic convergence theorem needs an invertible Hessian at the minimiser, and $f''(0) = 0$ here: the hypothesis fails.

check

True or false: gradient descent with exact line search on a quadratic in two variables moves along perpendicular directions in consecutive steps. Why?

Solution

True. The exact step minimises $\phi(\eta) = f(x_k - \eta g_k)$ , so $\phi'(\eta_k) = -\nabla f(x_{k+1})\cdot g_k = 0$ : the new gradient is perpendicular to the old one. (This holds for any smooth $f$ , not only quadratics.)

usenew (least squares by gradient descent)

For $f(w) = \frac12\lVert Xw - y\rVert^2$ with $X^{\mathsf T}X = \begin{bmatrix}5 & 3\\ 3 & 5\end{bmatrix}$ , find the best fixed learning rate and the number of steps to reduce the error by a factor $10^{-3}$ .

Solution

The eigenvalues of $X^{\mathsf T}X$ are $2$ and $8$ , so $\mu = 2$ , $L = 8$ , $\kappa = 4$ . The best step is $\eta^* = \frac2{10} = 0.2$ with rate $\rho^* = \frac35 = 0.6$ . We need $0.6^k \le 10^{-3}$ : $k \ge \dfrac{\ln1000}{\ln(1/0.6)} = \dfrac{6.908}{0.5108} = 13.5$ , so $14$ steps.

usenew (momentum numbers)

For $\mu = 1$ and $L = 400$ , give the optimal heavy-ball parameters and the rate, and compare the predicted number of steps with plain descent for a factor $10^{-6}$ .

Solution

$\kappa = 400$ , $\sqrt\kappa = 20$ . $\beta = \big(\frac{19}{21}\big)^2 = 0.8186$ , $\eta = \frac4{(20 + 1)^2} = \frac4{441} = 0.00907$ (with $\sqrt L + \sqrt\mu = 21$ ), rate $\frac{19}{21} = 0.9048$ . Steps: $\frac{\ln10^6}{\ln(21/19)} = \frac{13.82}{0.1001} \approx 138$ . Plain descent: rate $\frac{399}{401}$ , steps $\frac{13.82}{\ln(401/399)} = \frac{13.82}{0.005000} \approx 2763$ . The ratio is about $20 = \sqrt\kappa$ .

usenew (SGD)

An SGD run with $\eta = 0.05$ on a one-parameter least-squares loss has per-example variance $\sigma^2 = 9$ . Find the long-run variance of the iterate, and the batch size that brings it below $0.01$ .

Solution

Floor $= \dfrac{\eta\sigma^2}{2 - \eta} = \dfrac{0.45}{1.95} = 0.2308$ . With batch size $B$ : $\dfrac{\eta\sigma^2/B}{2 - \eta} = \dfrac{0.2308}{B} \lt 0.01$ needs $B \gt 23.1$ , so $B = 24$ .

usenew (Newton in two variables)

Do one Newton step for $f(x, y) = x^2 + xy + 2y^2 - 3x$ from $(0, 0)$ and verify that it lands on the minimiser.

Solution

$\nabla f = (2x + y - 3,\ x + 4y)$ , which at the origin is $(-3, 0)$ . $H = \begin{bmatrix}2 & 1\\ 1 & 4\end{bmatrix}$ (constant, positive definite: trace $6$ , determinant $7$ ). Solve $Hs = (3, 0)$ : from the second row $s_1 = -4s_2$ ; the first gives $-8s_2 + s_2 = 3$ , $s_2 = -\frac37$ , $s_1 = \frac{12}7$ . So the new point is $\big(\frac{12}7, -\frac37\big)$ . Check: $\nabla f = \big(\frac{24}7 - \frac37 - 3,\ \frac{12}7 - \frac{12}7\big) = (0, 0)$ . ✓ A quadratic is solved in one step.

prove

Prove that gradient descent with $\eta \le 1/L$ on an $L$ -smooth function satisfies $\sum_{k = 0}^{K-1}\lVert\nabla f(x_k)\rVert^2 \le \dfrac{2\big(f(x_0) - f_{\min}\big)}{\eta}$ if $f$ is bounded below by $f_{\min}$ , so that $\min_k\lVert\nabla f(x_k)\rVert \to 0$ at the rate $K^{-1/2}$ (no convexity needed).

Solution

Plan. Sum the descent lemma over the steps.

By the descent lemma, $f(x_{k+1}) \le f(x_k) - \frac\eta2\lVert\nabla f(x_k)\rVert^2$ . Summing from $k = 0$ to $K - 1$ telescopes: $f(x_K) \le f(x_0) - \frac\eta2\sum_{k = 0}^{K-1}\lVert\nabla f(x_k)\rVert^2$ . Since $f(x_K) \ge f_{\min}$ , $\sum_{k = 0}^{K-1}\lVert\nabla f(x_k)\rVert^2 \le \frac2\eta\big(f(x_0) - f_{\min}\big)$ . The smallest of the $K$ terms is at most the average, so $\min_k\lVert\nabla f(x_k)\rVert^2 \le \frac{2(f(x_0) - f_{\min})}{\eta K}$ . $\blacksquare$ For a non-convex loss this is the strongest general statement: the gradient becomes small, though the point need not be near a minimum.

provenew (unbiasedness)

Show that the minibatch gradient $g = \frac1B\sum_{j = 1}^B\nabla\ell_{i_j}(w)$ with the $i_j$ drawn independently and uniformly from $\{1, \dots, N\}$ has average $\nabla\ell(w)$ , and in the one-dimensional model of the SGD theorem has variance $\sigma^2/B$ .

Solution

Plan. Linearity of averages for the mean; independence for the variance.

For one index, the average of $\nabla\ell_i(w)$ over the uniform choice is $\frac1N\sum_i\nabla\ell_i(w) = \nabla\ell(w)$ , since the gradient of a sum is the sum of the gradients. A minibatch averages $B$ such terms, so its average is $\nabla\ell(w)$ . In the model $\nabla\ell_i(w) = w - a_i$ , so the minibatch gradient is $w - \bar a_B$ where $\bar a_B$ is the average of $B$ independent samples with variance $\sigma^2$ each; the variance of an average of $B$ independent terms is $\sigma^2/B$ . $\blacksquare$

provenew (the stability limit is sharp)

Show that for $f = \frac12\lambda x^2$ with $\lambda \gt 0$ in one variable, gradient descent converges for every start iff $0 \lt \eta \lt 2/\lambda$ , and for $\eta = 2/\lambda$ it oscillates forever without converging (unless started at $0$ ).

Solution

Plan. The update is multiplication by a constant.

$x_{k+1} = x_k - \eta\lambda x_k = (1 - \eta\lambda)x_k$ , so $x_k = (1 - \eta\lambda)^kx_0$ . This tends to $0$ for all $x_0$ iff $|1 - \eta\lambda| \lt 1$ , i.e. $0 \lt \eta\lambda \lt 2$ . At $\eta\lambda = 2$ the multiplier is $-1$ and $x_k = (-1)^kx_0$ does not converge for $x_0 \ne 0$ . For $\eta\lambda \gt 2$ the multiplier is below $-1$ and $|x_k|$ grows. $\blacksquare$

provenew (Newton is affine invariant)

Let $A$ be an invertible $n \times n$ matrix and $g(z) = f(Az)$ . Show that the Newton iterates for $g$ started at $z_0$ are $z_k = A^{-1}x_k$ , where $x_k$ are the Newton iterates for $f$ started at $x_0 = Az_0$ , whereas gradient descent does not have this property.

Solution

Plan. Use the chain rule (5.4) for the gradient and the Hessian of $g$ .

$\nabla g(z) = A^{\mathsf T}\nabla f(Az)$ and $H_g(z) = A^{\mathsf T}H_f(Az)A$ . The Newton step for $g$ at $z$ is $-H_g^{-1}\nabla g = -A^{-1}H_f^{-1}A^{-\mathsf T}A^{\mathsf T}\nabla f(Az) = -A^{-1}H_f(Az)^{-1}\nabla f(Az)$ , so $Az_{k+1} = Az_k - H_f(Az_k)^{-1}\nabla f(Az_k)$ , which is the Newton step for $f$ at $x_k = Az_k$ . By induction $x_k = Az_k$ . The gradient descent step for $g$ is $-\eta A^{\mathsf T}\nabla f(Az)$ , whose image under $A$ is $-\eta AA^{\mathsf T}\nabla f$ , not $-\eta\nabla f$ unless $AA^{\mathsf T} = I$ : changing variables changes the conditioning, and gradient descent notices. $\blacksquare$ (Newton’s method does not care about the condition number; this is why preconditioning is an approximation to it.)

stretchnew (momentum on a one-dimensional quadratic)

For $f(x) = \frac12\lambda x^2$ the heavy-ball recursion is $x_{k+1} = (1 + \beta - \eta\lambda)x_k - \beta x_{k-1}$ . Take $\beta = 0.5$ , $\eta\lambda = 0.5$ and compute $x_2, x_3, x_4$ from $x_0 = x_1 = 1$ ; then show that the rate is $\sqrt\beta$ when $|1 + \beta - \eta\lambda| \le 2\sqrt\beta$ .

Solution

Here $1 + \beta - \eta\lambda = 1$ , so $x_{k+1} = x_k - 0.5x_{k-1}$ . $x_2 = 1 - 0.5 = 0.5$ , $x_3 = 0.5 - 0.5 = 0$ , $x_4 = 0 - 0.25 = -0.25$ . The condition $|1| \le 2\sqrt{0.5} = 1.414$ holds, so the roots of $z^2 - z + 0.5 = 0$ are $\frac{1 \pm i}2$ , of absolute value $\sqrt{0.5} = 0.707 = \sqrt\beta$ : the sequence rotates and decays, which is the oscillation visible in the numbers (positive, zero, negative). In general for $|\tau| \le 2\sqrt\beta$ the roots of $z^2 - \tau z + \beta$ are complex conjugates (or equal) with product $\beta$ , hence each has absolute value $\sqrt\beta$ ; the error is a combination of $z^k$ and $\bar z^k$ , so it decays like $(\sqrt\beta)^k$ (times $1 + k$ in the double root case).

stretchnew (the contraction view)

For $f(x) = \frac12x^2 + \ln(1 + e^{x})$ in one variable, show that gradient descent converges for $\eta = 1$ using the contraction corollary, and find the contraction constant.

Solution

$f'(x) = x + \dfrac{e^x}{1 + e^x}$ and $f''(x) = 1 + \dfrac{e^x}{(1 + e^x)^2} \in \big[1, \frac54\big]$ (since $\frac{e^x}{(1 + e^x)^2} \le \frac14$ , with the maximum at $x = 0$ ). So $\mu = 1$ and $L = \frac54$ , and $\eta = 1 \lt 2/L = 1.6$ . The map $T(x) = x - f'(x)$ has $T'(x) = 1 - f''(x) = -\dfrac{e^x}{(1 + e^x)^2} \in \big[-\frac14, 0\big]$ , so $|T'| \le \frac14$ : $T$ is a contraction with constant $\rho = \frac14 = \max(|1 - \mu|, |1 - L|) = \max(0, \frac14)$ . By Banach, gradient descent converges to the unique solution of $f'(x) = 0$ (here $x + \frac{e^x}{1 + e^x} = 0$ , near $-0.4$ ), with error at most $\frac{(1/4)^k}{3/4}|x_1 - x_0|$ .

Summary

Gradient descent $x_{k+1} = x_k - \eta\nabla f(x_k)$ follows the steepest descent direction of 5.2. For an $L$ -smooth function a step $\eta \le 1/L$ lowers $f$ by at least $\frac\eta2\lVert\nabla f\rVert^2$ (descent lemma).
On a quadratic the error is multiplied by $I - \eta A$ : it converges iff $0 \lt \eta \lt 2/L$ , with rate $\max(|1 - \eta\mu|, |1 - \eta L|)$ , best at $\eta = \frac2{L + \mu}$ where the rate is $\frac{\kappa - 1}{\kappa + 1}$ . The condition number $\kappa$ of 3.5 sets the number of steps (about $\frac\kappa2\ln\frac1\varepsilon$ ), and the path zig-zags because the steep coordinate flips sign at every step. For least squares $\kappa = \kappa(X)^2$ .
For any $C^2$ function with Hessian eigenvalues in $[\mu, L]$ , $x \mapsto x - \eta\nabla f(x)$ is a contraction (Banach, 4.2). For convex smooth $f$ , $f(x_k) - f^* \le \frac{L\lVert x_0 - x^*\rVert^2}{2k}$ .
Momentum $x_{k+1} = x_k - \eta\nabla f(x_k) + \beta(x_k - x_{k-1})$ with the tuned $\eta, \beta$ improves the rate to $\frac{\sqrt\kappa - 1}{\sqrt\kappa + 1}$ : $\sqrt\kappa$ steps instead of $\kappa$ .
Stochastic gradient descent uses one example (or a minibatch) per step; the step is unbiased, the iterate has a noise floor proportional to $\eta\sigma^2/B$ , and steps must decrease ( $\sum\eta_k = \infty$ , $\sum\eta_k^2 \lt \infty$ ; Robbins and Monro, 1951) to converge exactly.
Newton’s method $x - H^{-1}\nabla f$ minimises the quadratic model, solves quadratics in one step, is invariant to conditioning, and converges quadratically near a minimiser with invertible Hessian; but it costs $n^2$ memory and $n^3$ time, is not monotone, and is attracted to saddles and maxima.
Training a network is gradient descent on a non-convex loss with the gradient from backpropagation (3.8): results depend on the start (one seed in ten failed in the lab).
In high dimension, saddle points are far more numerous than minima; near a saddle with weak negative curvature $\gamma$ , descent leaves after about $\frac{\ln(1/|y_0|)}{\eta\gamma}$ steps, so a flat loss curve is not proof of convergence.

Watch and read

The same ideas explained by someone else. Use these when an explanation here does not land.

5.10 Gradient descent, momentum and Newton's method

The picture: a ball rolling down a valley

The idea, in words

Definitions

Gradient descent

Learning rate

Loss function

One step lowers the function: the descent lemma

Descent lemma

Gradient descent on a quadratic: the eigenvalues decide

Gradient descent on a quadratic

Gradient descent is a contraction

Exact line search, and why the zig-zag is built in

Convex functions: progress without a second derivative

Gradient descent on a convex smooth function

Momentum

Momentum

Momentum on a quadratic

Stochastic gradient descent

Stochastic gradient descent

Noise floor of stochastic gradient descent

Newton’s method in many variables

Newton's method

Quadratic convergence of Newton's method

Training a small neural network

Saddle points in high dimension

Gradient descent near a saddle

Where it is used

Traps

Problems

Summary

Watch and read

Where this sits

Comes from

Leads to

Same idea elsewhere