Articulated Shape Matching Using Laplacian Eigenfunctions and Unsupervised Point Registration

In plain words

Picture a 3-D person or a wooden mannequin as a cloud of small cubes (voxels). Two scans of the same figure in different poses give two clouds that look very different, yet the neighbourhood structure inside each cloud changes little, except where a limb touches another part of the body. The method turns each cloud into a graph of neighbouring cubes and uses the graph's smooth "vibration modes" (eigenvectors) as a coordinate system that stays stable when limbs move. Matching two shapes then becomes lining up those coordinates.

The hard part is deciding which modes correspond and in which order. The usual rule of sorting modes by their eigenvalue is unreliable here. The authors instead describe each mode by a histogram of its values, which does not depend on how the cubes are numbered, and pair up the modes whose histograms agree. That rough alignment then starts a soft clustering step. It refines the match point by point and sets aside cubes that have no partner as outliers. The approach is aimed at large point sets and at poses where the body's connectivity changes.

TL;DR

• Input: two voxel sets of 3-D articulated shapes (for example one body in two poses, or a mannequin and a person) §5.
• Idea: Laplacian eigenvector embeddings of local graphs are congruent up to an orthogonal map, so dense matching reduces to point registration under an orthogonal transform §2.3.
• Initialisation: eigenvector histograms pick the pairing and signs; then EM with a uniform outlier component refines the registration §3 §4.
• Result: qualitative match figures and an EM iteration count of 4 to 9 on the examples; no ground-truth error and no quantitative baseline §5 tab. 1.
• Not released: no official code was found (see Paper vs code).

Problem

Classical spectral matching of weighted graphs uses the order of eigenvalues to align the eigenspaces. The authors state that this breaks down for large data sets or noisy data §1. Embedding methods also do not extend easily to articulated shapes and are sensitive to noise, outliers and missing or deformed parts §1. The setting assumed is voxel sets with a local, not fully connected, graph, so that a self-contact does not change every pairwise distance §5.1. Under exact local isometry the problem becomes an assignment under an orthogonal transform. In practice the two shapes differ and have different point counts, so the task must be posed as a maximum sub-graph isomorphism §2.3.

Method

Graph and Laplacian. Each shape is a set of voxel centres X_i in R^3. Edges join each voxel to its local neighbours, with weights from a Gaussian of the local shortest-path distance. The normalised matrix L = D^{-1/2} W D^{-1/2} is used, with D the diagonal degree matrix §2.1 eq. 1.

W_{ij} = \exp\left(-d^2(i,j)/\nu^2\right) \text{ if } i \in \mathcal{N}(j), \quad W_{ij} = 0 \text{ otherwise}, \quad W_{ii} = 0 eq. 1

Spectral reduction. The paper embeds each shape in K dimensions with eigenvectors of L (the eigenvalue 0 and its eigenvector are left out). Row i of the eigenvector matrix U^K gives the embedded coordinates x_i of voxel i §2.3.

Assignment under an orthogonal map. With orthogonal R^K and assignment ω, the alignment is written as the transpose of U_x^K equal to R^K times the transpose of U_y^K, with columns paired by ω §2.3 eq. 6. The paper writes the initial solution as R^K = S^K P^K, where S^K holds signs and P^K a permutation, and it allows R^K to be any K×K orthogonal matrix, including reflections §2.3.

Step 1: eigenfunction alignment. For each shape, take the leading eigenpairs (the paper starts from about 20 to 25 per shape §3). Compute a histogram per eigenvector. For each pair of eigenvectors (one from each shape), keep the smaller dissimilarity of the two sign choices. These costs fill a K×K matrix, and the Hungarian algorithm gives the best one-to-one pairing and the signs. Only the best-aligned pairs are kept §3.

Step 2: EM point registration. Observations are the embedded points of shape X. Centres are the embedded points of shape Y, mapped by R. Each centre is a Gaussian cluster with shared covariance Σ, and a uniform outlier cluster covers the volume V that contains the embedded shapes. Priors are equal across centres and π_out = 1 − M π_in §4. The E-step computes posteriors α_nm eq. 9. The M-step forms, for each centre, a posterior-weighted mean of observations eq. 11, updates R in closed form using an extension of the Arun et al. method to orthogonal matrices of any dimension §4, and updates Σ in the classical way §4. The paper states the EM objective as the negative expectation eq. 10.

Initialisation. R comes from Step 1. Σ = σ² I with σ of the order of 15 voxels. π_in is a small proportion of V §4 §5.

Key equations, as typeset in the paper:

\ln P(x_1,\ldots,x_N) = \sum_{n=1}^{N} \ln \sum_{m=1}^{M+1} \pi_m P(x_n \mid z_n = m) eq. 7

\alpha_{nm} = \frac{\exp(-d_\Sigma(x_n, R y_m))}{\sum_{i=1}^{M} \exp(-d_\Sigma(x_n, R y_i) + \kappa)}, \quad \kappa = \frac{(2\pi)^{K/2} (\det \Sigma)^{1/2} \pi_{out}}{V \pi_{in}} eq. 9

E(R \mid R^{(q)}) = \sum_{m=1}^{M} \xi_m \left( d_\Sigma(\bar{x}_m, R y_m) + \ln \det \Sigma \right) eq. 10

\bar{x}_m = \frac{\sum_{i=1}^{N} \alpha_{im} x_i}{\xi_m}, \qquad \xi_m = \sum_{i=1}^{N} \alpha_{im} eq. 11

Here d_Σ(a, b) = (a − b)^T Σ^{-1} (a − b) [footnote in §4].

Pseudocode for the two stages (the paper has no algorithm environment, so the citations point to the sections):

# Stage 1: eigenfunction alignment  [§2.1] [§2.3] [§3]
# Input:  voxel centres X (N, 3), Y (M, 3); embedding dimension K
# Output: initial orthogonal R (K, K), per-point posteriors (later)
def initial_alignment(X, Y, K):
    # Gaussian weights on local shortest-path distances  [§2.1] [eq. 1]
    W_x = local_gaussian_weights(X)     # nu and neighbourhood not stated
    W_y = local_gaussian_weights(Y)
    # Normalised affinity D^-1/2 W D^-1/2; which end of the spectrum  [§2.1]
    L_x = normalise_affinity(W_x)       # see Pitfalls for eigen convention
    L_y = normalise_affinity(W_y)
    # Leading eigenpairs per shape; the trivial pair is dropped  [§2.3] [§3]
    U_x = leading_eigenvectors(L_x, N_e)  # N_e = 20 to 25 in the example
    U_y = leading_eigenvectors(L_y, N_e)
    # Histogram signature of each eigenvector  [§3]
    # Bin count and histogram type are not stated.
    sig_x = [histogram(U_x[:, k]) for k in range(N_e)]
    sig_y = [histogram(U_y[:, l]) for l in range(N_e)]
    # Cost of pairing x-mode k with y-mode l, cheaper sign kept  [§3]
    # Dissimilarity measure C is not stated.
    for k in range(N_e):
        for l in range(N_e):
            A[k, l] = min(dissim(sig_x[k], sig_y[l]),
                          dissim(sig_x[k], histogram(-U_y[:, l])))
    # Hungarian assignment gives pairing P and signs S  [§3]
    P, S = hungarian(A)
    # Keep the K best-scoring pairs; the rule is not stated  [§3]
    R = S[:K, None] * P[:K, :K]         # R = S P, as an initial guess [§2.3]
    return R

# Stage 2: EM robust point registration  [§4] [eq. 9] [eq. 11]
# Input:  observations X (N, K); centres Y (M, K); init R, Sigma
# Output: refined R (K, K), Sigma, posteriors alpha (N, M)
def register_em(X, Y, R, Sigma, pi_in, V):
    # Outlier prior from the equal-prior constraint  [§4]
    pi_out = 1 - M * pi_in
    for step in range(n_iter):          # stop rule is not stated  [§5]
        # Outlier term kappa; it depends on det(Sigma)  [§4]
        kappa = ((2 * pi) ** (K / 2) * sqrt(det(Sigma))
                 * pi_out / (V * pi_in))
        # E-step: posterior of each point for each centre  [§4] [eq. 9]
        # The printed denominator and the d_Sigma factor are ambiguous;
        # see Pitfalls.
        alpha = posteriors(X, R @ Y.T, Sigma, kappa)
        # Posterior-weighted mean observation per centre  [§4] [eq. 11]
        xi = alpha.sum(axis=0)
        x_bar = (alpha.T @ X) / xi[:, None]
        # M-step: orthogonal R, closed form not written out  [§4]
        R = update_orthogonal(x_bar, Y, xi, Sigma)
        # M-step: covariance, "classically" estimated (not written out)  [§4]
        Sigma = update_covariance(X, Y, R, alpha, x_bar, xi)
    return R, Sigma, alpha

Hyperparameters

Results

Author-reported. Setup: three sets of poses of 3-D articulated objects (a wooden mannequin, a person, and a hand), captured with six calibrated cameras and a voxel reconstruction package. Hardware is not stated. The paper reports figures and a summary table, not error metrics against ground truth §5.

The paper's qualitative comparison is against using eigenvalue order alone. The authors write that EM would then "easily fall into a local minimum" in the mannequin case, but they show no numbers for it §5. They also report that most points are matched within a few iterations §5 and that the whole procedure takes "only a few seconds" for large voxel sets §5.1. No dataset release is stated.

Paper vs code

No official code.

Search, in order: (1) the arXiv abstract and metadata, and the TeX source, which have no code or project-page link; (2) GitHub repository search through the API for "articulated shape matching" and "robust point registration EM eigenfunction", both with zero results; a user search for Radu Horaud, zero results; the perception-inrialpes organisation page, not found. A general web search was not available in this run. No repository is cited.

Pitfalls

• Eigen convention conflicts with the text. The text says to drop eigenvalue 0 and use the K smallest eigenvalues §2.3. The matrix as written, D^{-1/2} W D^{-1/2}, has its trivial eigenvalue at 1, the largest value, so "smallest" and "leave out 0" do not match the matrix defined in §2.1. This is our reading of the formula; the paper does not resolve it. An implementer must choose either the leading eigenvectors of this affinity or the smallest eigenvectors of I − this affinity, then check that the trivial eigenvector is removed.
• Sign and eigenvector normalisation. Histograms of u and −u differ, so the sign choice drives the pairing §3. Normalisation of eigenvectors (unit norm or other) is not stated, and histograms depend on it.
• Histogram distance and bins are unstated. The dissimilarity C and the number of bins are not given §3. The result of Step 1 depends on both.
• Pair selection rule is unstated. The paper keeps "only the pairs of eigenfunctions with the best alignment score" §3, but it gives no threshold or count. The text uses K = 5 for the mannequin, while the Fig. 1 caption uses 7 for a different pair §5 fig. 1.
• Printed posterior is ambiguous. In eq. 9 the outlier term κ appears inside each of the M terms of the denominator, and there is no separate outlier term there. The mixture defined in §4 suggests that the outlier term belongs once in the normaliser. We read the printed form as a typesetting problem and recommend deriving the posterior from the mixture. The paper does not discuss it.
• Mahalanobis factor. The footnote defines d_Σ without the usual ½ factor §4. A Gaussian density has exp(−½ d_Σ), and the κ given in §4 is consistent with that normal density. The posterior in eq. 9 uses exp(−d_Σ), so the two are inconsistent by a factor of 2 in the exponent. Pick one convention and use it everywhere.
• Units of σ. The initial covariance is given as about 15 voxels §5, but the points are eigenvector coordinates, not voxel positions. The paper does not state how the scales relate, so the initial σ can only be set relative to the embedding scale. Check this against the embedding spread.
• Closed-form R is not written. The paper says R is found in closed form by extending Arun et al. to orthogonal matrices of arbitrary dimension §4. The paper allows R to include reflections §2.3 and states no determinant constraint, and it does not show how the weights ξ_m or the full covariance Σ enter the update. An implementer has to derive the update and decide how to handle reflections.
• Σ update is not written. The covariance is "easily estimated... as is classically done" §4. An implementer has to derive it from eq. 10.
• κ changes each step. κ depends on det Σ, so it must be recomputed in every E-step §4.
• Per-point labels are not defined. EM outputs one weighted mean per centre eq. 11. The paper says the posteriors classify points into inliers and outliers §4, but no decision rule (argmax, threshold) is stated. The one-to-one correspondence is taken from the posteriors, not given explicitly.
• Outlier volume V and π_in are unstated. Both set the outlier term. The paper gives no value, and does not say how V is measured from the embedded points.
• Stop criterion is unstated. Only iteration counts are reported tab. 1 §5.
• Proposition is not proven in the text. The equivalence between graph isomorphism and orthogonal point registration is stated as a proposition §2.3. The active text gives no formal proof, and the argument depends on a local-isomorphism assumption.
• Cost and hidden steps. The paper uses sparse Laplacians with ARPACK for the eigendecomposition §5.1. Histogram costs for 20 to 25 eigenvectors per shape, times two signs, are not discussed. Only "a few seconds" is reported, without hardware.
• No quantitative evaluation. There is no ground-truth error, no baseline, and no ablation of the outlier class or of the eigenfunction step §5. The data sets were captured by the authors and no release is stated.

Lineage

• Builds on: Umeyama's spectral graph matching §1 §2.2, which needs ordered distinct eigenvalues. Laplacian eigenmaps (Belkin and Niyogi) for the embedding and its local geometry §1 §2.3. Probabilistic point matching with Gaussian mixtures and robust weights §1. Arun et al. for the closed-form alignment, extended here §4.
• Compared against: eigenvalue-ordered alignment, used qualitatively in §5. Iterative closest point is cited as used with Laplacian embeddings of meshes by Jain and Zhang §1.
• Superseded by: not known from content in this wiki.