Limited-memory BFGS (L-BFGS)

The BFGS algorithm, an optimization method, improves upon approaches like Newton's method by avoiding direct Hessian computation and inversion. Despite this advancement, BFGS still requires storing and updating an approximation of the Hessian (or its inverse), typically as a dense $d \times d$ matrix. For models with millions or even billions of parameters, storing such a matrix is completely infeasible due to memory constraints ($O(d^2)$ memory). This is where the Limited-memory BFGS (L-BFGS) algorithm becomes essential.

The Core Idea: Trading Memory for Computation

L-BFGS cleverly sidesteps the $O(d^2)$ memory requirement of standard BFGS. Instead of storing the dense $d \times d$ Hessian approximation matrix $H_k$, L-BFGS only stores a small, fixed number, say $m$, of the most recent correction pairs $(s_i, y_i)$. Recall that:

• $s_k = x_{k+1} - x_k$ (change in parameters)
• $y_k = \nabla f(x_{k+1}) - \nabla f(x_k)$ (change in gradients)

Typically, $m$ is a small integer, often between 5 and 20, regardless of the problem dimension $d$. This reduces the memory requirement drastically to $O(md)$, which is linear in the number of parameters $d$ and thus manageable even for very large models.

The trade-off is that instead of having direct access to the matrix $H_k$, L-BFGS computes the matrix-vector product $H_k \nabla f_k$ (which is needed to determine the search direction $p_k = -H_k \nabla f_k$) implicitly using only the stored $(s_i, y_i)$ pairs. This computation involves a procedure often called the "L-BFGS two-loop recursion."

The L-BFGS Two-Loop Recursion

The heart of L-BFGS is how it reconstructs the search direction $p_k = -H_k \nabla f_k$ without ever forming the matrix $H_k$. It uses the BFGS update formula, which can be expressed recursively. The algorithm computes the product $H_k \nabla f_k$ using the $m$ most recent pairs $\{ (s_i, y_i) \}_{i=k-m}^{k-1}$.

The process involves:

1. Initialization: Start with an initial Hessian approximation $H_k^0$. A common choice is $H_k^0 = \gamma_k I$, where $I$ is the identity matrix and $\gamma_k$ is a scaling factor chosen to approximate the magnitude of the true Hessian. A popular choice is: $$\gamma_k = \frac{s_{k-1}^T y_{k-1}}{y_{k-1}^T y_{k-1}}$$ This scaling helps ensure the search direction has an appropriate magnitude.
2. First Loop (Backward Pass): Iterate backward from $i = k-1$ down to $k-m$. In this loop, the algorithm computes intermediate values $\alpha_i$ and updates the gradient vector $q = \nabla f_k$.
   • $\rho_i = 1 / (y_i^T s_i)$
   • $\alpha_i = \rho_i s_i^T q$
   • $q = q - \alpha_i y_i$
3. Scaling: Multiply the intermediate vector $q$ by the initial Hessian approximation $H_k^0$: $z = H_k^0 q$.
4. Second Loop (Forward Pass): Iterate forward from $i = k-m$ up to $k-1$. This loop uses the stored $\alpha_i$ values and the $(s_i, y_i)$ pairs to reconstruct the final direction.
   • $\beta_i = \rho_i y_i^T z$
   • $z = z + s_i (\alpha_i - \beta_i)$
5. Result: The final vector $z$ is the desired product, $z = H_k \nabla f_k$. The search direction is then $p_k = -z$.

This two-loop procedure requires only vector operations (dot products, scalar-vector multiplications, vector additions), resulting in a computational cost of approximately $O(md)$ per iteration, significantly cheaper than the $O(d^2)$ cost of standard BFGS updates or the $O(d^3)$ cost of Newton's method.

L-BFGS Algorithm Outline

A typical iteration of the L-BFGS algorithm looks like this:

1. Choose starting point $x_0$.
2. Set memory size $m$.
3. For $k = 0, 1, 2, ...$ until convergence: a. Compute the gradient $\nabla f_k = \nabla f(x_k)$. b. Compute the search direction $p_k = -H_k \nabla f_k$ using the two-loop recursion with the latest $m$ pairs $(s_i, y_i)$ and an initial Hessian approximation $H_k^0$. c. Find a suitable step length $\alpha_k > 0$ using a line search procedure (e.g., backtracking line search satisfying the Wolfe conditions). d. Update the parameters: $x_{k+1} = x_k + \alpha_k p_k$. e. Compute the new gradient $\nabla f_{k+1} = \nabla f(x_{k+1})$. f. Calculate $s_k = x_{k+1} - x_k$ and $y_k = \nabla f_{k+1} - \nabla f_k$. g. Curvature Check: If $s_k^T y_k > \epsilon$ (for some small $\epsilon > 0$), store the pair $(s_k, y_k)$. If the number of stored pairs exceeds $m$, discard the oldest pair $(s_{k-m}, y_{k-m})$. h. Check convergence criteria (e.g., $||\nabla f_{k+1}|| < \text{tolerance}$).

Advantages of L-BFGS

• Memory Efficiency: Its primary advantage. $O(md)$ memory makes it suitable for problems with a substantial number of parameters where storing $O(d^2)$ is impossible.
• Computational Efficiency: $O(md)$ computation per iteration is much faster than standard BFGS or Newton's method for large $d$.
• Good Convergence: Often achieves superlinear convergence rates, typically converging much faster (in terms of iterations) than first-order methods like Gradient Descent or Conjugate Gradient on many deterministic problems.

Disadvantages and Trade-offs

• Hyperparameter $m$: The size of the memory buffer, $m$, affects performance. A larger $m$ stores more curvature information, potentially leading to better search directions and faster convergence (fewer iterations), but increases the memory usage and computational cost per iteration. Typical values are small (e.g., 5-20).
• Line Search: Still requires a line search algorithm to find the step size $\alpha_k$, which can sometimes be computationally intensive itself, requiring multiple function and gradient evaluations per iteration.
• Curvature Condition: The update relies on the condition $s_k^T y_k > 0$, which implies that the function's curvature is positive in the direction of $s_k$. In non-convex regions, this might not hold. Implementations typically skip the update if this condition isn't met.
• Batch Method: L-BFGS, like standard BFGS and Newton's method, is fundamentally designed for deterministic (batch) optimization where the true gradient is computed using the entire dataset. Adapting it effectively to the stochastic settings common in deep learning (using mini-batches) is challenging, though some variants exist.

L-BFGS in Practice

L-BFGS is a workhorse algorithm in many scientific computing and optimization domains. In machine learning, it's often effective for:

• Training models where computing the full gradient is feasible (e.g., logistic regression, linear SVMs on moderately sized datasets).
• Problems where high precision is required near the optimum.
• Sometimes used as a final tuning step after initial training with stochastic methods like Adam or SGD.

While adaptive first-order methods like Adam are generally preferred for training large deep neural networks due to their compatibility with mini-batching and ability to navigate complex non-convex landscapes, understanding L-BFGS provides valuable insight into how curvature information can be leveraged for optimization, even under memory constraints. It represents a powerful bridge between first-order simplicity and second-order efficiency. Implementations are readily available in libraries like SciPy (scipy.optimize.minimize(method='L-BFGS-B')) and are used internally by some ML frameworks.