Mathematics1901foundational10 min read

On Lines and Planes of Closest Fit to Systems of Points in Space

عن الخطوط والمستويات الأقرب ملاءمةً لمنظومات النقاط في الفضاء

Pearson, K. — Philosophical Magazine

The problem

In the sciences of the late 1800s, researchers measured many variables for each observation — skull dimensions, meteorological readings, biological traits — and ended up with clouds of points in high-dimensional space. There was no systematic method to find the line, plane, or hyperplane that best summarizes such a cloud. The ordinary least-squares regression of the time minimized vertical distances to a line, which assumed one variable was "the answer" and the others were inputs. But when all variables are equally important — when you simply want the best geometric summary of the data — vertical distances are wrong. A new geometric criterion was needed.

The contribution

Pearson proposed minimizing the sum of squared perpendicular () distances from each point to the fitting line or plane — treating all variables symmetrically. Using calculus alone (no matrix algebra yet), he proved that the best-fit line passes through the data's center of mass and points in the direction of maximum . The best-fit plane contains the two directions of greatest variance, and so on for higher dimensions. This is the first formal description of what we now call Analysis (). The key insight: the directions that maximize spread are exactly the directions that minimize reconstruction error.

The impact

PCA became one of the most widely used algorithms in all of science. In it serves as a preprocessing step for high-dimensional data, a visualization tool for complex datasets, and the conceptual ancestor of autoencoders and . Every modern dimensionality-reduction technique — from kernel PCA to t-SNE to UMAP — descends from or reacts against Pearson's original geometric question. Over a century later, PCA remains the first tool reached for when a dataset has too many dimensions.

Imagine a photographer trying to capture a sculpture. She can only take a single flat photograph — one angle. She walks around the sculpture looking for the angle that shows the most detail: the silhouette that is widest, where features don't overlap and hide each other.

That "best angle" is the first principal component. The second-best angle — the one that reveals the most additional detail, viewed at right angles to the first — is the second principal component. Pearson's paper is the mathematical answer to: "what angle should the photographer choose?"

The problem: how do you summarize a cloud of points?

Suppose you measure the height and weight of 200 people. You get a scatter of 200 dots in 2D. Now add age, blood pressure, cholesterol — suddenly you have 200 dots in 5D. The data lives in a space too large to see.

You want a summary: a single line or flat surface that passes through the data cloud and captures its overall shape. But how do you define "best fit"?

The conventional method — ordinary least-squares regression — finds the line that minimizes vertical distances. This treats one variable as the output and the others as inputs. It works for prediction, but not for summarizing a cloud where all variables are equally important: if you swap which axis is "vertical," you get a different line.

Pearson's insight was to change the ruler. Instead of vertical distances, measure the perpendicular (shortest) distance from each point to the line. This treats all variables symmetrically — no variable is special. The line that minimizes the sum of squared perpendicular distances is the true geometric backbone of the data.

Open in Lab
Drag the line to rotate it. Watch how the total squared perpendicular distance changes — the minimum is PC1.
The demo wakes as you arrive…

The dual insight: minimum distance = maximum spread

Pearson discovered something elegant: the line that minimizes perpendicular distances from the points is the same line that maximizes the variance (spread) of the projected points along it.

Why? By the Pythagorean theorem. Each point's squared distance from the data center splits exactly into two parts: the squared distance along the line (the ) plus the squared to the line. The total of all points' squared distances from the center is fixed — it's a property of the data, not of the line. So if the perpendicular part shrinks, the projection part must grow by exactly the same amount:

Total (fixed) = Variance along the line + Perpendicular error

Minimizing the error is therefore identical to maximizing the variance. These are not two separate criteria — they are two sides of the same coin, joined by Pythagoras.

∑i=1n∥xi−xˉ∥2=∑i=1nproji2+∑i=1nperpi2\sum_{i=1}^{n} \|\mathbf{x}_i - \bar{\mathbf{x}}\|^2 = \sum_{i=1}^{n} \text{proj}^2_i + \sum_{i=1}^{n} \text{perp}^2_i
The Pythagorean decomposition — the heart of PCA — Total squared distance from the mean (fixed) = variance of projections along the line + sum of squared perpendicular errors. Maximizing one side automatically minimizes the other.
Open in Lab
Watch the Pythagorean split in real time — as one bar grows, the other must shrink by exactly the same amount.
The demo wakes as you arrive…

The recipe: covariance matrix → eigenvectors → principal components

Pearson solved the problem with calculus. The modern formulation (due to Hotelling, 1933) restates the same solution in the cleaner language of linear algebra, and it proceeds in three steps:

Step 1 — Center the data. Subtract the mean of each variable so the cloud sits at the origin. This ensures the best-fit line passes through the center.

Step 2 — Compute the . This symmetric matrix Σ\Sigma encodes how every pair of variables co-varies. The diagonal entries are the individual variances; the off-diagonal entries are the covariances. It is a compact map of the cloud's shape.

Step 3 — Find the eigenvectors and eigenvalues of Σ\Sigma. Each points in a principal direction — an axis of the cloud's elliptical shape. Its tells you how much variance sits along that axis. The eigenvector with the largest eigenvalue is the first principal component (PC1); the next largest gives PC2; and so on. All eigenvectors are mutually orthogonal — perpendicular to each other — so the components form a clean, non-overlapping coordinate system for the data.

Σ=1n∑i=1n(xi−xˉ)(xi−xˉ)⊤\Sigma = \frac{1}{n} \sum_{i=1}^{n} (\mathbf{x}_i - \bar{\mathbf{x}})(\mathbf{x}_i - \bar{\mathbf{x}})^\top
The covariance matrix — a compact map of the data's shape — Each entry Σᵢⱼ measures how variables i and j move together. The eigenvectors of this matrix are the principal component directions; the eigenvalues are the variances along them.
Σv=λv\Sigma \mathbf{v} = \lambda \mathbf{v}
The eigenvalue equation — finding the axes of the ellipse — v is an eigenvector (a principal direction), λ is its eigenvalue (the variance along that direction). Sort by decreasing λ and you have the principal components ranked by importance.
Open in Lab
Drag points to reshape the cloud and watch the eigenvectors (arrows) and eigenvalues update live.
The demo wakes as you arrive…

Dimensionality reduction: keeping the important, discarding the noise

Once you have the principal components ranked by variance, the practical payoff arrives: keep the top kk components and discard the rest.

If the first two PCs capture 95% of the total variance, the remaining dimensions are mostly noise — discarding them barely changes the data's structure. You've gone from, say, 100 dimensions to 2, and you can now plot and see the data.

This is . The kk retained principal components form a new coordinate system — a lower-dimensional — onto which you project every data point. The projection minimizes reconstruction error (by Pearson's original criterion) and maximizes retained variance (by Hotelling's equivalent criterion). It is the best possible linear compression of the data into kk dimensions.

Open in Lab
Watch 3D data compress to 2D then 1D. The percentage bar shows how much variance is retained at each step.
The demo wakes as you arrive…

PCA step by step

Below is a complete walkthrough of PCA on a small 2D dataset. Click through the steps to see how raw data becomes principal components:

Open in Lab
Step through the full PCA pipeline on a live dataset.
The demo wakes as you arrive…

The same idea in code

PCA from scratch in NumPypython

Simplified to show the idea — not the real implementation.

import numpy as np

def pca(X, k):
    """Reduce n×d data matrix X to n×k using PCA.

    1. Center the data (subtract the mean of each column).
    2. Compute the covariance matrix.
    3. Find eigenvectors (principal directions) & eigenvalues (variances).
    4. Project onto the top-k directions.
    """
    # Step 1: center
    mean = X.mean(axis=0)
    X_centered = X - mean

    # Step 2: covariance matrix (d × d)
    cov = (X_centered.T @ X_centered) / (len(X) - 1)

    # Step 3: eigendecomposition — sort by descending eigenvalue
    eigenvalues, eigenvectors = np.linalg.eigh(cov)
    order = np.argsort(eigenvalues)[::-1]
    eigenvalues = eigenvalues[order]
    eigenvectors = eigenvectors[:, order]

    # Step 4: project onto top-k principal components
    W = eigenvectors[:, :k]          # d × k projection matrix
    X_reduced = X_centered @ W       # n × k  (the reduced data)

    # How much variance did we keep?
    explained = eigenvalues[:k].sum() / eigenvalues.sum()
    print(f"Retained {explained:.1%} of total variance with {k} components")

    return X_reduced, W, eigenvalues

# That's it. Pearson's geometric insight in 15 lines.
# scikit-learn's PCA does exactly this (via SVD for speed).

Strengths and limitations

PCA has endured for over a century because it is simple, fast, and optimal among linear methods. It has a (no iterative optimization), scales well, and its output is interpretable — each principal component is a weighted combination of the original variables, so you can ask "what drives this axis?"

But PCA is linear. It finds flat subspaces. If your data lies on a curved surface — a Swiss roll, a spiral, a — PCA will miss the structure entirely. This is why nonlinear successors like kernel PCA, t-SNE, UMAP, and deep autoencoders were invented. Each one asks Pearson's original question — "what subspace best represents these points?" — but allows the answer to be curved.

PCA also assumes that variance equals importance. In some problems, the most interesting structure is in the low-variance directions (rare events, subtle signals). PCA will discard exactly these. Despite these limits, PCA remains the default first step: it's fast, lossless to compute, and tells you how much linear structure the data contains before you reach for heavier tools.

Why it mattered

  1. 1901

    Pearson — Lines and Planes of Closest Fit

    Karl Pearson posed the geometric question and solved it with calculus. The first formal description of what is now called PCA.

  2. 1933

    Hotelling — modern PCA via eigendecomposition

    Harold Hotelling reformulated PCA using the covariance matrix and eigenvectors, providing the computational method still used today.

  3. 1960

    First PCA visualization

    Jolicoeur & Mosimann published the earliest 2D PCA embedding we know of, analyzing turtle shell measurements — the beginning of PCA as a visual exploration tool.

  4. 1998

    Kernel PCA — beyond straight lines

    Schölkopf, Smola & Müller extended PCA to nonlinear manifolds using kernel functions, capturing curved structures that linear PCA misses.

  5. 2006

    Hinton — deep autoencoders generalize PCA

    Hinton & Salakhutdinov showed deep autoencoders outperform PCA for dimensionality reduction — learning nonlinear principal surfaces instead of flat principal planes.

  6. 2008

    t-SNE — visualization beyond PCA

    Van der Maaten & Hinton introduced t-SNE for 2D visualization, preserving local neighborhood structure that PCA's global linear projection misses.

  7. 2018

    UMAP — fast nonlinear embeddings at scale

    McInnes, Healy & Melville combined topological theory with fast optimization, making nonlinear dimensionality reduction practical for millions of points.

  8. 2026

    PCA today — still the first step

    PCA remains the default preprocessing step in genomics, finance, NLP, and computer vision. Every autoencoder, every VAE, every dimensionality-reduction tool descends from Pearson's geometric question.

Pearson didn't just find a line through a cloud of points. He showed that data has an internal geometry — axes along which it naturally wants to be described — and that a simple geometric principle can discover them without supervision. Every time a machine learning pipeline compresses, visualizes, or denoises data, it is building on this century-old foundation.

CitationPearson, K.. On Lines and Planes of Closest Fit to Systems of Points in Space. Philosophical Magazine, 1901.

Terms in this paper