David Romascano^{1,2}, Jonathan Rafael-Patino^{1}, Muhamed Barakovic^{1}, Alessandro Daducci^{3}, Jean-Philippe Thiran^{1,4}, and Tim B. Dyrby^{2,5}

White matter diffusion MRI enables non-invasive estimation of the axon diameter distribution, which is of interest as it modulates communication speed and delays between brain regions, and changes during development and pathology. Distribution mapping is challenging: current methods simplify it by either estimating the mean diameter, imposing parametric distributions, or combining non-parametric approaches with Double Diffusion Encoding. We present a non-parametric framework based on a PGSE protocol. Simulations show robust reconstruction of unimodal and bimodal distributions. The method is sensitive to population specific changes within bimodal distributions, as long as the underlying populations are separated by a minimum distance.

Reconstruction framework:

We assume the dMRI signal can be expressed as follows:^{3}$$y=f_{ic}\sum_{i=1}^{N_r}\Psi(d_i)S_{cyl}(d_i,d_\parallel,\Omega^\star)+(1-f_{ic})S_{tensor}(d_\perp,d_\parallel,\Omega^\star)$$where $$$f_{ic}$$$ is the Intra-axonal Compartment (IC) volume fraction, $$$\Psi(d)$$$ is the ADD, $$$S_{cyl}$$$ is the dMRI signal for a cylinder^{5} of diameter $$$d_i$$$, parallel diffusivity $$$d_\parallel$$$ and PGSE protocol settings $$$\Omega^\star$$$, and $$$S_{tensor}$$$ is the dMRI signal for a zeppelin^{6} with perpendicular diffusivity $$$d_\perp$$$. The signal is acquired in 60 directions. We focus on ex-vivo imaging and IC signal only, and therefore fixed $$$d_\parallel=0.6\times10m.s^{-1}$$$ and $$$f_{ic}=1$$$. Equation 1 can be expressed as a linear formulation $$$y=Ax$$$, where x can be recovered from the convex inverse problem:^{3,7,8}$$argmin_{x\geq 0}||Ax-y||^2_2+\lambda||\Gamma x_{ic}||_2^2$$Ill-posedness can be reduced by decreasing the mutual coherence of $$$A$$$ (using DDE^{9} for example), by using regularization $$$||\Gamma x||_2^2$$$),^{7} and/or by extracting the mean diameter index $$$a'$$$ instead of $$$\Psi$$$.^{4,8} We propose to combine Laplacian regularization $$$\Gamma=L$$$,^{7} with a PGSE protocol designed to maximize sensitivity to a range of diameters according to our biophysical model.

Protocol design and ADD reconstructions:

The PGSE parameter space is $$$\Omega={G,\Delta,\delta}$$$. We performed a grid-search on $$$\Omega$$$^{10} in order to find a set of 20 shells $$$\Omega^\star$$$ that maximized the sensitivity $$$S'(d_i)$$$ to a set of 20 diameters: $$$S'(d_i)=||S_{cyl}(d_i+\epsilon,\Omega_i)-S_{cyl}(d_i-\epsilon,\Omega_i)||_2$$$. We bounded $$$\Omega$$$ to $$$G_{max}=550mT/m$$$, 2.0<$$$\delta\leq$$$70ms and $$$\Delta\geq\delta+6ms$$$.

We computed the IC signal for cylinders with diameters sampled from different gamma distributions from histological samples,^{4} and reconstructed $$$\Psi$$$ from 100 noisy realisations (Rician noise, SNR=30, $$$\lambda=0.2$$$).

We then created bimodal distributions containing two gaussians: a first fixed population $$$N_1$$$ ($$$\mu_1$$$=4.0$$$\mu$$$m) and a second moving population $$$N_2$$$ with increasing mean $$$\mu_2$$$. We extracted two estimated means $$$\hat{\mu_1}$$$ and $$$\hat{\mu_2}$$$ by fitting a gaussian mixture model to the estimated $$$\Psi$$$. A second experiment tested if the method is sensitive to a reduction of 50% in the amplitude of $$$N_2$$$.

The 20 selected shells maximize the sensitivity to diameters between 4.5$$$\mu$$$m and 10.0$$$\mu$$$m (Figure 1). By maximizing sensitivity to a set of diameters, differences between columns of $$$A$$$ are increased, reducing its mutual coherence.^{3}

Unimodal distributions

Unimodal distributions are reliably reconstructed for all simulations (Figure 2). The mean diameter index^{4} $$$a'$$$ is robustly estimated for values down to almost 1.0$$$\mu$$$m, although the expected lower bound has previously been shown to be around 2$$$\mu$$$m.^{11} The smaller lower bound could be due to prediction or extrapolation properties of the regularization.

Resolving two axonal populations

As shown in Figure 3, $$$N_1$$$ and $$$N_2$$$ should be separated by at least 4.0$$$\mu$$$m for their respective means to be recovered ($$$\mu_1$$$ and $$$\mu_2$$$ are not within the 1st and 3rd quartiles of $$$\hat{\mu_1}$$$ and $$$\hat{\mu_2}$$$ for smaller separations). This minimal separation depends on amplitude and variance of the populations. When the two populations are separated enough, $$$\mu_1$$$ and $$$\mu_2$$$ are well estimated and population specific changes can be detected (Figure 4). Interestingly, preliminary results indicate that removing the last 10 shells of $$$\Omega^\star$$$ compromised bimodal reconstructions while preserving robustness for unimodal distributions, showing that $$$\Omega^\star_{10-20}$$$ provides information for bigger diameters.

1. Y. Assaf, T. Blumenfeld-Katzir, Y. Yovel, and P. J. Basser. AxCaliber: a method for measuring axon diameter distribution from diffusion MRI. Magn. Reson. Med. 2008;59:1347.

2. D. Barazany, P. J. Basser, and Y. Assaf. In vivo measurement of axon diameter distribution in the corpus callosum of rat brain. Brain. 2009;132:1210.

3. D. Benjamini, M. E. Komlosh, L. A. Holtzclaw, U. Nevo, and P. J. Basser. White matter microstructure from non parametric axon diameter distribution mapping. Neuroimage. 2016;135:333.

4. D. C. Alexander, P. L. Hubbard, M. G. Hall, E. A. Moore, M. Ptito, G. J. M. Parker, et al. Orientationally invariant indices of axon diameter and density from diffusion MRI. Neuroimage. 2010;52:1374.

5. P. van Gelderen, D. DesPres, P. C. van Zijl, and C. T. Moonen. Evaluation of restricted diffusion in cylinders. Phosphocreatine in rabbit leg muscle. J. Magn. Reson. B. 1994;103:255.

6. E. Panagiotaki, T. Schneider, B. Siow, M. G. Hall, M. F. Lythgoe, and D. C. Alexander. Compartment models of the diffusion MR signal in brain white matter: A taxonomy and comparison. Neuroimage. 2012;59:2241.

7. K. G. Hollingsworth and M. L. Johns. Measurement of emulsion droplet sizes using PFG NMR and regularization methods. J. Colloid Interface Sci. 2003;258:383.

8. A. Daducci, E. J. Canales-Rodríguez, H. Zhang, T. B. Dyrby, D. C. Alexander, and J.-Ph. Thiran. Accelerated Microstructure Imaging via Convex Optimization (AMICO) from diffusion MRI data. Neuroimage. 2015;105:32.

9. D. Benjamini, Y. Katz, and U. Nevo. A proposed 2D framework for estimation of pore size distribution by double pulsed field gradient NMR. J. Chem. Phys. 2012;137:224201.

10. I. Drobnjak, H. Zhang, A. Ianus, E. Kaden, and D. C. Alexander. PGSE, OGSE, and Sensitivity to Axon Diameter in Diffusion MRI: Insight from a Simulation Study. Magn. Reson. Med. 2016;75:688.

11. T. Dyrby, L. Sogaard, M. Hall, M. Ptito, and D. C. Alexander. Contrast and stability of the axon diameter index from microstructure imaging with diffusion MRI. Magn. Reson. Med. 2013;70:711.

12. R. Colello and M. Schwab. A role for oligodendrocytes in the stabilization of optic axon numbers. J.Neurosci. 1994;14.

13. G. C. DeLuca, G. C. Ebers, and M. M. Esiri. Axonal loss in multiple sclerosis: a pathological survey of the corticospinal and sensory tracts. Brain. 2004;127:1009.

Figure 1: (Left) Positions of 20 selected shells Ω* in the PGSE parameter space. Shells are ordered from Ω*$$$_1$$$ (G=550mT/m, δ=22ms, Δ=28ms) to Ω*$$$_{20}$$$ (G=70mT/m, δ=17ms, Δ=52ms). (Right) sensitivity profile to different diameters for the corresponding shells. Stars indicate sensitivity for the diameter Ω*$$$_i$$$ was selected for. Each color corresponds to a different shell.

Figure 2: (Left) estimated diameter index^{4} $$$a'$$$ compared to the ground-truth index for different distributions. (Middle) mean and standard deviation of the distribution coefficients estimated with our method, and compared to the ground-truth, for the distribution with smallest diameter index (bellow the diameter lower bound). (Right) same as (middle) for the distribution with biggest diameter index.

Figure 3: Boxplots for the 100 estimated μ_{1} and μ_{2}. Each line shows two boxplots summarizing values of μ_{1} (left) and μ_{2} (right) extracted for a given separation between the two populations N_{1} and N_{2} (y-axis). Ground-truth values for μ_{1} and μ_{2} are shown by yellow stars (μ_{1}=4.0μm for all separation distances, while μ_{2}=[5.0μm, 6.0μm, 7.0μm, etc..]). The figure mistakenly shows good performance for distances of 1.0μm and 2.0μm, because N_{1} and N_{2} are actually collapsed into a single peak (wrong ψ) and then split into correct μ_{1} and μ_{2} by the mixture model used to extract the estimated means (see Figure 4, Left).

Mean and variation of ADD coefficients obtained for bimodal distributions. (Left) bimodal distribution with separation of 2μm. (Middle) bimodal distribution with separation of 4μm. (Right) bimodal distribution with separation of 4μm and decreased volume fraction for N_{2}. All plots have a small boxplot on the top summarizing the distribution of estimated a' and the corresponding ground truth (yellow star), illustrating the robustness of a' but also the gain of information provided by mapping the distribution instead of the mean diameter.