Quantification of individual remyelination during short-term disease course by synthetic magnetic resonance imaging

Abstract MRI is an important diagnostic tool for evaluation of myelin content in multiple sclerosis and other CNS diseases, being especially relevant for studies investigating remyelinating pharmacotherapies. In this study, we evaluated a new synthetic MRI–based myelin estimation in methylenetetrahydrofolate reductase deficiency as a treatable primary demyelinating disorder and compared this method with established diffusion tensor imaging in both methylenetetrahydrofolate reductase deficiency patients and healthy controls. This is the first synthetic MRI–based in vivo evaluation of treatment-associated remyelination. 1.5 T synthetic MRI and 3 T diffusion MRI were obtained from three methylenetetrahydrofolate reductase deficiency patients at baseline and 6 months after therapy initiation, as well as from age-matched healthy controls (diffusion tensor imaging: n = 14, synthetic MRI: n = 9). Global and regional synthetic MRI parameters (myelin volume fraction, proton density, and relaxation rates) were compared with diffusion metrics (fractional anisotropy, mean/radial/axial diffusivity) and related to healthy controls by calculating z-scores and z-deviation maps. Whole-brain myelin (% of intracranial volume) of the index patient was reduced to 6 versus 10% in healthy controls, which recovered to a nonetheless subnormal level of 6.6% following initiation of high-dosage betaine. Radial diffusivity was higher at baseline compared with healthy controls (1.34 versus 0.79 × 10−3 mm2/s), recovering at follow-up (1.19 × 10−3 mm2/s). The index patient’s lesion volume diminished by 58% under treatment. Regional analysis within lesion area and atlas-based regions revealed lower mean myelin volume fraction (12.7Baseline/14.71Follow-up%) and relaxation rates, higher proton density, as well as lower fractional anisotropy and higher radial diffusivity (1.08 × 10−3Baseline/0.94 × 10−3Follow-up) compared with healthy controls. The highest z-scores were observed for myelin volume fraction in the posterior thalamic radiation, with greater deviation from controls at baseline and reduced deviation at follow-up. Z-deviations of diffusion metrics were less pronounced for radial and mean diffusivity than for myelin volume fraction. Z-maps for myelin volume fraction of the index patient demonstrated high deviation within and beyond lesion areas, among others in the precentral and postcentral gyrus, as well as in the cerebellum, and partial remission of these alterations at follow-up, while radial diffusivity demonstrated more widespread deviations in supra- and infratentorial regions. Concordant changes of myelin volume fraction and radial diffusivity after treatment initiation, accompanied by dramatic clinical and paraclinical improvement, indicate the consistency of the methods, while myelin volume fraction seems to characterize remyelinated regions more specifically. Synthetic MRI–based myelin volume fraction provides myelin estimation consistent with changes of diffusion metrics to monitor short-term myelin changes on individual patient level.


DTI preprocessing
MRI data was converted into Nifti format using dcm2niix (Chris Rorden's dcm2niiX version v1.0.20200331) 1 . High resolution 3D T1 weighted images from 3 T MRI were preprocessed by alignment to anterior and posterior commissure, cropping (FSL's robustfov 2 ), intensity normalization and brain extraction (using Freesurfer's autorecon1 pipeline 3 ) to optimize subsequent registration. After preprocessing of the diffusion images including denoising, Gibb's ringing correction, and bias field correction, which was conducted by use of the software MRtrix 4 , non-linear registration of the preprocessed b0 image to the preprocessed T1 weighted image was performed to account for susceptibility induced distortions 5 . For that, the T1 weighted image was resampled to a voxel size of 2 mm to reduce a possible bias by deformation to a resolution of 1 mm. Registration was conducted with ANTS registration pipeline SyNQuick, and the resulting transformation was applied to all diffusion images using ANTSapplyTransforms 6 . These preprocessing steps were performed on both the baseline and the follow-up images, whereby an additional rigid-body registration step (FMRIB's Linear Registration Tool -FLIRT) for the patient's follow-up data was done to align it to the baseline scan.
The resulting DTI maps were further edited using the first steps of the standard Tractography Based Spatial Statistics (TBSS) 7 pipeline in FSL. This included an additional pre-processing script to remove brain-edge artifacts and to zero the end slices, as well as a non-linear registration step to align each individual DTI map to the FMRIB58_FA target in standard MNI152 space. Since the final resolution is 1x1x1 mm by default, the resulting maps were down sampled to a voxel size of 2 mm for conformance with the other patient data. For the healthy control group, all maps were averaged using fslmaths to obtain a custom group template.
To create z-maps in MNI space, it was necessary to add an additional registration step for the patients' DTI data, because the registration to the FMRIB58_FA target within TBSS was not sufficient to achieve a good alignment to HC template. Thus, patients' 3D T1 weighted images were transformed into 2 mm MNI space with linear and non-linear registration (FSL's FLIRT and FNIRT). The resulting warping field was applied to the DTI maps. Finally, z-maps in 2 mm MNI space were obtained by calculation of To create z-maps, patient's SyMRI images were transformed into MNI152 space. For that, synthetic T1 weighted images were brain extracted using FSL's BET and resampled to 2 mm using FLIRT. They were then aligned to a 2 mm MNI template by use of FLIRT and FNIRT.
The resulting warping fields were applied to synthetic MVF, PD, R1, R2 maps.
SyMRI images of the healthy control group were aligned to a 2 mm MNI template to create a study specific template. For that, 3D T1w images, which were additionally obtained from the 1.5T scanner protocol for the healthy controls * , and synthetic T1 images were resampled to 2 mm first. Second, synthetic T1w images were aligned to 3D T1w images (FLIRT, 6 DOF).
Third, the warping field from registration of 3D T1w images to the MNI template (FNIRT) was applied to the synthetic T1w images.
Finally, z-maps were again obtained by calculation of

ROI based analysis
For quantification of DTI metrics, we extracted means and standard deviations for FA, MD, RD, and AD values for the whole brain and within selected regions of interest (ROI) which were derived from the JHU-ICBM-labels atlas. We included posterior thalamic radiation, corticospinal tract, cerebral peduncle, posterior limb of internal capsule, superior corona radiata, and corpus callosum as pathologically affected regions. In order to align the ROIs to DTI maps, FLIRT and FNIRT was used to register the FMRIB58_FA template to the patient's baseline FA map and the resulting warping field was then applied to the JHU-ICBM-labels atlas. ROI mean values of the patient and the healthy control group were obtained using fslstats.
For the assessment of mean MVF, PD, R1, and R2 derived from Synthetic MR maps, the T1 weighted images were registered to the structural T1 template in MNI152 space (FLIRT and FNIRT), and the resulting transformation was inverted and applied to the ROI maps of the JHU-ICBM-labels. For the resulting ROIs in SyMRI space, mean values of patients and the healthy control group were obtained using fslstats.
In addition, we extracted values within the patient's lesion area affected by the symmetrical paraventricular leukoencephalopathy. For that, the lesion map was automatically generated by the Synthetic MR Software, manually edited using FSLeyes to improve its accordance with the lesion based on the FLAIR sequence. The resulting lesion map was aligned to DTI maps, which were previously registered to 2 mm 3D T1 images, by applying the transformation from linear * 192 sagittal slices, repetition time: 2.2 s, echo time: 2.88 ms, inversion time: 900 ms, acquisition matrix: 256x246, voxel size: 1x1x1mm 3 registration of synthetic T1 image to 2 mm 3D T1 image. To compare the affected area to normal references, we registered the lesion mask to MNI space by applying the warping field from previous registration to the MNI target. Lesion values were again extracted using fslstats.