Week 3: Image Registration

What is image registration and why do we do it?

Aligning two images (same scale, orientation, position) so they are comparable. Needed to make difference images/overlays, compare satellite images over time, and fuse modalities (e.g. PET + CT). Without it, everything differs just because of rotation/shift.

Fixed vs moving image?

Fixed = target/reference, never changes. Moving = the image that is transformed to match the fixed one. Either image can take either role.

Registration transforms what: the image or its points?

The POINTS (coordinates), not the image itself: P_fixed = T · P_moving. Warping then applies this to the pixel grid.

Rigid transformation: what does it include and what is preserved?

Translation + rotation only. All Euclidean distances (and angles) are preserved, so shape and size stay the same.

Affine transformation: what does it include and what is preserved?

Rigid + scaling + shearing. Straight lines stay straight and parallel lines stay parallel (rectangles can become parallelograms).

Projective transformation: what does it include and what is preserved?

Affine + tilting. Straight lines stay straight but parallel lines can converge to a vanishing point.

Nesting of transformation classes?

Rigid ⊂ Affine ⊂ Projective (each class contains the ones before it). Nonlinear is a separate, more flexible category (local deformation).

Translation matrix (2D)?

[1 0 tx ; 0 1 ty ; 0 0 1]

Scaling matrix (2D)?

[sx 0 0 ; 0 sy 0 ; 0 0 1]

Shearing matrix (2D)?

[1 shx 0 ; shy 1 0 ; 0 0 1]

Rotation matrix (2D)?

[cos α −sin α 0 ; sin α cos α 0 ; 0 0 1]. (The landmark slide uses +sin top-right: just the opposite angle convention.)

Tilt (projective) matrix? What makes it special?

[1 0 0 ; 0 1 0 ; E F 1]. It is the only one with a non-trivial bottom row (E, F), which makes parallel lines converge.

Why use homogeneous coordinates (point (x, y, 1) and a 3x3 matrix)?

A 2x2 matrix always keeps the origin fixed, so it cannot translate. Adding a constant 1 to the point lets the last column (tx, ty) be added in, so ALL transformations are one matrix multiplication and can be chained.

Why does the bottom row of an affine matrix stay [0 0 1]?

It keeps the last coordinate equal to 1 (0·x + 0·y + 1·1 = 1), so the output is again (x', y', 1) and can be fed into the next transformation.

3D version of homogeneous coordinates?

Point (x, y, z, 1) and a 4x4 matrix. Top-left 3x3 = rotation/scale/shear, last column = translation, bottom row [0 0 0 1].

Does the order of combined transformations matter?

Yes. Matrix multiplication is NOT commutative: T_rot · T_shear ≠ T_shear · T_rot. With matrices, the rightmost one is applied first.

How do you undo a transformation?

Use the inverse: p = T⁻¹ · p'. T⁻¹ can be found by Gaussian elimination: row-reduce [T | I] until it becomes [I | T⁻¹].

Inverse of a translation by (1, 3)?

A translation by (−1, −3). Example: p' = (3, 5) → p = (2, 2).

Exercise: translate p = (2, 2) by (3, 4)?

p' = (5, 6).

Exercise: rotate p = (2, 2) by α = −90° (slide convention)?

T = [0 1 0 ; −1 0 0 ; 0 0 1], so p' = (2, −2).

Exercise: scale p = (2, 2) by s = 3?

p' = (6, 6).

Exercise: shear p = (2, 2) with shx = 1.5, shy = 0?

x' = 2 + 1.5·2 = 5, so p' = (5, 2).

Flowchart: how to choose a registration method?

Highly similar appearance → intensity-based (NCC, sum of squared differences). Somewhat similar → histogram-based (mutual information). Not similar → clear intensity relation? yes → mutual information; no → landmarks (user-defined, scanner-defined, feature-based).

Intensity/histogram methods vs landmark methods: how is T found?

Intensity/histogram: OPTIMISE a similarity metric iteratively. Landmarks: SOLVE equations directly.

Standard steps of registration?

Acquire images → choose similarity/landmarks → find transformation (rigid/affine/nonrigid) → apply transformation (with interpolation, e.g. cubic).

Landmark-based registration: how many equations does one point pair give?

2 (one for x', one for y'): x' = a11·x + a12·y + a13 and y' = a21·x + a22·y + a23.

Affine in 2D: number of unknowns and minimum points?

6 unknowns (a11..a23), 2 equations per point → minimum 3 points, which must be NON-collinear (not on one line). Use more than 3 in practice (least squares).

Rigid in 2D: number of unknowns and minimum points?

3 unknowns (α, a13, a23) because a11 = cos α = a22 and a12 = sin α = −a21. 2 points are enough in principle. Non-linear, so solved with e.g. SVD (Kabsch/Procrustes).

Affine in 3D: unknowns and minimum points?

12 unknowns, 3 equations per point → at least 4 non-coplanar points.

Min. control-point pairs (MATLAB): nonreflective similarity / similarity / affine / projective?

2 / 3 / 3 / 4. Polynomial: 6 (order 2), 10 (order 3), 15 (order 4). Piecewise linear: 4. Local weighted mean: 6 (12 recommended).

Trick for solving rigid registration with 2 points?

Subtract the two points WITHIN the same image to get an arrow. Translation cancels (both points shift equally), so comparing the moving arrow with the fixed arrow gives the rotation α. Then plug α back into one point to get a13, a23 (t = p' − R·p).

Why is 'moving − fixed' for the same landmark not useful for finding the rotation?

p' − p = (R − I)·p + t still contains the translation and depends on where the point is. Nothing cancels.

Real-data rigid method (centroids + SVD)?

Subtract each point set's centroid (removes translation), find the best rotation by minimising squared error with SVD (Kabsch/Procrustes), then t = centroid_fixed − R · centroid_moving.

Automatic landmark pipeline?

Detection → description → matching. Landmarks can also be placed manually.

What makes a good landmark and why are corners good?

A point you can find again in the other image. Corners are pinned down in both directions (moving in any direction changes what you see). Flat region: can't tell where you are. Edge: can slide along it with no change.

What does FAST stand for and what is it?

Features from Accelerated Segment Test. An automatic CORNER detector.

FAST: what is the 'ring'?

For a candidate pixel p, a circle of 16 pixels around it (radius ≈ 3), numbered 1–16 on the slide.

FAST: the corner condition?

p is a corner if N consecutive ring pixels (usually 12; notes say 12–16) are ALL brighter than I(p) + t, i.e. I(N) > I(p) + t, OR ALL darker than I(p) − t, i.e. I(N) < I(p) − t.

FAST: what is t?

A threshold/margin so that small noise differences don't count as brighter or darker.

FAST: why does a straight edge fail but a corner pass?

Straight edge: about half the ring (≈8 pixels) is on each side, so no run of 12. Sharp corner/tip: most of the ring (≥12 in a row) differs from p. Flat area: whole ring is similar.

FAST: why is it 'accelerated'?

It first checks only 4 ring pixels (1, 5, 9, 13). If fewer than 3 clearly differ from p, p is rejected immediately, so most pixels need very few comparisons.

In the cameraman FAST slide, what are the red circles?

Detected corners (results), not the ring. The zoom shows two cases: p dark with brighter ring (I(N) > I(p) + t) and p bright with darker ring (I(N) < I(p) − t).

Why does registration error grow away from the landmarks? (formula)

A small angular error θ at a landmark displaces a point at distance r by d = r · θ. Further away → larger error.

How does a 1-pixel landmark error become an angle error?

θ ≈ ε / L (ε = landmark misplacement, L = distance between landmarks). Landmarks close together → large θ; spread far apart → small θ.

Best practice for placing landmarks?

Spread them widely and SURROUND the region of interest (outer edge), as many as practical. Don't cluster them inside/in the middle: error is very low locally but huge outside.

Effect of the number of landmarks?

More landmarks lower the error overall (individual placement errors average out).

Warping?

Applying a geometric transformation to an image.

Why is a 90° rotation easy but 45° a problem?

90°: pixels map exactly onto the grid, nothing lost or invented. 45°: edges cut through pixels and transformed coordinates become real numbers (e.g. (2,2) → (2.83, ≈0)); the grid has no pixel there, so you must estimate an intensity → interpolation.

Nearest-neighbour interpolation?

Take the value of the closest pixel. Fast, but jagged/noisy.

Linear interpolation formula?

f = (d1·f2 + d2·f1) / (d1 + d2). f1, f2 = known values at x1, x2; d1 = distance x1→x, d2 = distance x→x2. The closer point gets the larger weight (hence the 'crossed' weights).

Linear interpolation example: pixels at 2 and 3 with intensities 100 and 200, point at 2.8?

d1 = 0.8, d2 = 0.2 → f = (0.8·200 + 0.2·100)/1 = 180.

Bilinear interpolation: how?

Uses the 2x2 = 4 surrounding pixels. Interpolate between f11 and f12 (v1, v2) to get f(x1, y); same for f21, f22 to get f(x2, y); then interpolate between those two using h1, h2. Smoother but a bit blurry; fits straight lines.

Bicubic (cubic) interpolation?

Uses a 4x4 = 16 pixel neighbourhood and fits smooth curves. Smoothest; closest to the ideal result.

5x rotation by −18° vs one rotation by −90°: what does it show?

Each rotation interpolates again, so errors accumulate. Nearest neighbour: jagged and noisy. Bilinear: blurry. Cubic: closest to the single-step result.

Intensity-based registration loop?

Start with transform type + initial matrix → transform moving image (bilinear interpolation) → compare with fixed image using a metric → optimiser decides to stop; if not, it proposes a new transformation matrix and the loop repeats → registered image.

MSD formula and conditions?

MSD = (1/N) Σ (I_moving,i − I_fixed,i)². A dissimilarity metric to be MINIMISED. Only works for the same modality (intensities directly comparable).

When is mutual information used?

When intensities are not directly comparable (different modalities). It is MAXIMISED. (Slide calls it a 'dissimilarity metric to be maximised' but it is really a similarity measure.)

Mutual information: core idea?

Based on the joint histogram: when images are well aligned, intensity pairs form tight clusters. MI measures how much knowing one image's intensity predicts the other's. Needs only a CONSISTENT relationship, not the same intensities.

MI examples from the slide (A vs A, B, C, negative D)?

A vs A = 4.48 (sharp diagonal). A vs B (slightly misaligned) = 1.57 (blurrier). A vs C (strongly misaligned) = 0.90. A vs negative D = 4.48 (anti-diagonal, as high as a perfect match).

Mutual information formula?

I(A,B) = Σ p_AB(a,b) · log[ p_AB(a,b) / (p_A(a) · p_B(b)) ]

MI formula: what does the fraction inside the log mean?

Actual joint probability divided by the probability expected if A and B were independent (product of marginals). Ratio 1 → log 0 → no information; ratio > 1 → pair occurs more than chance → positive contribution. Each term is weighted by p_AB.

MI slide example: why does I(A,B) = I(A2,B) = 0.995?

Pairs (1,5)×5, (2,4)×2, (3,2)×2 → probabilities 5/9, 2/9, 2/9. Terms: 5/9·ln(1.8) + 2/9·ln(4.5) + 2/9·ln(4.5) ≈ 0.327 + 0.334 + 0.334 = 0.995. A2 is A with relabelled intensities but the same pair COUNTS, so MI is identical: MI only cares about the pairing pattern, not intensity values.

Optimiser: gradient descent formula?

T_(n+1) = T_n − α ∇F(T_n). α = step size, ∇ = gradient, F = cost function. You 'walk downhill' on the cost surface.

Gradient descent: effect of step size?

Optimal: reaches the minimum quickly. Too small: converges but slowly (more computation). Larger: bounces around the minimum. Too large: diverges.

Key MATLAB step-gradient-descent parameters?

MaximumStepLength (initial step, 0.0625), MinimumStepLength (convergence tolerance, 1e-5), MaximumIterations (100), RelaxationFactor (step shrink factor, 0.5), GradientMagnitudeTolerance (1e-4).

Three ways to quantify registration accuracy?

Point-to-point distance, Dice coefficient, Hausdorff distance.

Point-to-point error?

Euclidean distance ||p_f − p'_m|| between a landmark in the fixed image and the same landmark after transformation.

Dice coefficient?

Dice = 2|A ∩ B| / (|A| + |B|) for regions A, B in the two images. 1 = perfect overlap, 0 = no overlap.

Hausdorff distance?

For each point find its shortest distance to the other set (both directions: P→Q and Q→P). The Hausdorff distance is the LONGEST of these shortest distances = worst local mismatch. (Dice = overall overlap.)

Template matching?

Finding a known small image (template) inside a larger image.

Normalised cross-correlation (NCC): how does it work?

Slide template t over image f. At each position (u, v) compute the Pearson correlation between the template and the image window → γ(u, v). Best match = highest γ. Uses the window's own mean and standard deviation.

NCC is invariant to what?

Linear intensity changes: new = a·old + b with a > 0. Offset (+b: brighter/darker) and slope (a: contrast) are removed by subtracting the mean and dividing by the standard deviation. E.g. Y = 2X + 5 still gives r = 1.

NCC is NOT invariant to what?

Anything changing the spatial pattern: rotation, scaling, shearing/other geometric changes, nonlinear intensity changes (gamma, partial inversion), noise, occlusion, different pose.

What do NCC values +1, 0 and −1 mean?

+1 = perfect match (up to brightness/contrast). 0 = no linear relationship. −1 = perfect inverse relationship (photographic negative); this is the WORST match.

NCC vs mutual information for a contrast-reversed (negative) image?

NCC gives −1 (bad match). MI treats it as a perfect match (only needs a consistent relationship). So MI is better for multimodal registration.

Application example: aim and registration process?

Register mass spectrometry imaging (10 µm pixels) to a microscopy image (0.5 µm pixels). 1) Coarse: ~5 manual landmarks + affine. 2) Fine: pre-process both images (threshold/smooth H&E red channel, k-means on MSI), then automatic affine intensity-based registration with mutual information. 3) Evaluate error with laser landmarks in different places.

Application example: result of the error evaluation?

Average ≈ 8.37 px (SD 2.69), at 0.537 µm/pixel ≈ 4.5 µm ± 1.45 µm: below the 10 µm MSI pixel size. Allows mapping molecular signals onto single cells.