Least Squares Problems

The Problem¶

Given $A \in \mathbb{R}^{m \times n}$ with $m > n$ and $b \in \mathbb{R}^m$, the system $Ax = b$ typically has no solution.

Why? The vector $b$ usually doesn’t lie in the column space of $A$.

Goal: Find $\hat{x}$ that minimizes the residual:

Linear Regression Example¶

The most common application: fitting a model to data.

Example: Given $N$ data points $(t_i, y_i)$, fit a polynomial $p(t) = c_0 + c_1 t + c_2 t^2$:

With $N > 3$ data points, this system is overdetermined.

Geometric Interpretation¶

The least squares solution finds the point in $\text{R}(A)$ closest to $b$:

        b
       /|
      / |  residual r = b - Ax̂
     /  |  is perpendicular to R(A)
    Ax̂--+---- R(A) (column space)

Key insight: The residual $r = b - A\hat{x}$ is orthogonal to the column space of $A$.

The Normal Equations¶

From the orthogonality condition $r \perp \text{R}(A)$, we can derive:

Derivation: We need $(b - A\hat{x}) \perp Av$ for all $v \in \mathbb{R}^n$.

This requires $A^T(b - A\hat{x}) = 0$, giving $A^T A \hat{x} = A^T b$.

Why Are They Called “Normal” Equations?¶

The name comes from the fact that the residual is normal (perpendicular) to the column space—not because they’re “standard” equations.

Properties of $A^T A$¶

When $A$ has full column rank:

Minimization Viewpoint¶

The least squares problem is equivalent to minimizing:

Taking the gradient and setting to zero:

yields the normal equations.

Observation: This is a quadratic function with positive definite Hessian $2A^T A$, so there’s a unique global minimum.

Multiple Linear Regression¶

In statistics notation, the least squares problem for regression:

where:

• $X \in \mathbb{R}^{N \times (k+1)}$ is the design matrix (observations of explanatory variables)

• $y \in \mathbb{R}^N$ is the response vector

• $\beta \in \mathbb{R}^{k+1}$ are the regression coefficients

The columns of $X$ typically include a column of ones (for the intercept).

The Pseudoinverse¶

The matrix $(A^T A)^{-1} A^T$ is called the Moore-Penrose pseudoinverse $A^+$:

So $\hat{x} = A^+ b$.

Properties:

• $A^+ A = I_n$ (left inverse)

• $A A^+ \neq I_m$ in general (not a right inverse)

• $A A^+$ projects onto $\text{R}(A)$

Numerical Considerations¶

Example: Polynomial Fitting¶

import numpy as np

# Data points
t = np.array([0, 1, 2, 3, 4])
y = np.array([1.0, 2.1, 3.9, 8.2, 15.8])

# Design matrix for quadratic fit
X = np.column_stack([np.ones_like(t), t, t**2])

# Normal equations (don't do this!)
# beta_bad = np.linalg.solve(X.T @ X, X.T @ y)

# Better: use lstsq which uses QR internally
beta, residuals, rank, s = np.linalg.lstsq(X, y, rcond=None)

Residual Analysis¶

After finding $\hat{x}$, the residual is:

Properties:

• $\|r\|_2^2$ is the sum of squared errors (SSE)

• $r \perp \text{R}(A)$

• $A^T r = 0$

The residual measures how well the model fits the data.

Summary¶