M5-T4 — Normal Equations

Learning objectives

After this lesson, you should be able to:

• Derive AᵀA x̂ = Aᵀb from the gradient and from orthogonality.
• Prove that AᵀA is symmetric, positive semidefinite, and invertible exactly when A has independent columns.
• Solve normal equations for a line fit by hand and for larger problems with a Cholesky factorisation.
• State that cond(AᵀA) = cond(A)² and explain what this does to accuracy.

Introduction and motivation

Every first-year fitting formula (the slope and intercept of a regression line) comes from the normal equations, and so do the closed-form expressions inside many estimators. They are worth knowing exactly. But a robot fitting a plane to millions of points, or a polynomial to a trajectory, needs to know when they fail. The answer is concrete: they work well for mild problems and lose digits quickly for ill-conditioned ones, because forming AᵀA squares the sensitivity of the problem.

The intuition

The previous lesson showed that the best fit leaves a residual perpendicular to every column of A. Write that as an equation: each column dotted with the residual is zero, which is Aᵀ(b − Ax) = 0. Rearranged, this is a square system for the unknowns, whatever the number of data points: (AᵀA) x = Aᵀb. A thousand measurements of two unknowns collapse to a 2 × 2 system.

Geometric and physical interpretation

Geometric view. The entries of AᵀA are the dot products between the columns of A (a Gram matrix), and the entries of Aᵀb are the dot products of the columns with the data. The normal equations say: the combination of columns that best matches b is the one for which, taking dot products with each column, the prediction and the data agree. When the columns are nearly parallel the Gram matrix is nearly singular. Physical view. In a line fit, AᵀA = [[Σx², Σx], [Σx, n]] contains the sums of the input values; fitting is solving two coupled equations about those sums. If all x values are almost equal (the data barely span the input range), the slope cannot be determined: the system is nearly singular.

The mathematics

AᵀA

the n × n Gram matrix of the columns: entry (i, j) is a_i·a_j.

Aᵀb

the n-vector of dot products of the columns with the data.

κ(A)

the condition number of A (ratio of the largest to the smallest singular value, Module 6).

The normal equations

AᵀA x̂ = Aᵀ b          x̂ = (AᵀA)⁻¹ Aᵀ b   (columns of A independent)

What it means. A square n × n system whose solution minimises the sum of squared residuals. The matrix AᵀA is symmetric and, for independent columns, positive definite, which is why it can be solved by Cholesky factorisation.

When is AᵀA invertible?

AᵀA invertible   ⇔   columns of A independent   ⇔   rank A = n

What it means. If the columns are dependent there is a nonzero x with Ax = 0, so (AᵀA)x = 0 and AᵀA is singular. Then the least-squares answer is not unique.

The numerical price

κ(AᵀA) = κ(A)²

What it means. The sensitivity of the problem is squared: with κ(A) = 10⁶ the normal equations behave like a problem with sensitivity 10¹², and roughly twelve digits of the sixteen available in double precision are lost.

Normal equations: AᵀA x̂ = Aᵀb. Valid and elegant, but they square the condition number, so use QR for anything ill-conditioned.

Step-by-step derivation: two routes to the normal equations, and why AᵀA is invertible

Route 1: set the gradient to zero. The cost is f(x) = xᵀAᵀAx − 2bᵀAx + bᵀb.

1. The gradient of xᵀMx for symmetric M is 2Mx, and the gradient of cᵀx is c.
2. So ∇f = 2AᵀAx − 2Aᵀb.
3. Setting ∇f = 0 gives AᵀAx = Aᵀb. Since the Hessian 2AᵀA is positive semidefinite, this stationary point is a global minimum.

Route 2: orthogonality. The residual r = b − Ax̂ must satisfy Aᵀr = 0 (previous lesson). That is Aᵀ(b − Ax̂) = 0, i.e. AᵀAx̂ = Aᵀb. Same equations, no calculus.

Properties of AᵀA. It is symmetric: (AᵀA)ᵀ = AᵀA. It is positive semidefinite: xᵀAᵀAx = ‖Ax‖² ≥ 0. It is invertible if and only if the columns of A are independent:

1. Suppose AᵀAx = 0. Multiply on the left by xᵀ: xᵀAᵀAx = ‖Ax‖² = 0.
2. So Ax = 0. If the columns are independent this forces x = 0, hence AᵀA has trivial null space and is invertible.
3. Conversely, if some x ≠ 0 has Ax = 0 then AᵀAx = 0, so AᵀA is singular.

The condition number squares. If A has singular values σ₁ ≥ … ≥ σn, then AᵀA has eigenvalues σ₁², …, σn² (Module 6). Its condition number is σ₁²/σn² = κ(A)². Squaring, in floating point with 16 digits, turns κ = 10⁶ into κ = 10¹², leaving about four correct digits instead of ten.

A worked example: a line fit by hand, and a hard fit by computer

By hand. Fit y = m x + c to (1, 1), (2, 2), (3, 2).

1. Form A = [[1, 1], [2, 1], [3, 1]] and b = (1, 2, 2). Then AᵀA = [[1 + 4 + 9, 1 + 2 + 3], [6, 3]] = [[14, 6], [6, 3]] and Aᵀb = (1·1 + 2·2 + 3·2, 1 + 2 + 2) = (11, 5).
2. Solve [[14, 6], [6, 3]] (m, c) = (11, 5). The determinant is 14·3 − 6·6 = 6. Cramer: m = (11·3 − 6·5)/6 = 3/6 = 0.5 and c = (14·5 − 6·11)/6 = 4/6 = 0.6667.
3. Check the orthogonality: the residual (−0.1667, 0.3333, −0.1667) has zero sum and zero dot product with (1, 2, 3) ✓.
4. Conditioning. κ(A) = 6.79 and κ(AᵀA) = 46.1 = 6.79²: squaring costs about one digit here, harmless.

By computer, a hard problem. Fit a degree-9 polynomial to 30 points on [0, 1] with known coefficients 1/(k + 1). The design matrix has κ(A) = 3.5·10⁶, so κ(AᵀA) = 1.2·10¹³. Solving the normal equations returns coefficients with a relative error of about 3·10⁻⁴; a QR-based solver returns an error of about 7·10⁻¹¹, a million times better, from exactly the same data.

Robotics application

Closed-form estimators. The slope and intercept formulas of a regression line, the mean and covariance of a Gaussian fit, and the sums used in plane-fitting code are all solutions of normal equations. For a few well-conditioned parameters, such as scale and offset of a range sensor, they are perfectly good.

Gauss-Newton. Nonlinear least squares in SLAM and bundle adjustment solves normal equations JᵀJ δ = −Jᵀr at every iteration, often with sparse Cholesky factorisation, because the matrix is large and well structured. The SLAM track's nonlinear least squares lesson develops it. Damping (Levenberg-Marquardt) adds a regulariser exactly of the kind introduced in the lesson on regularisation.

Information form. In estimation, AᵀA is the information matrix: how much the data constrain each parameter. A nearly singular AᵀA means some parameter combination is barely observable. This is the sensor-design reading of the equations; see Gaussian distributions and covariance.

Where they fail. High-degree polynomial trajectory fits, calibrations with nearly redundant parameters and long-baseline range fits are the situations where a robot's code should use QR.

Engineering implications and limitations

Practical rules

• Use normal equations only for small, well-conditioned problems. Check κ(A) (or the ratio of the largest to smallest diagonal entry of AᵀA as a crude warning). If it exceeds about 10⁴ switch to QR.
• Solve with Cholesky, never invert. AᵀA is symmetric positive definite, so use Cholesky factorisation and two triangular solves; do not compute (AᵀA)⁻¹.
• Centre and scale your data. Subtract the mean and divide by a typical size before fitting. A polynomial in x on [1000, 1001] is catastrophically ill-conditioned; the same polynomial in x − 1000.5 is fine.
• If Cholesky fails, the problem is rank-deficient. A failed factorisation of AᵀA is a diagnosis (dependent columns), not a bug to patch with a tiny diagonal addition without understanding why.
• Prefer QR or SVD in production. They solve the same problem without squaring the condition number; the next lesson explains.

Limits

The normal equations are exact in exact arithmetic; the only problem is rounding. For well-conditioned data, they are fast (cost about m n² + n³/3) and entirely adequate. The decision is about conditioning, not about right and wrong.

Common misconceptions

Engineering takeaways

1. Setting the gradient to zero (or the residual orthogonal to the columns) gives the normal equations AᵀA x̂ = Aᵀb.
2. AᵀA is symmetric positive semidefinite, and invertible exactly when the columns of A are independent.
3. The system size is n × n regardless of the number of data points; solve it with Cholesky.
4. κ(AᵀA) = κ(A)²: the normal equations can lose half the available digits.
5. Use them for small, well-scaled problems; use QR for anything ill-conditioned.

References and further reading

Primary papers and authoritative documentation.

• MIT 18.06SC Linear Algebra (OCW): lecture on least squares and the normal equations
• Boyd and Vandenberghe, Introduction to Applied Linear Algebra (VMLS): the chapter on least squares: the normal equations; the site did not load when we last checked
• Mathematics for Machine Learning: the chapter on linear regression: maximum likelihood and normal equations
• UCSD ECE276A: Sensing and Estimation in Robotics: least squares in estimation and the information matrix

Knowledge check

Three conceptual questions. Write an answer, then reveal the explanation. Your answers stay in your browser.

What’s next

Normal equations are elegant but fragile. QR decomposition solves the same problem by working on A directly with orthogonal factors, and is the method to reach for.

Related content

Nothing below is required to finish this lesson.

• Before
  M5-T3: Geometric Interpretation of Residuals

  See the least-squares answer as a projection onto the column space and the residual as orthogonal to it.

• Lab
  M5-L1: LiDAR Plane Fitting

  Fit the ground plane to a noisy LiDAR scan, find out why a naive fit fails on a wall, and fix it with a stable method.