Whole-brain ex-vivo quantitative MRI of the cuprizone mouse model

Animal model

All experiments were performed under approval of the local King’s College London ethics committee and the UK Animals (Scientific Procedures) Act 1986, according to the Home Office Project License Number 70/8480 (held by Diana Cash). For a thorough description of the cuprizone model, the reader is referred to Skripuletz et al. (2011). Seven adult male C57BL/6J (Harlan, UK) mice were housed communally and freely fed powdered rodent chow for five weeks, while a further eight were fed powdered chow mixed with 0.2% cuprizone (Bis(cyclohexanone)oxaldihydrazone; Sigma Aldrich, Irvine, UK). All animals were weighed weekly. The five week time-point was chosen as heavy demyelination is expected to be present. They were then killed by transcardiac perfusion with ice-cold heparinized saline followed by 4% buffered paraformaldehyde (PFA). The heads were removed, stored in PFA for 24 h and then rehydrated in Phosphate Buffered Saline preserved with 0.05% sodium azide at 4° for a minimum of 30 days (Cahill et al., 2012).

After image acquisition one control and one cuprizone mouse were excluded from the MRI analysis as the cerebellum had been damaged during sample preparation, preventing a good registration in this region, however the rest of the brain could be used for histology. For quantitative histology, high quality sections could be obtained from all cuprizone animals but only five controls.

MRI acquisition

MR Images were acquired using a 7 Tesla pre-clinical MR system (Agilent Technologies). Samples were immersed in fluorinated liquid to reduce susceptibility artefacts (Galden; Solvay) and loaded four at a time (see Fig. 1) into a 39 mm diameter transmit-receive birdcage coil (Rapid GmbH). Three sets of MR images were acquired; a 3D Fast-Spin Echo for structural analysis, a mcDESPOT protocol, and a DTI protocol.

The 3D FSE image had a matrix size of 256 × 256 × 256 with isotropic 112.5 μm voxels, TE/TR = 40/3,000 ms, an echo-train length of 16, echo-spacing 6.67 ms, readout bandwidth 62.5 kHz, and scan time 3 h 25 min. The DTI scans were acquired using a four shot EPI sequence with a 192 × 128 matrix, 40 slices, 150 × 225 × 500 μm voxel size, TE/TR = 43/4,000 ms, 10 averages and scan time was 3 h 7 min. The diffusion preparation consisted of 30 diffusion directions with b = 2,000 and four b = 0 images, with δ/Δ = 4/16 ms.

The mcDESPOT protocol consisted of a SPoiled Gradient Recalled (SPGR) scan, a balanced Steady-State Free-Precession (bSSFP) scan, and an Actual Flip-angle Imaging (AFI) scan for B1 inhomogeneity correction. The SPGR and bSSFP images both had a 192 × 192 × 192 matrix size with isotropic 150 μm voxels. The SPGR images had TE/TR = 5.148/20 ms, with a read-out bandwidth of 25 kHz, and were acquired at 12 flip-angles (4, 5, 6, 8, 10, 12, 16, 20, 24, 26, 28, 30°), scan time 2 h 28 min, and a strong diffusion spoiling scheme. The bSSFP images were acquired with TE/TR = 3/6 ms, a band-width of 62.5 kHz, also 12 flip-angles (8, 9, 10, 12, 15, 20, 30, 40, 50, 55, 60, 65°) and scan time 2 h 57 min. The extreme flip-angles were chosen to be the optimal for the expected single-component T1 and T2 values (Wood, 2015). The AFI scan had a matrix size of 96 × 96 × 96, isotropic 300 μm voxels, TE/TR1/TR2 = 4.304/20/100 ms, bandwidth 15.625 kHz, a flip-angle of 55° and scan time 37 min (Yarnykh, 2007; Yarnykh, 2010). The acquisition protocol for mcDESPOT allows the calculation of standard relaxometry (T1&T2) maps from the same data.

MRI analysis

The MR Images were first converted to NIFTI format from the manufacturer’s proprietary format, and then were processed using a combination of FSL (Jenkinson et al., 2012), ANTs (Avants et al., 2011) and in-house C++ software utilizing the ITK library, available from https://github.com/spinicist/QUIT. The processing pipeline consisted of several steps. The following operations were carried out in the native space of the acquired images.

1. A Tukey filter was applied in k-space to remove high frequency noise (Tagge et al., 2016). This is a common reconstruction step on clinical scanners but is not implemented on our system.

2. The B1 map was calculated from the AFI scan (Yarnykh, 2007). The single-component T1 map was then calculated from the SPGR data using the B1 map to correct the flip-angles. Because we acquired multiple flip-angles (required for mcDESPOT), we used a non-linear Levenberg-Marquadt algorithm instead of the common linearization method to fit the data. An initial value of 1s was chosen for T₁ and no issues were observed with convergence to local minima.

3. T2 and off-resonance maps were calculated from the SSFP data, B1&T1 maps using the DESPOT2-FM method (Deoni, 2009). To improve the speed of the fitting procedure, instead of the original stochastic fitting method we used a bounded non-linear trust-region optimizer (Zhu et al., 1997). To ensure the global optimum was found four different starting points for off-resonance were ranging from −1/(2*TR) to +1/(2*TR) were tried. The starting point for T2 was set to T1/10. Yarnykh & Yuan (2004) lists values of T2/T1 ranging from 0.044 to 0.065 in various in-vivo and ex-vivo samples, but we observed fewer fitting failures using a larger initial value.

4. The mcDESPOT parameter maps, consisting of the MWF, Intra-Extra cellular Water Fraction (IEWF), Free Water Fraction (FWF), T1&T2 values for each fractional component, and the myelin water residence time (τ_M), were calculated using a Gaussian-prior stochastic Region Contraction method (GRC) (Deoni & Kolind, 2014). The mcDESPOT model is known to be difficult to fit, particularly close to banding artefacts in bSSFP data (Lankford & Does, 2012; Hurley & Alexander, 2014). In order to stabilize the fit, we applied some additional simple heuristics; the off-resonance value of each component was fixed to the value calculated from DESPOT2-FM and we weighted the residuals of the bSSFP data by 3sin²(ϕ + Ψ)/4, where ϕ is the phase-increment value and Ψ is the accrued phase due to off-resonance in each TR. This weighting scheme ensures that data in the banding artefacts of one phase-increment is ignored in preference to data from the other phase-increments. The fittings ranges for the mcDESPOT parameters were chosen by observing the single-component T1&T2 values in regions of white and grey matter and are given in Table 2.

5. FSL topup and eddy were used to remove distortion and eddy current artefacts in the raw diffusion data (Andersson, Skare & Ashburner, 2003; Andersson & Sotiropoulos, 2016). We did not acquire data with a reversed phase-encode direction as the manufacturer’s sequence does not have this option. Instead, we synthesized a standard spin-echo image from the T1&T2 maps using the TE/TR of the diffusion sequence. The resulting image had no distortion and diffusion weighting, and was given to topup as a b = 0 image, which could then undistort the acquired diffusion images. The DTI parameter maps were then calculated and consisted of Mean Diffusivity (MD), Axial Diffusivity (AD), Radial Diffusivity (RD) and Fractional Anisotropy (FA) using FSL dtifit.

The mcDESPOT processing produces ten separate parameter maps (ignoring the B₀ and B₁ parameter maps that correct for field inhomogeneities). However the MWF, IEWF and FWF are defined as fractions that must sum to one, and so are not independent parameters. Hence, of these only the MWF was used for statistical analysis. Of the remaining parameters the myelin water residence time τ_M could potentially be an indicator of myelin sheath integrity. However, as will be shown below the current mcDESPOT methodology cannot reliably fit this parameter, so we did not analyse it further.

The following procedure was then used to split the images into individual subjects and register them to a common space:

1. The structural scan was bias-field corrected (Tustison et al., 2010), thresholded and the four largest connected-components were identified. These components were used as masks to separate each subject from the others. The center-of-gravity of each image was then calculated. Because the subjects were scanned in a consistent orientation (with the base of the brain towards the center of the sample tube), we calculated a simple rotation around the tube axis and translation that moved the subjects to the center of the field of view and oriented them in a standard manner. A simple rigid registration was then performed between each subject and an atlas image (Dorr et al., 2008) to ensure all samples were approximately aligned. This process is illustrated in Fig. 1.

2. A template image was constructed using the 3D FSE and FA images from all subjects in the study (Avants et al., 2010). The resulting 3D FSE template was then non-linearly registered to the atlas image. Including the FA maps in the registration process improved the alignment of the external capsule between groups (see Discussion).

3. All subjects were non-linearly registered to the study templates using their FSE and FA images. Logarithmic Jacobian determinants were calculated from the inverse warp fields in standard space to estimate apparent volume change. The transforms from native to study template, and from study to standard space were concatenated and applied to all relaxometry and DTI parameter maps to align them to the template. These images were resampled to match the voxel size of the template using a Gaussian interpolator. The full-width half-maximum of the interpolator was set to 100 μm for the relaxometry data and 125 μm for the DTI, due to their differing acquisition voxel sizes.

4. A brain parenchyma mask was created from the atlas labels by excluding Cerebrospinal Fluid (CSF) regions. The inverse transforms from the atlas to the study template and from the study template to each subject were applied to calculate the brain volume of each subject.

A group analysis was then carried out on all relaxometry maps, DTI maps and the Jacobian deteriminant images with permutation tests and Threshold-Free Cluster Enhancement (TFCE) using FSL randomize (Smith & Nichols, 2009; Winkler et al., 2014). The brain volume estimates were included as a regressor of no interest in the design matrix when analysing the Jacobian determinants, but not for the parameter maps. Animal weights and total brain volumes were compared using a separate two sample t-Test assuming unequal variance.

Histology

After imaging the brains were removed from the skulls and cryoprotected by immersion in 30% sucrose for at least 72 h. A total of 20 μm thick sections were cut in 12 series on a cryostat, collected onto chrome-gelatin coated slides and stored at −20 °C. Immunohistochemistry was performed on two of the series with washes between each step. After rehydrating in Tris buffer (TBS) for 3 × 5 min endogenous peroxidase activity was blocked by applying 1% hydrogen peroxide (H₂O₂) in TBS for 30 min at room temperature (RT), followed by a non-specific binding block with 10% skimmed milk powder in TBS with 2% Triton-X (TBS-X) for 2 h at RT. Sections were then incubated in primary antibodies for either microglia (rabbit anti-Iba-1, 1:2,000, 019-19741, Alpha Laboratories) or Myelin Basic Protein (rat anti-MBP, 1:1,000, ab7349; Abcam) diluted in 5% skimmed milk powder in TBS-X overnight at 4 °C. Following three washes in TBS-X sections were incubated in either biotinylated goat anti-rabbit or anti-rat antibody diluted in 5% skimmed milk powder in TBS-X (1:1,000; BA-1000 and BA-9400, respectively; Vector Laboratories Ltd.) for 2 h at RT followed by incubation in avidin-horseradish peroxidase complex (Vectastain ABC Elite, PK-6000; Vector Laboratories) for 1 h.

Immunoreactivity was visualized by incubating sections in 0.05% diaminobenzidine and 0.01% H₂O₂ for up to five minutes with exact timing being determined by the depth of colour of the sections. Sections were then rinsed in TBS and dehydrated in increasing concentrations of industrial methylated spirits (IMS) followed by xylene before cover-slipping with DPX mounting medium (Sigma Aldrich, UK).

A third series was stained with Luxol Fast Blue (LFB) to visualize myelin. Sections were placed in a 50/50 IMS/Histoclear solution overnight to remove fat from the tissue and then hydrated in 95% IMS. Sections were then placed into LFB solution at 56 °C overnight, followed by rinses in 95% IMS and distilled water before differentiating in lithium carbonate for 30 s followed by 30 s in 70% IMS and a rinse in distilled water. Differentiation was then checked microscopically and differentiation steps repeated until GM was clear and WM well defined. Once completed, sections were dehydrated in increasing concentrations of IMS followed by xylene before cover-slipping with DPX mounting medium. For analysis of all three stains we used three sections approximately −1.58, −1.82 and −2.06 mm posterior from Bregma.

Iba1 stained sections were analysed with a Zeiss Axioskop2 MOT microscope and design-based optical fractionator probe in Stereoinvestigator software (v11.03.1, MBF Bioscience). An ROI was drawn over the corpus callosum, cortex above and hippocampus below covering an area of 20 mm² in each section using a PlanApo 4× objective. Microglia population estimates were obtained using systematic random sampling with a sampling grid of 300 × 300 μm (Gundersen coefficient of error < 0.1) and a counting frame of 50 × 50 μm. Section thickness was manually defined as 15 μm with a dissector height of 14 μm and guard zones of 0.5 μm at the top and bottom of each frame. All microglia falling with the bounds of the counting frame and not touching the exclusion boundaries were counted using a PlanApo 40× objective at the monitor. Cell density was calculated as the ratio of the population estimate to the volume of the ROI, estimated using the Cavalieri principle (West, Slomianka & Gundersen, 1991).

Sections stained for MBP and LFB were analysed using thresholding techniques. The MBP sections were first pseudocoloured in Aperio Imagescope (v12; Aperio Technologies Inc.) where staining in each pixel of the corpus callosum ROI was categorised into one of three levels and a snapshot at ×4 magnification was then analysed by ImageJ (v 1.50b, NIH) software. The analysis included removing of background and measuring the percent coverage of the appropriate peak corresponding to corpus callosum staining. The same method was applied to LFB but omitting the pseudocolouring step and using Otsu threshold parameters in ImageJ. All results were statistically analysed in Prism v6.07 (GraphPad Software Inc.) using unpaired T-tests.

Validation of mcDESPOT sensitivity to myelination

A major aim of this study was to provide a pre-clinical validation of the MWF as measured by mcDESPOT as sensitive and specific to myelination state. In this regard the study can be regarded as success, principally because in contrast to Thiessen et al. (2013) we found a non-zero MWF in healthy control animals, and then observed a decrease in MWF in cuprizone treated mice.

However, there are some important caveats to this apparent success. The first is that it is not immediately clear what advantages mcDESPOT brings over conventional relaxometry, given the extra acquisition time and extensive processing required to produce the parameter maps. As shown in Fig. 5, the regions of significant change detected in T1 are most constrained to WM, the regions of T2 change extend a long way into the GM, and MWF is somewhere in-between. No regions of significant change were detected in the MWF that were not also detected in either the T1 or T2 map. Due to the co-localisation of inflammation with demyelination in this model, it is difficult to disentangle the impact of these different mechanisms to the MR parameters.

Moreover, because the MWF is defined as a fraction of total water in a voxel, it is obvious that a change in the absolute amount of IE-water will by definition change the MWF, although there has been no change in the absolute amount of myelin water. This means that by definition the MWF can only be sensitive and not specific to myelination state. As currently formulated, due to the need to normalize intensities between the SPGR and bSSFP acquisitions, mcDESPOT cannot be adapted to image absolute myelin and IE-water content, so further work is clearly needed in this area.

The final caveat is that the mcDESPOT model is difficult to fit correctly, and this has been remarked upon elsewhere (Lankford & Does, 2012; Zhang et al., 2014; Bouhrara & Spencer, 2015; Bouhrara et al., 2015). We found a lower MWF than has been reported for in-vivo mcDESPOT human studies, where typical values in WM are over 20% rather than the 10–20% reported here. These lower values are close to those reported by MET2 studies in human in-vivo and ex-vivo studies. However, we found these values to be sensitive to the fitting ranges used, in particular for τ_M, the residence time of water in myelin. As described above, this parameter is meaningless in GM regions and simply converges to the center of any chosen fitting range.

In WM, where the residence time is well defined, we found values below 50 ms. This is significantly shorter than that reported for human studies. This could be attributed to species differences, or the process of PFA fixation, which is known to disrupt biological membranes and introduce holes into otherwise impermeable structures (Zhang et al., 2012). This will increase the rate at which water can move between the myelin and IE-water pools. Hence, we do not think the low residence times are artefactual, but they do conflict with both the widely used MET2 model and the results of Bouhrara et al. (2015), which do not include exchange because it has been assumed that exchange is slow in relation to T2. Given the high CoV for τ_M these values should be treated with some caution and further work is needed to improve the accuracy with which this parameter can be extracted from mcDESPOT data.

The cuprizone model

The above results show that the effect of cuprizone treatment is not limited to the splenium of the corpus callosum. To the best of our knowledge this is the first study to use MRI to show extensive cuprizone-induced changes throughout both the cerebrum and cerebellum. In our data the effect of demyelination was particularly striking in the arbor vitae, which contains the deep cerebellar nuclei surrounded by heavily myelinated tracts. Cerebellar demyelination and inflammation has been previously shown in histological examinations, although changes in the cerebellar cortex were less pronounced and delayed (Groebe et al., 2009; Skripuletz et al., 2010). A recent experiment where rats were administered cuprizone and then serially imaged suggested that demyelination started in the cerebellum and progressed forwards in the brain (Oakden et al., 2016).

We found no difference in total brain volume, despite confirming that cuprizone causes significant weight loss (Nystad et al., 2014). We did not test for behavioural changes in our model. It is well known that despite causing profound changes in the brain, cuprizone treatment causes no large behavioural changes and only subtle motor and cognitive deficits, which would require specialized equipment to test for (Skripuletz et al., 2011).

The effect of cuprizone in the principal areas of the arbor vitae, splenium and external capsule is very strong in our data. In the arbor vitae, the change in T1 is approximately 20%, and in T2 it approaches 50%. This means that despite our low group sizes the statistical analysis appears extremely robust, with many voxels meeting the p < 0.05 threshold even using the conservative FWE multiple comparisons correction. T2 but not T1 detects changes in GM regions such as the cortex and thalamus, which are known to also be affected by cuprizone (Gudi et al., 2009; Fjær et al., 2013; Goldberg et al., 2015). Even more widespread changes in T2 could be seen using FDR instead of FWE correction (available in the on-line dataset), whereas changes in T1 remain restricted to areas that are heavily myelinated in healthy animals.

We also detected changes in several quantitative parameters in the hippocampus, where myelination and inflammatory changes have been well characterised by histology (Skripuletz et al., 2011; Goldberg et al., 2015). This is again a region which most imaging studies have neglected to either analyze or report. Our results indicate the possibility that previous imaging studies to document the longitudinal profile of the cuprizone model have been missing important areas of the brain, and inclusion of these areas in future studies is necessary.

It was beyond the scope of this paper to perform a fully quantitative histological evaluation of all brain areas in which the MRI changes were detected. Many studies to date have shown extensive histological and immunohistochemical changes due to cuprizone treatment, which match both our qualitative (Figs. 6 and 7) and quantitative results (Fig. 8). We limited our myelination quantification with LFB and MBP to the key areas of rostral corpus callosum, and microglial assessment with Iba1 to the surrounding cortex and hippocampus. As expected, these showed a robust and significant decrease in myelination in cuprizone mice, and a profound increase in number of activated microglia. It is also clear from qualitative observations of histological slides that areas of demyelination are spread throughout the brain and particularly salient at the external capsule, hippocampus and arbor vitae. The affected WM areas are accompanied and surrounded by clearly increased Iba1 staining and visibly enlarged microglial cell bodies, indicating ongoing microglial hypertrophy and hyperplasia. The colocation of demyelination and inflammation in the cuprizone model means that it is difficult to rule out a contribution from inflammation to changes in the MR parameters.

We believe this study is the first to estimate brain volumes and apparent volume change in the cuprizone model. Localized volume changes (both increases and decreases) were widespread and extended to regions that showed no change in relaxometry or DTI parameters. Volume increases were located in the corpus callosum, dorsal hippocampus, arbor vitae and sensorimotor cortex. These were corroborated by histology as apparent enlargement and distortion of WM tracts, reflecting previously reported underlying axonal swelling, damage and cellular infiltration (Song et al., 2005; Xie et al., 2010). In addition, we also found bilateral decreases in the frontal, cingulate and retrosplenial cortices as well as the caudate putamen and several other subcortical GM areas.

The lack of change in FA (discussed further below) was beneficial for this study as this parameter could be used to improve the registration process. Initially, we followed standard practice and used only the 3D FSE image for template creation and registration. However, due to the strong effect of cuprizone the FSE contrast in several WM regions is reversed compared to control animals. This led to subtle misalignments between the two groups, principally of the external capsule, resulting in nonsensical changes in volume and FA outside WM tracts. Incorporating the FA maps into the registration correctly aligned the external capsule in the cuprizone group.

This demonstrates a difficulty in using current automated registration methods, as fundamentally they assume that all input images are equivalent and this is clearly not the case when gross pathology or even subtle changes in T1&T2 are present (Cousins et al., 2013). Although non-linear deformation algorithms include regularisation methods to prevent excessive volume changes (van Eede et al., 2013), this was not sufficient in our data to prevent mis-registrations. However, incorporating a parameter that did not demonstrate gross changes between groups (FA) stabilised the registration and yielded a high quality result.

Diffusion measurements in the cuprizone model

DTI has become an extremely widespread method of assessing WM health in clinical studies (Jones, Knösche & Turner, 2013), but interpretation of the diffusivity and FA parameters can be difficult. In particular, RD and FA are often assumed to be a marker of myelin “integrity” (Song et al., 2005; Wheeler-Kingshott & Cercignani, 2009; Janve et al., 2013). Our results indicate two issues; firstly that both RD and AD (and hence MD) are sensitive to demyelination, and that FA appears to be far less sensitive to changes in myelination than diffusivities.

There are conflicting literature results for diffusion changes in the acute phase of the cuprizone model, which are summarised in Table 3. Results from the chronic (greater than six weeks) phase of cuprizone treatment have been omitted, as these do not relate to our time point. However, the reported results also vary widely in diffusion protocols, importantly the number of diffusion directions and strengths. It is possible that earlier papers using only six directions could not correctly calculate the tensor correctly for tracts in all orientations. Some caution should be applied to interpretations comparing diffusion imaging from in-vivo and ex-vivo subjects, because as discussed by Zhang et al. (2012), diffusivity parameters may alter as a result of tissue fixation. Another contributory factor is different lengths of cuprizone treatment, as the levels of demyelination and inflammation in cuprizone undergo progressive change between four and six weeks of treatment (Skripuletz et al., 2011). In particular, at week four demyelination is incomplete whereas microglial activation is highest, but by week five demyelination is complete but microglial activation is decreasing. Hence, as discussed for relaxometry above, there are competing effects of inflammation and demyelination on diffusivity measurements which are hard to disentangle using this model. Recent work (Guglielmetti et al., 2016) supports this hypothesis, as they found an initial decrease of diffusivity after three weeks of cuprizone treatment, followed by an increase at six weeks.

The limited areas of change in FA detected in the cuprizone model were surprising. There is only one area of convincing change, a decrease directly at the splenium of the corpus callosum, with an increase directly above. Comparison to histology suggests that this region of the corpus callosum appears to distort, reducing the ordered nature of the axons. Smaller potential changes are located in the cingulate cortex and on the edges of the arbor vitae, but there are no widespread effects throughout WM tracts as for the other parameter maps. This suggests that FA is insensitive to the myelination state of WM, although this may be due to the increased variance in the FA measurements.

It was beyond the scope of this work to include an assessment of axonal damage via neurofilament histological staining. However recent work indicates that during the acute phase (less than six weeks) of cuprizone treatment axonal damage should be minimal (Schregel et al., 2012; Alme et al., 2015), and only starts to accumulate during the chronic phase (Lindner et al., 2009). This suggests that the changes we see in FA are due to changes in the extracellular matrix or organisation of axons, rather than outright damage.

Comparison of DTI and relaxometry data

Our results are broadly in line with those of Santis et al. (2014), who compared mcDESPOT and DTI in healthy human subjects. They found a low correlation between FA and measures of myelination. They also found a similar trend of parameter variance, with T1 having the best performance and FA & MWF the worst. They found that T2 measurements performed significantly worse than diffusivities, while we found it performed slightly better. We attribute this to the increased number of phase increments in our bSSFP acquisition, which we found markedly decreased banding artefacts.

We calculated Pearson’s correlation ratio between the quantitative parameters across the entire brain and all mice using AFNI (Cox, 2012). Slices through the splenium and cerebellum are shown above and below the diagonal, respectively in Fig. 10. Red colors indicate strong positive correlations, and blue indicates strong negative correlations. Perhaps unsurprisingly, the relaxometry parameters and diffusivities form two separate groups of highly correlated parameters.

T1&T2 are highly correlated in WM regions, and particularly so where cuprizone acts, but are less strongly correlated in GM. They are both negatively correlated to MWF in similar patterns. Correlations between diffusivities and T1&T2 are poor except in areas where cuprizone has a strong effect. Correlations between FA and the other parameters are poor. MD, AD and RD are highly correlated in both GM and WM, except for AD and RD in the cerebral peduncles, fimbria and optic tract. These WM tracts are known to be resistant to cuprizone damage at five weeks (Goldberg et al., 2015), and we observed no qualitative changes in our histology in these regions. Hence, the low correlations are likely a result of different WM organisation in these regions compared to the highly ordered corpus callosum.

This analysis underscores the severity of the cuprizone model at five weeks—in the primarily affected regions, the difference between model and control is so complete that correlations across the entire dataset approach one, and we expect the same would hold true for a histological analysis. It would be useful to better understand the contribution of myelin and inflammation to the different MR parameters, similarly to how Desmond et al. (2016) developed a model of the contribution of myelin and metals to T1&T2. In our opinion, due to the strength of the correlations the single time point presented here is insufficient to build such a model. However, microglial activation is known to precede demyelination in the cuprizone model by approximately a week (Skripuletz et al., 2011), and is required for demyelination to occur (Praet et al., 2015). Hence, a longitudinal study with several earlier time points should observe a differential build up of inflammation and demyelination, allowing the two effects to be separated.