Home | About | Journals | Submit | Contact Us | Français |

**|**HHS Author Manuscripts**|**PMC2979342

Formats

Article sections

- Abstract
- 1. Introduction
- 2. Diffeomorphic Shape Analysis of Sulci and Gyri
- 3. Intrinsic Statistical Model for Sulci
- 4. Cortical Surface Registration using Landmark Curves
- 5. Results
- 6. Conclusion
- References

Authors

Related links

Proc IEEE Comput Soc Conf Comput Vis Pattern Recognit. Author manuscript; available in PMC 2010 November 12.

Published in final edited form as:

Proc IEEE Comput Soc Conf Comput Vis Pattern Recognit. 2010 June 1; 13-18 June: 475–482.

doi: 10.1109/CVPR.2010.5540177PMCID: PMC2979342

NIHMSID: NIHMS245081

Laboratory of Neuro Imaging, University of California Los Angeles Los Angeles, California, 90095, USA

Shantanu H. Joshi: sjoshi/at/loni.ucla.edu; Ryan P. Cabeen: rcabeen/at/loni.ucla.edu; Anand A. Joshi: ajoshi/at/loni.ucla.edu; Roger P. Woods: rwoods/at/loni.ucla.edu; Katherine L. Narr: narr/at/loni.ucla.edu; Arthur W. Toga: toga/at/loni.ucla.edu

See other articles in PMC that cite the published article.

We present an intrinsic framework for constructing sulcal shape atlases on the human cortex. We propose the analysis of sulcal and gyral patterns by representing them by continuous open curves in ^{3}. The space of such curves, also termed as the shape manifold is equipped with a Riemannian *L*^{2} metric on the tangent space, and shows desirable properties while matching shapes of sulci. On account of the spherical nature of the shape space, geodesics between shapes can be computed analytically. Additionally, we also present an optimization approach that computes geodesics in the quotient space of shapes modulo rigid rotations and reparameterizations. We also integrate the elastic shape model into a surface registration framework for a population of 176 subjects, and show a considerable improvement in the constructed surface atlases.

The cerebral cortex in mammalian brains is largely responsible for memory, thought processes, and sensory perceptions, and plays a key role in social abilities, language development, and abstract reasoning and thought, especially in human brains. The cortex comprises of neurons and nerve fibers and appears gray in color. The boundary of this gray matter region, also referred to as the cortical surface, can be observed as a highly convoluted exterior with valleys (sulci) and ridges (gyri). Cortical surfaces among populations exhibit highly varied folding patterns and it is of great interest to neuroscientists as well as neurologists to recognize, classify, and understand them. A surface-based morphometric analysis of the cortex is shown to have wide reaching applicability for the purpose of mental disease detection, progression, as well as prediction and understanding of normal and abnormal developmental behaviors. Surface-based morphometry has three major ingredients: i) representation of the cortical surface, ii) registration and alignment of one or more surfaces for construction of atlases, and iii) statistical analysis of deformations or warps explaining the variability of surface features in a given population. Surface registration aims at determining correspondences between a pair of surfaces by aligning several homologous features with respect to each other. The correspondence can be achieved either automatically, for e.g. using both local and global features as in the case of Fischl et al. [5, 2], or in a semi-automatedmanner, using expertly delineated sulcal and gyral landmarks as in the case of Thompson et al. [14, 7]. These landmarks are usually 3D continuous space curves corresponding to the deepest regions of the valleys for sulci, and topmost regions of the ridges for the gyri. Typically it is difficult for automated algorithms to locate, and classify sulcal and gyral patterns exclusively based on local geometric features. Thus the main advantage of using explicit landmarks is the incorporation of expert anatomical knowledge that improves the consistency in matching of homologous features. This in turn potentially improves statistical power in the neighborhood of the landmarks. Additionally, increasing the number of consistent landmarks also improves the alignment accuracy, thereby allowing more control in the registration process.

Previously, landmark curves have mostly been used as boundary conditions for various cortical alignment approaches. Their shapes have largely been ignored during registration. The sulcal and gyral fissures are highly correlated with the folding patterns of the cortex. Recently, there has been a renewed interest in studying the cortical folding patterns for the purpose of increased understanding of the organizational processes leading to the development of the cortex [4, 10, 3, 12]. Furthermore, the shapes of such patterns are also linked with brain function. In this paper we focus on modeling the shape of sulci and gyri, and incorporating it during the process of cortical registration and atlas construction. In the past, various researchers have modeled the sulci and gyri using different representations. Leow et al. [9] have demonstrated multiple curve comparison strategies that include homothetic mapping of curves discretized with equally spaced points, and whole-curve matching by relaxing point-wise correspondences between pairs of curves. However, they project three-dimensional sulcal geometry onto a two-dimensional plane, and then use diffeomorphic mappings between the flat planes to find point correspondences. Tao et al. [13] represent sulci using landmark points on curves, and build a statistical model by using a Procrustes alignment of sulcal shapes. Vaillant et al. [15] represent cortical sulci by medial surfaces of cortical folds. They subsequently parameterize the surface into salient isoparametric active contours. Although, the advantage of this model is that it represents entire cortical folds, a limitation of this method is the use of unit speed parameterizations of active contours for constructing Procrustes shape averages for sulci. This implies that the average shape is only invariant to translation, rotation, and scaling. Furthermore, for the latter two approaches, the shapes are represented by finite features or landmarks, and thus are limited in characterization of the rich geometric detail that manifests in the cortical folds that give the sulci their shapes.

In this paper, unlike the previous approaches, we will represent the sulci and gyri by parameterized three-dimensional curves. Instead of analyzing the curves directly, we will deal with their functional representations that map curve instances onto a shape manifold. The main contributions of this paper are as follows:

- An inverse-consistent diffeomorphic framework for matching sulcal shapes.
- An intrinsic sulcal shape atlas based on the Riemannian metric on the shape manifold.
- Integration of sulcal curve diffeomorphisms in driving cortical surface registrations.

To our knowledge, the proposed framework of direct diffeomorphic three-dimensional sulcal curve mappings have not been used in cortical registration before. Additionally, we also show comparisons of shape averages of our approach to the standard Euclidean curve matching that is normally used in surface registration.

The remainder of this paper is organized as follows. Section 2 outlines the shape modeling scheme for sulci and gyri. It deals with the shape representation, and specifies a Riemannian metric on the tangent space of the shape manifold. Section 3 outlines the procedure for computing statistical shape averages for sulci and gyri for a given population. Finally Section 4 incorporates the sulcal shape model in cortical surface registration, followed by results and conclusion.

In this section, we describe the modeling scheme used to represent sulcal and gyral shape features. We will represent the cortical valleys (sulci), and the ridges (gyri) by open curves. However unlike previous approaches which have used landmarks for representing the sulcal and gyral features, we will use continuous functions of curves for representing shapes. For example, Fig. 1 shows the lateral, medial and the ventral views of a cortical surface along with manually highlighted sulci and gyri.

We use a vector valued map to parameterize the 3D sulci and gyri as follows. Let *β* be a 3D, arbitrarily parameterized, open curve, such that *β* : [0, 2π] → ^{3}. We represent the shape of the curve *β* by the function *q* : [0, 2π] → ^{3} as follows,

(1)

Here, *s* [0, 2*π*],
, and (·, ·)^{3} is taken to be the standard Euclidean inner product in ^{3}. The quantity *q*(*s*) is the square-root of the instantaneous speed, and the ratio
is the instantaneous direction along the curve. The original curve *β* can be recovered upto a translation, using
. In order to make the representation scale invariant, we will normalize the function *q* as
. We now denote
as the space of all unit-length, elastic curves. On account of scale invariance, the space S^{q} becomes an infinite-dimensional unit-sphere and represents all open elastic curves invariant to translation and uniform scaling. The tangent space of S^{q} is easy to define and is given as
. Here each *w _{i}* represents a tangent vector in the tangent space of S

(2)

Now given two shapes *q*_{1} and *q*_{2}, the translation and scale invariant shape distance between them is simply found by measuring the length of the geodesic connecting them on the sphere. We know that geodesics on a sphere are great circles and can be specified analytically. The geodesic on S^{q} between the two points *q*_{1}, *q*_{2} S^{q} along a unit direction *f* *T*_{q1}(S^{q}) towards *q*_{2} for a time *t* is given by

(3)

where *t* is a subscript for time. Then the distance between the two shapes *q*_{1} and *q*_{2} is given by

(4)

The quantify * _{t}* is also referred to as the velocity vector field on the geodesic path

In order to match curves elastically, in addition to translation and scaling, we consider the following reparameterizations and group actions on the curve that preserve its shape. A rigid rotation of a curve is considered as a group action by a 3 × 3 rotation matrix *O*_{3} *SO*(3) on *q*, and is defined as *O*_{3} · *q*(*s*) = *O*_{3}*q*(*s*), ∀*s* [0, 2π]. A curve traveled at arbitrary speeds is said to be reparameterized by a non-linear differentiable map *γ* (with a differentiable inverse) also referred to as a diffeomorphism. We define *D* = {*γ* : S^{1} → S^{1}} as the space of all orientation-preserving diffeomorphisms. Then the resulting variable speed parameterizations of the curve can be thought of as diffeomorphic group actions of *γ* *D* on the curve *q*. This group action is derived as follows. Let *q* be the representation of a curve *β*. Let *α* = *β*(*γ*) be a reparameterization of *β* by *γ*. Then the respective velocity vectors can be written as
. The reparameterization of *q* by *γ* is defined as a right action of the group *D* on the set *C* and written as
. Thus we are interested in constructing the shape space as a quotient space of S^{q}, modulo shape preserving transformations such as rigid rotations and reparameterizations.

Altogether, the set of curves affected by the group actions above, partition the space S^{q} into equivalent classes. We now define the elastic shape space as the quotient space *S* = S^{q}/(*SO*(3) × *D*). The problem of finding a geodesic between two shapes in *S* is same as finding the shortest path between the equivalent classes of the given pair of shapes. Since the actions of the re-parametrization groups on *C* constitute actions by isometries, this problem also amounts to minimizing the length of the geodesic path, such that

(5)

where *d* is given by the distance in Eq. 4. In order to optimize Eq.5, we recognize that for a fixed rotation *O*_{3}, the distance *d _{e}* can be obtained by finding the optimal reparameterization between

Having defined the framework for curve shape matching, we now construct a sulcal atlas of a large population of curves in the shape space. However, in order to do so, we need the notion of a shape average based on the sulcal and gyral curves. Owing to the nonlinearity of the shape space, the computation of an average shape is not straightforward. There are two well known approaches of computing statistical averages in such spaces. The extrinsic shape average is computed by projecting the elements of the shape space in the ambient linear space, where an Euclidean average is computed, and subsequently projected back to the shape space. On the other hand, the intrinsic average, also known as the Karcher mean ([1, 8]) is computed directly on the shape space, and makes use of distances and lengths that are defined strictly on the manifold. It uses the geodesics defined via the exponential map, and iteratively minimizes the average geodesic variance of the collection of shapes. The extrinsic average is simple and efficient to calculate, however it ignores the nonlinearity of the shape space, as well as depends upon the method used for embedding the manifold in the Euclidean space. As a simple example, a unit circle can be embedded inside ^{2} in several ways, and each embedding possibly leads to different values of extrinsicmeans on the circle. We will adopt the intrinsic approach by computing the Karcher mean for a given set of shapes. Algorithm 1 describes the procedure for computing the Karcher mean of a collection of shapes under the *q* representation.

Using Algorithm 1, we now construct a sulcal shape atlas for a population of 176 subjects. These subjects underwent high-resolution T1-weighted MRI scans with spatial resolution 1×1×1 *mm*^{3} on a Siemens 1.5T scanner. After preprocessing the raw data, and registering it stereotaxically to a standard atlas space, the cortical surfaces for these subjects were extracted using an automated algorithm [6]. For each of these subjects, a total of 28 landmark curves were manually traced. Figure 3 shows the original 28 landmark curves for each of the 176 subjects for the left hemisphere overlaid together. For the remainder of this paper, we will demonstrate the results on the left hemisphere of the brain to simplify things, although the same procedure can be repeated for the right hemisphere as well. Figure 3 also shows the intrinsic sulcal shape averages of the 28 landmark curves using Algorithm 1, as well as the respective extrinsic Euclidean averages for the same. While computing the extrinsic average, each curve for the same landmark type was mapped to its *q* representation, thus making it scale and translation invariant. The Euclidean average of all the *q* functions was then computed after a pairwise rotational alignment. Both the Karcher mean shapes as well as the Euclidean averages were then mapped back to the curve space in order to visualize them.

Lateral, frontal, and medial views of, Top row: 28 landmark sulci and gyri for 176 subjects, Middle row: Euclidean sulcal shape averages for each landmark, Bottom row: Karcher means for each landmark. (Best viewed in color).

It is observed that the intrinsic averages although smooth, have preserved important features along the landmarks, thus representing the average local shape geometry along the sulci and gyri. This also implies that the shape average has not only captured the salient geometric features, but has also reduced the shape variability in the population. In order to demonstrate this property, we plot the variance of the shape deformation for each landmark type as captured by the velocity vector along the geodesic path, both for Euclidean extrinsic, and Riemannian intrinsic averages. Thus for each of the 28 landmark average shapes for 176 subjects, * _{i}, i* = 1, … , 28, we plot the quantity
, where

In this section, we briefly introduce our scheme for registering cortical surfaces. There have been prominent approaches [14, 7] for cortical registration in the neuroimaging community. The method, although based on the same conceptual framework of elastic registration, provides a slightly differentmodel for computing the deformation, and there are improvements in the implementation that increase flexibility and efficiency.

A general outline of the process is as follows. The first stage is to establish homology of the curve data, which is accomplished by matching the shapes in a Riemannian shape space framework. Next, the surfaces and curves are conformally mapped to the sphere, establishing a common space where deformation will be defined. Following this, there is a rotational alignment of the surfaces and curves to account for the differences in spherical mapping orientations. Next, the spherical mean of the curves is found to define the atlas curves for the domain of the deformation. Finally, the elastic deformation of the atlas on the sphere is found, which constrained to map the atlas curves to each case’s set of curves. The surfaces are reparameterized by this elastic deformation. The resulting surfaces then have homologous coordinate systems, allowing local comparisons across the group.

We define a set of *N* surfaces, {*M*_{1}, …,*M _{N}*} where

The first step of the process is to establish a correspondence between the homologous landmark curves by computing mappings * _{ij}* : [0, 2π] [0, 2π] such that for curve

Next, the surfaces are mapped to the sphere, where the deformation will be computed. In general, a set of bijective mappings is found from each surface to the unit sphere, {*ϕ*_{1}, …, *ϕ _{N}*} where

(6)

We then optimally rotate the curves as well as the spherical maps as, * _{i}* =

Once the meshes and curves have been aligned on the sphere, the mean curves are computed. These mean curves will be used as the atlas curves in the surface warping. The Karcher mean on the sphere is found for each vertex of each curve, across the group. In this method, an initial guess is found by the normalized average of the points. For this point, the tangent space is defined by the gnomonic projection. A new mean is computed in the tangent space, and it is mapped back to the sphere and the process is repeated until the difference between the previous and the new mean is smaller than a threshold. We can express the curve atlas as the set , where is the karcher mean of .

For surface *i*, the deformation of the atlas is *ϕ _{i}* :

(7)

Computationally, the atlas mesh is defined on the sphere by tessellating the sphere with a subdivided octahedron, shown in Figure 5 [11]. This representation is advantageous for its multiscale and flattened representation. Sub-division of each triangle is performed by adding vertices at the midpoints of the edges, then adding new edges connecting them, creating four triangles from the original. This scheme is used for defining the representations of the mesh in multiscale algorithms. Furthermore, this subdivision process ensures that the mesh can be flattened to a rectangular and regularly sampled grid. This flattening can be imagined as follows. First, choose one of the vertices of the octahedron, this will be mapped to the center of the grid. Then cut the four far edges that do not contain the center vertex. These edges are duplicated and define the boundary of the grid, and the opposite vertex maps to the four corners of the grid. This is illustrated in Figure 6. This flattened representation allows for efficient interpolation, smoothing and finite differences operations on the grid [11].

The deformation is computed iteratively using finite differences with a multigrid method. The octahedral domain is particularly suitable for this because the representation lends to simple prolongation, restriction and smoothing multigrid operations. The solver accounts for the spherical topology of the domain by solving the equation in the flattening of the mesh at the vertex that is being updated, mapping neighboring vertex values relative to this flattening, as described above. Each mesh is resampled by the deformed atlas, establishing vertex homology between the meshes.

In this section, we demonstrate results of cortical surface registration with and without the incorporation of the elastic sulcal shape analysis. The data consists of cortical surfaces of 176 subjects collected and processed according to the procedure described in Section 3. The sulcal shape atlas is already constructed and shown in Fig. 3. We now compute geodesics between the average shape of the landmark, and the set of all sulci belonging to that landmark type, and reparameterize the set of sulci according to inverses of the resulting diffeomorphisms. We then follow the steps outlined in Sec. 4 in order to warp all the surfaces meshes to the atlas. Figure 7 shows the flattened representation of a surface atlas colored according to curvature of the surface. The figure shows flattened atlases for surfaces without any landmark curves, surfaces with landmark curves using Euclidean analysis, and finally surfaces with landmark curves using the diffeomorphic shape analysis. In the absence of landmark curves, the atlas shows smoothing of corresponding features and information present in the population, while in the presence of landmarks, the general landmark features reappear in the atlas. However the most interesting case is shown in Fig. 7 (c), where not only do the landmark curves manifest themselves as before, but they also exhibit rich local geometry that has been preserved due to the elastic diffeomorphisms. Additionally, Fig. 8 shows three the lateral, dorsal, and medial views of the reconstructed cortical surface from the flattened representations. It is observed that the surface with diffeomorphic shape analysis shows richer geometric detail than the typical Euclidean analysis.

Average flattened representations for atlases: (a) without landmark curves, (b) with landmark curves using Euclidean analysis, (c) with landmark curves using the diffeomorphic shape analysis.

We have presented a geometric method for shape analysis of sulcal and gyral features and demonstrated their use in cortical surface registration. Although, the method models shape averages of cortical sulci and gyri, the theory itself is quite general and is applicable to diverse surface correspondence problems that use landmarks to guide the warping. Thus surface registration accuracy is improved, not only by the inclusion of increasing number of landmarks, but largely by the proposed intrinsic, elastic shape analysis framework.

This research was partially supported by the National Institute of Health through the National Center for Research Resources (NCRR) for Center for Computational Biology (CCB), Grant U54 RR021813.

1. Bhattacharya R, Patrangenaru V. Nonparametic estimation of location and dispersion on Riemannian manifolds. Journal of Statistical Planning and Inference. 2002 Jan;108(1-2):23–35.

2. Dale AM, Fischl B, Sereno MI. Cortical surface-based analysis. I. segmentation and surface reconstruction. NeuroImage. 1999 Feb;9(2):179–94. [PubMed]

3. Duchesnay E, Cachia A, Roche A, Riviere D, Cointepas Y, Papadopoulos-Orfanos D, Zilbovicius M, Martinot J-L, Regis J, Mangin J-F. Classification based on cortical folding patterns. IEEE Tran Med Imaging. 2007 Jan;26(4):553–565. [PubMed]

4. Fischl B, Rajendran N, Busa E, Augustinack J, Hinds O, Yeo BTT, Mohlberg H, Amunts K, Zilles K. Cortical folding patterns and predicting cytoarchitecture. Cereb Cortex. 2008 Aug;18(8):1973–80. [PubMed]

5. Fischl B, Sereno MI, Dale AM. Cortical surface-based analysis. II: Inflation, flattening, and a surface-based coordinate system. NeuroImage. 1999 Feb;9(2):195–207. [PubMed]

6. Holmes C, MacDonald D, Sled J, Toga A, Evans A. Cortical peeling: CSF/grey/white matter boundaries visualized by nesting isosurfaces. Lecture notes in computer science. 1996 Jan;1131:99–104.

7. Joshi A, Shattuck D, Thompson P, Leahy R. Surface-constrained volumetric brain registration using harmonic mappings. IEEE Tran Med Imagng. 2007 Dec;26(12):1657–1669. [PubMed]

8. Le H. Locating Fréchet means with application to shape spaces. Advances in Applied Probability. 2001 Jan;33(2):324–338.

9. Leow A, Yu C, Lee S, Huang S, P H. Brain structural mapping using a novel hybrid implicit/explicit framework based on the level method. NeuroImage. 2005 Jan;24:910–927. [PubMed]

10. Lohmann G, von Cramon DY, Colchester ACF. Deep sulcal landmarks provide an organizing framework for human cortical folding. Cereb Cortex. 2008 Jan;18(6):1415–1420. [PubMed]

11. Losasso F, Hoppe H, Schaefer S, Warren J. Smooth geometry images. ACM SIGGRAPH Symposium on Geometry processing. 2003 Jan;43:138–145.

12. Mangin J, Riviere D, Cachia A, Duchesnay E, Cointepas Y, Papadopoulos-Orfanos D, Scifo P, Ochiai T, Brunelle F, Regis J. A framework to study the cortical folding patterns. NeuroImage. 2004 Jan;23:S129–S138. [PubMed]

13. Tao X, Han X, Rettmann M, Prince J. Statistical study on cortical sulci of human brains. Lecture notes in computer science. 2001 Jan;:475–487.

14. Thompson P, Toga A. A surface-based technique for warping 3-dimensional images of the brain. IEEE Tran Med Imagng. 1996 Jan;15(4):402–417. [PubMed]

15. Vaillant M, Davatzikos C. Finding parametric representations of the cortical sulci using an active contour model. Medical Image Analysis. 1997 Jan;1(4):295–315. [PubMed]

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's Canada Institute for Scientific and Technical Information in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |