Newton's Method

Why It Matters

---

Newton's Update Rule

The intuition is simple: gradient descent moves in the direction of steepest descent (-∇f), while Newton's method moves in a direction that accounts for the curvature of the function. By multiplying the gradient by H⁻¹, Newton's method effectively rescales each direction according to how curved the function is in that direction.

---

Quadratic Convergence

---

Hessian Matrix

---

Damped Newton Method

The typical damped Newton algorithm at iteration k:

1. Compute gradient g_k = ∇f(θ_k) and Hessian H_k = ∇²f(θ_k)
2. Solve H_k d_k = -g_k for the Newton direction d_k
3. If H_k is not positive definite, modify it (e.g., add λI to make it positive definite)
4. Perform line search to find step size α_k
5. Update θ_{k+1} = θ_k + α_k d_k

---

Quasi-Newton Methods

---

Gauss-Newton Method

---

Comparison with First-Order Methods

---

Python Implementation

import numpy as np
from scipy.optimize import minimize, least_squares

# --- Newton's Method (using scipy with Hessian) ---

def rosenbrock(x):
    """Rosenbrock function: f(x,y) = (1-x)^2 + 100(y-x^2)^2"""
    return (1 - x[0])**2 + 100 * (x[1] - x[0]**2)**2

def rosenbrock_grad(x):
    """Gradient of Rosenbrock function."""
    dfdx = -2*(1 - x[0]) - 400*x[0]*(x[1] - x[0]**2)
    dfdy = 200*(x[1] - x[0]**2)
    return np.array([dfdx, dfdy])

def rosenbrock_hessian(x):
    """Hessian of Rosenbrock function."""
    h11 = 2 - 400*(x[1] - 3*x[0]**2)
    h12 = -400*x[0]
    h21 = -400*x[0]
    h22 = 200
    return np.array([[h11, h12], [h21, h22]])

# Newton's method with exact Hessian
result_newton = minimize(
    rosenbrock,
    x0=[-1.0, 1.0],
    jac=rosenbrock_grad,
    hess=rosenbrock_hessian,
    method='Newton-CG'
)
print(f"Newton: x = {result_newton.x}, f(x) = {result_newton.fun:.2e}")

# --- BFGS (Quasi-Newton) ---

result_bfgs = minimize(
    rosenbrock,
    x0=[-1.0, 1.0],
    method='BFGS'
)
print(f"BFGS:   x = {result_bfgs.x}, f(x) = {result_bfgs.fun:.2e}")

# --- L-BFGS-B (Bounded Quasi-Newton) ---

result_lbfgs = minimize(
    rosenbrock,
    x0=[-1.0, 1.0],
    method='L-BFGS-B',
    bounds=[(-2, 2), (-2, 2)]
)
print(f"L-BFGS: x = {result_lbfgs.x}, f(x) = {result_lbfgs.fun:.2e}")

# --- Gauss-Newton for Nonlinear Least Squares ---

def model(x, t):
    """Model: y = a * exp(b * t)"""
    return x[0] * np.exp(x[1] * t)

def residuals(x, t, y):
    """Residuals for least squares."""
    return y - model(x, t)

# Data
t_data = np.array([1.0, 2.0, 3.0, 4.0, 5.0])
y_data = np.array([2.71, 7.39, 20.09, 54.60, 148.41])

# Levenberg-Marquardt (Gauss-Newton variant)
result_gn = least_squares(
    residuals,
    x0=[1.0, 0.5],
    args=(t_data, y_data),
    method='trf'
)
print(f"Gauss-Newton: a = {result_gn.x[0]:.4f}, b = {result_gn.x[1]:.4f}")

# --- Comparison of convergence rates ---

from scipy.optimize import minimize as minimize_comparison

def rosenbrock_only(x):
    return (1 - x[0])**2 + 100 * (x[1] - x[0]**2)**2

methods = ['BFGS', 'L-BFGS-B', 'Nelder-Mead', 'Powell']
for method in methods:
    res = minimize_comparison(
        rosenbrock_only,
        x0=[-1.0, 1.0],
        method=method,
        options={'maxiter': 1000}
    )
    print(f"{method:12s}: f(x) = {res.fun:.2e}, iters = {res.nit}")

---

Applications in AI/ML

Key applications include:

• Logistic regression: Newton's method (IRLS) converges in 5-6 iterations for most problems
• Gaussian process regression: Matrix inversion is equivalent to solving the normal equations
• Variational inference: Natural gradients for mean-field and full-rank approximations
• Reinforcement learning: Trust region policy optimization (TRPO) uses a second-order approximation
• Hyperparameter optimization: L-BFGS for tuning differentiable hyperparameters
• Computer vision: Bundle adjustment in SLAM uses Gauss-Newton/Levenberg-Marquardt
• Scientific computing: Nonlinear PDE solvers, curve fitting, parameter estimation

---

Common Mistakes

---

Interview Questions

---

Practice Problems

---

Quick Reference

---

Cross-References