The Normal Equation
Gradient
descent finds the best-fit line the way a hiker finds a valley floor in thick fog:
take a small step downhill, check the slope again, take another step, and repeat — hundreds or
thousands of times — until the steps stop making progress. It always gets there eventually, but it
crawls.
For linear regression, that crawl is often unnecessary. The
cost
function is a smooth, bowl-shaped surface with exactly one lowest point — no false
valleys to get stuck in. And a bowl has a very convenient property: at its single lowest point, the
slope is zero in every direction at once. Instead of hunting for that point step
by step, you can write down the equations that say "the slope is zero here" and solve them
directly. That direct, one-shot solution is the normal equation — no hiking, no
learning rate to tune, no iterations at all. You calculate once, and you're already standing at the
bottom.
The formula
Stack every training example's features into rows of a matrix, and every label into a matching
vector. Then one closed-form expression hands back the optimal weights in a single step.
For a feature matrix X (one row per example, one column per feature —
plus a column of ones for the bias) and a label vector \vec{y}, the
weight vector that minimizes the squared-error cost is
\vec{\theta} = (X^{\mathsf T} X)^{-1} X^{\mathsf T} \vec{y}.
- X — the feature matrix: every row is one training
example, every column one feature.
- \vec{y} — the label vector: the true value for
each example, in the same row order as X.
- \vec{\theta} — the solved weight vector: the
intercept and every feature's coefficient, all at once.
Every piece of that formula is ordinary linear algebra you've already met:
transpose
X to build X^{\mathsf T}, multiply the two
matrices together, take the
inverse
of the result, and finish with a
matrix-times-vector
multiplication against X^{\mathsf T}\vec{y}. There's no calculus left to
do by hand — it's all baked into the formula already.
Try to beat the formula
Adjust your line and compare your cost with the optimal cost the normal equation
achieves (the faint line is its answer). No matter how carefully you tune by hand, you can only
ever match it — never beat it. That faint line is the provably best straight-line fit, computed in
one shot from the formula above, not found by trial and error.
Worked example: solving it by hand
Here is a tiny dataset with one feature — four points, chosen to be small enough to solve on
paper:
x = [1,\ 2,\ 3,\ 4] \qquad y = [3,\ 5,\ 6,\ 8]
Step 1 — build X and \vec{y}.
Add a column of ones so the formula can solve for the intercept too:
X = \begin{pmatrix} 1 & 1 \\ 1 & 2 \\ 1 & 3 \\ 1 & 4 \end{pmatrix} \qquad \vec{y} = \begin{pmatrix} 3 \\ 5 \\ 6 \\ 8 \end{pmatrix}
Step 2 — compute X^{\mathsf T}X and
X^{\mathsf T}\vec{y}. Every entry is just a sum over the four
rows:
X^{\mathsf T}X = \begin{pmatrix} 4 & 10 \\ 10 & 30 \end{pmatrix} \qquad X^{\mathsf T}\vec{y} = \begin{pmatrix} 22 \\ 63 \end{pmatrix}
Step 3 — invert X^{\mathsf T}X. For a
2\times 2 matrix the determinant is
4 \times 30 - 10\times 10 = 20, so:
(X^{\mathsf T}X)^{-1} = \frac{1}{20}\begin{pmatrix} 30 & -10 \\ -10 & 4 \end{pmatrix}
Step 4 — multiply through.
\vec{\theta} = (X^{\mathsf T}X)^{-1}X^{\mathsf T}\vec{y} works out to
\theta_0 = 1.5 (the intercept) and \theta_1 = 1.6
(the slope) — in one calculation, no iteration required. Run the same arithmetic
below and check it matches:
function normalEquation(x: number[], y: number[]): { intercept: number; slope: number } {
const n = x.length;
const sumX = x.reduce((s, xi) => s + xi, 0);
const sumY = y.reduce((s, yi) => s + yi, 0);
const sumXX = x.reduce((s, xi) => s + xi * xi, 0);
const sumXY = x.reduce((s, xi, i) => s + xi * y[i], 0);
// XtX = [[n, sumX], [sumX, sumXX]], Xty = [sumY, sumXY]
const det = n * sumXX - sumX * sumX;
const intercept = (sumXX * sumY - sumX * sumXY) / det;
const slope = (n * sumXY - sumX * sumY) / det;
return { intercept, slope };
}
const x = [1, 2, 3, 4];
const y = [3, 5, 6, 8];
const theta = normalEquation(x, y);
console.log("intercept:", theta.intercept.toFixed(2));
console.log("slope:", theta.slope.toFixed(2));
Plug x=1,2,3,4 back into
1.5 + 1.6x and you get 3.1,\ 4.7,\ 6.3,\ 7.9
— close to, but not exactly, the original y values. That's expected: no
straight line fits four slightly-noisy points perfectly, and the normal equation guarantees this
is the line that minimizes the total squared error, not one that passes through every point. It's
the exact same answer gradient descent would eventually crawl its way to, given enough steps.
Normal equation or gradient descent?
Both methods find the same optimal line — they trade off differently, and the right choice depends
on how big your problem is.
- The normal equation is exact and one-shot: no learning rate to tune, no
worrying about convergence. Its cost is dominated by inverting a
k\times k matrix (k being the number of
features), and matrix inversion scales roughly with k^3 — triple the
features and the work goes up by roughly 3^3 = 27\times. That's fine
for tens or hundreds of features, brutal for thousands, and impossible outright if
X^{\mathsf T}X isn't invertible.
- Gradient descent is iterative and needs a learning rate, but each step is
cheap, and it scales gracefully to huge numbers of features or examples. It's also the
only option for models like
logistic
regression or neural networks, which have no tidy closed-form formula at all.
Rough rule of thumb: with a handful to a few thousand features, the normal equation is often
simplest and plenty fast. Beyond that — or for any model without a magic formula — gradient
descent (or one of its many modern variants) takes over.
-
X^{\mathsf T}X must be invertible. If two features
are linearly
dependent (say, one column is just another rescaled, like a price in dollars and the
same price in cents), or if there are more features than training examples, then
X^{\mathsf T}X becomes singular — it has no inverse, and the formula
simply cannot be computed. The usual fallback is gradient descent, or adding
regularization,
which nudges the matrix back into being invertible.
-
It stops being practical long before it becomes impossible. Inverting a matrix
gets dramatically slower as it grows — doubling the number of features doesn't just double the
work, it multiplies it several times over. With a modest handful of features this finishes
instantly; with thousands of features (common in real datasets — think one column per word in a
vocabulary) it can grind to a halt. This is a big reason gradient descent dominates in
large-scale machine learning today, even though each of its individual steps is "less exact"
than the normal equation's single leap.
The normal equation looks thoroughly modern — matrices, inverses, a tidy one-liner — but the idea
underneath it is over two centuries old. It's exactly the
Gauss-and-Legendre
least-squares method that predicted where the "lost" dwarf planet Ceres would reappear in 1801,
wearing modern matrix notation. Gauss didn't have X^{\mathsf T}X
notation or a computer to invert it with — he solved the same equations by hand, the long way — but
the underlying formula is identical to the one you just used above.
Ever added a "trendline" to a scatter chart in a basic spreadsheet program? For an ordinary,
small, one-feature dataset, that trendline is almost certainly computed with the normal equation —
it's a small, one-shot calculation that a spreadsheet can finish instantly, with no need for
anything as heavyweight as iterative gradient descent. Every time you've clicked "add trendline,"
there's a good chance a tiny matrix inversion just happened behind the scenes.