Lab U5: Gradient descent and tiny training

Unit: Unit 5, Optimization and training
Role: Required
Textbook sections: directional derivatives and gradient descent; eigenvalues, eigenvectors, and diagonalization; symmetric matrices and quadratic forms; Hessians and the second derivative test; least squares from gradients; fixed-hidden-layer training; applications and computation recap
Core path: gradient descent updates, learning-rate diagnosis, debugging the descent sign, two-dimensional quadratic losses, Hessian eigenvalues, least squares as gradient zero, residual orthogonality, fixed-hidden-layer design matrices, and rank-one updates
Optional extension idea: try different learning rates or replace the fixed feature columns in Part 8

Computational tools used in this lab

Before starting, review these parts of Appendix B, NumPy and SymPy Quick Reference for the Labs:

• Appendix B.2: NumPy arrays, vectors, matrices, and shapes
• Appendix B.4: Elementwise arithmetic versus linear algebra
• Appendix B.8: Constructing common arrays and design matrices
• Appendix B.9: Numerical checks and roundoff
• Appendix B.10: NumPy linear algebra commands
• Appendix B.12: Small Python patterns used in the labs

The goal is to interpret the mathematical computation, not to memorize every command.

Part 0. Setup and reading habits

Math reminder. Unit 5 computations usually involve one of these patterns:

• a gradient descent update,
• an eigenvalue computation for a Hessian,
• a least-squares solve,
• a residual orthogonality check,
• a design matrix for fixed features,
• an outer product update.

Predict before running. Which objects below are vectors? Which objects are matrices? What shape should H_demo @ u have? What shape should np.outer(g, h) have?

import numpy as np

np.set_printoptions(precision=3, suppress=True)

x = np.array([4.0])
u = np.array([3.0, 4.0])
H_demo = np.array([
    [4.0, 0.0],
    [0.0, 1.0],
])

g = np.array([1.0, -2.0, 3.0])
h = np.array([4.0, 0.0])

x.shape, u.shape, H_demo.shape, g.shape, h.shape, H_demo @ u, np.outer(g, h).shape

Run and compare. The objects x, u, g, and h are vectors. The object H_demo is a matrix. The product H_demo @ u is a matrix-vector product, so it has the same number of entries as the number of rows of H_demo. The command np.outer(g, h) creates a matrix with one row for each entry of g and one column for each entry of h.

Interpretation check. Before interpreting any code in this lab, first identify the mathematical objects and their shapes.

Part 1. One-dimensional gradient descent

Math reminder. For (f(x)=x^2), the derivative is (f'(x)=2x). Gradient descent with learning rate (\alpha) uses [ x_{k+1}=x_k-\alpha f'(x_k). ]

In the code, we store (x) as a one-entry NumPy array so that the same update syntax also works for vector inputs later.

Predict before running. Starting from (x_0=4), should the first few loss values (f(x_k)) decrease when (\alpha=0.2)?

def f1(x):
    return float(x @ x)

def grad_f1(x):
    return 2 * x

x = np.array([4.0])
alpha = 0.2
losses = []
iterates = []

for k in range(8):
    iterates.append(float(x[0]))
    losses.append(f1(x))
    x = x - alpha * grad_f1(x)

iterates, losses

Run and compare. The list iterates stores the successive values of (x_k). The list losses stores (f(x_k)).

Interpretation check. The line x = x - alpha * grad_f1(x) is the gradient descent update. The minus sign matters: the derivative points toward increasing (f), so the negative derivative points toward decreasing (f), at least locally.

Answer sketch. For this simple example, the loss decreases quickly for (\alpha=0.2). This is evidence of useful progress for this starting point and learning rate, but it is not a proof that every gradient descent run converges.

Part 2. Learning-rate comparison

Math reminder. The learning rate controls the step size. Too small can be slow. Too large can be unstable.

Predict before running. For the same function (f(x)=x^2), which learning rate should be slow, which should be useful, and which should be unstable: (0.05), (0.20), or (1.05)?

def run_gd(alpha, steps=4):
    x = np.array([4.0])
    losses = []
    for k in range(steps):
        losses.append(round(f1(x), 3))
        x = x - alpha * grad_f1(x)
    return losses

{alpha: run_gd(alpha) for alpha in [0.05, 0.20, 1.05]}

Run and compare. The learning rate (\alpha=0.05) makes slow but steady progress. The learning rate (\alpha=0.20) makes faster useful progress. The learning rate (\alpha=1.05) is unstable in this example.

Interpretation check. A loss table is a diagnostic. It tells what happened for the displayed iterates; it does not by itself prove that a global minimum has been found.

Common mistake. A larger learning rate is not always better.

Part 3. Debugging the sign

Math reminder. Gradient descent uses a minus sign: [ x_{k+1}=x_k-\alpha\nabla f(x_k). ] Changing the minus sign to a plus sign gives a gradient ascent step.

Predict before running. Which loss list should decrease: the one labeled "minus sign" or the one labeled "plus sign"?

def run_gd_with_sign(alpha, sign=-1, steps=5):
    x = np.array([4.0])
    losses = []
    for k in range(steps):
        losses.append(round(f1(x), 3))
        x = x + sign * alpha * grad_f1(x)
    return losses

{
    "minus sign": run_gd_with_sign(0.2, sign=-1),
    "plus sign": run_gd_with_sign(0.2, sign=1),
}

Run and compare. The "minus sign" version decreases the loss. The "plus sign" version increases the loss.

Interpretation check. The expression grad_f1(x) points in the direction of steepest increase for this one-dimensional function. Subtracting it moves in the opposite direction.

Debugging habit. When a gradient descent loop makes the loss increase immediately, first check the sign and the learning rate.

Part 4. Two-dimensional quadratic loss

Math reminder. A quadratic loss often has the form [ f(\mathbf{x})=\frac12\mathbf{x}^T H\mathbf{x}+\mathbf{b}^T\mathbf{x}, ] where (H) is a symmetric matrix. Its gradient is [ \nabla f(\mathbf{x})=H\mathbf{x}+\mathbf{b}. ]

Predict before running. What shape should grad_quad(x) have? Why must it have the same shape as x?

H = np.array([
    [4.0, 1.0],
    [1.0, 2.0],
])
b = np.array([-1.0, 2.0])

def f_quad(x):
    return float(0.5 * (x @ H @ x) + b @ x)

def grad_quad(x):
    return H @ x + b

x = np.array([2.0, -1.0])
alpha = 0.1
g = grad_quad(x)
x_next = x - alpha * g

f_quad(x), g, x_next, f_quad(x_next)

Run and compare. The gradient g has the same shape as x. This must happen because the update subtracts alpha * g from x.

Interpretation check. The code computes one gradient descent step for a two-variable loss. The value f_quad(x_next) can be compared with f_quad(x) to check whether this step lowered the loss.

Common mistake. The expression H @ x + b is not a Hessian. It is the gradient for this quadratic loss. The Hessian is the matrix H.

Part 5. Hessian eigenvalues and curvature

Math reminder. At a critical point, the Hessian controls the quadratic term: [ f(\mathbf{a}+\mathbf{h})-f(\mathbf{a}) \approx \frac12\mathbf{h}^T H_f(\mathbf{a})\mathbf{h}. ] For a symmetric Hessian, eigenvalue signs diagnose curvature directions.

The command np.linalg.eigvalsh is specialized for symmetric matrices. For a general square matrix, np.linalg.eig computes eigenvalues and eigenvectors.

Predict before running. Which matrix below should correspond to a local minimum test, and which should correspond to a saddle point test?

H_min = np.array([
    [4.0, 1.0],
    [1.0, 2.0],
])

H_saddle = np.array([
    [4.0, 0.0],
    [0.0, -1.0],
])

np.linalg.eigvalsh(H_min), np.linalg.eigvalsh(H_saddle)

Run and compare. The matrix H_min has positive eigenvalues, so it is positive definite. The matrix H_saddle has one positive and one negative eigenvalue, so it is indefinite.

Interpretation check. A positive definite Hessian at a critical point gives a local minimum. An indefinite Hessian at a critical point gives a saddle point.

Common mistake. The Hessian test classifies a critical point. Away from a critical point, the first-order gradient term still matters.

Part 6. Curvature and stable step sizes

Math reminder. For the quadratic loss [ f(\mathbf{x})=\frac12\mathbf{x}^T H\mathbf{x}, ] the gradient is (H\mathbf{x}). Large Hessian eigenvalues mean steep curvature directions. A learning rate that is too large can overshoot in those directions.

Predict before running. Which learning rate should be slow, which should be useful, and which should be unstable: (0.05), (0.40), or (0.60)?

H = np.diag([4.0, 1.0])

def f_curv(x):
    return float(0.5 * (x @ H @ x))

def grad_curv(x):
    return H @ x

def run_quadratic_gd(alpha, steps=6):
    x = np.array([4.0, 4.0])
    losses = []
    for k in range(steps):
        losses.append(round(f_curv(x), 3))
        x = x - alpha * grad_curv(x)
    return losses

{alpha: run_quadratic_gd(alpha) for alpha in [0.05, 0.40, 0.60]}

Run and compare. The learning rate (\alpha=0.05) is slow. The learning rate (\alpha=0.40) makes useful progress. The learning rate (\alpha=0.60) is unstable for this quadratic loss.

Interpretation check. The largest eigenvalue of the Hessian corresponds to the steepest curvature direction. A step size can be too large because of one steep direction even when another direction is less steep.

Part 7. Least squares as gradient zero

Math reminder. Unit 4 derived least squares using residual orthogonality. Unit 5 derives the same normal equations by minimizing [ h(\mathbf{x})=|A\mathbf{x}-\mathbf{b}|^2. ] The gradient is [ \nabla h(\mathbf{x})=2A^T(A\mathbf{x}-\mathbf{b}). ]

Predict before running. After least squares, should A.T @ r be close to zero? Does that mean the residual vector itself must be zero?

A = np.array([
    [1.0, 0.0],
    [1.0, 1.0],
    [1.0, 2.0],
])
b = np.array([1.0, 2.0, 2.0])

xhat = np.linalg.lstsq(A, b, rcond=None)[0]
r = b - A @ xhat
grad_at_xhat = 2 * (A.T @ (A @ xhat - b))
H_lstsq = 2 * (A.T @ A)

xhat, r, A.T @ r, grad_at_xhat, np.linalg.eigvalsh(H_lstsq)

Run and compare. The vector xhat is the least-squares coefficient vector. The vector r is the residual (\mathbf{b}-A\hat{\mathbf{x}}). The expression A.T @ r is close to the zero vector.

Interpretation check. Unit 4 reads A.T @ r as residual orthogonality. Unit 5 reads grad_at_xhat as the zero-gradient condition. These are the same normal equations written from two viewpoints.

Common mistake. The residual can be nonzero even when A.T @ r is zero. Least squares asks for the closest reachable vector, not necessarily an exact solution.

Part 8. Fixed-hidden-layer least squares

Math reminder. Let (\sigma(t)=\tanh(t)). A fixed-hidden-layer model can have the form [ N_{\mathbf{c}}(t)=c_0+c_1\sigma(t)+c_2\sigma(t-1)+c_3\sigma(t+1). ] The features are nonlinear functions of (t), but the model is linear in the trainable coefficient vector (\mathbf{c}).

Predict before running. Which objects below are fixed features? Which object is trained?

t = np.array([-2.0, -1.0, 0.0, 1.0, 2.0])
y = np.array([4.0, 1.0, 0.0, 1.0, 4.0])

A = np.column_stack([
    np.ones_like(t),
    np.tanh(t),
    np.tanh(t - 1),
    np.tanh(t + 1),
])

c = np.linalg.lstsq(A, y, rcond=None)[0]
yhat = A @ c
residual = y - yhat

A.shape, c, yhat, np.linalg.norm(residual)

Run and compare. The columns of A are fixed feature vectors. The vector c is trained by least squares. The predictions are yhat = A @ c.

Interpretation check. The model is nonlinear as a function of the input (t), because it uses (\tanh(t)), (\tanh(t-1)), and (\tanh(t+1)). It is linear as a function of (\mathbf{c}), so least squares applies.

Common mistake. “Nonlinear model” and “least-squares training” are not contradictory here. The model is nonlinear in (t), but linear in the coefficients being trained.

Part 9. Rank-one update shape

Math reminder. If (\mathbf{g}\in\mathbb{R}^m) and (\mathbf{h}\in\mathbb{R}^d), then [ \mathbf{g}\mathbf{h}^T ] is an (m\times d) matrix. In NumPy, np.outer(g, h) explicitly forms this outer product when g and h are stored as one-dimensional arrays.

Predict before running. What shape should G have? Why should its rank be at most one?

g = np.array([1.0, -2.0, 3.0])
h = np.array([4.0, 0.0])

G = np.outer(g, h)
W = np.ones((3, 2))
alpha = 0.1
W_new = W - alpha * G

G, G.shape, np.linalg.matrix_rank(G), W_new

Run and compare. The matrix G has shape (3, 2). It is the outer product (\mathbf{g}\mathbf{h}^T). Since every column of this matrix is a scalar multiple of (\mathbf{g}), its rank is at most one.

Interpretation check. The line W_new = W - alpha * G subtracts a scaled rank-one update from W.

Common mistake. For one-dimensional NumPy arrays, g @ h.T does not create the outer product. Use np.outer(g, h) when the goal is (\mathbf{g}\mathbf{h}^T).

Part 10. Review code-reading bank

These are short checks for discussion or self-review. They are meant to be answered without writing long programs.

1. Update rule. What mathematical update is represented by x = x - alpha * grad_f(x)?
2. Learning rate. Which symbol in (\mathbf{x}_{k+1}=\mathbf{x}_k-\alpha\nabla f(\mathbf{x}_k)) corresponds to alpha?
3. Shape check. Why must grad_f(x) have the same shape as x?
4. Sign check. What usually happens if the minus sign is changed to a plus sign?
5. Loss history. What does losses.append(f(x)) record?
6. Learning-rate diagnosis. What does an increasing loss table suggest?
7. Hessian object. In a quadratic loss, what object is usually represented by H?
8. Eigenvalue signs. What do all positive Hessian eigenvalues imply at a critical point?
9. Eigenvalue signs. What does it imply if a Hessian has both positive and negative eigenvalues at a critical point?
10. Zero eigenvalue. Why does a zero Hessian eigenvalue make the second derivative test inconclusive?
11. Least squares. What does np.linalg.lstsq(A, b, rcond=None)[0] return?
12. Residual. In r = b - A @ xhat, what is r?
13. Residual orthogonality. What condition is checked by A.T @ r?
14. Residual warning. Does A.T @ r being close to zero mean that r is close to zero?
15. Gradient viewpoint. How does Unit 5 get the normal equations?
16. Fixed features. In the fixed-hidden-layer example, what are the columns of A?
17. Trained vector. Which vector is trained in yhat = A @ c?
18. Linearity check. Is the fixed-hidden-layer model linear in (t), in (\mathbf{c}), both, or neither?
19. Outer product. What shape does np.outer(g, h) have if g.shape == (m,) and h.shape == (d,)?
20. Rank-one update. Why does (\mathbf{g}\mathbf{h}^T) have rank at most one?

Answer sketches

1. (\mathbf{x}_{k+1}=\mathbf{x}_k-\alpha\nabla f(\mathbf{x}_k)).
2. The learning rate (\alpha).
3. The update subtracts one vector from another.
4. It moves toward steepest increase instead of steepest decrease.
5. The loss value at the current iterate.
6. The learning rate may be too large, or the update may have the wrong sign.
7. A Hessian matrix.
8. The Hessian is positive definite, so the critical point is a local minimum.
9. The Hessian is indefinite, so the critical point is a saddle point.
10. The quadratic approximation has a flat direction, so second-order information alone does not decide.
11. A least-squares coefficient vector.
12. The residual vector (\mathbf{b}-A\hat{\mathbf{x}}).
13. Whether the residual is orthogonal to the columns of (A).
14. No. The residual can be nonzero and still orthogonal to the columns of (A).
15. By setting (\nabla|A\mathbf{x}-\mathbf{b}|^2=0).
16. The fixed feature vectors (1), (\tanh(t)), (\tanh(t-1)), and (\tanh(t+1)).
17. The coefficient vector (\mathbf{c}).
18. It is nonlinear in (t), but linear in (\mathbf{c}).
19. It has shape (m, d).
20. Every column is a scalar multiple of (\mathbf{g}), so the column space has dimension at most one.