Gradient Descent & Convergence

Overview & Motivation

In Convex Analysis we established the geometric and analytic infrastructure that makes optimization tractable: convex sets are closed under line segments, convex functions have a single global minimum, and smooth convex functions admit tangent-line lower bounds. That machinery tells us what to optimize but not how. This topic is about the how.

Gradient descent is the simplest first-order algorithm for unconstrained minimization, and arguably the most important. The idea is disarmingly simple: to minimize $f$, start at an initial guess $x_0$ and repeatedly take steps in the direction that decreases $f$ fastest — the negative gradient $-\nabla f(x_k)$. Each step replaces the function by its local linear approximation and moves to where that approximation is most favorable.

What makes gradient descent remarkable is not the algorithm itself but the convergence theory that governs it. With the right assumptions — smoothness, convexity, and sometimes strong convexity — we can prove exactly how fast the iterates approach the optimum, and these rates turn out to be tight: no first-order method can do better in the worst case (at least without acceleration).

What We Cover

1. The Gradient Descent Algorithm — the update rule and its geometric interpretation as steepest descent.
2. Smoothness & the Descent Lemma — how $L$-smoothness provides a quadratic upper bound that guarantees progress at each step.
3. Convergence for Convex Functions — the $O(1/k)$ sublinear rate, proven via telescoping.
4. Strong Convexity & Linear Convergence — the condition number $\kappa = L/\mu$ and the geometric rate $(1 - 1/\kappa)^k$.
5. Step Size Selection — fixed, exact line search, and Armijo backtracking.
6. Nesterov Accelerated Gradient — the momentum trick that achieves the optimal $O(1/k^2)$ rate.
7. Coordinate Descent — updating one variable at a time, and when that beats full-gradient steps.
8. Mirror Descent & Bregman Divergence — generalizing Euclidean geometry for constrained domains.
9. Stochastic Gradient Descent — noisy gradients, mini-batches, and decreasing step sizes.
10. Computational Notes — NumPy implementations and convergence diagnostics.

---

The Gradient Descent Algorithm

The gradient $\nabla f(x)$ points in the direction of steepest increase of $f$ at $x$. To decrease $f$, we walk in the opposite direction.

Geometrically, each iterate replaces $f$ by its first-order Taylor approximation $f(x_k) + \nabla f(x_k)^\top (x - x_k)$ and moves in the direction that decreases this linear model fastest. The step size $\eta$ controls how far we trust the linear approximation — too large and we overshoot; too small and progress is glacially slow.

For a quadratic $f(x) = \frac{1}{2}x^\top A x$ with $A$ symmetric positive definite, the gradient is $\nabla f(x) = Ax$ and the update simplifies to $x_{k+1} = (I - \eta A)x_k$. The convergence behavior of this linear map is entirely determined by the eigenvalues of $I - \eta A$: convergence requires $|1 - \eta \lambda_i| < 1$ for every eigenvalue $\lambda_i$, which means $\eta$ must lie in $(0, 2/\lambda_{\max})$.

---

Smoothness & the Descent Lemma

The key to analyzing gradient descent is smoothness — a Lipschitz condition on the gradient that provides a global quadratic upper bound on $f$.

Smoothness says the gradient cannot change too fast. This is exactly the condition that allows us to guarantee progress at each gradient step: if the gradient doesn’t change too fast between $x_k$ and $x_{k+1}$, then the decrease predicted by the linear approximation actually holds (up to a quadratic correction).

The Descent Lemma tells us that at each point $x$, the function $f$ is sandwiched between a linear lower bound (from convexity, if $f$ is convex) and a quadratic upper bound (from $L$-smoothness). The gradient step minimizes the upper bound exactly, and the decrease is proportional to $\|\nabla f(x)\|^2$.

---

Convergence for Convex Functions

Armed with the Descent Lemma, we can prove the first convergence result. For a convex, $L$-smooth function, gradient descent with step size $\eta = 1/L$ achieves an $O(1/k)$ rate — sublinear, but dimension-free.

The $O(1/k)$ rate tells us that to reach $\varepsilon$-accuracy, we need $k \sim L\|x_0 - x^*\|^2 / \varepsilon$ iterations. This depends on the initial distance to the optimum but — crucially — not on the dimension $n$.

---

Strong Convexity & Linear Convergence

The $O(1/k)$ rate is slow. If we strengthen our assumption from plain convexity to strong convexity, the rate improves dramatically from sublinear to linear (i.e., geometric).

Strong convexity says that $f$ curves upward at least as fast as a quadratic with parameter $\mu$. Together with $L$-smoothness, this gives a “quadratic sandwich”:

$f(x) + \nabla f(x)^\top(y-x) + \frac{\mu}{2}\|y-x\|^2 \leq f(y) \leq f(x) + \nabla f(x)^\top(y-x) + \frac{L}{2}\|y-x\|^2$

The linear rate means $\log(f(x_k) - f(x^*))$ decreases at a constant rate per iteration — the convergence plot is a straight line on a log scale. To reach $\varepsilon$-accuracy, we need $k \sim \kappa \log(1/\varepsilon)$ iterations. The dependence on $1/\varepsilon$ is logarithmic (much better than the $1/\varepsilon$ for the convex case), but the price is a linear dependence on $\kappa$. Ill-conditioned problems (large $\kappa$) converge slowly because the contours are highly elongated ellipses, and gradient descent zigzags across the narrow valley.

κ = 5Step size:1/L2/(L+μ)Too largePlayResetStep 0 / 110 | f(x) = 7.92e+0 | η = 0.2000

---

Step Size Selection

So far we’ve used the fixed step size $\eta = 1/L$, which requires knowing $L$ in advance. In practice, there are several strategies.

Exact line search

This strategy chooses the step size that minimizes $f$ along the gradient direction:

$\eta_k = \arg\min_{\eta > 0} f(x_k - \eta \nabla f(x_k))$

For a quadratic $f(x) = \frac{1}{2}x^\top A x$, this has a closed-form solution:

$\eta_k^* = \frac{\|\nabla f(x_k)\|^2}{\nabla f(x_k)^\top A \nabla f(x_k)}$

Exact line search always converges at least as fast as fixed step size but is usually impractical for non-quadratic objectives.

---

Nesterov Accelerated Gradient

Gradient descent with step size $1/L$ converges at rate $O(1/k)$ for convex functions. Can we do better? Nesterov showed in 1983 that the answer is yes — and that $O(1/k^2)$ is optimal among all first-order methods.

The critical difference from vanilla gradient descent: the gradient is evaluated at the extrapolated point $y_k$, not at the current iterate $x_k$. The extrapolation step $y_k$ adds a “momentum” term — it looks ahead in the direction the iterates have been moving, then computes the gradient at that look-ahead point.

We state this theorem without proof — the proof requires a carefully constructed Lyapunov function argument that we refer the reader to Nesterov’s original paper.

κ = 20Show momentum arrowsPlayResetStep 0 / 150

---

Coordinate Descent

Instead of computing the full gradient $\nabla f(x) \in \mathbb{R}^n$, coordinate descent updates one coordinate at a time. This is attractive when $n$ is large and computing a single partial derivative is cheap.

The characteristic visual signature of coordinate descent is the “staircase” pattern on contour plots: each step moves along a coordinate axis, creating axis-aligned zigzags. Compare this with gradient descent, which moves along the gradient direction (which is generally not axis-aligned).

---

Mirror Descent & Bregman Divergence

Gradient descent implicitly uses Euclidean distance to measure progress: the update $x_{k+1} = x_k - \eta \nabla f(x_k)$ can be written as

$x_{k+1} = \arg\min_x \left\{ \nabla f(x_k)^\top x + \frac{1}{2\eta}\|x - x_k\|^2 \right\}$

But Euclidean distance is not always the right geometry. When optimizing over the probability simplex $\Delta = \{x \geq 0 : \sum x_i = 1\}$, the Euclidean distance treats points near the boundary the same as interior points, which is wasteful — near a vertex, most directions leave the simplex. Mirror descent replaces $\|x - y\|^2$ with a Bregman divergence matched to the geometry of the domain.

When $\phi(x) = \frac{1}{2}\|x\|^2$, the Bregman divergence reduces to $D_\phi(x,y) = \frac{1}{2}\|x - y\|^2$ — ordinary Euclidean distance. When $\phi(x) = \sum_i x_i \log x_i$ (the negative entropy), the Bregman divergence is the KL divergence $D_\phi(x, y) = \sum_i x_i \log(x_i / y_i) - x_i + y_i$.

The exponentiated gradient update naturally maintains positivity and the simplex constraint without explicit projection. Mirror descent with negative entropy is the natural algorithm for the simplex: its iterates curve along the simplex boundary rather than cutting through the interior in straight lines.

The connection to Proximal Methods is deep: mirror descent is a special case of the proximal point algorithm with Bregman divergence replacing Euclidean distance. Proximal operators generalize the gradient step to handle non-smooth regularizers like the $\ell_1$ norm.

Objective:Optimum at e₁ (c = [0.1, 0.8, 0.3])Optimum at e₂ (c = [0.3, 0.1, 0.8])Optimum at e₃ (c = [0.8, 0.3, 0.1])Near center (c = [0.3, 0.3, 0.4])η = 2.0PlayResetStep 0 / 40

---

Stochastic Gradient Descent

In machine learning, the objective is often an average over a large dataset: $f(x) = \frac{1}{N}\sum_{i=1}^N f_i(x)$. Computing the full gradient $\nabla f(x) = \frac{1}{N}\sum_i \nabla f_i(x)$ requires $N$ gradient evaluations per step, which is expensive when $N$ is large. Stochastic gradient descent (SGD) replaces the full gradient with a noisy but unbiased estimate.

The step size schedule is crucial. With a fixed step size $\eta$, SGD converges to a neighborhood of the optimum but cannot reach it — the noise prevents the iterates from settling down. A decreasing schedule $\eta_k = O(1/k)$ allows convergence to the exact optimum, at the cost of slower progress.

---

Computational Notes

Here we collect practical implementations. The gradient_descent function below supports fixed step size and Armijo backtracking, and returns the full convergence history for diagnostic plotting.

def gradient_descent(f, grad_f, x0, eta=None, L=None, tol=1e-8, max_iter=1000,
                     line_search='fixed'):
    """
    General-purpose gradient descent.

    Parameters
    ----------
    f : callable — objective function
    grad_f : callable — gradient of f
    x0 : ndarray — initial point
    eta : float — fixed step size (required if line_search='fixed')
    L : float — Lipschitz constant; used for eta=1/L if eta not given
    tol : float — gradient norm tolerance
    max_iter : int — iteration budget
    line_search : str — 'fixed' or 'armijo'

    Returns
    -------
    dict with keys 'x', 'f_vals', 'grad_norms', 'trajectory', 'n_iter'
    """
    if eta is None and L is not None:
        eta = 1.0 / L

    x = x0.copy()
    f_vals = [f(x)]
    grad_norms = [np.linalg.norm(grad_f(x))]
    trajectory = [x.copy()]

    for k in range(max_iter):
        g = grad_f(x)
        if np.linalg.norm(g) < tol:
            break

        if line_search == 'armijo':
            step = armijo_backtracking(f, grad_f, x, -g)
        else:
            step = eta

        x = x - step * g
        f_vals.append(f(x))
        grad_norms.append(np.linalg.norm(grad_f(x)))
        trajectory.append(x.copy())

    return {
        'x': x,
        'f_vals': np.array(f_vals),
        'grad_norms': np.array(grad_norms),
        'trajectory': np.array(trajectory),
        'n_iter': len(f_vals) - 1
    }

The Armijo backtracking function:

def armijo_backtracking(f, grad_f, x, d, alpha0=1.0, beta=0.5, c=1e-4, max_iter=50):
    """Armijo backtracking line search."""
    alpha = alpha0
    fx = f(x)
    gx = grad_f(x)
    for _ in range(max_iter):
        if f(x + alpha * d) <= fx + c * alpha * (gx @ d):
            return alpha
        alpha *= beta
    return alpha

Nesterov accelerated gradient:

def nesterov_gd(A, x0, n_iter=100):
    """Nesterov accelerated gradient for f(x) = 0.5 x^T A x."""
    L = eigvalsh(A)[-1]
    eta = 1.0 / L
    x, y = x0.copy(), x0.copy()
    traj = [x0.copy()]
    t = 1.0
    for k in range(n_iter):
        x_new = y - eta * (A @ y)
        t_new = (1 + np.sqrt(1 + 4 * t**2)) / 2
        y = x_new + ((t - 1) / t_new) * (x_new - x)
        x, t = x_new, t_new
        traj.append(x.copy())
    return np.array(traj)

Cyclic coordinate descent:

def cyclic_cd(A, x0, n_iter=50):
    """Cyclic coordinate descent for f(x) = 0.5 x^T A x."""
    x = x0.copy()
    traj = [x.copy()]
    for _ in range(n_iter):
        for i in range(len(x0)):
            grad_i = A[i, :] @ x
            x[i] -= grad_i / A[i, i]
            traj.append(x.copy())
    return np.array(traj)

Mirror descent on the simplex (exponentiated gradient):

def mirror_descent_simplex(grad_f, x0, eta, n_iter=50):
    """Mirror descent with neg-entropy on the probability simplex."""
    x = x0.copy()
    traj = [x.copy()]
    for _ in range(n_iter):
        g = grad_f(x)
        x_new = x * np.exp(-eta * g)
        x_new /= x_new.sum()
        x = x_new
        traj.append(x.copy())
    return np.array(traj)

SGD with mini-batch averaging:

def sgd_minibatch(A, x0, eta_schedule, noise_std=1.0, batch_size=1,
                  n_iter=500, seed=42):
    """SGD with Gaussian gradient noise and mini-batch averaging."""
    rng = np.random.RandomState(seed)
    x = x0.copy()
    traj = [x.copy()]
    for k in range(n_iter):
        eta = eta_schedule(k)
        true_grad = A @ x
        noise = rng.randn(batch_size, len(x)) * noise_std
        x = x - eta * (true_grad + noise.mean(axis=0))
        traj.append(x.copy())
    return np.array(traj)

Convergence diagnostics. When running gradient descent in practice, monitor four quantities: (1) the objective value $f(x_k)$, which should decrease monotonically for deterministic GD; (2) the gradient norm $\|\nabla f(x_k)\|$, which indicates proximity to a stationary point; (3) the step size history (if using backtracking); and (4) the trajectory in parameter space (for low-dimensional problems).

---

Connections & Further Reading

Gradient descent sits at the center of the Optimization track, building on convex analysis and connecting forward to proximal methods and duality.

Topic	Connection
Convex Analysis	Every convergence proof here relies on convexity and smoothness established in Convex Analysis. The Descent Lemma is the $L$-smooth analogue of the first-order convexity condition; strong convexity provides the quadratic lower bound that yields linear convergence.
The Spectral Theorem	The condition number $\kappa = L/\mu$ for quadratics equals $\lambda_{\max}/\lambda_{\min}$ of the Hessian. Ill-conditioned Hessians (large eigenvalue spread) cause the zigzag behavior that slows gradient descent.
Singular Value Decomposition	For least-squares problems $\min \|Ax - b\|^2$, the condition number of $A^\top A$ governs GD convergence. The SVD reveals this as $\sigma_{\max}^2/\sigma_{\min}^2$, connecting linear algebra to optimization dynamics.
Proximal Methods	Mirror descent generalizes to proximal operators: $\text{prox}_{\eta f}(x) = \arg\min_y \{f(y) + \frac{1}{2\eta}\|y - x\|^2\}$. Proximal gradient descent handles composite objectives $f + g$ where $g$ is non-smooth but has a cheap proximal operator (e.g., $\ell_1$ regularization).
Lagrangian Duality & KKT	Constrained optimization extends the unconstrained theory here. The KKT conditions generalize “gradient equals zero” to include constraints, and dual gradient descent operates on the Lagrangian dual function.

The Optimization Track

Convex Analysis
    ├── Gradient Descent & Convergence (this topic)
    │       └── Proximal Methods
    └── Lagrangian Duality & KKT