RESEARCH NOTES

Procrustes Analysis: Aligning Shapes

Two maps show the same places. One has been shifted, rotated, and drawn at a different scale. Their coordinates may differ greatly even though the arrangement of places has stayed the same. Procrustes analysis makes the comparison explicit: which coordinate differences may be removed, and what should be compared afterwards?

It helps separate shape from placement. A deeper question follows: are the differences being removed actually irrelevant to the research question?

There is a little dark humor in the name. In Greek mythology, Procrustes stretched or cut down travelers to force them into beds that did not fit. Knowing the story makes the algorithm’s name harder to look at innocently. Fortunately, no cutting is involved here: the permitted adjustments have to be written down as transformation constraints.

First decide which differences should disappear↗

“Align the shapes” does not specify a unique algorithm. Whether to allow scaling, or to treat a mirror image as the same shape, must be decided first.

Table 1. Transformation choices and the information retained
Allowed transformations Differences removed Information retained
Translation Overall position Distances, scale, and orientation
Translation and rotation Position and overall orientation Distances, scale, and handedness
Also uniform scaling Overall size Relative distances and angles
Also reflection Handedness Shape differences beyond mirroring

Uniform scaling multiplies every direction by the same number. It does not stretch one axis independently or move individual points freely. If orientation, size, or handedness matters to the study, removing it discards relevant information.

The version below is explicit: align a pair of matched point sets, allowing translation, uniform scaling, rotation, and reflection. Zero scale is a degenerate case and should be flagged.

With several shapes and no single fixed reference, the problem becomes generalized Procrustes analysis: aligning the configurations jointly rather than fitting each separately to a predetermined standard. This page starts with a pair of matched point sets.

The objective: bring matched points closer↗

Let X,Y∈Rn×dX,Y\in\mathbb R^{n\times d}, with one point per row. XX is the reference and YY is the set to align. Row ii must refer to the same place, landmark, or observation in both. The correspondence is an input, not something this step discovers.

The fitting objective is:

min⁡s,R,t∥sYR+1t−X∥F2subject toRTR=I,s≥0,t∈R1×d.\begin{aligned} \min_{s,R,t}\quad&\left\|sYR+\mathbf1 t-X\right\|_F^2\\ \text{subject to}\quad&R^{\mathsf T}R=I,\quad s\geq0,\quad t\in\mathbb R^{1\times d}. \end{aligned}

RR is a d×dd\times d orthogonal matrix, which may include reflection; ss is a uniform scale factor; tt is a row vector for translation; and 1\mathbf1 is an n×1n\times1 column of ones. The squared Frobenius norm sums squared residuals over all matched coordinates.

“Fitting residual” is more precise here than a generic “distance between distributions”. The inputs are matched points; the procedure does not estimate a distance between arbitrary probability distributions.

From equations to pseudocode↗

First compute the mean points (averaging across rows) xˉ,yˉ∈R1×d\bar x,\bar y\in\mathbb R^{1\times d}, then move each set’s centre to the origin:

Xc=X−1xˉ,Yc=Y−1yˉ.X_c=X-\mathbf1\bar x,\qquad Y_c=Y-\mathbf1\bar y.

Take the singular value decomposition (SVD) of the cross-product:

YcTXc=UΣVT,R=UVT.Y_c^{\mathsf T}X_c=U\Sigma V^{\mathsf T},\qquad R=UV^{\mathsf T}.

Under the reflection-allowing constraint used here, the optimal scale and the translation back to the reference coordinates are:

s=tr⁡(Σ)∥Yc∥F2,t=xˉ−syˉR,Y^=sYR+1t.\begin{aligned} s&=\frac{\operatorname{tr}(\Sigma)}{\|Y_c\|_F^2},\\ t&=\bar x-s\bar yR,\\ \widehat Y&=sYR+\mathbf1 t. \end{aligned}

The logic is straightforward: remove location, find an orthogonal alignment and a uniform scale, then put the result at the reference location. This pseudocode keeps those steps without language-specific array or plotting boilerplate.

INPUT: matched n-by-d X, Y
REQUIRE: finite; n>=2; d>=1

muX <- mean_point(X)
muY <- mean_point(Y)
Xc <- centre(X, muX)
Yc <- centre(Y, muY)
REQUIRE: both have spread

M <- transpose(Yc) * Xc
U, sigma, Vt <- SVD(M)
R <- U * Vt
energy <- squared_norm(Yc)
s <- sum(sigma) / energy
t <- muX - s * muY * R
Yhat <- s * Y * R
add t to each row of Yhat
E <- squared_norm(Yhat - X)

IF s is near zero:
    flag degenerate fit
RETURN Yhat, R, s, t, E

Here, mean_point averages the point coordinates across rows, centre subtracts that mean from each row, and squared_norm sums squared matrix entries. Point ordering must be checked before the algorithm runs. If all points coincide, there is no usable spread for estimating scale. With rank-deficient data, the optimal orientation may also be non-unique; one returned matrix is not necessarily the only explanation.

What should be checked beyond a close fit?↗

A known transformation gives a useful check: translate, rotate, and uniformly scale a point set, then confirm that the procedure recovers the matched locations. Deliberately change their relationships or scramble the correspondence and inspect whether the residual remains misleadingly small. Symmetric configurations can hide ordering mistakes, so this cannot replace checking landmark identities.

The residual also depends on reference units, the number of points, and preprocessing. Because this fit freely rescales YY towards XX, its residual should not be assumed to be a symmetric distance metric. Fix normalization and scoring conventions before comparing many pairs.

For implementation, official interfaces are useful, but their objectives differ. SciPy’s orthogonal_procrustes fits only an orthogonal transform, without automatic translation or scaling. Its returned scale is a sum of singular values, not the ss above. spatial.procrustes first centres and standardizes the inputs, then fits the transformation; its outputs are in processed coordinates. Different conventions need not produce matching return values.

From shapes to representations: a choice inside similarity↗

The same issue arises with neural representations. Activations elicited by the same inputs in two networks can be viewed as two point sets, while the names and ordering of their coordinate axes need not have matching physical meanings. Williams and colleagues frame such comparisons through shape metrics and transformation invariance, showing how the allowed transformations enter the definition of equivalence.

Representations with different dimensions first need a specified mapping into a common space; the algorithm on this page assumes equal column counts. If the correspondence is meant to generalize to new inputs, the inputs used to fit the transformation should be separated from those used to evaluate it. A close fit on the original points first answers a fitting question.

This asks more of a comparison than finding a lower residual. It requires deciding what information should survive. Increasing the flexibility of the transformation can also remove differences that deserved an explanation.

Similarity after alignment describes a relationship under chosen transformation conventions. Establishing whether two systems use the same computational mechanism requires behavioral, interventional, or other mechanism-directed evidence. Procrustes analysis can organize the coordinates of a comparison; it cannot decide what the comparison means.