PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Med Image Comput Comput Assist Interv. Author manuscript; available in PMC 2010 March 31.
Published in final edited form as:
Med Image Comput Comput Assist Interv. 2007; 10(Pt 1): 186–193.
PMCID: PMC2847798
NIHMSID: NIHMS188685

Cortical Hemisphere Registration Via Large Deformation Diffeomorphic Metric Curve Mapping

Abstract

We present large deformation diffeomorphic metric curve mapping (LDDMM-Curve) for registering cortical hemispheres. We showed global cortical hemisphere matching and evaluated the mapping accuracy in five subregions of the cortex in fourteen MRI scans.

1 Introduction

Computational algorithms have been developed to compare anatomical shape and functions of the brain in subjects with and without neuropsychiatric disease using high-resolution magnetic resonance (MR) imaging datasets. The most popular approach for making such comparisons is to first spatially normalize brain structures and carry functions (e.g. cortical thickness, functional response) using atlas coordinates, and then to perform hypothesis testing in the atlas. We call this approach an extrinsic analysis since it depends on the relation between individual brains and the atlas. This is a powerful approach allowing us to study a large number of populations. Due to high variability of the brain anatomy across subjects, the spatial normalization becomes a crucial step because its accuracy directly influences statistical testing being inferred. Therefore, there has been great emphasis by groups on the study of the brain coordinates via vector mapping.

Cortical hemisphere registration driven by sulcal and gyral curves has been paid a great attention since sulci and gyri preserve the sulco-gyral pattern of the cortex. Specially, the brain is functionally partitioned by the sulci whose variability has been shown relative to clinical stages (e.g [1]). Compared to landmark points, the sulcal and gyral curves are one-dimensional higher manifolds and partially incorporate the geometry of the cortex. As comparing with image mapping based on intensity, the sulcal and gyral curves indirectly incorporate the intensity information since they lie at the boundary of the gray and white matters or gray matter and cerebrospinal fluid (CSF). And these curves also incorporate anatomical labeling information that is missing in the intensity-based image matching methods. As comparing with mappings of images and cortical surfaces, the vector mapping of the curves has less computation as well as itself has its own valuable clinical applications. However, due to high variability of brain structures across subjects, a desirable curve matching algorithm for brain registration has to be able to handle cases, for instance, the initial and ending points of two curves do not necessarily correspond to the same structures and one curve can be a part of the other.

Our group has focused on a diffeomorphic metric mapping of anatomical coordinates in the setting of Large Deformation Diffeomorphic Metric Mapping (LDDMM). The basic diffeomorphic metric mapping approach provides one-to-one, smooth with smooth inverse transformations acting on the ambient space. Thus, connected sets remain connected, disjoint sets remain disjoint, and the smoothness of features such as curves and surfaces is preserved. Specially, it places the set of anatomical shapes into a metric space for understanding anatomical structures across subjects. The LDDMM approach has been adapted to landmarks, images, and tensors [2,3,4].

We present the LDDMM-Curve mapping algorithm for registering cortical hemisphere. We show that such a mapping approach does not require exact correspondence information in the sense that paired curves on two anatomies can be represented by unequal number of points and one can be partial of the other or partially overlap with the other. We applied this approach to fourteen MRI datasets and evaluated them in five subregions of the cortex.

2 Large Deformation Diffeomorphic Metric Curve Mapping (LDDMM-Curve)

2.1 Vector-Valued Measure and Its Norm

Since a curve is a geometric object, it cannot be uniquely reconstructed based on the locations of a set of points. We consider that our curve embedded in R3 is a one-dimensional manifold in the sense that the local region of every point on the curve is equivalent to a line which can be uniquely defined by this point and the tangent vector at this location. In our curve mapping setting, we incorporate both the location information of points and their tangent vectors.

Assume curve C is discretized into a sequence of points x=(xi)i=1n. If xn = x1, then curve C is closed. We thus approximate curve C by a sequence of points and tangent vectors, C=(cx,i,τx,i)i=1n, where cx,i=xi+1+xi2 is the center of two sequent points and τx,i = xi+1xi is the tangent vector at cx,i. One can associate curve C with a specific measure μC given by sum of vector-valued Diracs, μC=i=1n1τx,iδcx,i such that μC acts on the smooth vector w in the form of (μCw)=i=1n1τx,iw. Its norm is defined as

μCW2=i=1n1τx,iδcx,iW2=i=1n1j=1n1kW(cx,i,cx,j)τx,iτx,j,
(1)

where kW is the kernel in a reproducing kernel Hilbert space W of smooth vectors and W* is the dual of W . The theoretical derivation of this norm can be found elsewhere [5]. Roughly, the more rounded the curve is, the higher energy it possesses, since flat shaped parts of the curve tend to vanish.

We assume that curves C and S are discretized by point sequences x=(xi)i=1n and y=(yj)j=1m, the sequence of transformed points z=(ϕ(xi))i=1n gives a discretization of the deformed curve ϕ · C. The closeness of ϕ · C and S can be quantified via the norm of the difference between their measures in Hilbert space W*, a matching functional D(ϕ · C, S), given by

D(ϕC,S)=i=1n1τz,iδcz,ij=1m1τy,jδcy,jW2,
(2)

which is explicitly

D(ϕC,S)=i=1n1j=1n1kW(cz,i,cz,j)τz,iτz,j2i=1n1j=1m1kW(cz,i,cy,j)τz,iτy,j+i=1m1j=1m1kW(cy,i,cy,j)τy,iτy,j.
(3)

The first and last terms are intrinsic energies of the two curves. The middle term gives penalty to mismatching between tangent vectors of S and those of ϕ(C). The explicit form of D(ϕ · C, S) suggests that it does not require equal number of points on the two curves. One can be a part of the other or partially overlapped with the other. Optimal matching between curves will be defined via minimization of a functional composed of the matching term D and optionally a regularization term to guarantee smoothness of the transformations. Next we need to propose a model for deformation maps ϕ.

2.2 Variational Formulation of the LDDMM-Curve Mapping

We assume that shapes can be generated one from the other via a flow of diffeomorphisms, solutions of ordinary differential equations ϕ.t=vt(ϕt),t[0,1] with ϕ0 = id the identity map, and associated vector fields vt, t [set membership] [0, 1]. We compute pairs C, S, such that there exists a diffeomorphism ϕ transforming one to the other ϕ1 · C = S at time t = 1. The metric distance between shapes is the length of the geodesic curves ϕt · C, t [set membership] [0, 1] through the shape space generated from C connecting to S. The metric between two shape C, S takes the form

ρ(C,S)2=infvt:ϕ.t=vt(ϕt),ϕ0=id01vtV2dtsuch thatϕ1C=S.
(4)

where vt [set membership] V , a Hilbert space of smooth vector fields with norm || · ||V. In practice, the metric ρ and the diffeomorphic correspondence ϕ = ϕ1 between the pair of curves (C, S) is calculated via a variational formulation of the “inexact matching problem”. We associate for each pair of curves (C, S), a norm-squared cost D(ϕ1 · C, S); then the variational problem requires minimization of the functional

J(vt)=infvt:ϕ.t=vt(ϕt),ϕ0=idγ01vtV2dt+p=1PD(ϕ1C(p),S(p)),
(5)

where D(ϕ1 · C(p), S(p) is given in (3) for the pth pair of curves {C(p), S(p)}. P indicates the number of paired curves considered in the mapping. For the simplicity of notation.

To ensure that the resulting solutions to (5) are diffeomorphisms, vt must belong to a space, V, of regular vector fields [6,7], for t [set membership] [0, 1] with 01vtVdt<. We model V as Hilbert space with an associated kernel function kV. Define the trajectories xi(t) := ϕt(xi) for i = 1 , . . . , n. Then the general solution to this variational problem can be written in the form of

vt(x)=i=1nkV(xi(t),x)αi(t),
(6)

where αi(t) are referred to as momentum vectors, which are analogous to the momentum in fluid mechanics. These momentum vectors can be computed from trajectories xi(t) by solving the system of linear equations

dxi(t)dt=j=1nkV(xj(t),xi(t))αj(t),i=1,,n.
(7)

2.3 Spline Interpolation of Deformation

Assume x to be a point on a cortical surface. We would like to apply the deformation on x based on v^t found from the LDDMM-Curve mapping. Starting from any arbitrary state v^t of the time-dependant vector fields, let us find the optimal vector fields vt that keep trajectories xi(t) = ϕt(xi), where xi is a point on curve C in (5), unchanged. This is a classical interpolation problem in a Hilbert space setting, which finds vector fields that minimize the V -norm and satisfy for each time t and each 1 ≤ in the constraint vt(xi(t))=x.i(t) where xi(t) are the fixed trajectories. The solution to this interpolation is expressed as a linear combination of spline vector fields involving the kernel operator kV as given in (6) and its associated trajectory can be computed by (7). The spline interpolation of deformation in space of V guarantee diffeomorphism.

3 Experiments

3.1 Subjects, MRI Acquisition, and Data Processing

Fourteen subjects were randomly selected from schizophrenia and bipolar disorder studies in the Division of Psychiatric Neuroimaging at Johns Hopkins University School of Medicine. MRI scans were acquired using 1.5 T scanner and MPRAGE sequence (repetition time = 13.40 ms, echo time = 4.6 ms, flip angle = 20°, number of acquisition = 1, matrix 256 × 256) with 1 mm3 isotropic resolution across the entire cranium.

The cortical surfaces at the boundary between the gray matter and white matter were generated by freesurfer [8]. Sulcal curves are well observed in the cortical surface. We delineated six sulcal curves and five brain outline curves on each cortical surface via a semi-automated method [9]. As shown in Figure 1, the six curves include three of them in the lateral view, the superior frontal sulcus, the central sulcus, the Sylvian fissure, and three in the medial view, the sulcus on the top of corpus callosum, the parieto-occipital fissure, and the calcarine sulcus. These six sulcal curves are selected because they are either the sulci segmenting the brain into functional distinct regions or major sulci in the brain lobe. The five brain outline curves are labeled from the inferior to the superior and from the posterior to the anterior as numbered in Figure 1(B). To define each curve, we manually selected a pair of starting and ending points and then the dynamic programming procedure automatically tracked the curve between the two points along the sulcus (for six sulcal curves) or the gyrus (for five brain outlines) [9]. During the LDDMM-Curve mapping, the eleven curves of one subject were chosen as template curves to map to target curves of other surfaces via the LDDMM-Curve algorithm. Finally, the optimal deformation was applied for transforming the cortical surfaces back to the template surface coordinates. In our experiments, kernels of kV(x, y) and kW(x, y) were radial with the form eyx2σ2 id, where id is a 3 × 3 identity matrix. σV and σW represent the kernel sizes of kV and kW, respectively. They were experimentally adjusted as σV = 100 and σW = 5.

Fig. 1
Panels (A,B) respectively show a surface in the lateral and medial views. Panel (C,D) shows another example. Six sulci and five brain outlines are labeled on each surface.

3.2 Examples

Figure 2 illustrates two cortical surfaces deformed to the template surface coordinates via the LDDMM-Curve mapping algorithm. The template surface is in yellow and the two surfaces to be deformed are in green. The top rows of panels (A,B) indicate that the green surfaces are roughly in the same orientation as the yellow template surface before the mapping. However, the superior region and occipital lobe on the green surfaces are clearly not aligned with those on the template. Moreover, the sulco-gyral patterns are mismatched between the green and yellow template surfaces. We illustrate the mapped green surfaces superimposed with the template surface in the bottom rows of panels (A,B). Clearly, the overall structure of the green surfaces are well deformed to the shape of the template.

Fig. 2
Panels (A,B) respectively show one example of mapping surfaces using curves. The template surface is in yellow, while the other surfaces are in green. The top row of each panel illustrates the template superimposed with the target before the mapping; ...

3.3 Evaluation of the Mapping

There are two ways to evaluate the spatial normalization procedure, including measurements quantifying the mismatching of two objects as a direct evaluation method [10] as well as shape classification and group analysis of functional responses, e.g. [11] as an indirect evaluation method. We validate our mapping results via the direct evaluation method. We calculate the cumulative distribution of distances between the deformed targets and the template to quantify their closeness. We call this as surface distance graph from surface T to surface S, defined as the percentage of vertices on a template surface T having the distance to a surface S less than dmm. Let vsi and vti respectively be vertices on surfaces S and T. The distance of vti to S is defined by dti=minvsiSvsivti, where || · || is the Euclidean distance in R3. The surface distance graph is the cumulative distribution of dti. Since this measurement is sensitive to the folding structure, we chose five subregions of the cortex containing major sulci to be evaluated, as shown in Figure 3(A,B).

Fig. 3
Panels (A,B) respectively show the lateral and medial views of the template surface with colored subregions of the cortex. Panels (C-G) show surface distance graphs for each subregion of the cortex with gray lines corresponding to the deformed surfaces ...

Panels (C-G) of Figure 3 shows the surface distance graphs for five subregions of the cortex that respectively contain the superior frontal sulcus, the central sulcus, the sylvian fissure, the calcarine sulcus, as well as the parieto-occipital fissure. The black curves are surface distance graphs corresponding to each individual original surface; the red curve is the average graph among them. Similarly, the gray curves are surface distance graphs corresponding to each individual deformed surface; the green curve is the average graph among them. These five panels show that the mean average surface distance graph after the mapping is shifted to the top of one before the mapping, which indicates that these five structures were driven to the place close to their corresponding locations on the template. Compared to the spread of the surface distance graphs among the original surface, the graphs become much close to each other after the mapping, which suggests the removal of anatomical variation across subjects.

4 Conclusion

We present the LDDMM-Curve mapping algorithm in the discrete case to mapping multiple pairs of curves for cortical hemisphere registration. The relevant approach has been developed for dealing with geometric objects of surfaces in the LDDMM setting [12]. Both the LDDMM curve and surface mappings represent geometric curves or surfaces as vector-valued measure whose norm serves as a matching functional in the LDDMM variational formulation. Both of the approaches are powerful tools for brain registration. However, one deals with curve configuration and the other handles surface configuration. Compared to the surface mapping approach, the curve mapping has its own advantage of less computation. And also the selection of curves is often supervised and provides additional labeling information. The choice of the mapping algorithms depends on the interest of applications. For instance, to differentiate functional activations on the two banks of a sulcus, the LDDMM-Curve mapping would be a better choice since the LDDMM-Surface mapping does not incorporate anatomical label information.

Several curve mapping algorithms exist for the brain registration [13,1], which require the intermediate spherical step of deforming the cortex to a spherical representation. However, this intermediate step itself introduces a large distortion that does not consistently appear across subjects. Such mapping induced distortion can result in significant loss of power in extraction of structural information. In addition to that, the deformation may not be directly applied to the natural image space. However, our approach does not require an intermediate spherical representation of the brain and directly works on the natural coordinates of curves, which potentially provide better matching.

Acknowledgments

The work reported here was supported by grants: NIH R01 MH064838, NIH R01 EB00975, NIH P20 M071616, NIH P41 RR15241, and NSF DMS 0456253. The authors would like to thank Dr. J. Tilak Ratnanather and Dr. Barta for providing the data.

References

1. Thompson PM, Schwartz C, Lin RT, Khan AA, Toga AW. Three–dimensional statistical analysis of sulcal variability in the human brain. J. Neurosci. 1996;16(13):4261–4274. [PubMed]
2. Joshi S, Miller MI. Landmark matching via large deformation diffeomorphisms. IEEE Trans. Image Processing. 2000;9(8):1357–1370. [PubMed]
3. Beg MF, Miller MI, Trouvé A, Younes L. Computing large deformation metric mappings via geodesic flows of diffeomorphisms. International Journal of Computer Vision. 2005;61(2):139–157.
4. Cao Y, Miller M, Winslow R, Younes L. Large deformation diffeomorphic metric mapping of vector fields. IEEE Trans. Med. Imag. 2005;24:1216–1230. [PubMed]
5. Glaunès J, Trouvé A, Younes L. In: Modeling planar shape variation via hamiltonian flows of curves. Statistics and Analysis of Shapes. Krim H, Yezzi A, editors. Birkhauser; 2006.
6. Trouvé A. An infinite dimensional group approach for physics based models. Technical report electronically. 1995. available at http://www.cis.jhu.edu.
7. Dupuis P, Grenander U, Miller MI. Variational problems on flows of diffeomorphisms for image matching. Quaterly of Applied Math. 1998;56:587–600.
8. Fischl B, Sereno MI, Dale AM. Cortical surface-based analysis II: inflation, flattening, and a surface-based coordinate system. NeuroImage. 1999;9:195–207. [PubMed]
9. Ratnanather JT, Barta PE, Honeycutt NA, Lee N, Morris NG, Dziorny AC, Hurdal MK, Pearlson GD, Miller MI. Dynamic programming generation of boundaries of local coordinatized submanifolds in the neocortex: application to the planum temporale. NeuroImage. 2003;20(1):359–377. [PubMed]
10. Vaillant M, Qiu A, Glaunès J, Miller MI. Diffeomorphic metric surface mapping in subregion of the superior temporal gyrus. NeuroImage. 2007 [PMC free article] [PubMed]
11. Miller MI, Beg MF, Ceritoglu C, Stark C. Increasing the power of functional maps of the medial temporal lobe by using large deformation diffeomorphic metric mapping. Proc. Natl. Acad. Sci. 2005;102:9685–9690. [PubMed]
12. Vaillant M, Glaunès J. Surface matching via currents. In: Christensen GE, Sonka M, editors. IPMI 2005. Vol. 3565. Springer; Heidelberg: 2005. pp. 381–392. LNCS. [PubMed]
13. Vaillant M, Davatzikos C. Hierarchical matching of cortical features for deformable brain image registration. In: Kuba A, Sámal M, Todd-Pokropek A, editors. IPMI 1999. Vol. 1613. Springer; Heidelberg: 1999. pp. 182–195. LNCS.