MRI Analysis of White Matter Myelin Water Content in Multiple Sclerosis: A Novel Approach Applied to Finding Correlates of Cortical Thinning

A novel lesion-mask free method based on a gamma mixture model was applied to myelin water fraction (MWF) maps to estimate the association between cortical thickness and myelin content, and how it differs between relapsing-remitting (RRMS) and secondary-progressive multiple sclerosis (SPMS) groups (135 and 23 patients, respectively). It was compared to an approach based on lesion masks. The gamma mixture distribution of whole brain, white matter (WM) MWF was characterized with three variables: the mode (most frequent value) m₁ of the gamma component shown to relate to lesion, the mode m₂ of the component shown to be associated with normal appearing (NA) WM, and the mixing ratio (λ) between the two distributions. The lesion-mask approach relied on the mean MWF within lesion and within NAWM. A multivariate regression analysis was carried out to find the best predictors of cortical thickness for each group and for each approach. The gamma-mixture method was shown to outperform the lesion-mask approach in terms of adjusted R², both for the RRMS and SPMS groups. The predictors of the final gamma-mixture models were found to be m₁ (β = 1.56, p < 0.005), λ (β = −0.30, p < 0.0005) and age (β = −0.0031, p < 0.005) for the RRMS group (adjusted R² = 0.16), and m₂ (β = 4.72, p < 0.0005) for the SPMS group (adjusted R² = 0.45). Further, a DICE coefficient analysis demonstrated that the lesion mask had more overlap to an ROI associated with m₁, than to an ROI associated with m₂ (p < 0.00001), and vice versa for the NAWM mask (p < 0.00001). These results suggest that during the relapsing phase, focal WM damage is associated with cortical thinning, yet in SPMS patients, global WM deterioration has a much stronger influence on secondary degeneration. Through these findings, we demonstrate the potential contribution of myelin loss on neuronal degeneration at different disease stages and the usefulness of our statistical reduction technique which is not affected by the typical bias associated with approaches based on lesion masks.

1. Introduction

The mechanisms of cortical atrophy in multiple sclerosis (MS) are poorly understood. MS is typically characterized by the appearance of white matter (WM) lesions on T2 FLAIR magnetic resonance imaging (MRI) images. Although lesions are the most visible radiological aspect of the disease, used for both diagnosis (Polman et al., 2011) and clinical trials (Polman et al., 2006; Kappos et al., 2007; Mikol et al., 2008), disability has been shown to have a stronger correlation with gray matter (GM) cortical atrophy compared to WM measurements (Ramasamy et al., 2009; Calabrese et al., 2010, 2011; Narayana et al., 2013). It is then of particular interest to explore the mechanisms driving cortical thinning.

Wallerian degeneration as a result of chronic demyelination has been implicated as one potential mechanism for cortical GM loss (Peterson et al., 2001; Vercellino et al., 2005; Geurts and Barkhof, 2008), however the extent that demyelination contributes to this process has yet to be fully elucidated. Interestingly, it has been shown that the rate of GM atrophy is markedly different between the remitting relapsing (RR) and secondary progressive (SP) phases of the disease (Fisher et al., 2008). Therefore, we hypothesize that the relationship between WM demyelination and cortical thinning is different among these two disease stages.

Quantifying WM demyelination, especially myelin loss has been challenging, and several MRI sequences have been utilized to obtain an imaging biomarker of myelin integrity in MS, including diffusion MRI (Rovaris et al., 2005), magnetization transfer ratio (MTR) imaging (Horsfield, 2005) and myelin water fraction (MWF) imaging (Laule et al., 2004). Although promising, the clinical utility of MWF maps has been impeded by a prohibitively long acquisition time of the associated T2 relaxometry MRI sequence and challenging T2 data analysis (Kolind et al., 2009). To improve the efficiency of T2 data acquisition, we optimized a 3D T2prep GRE sequence (Nguyen et al., 2012), which is able to achieve full brain coverage within 10 min at 3T and feasible to integrate into a clinical imaging protocol.

WM damage is traditionally quantified globally from the “lesion load”: the total volume of visible lesions on a dedicated MRI sequence such as FLAIR. This approach has limitations as the normal-appearing (NA) WM was shown to have pathological changes (Evangelou et al., 2000), and similarly for “dirty WM” (Seewann et al., 2009) which includes voxels with intermediate intensity between lesion and NAWM. Additionally, while advanced MRI sequences provide metrics to quantify local WM damage, they still rely on the creation of associated lesion masks. These suffer from two main issues: there is inter-rater variability in their drawing or editing, and they are time-consumming to create, thus not adapted to large cohorts (Egger et al., 2017). Automated methods for mask generation alleviate somewhat both these issues, via high reproducibility and fast computer automation, respectively. However the associated algorithms are often calibrated to generate an output in agreement with a given experimenter and therefore are inherently biased. As a result, automated methods of lesion mask creation have around the same variability between themselves as between experimenters manually drawing or editing masks.

The aim of this work was 2-fold. First, to evaluate the association between cortical thickness and WM myelin, quantified from the measured MWF distribution in whole brain WM, and assess if this relationship differed between disease stages, that is between RRMS and SPMS patients. Second, to propose for the first time a model of the distribution of cerebral white matter MWF with a statistical reduction technique which circumvents the bias introduced by the use of lesion masks.

2. Materials and Methods

2.1. Subjects

Patient clinical and MRI data utilized for this analysis was part of a larger MRI and clinical database collection study at the Judith Jaffe Multiple Sclerosis Center, Weill Cornell Medicine, Department of Neurology, New York, USA. The MS database study was approved by Weill Cornell Medicine institutional review board and all subjects were consented for participation during their routine clinical visits. Informed and written consent for data collection and publication was obtained from study participants.

One-hundred and fifty-seven patients with the diagnosis of multiple sclerosis (Polman et al., 2005), including RRMS and SPMS were selected from our database. Patients were selected from our ongoing data collection, and informed consent was obtained from each subject. All RRMS and SPMS patients with a clinical MRI sequence that included the T2prep 3D spiral sequence were included in the analysis. Patient characteristics (Table 1) and clinical data were obtained within 3 months of the individual's brain MRI. The following clinical data was collected: gender, age, disease duration (DD) from initial symptom, Expanded Disability Status Score (EDSS), disease subtype, and duration of treatment with disease modifying therapy (DMT).

TABLE 1

Table 1. RRMS and SPMS group demographics and disease characteristics.

2.2. MRI Data Acquisition

Each participant underwent T1, FLAIR, and T2 relaxometry sequences with a 3Tesla GE scanner (HDxt 16.0) using 8-channel phased-array coil. The T1 sequence was an axial 3D inversion recovery fast spoiled gradient recalled echo (FSPGR) resulting in T1 volumes with 1.2 mm isotropic resolution. The T2 FLAIR volumes had 1.2 × 0.6 × 0.6 mm³ resolution. The T2 relaxometry sequence was a whole-brain T2prep 3D spiral as previously described by Nguyen et al. (2012) and included the following parameters: axial FOV = 24 cm; matrix size = 192 × 192 (interpolated to 256 × 256); slice thickness = 5 mm; number of slices = 32; sequence TR (time between subsequent T2prep pulses) = 2.5 s; spiral TR = 8.1 milliseconds (ms); spiral TE = 0.5 ms; flip angle = 10°; readout bandwidth = ±125 kHz; number of spiral leaves per segment = 64; 15 nominal T2prep times = 0, 5 ms, 10–40 ms (10 ms step), 60–140 ms (20 ms step), 180–300 ms (40 ms step); scan time = 10 min. A modified BIR-4 adiabatic pulse (De Graaf and Nicolay, 1998; Jenista et al., 2013) was used in the T2prep module to improve T2 weighting accuracy against increased B0 and B1 field inhomogeneities at 3T.

2.3. T2 Relaxometry Post-processing

Recently, we reported a multi-voxel spatial regularization approach to analyze multiexponential T2 decay data, with the goal of obtaining usable myelin maps from noisy but fast acquisitions (Kumar et al., 2012). The size and stability of the challenging minimization problem has been greatly reduced by our subsequent fitting method, called “Spatially constrained multi-Gaussian” algorithm, described in Raj et al. (2014). Briefly, it employs a non-linear model to recognize only three distinct water compartments contributing to the T2 spectrum, s, in each voxel: myelin water modeled by a gaussian G(μ_MW, σ_MW), intra/extra-cellular water modeld by another gaussian G(μ_CW, σ_CW) and cerebral spinal fluid modeled by a delta function δ(μ_CSF). In each voxel, s is then given by a linear combination of these contributions: s = c_MW G(μ_MW, σ_MW) + c_CW G(μ_CW, σ_CW) + c_CSF δ(μ_CSF). The MWF is defined as the ratio between the portion of the signal coming from myelin water, and the total signal from all water compartments: $MWF=\frac{c_{MW}}{c_{MW}+c_{CW}+c_{CSF}}$. The resultant MWF maps were shown to provide good spatial resolution and visual lesion discrimination (Figure 1). Note that the lesion-mask free method presented in this work and described in Section 2.6 does not depend on the method used to compute MWF.

FIGURE 1

Figure 1. Examples of T2 FLAIR images (left) from two MS patients (A,B) and corresponding MWF maps (right) demonstrate good spatial resolution and sensitivity to visually T2 hyperintense lesions. The color bar represents the range of MWF values. Hyperintense lesions on T2 Flair images are found to have low MWF values as compared to the surrounding NAWM. MWF also appears sensitive to dirty WM, as shown by the low MWF displayed in dark/light blue within the central white matter of (B).

2.4. T1 Post-processing and Cortical Thickness

FreeSurfer (http://freesurfer.net) was used to segment T1 images into GM and WM masks, which were visually checked and manually edited for misclassification due to WM T1-hypointensities associated with lesions. This editing was limited to ensuring that no WM was included in the GM mask. It was different from the lesion-mask editing procedure detailed in the next section (and applied to FLAIR contrast) which involved separating NAWM from every single lesion present within the global WM mask.

After completion, a trained neurologist (SG) gave a final approval on the resulting GM and WM masks. Cortical thickness measurements were then directly obtained with FreeSurfer from the corrected masks. The 1st echo of the T2 relaxometry sequence was linearly coregistered to the subjects T1 with FSL (http://fsl.fmrib.ox.ac.uk/fsl/fslwiki). The inverse of this transform was then applied to the WM mask to get it into the same space as the MWF maps.

2.5. FLAIR Post-processing and Lesion Masks

The white matter hyperintensity lesion masks were created from the T2 FLAIR images aligned onto FreeSurfer volume with the boundary-based registration (Greve and Fischl, 2009). The skullstripped FLAIR images were bias-corrected and segmented in three tissue classes with FSL FAST (Zhang et al., 2001) based on the image intensities. The third tissue class included both CSF tissue and WM lesions. Only voxels of the the bias-corrected FLAIR image belonging to that class and to the Freesurfer WM mask were kept as WM lesion candidate voxels. The associated volume was further thresholded to select only voxels with intensity greater than the 30th percentile of the FSL “robust intensity range” defined as the range of intensities between the 2nd and 98th percentile. The thresholded maps were finally binarized to provide the lesion masks.

The resulting WM lesion masks were overlaid on FLAIR images and manually edited to ensure an optimal match with the lesion hyper-intensities seen on the underlaid FLAIR contrast. A trained neurologist gave a final approval on these edits. The NAWM masks were defined as the voxels of the WM mask not belonging to the WM lesion mask. As a result the lesion and NAWM masks were mutually exclusive. The transform obtained described in the T1 post-processing section was used to obtain the lesion and NAWM masks in MWF space.

2.6. MWF Distribution

We analyzed MWF maps within WM masks for each patient using plots and descriptive statistics. The observed distributions of MWF values (at the patient level) appeared skewed and with one or two prominent peaks (Figures 2, 3). One of the most common skewed distribution is the gamma distribution, used in fields ranging from economics (Chatterjee and Chakrabarti, 2007; Chotikapanich and Griffiths, 2008) to meteorology (Smith, 2003; Yoo et al., 2005) as either a single distribution (Smith, 2003; Chatterjee and Chakrabarti, 2007) or mixture of two distributions (Yoo et al., 2005; Chotikapanich and Griffiths, 2008). The skewed distribution of our data also presented one or two peaks and exploratory analysis demonstrated a tight visual fit from a gamma model. These elements provided ground for describing the data distributions as either a gamma random variable or a mixture of two gamma-distributed random variables. In other words, we hypothesized that for each patient i, the MWF distribution p_i(MWF) followed a mixture of two gamma distributions of the form p_i(MWF) = λ_iΓ(MWF; α_{1, i}, β_{1, i})+(1−λ_i)Γ(MWF; α_{2, i}, β_{2, i}) where Γ denotes a gamma distribution, λ_i is the proportion of mixing between the two distributions, and the “shape” and “scale” pairs (α_{1, i}, β_{1, i}) and (α_{2, i}, β_{2, i}) are the parameters describing the shape and scale of each mixture component (Figure 2). In contrast, using simply the mean MWF would have lost the information contained in this distribution, especially the presence of two peaks.

FIGURE 2

Figure 2. Illustration of the same gamma mixture distribution with two different values of λ. m₁ and m₂ are the modes (most frequent values) of the low and high MWF components of the gamma mixture, respectively. λ is the gamma mixture mixing ratio and represents the proportion of low MWF to high MWF (when λ is equal to 0.5, the two components have the same contribution in the gamma mixture). To compare the relationship between lesion/NAWM masks and each gamma mixture component we defined ROI[m₁] and ROI[m₂] as all voxels with MWF value less than m₁ (in red) and greater than m₂ (in green), respectively. We assumed the lesion mask had more overlap with ROI[m₁] than ROI[m₂], and the NAWM mask more overlap with ROI[m₂] than ROI[m₁]. (Left) Distribution with a λ value less than 0.5 (exact parameters: m₁ = 0.06, m₂ = 0.17, λ = 0.25). (Right) Distribution with a λ value greater than 0.5 (exact parameters: m₁ = 0.06, m₂ = 0.17, λ = 0.75).

FIGURE 3

Figure 3. Gamma mixture fits of 3 RRMS subjects (left column) and 3 SPMS subjects (right column). Histograms are only represented for illustrative purposes and were not used during the fitting procedure.

The shape parameter controls the skewness of the distribution while the scale controls how spread the distribution is. The position of the peak in each component of the distribution, i.e., the distribution modes, are given by m₁ = (α₁−1) × β₁ and m₂ = (α₂−1) × β₂. Another reason we chose the gamma distribution as opposed to another one, is that the parameters are easier to interpret, and directly relate to explicit quantities such as the mean and modes of the distribution as just described. For the single gamma distribution, parameters were estimated with the maximum-likelihood (ML) method, while for the gamma mixture distribution they were estimated with the expectation maximization (EM) method (Dempster et al., 1977). A penalized likelihood ratio test was implemented to determine the goodness of fit of the mixture gamma compared to the single gamma model (Wong and Li, 2014) and choose between the two.

2.7. Interpretation of the Gamma Mixture Parameters in Terms of Lesion and NAWM

Our assumption was that the lesion mask was more associated to the low-MWF component at position m₁, compared to the high-MWF component at position m₂. Similarly, we expected the NAWM mask to be more associated to the high-MWF component than to the low-MWF component. To test the first assumption, we computed a metric of similarity between the lesion mask and a region of interest ROI[m₁] associated with the low-MWF component at m₁, and between the lesion mask and a region ROI[m₂] associated with the high-MWF component at m₂. We then tested if the metric of similarity was significantly higher between the lesion mask and ROI[m₁] than between the lesion mask and ROI[m₂]. The metric of similarity was chosen to be the Dice coefficient, commonly use in neuroimaging to compute the overlap of two ROIs (see e.g., Dayan et al., 2015). It is defined for two regions A and B as $Dice(A,B)=\frac{2\times |A\cap B|}{\left|A\right|+\left|B\right|}$ where |X| stands for the number of voxels in region X. ROI[m₁] was defined as the region with voxels having MWF value below m₁, and ROI[m₂] was defined as the region with voxels having MWF value above m₂ (Figure 2). We therefore tested if Dice(lesion mask, ROI[m₁]) > Dice(lesion mask, ROI[m₂]), using a Welch two sample t-test. Similarly, to test our second assumption that the NAWM mask is more associated to the high-MWF component than to the low-MWF component, we used a Welch two sample t-test to test that Dice(NAWM mask, ROI[m₂]) > Dice(NAWM mask, ROI[m₁]).

2.8. Group Analysis

To assess the association between the distribution of MWF for our cohort of MS patients and cortical thickness, we used a multivariate linear regression model. In order to compare our statistical reduction technique with a classical method based on lesion masks, we used two sets of possible predictors. For our gamma-mixture approach we considered as predictor variables the MWF distribution mixing factor λ, and distribution modes m₁ and m₂. For the standard approach based on lesion masks, we used mean MWF in the lesion mask validated by the trained neurologist, and mean MWF in the mutually exclusive NAWM mask. We also considered for both sets of predictors age, gender, DMT and DD. The best model for each approach and each group was chosen according to the procedure detailed in the Appendix.

One aim was to compare the best model obtained with predictors derived from our statistical reduction method and the best model obtained with predictors derived from the classical approach based on lesion masks. This was carried out in each group by comparing the adjusted R² of these two models. Another aim was to examine if the final models were different in the RRMS and SPMS group, and interpret the models in terms of their predictors. To help in this endeavor, we assessed the contribution of each coefficient to R² by calculating the average contribution of that coefficient when added to all possible models excluding that coefficient, as in Lindeman et al. (1980).

All p-values were two-sided with statistical significance evaluated at the 0.05 alpha level. Statistical analysis was conducted using R: A Language and Environment for Statistical Computing, R Development Core Team, Vienna, Austria, 2011, http://www.R-project.org.

3. Results

3.1. Multiple Metrics Measured in a Large Cohort of RRMS and SPMS Patients

In all patients, the T2 spiral sequence was easily placed within the clinical protocol given the feasible acquisition time of only 10 min. The combination of T2 spiral acquisition and multi-Gaussian post-processing algorithm produced whole brain, high-resolution MWF white matter maps with good lesion detection (Figure 1).

Cortical thickness measures and other patient characteristics in each of the RRMS and SPMS groups are summarized in Table 1. Univariate t-tests (Welch Two Sample t-tests) demonstrated significant differences in age, cortical thickness and mean whole-brain WM MWF between the RRMS and SPMS groups.

3.2. Modeling Whole Brain White Matter MWF Distribution

Lesion centric data analysis based upon semi-automated lesion mapping programs introduces a bias based upon the reviewer's determination of perceived regions of pathology. We know that regions considered to be dirty WM on T2 imaging have less myelin (Moore et al., 2008; Seewann et al., 2009), therefore we considered an unbiased approach to model each individual's whole brain white matter MWF data as either a gamma random variable or a mixture of two gamma-distributed random variables. The likelihood ratio test provided p-values less than 0.00001 in all patients, thus there was statistical evidence to conclude that a gamma mixture distribution was preferred over the single gamma distribution. As a result the gamma mixture model was chosen for all patients. Further, a univariate t-test demonstrated a significant difference in mean whole-brain WM MWF between the RRMS and SPMS groups (Table 1), suggesting a possible differential distribution of MWF in each of these groups. Figure 3 includes examples of 6 individual patients; the figure shows a histogram of the MWF raw data and the smooth green curve represents the fitted mixture distribution. Different distribution shapes according to MS groups are visible on that figure, with a trend to higher m₁ and m₂, and lower λ, in RRMS patients compared to SPMS patients. Note that the method used to carry out the gamma mixture fitting, expectation maximization, did not rely on the binning of the data or its representation as an histogram. Histograms have been created exclusively for visualization of the fitting procedure outcome.

3.3. Relationship between Gamma Mixture Components and Lesion/NAWM Masks

To interpret the gamma-mixture model parameters in terms of the underlying pathology, we tested that the lesion mask had a higher Dice coefficient with ROI[m₁], defined from the first gamma component, than with ROI[m₂], defined from the second component, in both RRMS and SPMS groups. The t-test demonstrated it was significantly higher in each group (p < 0.00001, Figure 4, top).

FIGURE 4

Figure 4. Comparison of dice coefficient between ROIs from gamma mixture components and lesion/NAWM masks in (left) RRMS and (right) SPMS groups. (Top) Comparison between lesion mask and ROI[m₁], and lesion mask and ROI[m₂] in (left) RRMS and (right) SPMS groups. The difference was shown to be significant in both RRMS (p < 0.00001) and SPMS (p < 0.00001) groups. (Bottom) Comparison between NAWM mask and ROI[m₁], and NAWM mask and ROI[m₂] in (left) RRMS and (right) SPMS groups. The difference was shown to be significant in both RRMS (p < 0.00001) and SPMS (p < 0.00001) groups. *p < 0.01, **p < 0.001, ***p < 0.0001, ****p < 0.00001.

Similarly we tested that the NAWM had a higher Dice coefficient with ROI[m₂] than with ROI[m₁] in both groups, and it was also significantly higher in each group (p < 0.00001, Figure 4, bottom).

This provided support for interpretating m₁ as associated to lesioned WM and m₂ to NAWM. An illustration of lesion and NAWM masks, as well as ROI[m₁] and ROI[m₂], is provided in Figure 5. It can be seen that ROI[m₁] tend to overlap with the lesion mask while also including dirty WM, while ROI[m₂] is mostly included within NAWM and tend to exclude lesion and dirty WM.

FIGURE 5

Figure 5. Illustration of lesion and NAWM masks, and their relationships with the ROIs defined from the gamma mixture, ROI[m₁] and ROI[m₂], in one (A) RRMS and one (B) SPMS patients. (Top) FLAIR images with hyperintensities corresponding to lesions. (Middle) Lesion mask (red) and ROI[m₁] (blue). (Bottom) NAWM mask (green) and ROI[m₂] (yellow). ROI[m₁] and ROI[m₂] have been created from the gamma mixture so that not to overlap one with the other (cf. Figure 2), and as a result they only represent part (approximately half) of each gamma component. Despite this limitation, it can be clearly seen that the lesion mask is more similar to ROI[m₁] than to ROI[m₂], and that the NAWM mask is more similar to ROI[m₂] than to ROI[m₁], as corroborated by the DICE coefficient calculations shown in Figure 4.

3.4. Predictors of Cortical Thickness for SPMS and RRMS Patients

We reasoned that the degree of cortical atrophy may directly correlate with whole brain myelin content. More precisely, we wanted to test the hypothesis that MWF distribution was associated with cortical thickness, and that a differential association existed between MS groups. We also investigated how a fully automated model performed to accomplish this objective compared to a classical approach based on manually edited lesion masks.

For the gamma-mixture approach, we utilized three metrics from the MWF distribution (λ, m₁, m₂) to explore the relationship of myelin and cortical thickness among RRMS and SPMS patients. For the standard approach we used the mean MWF in the lesion mask and the mean MWF in the NAWM mask. In both cases we controlled for other demographics and patient predictors, (age, gender, DD, DMT) and applied an exhaustive model search coupled to 10-fold cross-validation to select the best model for each approach and for each RRMS and SPMS group (cf. Appendix). For the classical approach, the final parameters kept in the model were {age, mean lesion MWF} with a resulting adjusted R² of 0.13 for the RRMS group, and were solely {mean NAWM MWF} with a resulting adjusted R² of 0.39 for the SPMS group. For the statistical reduction approach, the final parameters were {λ, m₁, age} with an adjusted R² of 0.16 for the RRMS group, and were {m₂} with an adjusted R² of 0.45 for the SPMS group. The greater adjusted R² for the gamma-mixture approach demonstrated that the automated method outperformed the lesion mask approach to explain the variance in cortical thickness in either patient group. Further, comparing the final models in each approach, {mean NAWM MWF} and {m₂} for the SPMS group, and {age, mean lesion MWF} and {λ, m₁, age} for the RRMS group, brought additional evidence of the surrogate relationship between NAWM and the second gamma-mixture component, and lesion and the first gamma-mixture component.

Analyzing in more details the statistical reduction model, we found that for the RRMS group, the regression coefficients of the predictor variables were −0.0031 (p < 0.005) for age (in years), 1.56 (p < 0.005) for m₁ (in MWF, a ratio) and −0.30 (p < 0.0005) for λ (a ratio). One can interpret these coefficients as showing that on average a cortical thickness decrease of 0.05 mm is associated with an age increase of 16 years (all other variables held constant), an m₁ decrease of 0.032 (all other variables held constant) or a λ increase of 0.17 (all other variables held constant). As a reminder, m₁ and m₂ are the two modes of the gamma mixture distribution, corresponding to the low and high MWF components of the distribution, determined to relate to lesion and NAWM, respectively, while λ is the mixing ratio determining the proportion of each component in the mixture. A high λ represents a higher contribution of low WM MWF to the gamma mixture while a low λ represents a higher contribution of high WM MWF (Figure 2).

The contribution of each coefficient to the absolute R² is shown in Figure 6. It can be seen that in decreasing order (in terms of contribution to R²) the coefficients are λ, age and m₁, all having a relative importance greater than 20%.

FIGURE 6

Figure 6. Relative contribution to total R² of the cortical thickness regression coefficients for the RRMS group.

For the SPMS group, the regression coefficient was 4.72 (p < 0.0005) for m₂, i.e., on average, a cortical thickness decrease of 0.05 mm is associated with an m₂ decrease of 0.011. The association between cortical thickness and the position of the gamma mixture second mode m₂ is illustrated in Figure 7.

FIGURE 7

Figure 7. Regression of cortical thickness from the position of the gamma mixture second mode m₂ in the SPMS group.

These results suggest a marked difference in the association between cortical thickness and WM MWF distribution in RRMS and SPMS patients.

4. Discussion

Progressive disability in MS is related to cortical thinning (Calabrese et al., 2007; Ceccarelli et al., 2008; Ramasamy et al., 2009; Roosendaal et al., 2011; Batista et al., 2012) and thus important insights into MS could be obtained from understanding the mechanisms driving this cortical loss. The cause of cortical GM atrophy, and whole brain GM damage in general, in multiple sclerosis is currently unclear. Although correlation between WM insult (in terms of lesion load) and GM atrophy has been shown in a number of studies, other investigations support other independent processes, such as inflammation. We hypothesized that these mechanisms differ according to the stage of the disease, notably WM pathology which is the focus of the present study. The statistical analysis undertaken allowed to explore the differential WM correlates of cortical thinning in RRMS and SPMS patients while avoiding the bias introduced by lesion masks.

4.1. Gamma-Mixture Fitting Applied to MWF

Multi-component T2 relaxometry (Koenig et al., 1990; MacKay et al., 1994, 2006; Whittall et al., 1997; Laule et al., 2004, 2007) is a promising noninvasive method to measure MWF. MWF has been shown to correlate highly with histological myelin measurement in rodents (Gareau et al., 2000; Kozlowski et al., 2008; McCreary et al., 2009) and myelin content in post-mortem MS samples, (Moore et al., 2000, 2008; Laule et al., 2006, 2008, 2011) but there have been limited studies translating this technique to study the association with MS disease variables. Through the clinical feasibility of our 10-min whole brain T2prep 3D spiral sequence, we were able to apply this technology to a large cohort of multiple sclerosis patients. Note that although we chose MWF for this study, any surrogate marker of myelin could be used with our novel statistical reduction technique, such as magnetization transfer ratio (Horsfield, 2005) or diffusion MRI metrics (Rovaris et al., 2005).

We illustrated the use of our whole brain white matter gamma-mixture fitting on some of the largest MWF quantitative analyses conducted on multiple sclerosis patients. Another hurdle when investigating WM damage is localizing abnormal WM. Semi-automated delineation of WM lesions has been the gold standard for this purpose although it requires expert knowledge, is particularly time consuming, and biased in nature. We proposed an alternative solution by examining all WM voxels without the use of user-defined lesion masks. Whole brain approaches have been carried out previously with histogram analyses, both from MTR and diffusion MRI data (see for example Van Buchem et al., 1996; Cercignani et al., 2001, respectively). In the current study we chose a whole WM approach which, contrarily to the aforementioned studies, preserves an interpretation of the model parameters in terms of damaged and NAWM. The visible skewness of the data and the possible presence of two peaks led to choose between a single gamma and a gamma mixture, with the latter demonstrated to provide the best fit in all subjects. The position of the peaks of the first and second components were shown to relate to the amount of pathological and NAWM, respectively, as shown by the significant difference in DICE coefficients. Beyond lesional WM, the first component is assumed to also include dirty WM with MWF values lower than the high MWF component. This results in larger ROI[m₁] volume compared to lesion masks, and would explain why the absolute value of Dice(lesion mask, ROI[m₁]) is not as high as between Dice(NAWM, ROI[m₂]): the relatively limited size |lesion mask| of the lesion volume constrains the Dice numerator |lesion mask∩ROI[m₁]|, and thus the Dice coefficient itself, not to be large. Explicitely and automatically distinguishing between NAWM, with high MWF, and lesional and dirty WM, with low MWF, is a strength of the gamma mixture model.

The process of biomarker development requires qualification, i.e., linking the biomarker to the biological process and clinical endpoints. This work aimed at qualifying MWF by demonstrating its relevance to explain some of the variance in patient cortical thickness and to differentiate between disease subtypes, namely RRMS and SPMS. Through our statistical approach, three variables describing the distribution of each patient's whole brain MWF (m₁, m₂, λ) were used as independent variables for the multi-variate analysis. Importantly, the distributions were created for each individual patient without a prior classification of the data; specifically we did not identify lesion or NAWM. This approach removes a potential bias introduced into lesion identification through manual delineation or correction if lesion-detection is partly automated. We compared it with models based on mean lesional and NAWM calculated from lesion masks, and found the adjusted R² of the gamma-mixture models to be higher. This demonstrated the superiority of our statistical reduction technique despite bypassing the time-consumming step of expert drawing of lesion masks. Importantly, the modeled mixture parameters could be interpreted in terms of the underlying WM aspects. Variation in its lowest mode m₁ corresponds to changes in abnormal WM, variation in m₂ corresponds to changes in NAWM and λ indicates the proportion of abnormal to NAWM.

4.2. Differential Relationship of Cortical Thinning between RRMS and SPMS Groups as Assessed from the Gamma Mixture Model

The multivariate regression models for RRMS and SPMS demonstrated that the MWF model parameters (m₁, m₂, λ) were significant in explaining the variance in cortical thickness, but formed complementary subsets of predictors depending on the group at stake: λ and m₁ (together with age) were significant predictors for RRMS patients while only m₂ was a significant predictor for SPMS patients. This suggests that different processes are ongoing, or similar processes but with different magnitudes, during these two stages of the disease. The variables m₁ and λ are associated with low MWF, assumed to represent pathological WM, and the fact we did not find them to be significant predictors of cortical GM atrophy in the SPMS group contrarily to the RRMS group is an interesting finding. This is in line with the recent work of Steenwijk et al. (2014) and Steenwijk et al. (2015) who similarly reported a differential relationship of WM damage and cortical loss among RRMS and SPMS patients. In this study traditional T2 hyperintense lesion load characterized WM pathology. It was shown to be a significant predictor of cortical thickness in the RRMS group but not in the SPMS group, for which no WM correlate of cortical thinning was found (Steenwijk et al., 2014). The advantage of our study was the utilization of an MRI metric more specific to myelin damage, namely MWF, with which we could more accurately characterize this relationship. It notably led us to find a surrogate of NAWM, m₂, highly correlated with cortical thinning in the SPMS group, and provided strong support for the existence of a differential relationship of cortical thinning between RRMS and SPMS groups.

The moderate adjusted R² (R² = 0.16) obtained for the RRMS group to explain cortical thickness suggests that an increase of abnormal WM in terms of relative amount and pathology (as indicated by λ and m₁, respectively) is associated to cortical thinning, although in a limited fashion. More precisely, a change in m₁ is assumed to result from focal demyelination, likely predominantly found within lesions and dirty WM. These results suggest that a small contribution of secondary degeneration is related to focal WM myelin injury and other sources not included in our model such as normal aging and primary disease related GM pathology (i.e., GM lesions Seewann et al., 2011), may drive the decrease in cortical thickness. The comparatively high R² calculated in the SPMS group with only a predictor corresponding to NAWM suggests that cortical thinning has little dependence on focal WM damage (of which m₁ is a surrogate marker). In contrast to RRMS, WM with higher MWF, associated to m₂ and often reprensenting the larger volume of WM, had a profound effect. It may be related to more diffuse inflammation at this stage of the disease and consequently to more global white matter damage (Kutzelnigg et al., 2005). Interestingly, as in Steenwijk et al. (2014), we did not find any influence of age on cortical thinning in SPMS patients, which suggests that age-related neuronal loss is likely outweighed by a shift of NAWM toward lower MWF, possibly illustrating overwhelming global white matter damage.

4.3. Limitations

The causality relationship between WM demyelination and cortical thinning is restricted to an assumption in this work, and a longitudinal study would be required to validate it. Additionally we investigated global GM atrophy, and did not refine the analysis in terms of individual GM regions. Both of these aspects will be tackled in future research using the profilometry framework we introduced recently in Dayan et al. (2016). It would be particularly useful for these two aims, as it can identify WM tract connecting a given pair of GM regions, and assess how MWF changes with time along these tracts. This approach would be advantageous not only in diseases involving whole brain WM, such as in MS, but also in disorders known to affect markedly specific WM tracts, such as the cerebellar peduncles in ataxia (Dayan et al., 2017) and the corticospinal tracts in amyotrophic lateral sclerosis (Turner et al., 2009). It will be combined with Kuceyeski et al. (2015) to explore the association between GM/WM damage and clinical scores, which was not done in this work. Another limitation is the lack of a cohort of normal controls to further interpret the model parameters. An effort dedicated to recruiting healthy subjects has been started to remedy this impediment. Finally, although we circumvented the creation of lesion masks and its associated bias, part of the processing pipeline still required some minimal editing of the border between GM and WM masks which could be affected by lesions.

5. Conclusion

The relationship between cortical GM loss and WM damage in MS is important to elucidate. We presented in this work a whole brain WM model to investigate this relationship at the RRMS and SPMS stage of the disease avoiding the use of lesion masks and their associated bias. Our results represent one of the first quantitative studies characterizing the distribution of whole brain WM MWF within a large cohort of MS patients and relating these measures to cortical atrophy.

Our findings supports the idea that secondary degeneration could play a role in cortical thinning during the first phase of MS as a result of focal demyelination, in parallel with other processes such as GM lesions. We also showed that a different mechanism affected SPMS patients, where high MWF WM presented a strong correlation with cortical GM loss, possibly driving this atrophy. A longitudinal analysis would be required to ascertain that cortical thinning appears after WM damage.

The data reduction technique presented was shown to unravel a variety of different relationships between WM and GM while collapsing the high dimensionality of MWF maps to a few parameters. As it does not rely on lesion masks, it avoids the associated experimenter bias and does not suffer from the potential time-consumming operations required to create such masks. These are important advantages which should make this technique suitable for other MS studies.

Our study suggests that MWF, as an indirect biomarker of myelin, is a promising tool to study the influence of myelin change on cortical atrophy. The rate of cortical loss has been suggested to relate to the rate of disease progression and our results provide further evidence that inducing remyelination would also function as a neuroprotective therapeutic strategy (Piaton et al., 2010). Coupled to our gamma distribution model, this imaging biomarker could help adapt therapeutics according to the different disease stages. In particular, intervention early on would ideally prevent conversion to the progressive stage of the disease. To this aim, future work will include examining the variability of these results across MWF calculation algorithms and MRI scanners.

Author Contributions

MD participated in the design, statistical modeling, data analysis, interpretation and writing of the work. SH participated in the design, statistical modeling, data analysis, interpretation and editing of the manuscript. KF, SP, EL, and EM participated in the data analysis and revising of the manuscript. TV participated in design and revising of the manuscript. TN and AR participated in the design, data analysis and revising of the manuscript. SG contributed as senior author and participated in the design, data analysis and editing of the manuscript.

Conflict of Interest Statement

TV is a speaker for Teva Neuroscience and Mallinckrodt, has received honoraria for advising Genzyme, and has received grant support from the National MS society, Biogen, and Mallinckrodt. SG has received honoraria for advising Genentech and Teva Neuroscience and has received grant support from the National MS society, Biogen, Novartis Pharmaceuticals, Mallinckrodt, Genzyme, and EMD Serono.

The other authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work was supported by the NMSS society grant RG 466-A-2 and the CTSC grant UL1 TR000456-06. The authors would like to acknowledge the help of Dr. Nancy Nealon and Dr. Jai Permumal at the Department of Neurology of Weill Cornell Medicine, New York, NY. The authors would like to thank the editor and two reviewers for their valuable comments and suggestions which helped improve and clarify this manuscript.

References

Batista, S., Zivadinov, R., Hoogs, M., Bergsland, N., Heininen-Brown, M., Dwyer, M., et al. (2012). Basal ganglia, thalamus and neocortical atrophy predicting slowed cognitive processing in multiple sclerosis. J. Neurol. 259, 139–146. doi: 10.1007/s00415-011-6147-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Calabrese, M., Atzori, M., Bernardi, V., Morra, A., Romualdi, C., Rinaldi, L., et al. (2007). Cortical atrophy is relevant in multiple sclerosis at clinical onset. J. Neurol. 254, 1212–1220. doi: 10.1007/s00415-006-0503-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Calabrese, M., Rinaldi, F., Grossi, P., and Gallo, P. (2011). Cortical pathology and cognitive impairment in multiple sclerosis. Exp. Rev. Neurotherapeut. 11, 425–432. doi: 10.1586/ern.10.155

PubMed Abstract | CrossRef Full Text | Google Scholar

Calabrese, M., Rinaldi, F., Mattisi, I., Grossi, P., Favaretto, A., Atzori, M., et al. (2010). Widespread cortical thinning characterizes patients with ms with mild cognitive impairment. Neurology 74, 321–328. doi: 10.1212/WNL.0b013e3181cbcd03

PubMed Abstract | CrossRef Full Text | Google Scholar

Ceccarelli, A., Rocca, M., Pagani, E., Colombo, B., Martinelli, V., Comi, G., et al. (2008). A voxel-based morphometry study of grey matter loss in ms patients with different clinical phenotypes. Neuroimage 42, 315–322. doi: 10.1016/j.neuroimage.2008.04.173

PubMed Abstract | CrossRef Full Text | Google Scholar

Cercignani, M., Inglese, M., Pagani, E., Comi, G., and Filippi, M. (2001). Mean diffusivity and fractional anisotropy histograms of patients with multiple sclerosis. Am. J. Neuroradiol. 22, 952–958. Available online at: http://www.ajnr.org/content/22/5/952.full

PubMed Abstract | Google Scholar

Chatterjee, A., and Chakrabarti, B. K. (2007). Kinetic exchange models for income and wealth distributions. Eur. Phys. J. B 60, 135–149. doi: 10.1140/epjb/e2007-00343-8

CrossRef Full Text | Google Scholar

Chotikapanich, D., and Griffiths, W. E. (2008). Estimating Income Distributions Using a Mixture of Gamma Densities. New York, NY: Springer New York.

Google Scholar

Dayan, M., Kreutzer, S., and Clark, C. A. (2015). Tractography of the optic radiation: a repeatability and reproducibility study. NMR Biomed. 28, 423–431. doi: 10.1002/nbm.3266

PubMed Abstract | CrossRef Full Text | Google Scholar

Dayan, M., Monohan, E., Pandya, S., Kuceyeski, A., Nguyen, T. D., Raj, A., et al. (2016). Profilometry: a new statistical framework for the characterization of white matter pathways, with application to multiple sclerosis. Hum. Brain Mapp. 37, 989–1004. doi: 10.1002/hbm.23082

PubMed Abstract | CrossRef Full Text | Google Scholar

Dayan, M., Olivito, G., Molinari, M., Cercignani, M., Bozzali, M., and Leggio, M. (2017). Impact of cerebellar atrophy on cortical gray matter and cerebellar peduncles as assessed by voxel-based morphometry and high angular resolution diffusion imaging. Funct. Neurol. 31, 239–248. doi: 10.11138/FNeur/2016.31.4.239

PubMed Abstract | CrossRef Full Text | Google Scholar

De Graaf, R., and Nicolay, K. (1998). Adiabatic water suppression using frequency selective excitation. Magn. Reson. Med. 40, 690–696. doi: 10.1002/mrm.1910400508

PubMed Abstract | CrossRef Full Text | Google Scholar

Dempster, A. P., Laird, N. M., and Rubin, D. B. (1977). Maximum likelihood from incomplete data via the em algorithm. J. R. Stat. Soc. 39, 1–38.

Google Scholar

Egger, C., Opfer, R., Wang, C., Kepp, T., Sormani, M. P., Spies, L., et al. (2017). MRI FLAIR lesion segmentation in multiple sclerosis: does automated segmentation hold up with manual annotation? Neuroimage 13, 264–270. doi: 10.1016/j.nicl.2016.11.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Evangelou, N., Esiri, M. M., Smith, S., Palace, J., and Matthews, P. M. (2000). Quantitative pathological evidence for axonal loss in normal appearing white matter in multiple sclerosis. Ann. Neurol. 47, 391–395. doi: 10.1002/1531-8249(200003)47:3<391::AID-ANA20>3.0.CO;2-J

PubMed Abstract | CrossRef Full Text | Google Scholar

Fisher, E., Lee, J. C., Nakamura, K., and Rudick, R. A. (2008). Gray matter atrophy in multiple sclerosis: a longitudinal study. Ann. Neurol. 64, 255–265. doi: 10.1002/ana.21436

PubMed Abstract | CrossRef Full Text | Google Scholar

Gareau, P. J., Rutt, B. K., Karlik, S. J., and Mitchell, J. R. (2000). Magnetization transfer and multicomponent t2 relaxation measurements with histopathologic correlation in an experimental model of ms. J. Magn. Reson. Imaging 11, 586–595. doi: 10.1002/1522-2586(200006)11:6<586::AID-JMRI3>3.0.CO;2-V

PubMed Abstract | CrossRef Full Text | Google Scholar

Geurts, J., and Barkhof, F. (2008). Grey matter pathology in multiple sclerosis. Lancet Neurol. 7, 841–851. doi: 10.1016/S1474-4422(08)70191-1

CrossRef Full Text | Google Scholar

Greve, D. N., and Fischl, B. (2009). Accurate and robust brain image alignment using boundary-based registration. Neuroimage 48, 63–72. doi: 10.1016/j.neuroimage.2009.06.060

PubMed Abstract | CrossRef Full Text | Google Scholar

Horsfield, M. (2005). Magnetization transfer imaging in multiple sclerosis. J. Neuroimaging 15(4 Suppl.), 58S–67S. doi: 10.1177/1051228405282242

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenista, E., Rehwald, W., Chen, E.-L., Kim, H., Klem, I., Parker, M., et al. (2013). Motion and flow insensitive adiabatic t2-preparation module for cardiac mr imaging at 3 tesla. Magn. Reson. Med. 70, 1360–1368. doi: 10.1002/mrm.24564

PubMed Abstract | CrossRef Full Text | Google Scholar

Kappos, L., Freedman, M., Polman, C., Edan, G., Hartung, H.-P., Miller, D., et al. (2007). Effect of early versus delayed interferon beta-1b treatment on disability after a first clinical event suggestive of multiple sclerosis: a 3-year follow-up analysis of the benefit study. Lancet 370, 389–397. doi: 10.1016/S0140-6736(07)61194-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Koenig, S. H., Brown, R. D. R., Spiller, M., and Lundbom, N. (1990). Relaxometry of brain: why white matter appears bright in MRI. Magn. Reson. Med. 14, 482–495. doi: 10.1002/mrm.1910140306

PubMed Abstract | CrossRef Full Text | Google Scholar

Kolind, S., Madler, B., Fischer, S., Li, D. K., and MacKay, A. (2009). Myelin water imaging: implementation and development at 3.0t and comparison to 1.5t measurements. Magn. Reson. Med. 62, 106–115. doi: 10.1002/mrm.21966

PubMed Abstract | CrossRef Full Text | Google Scholar

Kozlowski, P., Liu, J., Yung, A. C., and Tetzlaff, W. (2008). High-resolution myelin water measurements in rat spinal cord. Magn. Reson. Med. 59, 796–802. doi: 10.1002/mrm.21527

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuceyeski, A., Vargas, W., Dayan, M., Monohan, E., Blackwell, C., Raj, A., et al. (2015). Modeling the relationship among gray matter atrophy, abnormalities in connecting white matter, and cognitive performance in early multiple sclerosis. Am. J. Neuroradiol. 36, 702–709. doi: 10.3174/ajnr.A4165

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumar, D., Nguyen, T., Gauthier, S., and Raj, A. (2012). Bayesian algorithm using spatial priors for multiexponential t2 relaxometry from multiecho spin echo MRI. Magn. Reson. Med. 68, 1536–1543. doi: 10.1002/mrm.24170

CrossRef Full Text | Google Scholar

Kutzelnigg, A., Lucchinetti, C. F., Stadelmann, C., Bruck, W., Rauschka, H., Bergmann, M., et al. (2005). Cortical demyelination and diffuse white matter injury in multiple sclerosis. Brain 128(Pt 11), 2705–2712. doi: 10.1093/brain/awh641

PubMed Abstract | CrossRef Full Text | Google Scholar

Laule, C., Kozlowski, P., Leung, E., Li, D. K., Mackay, A. L., and Moore, G. R. (2008). Myelin water imaging of multiple sclerosis at 7 T: correlations with histopathology. Neuroimage 40, 1575–1580. doi: 10.1016/j.neuroimage.2007.12.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Laule, C., Leung, E., Li, D. K. B., Troboulsee, A. L., Paty, D. W., MacKay, A. L., et al. (2006). Myelin water imaging in multiple sclerosis: quantitative correlations with histopathology. Multiple Sclerosis 12, 747–753. doi: 10.1177/1352458506070928

PubMed Abstract | CrossRef Full Text | Google Scholar

Laule, C., Vavasour, I. M., Kolind, S. H., Traboulsee, A. L., Moore, G. R., Li, D. K., et al. (2007). Long t2 water in multiple sclerosis: what else can we learn from multi-echo t2 relaxation? J. Neurol. 254, 1579–1587. doi: 10.1007/s00415-007-0595-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Laule, C., Vavasour, I. M., Leung, E., Li, D. K., Kozlowski, P., Traboulsee, A. L., et al. (2011). Pathological basis of diffusely abnormal white matter: insights from magnetic resonance imaging and histology. Multiple Sclerosis 17, 144–150. doi: 10.1177/1352458510384008

PubMed Abstract | CrossRef Full Text | Google Scholar

Laule, C., Vavasour, I. M., Moore, G. R., Oger, J., Li, D. K., Paty, D. W., et al. (2004). Water content and myelin water fraction in multiple sclerosis. a t2 relaxation study. J. Neurol. 251, 284–293. doi: 10.1007/s00415-004-0306-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Lindeman, R., Merenda, P., and Gold, R. (1980). Introduction to Bivariate and Multivariate Analysis. Glenview, IL: Scott, Foresman.

Google Scholar

MacKay, A., Laule, C., Vavasour, I., Bjarnason, T., Kolind, S., and Madler, B. (2006). Insights into brain microstructure from the t-2 distribution. Magn. Reson. Imaging 24, 515–525. doi: 10.1016/j.mri.2005.12.037

CrossRef Full Text | Google Scholar

MacKay, A., Whittall, K., Adler, J., Li, D., Paty, D., and Graeb, D. (1994). In vivo visualization of myelin water in brain by magnetic resonance. Magn. Reson. Med. 31, 673–677. doi: 10.1002/mrm.1910310614

PubMed Abstract | CrossRef Full Text | Google Scholar

McCreary, C. R., Bjarnason, T. A., Skihar, V., Mitchell, J. R., Yong, V. W., and Dunn, J. F. (2009). Multiexponential t2 and magnetization transfer mri of demyelination and remyelination in murine spinal cord. Neuroimage 45, 1173–1182. doi: 10.1016/j.neuroimage.2008.12.071

PubMed Abstract | CrossRef Full Text | Google Scholar

Mikol, D., Barkhof, F., Chang, P., Coyle, P., Jeffery, D., Schwid, S., et al. (2008). Comparison of subcutaneous interferon beta-1a with glatiramer acetate in patients with relapsing multiple sclerosis (the rebif vs glatiramer acetate in relapsing ms disease [regard] study): a multicentre, randomised, parallel, open-label trial. Lancet Neurol. 7, 903–914. doi: 10.1016/S1474-4422(08)70200-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Moore, G. R., Laule, C., Mackay, A., Leung, E., Li, D. K., Zhao, G., et al. (2008). Dirty-appearing white matter in multiple sclerosis: preliminary observations of myelin phospholipid and axonal loss. J. Neurol. 255, 1802–1811. discussion: 1812. doi: 10.1007/s00415-008-0002-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Moore, G. R., Leung, E., MacKay, A. L., Vavasour, I. M., Whittall, K. P., Cover, K. S., et al. (2000). A pathology-MRI study of the short-t2 component in formalin-fixed multiple sclerosis brain. Neurology 55, 1506–1510. doi: 10.1212/WNL.55.10.1506

PubMed Abstract | CrossRef Full Text | Google Scholar

Narayana, P. A., Govindarajan, K. A., Goel, P., Datta, S., Lincoln, J. A., Cofield, S. S., et al. (2013). Regional cortical thickness in relapsing remitting multiple sclerosis: a multi-center study. Neuroimage 2, 120–131. doi: 10.1016/j.nicl.2012.11.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Nguyen, T. D., Wisnieff, C., Cooper, M. A., Kumar, D., Raj, A., Spincemaille, P., et al. (2012). T2 prep three-dimensional spiral imaging with efficient whole brain coverage for myelin water quantification at 1.5 tesla. Magn. Reson. Med. 67, 614–621. doi: 10.1002/mrm.24128

PubMed Abstract | CrossRef Full Text | Google Scholar

Peterson, J., Bö, L., Mörk, S., Chang, A., and Trapp, B. D. (2001). Transected neurites, apoptotic neurons, and reduced inflammation in cortical multiple sclerosis lesions. Ann. Neurol. 50, 389–400. doi: 10.1002/ana.1123

PubMed Abstract | CrossRef Full Text | Google Scholar

Piaton, G., Gould, R. M., and Lubetzki, C. (2010). Axon-oligodendrocyte interactions during developmental myelination, demyelination and repair. J. Neurochem. 114, 1243–1260. doi: 10.1111/j.1471-4159.2010.06831.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Polman, C., O'Connor, P., Havrdova, E., Hutchinson, M., Kappos, L., Miller, D., et al. (2006). A randomized, placebo-controlled trial of natalizumab for relapsing multiple sclerosis. N. Engl. J. Med. 354, 899–910. doi: 10.1056/NEJMoa044397

PubMed Abstract | CrossRef Full Text | Google Scholar

Polman, C., Reingold, S., Banwell, B., Clanet, M., Cohen, J., Filippi, M., et al. (2011). Diagnostic criteria for multiple sclerosis: 2010 revisions to the mcdonald criteria. Ann. Neurol. 69, 292–302. doi: 10.1002/ana.22366

PubMed Abstract | CrossRef Full Text | Google Scholar

Polman, C. H., Reingold, S. C., Edan, G., Filippi, M., Hartung, H. P., Kappos, L., et al. (2005). Diagnostic criteria for multiple sclerosis: 2005 revisions to the “mcdonald criteria.” Ann. Neurol. 58, 840–846. doi: 10.1002/ana.20703

PubMed Abstract | CrossRef Full Text | Google Scholar

Raj, A., Pandya, S., Shen, X., LoCastro, E., Nguyen, T. D., and Gauthier, S. A. (2014). Multi-compartment t2 relaxometry using a spatially constrained multi-gaussian model. PLoS ONE 9:e98391. doi: 10.1371/journal.pone.0098391

PubMed Abstract | CrossRef Full Text | Google Scholar

Ramasamy, D., Benedict, R., Cox, J., Fritz, D., Abdelrahman, N., Hussein, S., et al. (2009). Extent of cerebellum, subcortical and cortical atrophy in patients with ms. A case-control study. J. Neurol. Sci. 282, 47–54. doi: 10.1016/j.jns.2008.12.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Roosendaal, S., Bendfeldt, K., Vrenken, H., Polman, C., Borgwardt, S., Radue, E., et al. (2011). Grey matter volume in a large cohort of ms patients: relation to MRI parameters and disability. Multiple Sclerosis J. 17, 1098–1106. doi: 10.1177/1352458511404916

PubMed Abstract | CrossRef Full Text | Google Scholar

Rovaris, M., Gass, A., Bammer, R., Hickman, S., Ciccarelli, O., Miller, D., et al. (2005). Diffusion MRI in multiple sclerosis. Neurology 65, 1526–1532. doi: 10.1212/01.wnl.0000184471.83948.e0

PubMed Abstract | CrossRef Full Text | Google Scholar

Seewann, A., Vrenken, H., Kooi, E.-J., Van Der Valk, P., Knol, D., Polman, C., et al. (2011). Imaging the tip of the iceberg: visualization of cortical lesions in multiple sclerosis. Multiple Sclerosis J. 17, 1202–1210. doi: 10.1177/1352458511406575

PubMed Abstract | CrossRef Full Text | Google Scholar

Seewann, A., Vrenken, H., van der Valk, P., Blezer, E. L., Knol, D. L., Castelijns, J. A., et al. (2009). Diffusely abnormal white matter in chronic multiple sclerosis: imaging and histopathologic analysis. Arch. Neurol. 66, 601–609. doi: 10.1001/archneurol.2009.57

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, P. L. (2003). Raindrop size distributions: exponential or gammadoes the difference matter? J. Appl. Meteorol. 42, 1031–1034. doi: 10.1175/1520-0450(2003)042<1031:RSDEOG>2.0.CO;2

CrossRef Full Text | Google Scholar

Steenwijk, M., Daams, M., Pouwels, P., Balk, L., Tewarie, P., Killestein, J., et al. (2014). What explains gray matter atrophy in long-standing multiple sclerosis? Radiology 272, 832–842. doi: 10.1148/radiol.14132708

PubMed Abstract | CrossRef Full Text | Google Scholar

Steenwijk, M., Daams, M., Pouwels, P., Balk, L., Tewarie, P., Geurts, J., et al. (2015). Unraveling the relationship between regional gray matter atrophy and pathology in connected white matter tracts in long-standing multiple sclerosis. Hum. Brain Mapp. 36, 1796–1807. doi: 10.1002/hbm.22738

PubMed Abstract | CrossRef Full Text | Google Scholar

Turner, M. R., Kiernan, M. C., Leigh, P. N., and Talbot, K. (2009). Biomarkers in amyotrophic lateral sclerosis. Lancet Neurol. 8, 94–109. doi: 10.1016/S1474-4422(08)70293-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Buchem, M., McGowan, J., Kolson, D., Polansky, M., and Grossman, R. (1996). Quantitative volumetric magnetization transfer analysis in multiple sclerosis: estimation of macroscopic and microscopic disease burden. Magn. Reson. Med. 36, 632–636. doi: 10.1002/mrm.1910360420

PubMed Abstract | CrossRef Full Text | Google Scholar

Vercellino, M., Plano, F., Votta, B., Mutani, R., Giordana, M., and Cavalla, P. (2005). Grey matter pathology in multiple sclerosis. J. Neuropathol. Exp. Neurol. 64, 1101–1107. doi: 10.1097/01.jnen.0000190067.20935.42

PubMed Abstract | CrossRef Full Text | Google Scholar

Whittall, K. P., MacKay, A. L., Graeb, D. A., Nugent, R. A., Li, D. K., and Paty, D. W. (1997). In vivo measurement of t2 distributions and water contents in normal human brain. Magn. Reson. Med. 37, 34–43. doi: 10.1002/mrm.1910370107

PubMed Abstract | CrossRef Full Text | Google Scholar

Wong, T., and Li, W. (2014). Test for homogeneity in gamma mixture models using likelihood ratio. Comput. Stat. Data Anal. 70, 127–137. doi: 10.1016/j.csda.2013.09.001

CrossRef Full Text | Google Scholar

Yoo, C., Jung, K.-S., and Kim, T.-W. (2005). Rainfall frequency analysis using a mixed gamma distribution: evaluation of the global warming effect on daily rainfall. Hydrol. Process. 19, 3851–3861. doi: 10.1002/hyp.5985

CrossRef Full Text | Google Scholar

Zhang, Y., Brady, M., and Smith, S. (2001). Segmentation of brain mr images through a hidden markov random field model and the expectation-maximization algorithm. IEEE Trans. Med. Imaging 20, 45–57. doi: 10.1109/42.906424

PubMed Abstract | CrossRef Full Text | Google Scholar

Appendix

Selection of the Best Multivariate Linear Regression Model

For each approach and each group of MS patients, we had to choose the best model with respect to the N predictors available. This was implemented in two steps.

First, given a group and a number K of predictors, we did an exhaustive search to choose the best K-predictor model. This was done by comparing all possible K-predictor models, and selecting the one with the minimum residual sum of squares (RSS).

Second, among the best 1-predictor, 2-predictor, …, N-predictor models, we used 10-fold cross-validation to select the one with the minimum test-RSS. Briefly, the data is split in 10 subsets (“folds”), and selecting 1-fold in turn, the model is fitted on the remaining 9 folds before being tested on the selected fold. The test-RSS is the average RSS obtained from this iterative procedure on the 10 folds: test-RSS = mean(test-RSS_fold1, test-RSS_fold2, …test-RSS_fold10).