Search tips
Search criteria 


Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
Neuroimage. Author manuscript; available in PMC 2010 May 17.
Published in final edited form as:
PMCID: PMC2871383

Morphological Appearance Manifolds in Computational Anatomy: Groupwise Registration and Morphological Analysis


Existing approaches to computational anatomy assume that a perfectly conforming diffeomorphism applied to an anatomy of interest captures its morphological characteristics relative to a template. However, the amount of biological variability in a groupwise analysis renders this task practically impossible, due to the nonexistence of a single template that matches all anatomies in an ensemble, even if such a template is constructed by group averaging procedures. Consequently, anatomical characteristics not captured by the transformation, and which are left out in the residual image, are lost permanently from subsequent analysis, if only properties of the transformation are examined.

This paper extends our recent work [37] on characterizing subtle morphological variations via a lossless morphological descriptor that takes the residual into account along with the transformation. Since there are infinitely many [transformation, residual] pairs that reconstruct a given anatomy, we treat them as members of an Anatomical Equivalence Class (AEC), thereby forming a manifold embedded in the space spanned by [transformation,residual]. This paper develops a unique and optimal representation of each anatomy that is estimated from the corresponding AECs by solving a global optimization problem. This effectively determines the optimal template and transformation parameters for each individual anatomy, and eliminates respective confounding variation in the data. It, therefore, constitutes the second novelty, in that it represents a group-wise optimal registration strategy that individually adjusts the template and the smoothness of the transformation according to each anatomy. Experimental results support the superiority of our morphological analysis framework over conventional analysis, and demonstrate better diagnostic accuracy.

Keywords: Computational anatomy, groupwise registration, deformation based morphometry, voxel-based morphometry, morphological appearance manifolds

1. Introduction

Computational anatomy typically characterizes anatomical differences between a subject S and a template T by mapping the template space ΩT to the subject space ΩS through a diffeomorphic shape transformation h : ΩT → ΩS, x [mapsto] h(x). The resulting transformation h carries information about morphological differences between a subject and a template. Various descriptors may then be derived from h for quantifying these morphological characteristics [22, 3, 2, 47, 21, 19, 16, 34].

Such approaches date back to D’Arcy Thompson [46], who, in 1917, studied differences among species by measuring deformations of the coordinate grids from images of one species to those of another. [27] later capitalized on this idea to propose the deformable template based morphometrics, which were further extended by the landmark-based approach [12]. Later work [17, 38] provided a foundation for numerous approaches [28, 40, 39, 26, 36] based on the diffeomorphic transformations.

Deformation based morphometry (DBM), for instance, establishes group differences based on the local deformation, and displacement, [3, 25, 30, 20, 15], or the divergence of the displacement of various anatomical structures [45]. A feature of particular interest in this case is the Jacobian determinant (JD), which identifies regional volumetric changes [22, 1, 19]. Tensor based morphometry (TBM) [47, 34, 43] utilizes the tensor information for capturing local displacement. Variants of TBM include voxel compression mapping [24], which describes tissue loss rates over time.

Another class of methods known as voxel-based morphometry (VBM) [2, 21, 6, 13, 16] factors out global differences via a relatively low dimensional transformation, before analyzing them for anatomical differences. VBM is, therefore, considered as complementary to DBM or TBM, since the former utilizes the information not represented by the transformation. RAVENS maps [21] addressed some of the limitations of both approaches, in a mass preserving framework, by combining the corresponding complimentary features. An important characteristic of the resulting maps was their ability to retain the total amount of a tissue under a shape transformation in any arbitrary region by accordingly increasing or decreasing the tissue density.

Irrespective of the features employed in various types of analyses, their accuracy largely depends on the ability of finding a spatial transformation that allows perfect warping of anatomical structures. Such transformations are typically driven by image similarity measures, in the so-called intensity-driven methods [17, 55, 54, 1], either by directly employing intensity differences or via mutual information. Topological relationships among anatomical structures are maintained by imposing smoothness constraints either via physical models [17, 8] or directly on the deformation field. Although these methods have provided remarkable results, their ability to establish anatomically meaningful correspondences is not always certain, since image similarity does not necessarily imply anatomical correspondence. Alternative methods [48, 53] are based on feature identification and matching. Hybrid methods, such as [42], incorporate geometrically significant features to identify anatomically more consistent transformations.

Despite remarkable success of the above methods, their fundamental shortcomings result from inherent complexity of the problem. First, anatomical correspondence may not be uniquely determined from intensity-based image attributes, which drive these algorithms. Second, exact anatomical correspondence may not exist at all due to anatomical variability across subjects. For instance, it may not be possible to perfectly warp a single-folded cingulate sulcus to a double-folded cingulate sulcus through a biologically meaningful transformation. As a result, the choice of a template plays an important role in the accuracy of the analysis. Anatomies closer to the template are well represented by a diffeomorphism. However, large differences between an individual and the template lead to the residual information that the transformation does not capture. In other words, a priori fixing a single template for all subjects biases the analysis [13, 11]; this has been a fundamental assumption in computational anatomy. Although one may argue that this residual may be eliminated by a very aggressive registration, this may create biologically implausible correspondences, as illustrated in Fig. 1. Aggressive warping may also lead to noisy deformation fields that are unsuitable for subsequent statistical analysis. To avoid these situations, typically some level of regularity is always desirable in the transformation. Both parameters, i.e., the choice of a template as well as the level of regularity, are often arbitrarily or empirically chosen. More importantly, they are not optimized for each anatomy; individually this leads to unwanted confounding variation in the morphological measurements, which may reduce our ability to detect subtle abnormalities.

Figure 1
Effect of aggressive registration: (a) Subject; (b) Template; (c) Subject warped through a viscous fluid algorithm. Aggressively registering a bifolded sulcus to a single-folded sulcus creates very thin needle like structures, which are atypical of a ...

Fig. 2 shows a characteristic example of residuals that persist after registration of brain images. MRI scans of 31 human brains were spatially normalized to the template given in Fig. 2(a) using deformable registration of [42], with the smoothness parameter adjusted to avoid the creation of biologically incorrect warpings. The average of spatially normalized brains is given in Fig. 2(b), whose clarity indicates relatively good registration. A typical subject (Fig. 2(c)), on the other hand, when warped to the template (Fig. 2(d)) still exhibits significant residual as shown in Fig. 2(e). Depending on the template, this residual may be large enough to offset disease specific atrophy in a brain, and may easily confound subtle anatomical variations. The question under consideration is, therefore, how to carry out accurate statistical analysis in such situations.

Figure 2
(a) Template; (b) Mean of 31 spatially normalized brains; (c) A representative subject; (d) Spatial normalization of (b) using HAMMER; (e) Corresponding residual. While crispness of the mean brain indicates reasonably good anatomical correspondence, there ...

Some approaches have recently been proposed that utilize the mean brain as a template [23, 4]. In most practical cases, considerable differences still persist between some samples and the mean brain, and consequently, the residual is never negligible. A very promising approach in this situation is the groupwise registration [9, 50, 10], which solves the problem to a certain extent, in the sense that instead of minimizing individual dissimilarity measures it minimizes combined cost. Bhatia et al. [9], for instance, implicitly find a common coordinate system by constraining the sum of all deformations from itself to each subject to be zero. Joshi et al. [31] compute the most representative template image through a combined cost functional, by defining a metric on the diffeomorphism group of spatial deformations. Such groupwise registration based representations are, therefore, more consistent across the samples, even though they still fail to eliminate residual.

A very elegant approach that addresses residual image differences after a diffeomorphism is the so-called metamorphosis [49], which constructs paths that simultaneously spatially transform and change the intensity of the template, so that it matches a target. In this paper, we restrict our attention to a specific subset of such paths, i.e., those which reconstruct a given anatomy, and we attempt to learn the respective submanifold by varying parameters of the transformation and the template. We start with a framework that constructs a complete morphological descriptor1 (CMD) of the form [Transformation, Residual] for a given anatomy. Instead of fixing registration parameters arbitrarily or in some ad hoc manner, we build an entire class of such anatomically equivalent CMDs by varying templates and regularization levels, which we refer to as Anatomical Equivalence Class (AEC). This leads to a manifold representation of each anatomy that reflects variations due to the choice of the template and parameters. In other words, a given anatomy is not represented by a single measurement, but by infinitely many different, yet in some sense equivalent, representations. We find an optimal morphological signature (OMS) for each individual anatomy by allowing each anatomy to traverse its own manifold; the optimal positions on their respective manifolds are determined from all anatomies simultaneously, by optimizing certain criteria. Since this optimization is applied jointly to all images in the group, it effectively defines a group-wise registration procedure. Our results indicate that this morphological representation provides generally better and more robust detection of group morphological differences than the transformation alone. Moreover, it is largely invariant to template selection, and the amount of regularization in the transformation, thus avoiding their arbitrary selection.

The idea of residual based analysis was first proposed in [37], where individual anatomies were mapped to a global PCA subspace, and projections on the orthogonal subspace were used to remove confounding variations that result merely from varying the transformation parameters. A major limitation of this approach was its inability to construct a manifold from each individual, by jointly fitting a single hyperplane to all individuals together. The current paper builds on the work presented in [5], where we approximate each individual AEC by a hyperplane. In addition, [37] considered only the regularization parameter, whereas similar to [5] in this paper, we combine all confounding factors, including the effect of the choice of the template in a general framework that leads to optimal representation for morphological analysis.

The concept of an appearance manifold has gained a great deal of attention in the past 5 years in the computer vision community [33, 18, 41, 35], albeit in a different context. In particular, object recognition depends largely on lighting and pose, which are variable and are considered confounding parameters. Image appearance manifolds are often constructed by varying the lighting and pose parameters, and learning the resulting image variations. Morphological appearance manifolds herein are constructed in an analogous way: parameters such as regularization constants and templates are varied for each anatomy, thereby allowing us to construct an appearance manifold of that anatomy, and to factor out variations that do not reflect underlying biological characteristics.

The paper is organized as follows. We start with problem formulation and morphometric analysis framework given next. Proposed optimal descriptor will be presented in the later part of Section 2. Section 3 presents and discusses the experimental results on synthetic 2D and real 3D volumetric datasets with simulated atrophy. We conclude in Section 4 with a discussion of results and future directions.

2. Morphological Descriptor Framework

2.1. Motivation

The fundamental principle of computation anatomy is that differences between various individuals, Si, i = 1,…,L are characterized relative to a template T. This requires mapping the template space ΩT to individual subject space ΩSi through a diffeomorphism hi [set membership] HSi : ΩT → ΩSi, x [mapsto] h(x), where HSi is the set of all optimal diffeomorphic transformations that maximize some similarity criterion between T and normalized subject SiT. The resulting transformations hi, i = 1,…,L then carry morphological information of different individuals relative to the template [22, 38, 3, 30, 2, 47, 21, 19, 16, 34]. Ideally, the warped subject and the template should be identical, i.e., T(x) − Si(hi(x)) = 0, ∀x [set membership] ΩT. More importantly, hi must define anatomical correspondences for all individuals in the ensemble. Simultaneously achieving both goals is usually not possible, and residual error follows:


We, therefore, view transformation jointly with the residual for a complete representation, hii:=(hi,Rhi), which we refer to as CMD. Note that the subscript of x2133 denotes the transformation with superscript representing the individual. As a result, no information is lost, even if the transformation fails to capture some shape characteristics of the subject.

2.2. Anatomical Equivalence Class Framework

The CMD, hS, depends not only on the underlying anatomy S but also on transformation parameters, which collectively are denoted by a vector θ. An entire family of CMDs may be generated by varying θ [set membership] Dθ, where Dθ represents the domain of θ, thereby, defining a parametric manifold, should θ be a continuous variable. In short, hS is generated by combining vectorial forms of h and Rh.

Definition 1

Two CMDs (hp, Rhp) and (hq, Rhq) [set membership] X are anatomically equivalent (~) if for a given template T, they reconstruct the same anatomy S:


where X is the space spanned by all such CMDs, and y [set membership] ΩS.

This non-uniqueness of representation will be removed in Sections 2.3 and 2.4. We will discuss about specific examples of θ later; first we develop the framework in a general setting.

Assumption 1

For a given anatomy S, the class of anatomically equivalent CMDs, referred to as anatomically equivalent class (AEC), is generated by varying transformation parameters θ [set membership] Dθ:


forms a smooth manifold QS embedded in the subspace X [subset, dbl equals] Rn spanned by CMDs, n being twice the cardinality of discretized ΩT.

Note that all CMDs, each represented as a set of discretized components {hθS,xΩT} in Eq. (3), satisfy the constraint S(hθ(x)) = T(x) − Rhθ(x), ∀x [set membership] ΩT, which according to Definition 1 makes them anatomically equivalent. It is known that image appearance manifolds constructed by varying such parameters are continuous but not differentiable, if the images contain edges [52]. They can, however, be approximated by smooth manifolds, either by fitting parameterized manifolds to them, or by smoothing the images from which they are constructed [51]. Herein, we take the former approach, by approximating the manifold by a hyperplane within which it is contained. Since the dependence of hθS on diffeomorphism h is through the transformation parameters θ, in the rest of the text we drop the subscript to represent CMDs as x2133S(θ).

2.3. Anatomical Comparisons based on AEC Manifolds

Although the resulting AEC represents the entire range of variability in Rhθ, it is not always clear what θ should be selected for analysis. For instance, a small residual does not necessarily correspond to the best CMD, as we will experimentally show later, although minimizing residuals is a common target in most of deformable registration algorithms. It is, therefore, imperative to find an appropriate θ for each individual that is optimized to the underlying anatomy according to a certain criterion J, leading to optimal parameters Θ = (θ1,…,θL), where L is the number of subjects.

In order to define our criterion for optimality of Θ, we first consider the simplest case of two subjects. To find the distance between two anatomies represented by their respective manifolds, one may define J as the minimum separation between their corresponding manifolds. Physically this amounts to inter-orbital distance and is computed by moving along the manifolds such that the distance between corresponding points is minimized.


where d represents Euclidean distance defined on X, and FSA and FSB respectively denote the spaces of diffeomorphisms for anatomies A and B. While this works for two anatomies, comparing more than two anatomies becomes problematic as shown in Fig. 3, where the optimal representations for comparison between SA and SB differ from those for comparison between SB and SC.

Figure 3
Manifold structure and intersubject comparisons based on Euclidean distance.

To generalize the cost functional, we notice that for two subjects, optimization allows sliding along the respective manifolds to find two representations that yield minimum distance. These representations best highlight differences between these anatomies, since together they eliminate confounding effects of θ. For L anatomies, we allow their representations to slide along respective AECs to minimize the sum of pairwise squared distances of all individuals. The objective functional, therefore, becomes:


and the optimization is constrained to respective manifolds such that we in effect find optimal parameters as:


where x2133i(θi) is the CMD of subject i for Θ = (θ1,…,θL). Θ* represents the optimal selection of parameters, whose values provide the optimal θ for each individual, which differs, in general, across individuals.

It is trivial to show that the criterion of Eq. (6) minimizes the variance of morphological descriptors over entire ensemble with respect to confounding factors, leading to:




represents the mean descriptor. The resulting OMS, i(θi), corresponding to optimal parameters, θi, for each individual, and, therefore, removes arbitrariness due to these parameters, as illustrated schematically in Fig. 4.

Figure 4
OMS versus CMD: (a) Randomly selecting CMDs (random parameter selection) reduces inter-group separation. Dots marked with arrows represent group means; (b) OMS results in optimal separation between the two groups.

2.4. Optimal Morphological Signature

For simplicity and tractability, we approximate AEC manifolds with hyperplanes to solve the optimization problem of Eq. (7), as illustrated in Fig. 5.

Figure 5
Approximation of AEC manifolds with hyperplanes. τ and λ may be two of the confounding factors invariance to which is sought, as explained later in Section 2.5.

Each manifold is first independently represented in terms of its principal directions, computed through principal component analysis (PCA). If {Vj(i),j=1,,n} represent principal directions of subspace Xi in which the manifold QSi of subject Si is embedded, and Mi denotes corresponding subject mean, then the linear hyperplane approximating the corresponding AEC manifold is given by:


where αij[αijmin,αijmax],j=1,,n capture transformation dependent variability originally represented by θi. Basically, αij identifies the component of a complete descriptor x2133i(θ) along direction Vj(i), and by varying it in the interval [αijmin,αijmax], we allow sliding along the manifold. The bounds αijmin, and αijmax define the extents of the manifold, and may be computed from corresponding principal modes. In the experiments of this paper, these bounds were selected so as to allow moving twice the standard deviation along each principal direction.

The objective function of Eq. (7), therefore, becomes:




is the mean CMD across subjects for current correspondence αi = (αi1,…,αin).

Solution to this constrained problem is an algorithm that allows moving along individual hyperplanes, to minimize the objective function. As shown in Fig. 5, at each optimization iteration, an update of x2133i(αi) is computed, which yields the current floating mean M(A). The procedure is repeated until the minimum of Eq. (9) is attained. Analytically it leads to the following solution subject to constraints given above:


When combined with Eq. (8), optimal αi yields OMS, x2133i*, which is then used for subsequent analysis. It provides optimal combination of transformation and residual by finding optimal selection of transformation parameters θ.

2.5. Robustness to Template and Regularization Parameters

In this section, we particularize the above formulation to the problem at hand, by noting that the residual is mainly a consequence of two parameters, namely the template and the regularization of the diffeomorphism, as mentioned in Section 1. The AEC of a given anatomy, S, is generated by varying these two parameters λ [set membership] [proportional, variant]+ and τ [set membership] T with T representing the set of all possible templates. Since analysis eventually has to be carried out in a common space, we ultimately bring all warped anatomies to a common template space ΩT, as illustrated in Fig. 6. However, variations caused by the selection of different templates have already been represented by the intermediate templates τ (Fig. 6).

Figure 6
Constructing AECs: Each subjects is normalized to ΩT via intermediate templates at different smoothness levels of the warping transformations.

Suppose that for a given τ, F represents the set of all diffeomorphisms that warp2 S to τ to yield Sτ whereas GSτT denotes the set of diffeomorphisms that warp Sτ to T:


Then, the set of transformations that take S to T is given by:


Each element of F, GSτT, or x2130S, is a diffeomorphism that minimizes, in a standard deformable registration, an energy functional of the form [Data term] + λ [Smoothness term]. The regularity parameter λ, therefore, controls the level of flexibility of the transformation. For viscous fluid, for instance, it corresponds to the viscosity, whereas in elastic registration, it might represent Young’s modulus. For consistency, we reparameterize λ, such that small values correspond to more aggressive transformations, and, therefore, small residuals, in which case more information is carried by the transformation. Large λs, on the other hand, fuse more information in the residual at the cost of relatively smooth transformation.

The second parameter, τ provides a milestone between S and T to guide the registration algorithm. For a given λ, individual anatomies are first normalized to a stepping-stone template τ [set membership] T through fλ,τ [set membership] F as shown in Fig. 6. Intermediate results are then warped to T through gλ,τ [set membership] GSτT, which brings all anatomies to a single reference space, ΩT. We have found that appropriate selection of τ typically facilitates registration, and anatomies, which are otherwise not well represented by T, are captured effectively through x2130S. This is illustrated in Fig. 7, where indirect warping performs better than direct warping.

Figure 7
Intermediate templates aid registration: (a) Template of Fig 2(a); (b) Direct warping of Fig. 2(c); (c) Warping via an intermediate template.
An external file that holds a picture, illustration, etc.
Object name is nihms165523u1.jpg

2.6. Feature Extraction

The previously discussed framework was derived in a very general setting, with its applicability going beyond medical image analysis. In this paper, we leave out discussion on the application of this approach to other areas of computer vision, and instead focus on how to formulate features for application to biomedical image analysis.

So far we have used a general notation x2133h for CMDs, which is suitable for groupwise registration. However, several avenues open if one concentrates on the problem at hand. For instance, if one is interested in characterizing volumetric deficits for diagnosis of Alzheimer’s disease, Jacobian determinant (JD), Jh, of the transformation h is a feature of interest. For modeling brain development patterns, one may utilize some measure of diffusivity, such as fractional anisotropy [7]. For consistent labeling, intensity values provide a good feature of interest.

We are primarily interested in the amount of growth/shinakage of the brain tissue. The former is typically due to brain development and the latter to various neurodegenerative pathologies that cause tissue death. A natural feature of interest in this case is based on JD with the CMD defined as x2133h := (Jh, Rh). Combining Jh and Rh features in one descriptor is problematic since both represent different physical quantities. We propose two approaches to solve this problem.

2.6.1. Weighted Combination of Features

One may readily suggest weighting individual features with their z-scores. In this paper, we adopt a more natural combination that takes into account physical nature of the two quantities. It first involves the segmentation of warped anatomy, S(h(x)), into c tissues. For example, skull-stripped brain MR images are often segmented with c = 3 into grey matter (GM), white matter (WM), and cerebrospinal fluid (CSF). Pre-processed cardiac MR images are also typically segmented with c = 3 into muscle, blood, and fat. Based on segmented images one may readily extract tissues of interest. Suppose Ik(h(x)) and Kk(x) respectively denote indicator functions for k = 1,…,c tissues in the individual image S and the template T. Then corresponding residuals are defined as:


Since JD represents the degree of growth or atrophy, a similar quantity is derived by weighting tissue based residuals Rhk to have values of { −1, 0, +1}. With such weighting +1 indicates unit increase in volume, and −1 represents a unit decrease. Total tissue loss or gain as reflected by the residual may be computed by integrating the residual over the region of interest, which is tantamount to counting up or down the number of non-zero voxels (assuming that a voxel has a unit positive or negative volume), where the sign indicates volumetric growth or deficit. As a result, if the Jh shows an atrophy at x [set membership] ΩT, actual volumetric deficit at x may be lower or higher depending on the corresponding Rhk.

2.6.2. Tissue Density Maps

h and Rh may also be combined to form a morphometric descriptor based on residual and shape-based features, referred to as tissue density map (TDM) (RAVENS maps in [21]) that relates to the spatial distribution of anatomical tissues in the target anatomy:

Definition 2

A tissue density map (TDM), Dk : ΩT → ΩT, for k = 1,…,c tissues is defined as:


TDMs may be computed in a mass-preserving framework of [21] by warranting that the amount of tissues in S are preserved exactly under the transformation.

Theorem 1

TDM defined by Eq. (13) is a mass preserving map [37].


For any region R [subset or is implied by] ΩT defined in the template space,


As a consequence of this theorem, regions that contract under h−1 demonstrate an increase in density. Since a relatively larger anatomical region contracts to a smaller one, its density increases if the tissue mass is to be preserved, as illustrated in Fig. 8.

Figure 8
Schematic illustration of tissue density maps (TDMs). Left column shows subjects of varying sizes, which are warped to the template given in middle column. When expansion takes place during warping of an individual shape to the template, the tissue density ...

In further discussion, CMDs derived from these two types of features will be discussed.

3. Experiments

In this section, we present extensive experimental results to support the hypothesis that residual carries significant amount of information for identifying group differences and that OMS yields superior performance by maintaining group separation between normal and pathologic anatomies in a voxel-wise statistical comparison framework.

For comparison, we perform two types of tests on CMDs in addition to OMSs: (1) t-tests on individual log Jhλ,τ and Rhλ,τ components, and (2) T2 test on (log Jhλ,τ, Rhλ,τ) descriptor, to compute p values. This assigns a level of significance to each feature in terms of differences between healthy and pathologic anatomies. For CMDs, we randomly selected intermediate templates for each subject before conducting tests for all smoothness levels λ. Note that residual was smoothed with a Gaussian filter with various selections of smoothness parameter σ prior to statistical tests mainly due to two reasons. First, it ensures the Gaussianity of the smoothed residual. Second, since the residuals appear only on tissue boundaries, even if tissue atrophy is in the interior of the structure, smoothing produces a more spatially uniform residual. JD, on the other hand, was not smoothed for the T2 test due to its inherent smoothness.

Three datasets are considered comprising 2D toy data and 3D simulated cross-sectional and longitudinal MRI data.

It should also be mentioned that the proposed approach does not depend on the choice of the registration algorithm. We, therefore, utilized viscous fluid based registration for 2D dataset, whereas HAMMER for the 3D case. As discussed later, both algorithms yielded consistent results. In order to simulate the effect of variation in the regularization λ in transformation (shown in Fig. 6), we first registered all individuals with “maximal” flexibility (lowest λ) to multiple templates. These highly aggressive diffeomorphisms were then gradually smoothed to generate transformations hλ,τ with various levels of regularity. In addition, we smooth residual slightly (a Gaussian kernel of standard deviation of 0.5) before optimization in order to ensure regularity of the AEC manifolds.

3.1. 2D Synthetic Dataset

A 2D dataset of 60 shapes was generated by introducing random variability in 12 manually created templates given in Fig. 9. First, the shape of each template is represented by a number of control points, which are then randomly perturbed to simulate the variability across subjects resembling anatomical differences encountered in the gray matter folds of the human brain. Thinning was introduced (5%) in center one-third of the fold of 30 subjects to simulate atrophy (patient data). All subjects were spatially normalized to T12 via T1,…,T11, for smoothness levels of λ = 0,…,42 to construct individual AECs.

Figure 9
2D synthetic dataset – Templates T1,…,T12 simulating gray matter folds.

Effect of transformation parameters on registration accuracy is illustrated in Figs. 10 and and11.11. Note how registration accuracy deteriorates with inappropriate selection of λ and T. Anatomies that are topologically similar to the template are represented well, whereas others require an intermediate level of regularity to avoid the creation of biologically incorrect structures.

Figure 10
2D synthetic dataset – Dependence of residual on regularization λ: : See how registration deteriorates by changing the λ. (a) Subject; (b) Template; (c) Subject warped with the most aggressive transformation (small λ); ...
Figure 11
2D synthetic dataset – Dependence of residual on template T: See how registration deteriorates by changing the template. (a) Subject; (b) Template T12; (c) Subject warped to T12; (d) Template T1; (e) Subject warped to T1.

3.1.1. Voxel-wise Statistical Tests

First we compare the level of significance of group differences based on CMDs and OMSs by computing p-value maps for various values of σ (and λ for tests on CMDs). Minimum of p-value maps is plotted in Fig. 12 as a function of σ to indicate the best achievable performance for the two descriptors.

Figure 12
2D synthetic dataset – The effect of regularization parameter λ and smoothing σ on the performance of morphological descriptors without optimization, as depicted by minimum p-value plots based on the t-test for capturing significant ...

It may be observed from results based on CMDs that residual achieves considerably lower p values as compared to JD. The significance of both log JD (LJD) and the residual increases with σ and λ up to a point after which it starts degrading. Similarly, T2 test also shows best performance for intermediate values of λ (Fig. 13(a)), which means that an overly aggressive transformation tends to contaminate statistical analysis. This agrees with [14], which argued the need of moderate regularization in the transformation. These observations are also in accordance with our hypothesis that the residual carries anatomical information that is complementary to, and perhaps is more important than, the transformation. After optimization, the absolute minimum p value improves from 10−9 at λ = 23 for CMD to 10−10 for OMS (corresponding to σ = 13) as shown in Fig. 13(a).

Figure 13
2D synthetic dataset – OMS versus CMDs for capturing significant group differences based on Hotelling’s T2-tests: (a) Minimum p values; (b) Mean of p-values after thresholding them to 10−2. Optimization yields better performance ...

It should be noted that the minimum p value plots given in Figs. 12 and 13(a) do not completely represent spatial distribution. For classification, it is essential to use a descriptor that leads to small p-values in the entire ROI. For OMS, p values were consistently found to be small in the entire ROI, whereas for CMDs, it is typically an isolated voxel that yields low p-value. A comparison of mean p values in regions with p ≤ 10−2 is given in Fig. 13(b), which demonstrates that OMS (p = 10−4.5) slightly outperforms the best possible CMD (p = 10−4, λ = 39). It may, however, be immediately inferred from the figure that the OMS significantly outperforms the CMDs in general, i.e., over the entire range of λ values.

p-significant ROIs (p = 10−2) were computed for OMS and CMD corresponding to the parameter selection that yielded best performance (λ = 23 for CMD and σ = 13) as shown in Fig. 14. Note how OMS helps in precisely localizing atrophy, which is in accordance with the objective function of Eq. (9). On the other hand, CMDs fail to accurately localize atrophy, with a considerably large number of false positives.

Figure 14
2D synthetic dataset – T2-test based log10 p-value maps for 2D simulated data at σ = 13 thresholded to p ≤ 10−2: (a) CMD λ = 23; (b) OMS.

3.1.2. Invariance to Template Selection

In order to evaluate the robustness of OMS-based statistical analysis, we vary the template T to which all anatomies are finally warped (see Fig. 6). Three different choices of T were considered (T1, T2, T3) to set up three optimization problems. For each Ti, Hotelling’s T2 test was conducted on optimized as well as unoptimized (log10 Jhλ,τ, Rhλ,τ) descriptors. We specifically considered the best parametric selections for CMDs and OMS, which correspond to σ = 13 and λ = 23 as indicated by Fig. 13(b). p-value maps thus computed after appropriate thresholding (p = 10−4) are given in Fig. 15. Figs. 15(a)-(c) indicate that statistical test on CMDs lead to somewhat different regions of significant differences for different template selections. On the other hand, statistical analysis based on OMS is robust to variations in the choice of the template, and better agrees with the true underlying atrophy. In addition, OMS results in low p-values throughout the region of atrophy (p = 10−9 ~ 10−10), whereas those for CMDs are mostly in the range p = 10−5 ~ 10−7.

Figure 15
2D synthetic dataset – Hotelling’s T2 test on CMD and OMS corresponding to the best case of Fig. 13, i.e., σ = 13, and λ = 23. Note that p-value maps are thresholded to p ≤ 10−4. Results are shown on a log ...

3.1.3. Classification

In this section, we further test the performance of the CMD and OMS through pattern classification. For comparison, both OMS and CMD were considered corresponding to their best parameter selections. As indicated by the p-value maps of Figs. 14 and and19,19, not all voxels in morphological descriptors are discriminating. We, therefore, employ Hotelling’s T2 test to rank all features according to their associated p-values, and then select a subset of features with p < 10−2. Note that the resulting feature vectors, in general, form a low dimensional embedding in a high dimensional space. It is, therefore, necessary to find the embedding or to carry out dimension reduction to avoid the curse of dimensionality [29]. Since our proposed method is distance based, we used isomap [44] to find the embedding thereby preserving the distance structure. The neighborhood parameter was varied from 4 to 15 (one fourth of the dataset), and the dimensionality of the embedding was consistently found to be 3 as indicated by the elbow of the residual variance curve. Consequently, a support vector machine based classifier is learned to partition the embedding through corresponding low dimensional feature vectors.

Figure 19
Atrophy maps as captured by CMD and OMS for the 3D dataset – T2 test based p-value maps corresponding to the best results for each descriptor thresholded to p ≤ 10−5: (a) CMD λ = 7, σ = 4; (b) OMS with σ ...

In order to account for nonlinearity of the class boundary, we utilized a radial basis function (RBF) kernel. SVM classifier was trained and tested for both CMD and OMS through 5-fold cross validation, where classification rates were found to increase from 73% for CMD to 91% for OMS. This improvement clearly demonstrates the superiority of the proposed method over traditional approaches.

3.2. 3D MRI Sectional Data with Simulated Atrophy

The second dataset consisted of real MRI scans of 31 subjects. To simulate patient data, 10% atrophy was introduced in 15 randomly selected subjects in a spherical region as shown in Fig. 16 using the simulator of [32]. Five intermediate templates were selected for spatial normalization to generate AECs for all subjects with smoothing levels of λ = 0,…,7.

Figure 16
3D dataset: (a) A subject without atrophy; (b) With 10% simulated atrophy.

Minimum of p-value maps were computed for all values of σ (and λ for tests on CMDs). Minimum log10 p plots given in Fig. 17 show that residual achieves considerably lower p values as compared to LJD, again indicating the significance of residual for capturing group differences. Best performance is achieved at high λ (λ = 7). The dependence of (Jhλ,τ, Rhλ,τ) on λ is eliminated through optimization as indicated in Fig. 18 by OMSs. However, no improvement in minimum p value was achieved through optimization. In this particular example, CMDs perform slightly better (p = 10−9.75) than OMSs (p = 10−10.5). However, OMS appears superior to CMD due to its relative insensitivity to σ. The variation in p-values for OMS is in the range σ = 2−6 is 10−0.75, whereas that for CMDs is 103.5, which highlights that CMDs are much more affected by σ. Small variations in σ, therefore, may considerably degrade CMD-based analysis. OMS is, hence, not only more robust due to better dynamic range, but it also maintains the separation between the two groups as indicated by very low p-values. Importantly, the OMS is clearly better than the CMD if one considers the range of λ, and not the best value of λ. This is important since in practice we do not know the best value of λ and if we were to select the λ that yields the lowest p-values, we would introduce significant multiple comparison problems.

Figure 17
3D dataset – The effect of regularization parameter λ and smoothing σ on the performance of morphological descriptors without optimization, as depicted by minimum p-value plots based on the t-test for capturing significant group ...
Figure 18
3D dataset – Hotelling’s T2-test based minimum p-value plots for CMDs and OMS. The performance of CMDs is highly dependent on λ. Optimization, on the other hand, removes this dependency as evident from the largely stable curve ...

When p was thresholded to find regions with values 10−2, CMDs resulted in false negatives as shown in Fig. 19(a). On the other hand, OMS helped in precisely localizing atrophy Fig. 19(b), which is in accordance with the objective function of Eq. (9).

3.3. 3D MRI Simulated Longitudinal Data (Serial Scans)

The third set of experiments evaluated the concept of building AECs from a different perspective, by generating simulated longitudinal aging profile in an individual MRI scan. 50% atrophy was introduced in three different regions of the brain: (1) posterior cingulate; (2) hippocampus; and (3) superior temporal gyrus over 12 time points (simulating a period of 12 years). Such datasets are quite common in practical studies for modeling normal decay of GM versus disease specific atrophy. The reason for performing this experiment can be appreciated, if one considers Fig. 4, where each of the manifolds now correspond to the same anatomy measured at different time points. The true longitudinal change may then be better estimated if we measure distances across manifolds, rather than distances of individual measurements, which can artificially appear larger than true longitudinal variation.

GM TDM were computed by warping these 12 images to a common template via multiple (five) “intermediate” templates to construct anatomical manifolds. Each manifold, thus, accounted for the morphological description of a particular time point, simulating the situation where an anatomy evolves through a series of manifolds as an individual progresses in age. Factorization of variability along these manifolds is, therefore, necessary to find optimal descriptors that retain only temporal variation.

We compare the proposed approach with traditional TDM based analysis (most aggressive registration with direct warping), by evaluating mean temporal profiles of TDMs in regions with atrophy (Fig. 20) as well as regions without atrophy (Fig. 21). Results highlight the limitations of traditional analysis, which suffers from random fluctuations in TDMs resulting from arbitrariness due to transformation parameters. Consequently, they fail to correctly characterize temporal profiles. For instance, in Fig. 21, a linear regression model highlights atrophy, where in reality no atrophy was present. OMS, on the other hand, helps in minimizing this arbitrariness, to accurately account for underlying atrophy.

Figure 20
Longitudinal dataset – Comparison between most conforming CMD (traditional approach) and OMS in terms of mean GM TDM in the ROI in the presence of atrophy: (a) Posterior Cingulate; (b) Hippocampus; (c) Superior temporal gyrus. Dotted lines represent ...
Figure 21
Longitudinal dataset – Comparison between CMD and OMS in terms of mean TDM in the ROI when no atrophy is present. Note that random fluctuations in traditional approach may contribute to incorrect estimation of the temporal profile. Dotted lines ...

In order to understand the effect of optimization, regression maps were computed. For instance, that for traditional descriptor (Fig. 22(a)) indicates local tissue growth in the highlighted area, which is in fact a side effect of registration errors. It should be noted that a local growth in such cases contaminates the amount of atrophy modeled by traditional descriptors, as observed in Figs. 20 and and2121.

Figure 22
Longitudinal dataset – Comparison between CMD (traditional approach) and OMS in terms of rate of atrophy: (a) CMD based regression map for TDM; (b) OMS based regression map for TDM. Note how OMS was able to eliminate growth falsely picked up by ...

These effects were minimized by OMS based analysis, as indicated by Fig. 22(b), which correctly highlights atrophy in GM.

4. Conclusions

In this paper, we proposed a fundamentally novel framework for morphological analysis with three major contributions. First, the transformation, normalizing an anatomy to a common template space, was combined with the residual for a complete representation of the anatomy. Second, each anatomy was represented through an AEC manifold, reflecting variability due to different templates and regularization parameters. Third, an unsupervised approach was consequently developed for factoring out this unwanted variation due to these parameters, thereby measuring the true inter-subject differences. Such an approach ensures that optimal descriptors do not depend on group associations. Moreover, it yields similar morphological descriptors if the underlying anatomies are similar, irrespective of their group membership.

In the experiments presented in the paper, we focussed on atrophy based analysis, and validated the proposed approach with 2D synthetic and 3D real datasets with simulated cross-sectional as well as longitudinal atrophy. Two features were, therefore, considered, namely (LJD, residual), and TDM. Several improvements on traditional descriptors were readily observed. Residual was consistently found to be highly significant for characterizing group differences, as indicated by considerably low p values. T2 tests on CMDs confirmed our hypothesis that best performance is achieved for not so aggressive registration. Decreasing the level of regularity from very aggressive registration improves performance up to a certain extent after which it starts deteriorating. This suggests the importance of optimal selection of transformation parameters, which actually depends on the individual anatomy.

Optimization was found to introduce several improvements in the performance. First, it resulted in significantly low p-values indicating its ability to detect significant group differences. Second, the p-values were found to be consistently low in the entire ROI of the true atrophy, whereas the CMD results in very low p-values for only a few of the voxels. Consequently, the CMDs were found less discriminating for classification purposes. This was confirmed by the classification of morphological descriptors into healthy and pathologic anatomies through nonlinear support vector machines. As a consequence of optimization, classification rate was found to improve by 18%.

One of the most important aspects of optimization is its ability to precisely highlight the regions of significant differences. The CMD, however, resulted in relatively poorer localization of such regions, with considerable false positives and false negatives, thus lowering the specificity and sensitivity of the analysis. Invariance to the choice of the template was also readily observed for the OMS, in contrast to the unoptimized representation. It should be noted that the best performance of the CMD in terms of the lowest p-value requires intermediate values of λ, which are not exactly known a priori. On could potentially choose the λ that yields the most significant group differences. However, such an approach by construction amplifies the multiple comparison problem, which is prevalent in voxel-based statistical analysis of medical images. In addition, p-values were found to be less sensitive to σ after optimization, yielding a better dynamic range.

For the longitudinal dataset, the proposed approach was able to detect group differences where traditional analysis failed. Optimization helped in minimizing random fluctuations in the temporal profiles for more accurate characterization. Application of this approach for measuring longitudinal changes in serial scans is particularly interesting and encouraging. The problem of “jittery” measurements from serial scans is very important when evaluating subtle changes due to disease progression or response to treatment. These random variations that are unrelated to the true underlying morphological changes significantly reduce sensitivity and specificity of these measurements. Our approach effectively removed the confounding variations that are due to the template and smoothness parameter selection, and allowed us to obtain considerably more stable estimates of longitudinal change.

Future work includes application to real datasets for computer aided diagnosis as well as nonlinear modeling of AEC manifolds.

Table 1
List of abbreviations introduced in the paper


1The term “complete” is used, since this morphological descriptor does not discard any image information via consideration of the residual.

2Warp after the convergence of the registration algorithm. Note that for different λs, the algorithm will converge to different transformations.


1. Ashburner J, Friston K. Nonlinear spatial normalization using basis functions. Human Brain Mapping. 1999;7(4):254–266. [PubMed]
2. Ashburner J, Friston KJ. Voxel-based morphometry – the methods. NeuroImage. 2000;11(6):805–821. [PubMed]
3. Ashburner J, Hutton C, Frackowiak R, Johnsrude I, Price C, Friston K. Identifying global anatomical differences: deformation-based morphometry. Human Brain Mapping. 1998;6(6):348–357. [PubMed]
4. Avants BB, Epstein CL, Gee JC. Geodesic image normalization and temporal parameterization in the space of diffeomorphisms. MIAR. 2006:9–16.
5. Baloch S, Verma R, Davatzikos C. An anatomical equivalence class based joint transformation-residual descriptor for morphological analysis. IPMI. 2007:594–606. [PubMed]
6. Baron J, Chetelat G, et al. In vivo mapping of gray matter loss with voxel-based morphometry in mild Alzheimer’s disease. Neuroimage. 2001;14(2):298–309. [PubMed]
7. Basser P, Pierpaoli C. Microstructural and physiological features of tissues elucidated by quantitative-diffusion-tensor MRI. Journal of Magnetic Resonance, Series B. 1996;111:209–219. [PubMed]
8. Beg MF, Miller MI, Trouvé A, Younes L. Computing large deformation metric mappings via geodesic flows of diffeomorphisms. Int J Comput Vision. 2005;61(2):139–157.
9. Bhatia K, Hajnal J, Puri B, Edwards A, Rueckert D. Consistent groupwise non-rigid registration for atlas construction. ISBI’04. 2004:908–911.
10. Bhatia KK, Aljabar P, Boardman JP, Srinivasan L, Murgasova M, Counsell SJ, Rutherford MA, Hajnal JV, Edwards AD, Rueckert D. Groupwise combined segmentation and registration for atlas construction. MICCAI (1) 2007:532–540. [PubMed]
11. Blezek DJ, Miller JV. Atlas stratification. In: Larsen R, Nielsen M, Sporring J, editors. MICCAI (1). Vol. 4190 of Lecture Notes in Computer Science. Springer; 2006. pp. 712–719.
12. Bookstein F. Principal warps: Thin-plate splines and the decomposition of deformations. IEEE Trans on Pattern Analysis and Machine Intelligence. 1989;11(6):567–585.
13. Bookstein F. Voxel-based morphometry should not be used with imperfectly registered images. Neuroimage. 2001;14(6):1454–1462. [PubMed]
14. Cachier P. MICCAI ‘01: Proceedings of the 4th International Conference on Medical Image Computing and Computer-Assisted Intervention. Springer-Verlag; London, UK: 2001. How to trade off between regularization and image similarity in non-rigid registration? pp. 1285–1286.
15. Cao J, Worsley K. The geometry of the Hotelling’s T2 random field with applications to the detection of shape changes. Annals of Statistics. 1999;27(3):925–942.
16. Chetelat G, Desgranges B, et al. Mapping gray matter loss with voxel-based morphometry in mild cognitive impairment. Neuroreport. 2002;13(15):1939–1943. [PubMed]
17. Christensen G, Rabbit R, Miller M. Proc Conference on Information Sciences and Systems. Johns Hopkins University; 1993. A deformable neuroanatomy textbook based on viscous fluid mechanics; pp. 211–216.
18. Christoudias CM, Morency L-P, Darrell T. Light field appearance manifolds. ECCV. 2004:481–493.
19. Chung M, Worsley K, Paus T, Cherif C, Collins D, Giedd J, Rapoport J, Evans A. A unified statistical approach to deformation-based morphometry. NeuroImage. 2001;14(3):595–600. [PubMed]
20. Collins D, Paus T, Zijdenbos A, Worsley K, Blumenthal J, Giedd J, Rapoport J, Evans A. Age related changes in the shape of temporal and frontal lobes: an MRI study of children and adolescents. Soc Neurosci Abstr. 1998;24
21. Davatzikos C, Genc A, Xu D, Resnick S. Voxel-based morphometry using RAVENS maps: methods and validation using simulated longitudinal atrophy. Neuroimage. 2001;14:1361–1369. [PubMed]
22. Davatzikos C, Vaillant M, Resnick S, Prince J, Letovsky S, Bryan R. A computerized approach for morphological analysis of the corpus callosum. Journal of Comp Assisted Tomography. 1996;20(1):88–97. [PubMed]
23. Davis B, Lorenzen P, Joshi S. Large deformation minimum mean squared error template estimation for computational anatomy. ISBI’04. 2004:173–176.
24. Fox N, Crum W, Scahill R, Stevens J, Janssen J, Rossor M. Imaging of onset and progression of Alzheimers disease with voxel compression mapping of serial magnetic resonance images. Lancet. 2001;358:201–205. [PubMed]
25. Gaser C, Volz H, Kiebel S, Riehemann S, Sauer H. Detecting structural changes in whole brain based on nonlinear deformationsapplication to schizophrenia research. Neuroimage. 1999;10(2):107–113. [PubMed]
26. Geng X, Kumar D, Christensen GE. Transitive inverse-consistent manifold registration. IPMI. 2005:468–479. [PubMed]
27. Grenander U. Tech rep. Brown University; 1983. Tutorial in pattern theory: a technical report.
28. Grenander U, Miller MI. Computational anatomy: An emerging discipline. Quarterly of Applied Mathematics. 1998;56:617–694.
29. Hastie T, Tibshirani R, Friedman JH. The Elements of Statistical Learning. Springer; 2001.
30. Joshi S. Ph.D. thesis. Washington University; St Louis: 1998. Large deformation diffeomorphisms and gaussian random fields for statistical characterization of brain sub-manifolds.
31. Joshi S, Davis B, Jomier M, Gerig G. Unbiased diffeomorphic atlas construction for computational anatomy. NeuroImage. 2004 Jan;23(Suppl 1):S151–60. [PubMed]
32. Karachali B, Davatzikos C. Simulation of tissue atrophy using a topology preserving transformation model. IEEE Trans on Medial Imaging. 2006;25(5):649–652. [PubMed]
33. Lee K, Ho J, Yang M, Kriegman D. Video-based face recognition using probabilistic appearance manifolds. IEEE Conference on Computer Vision and Pattern Recognition. 2003:I, 313–320.
34. Leow A, Klunder A, Jack C, Toga A, Dale A, Bernstein M, Britson P, Gunter J, Ward C, Whitwell J, Borowski B, Fleisher A, Fox N, Harvey D, Kornak J, Schuff N, Studholme C, Alexander G, Weiner M, Thompson P. Longitudinal stability of MRI for mapping brain change using tensor-based morphometry. Neuroimage. 2006;31(2):627–640. [PMC free article] [PubMed]
35. Lina, TAKAHASHI T, IDE I, MURASE H. Construction of appearance manifold with embedded view-dependent covariance matrix for 3D object recognition. IEICE TRANSACTIONS on Information and Systems. 2008;E91-D(4):1091–1100.
36. Lorenzen P, Prastawa M, Davis B, Gerig G, Bullitt E, Joshi S. Multi-modal image set registration and atlas formation. Medical Image Analysis. 2006;10(3):440–451. [PMC free article] [PubMed]
37. Makrogiannis S, Verma R, Davatzikos C. Anatomical equivalence class: A computational anatomy framework using a lossless shape descriptor. IEEE Trans on Biomedical Imaging. 2007;26(4):619–631. [PubMed]
38. Miller M, Banerjee A, Christensen G, Joshi S, et al. Statistical methods in computational anatomy. Statistical Methods in Medical Research. 1997;6:267–299. [PubMed]
39. Miller M, Trouvé A, Younes L. On the metrics and Euler-Lagrange equations of computational anatomy. Annual Review of Biomedical Engineering. 2002;4:375–405. [PubMed]
40. Miller M, Younes L. Group actions, homeomorphisms, and matching: a general framework. International Journal of Computer Vision. 2001;41(1):61–84.
41. Shan C, Gong S, McOwan P. Appearance manifold of facial expression. ECCV Workshop on Computer Vision in Human-Computer Interaction. 2005:221.
42. Shen D, Davatzikos C. HAMMER: Hierarchical attribute matching mechanism for elastic registration. IEEE Transactions on Medical Imaging. 2002;21(11):1421–1439. [PubMed]
43. Studholme C, Cardenas V. Population based analysis of directional information in serial deformation tensor morphometry. MICCAI (2) 2007:311–318. [PMC free article] [PubMed]
44. Tenenbaum J, de Silva V, Langford J. A global geometric framework for nonlinear dimensionality reduction. Science. 2000;290:2319–2323. [PubMed]
45. Thirion J, Calmon G. Deformation analysis to detect quantify active lesions in 3D medical image sequences. IEEE Trans Med Imag. 1999;18:429–441. [PubMed]
46. Thompson D. On growth and form. Cambridge University Press; 1917.
47. Thompson P, Giedd J, Woods R, MacDonald D, Evans A, Toga A. Growth patterns in the developing human brain detected using continuum-mechanical tensor mapping. Nature. 2000;404(6774):190–193. [PubMed]
48. Thompson P, Toga A. A surface-based technique for warping three-dimensional images of the brain. IEEE Trans on Med Imaging. 1996;15:402–417. [PubMed]
49. Trouvé A, Younes L. Metamorphoses through lie group action. Foundations of Computational Mathematics. 2005;5(2):173–198.
50. Twining C, Cootes T, Marsland S, Petrovic V, Schestowitz R, Taylor C. A unified information-theoretic approach to groupwise non-rigid registration and model building. IPMI’05. 2005:1–14. [PubMed]
51. Wakin MB. Ph.D. thesis. Rice University; Houston: 2007. The geometry of low-dimensional signal models.
52. Wakin MB, Donoho DL, Choi H, Baraniuk RG. High-resolution navigation on non-differentiable image manifolds. ICASSP’05. 2005;5:1073–1076.
53. Wang Y, Peterson B, Staib L. 3D brain surface matching based on geodesics and local geometry. Computer Vision and Image Understanding. 2003;89:252–271.
54. Woods R, Dapretto M, Sicotte N, Toga A, Mazziotta J. Creation and use of a Talairach-compatible atlas for accurate, automated, nonlinear intersubject registration, and analysis of functional imaging data. Hum Brain Mapp. 1999;8(23):73–79. [PubMed]
55. Woods R, Mazziotta J, Cherry S. MRI-PET registration with automated algorithm. J Comput Assist Tomogr. 1993;17:536–546. [PubMed]