3D tract-specific local and global analysis of white matter integrity in Alzheimer's disease.

Alzheimer's disease (AD) is a chronic neurodegenerative disease characterized by progressive decline in memory and other aspects of cognitive function. Diffusion-weighted imaging (DWI) offers a non-invasive approach to delineate the effects of AD on white matter (WM) integrity. Previous studies calculated either some summary statistics over regions of interest (ROI analysis) or some statistics along mean skeleton lines (Tract Based Spatial Statistic [TBSS]), so they cannot quantify subtle local WM alterations along major tracts. Here, a comprehensive WM analysis framework to map disease effects on 3D tracts both locally and globally, based on a study of 200 subjects: 49 healthy elderly normal controls, 110 with mild cognitive impairment, and 41 AD patients has been presented. 18 major WM tracts were extracted with our automated clustering algorithm-autoMATE (automated Multi-Atlas Tract Extraction); we then extracted multiple DWI-derived parameters of WM integrity along the WM tracts across all subjects. A novel statistical functional analysis method-FADTTS (Functional Analysis for Diffusion Tensor Tract Statistics) was applied to quantify degenerative patterns along WM tracts across different stages of AD. Gradually increasing WM alterations were found in all tracts in successive stages of AD. Among all 18 WM tracts, the fornix was most adversely affected. Among all the parameters, mean diffusivity (MD) was the most sensitive to WM alterations in AD. This study provides a systematic workflow to examine WM integrity across automatically computed 3D tracts in AD and may be useful in studying other neurological and psychiatric disorders. Hum Brain Mapp 38:1191-1207, 2017.

INTRODUCTION

Alzheimer's disease (AD) is the most common form of dementia and accounts for 60%–70% of dementia cases [Burns and Iliffe, 2009]. AD is a neurodegenerative disorder characterized by progressive impairment in cognitive functions, including memory, affection, and executive function. Worldwide, around 47 million people have AD or related dementias. The annual global cost of AD is now estimated to be $818 billion—equivalent to 1.1% of the entire world's gross domestic product (Alzheimer's Disease International, http://www.alz.co.uk/research/statistics). AD not only brings huge distress to patients and their families, but also exerts extraordinary financial burden worldwide.

Neuroimaging has a key role in differential diagnosis of AD, and may also be used in assessing treatment efficacy. Magnetic resonance imaging (MRI) can identify biomarkers of disease progression and is non‐invasive and widely accessible. In particular, structural MRI (sMRI) has been widely used to study brain aging and AD. Specifically, sMRI can quantify widespread volumetric atrophy in AD, including changes in cortical volume and thickness [Thompson et al., 2001, 2007] and changes in individual regions of interest (ROIs), such as the medial temporal lobe and hippocampus [Apostolova et al., 2010; Morra et al., 2009; Schuff et al., 2009], precuneus [Karas et al., 2007], and frontal cortex [Frisoni et al., 2002]. However, standard anatomical MRI is less sensitive to white matter (WM) changes due to lack of discernible anatomical features within the WM.

Diffusion‐weighted imaging (DWI) [Basser et al., 1994]—a variant of MRI—measures water diffusion in brain tissue and captures microscopic WM alterations not detectable with sMRI. DWI is also non‐invasive and can be used to examine WM integrity and abnormalities in various brain disorders, such as Parkinson's disease [Zhang et al., 2011], multiple sclerosis [Filippi et al., 2001], autism [Jin et al., 2015b, 2015c], AD brain trauma [Dennis et al., 2015a, 2015c], and dementia [Daianu et al., 2016; Li et al., 2013]. So far, three major analytical methods have been applied to DWI, including standard ROI analysis [Hanyu et al., 1999], voxel‐based analysis (VBA) [Rose et al., 2008], and tract‐based spatial statistics (TBSS) [Smith et al., 2006], to quantify WM degeneration. For AD, ROI analysis often calculates some summary statistics (e.g., mean) in pre‐defined ROIs. They reveal degeneration in the corpus callosum [Hanyu et al., 1999], posterior cingulate [Kantarci et al., 2001], and WM of the frontal, temporal, parietal, and occipital lobes [Bozzali et al., 2002; Kantarci et al., 2001]. However, major drawbacks of ROI analysis include being a manual process; difficulty in identifying meaningful ROI boundaries, and especially WM fiber tracts, making statistical results from ROI analysis less stable; a partial volume effect in large or elongated ROIs (such as the fornix) that can lead to artifacts. Moreover, a stringent, but fundamental assumption of the ROI analysis method is that the diffusion properties in all voxels of the same ROI are essentially homogeneous, which is not necessarily true for DWI, or in a given disease.

Another approach is to perform voxel‐based analysis (VBA) after all subjects’ scans are aligned into a common atlas space; diffusion measures are then statistically compared on a voxel‐to‐voxel basis. VBA is automated and allows whole‐brain analysis. It can reveal more specific locations of effects. In AD, VBA has revealed WM degeneration in the frontal, temporal and occipital lobes, corpus callosum, superior longitudinal fasciculus, internal capsule, cingulum, and fornix [Rose et al., 2008; Teipel et al., 2007; Xie et al., 2006]. However, two major drawbacks of VBA include the issue of aligning homologous regions across subjects, and the effects of smoothing on statistical results [Jones et al., 2005]. Existing registration methods do not explicitly model the underlying fiber architecture of WM structure, including the neural systems and circuits affected.

More recently, tract‐based spatial statistics (TBSS) has been proposed to address some limitations of VBA [Smith et al., 2006], although it also has limitations [Schwarz et al., 2014]. It projects all the fractional anisotropy (FA) maps of individual subjects onto a common mean FA tract skeleton and performs voxelwise analysis along the skeleton. In AD, TBSS analyses reveal numerous regions of lower FA, including the parahippocampal WM, cingulum, inferior and superior longitudinal fasciculus, corpus callosum, and fornix [Liu et al., 2011]. Nonetheless, TBSS is limited in its representation of the neural pathways and systems in the WM. Voxels more distant from tract centers contribute with lower weight to the average value projected on the skeleton, and changes in such locations are harder to detect [Zalesky, 2011]. In addition, effects can be artificially split over multiple locations, making results harder to interpret if findings at one location may be driven by voxels elsewhere [Zalesky, 2011].

In summary, the previous WM analysis techniques focused on fitting predictive models, such as general linear model (GLM), to spatially averaged diffusion parameters on a tract, that is, ROI analysis (based on scalar values), VBA (statistics at each voxel in the image), or TBSS (statistics based on projections of neighboring voxels to the skeleton). However, none of these methods can provide a 3D detailed profile that reflects local alterations in a particular tract associated with AD. Here, we use two novel techniques to develop a tract‐specific analysis (TSA) framework to achieve this goal. First, our newly developed autoMATE (“automated Multi‐Atlas Tract Extraction”) pipeline [Jin et al., 2013, 2014] framework is operator‐independent and automatically labels fiber tracts from whole‐brain tractography using multi‐atlas label fusion. Instead of relying on an average statistic for the entire tract, we use point‐to‐point fiber matching to compare and visualize statistics from corresponding fiber points in 3D. Second, in terms of statistical analysis, we implement a new statistical analysis method, called FADTTS (“Functional Analysis for Diffusion Tensor Tract Statistics”), to associate diffusion parameters, such as FA and mean diffusivity (MD), along fiber tracts with a set of covariates of interest, such as age, sex, and diagnostic status [Zhu et al., 2010, 2011]. FADTTS outpowers traditional GLM for identifying more significant areas. It explicitly accounts for the spatial correlation and smoothness between neighboring points along major fibers and does not treat them as isolated points as GLM does—that is, diffusion parameters are treated as a function of fiber point locations. With our proposed techniques, we can show 3D profiles of differences in diffusion parameters on the tracts themselves—instead of only scalar mean values or a skeleton line; therefore, differences can be localized to a specific region of a tract and detailed changes can be revealed. Meanwhile, the global effects of AD on the tracts can also be delineated with FADTTS.

Furthermore, unlike prior studies that include only a few tracts per cohort [Zhang et al., 2009], we perform a comprehensive analysis of 18 major anatomically well‐defined WM tracts. These 18 WM tracts include representative association tracts, projection tracts, commissural tracts between cortices, and tracts in the limbic system. Studying these 18 tracts may offer a more complete landscape of how AD affects WM in the brain. Our study also includes a relatively large cohort of 200 subjects in three diagnostic groups (healthy elderly normal controls—NC, mild cognitive impairment—MCI, and AD).

MATERIALS AND METHODS

Subjects and Image Acquisition

DWI data were downloaded from the Alzheimer's Disease Neuroimaging Initiative (ADNI) database (http://adni.loni.usc.edu). As a multi‐site longitudinal study, ADNI's goal is to detect AD at the earliest possible stage, and identify ways to track the disease through biomarkers. The DWI data started to be collected during the “ADNI GO” phase in 2010 and continued with ADNI2 in 2011. We analyzed the baseline DWI data from a total of 200 participants, from three clinical groups: 49 NC (23M/26F; mean age standard deviation: 74.2 years 6.2 years), 111 people with MCI (70M/41F; 72.3 7.3), and 40 AD patients (23M/17F; 75.5 9.0).

DWI data were acquired with a 3D echo‐planar imaging (EPI) sequence from 3‐Tesla GE Medical Systems scanners. Acquisition parameters were: voxel size 2.7 × 2.7 × 2.7 mm3, image dimension 128 × 128 × 59, scan time = 9 min. About 46 separate volumes were acquired for each scan: 5 b 0 volumes with no diffusion sensitization and 41 diffusion‐weighted volumes (b = 1,000 s/mm2).

Image Preprocessing and Tractography

The DWI images were skull‐stripped using the Brain Extraction Tool (BET) from FSL (http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/) [Smith, 2002]. Those skull‐stripped images underwent eddy‐current correction using FSL to adjust for distortion from the EPI acquisition sequence. FA, MD (mm2/s), axial diffusivity (AxD) (mm2/s), and radial diffusivity (RD) (mm2/s) images were then computed from the DWI volumes using the diffusion tensor model.

Whole‐brain tractography was performed with Camino (http://cmic.cs.ucl.ac.uk/camino/), an open‐source software toolkit for DWI processing. Arguably, diffusion tensor imaging (DTI) can accurately track the orientation of a single fiber population within a voxel, but q‐ball reconstruction based on high angular resolution diffusion imaging (HARDI) provides orientation diffusion functions (ODFs) that can describe the orientations of multiple fiber populations [Tuch, 2004]. Similarly, probabilistic tractography is based on the ODF rather than principal direction in deterministic tractography, resulting in more accurate estimates of tract directions at the fiber crossing regions [Descoteaux et al., 2009]. Here, the spherical harmonic (SH) representation of the q‐ball ODF was reconstructed and a maximum of 3 local ODF peaks were detected at each voxel. A probabilistic tracking method—the Probabilistic Index of Connectivity method (PICo) [Parker et al., 2003]—was used to generate whole‐brain tractography. Seed points were selected at voxels where FA values were greater than 0.3. In PICo, Monte Carlo simulation generated streamlines that started from the selected seed voxels and followed the voxel‐wise probability density function profile estimated from the local ODF maxima with 4th‐order Runge‐Kutta interpolation at a step size of 1 mm. The maximum fiber turning angle was set to 45°/voxel. Tracking stopped at any voxel with an FA value less than 0.2. The performance of this parameter set has been validated in our previous study [Jin et al., 2014].

Label Fusion Based Tract Clustering

To extract fiber tracts from the whole‐brain tractography, we used a clustering algorithm‐autoMATE‐proposed in our previous works [Jin et al., 2013, 2014]. To guide the clustering, 5 representative WM tract atlases that consisted of 18 major WM tracts were manually constructed from the NC group (3 male and 2 female participants). We first non‐linearly registered the FA map of the single‐subject “Eve” template [Oishi et al., 2009] to the FA maps of those five subjects. The 130 parcellated ROIs of the Eve template were warped to those five subjects with the deformation fields generated by registration. The 18 WM tracts were then extracted with the ROIs based on a look‐up table [Zhang et al., 2010]. For example, the fibers that connected both occipital lobes were considered the occipital sub‐section of the corpus callosum. All the final WM tracts were manually edited to remove false positive streamlines. The names of the 18 WM tracts and their abbreviations are listed in Table 1.

Table 1

The names of the 18 WM tracts extracted and their abbreviations

No.	Tract name	Abbreviation
1	Left arcuate fasciculus	L‐ARC
2	Left anterior thalamic radiation	L‐ATR
3	Right anterior thalamic radiation	R‐ATR
4	Left cingulum	L‐CGC
5	Right cingulum	R‐CGC
6	Left corticospinal tract	L‐CST
7	Right corticospinal tract	R‐CST
8	Left inferior fronto‐occipital fasciculus	L‐IFO
9	Right inferior fronto‐occipital fasciculus	R‐IFO
10	Left inferior longitudinal fasciculus	L‐ILF
11	Right inferior longitudinal fasciculus	R‐ILF
12	Fornix	FNX
13	Corpus callosum, frontal	CC‐FNR
14	Corpus callosum, pre‐central gyrus	CC‐PRCG
15	Corpus callosum, post‐central gyrus	CC‐POCG
16	Corpus callosum, parietal	CC‐PAR
17	Corpus callosum, temporal	CC‐TEM
18	Corpus callosum, occipital	CC‐OCC

Similarly, for each new subject, we warped the 130 ROIs from the Eve atlas to the subject's space. The 18 WM tracts that traversed the corresponding ROIs were extracted. Then, we registered the FA map of each manually constructed atlas to the FA map of the subject and warped its WM tracts to the subject's space with the deformation fields. To refine the result, for each particular WM tract, we only kept the fibers whose Hausdorff distance to at least one of the warped atlas fibers was within a pre‐defined threshold (15 mm). For any pair of fibers and , the Hausdorff distance is defined as: , where is the Euclidean norm and the ordered pair ( , ) indicates an asymmetric distance from to , where and are the coordinate points along the fiber γ and the fiber γ [Gerig et al., 2004]. After the filtering, only the fibers whose locations and shapes were similar to the atlases were left.

For each of the 5 atlases, we obtained a set of the candidate fibers for each new subject for a particular tract. Each subject then had 5 sets of candidate fibers. Due to the variability of individual atlases, we used a label fusion scheme to combine the 5 sets. We ranked the “candidates” based on the average distances from the candidate fibers to the specific tracts of the 5 atlases. Only those fibers that ranked in the top percentage (varied per each tract) were kept [Jin et al., 2014]. Unlike using a single atlas, the label fusion scheme that combines multiple atlases can tune the final result by removing outliers in one set of candidate fibers or adding back the missing fibers in another set.

Fiber Matching

After all WM tracts were extracted, we established a fiber correspondence for each tract across the population to perform group studies. In our fiber matching scheme [Jin et al., 2014], we first selected a representative sample for each of the 18 WM tracts from the NC group (the degeneration of their WM was the least). To ease computation, we resampled the fiber points of each tract to 20 points per fiber in equal arc length distance. For each subject, the resampled representative set of tracts was warped to its space with registration. For each point on a warped representative tract, we searched its neighborhood of 10 mm. If there were fibers intersecting the neighborhood, the closest projection point on those fibers (defined as the point that had the shortest distance to the representative point) would be taken as the corresponding point of that representative point in that subject. If there were no intersecting fibers, which was very likely in some subjects due to their age and clinical condition, the warped representative point itself was considered as the corresponding point. Finally, the FA, MD, AxD, or RD values at those corresponding points over the population were interpolated as the measures used for group studies.

Statistical Analysis

We performed statistical analysis using the 3D FA/MD/AxD/RD profiles for each of the 18 WM tracts over the three groups (NC, MCI, and AD), respectively. Here, we used FADTTS, a varying coefficient model that delineates the variability structure of diffusion parameters along the WM fiber tract associated with a set of covariates, such as diagnostic status [Zhu et al., 2010, 2011]. This method can outperform the GLM in terms of both statistical power and prediction accuracy [Zhu et al., 2011]. Mathematically, we treated the diffusion parameters (e.g., FA and MD) along fiber tracts as mathematical functions of the fiber positions. For a particular WM tract in the th subject, denotes the arc length for the th point on the fiber relative to one of the end points and is its associated FA, MD, AxD, or RD value, where ; N, and . Here, is the total number of fiber points, is the subject number, and is the arc length of the fiber. Then the varying coefficient model can be described as follows: where is a vector of the covariates of interest for the subject and included age, sex, and diagnostic status, while represents the measurement error. For diagnostic status, a dummy variable is adopted to indicate whether the subject belongs to one specific diagnostic group, that is, NC, MCI, or AD. and are coefficient functions at arc length , where describes the average curve across all subjects and , , and are the coefficients of age, sex, and diagnostic status, respectively. Moreover, characterizes both individual curve variations from and the spatial correlation between and for different and . That is, measures both subject‐specific variability and location‐specific variability.

In this study, we assessed the extent of WM fiber degeneration across different clinical groups associated with diagnostic status. We investigated the problem at two levels—both locally and globally. Locally, we performed a groupwise comparison between three groups at each fiber point for a particular WM tract, while at the global level, we conducted the group comparison on each fiber. At the local level, we separated the cases by the direction of change of the parameters for the compared group relative to the baseline group in each comparison. Such questions can be formulated into the following two sets of hypotheses:

Equations (2) and (3) describe the local level comparison, whereas Eq. (4) is the global level comparison. For the th fiber, we then defined the local test statistic and the global test statistic [Zhu et al., 2011] as: where and are the estimates of and , respectively, where is the variance of η( . Zhu et al. [2011] have shown that those statistics have a desirable asymptotic distribution, namely, a weighted distribution. We then used a wild bootstrap method to effectively compute the ‐values of those statistics [Zhu et al., 2012]. Unlike the traditional bootstrap, we kept the independent variables (age, sex, and diagnostic status) at their original sample values, but generated the response variable (FA/MD/AxD/RD) based on the residual calculated from the fitted model under the null hypothesis.

RESULTS

Tract Visualization

Figure 1 shows the tract labeling results from the whole‐brain tractography, where the 18 WM tracts of a representative subject from each group (NC, MCI, and AD) were extracted. WM alterations can be noted in the gradually decreased number of fibers from NC to MCI to AD in the genu of CC in the top view.

Figure 1

Representative clustering results for 18 WM tracts in individuals from the three groups (NC, MCI, and AD). Left side, back, and top views of tract overlays are displayed. [Color figure can be viewed at http://wileyonlinelibrary.com]

Local Tract Profiles

At the local level, we tested Eqs. (2) and (3) and detected differences in FA, MD, AxD, and RD, associated with AD at each individual fiber point across groups. Figures 2, 3, 4, 5 show the 3D local profiles of these 4 parameters of the 18 WM tracts in the group comparisons for MCI versus NC (NC—baseline), AD versus MCI (MCI—baseline), and AD versus NC (NC—baseline), respectively. The – P‐values correspond to the color bars. Redder colors indicate greater differences. We corrected for multiple comparisons across all points on each tract by using the false discovery rate (FDR) [Benjamini et al., 2001]. Group differences were considered statistically significant at points with –log values > 1.3 (FDR corrected P = 0.05). We marked the significant areas with positive change in the parameter with the “+” sign and negative change with the “–” sign. For example, Figure 2 presents both the increased (“+”) and decreased (“–”) FA profiles of tracts for all three comparisons.

Figure 2

The 3D FDR corrected color maps of the 18 WM tracts reveal differences in FA values in pairwise comparisons of the three diagnostic groups. [Color figure can be viewed at http://wileyonlinelibrary.com]

Figure 3

The 3D FDR corrected color maps of the 18 WM tracts reveal differences in MD values in pairwise comparisons of the three diagnostic groups. [Color figure can be viewed at http://wileyonlinelibrary.com]

Figure 4

The 3D FDR corrected color maps of the 18 WM tracts reveal differences in AxD values in pairwise comparisons of the three diagnostic groups. [Color figure can be viewed at http://wileyonlinelibrary.com]

Figure 5

The 3D FDR corrected color maps of the 18 WM tracts reveal differences in RD values in pairwise comparisons of the three diagnostic groups. [Color figure can be viewed at http://wileyonlinelibrary.com]

All 18 tracts were affected by AD but at different segments of the tracts. As to the FA profile (Fig. 2), in MCI versus NC, only small differences in FA were detectable in most tracts. As the disease progressed, alterations became widespread in certain regions. For example, in projections of the bilateral ATRs, all the segments of the CC, the posterior portion of the bilateral CGCs, the crossing regions of the bilateral CST with the CC, and the posterior part of the bilateral IFOs and the L‐ILF. The entire FNX was highly affected by AD. Multiple regions in the fimbria and body of the FNX showed degeneration even at the MCI stage and it spread to the entire tract at the AD stage. FA values decreased in most of the affected regions from NC to MCI to AD, while FA values of certain fiber crossing regions, such as the crossing regions of the bilateral CSTs with the CC, the middle section of the L‐ARC, and the posterior portion of the bilateral IFOs actually became higher. The 3D MD profiles (Fig. 3) showed greater differences than those for FA in most tracts, especially when comparing AD versus NC, including the L‐ARC, the bilateral ATRs, most of the CC segments, the posterior CGCs, the FNX, and the bilateral ILFs. In contrast, MD values increased in those significant regions. The AxD (Fig. 4) and RD (Fig. 5) profiles showed similar alterations in those tracts as those of MD when the disease progressed; however, the differences were not as obvious between groups as those in MD.

Figure 6a–c shows the percentage of points across the 18 tracts whose FA/MD/AxD/RD values showed statistically significant differences (FDR corrected P < 0.05). Comparing MCI with NC (Fig. 6a), FNX (54%), CC‐TEM (27%), and CC‐POCG (25%) were the top three ranked tracts with the highest percentages of decreased FA; FNX (72%), L‐CST (64%), and L‐ARC (58%) were the top three for the increased MD. Similarly, the same tracts, that is, FNX, L‐CST, and L‐ARC, were again the top three for the increase of AxD, yet with lower percentages, that is, 50%, 45%, and 43%, respectively. As to the elevated RD, FNX still showed the most change (68%), while L‐IFO (28%) and CC‐PRCG (22%) were ranked the second and the third, respectively.

Figure 6

The percentages of points in the 18 WM tracts with statistically significantly lower FA and elevated MD/AxD/RD between (a) MCI versus NC, (b) AD versus MCI, and (c) AD versus NC (the latter group is the baseline), with FADTTS (FDR corrected P < 0.05). The tracts are ranked in order from the highest percentage to the lowest for MD. [Color figure can be viewed at http://wileyonlinelibrary.com]

Overall, the percentages of significant differences were lower in the comparison of MCI versus AD (Fig. 6b). L‐ILF (32%), CC‐OCC (28%), and CC‐POCG (26%) were ranked the top three tracts for the negative change of FA, whereas FNX (38%), CC‐POCG (35%), and CC‐FRN (33%) were listed as the top three for the positive change of MD. For AxD, L‐ARC (49%), FNX (29%), and R‐CST (25%) topped the list for the elevated values, while CC‐POCG (27%), FNX (24%), and L‐ATR (24%) did in RD, respectively.

It is not surprising that the comparison of AD versus NC showed the largest percentages of significant group differences (Fig. 6c). Specifically, FNX (70%), CC‐OCC (54%), and CC‐POCG (47%) were ranked the top three for lower FA values, whereas FNX (91%), L‐ARC (84%), and CC‐PAR (62%) showed the higher MD values. AxD and RD showed overall less change than MD with the most positive changes at FNX (89%), L‐ARC (81%), and L‐CST (50%) for AxD and FNX (90%), CC‐POCG (54%), and CC‐PRCG (51%) for RD, respectively. Again, the spatial extent of significant differences in the MD profiles was higher than the other three profiles for most tracts.

The largest difference in spatial extent of significant group differences between the MD and FA results was 54% in the L‐ARC between MCI versus NC, 24% in the FNX between AD versus MCI, and 76% in the L‐ARC between AD versus NC, respectively. There were no detectable differences in FA in the R‐ILF between AD versus NC, whereas 47% of the tract showed the increase in MD. The only exception was the CC‐OCC whose percentages of detected differences in FA were similar to those based on MD. Compared with MD, AxD and RD showed less alteration for most of the tracts across the three group comparisons, except in L‐ARC between AD versus MCI, AxD detected larger extent of difference than MD (49% vs. 31%).

The profiles of other directions, such as higher FA or lower MD, AxD, and RD from NC to MCI to AD showed little change in most of the tracts, except that the middle sections of bilateral CSTs and L‐ARC and the posterior sections of bilateral IFOs showed higher FA between MCI versus AD and between NC and AD (Fig. 7a–c). We also saw some small isolated regions with higher FA or lower MD values, which may be due to noise or sample size, and may not be maintained as the cohort becomes larger.

Figure 7

The percentages of points in the 18 WM tracts with statistically significantly elevated FA and lower MD/AxD/RD between (a) MCI versus NC, (b) AD versus MCI, and (c) AD versus NC (the latter group is the baseline), with FADTTS (FDR corrected P < 0.05). The tracts are ranked in order from the highest percentage to the lowest for FA. [Color figure can be viewed at http://wileyonlinelibrary.com]

Global Tract Profiles

At the global level, we tested Eq. (4) and estimated the differences in WM integrity on individual fibers between groups. Figure 8a–c lists the percentages of fibers in each of the 18 tracts that showed significant differences in the group‐wise comparisons for FA, MD, AxD, and RD, respectively (FDR corrected P < 0.05). In general, the percentages of detected differences were higher at the fiber level than the individual point level for the 4 parameters. For example, the percentages increased from 44% to 96% in the L‐ATR and from 54% to 81% in the R‐ATR for MD between AD versus NC, respectively. However, there were a few exceptions. For instance, the percentages decreased from 6% to 1% in the L‐CGC and from 5% to 0.2% in the R‐CGC for FA between MCI versus NC, respectively. Compared with FA, MD was the more sensitive diffusion parameter in terms of difference detection ‐ the percentages were much higher for most of the tracts at the global level as well. However, AxD and RD were very close to MD in terms of sensitivity at this level, compared with the differences at the local level.

Figure 8

The percentages of fibers in the 18 WM tracts with statistically significant FA/MD/AxD/RD differences between (a) MCI versus NC, (b) AD versus MCI, and (c) AD versus NC, with FADTTS (FDR corrected P < 0.05). The tracts are ranked in order from the highest percentage to the lowest for MD. [Color figure can be viewed at http://wileyonlinelibrary.com]

Table 2 lists the global FADTTS statistics that averaged all the local statistics calculated with Eq. (5) throughout the tract for all the 18 tracts. To compare statistics across the tracts, the local statistic at each point of a fiber for a particular tract was divided by the median value on that fiber. The scaled mean statistics can be taken as another measure to represent the effects of AD on the tracts. The larger the global FADTTS statistic the more prevalent differences in the tract. For example, between AD versus NC, the top three tracts that had the greatest global FADTTS statistics were FNX (12.3), CC‐PRCG (7.6), and CC‐POCG (7.4) for FA and FNX (19.5), L‐ARC (12.1), and CC‐PAR (9.0) for MD, respectively. These results were in line with the percentages of significant local differences for each tract shown in Figure 6a–c. We also noticed that the mean FADTTS statistics for MD were mostly larger than those for FA, AxD, and RD.

Table 2

The scaled global FADTTS statistics of the 18 WM tracts for FA/MD/AxD/RD

FADTTS	MCI vs. NC				AD vs. MCI				AD vs. NC
	FA	MD	AxD	RD	FA	MD	AxD	RD	FA	MD	AxD	RD
L‐ARC	1.0	6.6	5.3	2.7	2.9	4.1	6.4	2.2	2.8	12.1	16.9	3.2
L‐ATR	2.0	2.7	2.1	2.9	1.9	3.8	3.3	3.8	3.1	7.3	6.1	7.4
R‐ATR	2.1	2.8	2.1	2.7	1.7	3.9	3.6	3.7	2.7	5.5	4.3	5.6
CC‐FRN	2.2	2.2	2.0	2.2	2.2	3.7	2.5	3.6	4.6	5.8	4.1	6.1
CC‐PRCG	3.3	5.9	4.9	6.2	6.8	4.0	3.6	3.6	7.6	7.9	6.4	8.1
CC‐POCG	4.6	3.2	2.6	3.7	3.3	4.2	3.1	4.6	7.4	8.0	5.1	9.5
CC‐PAR	2.5	3.4	1.8	3.6	2.3	3.3	3.1	3.1	4.6	9.0	5.7	8.3
CC‐TEM	3.2	4.5	3.9	4.1	1.4	1.7	2.0	1.5	4.3	7.1	7.6	5.8
CC‐OCC	2.6	2.7	1.6	3.1	2.6	2.3	1.6	2.6	5.8	5.5	3.1	6.5
L‐CGC	1.2	4.1	3.1	2.5	1.2	3.1	2.4	2.2	1.5	7.4	5.4	4.2
R‐CGC	1.2	3.1	4.0	1.6	3.1	3.4	1.6	4.0	2.6	5.3	3.4	4.3
L‐CST	3.6	8.9	8.2	5.3	3.9	2.8	4.0	2.6	6.0	8.2	13.6	4.1
R‐CST	3.1	5.2	3.7	4.4	4.4	6.8	5.0	3.1	4.8	6.7	8.8	3.5
FNX	7.7	10.4	8.3	10.8	3.1	5.1	5.3	4.7	12.3	19.5	16.4	19.7
L‐IFO	6.4	4.2	5.6	5.4	4.0	2.4	2.2	3.0	6.9	6.4	7.1	10.2
R‐IFO	3.1	3.5	3.8	3.2	4.2	3.8	3.5	3.9	5.6	4.5	4.4	5.4
L‐ILF	1.4	3.1	3.1	2.5	3.7	2.7	2.4	2.7	4.5	6.2	2.4	5.9
R‐ILF	1.1	2.6	1.6	1.2	1.8	2.1	1.4	1.7	1.1	5.4	2.9	2.8

DISCUSSION

DWI measures water diffusion in WM fiber tracts and is now widely used to investigate the integrity of WM structures in vivo. With DWI, the traditional ROI analysis sums the values of diffusion parameters such as FA and MD of a specific region inside the WM tracts, which are typically “seeded” by an operator; in contrast, the ROIs of TSA are the fiber tracts automatically extracted and clustered in our approach; diffusion parameters are measured along the entire tracts. Compared with VBA, TSA relies on the underlying fiber architecture formulated by tractography. It can better isolate individual anatomical structures and avoid contamination resulting from mis‐registration and the arbitrary choice of smoothing incurred in VBA.

Our 3D TSA framework, which we refer to as “autoMATE,” further enhances traditional TSA by providing 3D profiles of diffusion parameters on automatically clustered WM tracts. The autoMATE framework extracts WM tracts automatically instead of requiring manual ROI placement as is typical in traditional TSA. To construct 3D profiles, we proposed a fiber matching scheme that matches the corresponding points on the tracts across the population for 3D point‐wise and fiber‐wise comparisons, while traditional TSA only compares the mean values of the extracted tracts between groups. TBSS [Smith et al., 2006] is a popular method that projects all subjects’ normalized FA data to a mean FA skeleton and then applies voxel‐wise statistics. It aligns FA maps from multiple subjects and avoids the arbitrariness of making a spatial smoothing choice when creating the mean skeleton. However, it loses a large amount of local information by reducing the entire tract to its skeleton. On the other hand, our TSA framework extracts the entire tracts and the 3D functional analysis can reveal detailed local alterations in WM. Although these two methods cannot be compared directly as their dimensionalities and domains of definition are different, our method can provide supplementary information on WM integrity that TBSS cannot provide. The two methods are illustrated in Figure 9. We can see that TBSS provided only approximate alteration locations on a skeleton at the whole brain scale (Fig. 9a, [Acosta‐Cabronero et al., 2010]), while our 3D TSA framework can reveal detailed local alterations in each tract (Fig. 9b).

Figure 9

(a) An illustration of the TBSS results for reduced FA in AD patients [Acosta‐Cabronero et al., 2010]. (b) An illustration of our 3D TSA framework results on three WM tracts: anterior thalamic radiation (ATR), corpus callosum (CC), and cingulum (CGC) in AD. [Color figure can be viewed at http://wileyonlinelibrary.com]

Furthermore, in most of the previous AD studies of WM degeneration, the authors focused on only one or a few tracts of interest such as the cingulate, inferior occipitofrontal fasciculus, and corpus callosum [Preti et al., 2012; Zhang et al., 2009]. In this study, we examined the systematic effects of AD on 18 major WM tracts and provided a comprehensive overview of the effects of AD on the integrity of WM. We also analyzed a relatively large cohort of 200 subjects to increase the statistical power of our findings. Sample sizes in prior studies were generally much smaller, usually with only tens of subjects, for example, 76 in [Preti et al., 2012] and 18 in [Zhang et al., 2009]. We conducted three‐group pairwise comparisons (NC, MCI, and AD) to show gradual WM degeneration as the disease progressed.

At present, the exact pathological mechanisms leading to WM degeneration in AD are not fully understood. Prior studies have reported widespread cortical GM loss in AD patients [Frisoni et al., 2007; Thompson et al., 2003]. Microstructural WM degeneration may be a secondary degeneration arising from GM pathology. It is an anterograde process, resulting from the degeneration within neurons, leading to loss of axons and myelin [Englund, 1998]. Oligodendrocyte death and reactive gliosis, considered as direct WM degradation rather than a consequence of primary GM damage, have also been reported [Sjöbeck and Englund, 2003]. From the neuroimaging perspective, FA, MD, AxD, and RD are the most commonly derived quantitative measures from DWI as the indicators for WM structural alterations. FA is a scalar value between 0 and 1 that measures anisotropic water diffusion, which can reflect the degree of the directionality of cellular microstructures inside WM tracts [Chua et al., 2008]. Anisotropy increases as FA increases from 0 to 1. MD is a scalar that measures the mean water diffusion at each voxel. Because abnormalities in the cellular microstructures of WM decrease the barriers to free diffusion of water, MD increases and diffusion becomes more isotropic, resulting in lower FA [Stebbins and Murphy, 2009]. AxD is the principal eigenvalue of diffusion tensor, while RD is the average of the second and the third eigenvalues. It is believed that AxD is related to axon degeneration and RD seems to be modulated by myelin [Song et al., 2002]. Therefore, the changes of these 4 parameters between normal aging people and MCI or AD patients can reveal the important information about tract integrity in MCI and AD patients. It is also noticeable from Figures 2 to 7 that compared to the other 3 parameters, that is, FA, AxD, and RD, MD may be a more sensitive indicator for catching WM degeneration for most of the tracts, which was also suggested in earlier studies [Nir et al., 2013].

Figures 6 and 7 indicate widespread WM degeneration caused by AD across all the 18 major WM tracts examined, based on the percentages of points or fibers whose FA/MD/AxD/RD shows statistically significant change between groups with FADTTS. In the limbic network, lower FA and higher MD/AxD/RD were found in the FNX and the CGC. The FNX is most affected by AD in terms of diffusion parameters among all the 18 WM tracts. The FNX is widely associated with AD and is the main WM pathway emerging from the hippocampus, a region well known for its atrophy in AD [Apostolova et al., 2010]. The mean FA and MD of the FNX show a significant correlation with hippocampal volumes and may be useful in predicting disease progression [Mielke et al., 2012]. These measurements correlate with clinical scores [Jin et al., 2015a]. For the bilateral CGCs, compared to NC, the MCI group shows little change in terms of FA. As the disease progresses, the middle and posterior portions demonstrate greater alterations and the trend of alterations becomes more obvious in the MD/AxD/RD maps. Prior studies reported a close connection between hippocampal volume reduction and the FA/MD values of the posterior portions of the CGC in AD [Nakata et al., 2009; Zhang et al., 2007], which may lead to secondary degeneration in the tract.

We also examined three long cortico‐cortical association tracts including the ARC, the IFO, and the ILF. They all show regions of lower FA and higher MD/AxD/RD in AD versus NC. The L‐ARC connects Broca's area, located in the left posterior frontal lobe, as well as Wernicke's area, located in the left posterior temporal lobe. Due to high biological variation of ARC and the limitation of current DWI technology [Catani et al., 2007; Guevara et al., 2012], the right ARC may be too small to be fully extracted in some of the subjects. Therefore, we did not include it in our study as it would make a group comparison difficult. For the left ARC, FA seems to decrease more significantly near the temporal lobe region, while MD/AxD increases along the entire tract as the disease progresses. As both Broca's area and Wernicke's area are associated with speech processing, degeneration of the L‐ARC may contribute to language impairment in AD [Ferris and Farlow, 2013]. The IFO connects the orbitofrontal areas and the occipital and posterior temporal lobes, whereas the ILF connects the anterior temporal lobe and the occipital lobe. They have been suggested to be related to object recognition, visual emotional processing, and visual memory [Chanraud et al., 2010]. Greater differences (redder areas in the FA/MD maps in Figs. 2 and 3) occur in the regions closer to the occipital lobes in both the ILFs and the IFOs and the middle portions of the IFOs, which may suggest cognitive decline in AD pathology [Kitamura et al., 2013]. In our findings, the IFO and the ILF on the left side seem to show more extensive differences than those on the right side.

The CC is the largest commissural tract with extensive connections between the two hemispheres. Therefore, it is not surprising to see its progressive atrophy in AD [Teipel et al., 2002]. Prior studies have reported that the genu [Xie et al., 2006] or splenium [Rose et al., 2000] of the CC are adversely affected by AD. Here, our findings demonstrate that the affected regions are widespread across the entire CC at each of the six segments we defined. The alterations are not so obvious between MCI versus NC, but they become very pronounced between AD versus NC, especially in the MD map. In fact, the spatial extent of group differences in the CC remains the top of the list in Figure 6 among the 18 WM tracts that we investigated.

We examined two projection tracts, the ATR and the CST. The ATR consists of fibers connecting between the thalamus and the frontal cortex. The ATR showed significantly lower WM integrity in the distal projections closer to the frontal cortex in AD, while the degree was much less in both NC and MCI. Alterations in the ATR in AD patients may lead to further cognitive impairment, as the disease progresses. Prior studies also examined poor visual evaluation and reduced anisotropy of the ATR in AD patients [Niida et al., 2013]. As the CST is involved in control of motor function of the body, only a few studies have analyzed the tract so far [Douaud et al., 2011; Teipel et al., 2007]. Here, our results show significant differences in the crossing fiber regions between the CC and the CST [Wang et al., 2016] with both higher FA and higher MD. The higher FA suggests a more linear shape of diffusion tensor, perhaps indicating loss of crossing fibers in that region. We also notice the same phenomenon in the middle segment of the L‐ARC, and posterior portion of bilateral IFO, where crossing fibers are easily seen. Probabilistic whole‐brain tractography with PICo [Parker et al., 2003] was used to explicitly account for fiber crossings in those regions. The region of the CST near the cerebral peduncle shows lower FA and higher MD/AxD/RD, suggesting disruption associated with AD in that region as well.

There are several limitations of the proposed TSA approach. It is only useful when the tract of interest is already known and can be anatomically identified. Therefore, our analysis is limited to only those major WM tracts and misses some small tracts, such as short cortico‐cortical association fibers [Zhang et al., 2010]. The reason is that they are not anatomically well‐defined and not every subject has the “same” U‐fibers. In the future work, we will continue adding newly found tracts to our WM atlases to perform more comprehensive analysis. Another limitation of TSA is that the reconstructed tracts may fail to follow the true fiber tract pathways due to limitations of tractography, which include the choices of parameters such as turning angles, and stopping criteria [Dennis et al., 2015b; Zhan et al., 2013, 2015]. In future, more advanced tractography algorithms will be incorporated into our pipeline to increase the accuracy or predictive power of the TSA analysis.

CONCLUSION

We have presented a comprehensive WM analytical framework to conduct TSA of WM in AD. Unlike prior studies, our method allows a systematic assessment of diffusion parameters along 3D tract profiles associated with different clinical and genetic variables of interest (e.g., diagnosis and genetic markers). AutoMATE, automatically extracts anatomically meaningful WM tracts, and establishes a correspondence among fiber points across subjects in a specific study. Our statistical analysis method, FADTTS, captures the variability and associations in diffusion properties along the tracts associated with diagnostic status. We applied the framework to a large cohort of 200 subjects and performed a comprehensive data analysis on 18 WM tracts. Our framework will be further extended to study WM integrity in other neurological and psychiatric diseases.