|Home | About | Journals | Submit | Contact Us | Français|
Due to its unique sensitivity to tissue microstructure, diffusion-weighted magnetic resonance imaging (MRI) has found many applications in clinical and fundamental science. With few exceptions, a more precise correspondence between physiological or biophysical properties and the obtained diffusion parameters remain uncertain due to lack of specificity. In this work, we address this problem by comparing diffusion parameters of a recently introduced model for water diffusion in brain matter to light microscopy and quantitative electron microscopy. Specifically, we compare diffusion model predictions of neurite density in rats to optical myelin staining intensity and stereological estimation of neurite volume fraction using electron microscopy. We find that the diffusion model describes data better and that its parameters show stronger correlation with optical and electron microscopy, and thus reflect myelinated neurite density better than the more frequently used diffusion tensor imaging (DTI) and cumulant expansion methods. Furthermore, the estimated neurite orientations capture dendritic architecture more faithfully than DTI diffusion ellipsoids.
The utility of diffusion weighted imaging, in particular in acute stroke (Hjort et al., 2005; Moseley et al., 1990), continues to drive research in the underlying biophysical mechanisms. During a typical diffusion experiment, water molecules in living tissue probe length scales on the order of 5 to 20 µm, making diffusion sensitive to a plethora of microstructural and physiological detail. The exquisite sensitivity of diffusion weighted MRI to tissue microstructure fuels the hope that the MRI scanner may eventually serve as an in vivo microscope of the brain and body. This would allow scientists and clinicians to obtain information not only about fiber tract trajectories (Ciccarelli et al., 2008), but also more detailed traits of tissue morphology such as, e.g., cell and fiber density, compartment sizes, and various cytoarchitectural attributes. Many studies have shown that parameters routinely obtained from diffusion weighted MRI, e.g. fractional anisotropy and mean diffusivity, act as biomarkers for different tissue properties (Basser, 1995; Basser, 1997; Gupta et al., 2008; Kozlowski et al., 2008; Tofts, 2003). Such MRI-derived physiological stains are very useful, but often suffer the fundamental disadvantage of lacking specificity—thus, more direct measures of histological characteristics would be advantageous. In order to achieve this goal, careful biophysical modeling of diffusion and independent validation using existing histological methods is a prerequisite. An early example of this approach demonstrated that membrane permeability and cellular surface to volume fraction could be obtained from modeling diffusion in packed erythrocytes (Latour et al., 1994). More recently, Yablonskiy et al. modeled diffusion in the airways of the lung, and showed that quantitative geometric information about the alveoli could be extracted from diffusion MRI (Sukstanskii and Yablonskiy, 2008; Woods et al., 2006; Yablonskiy et al., 2002). Two groups have recently obtained good agreement to histology in measuring axon diameters with q-space imaging, (Bar-Shir and Cohen, 2008; Ong et al., 2008). However, the general application of q-space methods for inferring microstructure requires careful interpretation, for example in the presence of permeable barriers (Sukstanskii et al., 2004). Assaf et al. (2008) accomplished good agreement between histological measurements of axonal diameter distributions and diffusion MRI by utilizing their modeling framework “AxCaliber.”
Recently, we introduced a biophysical model of diffusion in the brain (Jespersen et al., 2007). The model was based on a biophysical description of brain microstructure, and was shown to describe high b-factor diffusion weighted MRI experiments of a baboon brain better than alternative models. The fundamental assumption of the model is that diffusion can be described in terms of two non-exchanging components. One component is associated with diffusion in cylindrically symmetric structures, such as cell processes with exchange of water being sufficiently slow to be considered impermeable on the time scale of the diffusion experiment. Dendrites and axons, collectively termed neurites, were assumed to fulfill these assumptions. The net signal from this component then arises as a sum of the signal from all neurites weighted by an orientation distribution function, i.e. a probability density function specifying the number of neurites in every direction. The second component of the diffusion signal accounts for diffusion everywhere else, in particular in cell bodies, extracellular space, and glia cells. Here diffusion is hindered and described in terms of Gaussian isotropic diffusion with an effective diffusion constant. Several cytoarchitectural parameters can be extracted from this framework, one of the more interesting being the neurite density. Noninvasive, and eventually in vivo access to this parameter could be useful in several research areas, e.g. brain mapping (Schleicher et al., 1999), postnatal ontogeny (Amunts et al., 1997), aging (Stark and Pakkenberg, 2004), comparative neuroanatomy (Sherwood et al., 2004; Zilles et al., 1986), multiple sclerosis (Peterson et al., 2001; Trapp et al., 1998), schizophrenia (Colon, 1972), Alzheimer's disease (Stark et al., 2005), and alcoholism (Tang et al., 2004).
The purpose of the present work was to compare model parameters to measurements obtained from light microscopy and quantitative electron microscopy (stereology) in order to evaluate the interpretation of the model parameters. In particular, MR diffusion measurements of neurite density obtained by applying the diffusion model were compared to measures of neurite density and cell density at the light and electron microscopic level. The results were contrasted to an identical analysis using parameters from the DTI and cumulant expansions. Furthermore, we assessed the accuracy of the estimated neurite orientation functions qualitatively.
We analyzed the data in terms of the model presented previously (Jespersen et al., 2007), as well as an extension to account more directly for anisotropic diffusion in the extra-cellular space, and compared to two other diffusion models in common use, the diffusion tensor model and a fourth order cumulant expansion. In the model introduced in (Jespersen et al., 2007) (Model I), the signal S consisted of two non-exchanging components, one with cylindrical symmetry describing diffusion in neurites (dendrites and axons), and a hindered component accounting for diffusion in the remaining compartments, among which fast exchange was assumed:
Here, q = γδg is the diffusion wave vector, γ the gyromagnetic ratio, δ the duration of the diffusion gradients, g the diffusion gradient, Δ the time between the leading edges of the two diffusion gradients, S0 the signal amplitude at q = 0, Sc the signal component arising from molecules in a collection of cylindrically symmetric structures, v the associated water volume fraction, and Sh the signal from the hindered component modelling diffusion in the extra-cylindrical space. In (Jespersen et al., 2007), the model was aiming in particular at describing diffusion in gray matter, and the extra-cylindrical component Sh was assumed to be well approximated by isotropic Gaussian diffusion Sh = exp(-bDeff), where b = q2(Δ – δ/3) and Deff is the effective diffusion constant. Here, we extend the model to explicitly incorporate anisotropy in the extra-cylindrical compartment (Model II) by generalizing the Gaussian diffusion constant of the isotropic compartment to a diffusion tensor. Thus, for model II we use Eq. (1.1) with
where summation over repeated indices is implied (Einstein summation convention), and D is the diffusion tensor in the extra-cylindrical space. When discussing the results for the model parameters, we will refer to the eigenvalues of diffusion tensors as D1, D2, and D3 with D1≥D2≥D3. Models I and II both make use of an expansion of the orientation distribution function of neurites in spherical harmonics (Laplace expansion). Employing the Gaussian approximation for diffusion in each neurite, the signal contribution Sc can be viewed as the integral of a cylindrically symmetric diffusion tensor model over a distribution f(θ,) of principal directions. Here we terminate the expansion of f at L=4 (Jespersen et al., 2007) such that
where flm are the expansion coefficients of the orientation distribution function to be determined from the data, DL and DT are the longitudinal (parallel to the neurite) and transverse (perpendicular to the neurite) diffusion constants for water inside the neurites, and
In addition, a standard DTI fit and a fourth-order cumulant tensor expansion (CumExp) were also obtained for comparison. The fourth-order cumulant expansion (Risken, 1984), also known as the kurtosis model (Jensen et al., 2005; Lu et al., 2006), is described here by (see also (Liu et al., 2004))
W is a rank 4 fully symmetric tensor with 15 independent parameters, fulfilling
in the narrow pulse approximation (Callaghan, 1991), where R is the spin displacement and i,j,k,l=x,y,z label its Cartesian components. We examine, in particular, the directionally averaged W
as well as mean kurtosis (Jensen et al., 2005; Lu et al., 2006). It is important to note that the standard DTI model is also a cumulant expansion terminated at second instead of fourth order. The cumulant expansion is a mathematical representation of the Fourier transform of a probability density function (the characteristic function, corresponding here to the MR signal) valid under quite general circumstances (Risken, 1984). As such, it does not take into account specific tissue microstructure, and despite showing great potential as a clinical biomarker (Cheung et al., 2009; Falangola et al., 2008; Latt et al., 2008), its parameters can be difficult to interpret in terms of biological properties, although attempts are being made (Frohlich et al., 2006; Frohlich et al., 2008). This is in contrast to models I and II, which are biophysical models explicitly derived from detailed microstructural descriptions of brain tissue and formulated in terms of biophysical tissue parameters. Of course, the cumulant expansion is not expected to adequately describe the diffusion attenuated MR signal in the whole range b-values. Most likely, it is sufficient when bD is on the order of or less than unity, but a detailed discussion of this question is beyond the scope of this work. Suffice it to say that with the diffusion weighting employed here, we are probing a regime covering up to bD≈5 as judged by the maximum b-value and the trace of the diffusion tensor.
Three male Wistar rats (270–310 g) were exsanguinated during intraaortic perfusion fixation with isotonic saline containing heparin (10 IU/mL), followed by 4% paraformaldehyde in phosphate-buffered saline (PBS) (pH 7.4). The rat brains were removed and immersion-fixed in a fresh 4% paraformaldehyde solution at room temperature for 2 weeks. Each brain was divided into two hemispheres, one for stereology and the other one for histology and MRI.
Prior to MRI, the brains were washed for 2 days in a phosphate-buffered saline solution (PBS) pH=7.4 in order to remove formalin and minimize concomitant T2 degradation. Each specimen was then placed in a 10-mm standard NMR glass tube, the axial direction of the cerebrum parallel to the tube axis, and positioned in a 16.4 T (700 MHz for 1H) vertical-wide-bore Bruker Avance II spectrometer (Bruker BioSpin GmbH, Rheinstetten, Germany), equipped with a gradient system capable of up to 300 Gauss/cm. All experiments were performed at 21 °C A standard spin echo Stejskal–Tanner diffusion weighted sequence was used to acquire a total of 144 unique diffusion directions chosen from a 144 point spherical 16-design (Hardin and Sloane, 1996). A total of 16 linearly spaced shells were used from b=0 to b=15000 s/mm2, i.e. 9 directions on each shell (but note that the 9 directions change between shells). The remaining diffusion and imaging parameters were as follows: TR=3 s, TE=14.7 ms, data matrix 128×128, field of view 12.8 mm×12.8 mm, slice thickness 0.5 mm, and Δ/δ=8/2 ms. Each data set therefore represents more than 15 h of scan time.
The 144 diffusion weighted images were fit to each of the four models on a pixelwise fashion using the nonlinear least squares Levenberg–Marquardt algorithm as implemented in Matlab (The MathWorks, Inc., Natick, MA). We analyzed the reproducibility of the chi-square (χ2) minimum by plotting the mean and standard deviation of the parameters of interest obtained from fits over 50 pixels uniformly distributed over one of the rat brains using 50 different randomly chosen starting values for each pixel. Here , where Yi is ith datapoint corresponding to qi, and σ is the standard deviation of the noise in the real channel estimated from the signal average in an empty region of the image. If there were a substantial redundancy of parameters, we would expect a large number of local chi-square minima, resulting in turn in a high degree of uncertainty in the parameter estimations. For Model I, variation in the estimate of the neurite density was quite substantial, but this was fully alleviated by constraining the transverse diffusivity DT to 0, a reasonable approximation when the diffusion length, which is on the order of 5 µm in our experiments, is large compared to the diameter of the neurite. For myelinated rat axons with a typical outer diameter r of 0.3 µm (found with stereology, see Methods and also (Yang et al., 2008)), the effective transverse diffusivity is bounded from above by approximately DT~r2/4Δ≈0.003 µm2/ms, and this is expected to be significantly smaller than DL. Furthermore, bmaxD-<0.05, and for these reasons, setting DT = 0 seems to be a reasonable way to avoid the instability in the fit.
Due to the constraints that the orientation distribution function is normalized (fixing f00), real (flm=(−1)mfl−m) and has antipodal symmetry (flm=0 for odd l), the number of free parameters in the spherical harmonics expansion of order L (L even) is 5+9+…+ (2L +1) = L(L + 3) / 2. Thus the resulting number of parameters in the models are 18 (model I), 23 (model II), 7 (the DTI model), and 22 (CumExp).
Following MRI, the brains were re-immersed in formalin for weeks before cryo-protection with 30% sucrose, CO2-freezing and cryostat cutting into 40 µm coronal sections. They were then exposed to two types of histochemical stains, a standard Nissl-stain for neuronal cell bodies and an autometallographic (AMG) myelin stain (Larsen et al., 2003). The optical staining intensity of AMG reflects the volume fraction of myelin, and is therefore expected to be correlated with the volume fraction of axons, whereas the optical staining intensity of Nissl stain reflects the volume fraction of cell bodies, and is expected to be inversely correlated with the volume fraction of axons. The sections were digitally photographed under constant lighting conditions with an approximate resolution of 5.6 µm per pixel. Before comparison to the MRI images, the histological images were smoothed to yield an MRI comparable resolution of approximately 100 µm isotropic. Approximately corresponding slices from histology and MRI were identified from unique anatomical landmarks.
For use in comparing MRI model parameters with optical staining intensities, we defined a number of intensity-based ROIs: the v parameter maps were segmented linearly into five intensity levels, and one or two ROIs containing 24 pixels were chosen in each as shown in Fig. 1A. For later comparison, we also outlined 8 anatomical ROIs, defined in Fig. 1B: hindlimb cortex, parietal cortex, hilus hippocampus, stratum moleculare fasciae dentatae, stratum radiatum hippocampus, stratum moleculare hippocampus, commisura anterior, and striatum. These regions were selected to reflect a wide spectrum of dendrite and myelinated axon densities, and used to independently evaluate the relations obtained from the intensity based ROIs.
The other rat brain hemisphere was cut into ~2 mm coronal slices and corresponding slices from MRI were identified using the above-mentioned unique anatomical landmarks. From each anatomically defined region (striatum, commisura anterior, parietal cortex, hindlimb cortex, hilus fasciae dentatae, and stratum radiatum hippocampus) ~2×2×2 mm cubes were selected in the X-, Y-, and Z-direction. These cubes were embedded in TAAB 812 Epon (TAAB, Berkshire, England) with the EMS LYNX Automated Tissue Processor. Ultrathin sections (60–90 nm) were cut on a RMC MT6000 ultramicrotome from each regional sample in three directions and collected on pioloform-coated slot grids. Sections were observed in a Philips CM 10 electron microscope. Electron micrographs were taken with a digital camera (Kodak megaplus 1.6i) at an initial magnification of 13,500× and digitally enlarged to a final magnification of 32,750×. Counting fields were positioned in the selected region by a systematic, uniform random technique (Nyengaard, 1999). The micrographs were saved and later analyzed using the iTEM software (Olympus Soft Imaging Solutions GmbH).
The number of axonal profiles per area was estimated by the 2D counting frame. The diameter of the sampled axons was measured as the biggest diameter perpendicular to the longest axis. The volume fraction of axons, dendrites, neuronal somas, extracellular volume, and vessels were estimated by point counting.
In Fig. 2, images of Akaike's information criterion (AIC) are shown for models I, II, DTI, and CumExp. AIC provides a means of comparing models with differing numbers of free parameters by penalizing chi-square by two times the number p of free parameters in the model: AIC=χ2 + 2p (Akaike, 1974). From the AIC images, we conclude that DTI accounts only poorly for the variation of the diffusion signal with high diffusion weighting, as expected when bD1. The most faithful description of the diffusion signal by the DTI model occurs in areas in which almost no myelinated axons are present, e.g. the hippocampal stratum radiatum and stratum moleculare, and the hindlimb and parietal cortices, cf. Fig. 1. This trend is shared by all models, and is not unexpected given the very high b-values probed here. The cumulant expansion goodness of fit is of slightly poorer quality than models I and II. These models (I and II), generally attain the lowest values of the AIC, indicating that they describe data better than the other models. Nevertheless, the corresponding AIC images are clearly not pure noise, demonstrating that not all information in the diffusion signal has been captured by the models. Model II has a lower AIC in white matter as compared to model I as intended, but this comes at the expense of a significantly larger number of parameters (23 versus 18), and a correspondingly higher computational demand. A common trend of all the models is a relatively poor performance in the fimbria and the optical tract. This is particularly conspicuous in model I, whereas model II and the cumulant expansion appear to perform slightly better. We believe the reason has to do with the fact model I was designed for gray matter and is therefore not optimal in highly anisotropic white matter fiber bundles. In fact, the fimbria and the optical tract are highly uni-directional fiber bundles: even the corpus callosum is presumably less anisotropic due to splay of individual fibers in the plane perpendicular to the image, and indeed, the models appear to fare better in the corpus callosum than in the fimbria. A careful look on the FA maps of the DTI model (see below) is consistent with this interpretation: the average of FA in the fimbria is approximately 0.77 whereas the average FA of the corpus callosum is 0.70. The impression of the AIC being higher in the fimbria is further accentuated by partial volume effects in the AIC images: in the fimbria, several pixels are located in the interior of the bundle, whereas the AIC in the corpus callosum may be slightly underestimated due to contributions from tissue in the immediate vicinity, where the model shows a better fit.
A subset of the resulting parameter maps are illustrated in Fig. 3: for model I, we show in Fig. 3a the isotropic effective diffusivity Deff, and longitudinal diffusivity DL; for model II (Fig. 3b), longitudinal diffusivity DL as well as the mean diffusivity and fractional anisotropy FA, of the diffusion tensor in the extra-cylindrical space; for the DTI model (Fig. 3c) and CumExp (Fig. 3d), the ADC and FA are shown, in addition to the directionally averaged fourth order cumulant W for CumExp, as defined in Eq. (1.7).
The effective diffusion constant Deff of model I mainly reflects contrast between gray matter and white matter. In the gray matter areas such as in the cortex, typical values of Deff are around 0.4 µm2/ms, whereas typical values in the fiber tracts are lower, approximately 0.1 µm2/ms. Such a low value in the fiber tracts is presumably an artifact of the model failing to account properly for highly anisotropic extracellular diffusion—indeed, the mean diffusivity of model II (see below) has a higher value in the corpus callosum than in the surrounding gray matter. Areas with high cell body density, in particular hilus fascia dentatae and the pyramidal cell layers in the hippocampus have a higher than average Deff, and this observation applies equally to all other diffusivity constants in all four models. The longitudinal diffusivity is somewhat overestimated, and is generally at the level or above the diffusivity of pure water at the same temperature (D≈2.5 µm2/ms). This is not plausible physically, and we ascribe it to the inherent approximations of the model. The contrast in this parameter map is low, in agreement with the well-established constancy of the ADC across brain regions, across individuals, and even across species (see e.g. (Pierpaoli et al., 1996) and references therein).
The impression of the MD maps from model II is quite similar to the effective diffusivity from model I, only it is substantially higher than Deff in white matter fiber bundles such as the corpus callosum. The longitudinal diffusivity is also comparable to that obtained from model I, except it is generally lower and even markedly lower at some locations along the corpus callosum and in the cortex. The fractional anisotropy of the extra-cylindrical component shows very high degree of white/gray matter contrast as expected.
The ADC map from the standard DTI fit resembles the MD map from model II, but even more the effective diffusivity Deff of model I. The FA map is similar to model II, although with generally lower values and its variation appears more smooth, e.g., in the cortex and in the thalamus.
The ADC map of the kurtosis model is strikingly flat, as the other diffusivities described so far, and deviates substantially from its average value of around 0.4 µm2/ms only in small areas in the thalamus and in the vicinity of white matter fiber tracts where it is lower, and in hippocampal cell layers where is higher. The qualitative behavior of the contrast variation is similar to Deff. The corresponding FA map is only barely distinguishable visually from the FA map of the DTI model. The orientationally averaged kurtosis W varies between approximately 0.6 and 1.2 µm4, and is large in particular in areas of incoherently oriented myelinated axons, such as the parietal cortex and thalamus.
In Fig. 4, all diffusion tensor eigenvalues are shown together for comparison: on the top row are eigenvalues of the extracellular diffusion tensor of model II, the middle row shows eigenvalues from the DTI tensor, and on the bottom row are eigenvalues from the diffusion tensor part of the fourth order cumulant expansion. The three models have quite similar D3 maps, with a few differences in the thalamus, for example. The differences between the D2 maps are slightly more obvious, whereas the D1 maps clearly differ in both white and gray matter. The images of the two minor diffusion eigenvalues of model II qualitatively resemble Deff more than MD, suggesting that Deff in model I is biased towards transverse extracellular diffusion. This may, in turn, be related to the observed overestimation of DL particularly pronounced in model I. The major diffusion eigenvalues in white matter are generally higher for model II than for the other two models, explaining the differences in the ADC and MD maps noted earlier. That is, D1 of model II in the corpus callosum is sufficiently high to render the mean diffusivity there higher than average, which is not the case for the other two models. Keeping in mind the fact that DTI and CumExp diffusivities reflect properties of intra- and extracellular diffusion, these differences indicate that intracellular parallel diffusion is slower than extracellular parallel diffusion in fixated rat white matter tissue.
In Fig. 5, we demonstrate the qualitative correspondence between the slices from light microscopy and the MRI parameter maps. Here, the neurite density v calculated from models I and II is shown alongside the AMG histology image. The neurite density from models I and II both have values in the entire range between 0 and 1, and generally display similar features, although the gray/white matter contrast of model I is higher than that of model II. Neurite density is generally very high in the major white matter fiber tracts such as the corpus callosum, and lower in gray matter areas such as the cortex. In fact, it seems to very accurately follow the contrast in the AMG staining, reflecting the density of myelin. This holds true not only in the major white matter fiber tracts as mentioned above, but also for the myelinated axons in e.g. central parts of the parietal cortex. Conversely, in anatomical regions containing mainly non-myelinated axons and dendrites such as the stratum radiatum of the hippocampus, v is very small. These observations will be substantiated shortly with the quantitative analysis. Therefore, with the diffusion parameter settings in the current experiment, the neurite density of models I and II is seen to be mainly sensitive to myelinated axon density. However, it is important to emphasize that even though the dendritic contribution to the neurite density appears to be underestimated, it is still being captured by the models, as evident in the fiber orientation distribution maps to be discussed shortly.
A table of Pearson's linear correlation coefficients between parameter and light microscopy values from Fig. 1a intensity based ROIs is found in Table 1, showing that several model parameters display some degree of correlation to the optical staining intensity of the AMG. We interpret this to reflect a general sensitivity of diffusion weighted MRI to myeloarchitecture, being as it is a significant barrier to diffusing water molecules. Nevertheless, the neurite density of model I shows the strongest correlation, as seen in the table. Model II has a somewhat lower correlation coefficient. This could be due to chance fluctuations, but may also be due to the fact that the entire anisotropic contribution must be numerically accounted for by the neurite component in model I, whereas it may be partly shared between the two components, extracellular DTI and cylindrical neurites, in model II. Constraining model II to having a cylindrically symmetric diffusion tensor in the extra cylindrical space increases this correlation coefficient to 0.85. In the other two models (DTI and CumExp) as well as model II, the smallest diffusion tensor eigenvalue D3 is the model parameter with the highest correlation. In the case of model II, this correlation is natural since diffusing water molecules in extracellular space follow more tortuous paths perpendicular to closely packed fiber bundles in white matter than in gray matter. DTI is a phenomenological approach, and it is not surprising that it has parameters that reflect microstructure in some way—however, the diffusion tensor in the DTI model reflects properties of intracellular as well as extracellular diffusion, and it is difficult to predict a priori which parameters correlate with a given biophysical trait. This is in contrast to models I and II, which were developed by taking known brain microstructure into account. Note also that the fitting quality was substantially poorer for DTI than models I and II (we are well into the regime of non-monoexponential signal decay), and it is not clear how the residual information in the signal from high b-values affects the estimate of e.g. D3. We note also that DTI parameters obtained from a fit constrained to lower b-values such that bD<1 show substantially reduced correlation with staining. In model II, DTI and CumExp, FA is highly negatively correlated to the Nissl staining intensity also, presumably a simple consequence of less cell bodies appearing in the white matter where anisotropy is higher.
Scatter plots of the main model parameters versus AMG histology intensity values are shown in Fig. 6 along with the best linear fit and 95% confidence lines. The statements made above are also appreciated from the 95% confidence lines in Fig. 6, illustrating that the prediction of staining intensity from model I neurite density (Fig. 6a) is likely to be more accurate than any of the parameters with the next highest correlations, model II v and D3 of DTI and CumExp, since more points fall near to or within the 95% confidence region. D3 from DTI and CumExp is a weighted contribution of intra- and extracellular diffusion, and can therefore be large not only in regions with few myelinated axons, but also in densely packed fiber bundles of large diameter axons, possibly explaining the relatively large scatter in the data for large D3. The DTI and CumExp fractional anisotropy and histology deviate from a linear relationship in particular in areas of low FA (data not shown), which may correspond to a relatively wide range of AMG staining intensities. This may happen because low fractional anisotropy is found not only in gray matter where myelinated axon density is also low, but also in white matter regions of e.g. crossing fibers where the myelinated axon density is high.
In Fig. 7, the obtained linear fits relating model parameters to the AMG staining intensity are used to independently evaluate the predictive ability in all three of the hemispheres in the anatomically defined regions shown in Fig. 1B. The best approximations to the AMG optical staining intensity are afforded by the neurite density v of model I, indicating a more direct correspondence between the two measures. Model I v offers the closest match for half of the ROIs, whereas the best estimates in the remaining ROIs are divided among the parameters from the other three models. Nevertheless, there is a general overestimation of the staining intensity by the neurite density v of model I, the only exceptions being the parietal cortices of two animals and the anterior commissure of one animal.
These results indicate that with the experimental settings used here, our measure of neurite density v mainly reflects the density of myelinated axons, whereas the contribution from dendrites, e.g., in the stratum radiatum of the hippocampus is underestimated. The two components in the models are mainly distinguished by the underlying symmetry, and as such, v is essentially a “cylinder density,” proportional to the number of water molecules experiencing diffusion in an approximately cylindrical environment (i.e., can be incorporated in Sc) during the time course of the experiment. All processes, axons, dendrites, and glial processes, can presumably be locally approximated by cylinders, but this geometry restricts diffusion only if the permeability of the corresponding membranes is sufficiently low to prevent most water molecules from escaping. The relevant timescale is the diffusion time, and the reason that myelinated axons dominate v could be that the exchange time of typical dendrites is on the order of or less than the 10 ms employed here. Conversely, the neurite volume fraction of model I may become overestimated in white matter fiber tracts, since diffusing water molecules in the extracellular space arguably experience cylindrical symmetry. Other properties affect the quality of the models, for example branching, bending and other deviations from cylindrical geometry on distances smaller than the diffusion length (approximately 5 µm). Finally, it should be kept in mind that v is an effective density, and is, in particular, weighted by compartmental T2 relaxation (Chin et al., 2002). Despite these complications, the cylinder density of model I presents the highest degree of correlation to myelin staining patterns, as well as providing the best fit to the MR diffusion data (together with model II).
Stereology confirms this picture: Stereological measurements from electron microscopy of myelinated axonal density, cell soma density and extracellular volume (ECV) are compared to main model parameters in Table 2. The best correlation to the myelinated axonal volume density is achieved by v of model I (and D2 of model II), but again many parameters have a significant correlation to the myelinated axonal density, although not as high as in the case of the cylinder density. In addition, fractional anisotropy of DTI and CumExp is highly correlated also to cell soma density. Interestingly, only W from CumExp (data not shown) shows a large and significant correlation to extracellular volume. Finally, we examined the consistency of our two approaches by correlating stereology to optical AMG intensity. A large and significant correlation (ρ=0.9390, p<0.002) was found, indicating that the optical staining intensity indeed reflects myelinated neurite volume fraction.
The cylinder orientation distribution functions (ODF) in two different regions (as outlined in Fig. 8) for model I and diffusion ellipsoids from the DTI model (based on the same data) are shown in Fig. 9 and Fig. 10 overlaid on an S0 (i.e. T2 weighted) image. The surfaces of the ODFs are scaled to a common maximum with the color coding reflecting the original magnitude (red=high, blue=low). The orientation distribution function for model I demonstrates that the white matter directions are captured faithfully, see for example the corpus callosum/deep white matter in Fig. 9 where the ODFs are highly oriented. More importantly, in the stratum radiatum of the hippocampus (Fig. 9) where the primary stem-dendrites radiate perpendicular to the pyramidal cell-layer, the estimated neurite orientations closely follow the anatomically expected pattern, despite the apparent underestimation of dendrite density noted earlier. The same observations hold for the supragranular layers of the parietal cortex (Fig. 10), where an equally strong radial pattern is seen, in agreement with previous reports (Kroenke et al., 2006; Neil et al., 1998). The diffusion ellipsoids reveal the same tendencies, but much less conspicuously due to the more isotropic appearance of the ellipsoids in gray matter. A DTI fit restricted to lower b-values (b≤3000 s/mm2) only does not change this observation.
In this work, we have compared output from various diffusion models to a number of different measures obtained from light microscopy and stereological methods applied to electron microscopy. Models I and II are two variations of a biophysical diffusion model introduced in (Jespersen et al., 2007), whereas DTI and kurtosis models essentially are cumulant expansions of the spin phase distribution (i.e., the MR signal) terminated at second (DTI) or fourth order (CumExp). The cumulant expansion is a generally valid mathematical construction applicable to a wide class of probability distributions (Risken, 1984), and even though it has great value as a clinical biomarker (Cheung et al., 2009; Falangola et al., 2008; Latt et al., 2008), it is difficult to relate the parameters directly to biological properties, although attempts are being made (Frohlich et al., 2006; Frohlich et al., 2008). This is in contrast to models I and II, the parameters of which are derived from a biological model of brain cytoarchitecture.
In summary, we found that neurite density of our diffusion models showed a very high degree of correlation to the optical staining intensity of the AMG based myelin stain and stereological estimates of axonal volume fraction, indicating that the model allows an accurate estimation of myelinated neurite density in both gray and white matter. In addition, the model was able to describe dendrite architecture, e.g., in the hippocampus and cortex in agreement with known anatomy. In future work, we will examine the effect of varying the diffusion time as this could be conceived to affect the constitution of the cylindrical compartment, depending on the exchange time of its various components. In the same vein, we will explore the significance of the chosen range of diffusion weighting, and incorporate a distribution of cylinder radii and transverse relaxation constants in the modeling (Hansen and Vestergaard-Poulsen, 2006). Modeling the effect of cylinder bending and branching is a focus of current efforts.
The authors wish to thank the Danish National Research foundation (CFIN and inSpin), and the Danish Council for Strategic Research, for funding. SNJ thanks Martin Carlsen and Daniel Otykier for technical assistance. The histological lab assistance of Ms. Dorete Jensen is greatly acknowledged. The authors thank two anonymous reviewers whose comments led to improvements of the paper.