|Home | About | Journals | Submit | Contact Us | Français|
Diffusional kurtosis imaging (DKI) is a clinically feasible extension of diffusion tensor imaging that probes restricted water diffusion in biological tissues using magnetic resonance imaging. Here we provide a physically meaningful interpretation of DKI metrics in white matter regions consisting of more or less parallel aligned fiber bundles by modeling the tissue as two non-exchanging compartments, the intra-axonal space and extra-axonal space. For the b-values typically used in DKI, the diffusion in each compartment is assumed to be anisotropic Gaussian and characterized by a diffusion tensor. The principal parameters of interest for the model include the intra- and extra-axonal diffusion tensors, the axonal water fraction and the tortuosity of the extra-axonal space. A key feature is that these can be determined directly from the diffusion metrics conventionally obtained with DKI. For three healthy young adults, the model parameters are estimated from the DKI metrics and shown to be consistent with literature values. In addition, as a partial validation of this DKI-based approach, we demonstrate good agreement between the DKI-derived axonal water fraction and the slow diffusion water fraction obtained from standard biexponential fitting to high b-value diffusion data. Combining the proposed WM model with DKI provides a convenient method for the clinical assessment of white matter in health and disease and could potentially provide important information on neurodegenerative disorders.
Diffusion weighted imaging (DWI) is a widely applied and clinically important MRI method used to measure the micron-scale displacement of water molecules in the brain. Diffusion on this length scale is very sensitive to the microstructure of neural tissue, being strongly affected by the number, orientation and permeability of barriers (e.g. myelin) and the presence of various cell types and organelles (e.g. neurons, dendrites, axons, neurofilaments and microtubules) (Beaulieu, 2002). Moreover, as the tissue microarchitecture is closely associated with function, DWI offers a unique and very powerful method to study brain pathology.
By far the most widely applied DWI technique to date is diffusion tensor imaging (DTI), in which the apparent diffusion tensor is estimated from the measurement of the apparent diffusion coefficient (ADC) along multiple directions (Basser et al., 1994). Several rotationally invariant diffusion metrics can be extracted from a DTI-analysis, including the mean diffusivity (MD) and the fractional anisotropy (FA), which are both popular markers of white matter (WM) integrity (Pierpaoli and Basser, 1996). In addition, DTI is also a commonly used method for fiber tractography, i.e. reconstructing the pathways of major WM fiber tracts through the brain (Basser et al., 2000). Although DTI is an important technique for investigating mechanisms of health and disease in brain WM (Thomason and Thompson, 2011), among its limitations are the inability of DTI-based fiber tractography to resolve fiber crossings, and the lack of specificity to histological features.
While DWI has the potential to fully characterize the water diffusion properties of the brain, it is well recognized that DTI yields only a fraction of the information potentially accessible with DWI, which is mainly due to the fact that DTI is based upon a Gaussian approximation of the diffusion displacement probability function. Non-Gaussian diffusion is readily observed in the brain when applying diffusion gradients such that the corresponding b-value (diffusion weighting) is significantly higher than the typical DTI b-value of 1000 s/mm2 (Assaf and Cohen, 1998; Niendorf et al., 1996). The non-Gaussian diffusion effects in the brain are believed to arise from diffusion restricted by barriers, such as cell membranes and organelles, as well as the presence of distinct water compartments with differing diffusivities.
Several techniques for assessing non-Gaussian diffusion have been developed (Alexander et al., 2002; Jensen and Helpern, 2010; Liu et al., 2004; Maier et al., 2004; Tuch, 2004; Wedeen et al., 2005). Among them, diffusional kurtosis imaging (DKI) has been proposed as a minimal extension of DTI that enables the quantification of non-Gaussian diffusion through the estimation of the diffusional kurtosis, a quantitative measure of the non-Gaussianity of the diffusion process (Jensen et al., 2010; Jensen et al., 2005; Lu et al., 2006). A typical DKI-protocol for brain requires a maximum b-value of 2000 s/mm2 and DWI measurements along a minimum of 15 different directions (Tabesh et al., 2010). Quantitative rotationally invariant diffusion metrics can be extracted from the DKI-analysis, such as the mean kurtosis (MK), radial kurtosis and axial kurtosis, that are of potential interest to the study of white and gray matter integrity. So far, DKI has shown promising preliminary results for several brain diseases including stroke (Jensen et al., 2010), attention-deficit hyperactivity disorder (ADHD) (Helpern et al., 2010), the staging of glioblastomas (Raab et al., 2010), as well as normal aging (Falangola et al., 2008). Additionally, DKI is potentially useful in tractography for resolving crossing fibers (Lazar et al., 2008). However, similar to DTI, DKI metrics of non-Gaussianity are pure diffusion measures and lack microstructural and pathological specificity. Furthermore, a clear explanation for the microscopic origin of the diffusional kurtosis in WM has not been previously given.
The extraction of cell properties and histological details of WM necessarily relies on biophysical modeling of the DWI signal, and on the subsequent interpretation of the model parameters in terms of intrinsic tissue properties. The most basic model used to analyze high b-value data is the biexponential model that is based on the assumption of two non-exchanging compartments: one exhibiting fast diffusion, and the other slow diffusion. Despite good fits of the DWI signal, the original attempt to assign the two compartments to the intra- and extracellular space is still a subject of a debate (Assaf and Cohen, 1998; Clark and Le Bihan, 2000; Kiselev and Il'yasov, 2007; Maier et al., 2004; Mulkern et al., 2000; Niendorf et al., 1996). However, for the WM, assigning the highly restricted water diffusion inside the axons to the slow compartment and the less hindered diffusion in the extra-axonal space to the fast compartment has been justified experimentally (Assaf and Basser, 2005; Assaf and Cohen, 2000; Assaf et al., 2004) and theoretically (Fieremans et al., 2010b).
A number of advanced morphology-based models have been proposed to interpret DWI in brain WM. As an early and comprehensive model, Stanisz et al. represented bovine optic nerve tissue by three compartments formed by spherical glial cells, prolate ellipsoidal axons and the extracellular space. By using this analytical model, the compartment parameters, such as volume fractions, compartment size, membrane permeability and diffusivity, could be estimated for fixed tissue (Stanisz et al., 1997). The less elaborate CHARMED model (Assaf and Basser, 2005; Assaf et al., 2004) assumes two types of diffusion in the brain: restricted diffusion inside impermeable cylindrical axons and hindered diffusion in the extra-axonal space, allowing estimation of the compartment volume fractions and diffusivities for the human brain in vivo. In the framework of “Axcaliber”, the CHARMED model was further developed to extract the axonal diameter distribution, which was evaluated ex vivo on pig spinal cord (Assaf et al., 2008) and in vivo in the corpus callosum of a rat (Barazany et al., 2009). A similar model of two non-exchanging compartments has been developed by Jespersen et al, wherein the restricted diffusion component arises from an angular distribution of narrow cylinders (representing the axons and dendrites), allowing one to estimate the compartment volume fractions and diffusivities, as well as the intra-voxel distribution of fiber orientations, as demonstrated in fixed brain tissue of the rat and baboon (Jespersen et al., 2010; Jespersen et al., 2007). Recently, Alexander et al. proposed a four-compartment brain WM model that allows the axon diameter and density to be derived, as illustrated in fixed monkey brain and in vivo human brain (Alexander et al., 2010).
To extract all the features of the models summarized above, DWI data are needed for several high b-values (i.e., b ≥ 3000 s/mm2), multiple diffusion gradient directions and/or different diffusion times, which necessitates long scan times and limits the applicability of these models for most clinical studies. Alternatively, DKI is a clinically feasible technique with acquisition times only a few minutes longer than conventional DTI. However, as DKI metrics of non-Gaussianity are model-independent, they must be augmented with a tissue model to help interpret the physical meaning of any changes associated with disease processes.
In this work, we focus on WM regions consisting of more or less parallel aligned fiber bundles and propose a model of diffusion in the WM that is suitable for DKI analysis and provide a more meaningful physical interpretation of DKI diffusion metrics in WM. We first introduce the WM diffusion model of two non-exchanging compartments: the intra-axonal space, consisting of impermeable cylindrical axons (IAS), and the extra-axonal space (EAS). Next, we demonstrate how the diffusion in each compartment appears to be Gaussian for the b-values typically used in DKI and hence can be described by compartment-specific diffusion tensors. Combining this model with DKI provides analytical expressions for the intra- and extra-axonal diffusion tensors, and allows for quantification of the axonal water fraction (AWF) and of the tortuosity of the EAS. We use then the newly proposed model to characterize human brain WM in vivo and discuss the biological significance of the tissue parameters as derived using DKI. Finally, we compare the AWF obtained from DKI-analysis to the slow diffusion fraction obtained from conventional biexponential fitting to high b-value diffusion data.
In this section, we first describe the WM model, discuss the validity of its assumptions, and outline its relation to DKI. We then specify the datasets, imaging protocols and data processing needed to evaluate the derived WM parameters.
In this work, we make the following two main assumptions about in vivo diffusion in brain WM:
Based on these assumptions, both the IAS and EAS are modeled by the corresponding compartmental diffusion tensors Da and De. The axonal water fraction (AWF), denoted by the symbol f, is the fraction of MRI visible water in the axons relative to the total visible water signal. The DWI signal intensity, S, in the direction n, as a function of the diffusion weighting, b, in such a system is then described by
Below, we first justify Eq. (1) and then outline the procedure to extract the WM parameters based on the DKI representation of the DWI signal.
We substantiate here the assumption made in Eq. (1) of Gaussian diffusion in the EAS and the IAS. For the EAS, the mean square displacement of water molecules during a typical DWI experiment (, where De is the free extra-axonal diffusivity ) is much greater than the correlation length of the axonal packing (~1 μm), so that the tortuosity asymptote is reached (Fieremans et al., 2010b). Hence, it is reasonable to assume that the diffusion profile in the EAS at long times, while anisotropic, is Gaussian in every direction and can be represented by the diffusion tensor De.
The diffusion in the IAS is highly restricted and in general non-Gaussian. However, the fact that there is no observed time dependence of the diffusion metrics at clinically available diffusion times (Clark et al., 2001) indicates that any effects of non-Gaussian diffusion due to barriers become negligible in the direction along the axons at long diffusion times. It is then reasonable to assume Gaussian diffusion in the intra-axonal compartment for a voxel consisting of perfectly aligned axons with zero radius as the transverse diffusivity is zero in this limit (Fieremans et al., 2010b). Hence, the DWI signal S becomes trivially Gaussian in all directions and is described by
where Da is the free intra-axonal diffusivity and θ is the angle with respect to the axon axis. The zero radius approximation is valid for long diffusion times t >> R2/Da, where R is the axon radius. For typical values, R ≈ 1 μm and Da ≈ 1 μm2/ms, we find t >> 1 ms, which is well satisfied for clinical DWI experiments, where t ~ 50ms.
When the axons are not perfectly aligned within a voxel, Gaussianity becomes dependent on both the direction distribution of the axons and the b-value. A general criterion for the diffusion to appear Gaussian in a given direction can be derived based on the cumulant expansion (Jensen et al., 2005; Kiselev, 2010),
where S is the DWI signal and D, K are the apparent diffusion coefficient and diffusional kurtosis in that direction.
Whereas for the EAS and the case of perfectly aligned axons described earlier, the higher order cumulants vanish due to zero K or D, this is not necessarily true for any axonal geometry. To a good approximation, a compartment will appear Gaussian if the b-value
so that the 2nd and (presumably) higher order terms in b become negligible in Eq. (3). As a practical matter, we adopt the criterion
as the condition for being effectively Gaussian. This criterion can be used as a guide to whether the diffusion appears Gaussian in a specific axonal geometry for a given b-value.
As an example, consider a voxel that contains two crossing fiber bundles with equal volume fractions, oriented at polar angles (θA,θB) relative to a particular diffusion direction of interest, as illustrated in Fig. 1 (a). For such a system, D and K for the axonal compartment can be analytically derived (see Appendix A) and the b-value condition for Gaussianity (5) can then be written as
where Da is the free intra-axonal diffusivity. This condition is evaluated numerically in Fig. 1 (b), showing a substantial region of approximate Gaussianity in (θA,θB) space. In particular, for two axonal bundles intersecting at an angle ≤ 30 degrees, this axonal geometry will look essentially Gaussian for a typical DKI protocol (i.e., a maximum b-value of 2000 s/mm2) in all directions. However, this axonal geometry may not necessarily be accurately modeled as a Gaussian compartment for higher b-values.
In the more general case of a voxel containing axon bundles that intersect at larger angles or for randomly oriented axons in a plane, the diffusion in such geometries will not be Gaussian anymore in all directions. Yet, the derivation in Appendix A illustrates that an approximately coplanar axonal geometry with an out of plane angular spread of up to ± 45° will, for a typical DKI protocol, be effectively Gaussian for diffusion in the direction perpendicular to the plane.
Directly fitting the DWI signal for many directions simultaneously to Eq. (1) is a hard multi-parameter non-linear problem, whose solution is numerically challenging and may in practice be unstable to noise. An alternative, based on the biexponential fitting for each diffusion direction separately (Maier et al., 2004), has been shown to require high b-values which are typically above those utilized clinically (Kiselev and Il'yasov, 2007). Here we suggest utilizing the DKI metrics, which are determined by a straightforward linear fitting procedure, together with a set of relationships that connect these to the AWF and the tensors Da and De. In what follows next, we provide analytical expressions for the AWF, Da and De based on the diffusion and kurtosis tensor. A detailed derivation is given in Appendix B.
On a voxel-by-voxel basis, the AWF is estimated by
where Kmax is the maximum kurtosis over all diffusion directions (based on the kurtosis tensor). Eq. (7) assumes there is a direction for which the IAS diffusivity Da,i=0, while for the other directions the diffusion in the IAS appears Gaussian, according to criterion (5). This approximation for the AWF has a broad applicability, e.g. the derivation in Appendix A illustrates that in a quasi-coplanar axonal geometry the kurtosis as measured in the direction perpendicular to the fiber plane is typically the maximum kurtosis and thus indicates the water fraction of the restricted compartment.
When the assumption of above (Da,i=0) is not precisely met (e.g. imperfect alignment of the axons or due to effects of finite axonal radius relative to the diffusion length), Eq. (7) becomes a lower bound for the AWF. In that case, a more accurate estimate of the AWF is given by
where D1, D2, D3 are the diffusion coefficients and K1, K2, K3 the corresponding kurtoses along the axis directions of a chosen reference frame. We note that the term is rotationally invariant as discussed in Appendix B. Da is the axonal diffusion coefficient with Da,min being a lower bound for Da, defined by:
Da can then be estimated by selecting the maximum Da,min over all voxels within a specified region of interest (ROI). This estimate will be exact if the ROI contains at least one voxel wherein Da,i=0 in a certain direction. This procedure assumes that Da is the same for all the voxels within the ROI.
Alternatively, the AWF can be found from a conventional biexponential fit to high b-value diffusion data, as in many previous studies (Inglis et al., 2001; Maier et al., 2004; Mulkern et al., 2000; Niendorf et al., 1996), which will be performed in this study to evaluate our proposed WM model parameter estimates based on DKI metrics. The biexponential fitting approach requires, in practice, much higher maximum b-values than DKI (7000 s/mm2 rather than 2000 s/mm2), and hence, according to criterion (5), is only applicable in regions with nearly aligned fiber directions.
With an estimate for the AWF, given by Eq. (7) or (8), the IAS and EAS compartment diffusion coefficients in a given direction, Da,i and De,i, can be derived from the diffusion coefficient Di and kurtosis Ki in that direction by:
We note that Eqs. (10) and (11) have the underlying assumption that Da,i ≤ De,i, which is justified by the arguments described in Appendix B and Fig. B.1(a-b). It also follows from Appendix B and Fig. B.1(c) that Eqs. (10) and (11) are true in any direction. Hence, by choosing 6 or more independent directions, we can then reconstruct (using the standard DTI method) the full diffusion tensors Da and De.
By applying this WM model to the DKI metrics, several compartment-specific metrics can be derived, in addition to the usual DTI/DKI parameters. The same rotationally invariant measures as derived in standard DTI can now be obtained specifically for each compartment tensor. For this study, we focus on the following WM metrics of interest:
All of these WM metrics assume that Eq. (1) is valid in the considered voxel, as discussed earlier in section 2.1.2. This implies that Eq. (12) is a good approximation for the intra-axonal along axis diffusivity in voxels that contain fiber bundles with an angular spread of less than 30 degrees. In addition, the definition of the tortuosity (Eq. 15) is only meaningful in the specific case of a voxel containing a single fiber direction so that De, (Eq. 13) represents then the EAS diffusivity in the direction along the fibers and De, (Eq. 14) the EAS diffusivity perpendicular to the fibers.
In what follows, we will estimate and evaluate each of these WM parameters in healthy volunteers. For the AWF, we will compare the estimated values with the values obtained from a biexponential fit.
Test data sets for evaluating the model come from three healthy human volunteers. Subject 1 is a 27 year old female, subject 2 is a 28 year old female, subject 3 is a 28 year old male. The subjects were scanned with informed consent obtained as approved by our Institutional Review Board.
MRI images were acquired using a 3T wide-bore Siemens Verio system (Siemens Medical Solutions, Erlangen, Germany) with a transmission body coil and a 12-channel head coil for reception.
Whole brain T1-weighted magnetization-prepared rapid gradient echo (MPRAGE) images were acquired with TR = 2200 ms, TI = 1100 ms, TE = 2.3 ms, matrix = 256 x 256, FOV = 256 × 256 mm2, slices per slab = 192, slice thickness = 1 mm, band width = 260 Hz/pixel in a total time of 4:30 min.
Diffusion-weighted images were acquired along 30 gradient directions for b = 0, 1000, 2000 s/mm2 with a twice-refocused spin-echo (TRSE) echo planar imaging sequence (Reese et al., 2003) with TR = 8700 ms, TE = 96 ms (corresponding diffusion time ~ 50 ms), matrix = 82 × 82, FOV = 222 × 222 mm2, 40 slices, slice thickness = 2.7 mm, no gap, NEX = 11 for b = 0, NEX = 2 for b = 1000, 2000 s/mm2, band width = 1355 Hz/pixel in a total time of 20 min.
Additional DWIs were acquired along one direction for 16 b-values (0, 100, 500 to 7000 in increments of 500 s/mm2). The gradient direction was the same as the direction of the slice selection, which was chosen perpendicular to the anterior commissure-posterior commissure (AC-PC) plane. A standard Stejskal-Tanner (Tanner and Stejskal, 1968) sequence was used instead of the TRSE diffusion preparation to minimize the echo time (TE = 132 ms) and repetition time (TR = 12,500 ms) in order to obtain a higher SNR and shorter scan time, respectively. Other imaging parameters were: matrix = 104 × 104, FOV = 280 × 280 mm2, 40 slices, slice thickness = 2.7 mm, no gap, NEX = 8 for all b-values, band width = 1502 Hz/pixel. The total acquisition time was 36 min.
Data preprocessing with SPM8 (Statistical Parametric Mapping, Wellcome Department of Imaging Neuroscience, University College London, UK) included 3D motion correction by aligning all DWIs to the first b = 0- image and removing cerebrospinal fluid (CSF) from the DWIs by applying a binary mask derived from thresholding the CSF probability map (created with SPM from segmentation of the b = 0 image) at a cut-off probability value of ≤ 0.2.
The co-registered and masked DWI images acquired in 30 directions for b = 0, 1000, 2000 s/mm2 were then further processed using in-house software (Diffusional Kurtosis Estimator (DKE) (Tabesh et al., 2010)) running in Matlab and the diffusion and kurtosis tensors were calculated on a voxel-by-voxel basis using a weighted linear least squares (LLS) fitting algorithm. The AWF was derived based on Eq. (7) where Kmax is taken as the maximum value of calculated kurtosis-values along 10,000 randomly chosen directions based on the kurtosis tensor. The generated f-map, together with the diffusion and kurtosis tensor element maps were then used to derive parametric maps of the diffusion coefficients of the axonal and extra-axonal tensors Da and De along the same 30 directions that were used for the DWI acquisition based on Eqs. (10) and (11). Next, the eigenvalues of Da and De were obtained using the standard DTI method (Basser et al., 1994) and the specific WM parameters of interest Da, De,, De, and α were derived based on Eqs. (12)-(15).
The idealized WM model described above assumes that the axons are more or less parallel and becomes less appropriate in regions of complex fiber architecture such as fiber crossings of significant fanning. Additionally, to justify the use of a biexponential fit to the signal decay averaged over an ROI (see further in section 2.3.3), precise unidirectionality of the fibers in the specific ROI is required. Therefore we limit this study to voxels in which a single fiber orientation is expected corresponding to straight parallel fibers. We computed the coefficients of linearity , planarity , and sphericity , where λ1, λ2, and λ3 are the eigenvalues of the overall diffusion tensor D (λ1 ≤ λ2 ≤ λ3), that describe how close D is to the generic cases of line, plane and sphere (Westin et al., 1999). For further processing we included only voxels that fulfilled
Fiber bundles oriented in a plane perpendicular to the slice direction were identified for further analysis. ROIs were drawn on the color-coded FA map considering only voxels that are colored in red (left-right direction) and/or green (anterior-posterior direction). In addition, ROIs for which the mean values of cL, cP and cS did not fulfill criterion (16) were excluded. ROIs meeting the criteria described above were the genu, midbody and splenium of the corpus callosum and selected parts of the left and right part of forceps major (fmajor) and the inferior fronto-occipital fasciculus (ifo).
For each ROI and b-value (b = 0 – 7000 s/mm2), the signal of the DWIs acquired along the slice direction was averaged. The mean signal S(b) was fitted by a biexponential function
with Ds ≤ Df. The fitting procedure was performed in Matlab using the trust region algorithm with the robust option set to bisquare, enabling a fit that minimizes the summed square of residuals, and down-weights outliers using bi-square weights. The fitted slow and fast diffusion coefficients, Ds and Df, are then compared to the compartment diffusivities along the slice direction, Da,slice and De,slice. Similarly, the fraction fbiexp corresponding to Ds is compared to the AWF of our DKI WM model. The AWFs in the voxels of each ROI were estimated by either Eq. (7), or by Eq. (8) where Da is selected as the maximum Da,min over the whole ROI. To reduce the effect of noise, the maximum Da,min was approximated as max(Da,min) ≈ mean(Da,min) + SD(Da,min), where SD is the standard deviation over the whole ROI.
This section shows the results of testing the new WM DKI model; first we evaluate the parametric maps of three healthy young adults, then we show the comparison between the DKI-derived AWF and the slow diffusion fraction as obtained from biexponential fitting in WM regions where the biexponential model is physically justified.
The low b-value diffusion data acquired over multiple directions were used for testing the DKI WM model. Parametric maps are shown in Fig. 2 for the WM metrics of interest (the AWF based on Eq. (7); the axonal diffusivity, Da; the axial EAS diffusivity, De,; the radial EAS diffusivity, De,; and the tortuosity of the EAS, α) for an axial slice through the genu and splenium of the corpus callosum of one subject. The AWF is highest in the splenium of the corpus callosum. The IAS diffusivity Da is found to be smaller than the axial EAS diffusivity De,, and both are observed to be the highest in the corpus callosum (~1.2 μm2/ms and ~ 2.5 μm2/ms respectively), while the radial EAS diffusivity De, is found to be slightly lower in those regions. Consequently, the tortuosity of the EAS, as defined by Eqn. (15), is the highest in these regions of the corpus callosum.
The WM parameters are shown to be very reproducible between subjects. The histograms of the WM parameters of the same subject as in Fig. 2 are shown in Fig. 3. The parametric maps and histograms of the other two subjects appear to be very similar, which is demonstrated in Fig. 4, plotting the mean values and standard deviation of the WM parameters over all selected WM-voxels for each subject. Over all subjects, average values and standard deviations were found for the AWF, f = 0.49 ± 0.07, the axonal diffusivity, Da = 0.99 ± 0.18, the axial EAS diffusivity, De, = 2.26 ± 0.31, the radial EAS diffusivity, De, = 0.87 ± 0.16, and the tortuosity of the EAS, α = 2.75 ± 1.13.
For all the results shown in Figs. 2--4,4, only voxels were selected with a single fiber orientation according to criterion (16), which were respectively 19%, 18% and 24% of all WM voxels for the three subjects. In addition, a small number of voxels in the corpus callosum near the ventricles (< 0.3% of all WM voxels) was excluded from the analysis in which either D3 is artifactually small (< 0.08 μm2/ms) or De, is artifactually high (> 3 μm2/ms = diffusion coefficient of water at 37°C). To derive these ratios, a mask of all WM voxels was created based on the segmented MPRAGE-images (cut-off probability of ≤ 0.99) that were co-registered and resliced to the b = 0 – image.
The additional high b-value diffusion data acquired in the slice direction were used to compare the DKI WM model to the biexponential model. ROIs were chosen such that the voxels were in fiber bundles perpendicular to the diffusion gradient, i.e. the corpus callosum, the forceps major and parts of the inferior fronto-occipital fasciculus. The biexponential model, given by Eq. (17), fits the ROI-averaged b-dependent signal up to 7000 s/mm2 very well for all ROIs (R2 ≥ 0.997), as illustrated in Fig. 5 for one subject.
The theory predicts that the fitted slow and fast diffusion coefficient, Ds and Df, correspond to the IAS and EAS compartment diffusivities in the slice direction, Da,slice and De,slice, respectively. The fitted slow diffusion coefficient Ds is on average 0.06 ± 0.02 μm2/ms, whereas the fast diffusion coefficient Df is on average 0.75 ± 0.17 μm2/ms over all selected ROIs. As a comparison, Da,slice is on average 0.08 ± 0.02 μm2/ms and De,slice is on average 0.86 ± 0.15 μm2/ms over the same ROIs as used in the biexponential analysis.
Similarly, the fraction fbiexp of the slow diffusion compartment in Eqn. (17) can also be compared to the AWF derived using our DKI WM model, approximated either by fKmax in Eq. (7) or by fDa in Eq. (8). To assess agreement between the parameters obtained from the biexponential model on the one hand and the DKI WM model on the other hand, Bland-Altman plots of the difference against their mean were created for each parameter set in Fig. 6 (Bland and Altman, 1986). The Bland-Altman plot for (Ds, Da,slice) is plotted in Fig. 6(a), showing a mean bias of -0.03 μm2/ms and lower and upper limits of agreement of - 0.12 μm2/ms and 0.05 μm2/ms. The Bland-Altman plot for (Df, De,slice) is plotted in Fig. 6(b), showing a mean bias of 0.22 μm2/ms and lower and upper limits of agreement of -0.05 μm2/ms and 0.49 μm2/ms. The Bland-Altman plot for (fbiexp, fKmax) is plotted in Fig. 6(c), showing a mean bias of 0.05 and lower and upper limits of agreement of -0.02 and 0.12, and the Bland-Altman plot for (fbiexp, fDa) is plotted in Fig. 6(d), showing a mean bias of 0.01 and lower and upper limits of agreement of -0.04 and 0.07.
DKI has been introduced as a clinically feasible technique to study restricted diffusion in brain WM (Jensen and Helpern, 2010; Jensen et al., 2005). It is a fast and robust method that quantifies the diffusion coefficient and diffusional kurtosis, both physically well-defined, model independent, diffusion metrics. In this work, we provide a physical interpretation for the diffusional kurtosis by augmenting the DKI metrics with an idealized model for WM tissue, as described in detail in section 2.1.1. The WM is modeled here by two non-exchanging compartments, the IAS and the EAS, as captured by Eq. (1). Similar multi-compartment models based on the IAS and EAS are commonly used with high b-valuediffusion data in brain WM (Alexander et al., 2010; Assaf et al., 2004; Jespersen et al., 2007; Sen and Basser, 2005) and are shown to explain simulated and experimental data very well (Fieremans et al., 2010b; Panagiotaki et al., 2009). We also note that other bi- and multiexponential models have been proposed that take a form similar to Eq. (1), but model the diffusion originating purely from the IAS consisting of multiple fiber orientations (Anderson, 2005).
Fitting Eq.(1) directly to the raw diffusion data is a hard non-linear problem that has, to our knowledge, not been attempted and would likely be confounded by multiple local minima in the parameter space. We provide here a straightforward and practical method to solve Eq. (1) by utilizing the DKI metrics. In addition, our DKI-approach facilitates selecting the one solution as being the most physically reasonable, as the derivation in appendix B shows how the extra information contained in the kurtosis tensor allows one to resolve the sign ambiguity of Eq. (B3 - B4) and conclude that De,i ≥ Da,i.
Applying model (1) to in vivo human WM DKI data yields estimates of the AWF and the compartment-specific diffusion tensors, from which scalar parameters of potential interest can be derived, such as the radial and axial intra- and extra-axonal diffusivities and the tortuosity of the EAS. All these parameters values are consistent between subjects (Fig. 4) and agree with prior studies as discussed below in more detail. We also discuss how this approach can be useful in clinical applications, and finally we address some potential limitations.
The axial compartment diffusivities, Da and De,, provide estimates for the intrinsic intra- and extra-axonal diffusion coefficient of water in human WM. The axial diffusion is found to be higher in the EAS than in the IAS, which follows from Appendix B and is illustrated in Fig. B.1. This observation can be understood by the presence of cytoplasm and organelles in the IAS that plausibly slows down the diffusion. A similar experimental finding has been reported in fixed rat brain WM (Jespersen et al., 2010), in which a two-compartment model was used to interpret the high b-value diffusion data.
The best estimate for the “free” intra- and extra-axonal diffusion coefficients of water in human WM are probably the values of Da and De, in those voxels of the corpus callosum with the strongest fiber alignment. The IAS diffusion coefficient Da is then about 1.2 μm2/ms, resulting in 0.4 for the ratio of Da relative to the free diffusion constant of water at 37°C. Remarkably, Kroenke et al. (2004) obtained the same ratio for the in vivo measured parallel diffusion coefficient of NAA in the corpus callosum relative to its diffusion coefficient in dilute aqeous solution. That study also reports a higher ratio of 0.46 as measured in large voxels within rat brains. Differences might be due to imperfect fiber alignment within the voxel, and/or the fact that NAA diffuses much slower than water (0.36 μm2/ms versus 1.2 μm2/ms). As a result, NAA has a shorter diffusion length (6 μm) than water (11 μm), and its diffusion is thus less hindered by intra-axonal restrictions such as axonal varicosities, that have an overall mean spacing of 5.2 μm (Shepherd and Raastad, 2003). The IAS diffusivity in the direction perpendicular to the fibers, Da,slice, is found to be almost zero, which supports the assumption that the diffusion in the IAS is almost fully restricted.
The axial EAS diffusion coefficient in the corpus callosum is about 2.5 μm2/ms, which is 17% lower than the free diffusion constant of water at 37°C (3 μm2/ms), indicating that the water molecules are only weakly hindered by membranes in the EAS in the direction along the fibers. The radial EAS diffusivity, De,, is found to be much smaller than De,, which supports the EAS diffusion in the direction perpendicular to the fibers being strongly restricted.
The EAS tortuosity provides an indirect measure of the myelinated axonal fraction (including the myelin), as it increases for decreasing EAS volume fraction (Fieremans et al., 2008). The higher tortuosity values observed in the corpus callosum (Fig. 2) can be explained by the presence of myelin and a high axonal density. The mean value α = 2.79 we found in this study for WM is close to literature values. As a comparison, a mean tortuosity value α = 3.1 has been measured in brain extracellular space (Kroenke and Neil, 2004) using diffusion MRI. Similarly, a mean tortuosity value α = 2.89 has been found for the corpus callosum of adult rat brain using the standard iontophoresis method with tetramethylammonium (TMA) as probing molecule (Syková and Nicholson, 2008) (this value is derived from the reported value of 1.7 for the tortuosity, operationally defined as , and squared according to our definition Eq. (15)).
The AWF provides a measure of the IAS water volume relative to the EAS water volume and neglects the myelin water that practically does not contribute to the DWI signal due to fast transverse relaxation. This parameter should not be confused with the total cellular volume fraction of approximately 0.8, as the EAS also contains glial cells. The highest AWF values are observed in the corpus callosum, as shown in Fig. 2 and Fig. 6, as expected. The mean value f = 0.49 we found in this study for WM is in good agreement with normalized water fractions as estimated using multiexponential analysis of T1- and T2-relaxation measurements in human WM (Hwang et al., 2010; Lancaster et al., 2003; Lancaster et al., 2005) and rat trigeminal nerve (Does and Gore, 2002). Our results confirm that in the WM there is no real disagreement between the fast and slow diffusion fractions from the biexponential fit and the physical water fractions of the EAS and IAS (Assaf and Basser, 2005; Maier et al., 2004).
When neglecting the contribution of the myelin volume, the measured AWF can be compared to histological volume fractions: Tang et al. (1997) performed stereology on human WM and found a mean value of 0.33 for the volume fraction of myelinated fibers. This value agrees well with the earlier reported mean AWF-value of 0.3 over all WM voxels with FA > 0.25 (Fieremans et al., 2010a). The mean AWF-value of 0.49 we report here is higher, probably because we limited this study to those WM-voxels (according to Eq. (16)) that have the highest fiber density. For example, the AWF-values we find in the corpus callosum (Fig. 6) compare favorably to the reported value of about 0.7 for the axonal volume density of the commisura anterior of rat brain based on AMG staining (Jespersen et al., 2010).
The values we find for the AWF using DKI data are also similar to the values found in other advanced multi-compartment models that typically require analysis of high b-value diffusion data (Alexander et al., 2010; Jespersen et al., 2010; Jespersen et al., 2007; Panagiotaki et al., 2009). We compared the DKI WM model here against the biexponential model using high b-value data. Assuming that the fraction of the slow diffusion component corresponds to the axonal compartment, we find a reasonable agreement between this fraction and the AWF. Although the agreement seems slightly better for the AWF based on an ROI estimate of Da (Eq. (8)) in Fig. 6 (d)), the voxel-based approach (Eq. (7)) in Fig. 6 (c)) is more objective and easier to implement, while also providing a fairly good estimate. These preliminary results suggest that approximately the same information about compartment diffusivities and the AWF can be derived from DKI-data as from high b-value diffusion data using biexponential or more complicated fitting. The advantage of our approach is that the DKI data-acquisition time is much shorter because the DKI WM-analysis requires only low b-values (up to 2000 s/mm2) and puts less demand on the hardware. Thanks to the relatively low b-value range, DKI also provides an adequate signal-to-noise (SNR) making the technique useful in a clinical setting.
Another potential advantage of the DKI WM model over other advanced high b diffusion models lies in the fact that because of the smaller b-range, the two-compartment model has a broader applicability in the case of voxels containing non-collinear axons. For the data analysis and ROI selection, we applied a mask that selected voxels with a homogeneous axon orientation (Eq. (16)). Whereas this mask is most likely a prerequisite for applying biexponential fitting using high b-value diffusion data, the derivation in section 2.1.1 (together with Fig. 1 and Appendix A) suggests that the DKI WM model might be valid in significantly more voxels than those selected by the WM mask. Future work will investigate in more detail the applicability of the DKI WM model in less highly oriented voxels.
This model offers new perspectives for the clinical assessment of WM tissue. DKI is relatively easy to implement on clinical scanners and a typical DKI protocol requires 7 to 20 minutes for full brain coverage depending on the hardware and imaging parameters (Jensen and Helpern, 2010; Tabesh et al., 2010). The current DKI WM model then allows for the estimation of various microstructural parameters from the DKI data, such as the compartment diffusivities and the water fractions of each compartment. The axonal water fraction is formally determined by the maximum kurtosis over all directions (Eq. (7)), but can in practice often be approximated by the radial kurtosis. The axial kurtosis reflects then the diffusional heterogeneity between the IAS and the EAS.
These parameters may help elucidate the meaning of changes in DKI metrics as observed in clinical studies. As an example, Jensen et al. (2010) reported substantial increases in the mean kurtosis within ischemic lesions for stroke patients. More specifically, for the lesions with strongly oriented axon bundles, the increase in the axial kurtosis was much greater than that in the radial kurtosis. Using our approach, it follows then naturally that this change is due to a large decrease of the intra-axonal diffusivity, as also suggested in their study. The AWF and the tortuosity are likely to be sensitive to changes in the number of axons, the myelin volume and axonal geometry and could potentially provide important information for assessing, e.g. multiple sclerosis (Inglese and Bester, 2010; Warlop et al., 2009; Warlop et al., 2008), Alzheimer's disease (Bartzokis), and other neuro-psychiatric and/or neuro-degenerative disorders that may be related to myelin dysfunction (Nave, 2010). More advanced analysis methods such as tract-based spatial statistics (TBSS) (Smith et al., 2006) and tractometry (Bells et al., 2011) seem particularly well suited for our model, as these naturally restrict the analysis to WM voxels consisting of strongly aligned fibers. Future work will focus on the clinical applications of the DKI WM model.
The present WM model is highly idealized, with some aspects of water diffusion in WM not explicitly included. One such limitation is that the CSF is not considered as a separate compartment. A possible effect of different transverse relaxation times for the intra- and extra-axonal compartments is also neglected, that may influence the estimates of the AWF and the tortuosity (Frøhlich et al., 2008). Furthermore, the DKI WM model is based on relatively low b-values diffusion data in a small number of diffusion directions, which inherently limits the amount of information that can be extracted. This is in contrast to other advanced diffusion models, based on high b-value and/or many gradient directions, aiming to characterize the axon diameter, diameter distribution or fiber orientation distribution (Alexander et al., 2010; Assaf et al., 2008; Barazany et al., 2009; Jespersen et al., 2010). The limited information available in a DKI-dataset, due to its low b-value and limited angular information, makes it impossible to incorporate all of these as well as some other features into our model. Nonetheless, our proposed DKI WM model seems to be a useful description of the actual biophysical structure that yields realistic values for all microscopic parameter values included in the model. Important advantages of our approach are both its simplicity and the practical benefits of only requiring relatively low b-value diffusion data.
We have proposed an idealized two-compartment no exchange diffusion model of white matter suitable for analysis with diffusional kurtosis imaging (DKI) diffusion metrics. Based on this model, standard DKI metrics can be used to estimate the intra- and extra-axonal diffusivities, the axonal water fraction, and the tortuosity of the extra-axonal geometry. Values for these parameters obtained in healthy young adults agree well with those of prior studies. The axonal water fraction and the tortuosity provide information related to axonal and myelin density, which may be useful in assessing myelin-associated neuropathologies. Since a DKI dataset can be acquired within a few minutes, this approach can be applied for the quantitative assessment of white matter tissue properties in a clinical setting.
We thank Caixa Hu, Dmitry Novikov and Ali Tabesh for technical assistance and discussions. This work was supported by National Institutes of Health research grants (1R01AG027852 and 1R01EB007656 (to J.A.H)), the Litwin Foundation (to J.A.H.) and the King Baudouin Foundation (Henri Benedictus Fund (to E.F.)).
The bundles of axons are regarded as Gaussian compartments in our DKI WM model. We propose criterion (5) to test the validity of this assumption. Here we apply this criterion for the two axonal compartment geometries discussed in section 2.1.1. For that, we derive expressions for the compartmental D and K in these systems.
The first example discussed in section 2.1.2 is a voxel consisting of two crossing fiber bundles, oriented at polar angles (θA, θB) relative to a particular diffusion direction of interest, as illustrated in Fig. 1 (a). When idealizing the axons as having zero radius, the total axonal diffusion coefficient in such a system is given by
where we have assumed for simplicity the two bundles have the same water fractions and Da is the free intra-axonal diffusivity. The kurtosis for a given direction is derived using the standard formula for the multiple compartment kurtosis (Jensen and Helpern, 2010):
Next, we consider the general case of axons that are oriented according to a direction distribution of F(n), where n is a unit vector parallel to a particular axon. The distribution of directions is normalized so that
where dΩn is a solid angle element and the integral is taken over all possible directions. The diffusion coefficient in a direction n is then given by
and, similarly, the diffusional kurtosis is given by
Now let us consider a spherical coordinate system with n corresponding to the z-axis. Then we have
with (θ, ϕ) being the usual spherical angles. Note that G(θ) is normalized so that
Now we consider the 2nd example of axons in a nearly coplanar configuration with the diffusion gradient applied perpendicular to the axons. Let us assume that the axons are uniformly distributed within a range of angles defined by
independent of θ0, which using criterion (5), leads to the b-value condition
So even for a π/2-range (θ0 = π/4), taking Da = 1 μm2, we find b ≤ 2500 s/mm2, which is adequate for DKI.
Here we outline the derivations of some of the basic results presented in the section 2.1.3. To develop the DKI WM model, we choose an arbitrary reference frame where D1, D2, D3 are the measured diffusion coefficients along the axis directions 1, 2 and 3 and K1, K2, K3 the corresponding kurtoses. The diagonal elements of the compartmental diffusion tensors Da and De in Eq. (1) in this basis are Da,1, Da,2, Da,3 and De,1, De,2, De,3 respectively. By modeling the system as two non-exchanging compartments, the overall diffusion coefficients, D1, D2, D3, can be described as a function of the model parameters:
for i = 1, 2 and 3, where f is the volume fraction of the IAS.
Similarly, the along axis kurtoses, K1, K2 and K3, are related to the model parameters by
for i = 1, 2 and 3.
for i = 1, 2 and 3.
For the small axon approximation, one can show that
where Da is the free longitudinal intra-axonal diffusion coefficient along the axis. In addition, the mean extra-axonal diffusivity is
We first resolve the sign ambiguities in Eqs. (B.3) and (B.4) in the reference frame in which the overall D is diagonal: Taking De,3 ≥ Da,3 makes physical sense since the 3rd direction (which by assumption has the smallest diffusion coefficient) should be roughly perpendicular to the axon direction. Thus, we take the upper sign for the 3rd direction. The signs for the 1st and 2nd directions are then determined by
which is, theoretically, 1 for the upper sign and –1 for the lower sign. In Eq. (B.7) Wijkl is the kurtosis tensor, here defined in the frame of reference where the overall D is diagonal. The derivation of Eq. (B.7) follows from the properties of the kurtosis tensor as described by (Jensen et al., 2005; Lu et al., 2006) but is omitted due to length considerations. Empirically ηi is found to be about 1 for most voxels in healthy young brain, as illustrated in the example in Fig. B.1 (a) and (b). Hence, De,i ≥ Da,i. in the eigenframe in which the overall D is diagonal.
where D̄ (D1 + D2 + D3)/3 the total mean diffusivity. Applying Eq. (B.2) leads to
Combining with Eq. (B.8) yields the following expression:
Note that the right side of Eq. (B.10) is independent of the choice of reference frame. The same must then be true for the left side. Now let us assume that we start in the eigenframe of the overall D for which ±i = + for i = 1, 2, 3 and imagine continuously rotating the reference frame so that the values for (D1, D2, D3) and (K1, K2, K3) also change continuously. For the left side, continuity implies that the sign ambiguity ±i can only flip for a reference frame in which Di and/or Ki vanish. For brain tissue, we know that Di does not vanish, as all the diffusion tensor eigenvalues are positive. Thus, we can only have a sign flip if there is a diffusion direction for which Ki vanishes. As a corollary, we see that if Kmin> 0, where Kmin is the kurtosis minimum for all directions, then if the upper signs apply in one frame of reference then they must also apply in all reference frames. Empirically, Kmin is found to be positive (> 0) for most voxels in healthy young brain white matter, as illustrated in the example in Fig. B.1 (c). Therefore, we take the upper in sign in Eqs. (B.3) and (B.4) henceforth, so that De,i ≥ Da,i in any frame of reference. With this preferred sign choice, Eqs. (B.3) and (B.4) can be rewritten as Eqs. (10) and (11) in section 18.104.22.168.
In order to derive expressions for f, we express first all the undetermined model parameters in Eqs. (B.1) and (B.2) via Da, the free longitudinal intra-axonal diffusion coefficient along the axis, and then try to estimate Da. According to our sign choice, solving Eq. (B.10) for f gives
Thus the AWF can be found from the DKI measurements plus an estimate for Da. What is left now is to estimate Da. The physical requirement Da,i ≥ 0, combined with Eq. (B.4), yields the condition
Applying this to Eq. (B.10) leads to
We then have the lower bound of
So the largest kurtosis value Ki = Kmax in the denominator gives the best lower bound and hence the best estimate for Da.
Now let us define
where we vary the right hand side over all voxels within a specified ROI. We may then reasonably make the approximation . This would be exact if the axons are perfectly aligned in any of the voxels.
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.