Optimization for Data Science · Lecture 7
Learning curvature with L-BFGS
Gradient descent uses a scalar step size. Accelerated gradient descent adds momentum. L-BFGS uses recent steps and gradient differences to learn a geometry. Compare all three on the same problem, from the same starting point.
1A narrow valley
Watch L-BFGS learn the curvature of a rotated quadratic. GD uses its optimal constant step size for this quadratic, and AGD uses the exact L and μ. L-BFGS uses only function values and gradients.
- Press Run at slow speed to see L-BFGS change direction after collecting curvature pairs. Click the contour plot to choose a new start.
- Press Show result, then switch the convergence plot from iterations to gradient evaluations. L-BFGS pays for every line-search trial.
- Set κ = 1: GD and AGD solve the quadratic in one step. Increase κ to see where learning curvature helps.
2A little memory, many parameters
The same condition number can hide very different spectra. Start with clustered eigenvalues, then spread them out and change the number of stored pairs. The starting point has equal objective energy in every eigendirection.
- Compare Five clusters and Spread out at the same dimension and condition number.
- Change the memory from m = 1 to m = 10 or 20. Reusing more curvature information can reduce the work needed.
- Increase the dimension while keeping m fixed. The curvature-pair storage grows as 2md numbers; a dense inverse-Hessian matrix needs d².
3Logistic regression on real data
Fit a classifier to the Wisconsin Diagnostic Breast Cancer dataset: 569 samples, 30 measurements of cell nuclei, and a benign/malignant label. Unlike the quadratics above, the curvature changes as the classifier learns. All three methods minimize the same regularized training loss from zero.
f(w) = (1/n) ∑ᵢ [log(1 + exᵢᵀw) − yᵢxᵢᵀw] + (λ/2) ‖w‖².
- Run a few iterations and inspect the held-out predictions. A point to the wrong side of ½ is misclassified.
- Press Show result and compare both iterations and full-training-set gradient evaluations.
- Reduce λ. The guaranteed strong-convexity constant becomes smaller; compare how the methods respond.
Data, preprocessing, and reference solution
A fixed, stratified 80/20 split separates training and held-out samples. Each feature is centered and scaled using the training set only. We append an intercept and penalize all 31 coefficients, including the intercept, so μ = λ is a valid global bound. The held-out set is used only for the displayed predictions and accuracy.
GD and AGD use a global smoothness bound L ≥ ‖X‖²/(4n) + λ. These global bounds can be conservative along the trajectory. The displayed L/μ is a guaranteed curvature ratio, not an exact condition number of every Hessian. L-BFGS uses the same two-loop recursion and Wolfe search as above.
The objective plot shows (f(wₖ) − fref)/(f(w₀) − fref). A separate high-accuracy Newton solve computes fref; its error is bounded by ‖∇f(wref)‖²/(2λ). This reference is used only to measure progress and stopping, never to choose an optimization step. Its setup cost and prediction diagnostics are excluded from the method counts.
Open the experiment to compute the reference solution.
Classification accuracy may stabilize well before the objective gap becomes tiny. Optimizing the training objective more accurately does not necessarily improve held-out accuracy.
Data: Wolberg, Mangasarian, Street, and Street (1993), Breast Cancer Wisconsin (Diagnostic), UCI Machine Learning Repository, licensed CC BY 4.0. The original numeric measurements are embedded in this file; sample identifiers are omitted, and labels are encoded as benign = 0, malignant = 1. Standardization and splitting are performed locally.
4Thousands of words, a little memory
Train a spam classifier on real SMS messages. Each word in the training vocabulary has its own coefficient, plus an intercept. The data are sparse: a message contains only a small fraction of the vocabulary. L-BFGS learns curvature while retaining only a few pairs of parameter and gradient differences.
Open this experiment to compute the comparison.
What limited memory saves
Curvature storage at 8 bytes per number (Float64).
| Representation | Storage |
|---|---|
| Full BFGS d × d dense matrix | — |
| L-BFGS m pairs of length d | — |
Full BFGS storage is an estimate; that matrix is never allocated here. Both methods also need working vectors. These figures exclude the dataset and plotting diagnostics.
Direction computation: O(d²) for a dense inverse-BFGS matrix versus O(md) for the L-BFGS two-loop recursion. Sparse data also make gradient evaluations inexpensive.
Gradient convergence
‖∇f(wₖ)‖₂ / ‖∇f(w₀)‖₂, logarithmic vertical scale.
Compute the comparison to see the work needed.
- Compare the storage figures with the 31-parameter example above. Here the dimension comes from the vocabulary.
- Change m from 1 to 10 or 30. More history costs more memory and may reduce the optimization work.
- Compare iterations and gradient evaluations, then inspect the final spam recall and F1 score. Faster optimization and better predictive performance are different questions.
Data and numerical details
We use the UCI SMS Spam Collection by Almeida and Hidalgo (2011), licensed CC BY 4.0. Repeated messages are removed before a fixed stratified 80/20 split. Vocabulary and inverse-document frequencies are learned from training messages only. Each message is represented by nonnegative TF-IDF word features, normalized to unit Euclidean norm when nonempty. The intercept is stored separately. Unseen held-out words are ignored.
All methods fit mean logistic loss plus λ‖w‖²/2, penalizing the intercept too. We use μ = λ and a valid global L ≥ λ + ‖X‖²/(4n), computed with sparse matrix-vector products. GD uses 2/(L + μ); AGD uses step 1/L and momentum (√L − √μ)/(√L + √μ). L-BFGS uses the same two-loop recursion and strong Wolfe conditions as the preceding experiments.
This experiment measures the relative gradient norm, so it needs no reference optimizer or dense Hessian. Computation runs in a background worker. Only scalar iteration traces and final coefficient vectors are retained; we do not store a d-dimensional parameter vector for every animation frame. The reported optimizer storage excludes these visualization traces.
Algorithms and comparison details
Quadratic experiments: f(w) = ½ wᵀAw, A ≻ 0, w⋆ = 0, f⋆ = 0, L = 1, μ = 1/κ.
The first three experiments start every method at the same point and stop on the relative objective gap (f⋆ = 0 for the quadratics) or the iteration budget. The SMS experiment instead uses a common relative-gradient-norm target. The plots use the actual iterates, with no smoothing or fitted convergence rates. Values below the bottom of the logarithmic plot (including zero) are placed at its lower edge; the numerical readouts retain their actual values. All convergence-plot horizontal axes are linear; the objective-gap vertical axis is logarithmic.
Gradient descent
wk+1 = wk − η∇f(wk)
The default η = 2/(L + μ) minimizes the worst contraction factor of a constant-step GD iteration on this quadratic. You can also select η = 1/L.
Accelerated gradient descent
wk+1 = zk − ∇f(zk)/L
zk+1 = wk+1 + β(wk+1 − wk)
z0 = w0 and β = (√κ − 1)/(√κ + 1). This is Nesterov's strongly convex variant, using exact curvature bounds on the quadratics and valid global bounds on logistic regression. Its objective values need not decrease at every step.
Limited-memory BFGS
pk = −Hk∇f(wk), wk+1 = wk + ηkpk
The two-loop recursion uses the last m accepted (s, y) pairs. Hk(0) = γkI, with γ = sᵀy/(yᵀy) from the newest pair; γ0 = 1. Strong Wolfe search uses c1 = 10⁻⁴ and c2 = 0.9, starting with η = 1.
Counting work. Gradient evaluations count calls used by the optimization algorithm, including all rejected and accepted L-BFGS line-search trials. Function evaluations include the objective values used to display convergence. Extra gradient norms computed only for diagnostics are excluded from the gradient-work axis. In the SMS experiment, accepted-point gradients are required for stopping and are included in the counts. Neither axis is a wall-clock benchmark: L-BFGS also does O(md) arithmetic per direction. The displayed storage comparison counts only curvature pairs versus a dense matrix, at 8 bytes per number; it excludes the problem, histories, and working vectors.
The problems. In two dimensions, A is a rotated diagonal matrix. The higher-dimensional example uses a fixed orthogonal Householder rotation of a diagonal matrix, evaluated without forming a dense matrix. Its five eigenvalue clusters have a small internal spread; the alternative is a logarithmically spaced spectrum. The construction is deterministic and reproducible.
These experiments show why L-BFGS can be very effective on deterministic, ill-conditioned problems. They do not assert that it always beats AGD, or that a fixed memory gives a general superlinear convergence guarantee.
References: Liu and Nocedal (1989); Bottou, Curtis, and Nocedal, §6.2.