Inverse-Consistent Surface Mapping With Laplace-Beltrami Eigen-Features

Yonggang Shi; Jonathan Morra; Paul Thompson; Arthur Toga

doi:10.1007/978-3-642-02498-6_39

. Author manuscript; available in PMC: 2010 Nov 2.

Published in final edited form as: Inf Process Med Imaging. 2009;21:467–478. doi: 10.1007/978-3-642-02498-6_39

Inverse-Consistent Surface Mapping With Laplace-Beltrami Eigen-Features

Yonggang Shi, Jonathan Morra, Paul Thompson, Arthur Toga

PMCID: PMC2970526 NIHMSID: NIHMS242513 PMID: 19694286

Abstract

We propose in this work a novel variational method for computing maps between surfaces by combining informative geometric features and regularizing forces including inverse consistency and harmonic energy. To tackle the ambiguity in defining homologous points on smooth surfaces, we design feature functions in the data term based on the Reeb graph of the Laplace-Beltrami eigenfunctions to quantitatively describe the global geometry of elongated anatomical structures. For inverse consistency and robustness, our method computes simultaneously the forward and backward map by iteratively solving partial differential equations (PDEs) on the surfaces. In our experiments, we successfully mapped 890 hippocampal surfaces and report statistically significant maps of atrophy rates between normal controls and patients with mild cognitive impairment (MCI) and Alzheimer’s disease (AD).

1 Introduction

Surface mapping is an important problem in medical image analysis with applications in population studies[1] and the creation of 3D shape prior models[2]. While significant progresses have been made with the development of various techniques[3–9], most of them focus heavily on the regularization part of the problem. The lack of informative data terms makes it still a challenging problem to resolve the ambiguity in defining homologous points on smooth surfaces. To this end we propose in this paper a novel variational framework to compute maps directly between surfaces. By combining data terms from Laplace-Beltrami eigen-features and regularizing forces from inverse consistency and harmonic energy, our method can generate robust and high quality maps for a class of subcortical structures with clinical significance.

Many interesting techniques have been proposed in previous work to solve various surface mapping problems. A popular approach is to first map the surface onto the sphere and then solve the registration problem in this canonical domain [5–7, 9]. By viewing the surface as a subset of ℝ³, successful image registration techniques were applied to compute surface maps with possible landmark constraints [3, 4, 8]. Landmark curves detected by shape context features were used to guide the mapping of hippocampal surfaces[10]. The medial model provides a compact surface representation and is also very useful to construct maps between surfaces[11–13].

We propose in this work a variational approach to automatically compute both the forward and backward maps between two surfaces with inverse consistency. An important feature of our method is that we use the Reeb graph of the Laplace-Beltrami eigenfunctions to define intrinsic features for the data fidelity term in our energy function. These features provide quantitative descriptions of the salient geometry shared by elongated structures such as the hippocampus and significantly remove the ambiguity in finding correspondences for such surfaces. Another interesting feature of our method is that we incorporate both inverse consistency[14] and harmonic energy [15] as the regularization term to improve the quality of maps. The harmonic energy encourages the smoothness in the maps and the inverse consistency term helps reduce area distortions. In our experiments, we perform extensive tests to demonstrate the quality of the maps generated by our method. We also demonstrate the robustness of our method with its successful mapping of 890 hippocampal surfaces and report statistically significant results.

The rest of the paper is organized as follows. In section 2, we develop the variational framework that combines data terms of intrinsic features and regularization terms including both inverse consistency and harmonic energy. We propose the Laplace-Beltrami eigen-features for the modeling of elongated structures such as hippocampus in section 3. Experimental results are presented in section 4 to demonstrate the quality and robustness of our mapping algorithm. Finally conclusions and future work are discussed in section 5.

2 A Variational Framework

In this section, we propose a variational framework for the direct mapping of surfaces. Let ( Inline graphic , g) and (, h) be two Riemann surfaces with their metric tensor g and h. To obtain inverse consistency, we compute two maps jointly in our method. We denote u₁: → as the map from to , and u₂: → the map from to . We also assume there are L feature functions defined on each surface and denote $ξ_{1}^{j} : M_{1} \to R$ and $ξ_{2}^{j} : M_{2} \to R (1 \leq j \leq L)$ as the j-th feature function on Inline graphic and , respectively.

In our variational framework, the maps are computed as the minimizer of the following energy function:

E = E_{D} + E_{I C} + E_{H} .

(1)

There are three terms in the energy function and all of them are symmetric with respect to Inline graphic and . We call the first energy E_D as the data fidelity term and it is defined as:

E_{D} = \sum_{j = 1}^{L} α_{D}^{j} [\int_{M_{1}} {(ξ_{1}^{j} - ξ_{2}^{j} \circ u_{1})}^{2} d M_{1} + \int_{M_{2}} {(ξ_{2}^{j} - ξ_{1}^{j} \circ u_{2})}^{2} d M_{2}] .

(2)

This energy penalizes the mismatch of feature functions induced by the two maps, and the parameters $α_{D}^{j}$ are used to assign proper weights for different feature functions.

The other two terms are for regularization, where E_IC encourages inverse consistency and E_H is the harmonic energy and it ensures the smoothness in the maps. Let I denote the identity function such that I(x) = x. We define the energy E_IC as:

E_{I C} = α_{I C} [\int_{M_{1}} {(I - u_{2} \circ u_{1})}^{2} d M_{1} + \int_{M_{2}} {(I - u_{1} \circ u_{2})}^{2} d M_{2}]

(3)

where α_IC is a regularization parameter. For inverse consistency, this energy penalizes the difference between the identity map I and the composition of the two maps u₂ ◦ u₁ and u₁ ◦ u₂. The harmonic energy term E_H is defined as [15, 16]:

E_{H} = α_{H} [\int_{M_{1}} {| | J_{u_{1}} | |}^{2} d M_{1} + \int_{M_{2}} {| | J_{u_{2}} | |}^{2} d M_{2}]

(4)

where J_u₁ and J_u₂ are the Jacobian of the two maps defined on the surfaces, and α_H is the regularization parameter for this term.

To minimize the energy function with respect to the maps u₁ and u₂, we derive their gradient flows and compute them iteratively via the solution of these partial differential equations(PDEs) on the two surfaces Inline graphic and . For the energy E_D, the gradient flows of u₁ and u₂ are:

\begin{array}{l} \frac{\partial E_{D}}{\partial u_{1}} = - 2 \sum_{j = 1}^{L} α_{D}^{j} (ξ_{1}^{j} - ξ_{2}^{j} \circ u_{1}) \nabla_{M_{2}} ξ_{2}^{j} (u_{1}) \\ \frac{\partial E_{D}}{\partial u_{2}} = - 2 \sum_{j = 1}^{L} α_{D}^{j} (ξ_{2}^{j} - ξ_{2}^{j} \circ u_{2}) \nabla_{M_{1}} ξ_{1}^{j} (u_{2}) \end{array}

(5)

where Inline graphic and denote the intrinsic gradient on the surfaces. For the energy E_IC, the gradient flows are:

\begin{array}{l} \frac{E_{I C}}{\partial u_{1}} = - 2 α_{I C} [J_{u_{2}}^{T} (I - u_{2} \circ u_{1}) + (u_{2}^{- 1} - u_{1}) ∣ J_{u_{2}^{- 1}} ∣] \\ \frac{E_{I C}}{\partial u_{2}} = - 2 α_{I C} [J_{u_{1}}^{T} (I - u_{1} \circ u_{2}) + (u_{1}^{- 1} - u_{2}) ∣ J_{u_{1}^{- 1}} ∣] \end{array}

(6)

where $u_{1}^{- 1}$ and $u_{2}^{- 1}$ are the inverse maps of u₁ and u₂, and the determinant of their Jacobian are denoted as $∣ J_{u_{2}^{- 1}} ∣$ and $∣ J_{u_{2}^{- 1}} ∣$ .

To derive the gradients of the harmonic energy, we denote (x¹, x²) and (y¹, y²) as the local coordinates of Inline graphic and . In the local coordinates, the maps can be represented as $u_{1} = (u_{1}^{1}, u_{1}^{2})$ and $u_{2} = (u_{2}^{1}, u_{2}^{2})$ . The gradient flow of u₁ can then be expressed in the form of Einstein’s summation as [15]:

\frac{\partial E_{H}}{\partial u_{1}^{r}} = - α_{H} [Δ_{M_{1}} u_{1}^{r} + g^{α β} Γ_{2, p q}^{r} \frac{\partial u_{1}^{p}}{\partial x^{α}} \frac{\partial u_{1}^{q}}{\partial x^{β}}] (r = 1, 2)

(7)

where Inline graphic is the Laplace-Beltrami operator on , g^αβ = (g_αβ)⁻¹, and $Γ_{2, p q}^{r}$ is the Christoffel symbol on the manifold . Similarly, the gradient flow of u₂ is:

\frac{\partial E_{H}}{\partial u_{2}^{s}} = - α_{H} [Δ_{M_{2}} u_{2}^{s} + h^{p q} Γ_{1, α β}^{s} \frac{\partial u_{2}^{α}}{\partial y^{p}} \frac{\partial u_{2}^{β}}{\partial y^{q}}] (s = 1, 2)

(8)

where Inline graphic is the Laplace-Beltrami operator on , h^pq = (h_pq)⁻¹, and $Γ_{1, α β}^{s}$ is the Christoffel symbol on the manifold .

By combining the above results, we have the gradient descent flows of u₁ and u₂ to minimize the energy:

{\begin{array}{l} \frac{\partial u_{1}}{\partial t} = - \frac{\partial E_{D}}{\partial u_{1}} - \frac{\partial E_{I C}}{\partial u_{1}} - \frac{\partial E_{H}}{\partial u_{1}} \\ \frac{\partial u_{2}}{\partial t} = - \frac{\partial E_{D}}{\partial u_{2}} - \frac{\partial E_{I C}}{\partial u_{2}} - \frac{\partial E_{H}}{\partial u_{2}} \end{array}

(9)

To numerically solve these two equations on Inline graphic and , we use the approach of solving PDEs on implicit surfaces[17, 18, 16] and represent and as a signed distance function φ and ψ, respectively. With the implicit representation, we have a very simple formulation for the gradient flow of the harmonic energy. Using the signed distance function, we can express gradient operators on surfaces in terms of conventional gradient operators in ℝ³. For example, the intrinsic gradient Inline graphic f of a function f on can be expressed as Π_∇_φ ∇f, where ∇f is the gradient of f in ℝ³ and Π_∇_φ = I − ∇φ∇φ^T is the projection operator. By substituting the implicit form of gradient operators into (9), we can write the gradient descent flow of u₁ and u₂ as:

{\begin{array}{l} \frac{\partial u_{1}}{\partial t} = 2 \sum_{j = 1}^{L} α_{D}^{j} (ξ_{1}^{j} - ξ_{2}^{j} \circ u_{1}) Π_{\nabla ψ} \nabla ξ_{2}^{j} (u_{1}) + α_{I C} Π_{\nabla ψ} [J_{u_{2}}^{T} (I - u_{2} \circ u_{1}) + (u_{2}^{- 1} - u_{1}) ∣ Π_{\nabla φ} J_{u_{2}^{- 1}} ∣] + α_{H} Π_{\nabla ψ} \nabla \cdot (Π_{\nabla φ} J_{u_{1}}^{T}) \\ \frac{\partial u_{2}}{\partial t} = 2 \sum_{j = 1}^{L} α_{D}^{j} (ξ_{2}^{j} - ξ_{1}^{j} \circ u_{2}) Π_{\nabla φ} \nabla ξ_{1}^{j} (u_{2}) + α_{I C} Π_{\nabla φ} [J_{u_{1}}^{T} (I - u_{1} \circ u_{2}) + (u_{1}^{- 1} - u_{2}) ∣ Π_{\nabla ψ} J_{u_{1}^{- 1}} ∣] + α_{H} Π_{\nabla φ} \nabla \cdot (Π_{\nabla ψ} J_{u_{2}}^{T}) \end{array}

(10)

For numerical efficiency, all computations are performed on narrow bands surrounding the surfaces. More details of the computational schemes can be found in [16].

3 Laplace-Beltrami Eigen-Features

To use the above variational framework, it is important to design proper features that can capture the common geometry across surfaces. This is generally a difficult problem and the solution is application dependent. In this section, we propose two feature functions to characterize the global geometry of elongated structures with neuroanatomical significance such as hippocampus and putamen. Both features are invariant to scale differences and natural pose variations. For a surface Inline graphic , the first feature function ξ₁: → ℝ characterizes the tail-to-head trend of elongated structures, and the second feature ξ₂: → ℝ describes the lateral profile of the surface. We next develop the algorithm to compute both features from the Laplace-Beltrami eigenfunction of Inline graphic .

For a surface Inline graphic , the eigenfunctions of its Laplace-Beltrami operator is de-fined as [19, 20]:

Δ_{M} f = - λ f

(11)

The spectrum of Inline graphic is discrete and the eigenvalues can be ordered as 0 = λ₀ ≤ λ₁ ≤ λ₂ ≤ ···. For λ_i, the corresponding eigenfunction of λ_i is denoted as f_i. For the first eigenvalue λ₀ = 0, the eigenfunction f₀ is constant, so it is not useful to describe the shape. Among the rest eigenfunctions, the second eigenfunction f₁ is particularly interesting for modeling the global characteristics of elongated subcortical structures such as the hippocampus. If we view f₁ as a map from Inline graphic to ℝ, it has the following property[21]:

f_{1} = arg min_{f ⊥ f_{0}, | | f | | = 1} \int_{M} {∣ \nabla_{M} f ∣}^{2} d M .

(12)

with

λ_{1} = \int_{M} {∣ \nabla_{M} f ∣}^{2} d M .

(13)

Thus this function f₁ can be viewed as the smoothest, non-constant map from Inline graphic to ℝ in the space orthogonal to f₀. To numerically compute the eigenfunction, we represent as a triangular mesh = (, ), where is the set of vertices and is the set of triangles. By using the weak form of (11) and the finite element method, we can compute the eigenfunction by solving a generalized matrix eigenvalue problem:

Q f = λ U f

(14)

where Q and U are matrices derived with the finite element method.

To show how the eigenfunction f₁ can model the global shape of elongated structures, we construct its Reeb graph to obtain an explicit representation [22]. For a function f₁ defined on a manifold Inline graphic , its Reeb graph is defined as the quotient space with its topology defined via the equivalent relation x ≃ y if f₁(x) = f₁(y) for x, y ∈ . To numerically construct the Reeb graph, we trace a set of level contours of the eigenfunction f₁ on the triangular mesh representation of Inline graphic . To ensure the level contours distribute evenly over the entire surface, we use an adaptive sampling scheme developed in []. These level contours are used as the nodes of the Reeb graph and the connectivity of these nodes are established according to the neighboring relations of level contours. As an example, we show in Fig. 1(a) the Reeb graph of a hippocampus. By representing each node of the Reeb graph as the centroid of the level contour, we can see this graph has a chain structure and provides a compact model of the essentially one dimensional, tail-to-head trend of the hippocampus.

Fig. 1 — The Reeb graph of a hippocampus. (a) The Reeb graph of the f₁. (b) The surface patches interpolating the level contours.

Based on the Reeb graph of f₁, we define the first feature function ξ₁. Because the eigenfunction is generally a Moss function [23], its Reeb graph has a tree structure for genus zero surfaces. For hippocampus, the Reeb graph of f₁ typically has a chain structure as shown in Fig. 1. In the case there are branches in the Reeb graph, we prune the smaller branches according to the size of the associated level contour to ensure the pruned graph has a chain structure. We order the nodes on this chain with the increase of the function f₁ and denote the level contour at each node as C_i(i = 1, ···, N). Each contour is digitized into K points C_i = [C_i_,1, C_i_,2,, ···, C_i_,_K]. Because these points are obtained from the vertices of the mesh, we have the following relation

C = A V

(15)

where C = [C₁, C₂, ···, C_N]^T is the set of all the points on level contours, and the matrix A represents the linear interpolation operation that generates the level contours. To quantitatively describe the tail-to-head trend of the shape, we first define the feature function ξ₁ on the level contours as ξ₁(C_i_,_k) = −1 + 2 * i/N for points on C_i. To define the function ξ₁ on the entire mesh, we solve the following regularized linear inverse problem:

{| | ξ_{1} (C) - A ξ_{1} (V) | |}^{2} + β ξ_{1} {(V)}^{T} Q ξ_{1} (V)

(16)

where ξ₁(C) and ξ₁( Inline graphic ) are vectors of the values of ξ₁ on the level contours and the vertices of the mesh, respectively. The matrix Q is the same as in (14) and the term ξ₁()^T Qξ₁() encourages smoothness of the function ξ₁. For this least square problem, we have the solution for ξ₁() as

ξ_{1} (V) = {(A^{T} A + β Q)}^{- 1} A^{T} ξ_{1} (C) .

(17)

We define the second feature function ξ₂ also based on the eigenfunction of the Laplace-Beltrami operator to characterize the lateral profile of elongated structures. For each level contour C_i, we generate a surface patch that approximates the minimal surface of C_i. This is achieved by first building a Delaunay triangulation of C_i_,_k (k = 1, ···, K) using the software triangle [24] and then applying Laplacian smoothing to this mesh to obtain a smooth patch interpolating the interior of the boundary C_i as shown in Fig. 1(b). For this surface patch, we compute the second eigenfunction of its Laplace-Beltrami operator and denote it as $γ_{2}^{i}$ . The Reeb graph of $γ_{2}^{i}$ is then computed with N level contours D_i ordered with the increase of the function $γ_{2}^{i}$ . For each level contour D_i, we assign a value −1 + 2 * i/N to describe its lateral position on the surface. The value of the feature function ξ₂ on C_i_,_k is defined using linear interpolation from the values of neighboring level contours. Once we define the feature function ξ₂ on the level contours, we can compute its value on the vertices of the entire mesh similarly as:

ξ_{2} (V) = {(A^{T} A + β Q)}^{- 1} A^{T} ξ_{2} (C)

(18)

where ξ₂(C) and ξ₂( Inline graphic ) are the vectors of values of ξ₂ on the level contours and the vertices, respectively.

As an illustration, we show the feature functions of a hippocampus in Fig: 2 (a) and (b) with the parameter β = 10, N = 100, K = 100, which clearly illustrates the power of the eigen-features in characterizing the relative locations of points on the surfaces.

4 Experimental Results

In this section we present experimental results to demonstrate our method. We will first illustrate our algorithm on the mapping of two types of surfaces: hippocampus and putamen. After that, we apply our method to a data set of 890 hippocampal surfaces and present statistically significant results.

4.1 Hippocampus Results

In the first experiment, we tested our algorithm on nine hippocampal surfaces. We chose one of them as the atlas and computed maps between this surface and the other eight surfaces. The parameters in the energy function are $α_{D}^{1} = 100, α_{D}^{2} = 50$ , α_IC = 6, α_H = 1. To start the iterative algorithm, we use the same approach that we developed in [] to automatically find an initial map. For all the surfaces, the mapping algorithm converged in less than 1500 iterations and the computational process took less than 10 minutes on a PC of 1.6 GHz and 1.5 GB memory. The mapping results between the atlas and the eight surfaces are visualized in Fig. 3 by projecting a texture pattern onto these surfaces with the correspondences established by the computed maps. While the surfaces vary quite significantly, we can see the corresponding parts are matched correctly and this shows the robustness of our method to structural variations across population.

Fig. 3 — A visualization of the hippocampal mapping results. The texture pattern on the atlas surface in the center is projected onto the other eight surfaces.

To illustrate the importance of the Laplace-Beltrami eigen-features and inverse consistency, we turned off these terms separately and measured quantitatively their impact on the quality of maps. In our experiments, we computed the standard deviation of |J_u₁| and |J_u₂| as an indicator of the quality of maps because it measures the area distortion resulting from the maps. As shown in Fig. 4(a) and (b), smaller area distortions have been achieved for all the surfaces by incorporating both the eigen-features and inverse consistency.

Fig. 4 — The effect of eigen-features and inverse consistency on the quality of maps. (a) The standard deviation of the Jacobian with ( $α_{D}^{1} = 100, α_{D}^{2} = 50$ ) and without ( $α_{D}^{1} = 0, α_{D}^{2} = 0$ ) the eigen-features. (b) The standard deviation of the Jacobian with (*α_IC* = 6) and without (*α_IC* = 0) inverse consistency.

4.2 Putamen Results

In the second experiment, we tested our method on a group of 30 putamen surfaces. By choosing one of the surfaces as the atlas, we computed the maps from the other 29 surfaces to the atlas with the same parameters as in the first experiment. Using the correspondences across the 30 putamen surfaces that have been established by the maps, we can build a shape prior model for the putamen by performing a principal component analysis (PCA) [2]. The mean shape of the putamen is shown at the top of Fig. 5(a), where a texture pattern is also generated for visualization. The eigenvalues computed from the PCA process are plotted at the bottom of Fig. 5(a) where we can see most of the variations are captured by the first few components. As an illustration of the shape prior model, we plotted the shapes generated with the first three principle components. The texture pattern on the mean shape is also plotted on each synthesized surface in Fig. 5(b) to visualize correspondences across the shapes. From the results we can see the shape prior model is able to generate valid shapes with a quite large range (±5σ_i) of coefficients for the principal components.

Fig. 5 — The result of PCA on the 30 putamen surfaces. (a) Top: the mean shape. Bottom: the eigenvalues of PCA. (b) Shapes generated by varying the coefficients of the first three principal components.

4.3 ADNI Results

In the third experiment, we apply our mapping algorithm to a data set of 890 left hippocampal surfaces. The hippocampi were segmented automatically with the algorithm in [25] from the screening and 12-month follow-up scans of 445 subjects from the Alzheimer’s Disease Neuroimaging Initiative(ADNI)[26]. Among the subjects, there are 136 normal controls (NC), 228 patients with mild cognitive impairment(MCI), and 81 patients with Alzheimer’s disease(AD). The goal of this experiment is to test the differences of normal aging, MCI, and AD in terms of the atrophy rates in the hippocampus.

As a first step in this experiment, we computed maps from all 890 surfaces to the atlas surface used in the first experiment as shown in the center of Fig. 3, and the same parameters as in the first experiment were used for all surfaces. In order to compute such a large amount of surface maps, we implemented our method on a grid of more than 1000 CPUs and this enables us to compute all the maps in parallel.

To perform statistical analysis, we first generated a regular mesh representation of the atlas surface with 2000 vertices. Using the maps computed above, we projected the triangular mesh of the atlas onto the 890 surfaces. As a result, all surfaces were represented with the same triangulation and their vertices have one-to-one correspondences. To quantify the atrophy rate of the hippocampus locally, we define a thickness measure at each vertex of the mapped surfaces as the distance from the vertex to the Reeb graph of the first Laplace-Beltrami eigenfunction as shown in Fig. 1(a). By viewing this Reeb graph as a medial core of the hippocampus, our thickness measure is very similar to the definition in [1] except that our definition is completely intrinsic and invariant to pose variations such as translation, rotation and reflection. For the hippocampal surfaces of a subject, the atrophy rate at each corresponding vertex then equals the change of the thickness over the 12-month period as a percentage of the thickness at the screening scan.

Given the atrophy rates at corresponding vertices of all subjects, we performed two group analyses. In the first group study, we compared the NC group with the MCI group. At each vertex, a one-tailed t-test is applied to test the hypothesis that the MCI group has a higher atrophy rates than the NC group. The resulting p-value map is plotted onto the mean shape of the NC group in Fig. 6(a) and (b) in the top and bottom view. This map shows the significance level across different regions of the hippocampus. To correct for multiple comparisons, we performed 10,000 permutation tests and an overall statistically significance p-value 0.0063 was obtained. In the second group analysis, we compared the NC group with the AD group. Similarly, a one-tailed t-test was applied to test the hypothesis that the AD group has a higher atrophy rate than the NC group. The p-value map of this group analysis is plotted in Fig. 6(c) and (d). We also applied 10,000 permutation tests to correct for multiple comparisons and the overall p-value is 0.001, which shows the map is statistically significant. By comparing the maps in Fig. 6(c)(d) and Fig. 6(a)(b), we can clearly see the expansion of the regions with higher atrophy rates from the MCI to the AD group.

Fig. 6 — The p-value map of atrophy rates. (a)(b) NC vs MCI. (c)(d) NC vs AD.

5 Conclusions and Future Work

We proposed a general framework for surface mapping with both eigen-feature-based data terms and regularizing terms for inverse consistency and smoothness. We successfully demonstrated its application in computing maps between elongated structures such as hippocampus and putamen. For future work, we will design new eigen-features and apply the current framework to other anatomical structures. We will also perform validations to compare our method with previous techniques.

References

1.Thompson PM, Hayashi KM, de Zubicaray GI, Janke AL, Rose SE, Semple J, Hong MS, Herman DH, Gravano D, Doddrell DM, Toga AW. Mapping hippocampal and ventricular change in Alzheimer disease. NeuroImage. 2004;22(4):1754–1766. doi: 10.1016/j.neuroimage.2004.03.040. [DOI] [PubMed] [Google Scholar]
2.Cootes T, Taylor C, Cooper D, Graham J. Active shape models-their training and application. Computer Vision and Image Understanding. 1995;61(1):38–59. [Google Scholar]
3.Christensen GE, Rabbitt RD, Miller MI. Deformable templates using large deformation kinematics. IEEE Trans Imag Process. 1996;5(10):1435–1447. doi: 10.1109/83.536892. [DOI] [PubMed] [Google Scholar]
4.Joshi S, Miller MI. Landmark matching via large deformation diffeomorphisms. IEEE Trans Imag Process. 2000;9(8):1357–1370. doi: 10.1109/83.855431. [DOI] [PubMed] [Google Scholar]
5.Gerig G, Styner M, Jones D, Weinberger D, Lieberman J. Shape analysis of brain ventricles using SPHARM. Proc Workshop on Mathematical Methods in Biomedical Image Analysis. 2001:171–178. [Google Scholar]
6.Davies RH, Twining CJ, Allen PD, Cootes TF, Taylor CJ. Shape discrimination in the hippocampus using an MDL model. Proc IPMI. 2003:38–50. doi: 10.1007/978-3-540-45087-0_4. [DOI] [PubMed] [Google Scholar]
7.Gu X, Wang Y, Chan TF, Thompson PM, Yau ST. Genus zero surface conformal mapping and its application to brain surface mapping. IEEE Trans Med Imag. 2004;23(8):949–958. doi: 10.1109/TMI.2004.831226. [DOI] [PubMed] [Google Scholar]
8.Wang L, Miller JP, Gado MH, McKeel DW, Rothermich M, Miller MI, Morris JC, Csernansky JG. Abnormalities of hippocampal surface structure in very mild dementia of the alzheimer type. NeuroImage. 2006;30(1):52–60. doi: 10.1016/j.neuroimage.2005.09.017. [DOI] [PMC free article] [PubMed] [Google Scholar]
9.Yeo B, Sabuncu M, Vercauteren T, Ayache N, Fischl B, Golland P. Spherical demons: Fast surface registration. Proc MICCAI. 2008:745–753. doi: 10.1007/978-3-540-85988-8_89. [DOI] [PMC free article] [PubMed] [Google Scholar]
10.Shi Y, Thompson PM, de Zubicaray G, Rose SE, Tu Z, Dinov I, Toga AW. Direct mapping of hippocampal surfaces with intrinsic shape context. NeuroImage. 2007;37(3):792–807. doi: 10.1016/j.neuroimage.2007.05.016. [DOI] [PMC free article] [PubMed] [Google Scholar]
11.Pizer SM, Fritsch DS, Yushkevich PA, Johnson VE, Chaney EL. Segmentation, registration, and measurement of shape variation via image object shape. IEEE Trans Med Imag. 1999;18(10):851–865. doi: 10.1109/42.811263. [DOI] [PubMed] [Google Scholar]
12.Styner M, Gerig G, Lieberman J, Jones D, Weinberger D. Statistical shape analysis of neuroanatomical structures based on medial models. Med Image Anal. 2003;7(3):207–220. doi: 10.1016/s1361-8415(02)00110-x. [DOI] [PubMed] [Google Scholar]
13.Yushkevich PA, Zhang H, Gee JC. Continuous medial representation for anatomical structures. IEEE Trans Med Imag. 2006;25(12):1547–1564. doi: 10.1109/tmi.2006.884634. [DOI] [PubMed] [Google Scholar]
14.Christensen G, Johnson H. Consistent image registration. IEEE Trans Med Imag. 2001;20(7):568–582. doi: 10.1109/42.932742. [DOI] [PubMed] [Google Scholar]
15.Eells J, Lemaire L. Two reports on harmonic maps. World Scientific. 1995 [Google Scholar]
16.Mémoli F, Sapiro G, Osher S. Solving variational problems and partial differential equations mapping into general target manifolds. Journal of Computational Physics. 2004;195(1):263–292. [Google Scholar]
17.Osher S, Sethian J. Fronts propagation with curvature-dependent speed: algorithms based on Hamilton-Jacobi formulations. Journal of computational physics. 1988;79(1):12–49. [Google Scholar]
18.Bertalmío M, Cheng L, Osher S, Sapiro G. Variational problems and partial differential equations on implicit surfaces. Journal of Computational Physics. 2001;174(2):759–780. [Google Scholar]
19.Reuter M, Wolter F, Peinecke N. Laplace-Beltrami spectra as Shape-DNA of surfaces and solids. Computer-Aided Design. 2006;38:342–366. [Google Scholar]
20.Qiu A, Bitouk D, Miller MI. Smooth functional and structural maps on the neocortex via orthonormal bases of the Laplace-Beltrami operator. IEEE Trans Med Imag. 2006;25(10):1296–1306. doi: 10.1109/tmi.2006.882143. [DOI] [PubMed] [Google Scholar]
21.Belkin M, Niyogi P. Laplacian eigenmaps for dimensionality reduction and data representation. Neural Computation. 2003;15(6):1373–1396. [Google Scholar]
22.Reeb G. Sur les points singuliers d’une forme de Pfaff completement integrable ou d’une fonction nemérique. Comptes Rendus Acad Sciences. 1946;222:847–849. [Google Scholar]
23.Uhlenbeck K. Generic properties of eigenfunctions. Amer J of Math. 1976;98(4):1059–1078. [Google Scholar]
24.Shewchuk JR. Delaunay refinement algorithms for triangular mesh generation. Comput Geom Theory & Applications. 2002;22(1–3):21–74. [Google Scholar]
25.Morra J, Tu Z, Apostolova L, Green A, Avedissian C, Madsen S, Parikshak N, Hua N, Toga A, Jack C, Weiner M, Weiner M, Thompson P the ADNI. Validation of a fully automated 3D hippocampal segmentation method using subjects with Alzheimer’s disease, mild cognitive impairment, and elderly controls. NeuroImage. doi: 10.1016/j.neuroimage.2008.07.003. in press. [DOI] [PMC free article] [PubMed] [Google Scholar]
26.Mueller S, Weiner M, Thal L, Petersen R, Jack C, Jagust W, Trojanowski J, Toga A, Beckett L. The Alzheimer’s disease neuroimaging initiative. Clin North Am. 2005;15:869–877. xi–xii. doi: 10.1016/j.nic.2005.09.008. [DOI] [PMC free article] [PubMed] [Google Scholar]

[R1] 1.Thompson PM, Hayashi KM, de Zubicaray GI, Janke AL, Rose SE, Semple J, Hong MS, Herman DH, Gravano D, Doddrell DM, Toga AW. Mapping hippocampal and ventricular change in Alzheimer disease. NeuroImage. 2004;22(4):1754–1766. doi: 10.1016/j.neuroimage.2004.03.040. [DOI] [PubMed] [Google Scholar]

[R2] 2.Cootes T, Taylor C, Cooper D, Graham J. Active shape models-their training and application. Computer Vision and Image Understanding. 1995;61(1):38–59. [Google Scholar]

[R3] 3.Christensen GE, Rabbitt RD, Miller MI. Deformable templates using large deformation kinematics. IEEE Trans Imag Process. 1996;5(10):1435–1447. doi: 10.1109/83.536892. [DOI] [PubMed] [Google Scholar]

[R4] 4.Joshi S, Miller MI. Landmark matching via large deformation diffeomorphisms. IEEE Trans Imag Process. 2000;9(8):1357–1370. doi: 10.1109/83.855431. [DOI] [PubMed] [Google Scholar]

[R5] 5.Gerig G, Styner M, Jones D, Weinberger D, Lieberman J. Shape analysis of brain ventricles using SPHARM. Proc Workshop on Mathematical Methods in Biomedical Image Analysis. 2001:171–178. [Google Scholar]

[R6] 6.Davies RH, Twining CJ, Allen PD, Cootes TF, Taylor CJ. Shape discrimination in the hippocampus using an MDL model. Proc IPMI. 2003:38–50. doi: 10.1007/978-3-540-45087-0_4. [DOI] [PubMed] [Google Scholar]

[R7] 7.Gu X, Wang Y, Chan TF, Thompson PM, Yau ST. Genus zero surface conformal mapping and its application to brain surface mapping. IEEE Trans Med Imag. 2004;23(8):949–958. doi: 10.1109/TMI.2004.831226. [DOI] [PubMed] [Google Scholar]

[R8] 8.Wang L, Miller JP, Gado MH, McKeel DW, Rothermich M, Miller MI, Morris JC, Csernansky JG. Abnormalities of hippocampal surface structure in very mild dementia of the alzheimer type. NeuroImage. 2006;30(1):52–60. doi: 10.1016/j.neuroimage.2005.09.017. [DOI] [PMC free article] [PubMed] [Google Scholar]

[R9] 9.Yeo B, Sabuncu M, Vercauteren T, Ayache N, Fischl B, Golland P. Spherical demons: Fast surface registration. Proc MICCAI. 2008:745–753. doi: 10.1007/978-3-540-85988-8_89. [DOI] [PMC free article] [PubMed] [Google Scholar]

[R10] 10.Shi Y, Thompson PM, de Zubicaray G, Rose SE, Tu Z, Dinov I, Toga AW. Direct mapping of hippocampal surfaces with intrinsic shape context. NeuroImage. 2007;37(3):792–807. doi: 10.1016/j.neuroimage.2007.05.016. [DOI] [PMC free article] [PubMed] [Google Scholar]

[R11] 11.Pizer SM, Fritsch DS, Yushkevich PA, Johnson VE, Chaney EL. Segmentation, registration, and measurement of shape variation via image object shape. IEEE Trans Med Imag. 1999;18(10):851–865. doi: 10.1109/42.811263. [DOI] [PubMed] [Google Scholar]

[R12] 12.Styner M, Gerig G, Lieberman J, Jones D, Weinberger D. Statistical shape analysis of neuroanatomical structures based on medial models. Med Image Anal. 2003;7(3):207–220. doi: 10.1016/s1361-8415(02)00110-x. [DOI] [PubMed] [Google Scholar]

[R13] 13.Yushkevich PA, Zhang H, Gee JC. Continuous medial representation for anatomical structures. IEEE Trans Med Imag. 2006;25(12):1547–1564. doi: 10.1109/tmi.2006.884634. [DOI] [PubMed] [Google Scholar]

[R14] 14.Christensen G, Johnson H. Consistent image registration. IEEE Trans Med Imag. 2001;20(7):568–582. doi: 10.1109/42.932742. [DOI] [PubMed] [Google Scholar]

[R15] 15.Eells J, Lemaire L. Two reports on harmonic maps. World Scientific. 1995 [Google Scholar]

[R16] 16.Mémoli F, Sapiro G, Osher S. Solving variational problems and partial differential equations mapping into general target manifolds. Journal of Computational Physics. 2004;195(1):263–292. [Google Scholar]

[R17] 17.Osher S, Sethian J. Fronts propagation with curvature-dependent speed: algorithms based on Hamilton-Jacobi formulations. Journal of computational physics. 1988;79(1):12–49. [Google Scholar]

[R18] 18.Bertalmío M, Cheng L, Osher S, Sapiro G. Variational problems and partial differential equations on implicit surfaces. Journal of Computational Physics. 2001;174(2):759–780. [Google Scholar]

[R19] 19.Reuter M, Wolter F, Peinecke N. Laplace-Beltrami spectra as Shape-DNA of surfaces and solids. Computer-Aided Design. 2006;38:342–366. [Google Scholar]

[R20] 20.Qiu A, Bitouk D, Miller MI. Smooth functional and structural maps on the neocortex via orthonormal bases of the Laplace-Beltrami operator. IEEE Trans Med Imag. 2006;25(10):1296–1306. doi: 10.1109/tmi.2006.882143. [DOI] [PubMed] [Google Scholar]

[R21] 21.Belkin M, Niyogi P. Laplacian eigenmaps for dimensionality reduction and data representation. Neural Computation. 2003;15(6):1373–1396. [Google Scholar]

[R22] 22.Reeb G. Sur les points singuliers d’une forme de Pfaff completement integrable ou d’une fonction nemérique. Comptes Rendus Acad Sciences. 1946;222:847–849. [Google Scholar]

[R23] 23.Uhlenbeck K. Generic properties of eigenfunctions. Amer J of Math. 1976;98(4):1059–1078. [Google Scholar]

[R24] 24.Shewchuk JR. Delaunay refinement algorithms for triangular mesh generation. Comput Geom Theory & Applications. 2002;22(1–3):21–74. [Google Scholar]

[R25] 25.Morra J, Tu Z, Apostolova L, Green A, Avedissian C, Madsen S, Parikshak N, Hua N, Toga A, Jack C, Weiner M, Weiner M, Thompson P the ADNI. Validation of a fully automated 3D hippocampal segmentation method using subjects with Alzheimer’s disease, mild cognitive impairment, and elderly controls. NeuroImage. doi: 10.1016/j.neuroimage.2008.07.003. in press. [DOI] [PMC free article] [PubMed] [Google Scholar]

[R26] 26.Mueller S, Weiner M, Thal L, Petersen R, Jack C, Jagust W, Trojanowski J, Toga A, Beckett L. The Alzheimer’s disease neuroimaging initiative. Clin North Am. 2005;15:869–877. xi–xii. doi: 10.1016/j.nic.2005.09.008. [DOI] [PMC free article] [PubMed] [Google Scholar]

PERMALINK

Inverse-Consistent Surface Mapping With Laplace-Beltrami Eigen-Features

Yonggang Shi

Jonathan Morra

Paul Thompson, PM

Arthur Toga

Abstract

1 Introduction

2 A Variational Framework

3 Laplace-Beltrami Eigen-Features

Fig. 1.

Fig. 2.