PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
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 August 15.
Published in final edited form as:
PMCID: PMC2857762
NIHMSID: NIHMS189090

Multi-contrast large deformation diffeomorphic metric mapping for diffusion tensor imaging

Abstract

Diffusion tensor imaging (DTI) can reveal detailed white matter anatomy and has the potential to detect abnormalities in specific white matter structures. Such detection and quantification are, however, not straightforward. The voxel-based analysis after image normalization is one of the most widely used methods for quantitative image analyses. To apply this approach to DTI, it is important to examine if structures in the white matter are well registered among subjects, which would be highly dependent on employed algorithms for normalization. In this paper, we evaluate the accuracy of normalization of DTI data using a highly elastic transformation algorithm, called large deformation diffeomorphic metric mapping. After simulation-based validation of the algorithm, DTI data from normal subjects were used to measure the registration accuracy. To examine the impact of morphological abnormalities on the accuracy, the algorithm was also tested using data from Alzheimer’s disease (AD) patients with severe brain atrophy. The accuracy level was measured by using manual landmark-based white matter matching and surface-based brain and ventricle matching as gold standard. To improve the accuracy level, cascading and multi-contrast approaches were developed. The accuracy level for the white matter was 1.88 ± 0.55 and 2.19 ± 0.84 mm for the measured locations in the controls and patients, respectively.

Keywords: Human, White matter, Magnetic resonance imaging, Diffusion tensor, Normalization, LDDMM

Introduction

Diffusion tensor imaging (DTI) can reveal the detailed anatomy of white matter tracts, based on diffusion anisotropy and fiber orientation information (Basser et al., 1994; Pierpaoli et al., 1996; Makris et al., 1997; Pajevic and Pierpaoli, 1999; Catani et al., 2002; Wakana et al., 2004; Mori et al., 2005). Theoretically, this enables us to study the effects of diseases on specific white matter tracts. However, quantification of the delineated white matter anatomy is not always straightforward, especially not for group studies when using population-based maps.

The image analysis based on image normalization is one of the most widely used methods for automated quantification of anatomical features revealed in MR images. The normalization method has been applied to morphological analyses of T1-weighted images (for review, see, e.g., Ashburner and Friston, 2000) and, more recently, to DTI studies (Alexander et al., 2001; Jones et al., 2002; Park et al., 2003; Xu et al., 2003; Zhang et al., 2006; Muller et al., 2007; Yushkevich et al., 2008). In this approach, acquired images are transformed (normalized) to a template, which makes it possible to perform subsequent voxel-by-voxel analyses. In this process, the normalization accuracy significantly influences the result (Bookstein, 2001). In our previous paper, we estimated the accuracy level for normalization-based analyses of DTI data using a simple linear transformation (Mori et al., 2008). By using a population-averaged template (ICBM-DTI-81), most of the core white matter regions could be registered within 3 mm for the normal adult population. However, this registration quality cannot be guaranteed for patient populations in which there are substantial anatomical abnormalities. Registration errors complicate the interpretation of population analyses. For example, any detected FA abnormality could be due to real FA differences caused by microstructural changes such as axonal loss or myelin abnormalities or due to regional morphological differences (e.g., enlarged ventricles) causing registration errors. To separate out the effects of contrast (such as FA) from those of morphological changes, the use of nonlinear transformation is the next logical step.

In this paper, we tested the a highly non-linear transformation method, called large deformation diffeomorphic metric mapping (LDDMM) (Miller et al., 1993a,b; Miller et al., 1997) to normalize DTI data. The LDDMM algorithm calculates diffeomorphic transformations between images of anatomical configurations by computing the geodesic flow in the space of diffeomorphisms between these images. This allows the study of anatomical structures as a metric space. The diffeomorphic transformations are invertible and are smooth functions with a smooth inverse. In these transformations, the disjoint shapes remain disjoint and there is no fusion of points because of the one-to-one property of diffeomorphisms. Connected shapes also remain connected because of the continuity property, and the smoothness of the object boundaries are preserved because of the smoothness property of the diffeomorphisms.

There are several unique steps required to apply a non-linear transformation method to DTI data. First, DTI can produce various types of image contrasts, so the first step is to choose which contrasts will be used to drive the registration. Second, to apply the resultant transformation to the tensor field, tensor reorientation must be performed. In this study, we first tested our tensor transformation algorithm using a simulated DTI phantom. The method was then tested on a normal elderly population and Alzheimer’s disease (AD) patients with severe brain atrophy. To ensure excellent registration accuracy for the entire brain, multi-contrast LDDMM was developed.

Methods and Materials

Subjects

Institutional Review Board approval was obtained for the study and written, informed consent, including HIPAA compliance, was obtained from all subjects. DTI data from 13 AD patients (75.4 ± 7.0 years old) and age-matched 18 normal adult control subjects (73.8 ± 7.6 years old) were used in this study.

MRI scans

The 3.0 T MR scanners (Gyroscan NT, Philips Medical Systems) were used. DTI data were acquired with a single-shot, echo-planar imaging (EPI) sequence with sensitivity encoding (SENSE), using a parallel-imaging factor of 2.5 (Pruessmann et al., 1999). The imaging matrix was 96 × 96 and zero-filled to 256 × 256 pixels. The field-of-view was 211 × 211 mm. Transverse sections of 2.2 mm thickness were acquired parallel to the anterior commissure–posterior commissure line. A total of 55–60 sections covered the entire hemisphere and the brainstem without gaps. Diffusion weighting was encoded along 32 independent orientations, and the b-value was 700 mm2/s. Five additional images with minimal diffusion weighting (b≈33 mm2/s) were also acquired (called b0 images hereafter). The scanning time per dataset was approximately 4 min, which was repeated twice to improve the SNR.

DTI data processing

The DTI datasets were transferred to a personal computer running a Windows platform and were processed using DtiStudio (www.MriStudio.org) (Jiang et al., 2006). Images were first realigned for co-registration and eddy current distortion correction, using the affine transformation of the Automatic Image Registration (AIR) package (Woods et al., 1998), in which the first minimally diffusion-weighted image was used as a template. The six elements of the diffusion tensor were calculated for each pixel using multivariate linear fitting. After diagonalization, three eigenvalues and eigenvectors were obtained. For the anisotropy map, fractional anisotropy (FA) was used (Pierpaoli and Basser, 1996). The eigenvector associated with the largest eigenvalue was used as an indicator of fiber orientation.

Normalization process

All normalization procedures, including linear and LDDMM transformation, were performed using an in-house program called Landmarker (www.MriStudio.org, Kennedy Krieger Institute and Johns Hopkins University, X. Li, H. Jiang, and S. Mori). As a target for brain normalization, we used a single-subject white matter atlas called JHU-DTI-MNI atlas (also known as the “Eve atlas”), which is built in into the Landmarker and contains T1, T2, and DTI-derived contrasts. The atlas is in the ICBM-152 coordinate system and has a matrix size of 181 × 217 × 181 (1 mm isotropic resolution).

For the initial atlas-subject registration, affine transformation was used using the least diffusion-weighted images (b0 images). The transformation matrix was then applied to the calculated diffusion tensor field, based on the method described by Alexander et al. (2001) and Xu et al. (2003). This process took less than 1 min with 2–3 GHz Xeon processors.

After the linear normalization, the images were submitted for LDDMM. In this study, we tested several registration approaches in terms of choices of images that drive the registration. These were: b0 image only; FA map only; and b0 + FA dual-contrast registration. Because of the substantial computer resources required by the multi-contrast LDDMM, the data were resampled to 96 × 128 × 96 (2 × 2 × 2 mm) before the calculation and the computational time was approximately 5 h using one of the 32 CPUs in a cluster computer. Because of this extensive computational requirement, the Landmarker adopts a centralized “remote” strategy. Namely, the software creates a data packet that is automatically transferred to a designated centralized computational resource with enough computational power. Once the computation is completed, the users can activate the data retrieval function of Landmarker, through which the results can be automatically retrieved. This remote approach can be implemented on local Linuxor Unix servers with proper configurations. The computation time is typically 6 h with a 3 MHz Zeon processor.

Details of LDDMM-based image registration

Single-contrast

The deformable template model of Granander (Granander and Miller, 1996) models the observed anatomical images, I, as an orbit under the group of transformations, G, acting on a template image, I0.

I=GI0={I0°ϕ1,ϕG}.
(1)

Utilizing this model, for any given two images, I0, I1: Ω[subset or is implied by]R3R, the LDDMM algorithm (Beg, 2003; Beg et al., 2005) calculates the diffeomorphic transformation, [var phi]: ΩΩ, registering the images such that I1=I0 ° [var phi]−1. Ω R3 is the 3D cube on which the data is defined. The optimal transformation, [var phi], is generated as the end point, ϕ=φ1ν, of the flow of smooth time-dependent vector field, νt[set membership]V, t[set membership] [0,1], with the ordinary differential equation,

dφtνdt=νt(φtν),t[0,1],
(2)

where ϕ0 is the identity transformation such that ϕ0 (x) = x,x[set membership]Ω. Then, the optimal transformation, [var phiv with circumflex], is calculated by integrating the vector field that is found by minimizing the following equation.

ν^=arg minν:dφtν/dt=νt(φtν)(01νtV2dt+1σ2I0°ϕ1I1L22).
(3)

The solution of Eq. (3) is ensured to be in the space of diffeomorphisms, by enforcing smoothness on the vector fields, ν[set membership]V. The required smoothness is enforced by defining the norm on V as | f | V = | Lf |L2. L is a differential operator defined as L=−α[nabla]2 + γ I3 × 3, where I3 × 3 is the identity operator and [nabla]2 is the Laplacian operator. | · | L2 is the L2 norm for the square integrable functions defined on Ω. The gradient of the cost in Eq. (3) is:

νEt=2νtK(2σ2Dφt,1νJt0(Jt0Jt1)),
(4)

where the notation φs,t=φs°φt1 is used. In Eq. (4), Jt0=I0°φt,0 and Jt1=I1°φt,1. |Df| is the determinant of the Jacobian matrix. K:L2(Ω,R3)→V is a compact self-adjoint operator, defined by a,bL2=Ka,bV, which satisfies K(LL)g=g for any smooth vector field g[set membership]V. The parameter σ provides weighting between data matching and smoothness regularization terms. In the LDDMM algorithm, Eq. (3) is solved with a gradient descent algorithm using Eq. (4). In Eq. (4), the effect of the operator, K=(LL)−1, is low-pass filtering. The parameters α and γ define the magnitude and frequency response of this filter. Fig. 1 demonstrates the effects of these parameters on the LDDMM-based image transformation. In this simulation, a simple circle object with a radius of 30 pixels (r=30) is transformed to a larger circle (r=40). The low α/γ ratio leads to less elastic transformation, similar to the linear normalization. As the ratio decreases, the transformation is more localized to the edge of the circle where the contrast difference between the two circles is concentrated. Fig. 1B shows the deformed rectangular grid for different α/γ ratios. Fig. 1C shows the norm of the displacement vector as a function of distance from the center of the circle for deformations with different α/γ ratios. As the ratio decreases, it can be seen that the transformation becomes more nonlinear and localized at the peripheral where the contrast difference resides.

Fig. 1
Simple simulations to demonstrate the effect of the α/γ ratio, in which a circle of the 30-pixel radius is transformed to the 40-pixel radius, (A) with varying α/γ ratios (B and C). In (B), transformations of a grid are ...

Multi-contrast

For any given two vector valued images, I0 = [I01, I02,…,I0C] and I1 = [I11, I12,…,I1C] with I0c, I1c: Ω[subset or is implied by]R3R,c=1,…,C, the diffeomorphic transformation,[var phi], registering these two images, such that I1 = I0 ° [var phi]−1 or [I11, I12,…,I1C]=[I01 °[var phi]−1, I02 °[var phi]−1,…,I0C °[var phi]−1], is assumed to be generated as the end point of the flow of the smooth time-dependent vector field, νt[sm epsilon]V, with the ordinary differential equation, φtν/t=νt(φtν),t[0,1]. The optimal transformation, [var phiv with circumflex] is calculated by integrating the vector field, which is found by minimizing the following equation.

ν^=arg minν:dφtν/dt=νt(φtν)(01νtV2dt+c=1c{1σc2I0c°ϕ1I1cL22}).
(5)

Here the index c denotes the contrast images. Again, the solution for this equation is ensured to be in the space of diffeomorphisms, by enforcing smoothness on the vector fields ν[set membership]V. The gradient descent algorithm can be used to solve this equation. The gradient of the cost in Eq. (5) is

νEt=2νti=1N(K(2σc2Dφt,1νJt0c(Jt0cJt1c))),c=1,,C.
(6)

In this equation, the parameters, σc, control the weighting of contrast-matching terms and smoothness regularization terms. In this study, we used equal weighting for FA and b0 images.

Tensor transformation

The LDDMM algorithm is a powerful technique to calculate diffeomorphic transformations between scalar-valued images. However, in its original form, the LDDMM cannot be directly used to register DTI images. With the extension of the LDDMM as a multi-contrast image-matching algorithm, we can seed various scalar-valued isotropy or anisotropy images obtained from DTI data as contrasts to the multi-contrast LDDMM algorithm, and we can calculate diffeomorphisms between tensor images. Deformation of DTI images using the calculated transformations also requires the reorientation of the tensor at each voxel. Let M: Ω[subset or is implied by] R3R3 × 3 be a tensor image and let [var phi] be a diffeomorphism registering the tensor image, M, to another tensor image calculated using multi-contrast LDDMM. To transform the tensor image, M, with the diffeomorphism, [var phi], we use the method in Alexander et al. (2001), which is based on the Gram-Schmidt orthonormalization. At each point x[set membership]Ω of the tensor image, the tensor, M (x), can be decomposed as, M(x)=λ1ν1ν1T+λ2ν2ν2T+λ3ν3ν3T, where λ1, λ2, λ3 with λ1 ≥ λ2 ≥ λ3 are eigenvalues and ν1, ν2, and ν3 are the corresponding unit eigenvectors. The transformation on M(x) with [var phi] can be defined as:

ϕM(x)=(λ1ν^1ν^1T+λ2ν^2ν^2T+λ3ν^3ν^3T)ϕ1(x)=Mϕ1(x)=[M11ϕ1(x)M12ϕ1(x)M13ϕ1(x)M12ϕ1(x)M22ϕ1(x)M23ϕ1(x)M13ϕ1(x)M23ϕ1(x)M33ϕ1(x)],where

ν^1=Dϕν1Dϕν1,ν^2=Dϕν2ν^1,Dϕν^2ν^1(Dϕν22ν^1,Dϕν22)1/2,ν^3=ν^1 × ν^2.

In these equations, D[var phi] is the Jacobian matrix of [var phi]. This tensor transformation strategy depends on the assumption that the tissue microstructure does not change during the transformation, which means that the tensor shape (3D ellipsoid shape) in each voxel does not change during deformation, i.e., the eigenvalues do not change, but the principal eigenvector changes according to the Jacobian matrix, D[var phi]. The plane generated by the eigenvectors, ν1, ν2, changes to the plane generated by D[var phi]v1, D[var phi]v2.

Application of the multi-contrast LDDMM to human brains

For all 13 AD subjects and 18 control subjects, we mapped each subject to the template using one of the following approaches:

  1. Only b0 images were used with the single-contrast LDDMM and the transformations, ϕib0, i = 1,…,31, were calculated (b0-LDDMM).
  2. Only FA images were used with the single-contrast LDDMM and the transformations, ϕiFA, i = 1,…,31, were calculated (FA-LDDMM).
  3. Both b0 and FA images were used with the multi-contrast LDDMM and the transformations, ϕib0FA, i = 1,…,31, were calculated (b0 + FA-LDDMM).

In these transformations, the subscripts indicate the subject number (13 AD and 18 control subjects) and the superscripts denote the type of contrast images used in mappings.

Measurement of registration quality

Landmark-based measurements

We used manual landmarks as the gold standard for the degree of matching of various white matter structures between the atlas and subject data after normalization. For the landmarks, we used a so-called “standard landmark set,” which consists of 237 landmarks at easy-to-identify white matter structures in the JHU-DTI-MNI atlas, as previously described (Mori et al., 2008). These landmarks were transferred to the normalized patient images and manually moved to corresponding structures. The average distance between the positions of the 237 landmarks in the template and the normalized images represents the normalization quality of the white matter.

To compare the registration results with different approaches (affine, b0-LDDMM, FA-LDDMM, and b0 + FA-LDDMM), we performed a one-sided Kolmogorov–Smirnov (KS) test on cumulative distributions of the pair-wise landmark distances, in which the number of landmarks (Y axis) within a given distance (X axis) is plotted. For the calculation of the cumulative distributions, the averaged distance function d:{1,2,…,M}→[0,G][subset or is implied by]R is defined on the set of M landmark points of the template image and this function has the maximum value G. The ECDF (empirical cumulative distribution function) of d is

ECDF(i)=n(i)M,i=1,,M,

where n(i) is the number of elements or the cardinality of the set S = {x: d (x) ≤d (i)}.

We hypothesized that there was no significant difference between two distributions against the alternative that the two distributions were significantly different. The empirical distribution estimate of KS statistics was obtained using the permutation-based resampling. The entire 237 landmarks were resampled as a whole without replacement for each distribution. Then the KS statistic was computed for each resampled dataset. By repeating this process ten thousand times, an empirical null distribution of the KS statistic was constructed, and the P-value was calculated as a percentage of the KS statistics greater than the KS value of the original dataset at 5% significance level.

Surface-based measurement

Landmark-based matching is effective for measuring the registration quality of the white matter. However, it is difficult to apply to cortical and ventricular matching because of the paucity of easy-to-identify landmarks, especially in the cortex. We defined the outer brain surfaces and the ventricle surface of the b0 images in the template and in each subject and created triangulated meshes. For the surface generation, the image voxels were regarded as values defined as points in a rectangular lattice and eight neighboring voxels formed the eight vertices of a cube. Each cube was decomposed into five tetrahedra. Then the intensity values in the image were used to find whether a given tetrahedron intersected the isosurface. These intersections formed a triangular mesh that presents the surface. Each surface can be represented as a set of vertex coordinates of these triangulated meshes (Fig. 2). The Euclidian distances of the closest vertices between the template and the normalized subject images were calculated. Then, the empirical cumulative distribution curves were calculated similar to the landmark-based measurement. To compare the results with different registration approaches, we again performed a one-sided Kolmogorov–Smirnov (KS) test on the distance cumulative distributions. Note that this surface-to-surface matching measurement does not necessarily measure cortex-to-cortex registration quality, which is beyond the scope of this paper.

Fig. 2
A flowchart for the surface-based measurement of registration accuracy. The brain and ventricle surfaces were defined by triangular meshes for the template and the normalized subject brains. Then distances between the closest triangle vertices of the ...

Results

Analysis of the impact of the α/γ ratio

The impact of the α/γ ratio on the transformation is demonstrated in Fig. 1. In Fig. 3, we tested the effect of the ratio using circular objects with a varying degree of sinusoidal waves (amplitude=5 pixels) as a template. The total number of the sinusoidal waves around the perimeter (r=40 pixels) is changed to 10, 20, 30, and 40, in which the wavelength is 25.13, 12.57, 8.38, and 6.28 pixels, respectively. A circle object with the same diameter was transformed to these targets with different α/γ ratios and the mismatch was measured by counting the number of pixels that were not aligned to the template before and after the transformation. This simulation provides an idea about the relationship between the elasticity of the transformation controlled by the α/γ ratio and the complexity of the target structures. The results indicate that the transformation provided by a ratio larger than 0.01 is not elastic enough and no improvement (=1) was observed for this type of highly non-linear shape difference. The matching quality improves rapidly as the ratio decreases. For the low-frequency modulation (f=10), no improvement is observed below 0.001, but as the frequency increases (f>10), higher elasticity (α/γ<0.001) is required.

Fig. 3
Simulation of transformation accuracy of LDDMM with a varying α/γ ratio. Circular objects with different degrees of sinusoidal perimeter modulation (frequency=10, 20, 30, and 40) are created from a circle. All objects have a radius of ...

Fig. 4 demonstrates the impact of the α/γ ratio on the transformation of human brains; a control subject and an AD patient with severe anatomical deformation are depicted. The color images are Jacobian maps of the transformation from template to subject showing how much change in local tissue volume was introduced by the transformation. For the control subject, the gross mis-registration is appreciable after affine transformation, which is mostly removed by LDDMM with α/γ = 0.01. The Jacobian map at α/γ = 0.01 clearly shows the degree of local (non-linear) expansion (Jacobian>1) or shrinkage (Jacobian<1) to match the brain shape of the atlas to that of the subject. As the α/γ ratio decreases, the more local and larger transformation is exerted, especially in cortical areas, as evidenced by the Jacobian maps, even though matching improvement is not immediately obvious visually. The difference between an α/γ of 0.0025 and 0.001 is minimum and is confined to the cortex. For the AD patient, the improvement is visually appreciable up to α/γ =0.005 (please notice how the ventricle shrinks and the corpus callosum enlarges, approaching those of the template).

Fig. 4
Transformation of a normal adult and a severely abnormal case (AD patient) with gradually decreasing α/γ ratios. For a visual clue, the outline of the white matter of the template is superimposed on the data. The color maps are Jacobian ...

In Fig. 4, the α/γ ratio is gradually decreased (called a “cascading” ratio, hereafter). This ensures that there is only a small amount of required transformation at each step up to 0.001. If such a highly elastic transformation is applied directly to linearly normalized images without proper initialization, the registration may converge to an apparently wrong solution, as demonstrated in Fig. 5, especially when the patient anatomy is dramatically different from the template. The cascading approach is more computationally expensive, but ensures more robust results. Based on these results, we adopted a three-step LDDMM with a decreasing α/γ of 0.01–0.005–0.0025 for all further experiments.

Fig. 5
Examples of transformation with small α/γ ratios (highly elastic transformation) with and without the cascading approach. Performing the highly elastic transformation directly on linearly registered images may lead to an inaccurate solution ...

Qualitative examination of the improvement by multi-contrast LDDMM

In Fig. 6, an example of LDDMM normalization is shown for an AD patient with severe brain atrophy. In this figure, results from the four different transformation approaches are compared: affine (before LDDMM); b0-LDDMM; FA-LDDMM; and b0 + FA-LDDMM. From visual inspection, the following points can be immediately appreciated. First, the affine transformation (Fig. 6E) is capable of matching the overall outer brain shape, but fails to register inner structures. Second, when the b0 image is used for LDDMM (Fig. 6C), not only the brain surface, but also the ventricle surface, is well registered. However, the white matter structures revealed by the FA map are not well-matched. The opposite results are obtained when an FA map is used (Fig. 6D), which carries rich anatomical information about the white matter, but the gray matter–CSF boundary is not well depicted. As a result, while the white matter matching looks excellent, there is a gross mismatch at the brain and ventricle surfaces. The entire brain has noticeably better registration accuracy when the FA and b0 images are used simultaneously to drive the registration (Fig. 6B).

Fig. 6
Demonstration of LDDMM using an FA image only, a b0 image only, and an FA + b0 dual contrast. The upper and bottom rows show b0 and FA images from the atlas (first column) and an AD patient (2nd–4th columns). The b0 images carry a strong contrast ...

Quantitative measurements of transformation accuracy

To measure the registration accuracy, we performed two types of the measurements using control and AD populations. First, we measured the registration accuracy of the white matter by placing 237 landmarks on major white matter structures in both the template and the patient image. Second, the brain and ventricle surfaces were defined using the b0 images and surface matching was measured.

The results of the landmark-based white matter measurements are shown in Fig. 7. After affine transformation, registration accuracy measured by average landmark distance was 3.51 ± 1.16 and 4.40 ± 1.67 mm for the normal and AD populations, respectively. When FA maps were used to drive the LDDMM, the registration was improved to 1.65 ± 0.47 and 1.84 ± 0.65 mm for the normal and AD populations, respectively. The test–retest reliability of the landmark placement was 1.58 ± 0.60, indicating that the matching quality approached the accuracy of this measurement. The b0 + FA-LDDMM achieved registration accuracy similar to the FA-LDDMM (1.88 ± 0.55 and 2.19 ± 0.84 mm for the control and AD groups, respectively), but b0-LDDMM showed hardly any improvement over the affine registration (3.12 ± 0.54 and 3.99 ± 0.84 mm). A statistically significant improvement in the registration accuracy was observed by using FA (P-value=0.01) or b0+ FA (P-value = 0.001) mappings for 5% significance level.

Fig. 7
Results of landmark-based accuracy measurements of LDDMM transformation of the white matter. The 237 landmarks were manually placed on the template and the subjects in the corresponding locations after three different types of LDDMM registration (FA only, ...

In Figs. 8A and B, results of the brain surface matching are shown. The FA-driven LDDMM leads to poor normalization quality. The b0 image with high contrast for the brain boundary leads to significant improvement in normalization compared to the affine transformation. The same trend is observed for the ventricle shape matching (Figs. 8C and D). Although significant improvement is found by the FA-LDDMM, the b0 contrast is necessary for better registration accuracy.

Fig. 8
Measurements of surface-matching accuracy for the brain (A and B) and the ventricle (C and D). Images (A) and (C) show 3D maps of the residual difference after affine, b0-LDDMM, FA-LDDMM, and b0 +FA-LDDMM transformation of the brain (A) and the ventricle ...

Discussion

Choice of template

Several studies have measured registration accuracy for white matter structures (Grachev et al.,1999; Ardekani et al., 2005) using T1-weighted images for brain normalization. These approaches, in which there is only limited contrast within white matter, have shown that most white matter structures can be accurately registered within 3 mm by using linear and non-linear registration tools that are widely available, such as AIR and SPM. In our previous report using DTI data, the young adult population was normalized to a population-averaged DTI template (ICBM-DTI-81) (Mori et al., 2008), and we also found that 90% of the landmarks were within 3 mm. In the present study, we measured approximately 5 mm/90% as the normalization quality by affine transformation (Fig. 7A). This discrepancy is most likely due to the choice of the template. The choice of a population-averaged atlas, such as ICBM-DTI-81 as a template, ensures maximum matching between the template and linearly normalized subject data. Therefore, if linear transformation is used, the population-averaged atlas is a preferable choice as a template.

In this study, we employed a single-subject “Eve” atlas (JHU-DTI-MNI). There are several important points we would like to make regarding the choice of this template. First, highly non-linear transformation methods may not work properly with a population-averaged template, such as ICBM-DTI-81, in which the anatomical definition is blurred due to averaging. This is not a substantial issue for linear normalization, which is mostly driven by a large contrast change at the outside boundary of the brain, but the blurred internal structures could easily confuse highly non-linear transformation. Our initial testing of LDDMM led to overly inflated white matter when a population-averaged map is used as a template. Second, if a single subject is used as a template, the choice of the template could be an important issue. If the normalization algorithm is perfect, the template simply serves as the origin of coordinates to measure anatomical variability and the location of the origin may not be important as long as we are interested in differences among groups. However, in reality, it is preferable that the template image is similar to subject images so that the algorithm is not trapped by improper local minima during the transformation process. This issue is also related to the question of whether we should use an age-matched template. Our present study cannot address this issue completely. However, it is encouraging that the white matter of both the normal adult and AD patients were registered to the single-subject JHU-DTI-MNI atlas with excellent accuracy.

Choice of non-linear transformation parameters

Unlike linear transformation, non-linear transformation results are strongly dependent on transformation parameters. For LDDMM, the α/γ ratio determines the elasticity. The higher the elasticity, the more complicated shape differences could be matched (Fig. 4). However, the elasticity also has a higher probability of being trapped to a biologically wrong solution (Fig. 5).

In LDDMM, the α/γ ratio determines the amount of elasticity. It is not straightforward to determine which ratio is correct. The simulations in Figs. 1, ,3,3, ,44 and and55 should provide some ideas about the impact of this ratio on the image transformation and the ability of LDDMM to match the fine detail of the anatomy. Our empirical observation is as follows: 1) a ratio larger than 0.01 leads to conservative transformation; 2) for normal brain anatomy, for which affine transformation can match most deep white matter structures within 5 mm, single LDDMM transformation with a ratio=0.005 should provide the registration accuracy reported in this paper; and 3) for severely abnormal anatomy, it is important to employ cascading LDDMM transformation with a decreasing α/γ ratio, as reported in this paper. The downside of the cascading approach is the longer computation time.

Normalization accuracy by LDDMM

As expected, the affine transformation led to poor registration for both white matter and the ventricle of AD patients. To improve the normalization quality, non-linear transformation is needed. The results of LDDMM normalization are strongly dependent on image contrasts (Fig. 6). When b0 images are used to drive LDDMM, it has little impact on white matter registration because of the lack of contrasts within the white matter. b0 images, on the other hand, can match the brain and ventricle surfaces much more precisely. The opposite is true when FA maps are used for LDDMM. Namely, white matter structures are better registered while the brain surfaces are poorly matched. By combining these two contrasts using the two-contrast approach, the registration quality improves within the measurement errors both for the young adult and AD patient groups. Because these two images have strong contrasts in mutually complementary locations, the equal weighting factor, σc, for FA and b0 maps, seem a reasonable choice. In the future, if more contrasts are added, non-equal weighting may need to be examined.

Even if the normalization is accurately performed, it is always advised to exercise careful interpretation of image analysis results after normalization, especially when the data contain severe atrophy and a large amount of deformation is necessary (and therefore more pixel interpolation and partial volume effects). For example, excessive pixel interpolation may lead to loss of FA which reflects tissue atrophy, but not axonal or myelin damage.

In this paper, we did not perform a comparison with other widely used normalization methods such as SPM. In a recent publication by Klein et al. (2009), highly elastic registration tools were compared with the default setting of SPM and as expected, the former tools in general delivered improved registration accuracy. In this paper, we demonstrated the impact of the parameters that control the transformation elasticity. Also demonstrated is the potential pitfall of employing too elastic parameters and importance of gradual increase of elasticity during the transformation (Fig. 5). These points exemplify difficulties of comparing different registration tools because the results may vary depending on employed parameters and how the tools are used. It is therefore important to know how each registration tool works and carefully evaluate the impact of parameters. The multi-contrast approach and the cascading elasticity control could be combined with other registration tools to improve their accuracy.

Scalar- and tensor-based normalization

There are two unique issues for DTI normalization. First, when we choose contrasts to drive LDDMM, it is possible to use tensor orientation information. The second issue is the transformation of the tensor field, which requires the reorientation of the first eigenvector. The former issue has not been discussed in this study, where we used scalar images (b0 and FA images) to drive LDDMM. By inspecting DTI-derived images, it becomes immediately apparent that off-diagonal tensor element and eigenvector maps carry far stronger contrasts for intra-white matter structures, and thus, they could be more effective to drive normalization. It would be interesting in future work to see further improvement of normalization quality by adding these contrasts to the multi-contrast LDDMM. In the past, several methods have been postulated to use these contrasts for the normalization of DTI (Park et al., 2003; Cao et al., 2005; Zhang et al., 2006). However, these contrasts are rotationally variant, and therefore, require constant re-calculation during deformation. As a result, the method by Cao et al. (2005), while it may further improve the registration accuracy, requires far longer computational time. It is, therefore, important to carefully evaluate the added benefit with respect to increased computational loads.

In conclusion, we evaluated an LDDMM-based normalization method by testing single- and two-contrast approaches using b0 and FA contrasts, and these were applied to young adults and AD patients. Based on manual landmark-based measurements as gold standard, we found that the two-contrast (FA+b0) approaches can register the entire brain with higher spatial accuracy. The impact of a parameter that controls the elasticity of the transformation (α/γ) was also demonstrated. Using our brain atlas with 1 mm spatial resolution, the α/γ ratio of 0.01 provides conservative non-linear transformation. For brains with normal anatomy, an α/γ = 0.005 should achieve a high degree of matching for the white matter structures. For brains with severe anatomical changes, the cascading approach, in a range of 0.01–0.001, is recommended.

Acknowledgments

The authors thank Mr. Joseph S. Gillen, Ms. Terri Brawner, Ms. Kathleen Kahl, and Ms. Ivana Kusevic for experimental assistance. This publication was made possible by grant numbers P41RR015241, U24RR021382, R01EB004130, MH071616, HL085343, and R01AG20012, and by the Johns Hopkins Alzheimer’s Disease Research Center (P50AG005146). The National Center for Research Resources (NCRR) is a component of the National Institutes of Health (NIH). Its contents are solely the responsibility of the authors and do not necessarily represent the official view of NCRR or NIH. This research was supported by NIH grants. The images acquired on Alzheimer’s patients were supported by a methods development grant from Glaxo-Smith-Kline. Equipment used in the study is manufactured by Philips. Dr. van Zijl is a paid lecturer for Philips Medical Systems. Dr. van Zijl is the inventor of technology that is licensed to Philips. This arrangement has been approved by Johns Hopkins University in accordance with its conflict of interest policies.

References

  • Alexander DC, Pierpaoli C, Basser PJ, Gee JC. Spatial transformations of diffusion tensor magnetic resonance images. IEEE Trans Med Imag. 2001;20:1131–1139. [PubMed]
  • Ardekani BA, Guckemus S, Bachman A, Hoptman MJ, Wojtaszek M, Nierenberg J. Quantitative comparison of algorithms for inter-subject registration of 3D volumetric brain MRI scans. J Neurosci Methods. 2005;142:67–76. [PubMed]
  • Ashburner J, Friston KJ. Voxel-based morphometry—the methods. NeuroImage. 2000;11:805–821. [PubMed]
  • Basser PJ, Mattiello J, LeBihan D. Estimation of the effective self-diffusion tensor from the NMR spin echo. J Magn Reson B. 1994;103:247–254. [PubMed]
  • Beg MF. PhD thesis. The Johns Hopkins University; Jul, 2003. Variational and computational methods for flows of diffeomorph-isms in image matching and growth in computational anatomy.
  • Beg MF, Miller MI, Trouve A, Younes L. Computing large deformation metric mappings via geodesic flows of diffeomorphisms. International Journal of Computer Vision. 2005 February;61:139–157.
  • Bookstein FL. “Voxel-based morphometry” should not be used with imperfectly registered images. NeuroImage. 2001;14:1454–1462. [PubMed]
  • Cao Y, Miller MI, Winslow RL, Younes L. Large deformation diffeomorphic metric mapping of vector fields. IEEE Trans Med Imag. 2005;24:1216–1230. [PubMed]
  • Catani M, Howard RJ, Pajevic S, Jones DK. Virtual in vivo interactive dissection of white matter fasciculi in the human brain. NeuroImage. 2002;17:77–94. [PubMed]
  • Grachev ID, Berdichevsky D, Rauch SL, Heckers S, Kennedy DN, Caviness VS, Alpert NM. A method for assessing the accuracy of intersubject registration of the human brain using anatomic landmarks. NeuroImage. 1999;9:250–268. [PubMed]
  • Granander U, Miller MI. Computational anatomy: an emerging discipline. Stat Comput Graphics Newsl. 1996;7:3–8.
  • Jiang H, van Zijl PC, Kim J, Pearlson GD, Mori S. DtiStudio: resource program for diffusion tensor computation and fiber bundle tracking. Comput Methods Programs Biomed. 2006;81:106–116. [PubMed]
  • Jones DK, Griffin LD, Alexander DC, Catani M, Horsfield MA, Howard R, Williams SC. Spatial normalization and averaging of diffusion tensor MRI data sets. NeuroImage. 2002;17:592–617. [PubMed]
  • Klein A, Andersson J, Ardekani BA, Ashburner J, Avants B, Chiang MC, Christensen GE, Collins L, Hellier P, Song JH, Jenkinson M, Lepage C, Rueckert D, Thompson P, Vercauteren T, Woods RP, Mann JJ, Parsey RV. Evaluation of 14 nonlinear deformation algorithms applied to human brain MRI registration. NeuroImage 2009 [PMC free article] [PubMed]
  • Makris N, Worth AJ, Sorensen AG, Papadimitriou GM, Reese TG, Wedeen VJ, Davis TL, Stakes JW, Caviness VS, Kaplan E, Rosen BR, Pandya DN, Kennedy DN. Morphometry of in vivo human white matter association pathways with diffusion weighted magnetic resonance imaging. Ann Neurol. 1997;42:951–962. [PubMed]
  • Miller MI, Christensen GE, Amit Y, Grenander U. Mathematical textbook of deformable neuroanatomies. Proc Natl Acad Sci. 1993a;90:11944–11948. [PubMed]
  • Miller MI, Christensen GE, Amit Y, Grenander U. Mathematical textbook of deformable neuroanatomies. Proc Natl Acad Sci U S A. 1993b;90:11944–11948. [PubMed]
  • Miller M, Banerjee A, Christensen G, Joshi S, Khaneja N, Grenander U, Matejic L. Statistical methods in computational anatomy. Stat Methods Med Res. 1997;6:267–299. [PubMed]
  • Mori S, Wakana S, Nagae-Poetscher LM, van Zijl PC. MRI Atlas of Human White Matter. Elsevier; Amsterdam, The Netherlands: 2005.
  • Mori S, Oishi K, Jiang H, Jiang L, Li X, Akhter K, Hua K, Faria AV, Mahmood A, Woods R, Toga AW, Pike GB, Neto PR, Evans A, Zhang J, Huang H, Miller MI, van Zijl P, Mazziotta J. Stereotaxic white matter atlas based on diffusion tensor imaging in an ICBM template. NeuroImage. 2008;40:570–582. [PMC free article] [PubMed]
  • Muller HP, Unrath A, Ludolph AC, Kassubek J. Preservation of diffusion tensor properties during spatial normalization by use of tensor imaging and fibre tracking on a normal brain database. Phys Med Biol. 2007;52:N99–N109. [PubMed]
  • Pajevic S, Pierpaoli C. Color schemes to represent the orientation of anisotropic tissues from diffusion tensor data: application to white matter fiber tract mapping in the human brain. Magn Reson Med. 1999;42:526–540. [PubMed]
  • Park HJ, Kubicki M, Shenton ME, Guimond A, McCarley RW, Maier SE, Kikinis R, Jolesz FA, Westin CF. Spatial normalization of diffusion tensor MRI using multiple channels. NeuroImage. 2003;20:1995–2009. [PMC free article] [PubMed]
  • Pierpaoli C, Basser PJ. Toward a quantitative assessment of diffusion anisotropy. Magn Reson Med. 1996;36:893–906. [PubMed]
  • Pierpaoli C, Jezzard P, Basser PJ, Barnett A, Di Chiro G. Diffusion tensor MR imaging of human brain. Radiology. 1996;201:637–648. [PubMed]
  • Pruessmann KP, Weiger M, Scheidegger MB, Boesiger P. SENSE: sensitivity encoding for fast MRI. Magn Reson Med. 1999;42:952–962. [PubMed]
  • Wakana S, Jiang H, Nagae-Poetscher LM, Van Zijl PC, Mori S. Fiber tract-based atlas of human white matter anatomy. Radiology. 2004;230:77–87. [PubMed]
  • Woods RP, Grafton ST, Holmes CJ, Cherry SR, Mazziotta JC. Automated image registration: I. General methods and intrasubject, intramodality validation. J Comput Assist Tomogr. 1998;22:139–152. [PubMed]
  • Xu D, Mori S, Shen D, van Zijl PC, Davatzikos C. Spatial normalization of diffusion tensor fields. Magn Reson Med. 2003;50:175–182. [PubMed]
  • Yushkevich PA, Zhang H, Simon TJ, Gee JC. Structure-specific statistical mapping of white matter tracts. NeuroImage. 2008;41:448–461. [PMC free article] [PubMed]
  • Zhang H, Yushkevich PA, Alexander DC, Gee JC. Deformable registration of diffusion tensor MR images with explicit orientation optimization. Med Image Anal. 2006;10:764–785. [PubMed]