Introduction

Myelin is important in the transmission of neural information. It maintains the integrity of neural fibers and enhances the speed of propagation of action potentials, which are essential for the proper function of the brain1,2. Measuring myelin in the brain by magnetic resonance imaging (MRI) is important for evaluating the development and aging of healthy humans3,4,5. It is also important for estimating the progression of degenerative6 or demyelinating diseases7. Conventional MRI is highly sensitive to tissue contrast, but generally unspecific to tissue properties such as myelin content. Furthermore, lengthy scanning time has hindered the routine clinical use of MRI to obtain myelin measurements. Recently, rapid simultaneous relaxometry based on a single pulse sequence has been developed8. It quantifies the longitudinal relaxation rate (R1), transverse relaxation rate (R2), proton density (PD), and local B1 field in approximately 6 minutes for full head coverage9. The estimated B1 field is used for correction of local variations in flip angle. It is possible to create contrast-weighted images (the technique is called ‘synthetic MRI’) such as T1-weighted (T1w), T2-weighted (T2w), and fluid-attenuated inversion recovery (FLAIR) images, based on the acquired quantitative values, thus obviating the need for acquiring these contrast-weighted images separately10. At the same time, automatic brain segmentation11 and myelin measurement12 are also possible using the acquired quantitative values. These can be done with a dedicated software called ‘SyMRI’ with post-processing time less than 1 minute9. Thus, myelin measurements can now be performed within the limits of clinically allowed scanning time. The myelin model infers myelin volume fraction (MVF) in a voxel based on the effect of myelin on intra- and extracellular water relaxation rates due to magnetization exchange12. The observed R1 and R2 rates of intra- and extracellular water increase in the vicinity of fast relaxing myelin water. On the other hand, the observable PD decreases because myelin water decays much faster than non-myelin water. The SyMRI myelin measurement has been validated on 12 human cadavers using Luxol Fast Blue staining of histological specimens13. A repeatability study has reported a very low error (coefficient of variation, 0.59% for 0.8 mm in-plane resolution) for whole-brain myelin volume calculated using SyMRI14. Myelin volume measured by SyMRI has been shown to depend on age in pediatric populations, especially in children under 4 years old, thus indicating a correlation of this method with the normal myelination process15,16. This method has also been used in studies investigating patients with multiple sclerosis (MS)17,18, Sturge-Weber syndrome19, and cerebral autosomal dominant arteriopathy with subcortical infarcts and leukoencephalopathy (CADASIL)20, with promising results. However, correlation of SyMRI myelin measurement with other MRI techniques sensitive to myelin has not been investigated so far.

There are several other techniques for myelin measurement, including myelin water imaging21,22, macromolecular tissue volume derived from normalized PD mapping23, and magnetization transfer (MT) imaging2. MT is a phenomenon where the proton spins bound to macromolecules, once excited by a radiofrequency pulse, transfer a part of their energy to the neighboring mobile proton spins24. MT imaging estimates the macromolecular proton pool size with ultra-short T2 relaxation by transfer of magnetization to the observable mobile water pool25. MT ratio (MTR) has been widely used based on this theory and shown to correlate well with histological myelin content26,27, but also with other properties such as R124. R1 also correlates strongly with myelin28, meaning that MTR and R1 work against each other and R1 mitigates the power of MTR as a measure of myelin. Further, R1 is also sensitive to iron, calcium content, and axon size29 and count30, thus making the relationship between MTR and actual myelin content nonlinear. MT saturation (MTsat) imaging was developed to improve MTR, by decoupling MTR from R131. MTsat shows higher contrast in the brain than MTR does31, and has been shown to correlate more with disability metrics than MTR in patients with MS32. MTsat has also been shown to correlate well with quantitative MT measures25, which reduces dependency of MT imaging on sequence parameters. However, quantitative MT imaging is time-consuming and the post-processing is still challenging.

T1w/T2w ratio is another approach for assessing myelin content in the cortical gray matter, originally developed to map myeloarchitecturally distinct cortical regions for parcellation of cerebral cortex, thus providing a connectivity measurement33,34. Pixel intensity on T1w and T2w images is assumed to be directly and inversely proportional to myelin contrast, respectively. Thus, ratio of these images is thought to accentuate the intrinsic contrast of myelin. Because intensity scaling of T1w and T2w images differ across scanners and acquisition protocols, Ganzetti et al.35 have suggested that calibration of their intensities prior to making their ratio can increase the reproducibility of T1w/T2w ratio. Although T1w/T2w ratio is not a direct index of myelin, it is still considered a proxy of myelin content36. While intracortical myelin content across different ages has been evaluated using this method36,37, myelination of white matter (WM) in neonatal brains has also been investigated using this method38,39. Further, the test-retest reliability of T1w/T2w ratio has been reported to be high40. Recent histological studies investigated T1w/T2w ratio in patients with MS, showing that T1w/T2w ratio was significantly different between myelinated and demyelinated cortex in MS patients41, and also significantly different in the cortex between early-stage MS and healthy controls42. Because T1w and T2w images are routinely acquired as part of brain MRI protocols, this technique does not increase scanning time. However, the specificity of T1w/T2w to actual myelin content has been doubted by recent studies40,43.

As mentioned above, there are several different methods to estimate myelin volume in the brain. However, investigation of correlation among different methods is scarce. Specifically, no study has examined the correlation of SyMRI as a myelin imaging tool with other methods. Therefore, the aim of this study was to compare SyMRI with two other putative myelin measurement techniques by investigating the correlation of SyMRI with MTsat and T1w/T2w ratio in WM and gray matter (GM).

Results

Scatterplots and Mean Values of MVFMTsat, MVFSyMRI, and MVFT1w/T2w

The calibration factors for MVFMTsat and MVFT1w/T2w were 8.40 and 14.5, respectively, so that their means in the WM equaled that of MVFSyMRI. The scatterplots of these three MVF metrics are shown in Fig. 1. Table 1 shows the mean and standard deviation (SD) of each MVF metric after calibration, and MTsat and T1w/T2w ratio before calibration in each tissue region, with the percentage of MVF in subcortical or cortical GM to that in WM. Because both MVFMTsat and MVFT1w/T2w were calibrated to MVFSyMRI, so that their mean values in the WM were equal, the mean values of WM for all these metrics were the same. The contrasts among WM and subcortical GM, and WM and cortical GM were significantly higher for MVFSyMRI and lower for MVFT1w/T2w than other MVF metrics (p < 0.001).

Figure 1
figure 1

Scatterplots showing correlations among MVFMTsat, MVFSyMRI, and MVFT1w/T2w. For WM, the correlation between MVFSyMRI and MVFMTsat is stronger than the correlation between MVFT1w/T2w and MVFSyMRI or MVFMTsat.

Table 1 MVFMTsat, MTsat, MVFSyMRI, MVFT1w/T2w, and T1w/T2w ratio in WM, subcortical GM, and cortical GM, with the percentage of MVF in subcortical or cortical GM to that in WM.

Correlation Coefficients among MVFMTsat, MVFSyMRI, and MVFT1w/T2w

Table 2 shows the Spearman’s ρ correlation coefficients with their 95% confidence intervals (CIs) among MVF metrics. Correlations were significant for all regions—alone or combined—among these metrics (p < 0.001). In the WM and subcortical GM, the correlation coefficient was the highest between MVFMTsat and MVFSyMRI (p < 0.001 in the WM and p = 0.005 in the subcortical GM). In the WM, MVFT1w/T2w showed only weak to moderate correlation with MVFMTsat or MVFSyMRI. In the cortical GM, the correlation coefficient was the highest between MVFSyMRI and MVFT1w/T2w (p < 0.001), with MVFMTsat vs. MVFT1w/T2w showing the lowest value (p = 0.011). In all regions combined, all these metrics showed strong correlations. Correlation coefficients of MVFMTsat vs. MVFSyMRI and MVFSyMRI vs. MVFT1w/T2w were comparable (p = 0.62) and higher than that of MVFMTsat vs. MVFT1w/T2w (p < 0.001) Table 3 shows the Spearman’s ρ correlation coefficients among MVF metrics in individual areas representative of 10 WM, 2 subcortical GM, and 4 cortical GM, and their mean values. Out of 10 WM ROIs, 8 showed significant correlations between MVFMTsat and MVFSyMRI. The 2 WM ROIs that did not show significant correlation were genu and splenium of corpus callosum, which showed the highest MVFSyMRI. Meanwhile, only 3 and 4 ROIs showed significant correlation between MVFMTsat and MVFT1w/T2w, and MVFSyMRI and MVFT1w/T2w, respectively. Both of the 2 subcortical GM ROIs showed significant correlations in all comparisons, with comparison between MVFMTsat and MVFSyMRI showing the highest and strong correlation coefficients. For all the 4 cortical GM ROIs, comparison among MVFSyMRI and MVFT1w/T2w revealed the highest and significant correlations, whereas only 1 ROI (precentral) showed significant correlation between MVFSyMRI and MVFMTsat, and no significant correlation was observed between MVFMTsat and MVFT1w/T2w.

Table 2 Correlation among MVFMTsat, MVFSyMRI, and MVFT1w/T2w for WM, subcortical GM, cortical GM, and all regions.
Table 3 Correlation among MVFMTsat, MVFSyMRI, and MVFT1w/T2w for 48 WM ROIs, and MVFMTsat, MVFSyMRI, and MVFT1w/T2w in each ROI.

Regression Analysis of MVFSyMRI and MVFT1w/T2w as a Function of MVFMTsat

Table 4 shows the values of the intercept and slope with their standard error in each region—alone or combined—for MVFSyMRI and MVFT1w/T2w as a function of MVFMTsat. In WM, cortical GM, and all regions combined, significant difference was detected between the slopes of MVFSyMRI and MVFT1w/T2w, with that of MVFSyMRI nearer to 1. In subcortical GM, slopes of MVFSyMRI and MVFT1w/T2w did not show statistical significance, and y-intercepts differed significantly with that of MVFT1w/T2w nearer to 0.

Table 4 Intercept and slope of MVFSyMRI and MVFT1w/T2w as a function of MVFMTsat for each region—alone or combined.

Discussion

In this study, we investigated the concurrent validity of SyMRI myelin measurement method by comparing SyMRI with MTsat and T1w/T2w ratio in WM and GM. As part of the study, we tried to estimate the absolute myelin partial volume in a voxel by these three methods. SyMRI directly estimates MVF of a voxel by bloch simulation. On the other hand, MTsat and T1w/T2w ratio require calibration to be used as quantitative measures of myelin content. Thus, we calibrated MTsat and T1w/T2w ratio for their means in the whole WM to be equal to that of MVFSyMRI, partly because calibration method does not affect correlation coefficient and contrast between WM and cortical or subcortical GM. In this study, the mean MVFSyMRI in the WM was 30.70%. This corresponds to the previously reported values (around 25–30%) of MVF in WM, investigated by histology2,44. This value also corresponds to the results of MVF investigated using SyMRI for WM of cadavers (30.98%)13 and normal-appearing WM of MS patients (32.88% and 30.96%)17,18. For GM, reports on investigation into MVF by histology are rather scarce and most were performed with optical density using Luxol Fast Blue stain, which could be used only in comparison with the values of other brain microstructures45. Previous studies that investigated volume fraction of myelin in the brain showed optical densities of subcortical and cortical GM to be around 49–67% and 9.8–36% that of WM, respectively13,46. In our study, MVFSyMRI corresponded to the results of these histological studies in cortical GM better than MVFMTsat and MVFT1w/T2w. For subcortical GM, MVFMTsat and MVFSyMRI were comparable and these showed better correspondence to previous histological study than MVFT1w/T2w. In terms of WM to GM contrast, we conclude that MVFSyMRI was the best fit to the results of previous histological studies among the metrics investigated in our study.

In our study, we investigated the correlation among three different metrics for myelin content. The aim was to show the concurrent validity of MVFSyMRI by MVFMTsat and MVFT1W/T2W. For WM, MVFSyMRI showed strong and higher correlation with MVFMTsat than MVFT1w/T2w. In regression analysis, the slope was closer to 1 for MVFSyMRI than MVFT1w/T2w as a function of MVFMTsat in WM. These results are in line with the study by Arshad et al.40. They investigated the correlation between T1w/T2w ratio and myelin water fraction in WM, and found that T1w/T2w ratio poorly correlated with myelin water fraction and correlated more with geometric mean of multi-echo T2 relaxation, which had been shown to correlate with axon diameter based on histology, rather than myelin content47. Another study also showed poor correlation between T1w/T2w and myelin water fraction43. Therefore, T1w/T2w ratio may not be a suitable candidate as a measure of myelin in WM. In cortical GM, these three MVF metrics showed moderate to strong correlations to each other, with MVFSyMRI and MVFT1w/T2w showing a higher correlation. However, we cannot determine which is the best measure for estimating myelin content in GM among these three metrics at this moment. Myeloarchitecture is different among cortical areas, and high-resolution T1w/T2w ratio has been widely used for cortical parcellation, especially in the Human Connectome Project, showing good results48. In a future study, comparison of these metrics for the ability of cortical parcellation should be investigated. However, recent histological study showed that T1w/T2w ratio in the cerebral cortex correlated well with dendrites, but not with myelin, even though the sample size was small (9 MS patients)42. There is a possibility that T1w/T2w ratio does not reflect actual myelin content in the brain. All regions in aggregate showed strong correlation coefficients in all comparisons (i.e. MVFMTsat vs. MVFSyMRI, MVFMTsat vs. MVFT1w/T2w, and MVFSyMRI vs. MVFT1w/T2w). This may be because subgroups with different microstructures were included in the analysis.

When we analyzed individual structures representative of WM, subcortical GM, and cortical GM, the correlation coefficients showed similar tendency to those shown for each segment as a whole. Of note, only genu and splenium of corpus callosum out of the 10 WM ROIs did not show significant correlation between MVFMTsat and MVFSyMRI, with these showing the highest MVFSyMRI. This may be because SyMRI does not assume nonphysiological MVF higher than 40%12, and disagreement may have occurred between SyMRI and MTsat with high values.

Determination of the precise relationship between MRI measures of myelin and actual MVF is especially important for calculating the g-ratio, which is the ratio of the inner and the outer diameter of a myelinated nerve fiber49. Calculation of the g-ratio by MRI can be performed with myelin imaging in combination with diffusion MRI, such as diffusion tensor imaging (DTI) and neurite orientation dispersion and density imaging (NODDI)49,50. Because diffusion MRI alone is not sufficient to estimate axon volume fraction49, precise measurement of myelin is necessary for correct g-ratio calculation. Furthermore, g-ratio could complement MVF measurements in understanding tissue microstructure, because MVF only cannot differentiate partial demyelination of neuronal fibers from loss of axons, with the remaining axons fully myelinated. Thus, g-ratio can provide a more complete picture of the microstructure, which is important for understanding plasticity of the normal brain51 and may also be important for the care of patients with MS in choosing immunotherapy or remyelination therapy25. Because we could not perform histological measurements of actual myelin content in this study, we calibrated MTsat and T1w/T2w ratio to MVFSyMRI. Even though we assumed zero-intercept upon calibration of MVFMTsat and MVFT1w/T2w to MVFSyMRI, we detected a non-zero intercept when linear regression was performed. This means that at least two of these MVF metrics are not perfectly specific to myelin content in the brain. Although it may be expected that MTsat is also sensitive to macromolecules other than myelin, the specificity of our MVF metrics to actual myelin content should be investigated more precisely in future histological studies. We should also be aware that scaling factors depend on the acquisition protocol and post-processing, and should be carefully determined for each investigation25.

Rapid relaxation of myelin water cannot be directly measured by the SyMRI sequence, but the presence of MVF can be inferred by its effect of magnetization exchange with the slower cellular relaxation, as well as the decrease in observed PD. This is an indirect measurement and may have some limitations when compared with a more direct approach, such as myelin water fraction, which estimates T2 distribution of water including myelin water by fitting multi-exponential T2 decay22 and has been shown to correlate well with histological myelin content in patients with MS52. However, for clinical use, the robustness and easy implementation may be more important. SyMRI myelin measurement has been shown to have good repeatability, which is important for longitudinal studies14. In addition to myelin measurements, any contrast-weighted image can also be generated by SyMRI53, thus obviating the need for further conventional scans.

There are several limitations in this study. First, the resolutions of the images were different between MVFSyMRI or T1w/T2w ratio (2D acquisition) and MTsat (3D acquisition). Even though the difference in resolution could introduce deviation in the quantification, this would have been offset by a large number of ROIs used in this study. However, the analyses of 2D and 3D images by consistent methods was a challenge in our study. Rather than co-registering these images, we registered ROIs in template space to 2D or 3D space for each subject. Co-registration may cause some mis-registration, which will result in inappropriate comparison of voxels derived from different tissues. When we applied the ROIs to each MVF map, we used partial volume maps of GM, WM, or both, with thresholding, to minimize partial volume effects. Second, T1-weighted images for T1w/T2w ratio were acquired by a spin-echo sequence, even though mostly gradient-echo sequences have been used for calculating T1w/T2w ratio33,35,36,40,48. Because T1w/T2w ratio is a semi-quantitative value, different acquisitions may introduce different contrasts. However, T1w/T2w ratio has been shown to give very similar overall results when acquired on different scanners with different sequences and different field strengths33,35. Third, the myelin measurement methods investigated in this study may show variable behaviors in diseased brains from healthy brains, not only due to demyelination but also due to edema, inflammation, iron accumulation, or atrophy. This should be investigated in future studies. For example, MTR seems to correlate with not only myelin but also with change in water content caused by inflammation or edema in patients with MS54. Even though we assumed a linear relationship for calibration of MVF values, this assumption may not hold true in diseased brains.

In summary, we compared MTsat, MVFSyMRI, and T1w/T2w ratio as quantitative measures of myelin in the brain. We calibrated MTsat and T1w/T2w in WM to be equal to MVFSyMRI in WM (MVFMTSat and MVFT1w/T2w). Correlation of these metrics in WM was strong and higher between MVFMTsat and MVFSyMRI than between MVFT1w/T2w and MVFMTsat or MVFSyMRI, indicating that MVFMTsat and MVFSyMRI are similarly suited to measure myelin in the WM, whereas MVFT1w/T2w may be less optimal. In GM, moderate to strong correlation was observed among these metrics. However, further studies performing cortical parcellation using these measures or investigating the correlation between each MVF metric and histology should be conducted before concluding which is the best measure for estimating myelin content in GM.

Materials and Methods

Study Participants

Twenty healthy volunteers (9 male and 11 female, mean age 55.3 years, age range 25–71 years) were included in this study. These subjects were screened by a questionnaire for neurological or psychological symptoms, or history of neurologic diseases. Acquired images were also screened for moderate-to-severe WM ischemic lesions (Fazekas grade 2 or more55), asymptomatic cerebral infarction, or regional brain atrophy.

Ethical issue

All data from the patients were obtained in accordance with the 2013 revised Helsinki Declaration of 1964. We provided participants with detailed information, and written informed consent was obtained from all participants. The Ethical Committee of Juntendo University Hospital approved the study.

MRI Acquisition Protocol for SyMRI

All subjects were scanned on a single 3T MRI scanner (MAGNETOM Prisma, Siemens Healthcare, Erlangen, Germany) using a 64-channel head coil. QRAPMASTER (an acronym derived from ‘quantification of relaxation times and proton density by multi-echo acquisition of a saturation-recovery by using turbo spin-echo readout’ for simultaneous tissue relaxometry) was performed for all subjects. QRAPMASTER is a two-dimensional (2D) axial multi-slice, multi-echo, and multi-saturation delay saturation-recovery turbo spin-echo acquisition method with which images are collected with different combinations of echo times (TEs) and saturation delay times. In our institution, combinations of 2 TEs and 4 delay times were used to make a matrix of 8 complex images that were then used to quantify longitudinal R1 relaxation and transverse R2 relaxation rates and PD by using SyMRI software 8.0 (SyntheticMR, Linköping, Sweden). The TEs were 22 and 99 ms, and the delay times were 170, 620, 1970, and 4220 ms. The repetition time (TR) was 4250 ms. The other parameters used for QRAPMASTER were as follows: field of view (FOV) 230 × 186 mm; matrix 320 × 260; echo-train length 10; bandwidth 150 Hz/pixel; parallel imaging acceleration factor 2; slice thickness/gap 4.0 mm/1.0 mm; 30 sections; and acquisition time 5 min 8 sec.

Processing of SyMRI Data

Based on the R1, R2, and PD values acquired by QRAPMASTER, myelin volume fraction (MVFSyMRI) was also calculated automatically on SyMRI software. This model for myelin measurement hypothesizes 4 compartments in the brain: myelin, cellular, free water, and excess parenchymal water partial volumes12. The model assumes that the relaxation behavior of each compartment contributes to the effective relaxation behavior of an acquisition voxel. The R1, R2, and PD values of free water and excess parenchymal water partial volumes were fixed to those of cerebrospinal fluid (CSF) (R1, 0.24 sec−1; R2, 0.87 sec−1; PD, 100%)8. The R2 of myelin partial volume was fixed to the literature value of 77 sec−156. Optimization of other model parameters were done by performing simulation by running Bloch equations for observable R1, R2, and PD properties in a spatially normalized and averaged brain from a group of healthy subjects12. In this model, the magnetization exchange rates between partial volume compartments are also considered. A lookup grid was made in R1-R2-PD space for all possible distributions (ranging from 0% to 100%) of the four partial volumes. The measured R1, R2, and PD values were projected onto the lookup grid, for estimating the MVFSyMRI in each voxel. Although other methods for myelin imaging require scaling factors to estimate MVF from measured macromolecular pool size or myelin water fraction, assuming linear proportionality2, we omitted this procedure because MVFSyMRI directly estimates the volume fraction of myelin in a voxel12.

Processing of T1w/T2w ratio

Synthetic T1w and T2w images were produced from QRAPMASTER data. Parameters used for T1w images were: TR 500 ms; and TE 10 ms. Parameters used for T2w images were: TR 4500 ms; and TE 100 ms. These T1w and T2w images were intrinsically aligned. Synthetic T1w and T2w images were skull-stripped using the intracranial mask generated by SyMRI software57. In conventional MRI, imperfection of B1 field affects T1w and T2w images, generating intensity non-uniformity in these images. It has been proposed that this non-uniformity should be corrected before the ratio of these images is calculated, because a ratio does not adequately cancel the intensity non-uniformity35. The QRAPMASTER sequence acquires the B1 field map and the acquired quantitative data are automatically corrected for local B1 field when processed by SyMRI software9. Because T1w and T2w images are non-quantitative, the intensity scaling may vary among different individuals, sequences, or scanners. To minimize the effect of intensity scaling, we applied an external linear calibration to these contrast-weighted images as proposed by Ganzetti et al.35, which would provide a more consistent range of T1w and T2w intensities even across different datasets. Two masks of anatomical structures external to the brain—one with high T1w signal intensity and low T2w signal intensity (temporalis muscle) and the other with opposite properties (eye)—were used for calibration. These regions were defined in the MNI152 space using the ICBM152 template images (http://www.bic.mni.mcgill.ca/ServicesAtlases/ICBM152NLin2009) and then warped to each subject’s space using the registration matrix described below in the ROI Analysis section. Distribution peaks (modes) of intensity values were recorded for these regions of interest (ROIs) in each subject. In ICBM152 template images, we recorded the modes as reference values for the eyes as following: 28.2 for T1w images and 99.9 for T2w images. For the temporalis muscle, the values were: 58.6 for T1w images and 21.1 for T2w images. The linear scaling of either T1w or T2w images was performed using the following equation35:

$${I}_{C}=[\frac{{E}_{R}-{M}_{R}}{{E}_{S}-{M}_{S}}]\times I+[\frac{{E}_{S}{M}_{R}-{E}_{R}{M}_{S}}{{E}_{S}-{M}_{S}}]$$
(1)

where I and IC represent the images before and after calibration. ES and MS are the mode intensity values of each subject’s eye and muscle masks, respectively, and ER and MR show the reference values in template images of eye and muscle masks, respectively. After calibrating the T1w and T2w images, their ratio was calculated to produce the T1w/T2w ratio images.

Acquisition and Processing of MTsat

Three three-dimensional (3D) multi-echo fast low-angle shot (FLASH) sequences were performed with predominant T1-, PD-, and MT-weighting for all subjects. For T1w images, TR/excitation flip angle α = 10 ms/13° were used; for PD- and MT-weighted images, 24 ms/4° were used. For MT-weighted images, excitation was preceded by an off-resonance Gaussian-shaped RF pulse (frequency offset from water resonance 1.2 kHz, pulse duration 9.984 ms, and nominal flip angle 500°). For the other parameters, the following was used: slice thickness 1.8 mm; 104 slices; FOV 224 × 224 mm; matrix 128 × 128, parallel imaging using GRAPPA factor 2 in phase-encoding direction; 7/8 partial Fourier acquisition in the partition direction; bandwidth 260 Hz/pixel; and total acquisition time 6 min 25 sec.

These three images were used to calculate the MTsat index31. First, the apparent longitudinal relaxation rate R1app was calculated as follows:

$${R}_{1{\rm{app}}}=\frac{1}{2}\frac{{S}_{{\rm{T}}1}{\alpha }_{{\rm{T}}1}/{{\rm{TR}}}_{{\rm{T}}1}-{S}_{{\rm{PD}}}{\alpha }_{{\rm{PD}}}/{{\rm{TR}}}_{{\rm{PD}}}}{{S}_{{\rm{PD}}}/{\alpha }_{{\rm{PD}}}-{S}_{{\rm{T}}1}/{\alpha }_{{\rm{T}}1}}$$
(2)

where ST1 and SPD denote signal intensities of T1w and PD-weighted images, respectively; TRT1 and TRPD denote TR of T1w and PD-weighted images, respectively; and αT1 and αPD denote excitation flip angles of T1w and PD-weighted images, respectively.

Secondly, the apparent signal amplitude Aapp was calculated as follows:

$${A}_{{\rm{app}}}={S}_{{\rm{PD}}}{S}_{{\rm{T}}1}\frac{{{\rm{TR}}}_{{\rm{PD}}}{\alpha }_{{\rm{T}}1}/{\alpha }_{{\rm{PD}}}-{{\rm{TR}}}_{{\rm{T}}1}{\alpha }_{{\rm{PD}}}/{\alpha }_{{\rm{T}}1}}{{S}_{{\rm{T}}1}{{\rm{TR}}}_{{\rm{PD}}}{\alpha }_{{\rm{T}}1}-{S}_{{\rm{PD}}}{{\rm{TR}}}_{{\rm{T}}1}{\alpha }_{{\rm{PD}}}}$$
(3)

Thirdly, the apparent MT saturation δapp was calculated as follows:

$${\delta }_{{\rm{app}}}=({A}_{{\rm{app}}}{\alpha }_{{\rm{MT}}}/{S}_{{\rm{MT}}}-1){R}_{1{\rm{app}}}{{\rm{TR}}}_{{\rm{MT}}}-{\alpha }_{{\rm{MT}}}^{2}/2$$
(4)

where SMT, TRMT, and αMT denote signal intensity, TR, and excitation flip angle of MT-weighted image, respectively.

The apparent MT saturation is inherently robust against differences in relaxation rates and inhomogeneities of RF transmit and receive field compared with conventional MTR imaging31,58. Furthermore, we also corrected for small residual higher-order dependencies of the MT saturation on the local RF transmit field to further improve spatial uniformity, as suggested by Weiskpof et al.59:

$$M{T}_{{\rm{sat}}}=\frac{{\delta }_{{\rm{app}}}(1-0.4)}{1-0.4R{F}_{{\rm{local}}}}$$
(5)

where RFlocal is the relative local flip angle α compared to the nominal flip angle. RFlocal was calculated by dual-angle method60. For this method, two additional B1 maps using echo-planar imaging with nominal 10° and 20° flip angles were acquired in short acquisition time (around 10 seconds each). The first image was acquired after excitation with a flip angle α and had a magnitude proportional to sin(α). The second image was acquired after excitation with a flip angle 2α and had a magnitude proportional to sin(2α). The ratio of the two acquisitions was formed giving:

$$\frac{\sin \,\alpha }{\sin \,2\alpha }=\frac{1}{2\,\cos \,\alpha }$$
(6)

from which the local flip angle α was calculated.

ROI Analysis

We used Johns Hopkins University (JHU) ICBM-DTI-81 WM labels atlas61,62 and the automated anatomical labeling (AAL) atlas63,64 to define WM and GM ROIs, respectively. The JHU ICBM-DTI-81 WM labels atlas comprised 48 WM ROIs; AAL comprised 116 ROIs including 12 subcortical GM ROIs. Even though MVFSyMRI and T1w/T2w ratio were in an identical space with the same resolution and slice thickness, MTsat had a different resolution and slice thickness. To ensure that ROIs were placed in the same anatomical position in these different spaces, we warped the above ROIs to each metric map.

For generating the warp field to convert ROIs in the template space to each subject’s space, we first used the FMRIB Software Library (FSL) linear and nonlinear image registration tool (FLIRT and FNIRT)65,66 to register synthetic T1w and 3D T1w images to the MNI152 template. The generated warp fields were saved and inverted so they could be applied to all ROIs, including the eye and temporalis muscle masks. Next, to remove the partial volume effects from other tissues, we segmented synthetic T1w and 3D T1w images into WM, GM, and CSF using FMRIB’s Automated Segmentation Tool (FAST)67. These segmented images of WM and GM were used as masks and applied to MVFSyMRI, T1w/T2w ratio, and MTsat. These tissue masks were thresholded at 0.95 to make sure that the masks contained WM or GM with a probability of 0.95 or higher. WM plus GM tissue masks were also made and thresholded at 0.95. For MVFSyMRI and T1w/T2w ratio, we used tissue masks based on the synthetic T1w images; for MTsat, we used tissue masks made from 3D T1-weighted images. For applying the ROIs from the JHU ICBM-DTI-81 WM labels atlas, we used MVFSyMRI, T1w/T2w ratio, and MTsat masked by WM tissue masks. For applying the ROIs from the AAL atlas to cortical GM, we used MVFSyMRI, T1w/T2w ratio, and MTsat masked by GM tissue masks. For applying the ROIs from the AAL atlas to subcortical GM (e.g., thalamus), we used MVFSyMRI, T1w/T2w ratio, and MTsat masked by GM plus WM tissue masks, because many parts of subcortical GM were segmented as WM by FAST. After warping, all ROIs were inspected for gross registration errors. Upon ROI analysis, the mean values were recorded for further analysis. Examples of ROI placement are shown in Fig. 2.

Figure 2
figure 2

Examples of ROI placement are shown for a 56-year-old female subject. (A and B) show transformed ROIs overlaid on 2D synthetic and 3D T1-weighted images in the same subject, respectively. Transformed ROIs for cortical GM and WM were masked by GM and WM partial volume maps thresholded at 0.95, respectively. For subcortical GM ROIs, GM plus WM partial volume maps thresholded at 0.95 were used for masking. For analysis, ROIs transformed to 2D synthetic T1-weighted images were applied to MVFSyMRI and T1w/T2w ratio, and ROIs transformed to 3D T1-weighted images were applied to MTsat.

Calibration of MVF

Even though SyMRI directly estimates MVF of a voxel, MTsat and T1w/T2w cannot be used as quantitative myelin markers as they are. For calibration of MTsat and T1w/T2w ratio to be used for quantifying myelin in the brain, we assumed a linear relationship between MVFSyMRI, MTsat, T1w/T2w ratio, and actual myelin content, as described previously for MTsat68. In the brain, not only myelin, but also other microstructures contribute to the values of MTsat and T1w/T2w ratio. However, if we assume a linear relationship between MTsat or T1w/T2w ratio and actual myelin content, MTsat or T1w/T2w ratio would also correlate linearly with non-myelin microstructures. Hence, the intercept of the regression line of actual myelin on MTsat or T1w/T2w would be near to zero. Since several studies have calibrated scaling factors of myelin sensitive metrics by healthy WM25,49,68, we also decided to calibrate MTsat and T1w/T2w ratio by values of WM. We determined the scaling factors of T1w/T2w ratio and MTsat by making the means of these values in all the 48 WM ROIs equal to the mean MVFSyMRI. We denoted calibrated MTsat and T1w/T2w ratio as MVFMTsat and MVFT1w/T2w, respectively. Maps of MVFMTsat, MVFSyMRI, and MVFT1w/T2w are shown in Fig. 3. After calibration, we performed ROI analysis again for MVFT1w/T2w and MVFMTsat as described in the previous section and mean values were recorded.

Figure 3
figure 3

MVFMTsat, MVFSyMRI, and MVFT1w/T2w maps of the same subject as Fig. 2 are shown. Because MVFMTsat and MVFT1w/T2w were calibrated for their mean in the whole WM to be equal to the mean MVFSyMRI, these maps look similar to each other in WM. On the contrary, these maps show great variability in GM, with MVFSyMRI showing the highest contrast between WM and GM, and MVFT1w/T2w showing the lowest contrast between WM and GM.

Statistical analysis

For MVF values, normality was tested with the Shapiro-Wilk test. All of the datasets were not normally distributed; therefore, we used the Steel-Dwass test, which is a nonparametric test for multiple comparisons, to compare the contrast among WM and cortical GM, and WM and subcortical GM for the three MVF metrics, and used Spearman’s rank order correlation coefficient to investigate the correlation among MVF metrics for WM, subcortical GM, and cortical GM. Spearman’s ρ correlation coefficients were classified by using the following definitions: 0–0.30, very weak; 0.30–0.50, weak; 0.50–0.70, moderate; 0.70–0.90, strong; and 0.90–1.00, very strong69. Comparison of correlation coefficients among MVFMTsat vs. MVFSyMRI, MVFMTsat vs. MVFT1w/T2w, and MVFSyMRI vs. MVFT1w/T2w were performed in WM, subcortical GM, and cortical GM. This was performed with the Z test for the equality of the two correlations after Fisher r-to-Z transformation70. In addition to analyzing each segment as a whole, we also performed correlation analysis in individual structures representative of WM (genu of corpus callosum, splenium of corpus callosum, anterior limb of internal capsule, posterior limb of internal capsule, anterior corona radiata, superior corona radiata, posterior corona radiata, posterior thalamic radiation, external capsule, and superior longitudinal fasciculus), subcortical GM (pallidum and thalamus), and cortical GM (precentral, postcentral, Heschl, and lingual). Other than corpus callosum, we used bilateral regions aggregately in the analysis. Simple linear regression analysis was performed on the MVFSyMRI and MVFT1w/T2w as a function of MVFMTsat. The regression lines for MVFSyMRI and MVFT1w/T2w were compared by analysis of covariance to determine if they were significantly different from each other in WM, subcortical GM, cortical GM, and all regions combined. All statistical analyses were performed with the software package R, version 3.2.1 (http://www.r-project.org/). A 2-sided p value < 0.05 was considered significant.

Data availabillity

The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.