Estimation of microstructural properties of white matter from multiple orientation GRE signal simulations of realistic models
Renaud Hedouin1, Kwok-Shing Chan1, Riccardo Metere1, and Jose P Marques1

1Donders institute, Radboud university, Nijmegen, Netherlands


This study presents the creation of 2D white matter models, based on real histologically derived axon shapes, with large range of microstructure parameters (FVF, g-ratio). These models are used to simulate the complex gradient echo signal evolution under different main magnetic field orientations for (amongst other parameters) varying magnetic susceptibility and water density in the myelin compartment. A deep learning network, trained from those data, shows its capacity to recover parameter microstructure properties as g-factor and susceptibility on test data.


The gradient echo (GRE) MRI signal evolution is affected both in the magnitude and phase depending on the magnetic susceptibility of its various compartments with respect to the main static field B0. These effects have been successfully simulated with the hollow cylinder model [1]. However, recent work using realistic models of white matter has shown that this simple representation is not accurate enough to simulate the complex GRE signal [2]. In this study, we extend that work and present methods to: (i) generate models of white matter microstructure with different fiber volume fraction (FVF) and g-ratio based on real axon shapes; (ii) use these models to simulate GRE signal; (iii) train a deep learning network and demonstrate its ability to recover relevant microstructure specific parameters.


White Matter Model creation: The segmentation of an electron microscopy image of the spinal cord of a dog [3] was performed using AxonSeg opensource software [4], to obtain a collection of ~600K myelinated axons. To create each model, 400 axons were randomly picked and located in a large 1000x1000 grid. Subsequently, an in-house developed axon packing algorithm derived from [5] was applied to achieve the desired FVF (see Fig 1). In addition, a method to change the g-ratio of the model was used (see Fig 1) resulting in artificial but realistic 2D representation of WM including 3 compartments (axon, myelin, and extra-axonal) with expected axon diameter distribution [6] and different microstructure related parameter (FVF, g-ratio).
The radial arrangement of phospholipids in the myelin sheaths results in an anisotropic component (Xa) in addition to the isotropic component (Xi) of susceptibility (see Fig 2). Its tensorial form was calculated using a method, similar to [2], better able to cope with the non-cylindrical axon geometries. The field perturbation was computed taking into account the B0 orientation with respect to the generated model [7] and then, the complex signal evolution along TE was derived (see Fig 3) assuming a different relaxation time (T2*) and proton density (ρ) specific to each compartment [8].
The complex GRE signal evolution in WM can be used to estimate fiber orientation when the sample is rotated with respect to B0 with more than 7 directions [1]. Yet, DWI is the state-of-the-art method to compute fiber orientation and here, we take advantage of the signal variations due to orientation to map microstructure properties of WM.
The signal was simulated with 9 different orientations that were concatenated into one single vector Stotal as follow (see Fig 4):

Stotal= [θ1, abs(S1), phase(S1), … , θ9, abs(S9), phase(S9)]

where θi is the angle of the fibers with B0, abs(Si) and angle(Si) are the normalized magnitude and phase at 12 different TE’s acquired/simulated with the sample at the i-th orientation with respect to B0.
This signal vector is defined by 6 parameters: FVF, g-ratio, T2_myelin, T2_IntraExtraAxonal, Myelin Water Concentration (ρ) and Xi. The remaining parameters were fixed (Xa = -0.1 ppm, Water concentration of intra and extra water = 1). The dictionary is composed of ~20M vectors, including 8 repetitions with different geometric models created with similar FVF and g-ratio and between 5 and 10 entries for each parameter. Deep learning was performed on this dictionary using Keras [7]. The neural network (2 dense layers, loss function : mse) was trained for multi-regression of the 6 parameters on 7 models and tested on the last model.


The network has demonstrated its ability to recover the parameters with a small error on the test dataset (see Fig 5), particularly when no noise was added neither during the training or on the evaluated signal. Even after adding 1% noise to the model, it was possible to successfully determine the susceptibility difference between myelin and the intra and extra axonal compartments as well as the g-ratio and the T2* of the free water compartment.


Future work will evaluate if these realistic 2D models are a good representation of 3D model where fiber dispersion has to be accounted for. We have shown that some microstructural properties are recoverable using multiple orientation GRE data combined with prior axonal orientation knowledge. Future work will be devoted to the optimization of the sample rotations needed to ensure that the dictionary is better able to differentiate parameters. New parameter spaces will be searched to ensure that FVF and myelin water density, that are currently highly correlated, can be better teased apart. This approaches and deep learning decoding will then be applied to ex-vivo multi-echo data where that rotational freedom exists.


This work is part of the research programme with project number FOM-N-31/16PR1056/RadboudUniversity, which is financed by the Netherlands Organisation for Scientific Research (NWO).


[1] Wharton, Samuel, and Richard Bowtell. "Gradient echo based fiber orientation mapping using R2* and frequency difference measurements." Neuroimage 83 (2013): 1011-1023.

[2] Xu, Tianyou, et al. "The effect of realistic geometries on the susceptibility‐weighted MR signal in white matter." Magnetic resonance in medicine 79.1 (2018): 489-500.

[3] Cohen-Adad, et al. (2018, October 16). White Matter Microscopy Database. https://doi.org/10.17605/OSF.IO/YP4QG

[4] Zaimi, Aldo, et al. "AxonSeg: open source software for axon and myelin segmentation and morphometric analysis." Frontiers in neuroinformatics 10 (2016): 37.

[5] Mingasson, Tom, et al. "AxonPacking: an open-source software to simulate arrangements of axons in white matter." Frontiers in neuroinformatics 11 (2017): 5.

[6] Pajevic, Sinisa, and Peter J. Basser. "An optimum principle predicts the distribution of axon diameters in normal white matter." PLoS One 8.1 (2013): e54095.

[7] Liu, Chunlei. "Susceptibility tensor imaging." Magnetic Resonance in Medicine: An Official Journal of the International Society for Magnetic Resonance in Medicine 63.6 (2010): 1471-1477.

[8] Wharton, Samuel, and Richard Bowtell. "Fiber orientation-dependent white matter contrast in gradient echo MRI." Proceedings of the National Academy of Sciences (2012): 201211075.

[9] Chollet, François, et al. “Keras”, https://keras.io/


Fig 1. Left: Axon packing illustration. 400 axons randomly picked are regularly placed on a 1000x1000 grid, the extra-axonal space is represented in blue and the green axons are surrounded by their yellow myelin sheaths. The axons are attracted to the center of the image and repulse each other to avoid overlap. The packing process occurs to achieve a high FVF value (0,85).

Right: WM Models. Axons are randomly removed from the packed area to reach an expected FVF. Then, keeping the same myelinated axon shapes,the mean g-ratio is modified by dilatation/erosion of the myelin to obtain a model with expected FVF and g-ratio

Fig 2. Left panel: Top left figure represent myelin phospholipid orientation of one axon. The 3 other images represent the simulated field with various B0 orientation (polar angle θ and the azimuthal angle Ф) and fixed susceptibility values (Xi = - 0,1 ppm, Xa = -0,1 ppm).

Right panel: Corresponding magnitude and phase of the signal with different θ and Ф values. As expected, there is an important signal variation as a function of θ, while the Ф has a smaller but non-negligible effect. This observation leads us to take the entire orientation of the B0 magnetic field (θ and Ф) into our model.

Fig 3. Each row represents a model with specific FVF and g-ratio values. First column: model representation. Second column: corresponding field under a perpendicular B0 orientation (θ =90°, ϕ=0°) with Xi = -0,2 ppm and Xa = -0,1 ppm. Two last columns : magnitude and phase of the signal for different Xi values with fixed Xa = -0,1ppm, T2_myel = 0,015 ms, T2_IntraExtraAxonal = 0,05 ms, Myelin water concentration = 0,5.

Fig4. Entire signal simulation of the same 12 TE with 9 different B0 orientations (represented at the top right) for several Xi values. The difference of signal behavior between orientations provide valuable additional information. Arrows on the right side reflect the orientations of the magnetic field while our fiber model was originally oriented along z axis.

Fig5. The table represents the rmse and the normalized rmse of the parameter estimation of noiseless data. The 4 plots show for FVF, g-ratio, T2_axon, and Xi, the RMSE dependence on the noise level (in %) of the test data. Red line shows the performance of the deep learning algorithm with noiseless training data; blue line the noise of the training dataset was matched to the test dataset; green line training dataset had 0.5% noise level. As expected, the same level of noise for training and test data is the best solution. Without prior knowledge, a small level of noise on training data prevents overfitting.

Proc. Intl. Soc. Mag. Reson. Med. 27 (2019)