Front. Neuroinform. Frontiers in Neuroinformatics Front. Neuroinform. 1662-5196 Frontiers Media S.A. 10.3389/fninf.2014.00028 Neuroscience Original Research Article Isotropic non-white matter partial volume effects in constrained spherical deconvolution Roine Timo 1 * Jeurissen Ben 1 Perrone Daniele 2 Aelterman Jan 2 Leemans Alexander 3 Philips Wilfried 2 Sijbers Jan 1 1iMinds-Vision Lab, Department of Physics, University of Antwerp Antwerp, Belgium 2Ghent University-iMinds/Image Processing and Interpretation Ghent, Belgium 3Image Sciences Institute, University Medical Center Utrecht Utrecht, Netherlands

Edited by: Daniele Marinazzo, University of Gent, Belgium

Reviewed by: Fang-Cheng Yeh, Carnegie Mellon University, USA; Shawna Farquharson, The Florey Institute of Neuroscience and Mental Health, Australia

*Correspondence: Timo Roine, iMinds-Vision Lab, Department of Physics, University of Antwerp, Universiteitsplein 1, Building N, 2610 Wilrijk, Antwerp, Belgium e-mail: timo.roine@uantwerpen.be

This article was submitted to the journal Frontiers in Neuroinformatics.

 28 03 2014 2014 8 28 06 12 2013 02 03 2014 Copyright © 2014 Roine, Jeurissen, Perrone, Aelterman, Leemans, Philips and Sijbers. 2014

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

Diffusion-weighted (DW) magnetic resonance imaging (MRI) is a non-invasive imaging method, which can be used to investigate neural tracts in the white matter (WM) of the brain. Significant partial volume effects (PVEs) are present in the DW signal due to relatively large voxel sizes. These PVEs can be caused by both non-WM tissue, such as gray matter (GM) and cerebrospinal fluid (CSF), and by multiple non-parallel WM fiber populations. High angular resolution diffusion imaging (HARDI) methods have been developed to correctly characterize complex WM fiber configurations, but to date, many of the HARDI methods do not account for non-WM PVEs. In this work, we investigated the isotropic PVEs caused by non-WM tissue in WM voxels on fiber orientations extracted with constrained spherical deconvolution (CSD). Experiments were performed on simulated and real DW-MRI data. In particular, simulations were performed to demonstrate the effects of varying the diffusion weightings, signal-to-noise ratios (SNRs), fiber configurations, and tissue fractions. Our results show that the presence of non-WM tissue signal causes a decrease in the precision of the detected fiber orientations and an increase in the detection of false peaks in CSD. We estimated 35–50% of WM voxels to be affected by non-WM PVEs. For HARDI sequences, which typically have a relatively high degree of diffusion weighting, these adverse effects are most pronounced in voxels with GM PVEs. The non-WM PVEs become severe with 50% GM volume for maximum spherical harmonics orders of 8 and below, and already with 25% GM volume for higher orders. In addition, a low diffusion weighting or SNR increases the effects. The non-WM PVEs may cause problems in connectomics, where reliable fiber tracking at the WM–GM interface is especially important. We suggest acquiring data with high diffusion-weighting 2500–3000 s/mm2, reasonable SNR (~30) and using lower SH orders in GM contaminated regions to minimize the non-WM PVEs in CSD.

 diffusion MRI fiber orientation partial volume effect constrained spherical deconvolution gray matter INTRODUCTION

Diffusion-weighted (DW) magnetic resonance imaging (MRI) is a non-invasive imaging method to investigate tissue microstructure via the measurement of the displacement of water molecules (Stejskal and Tanner, 1965; Jones, 2010). Diffusion in white matter (WM) neural tracts is anisotropic: it is larger parallel to the tract than in the perpendicular direction. In liquid, such as cerebrospinal fluid (CSF), diffusion is isotropic, i.e., equal in all directions. This diffusion property can be exploited to extract fiber orientations from DW data and investigate neural tracts in the brain WM using fiber tractography algorithms (Conturo et al., 1999; Basser et al., 2000; Mori and van Zijl, 2002; Jones, 2008; Tournier et al., 2010; Jeurissen et al., 2011).

The image resolution in DW-MRI is typically about 2–3 mm in all directions. Thus, significant partial volume effects (PVEs) are present in the measured signal (Alexander et al., 2001; Vos et al., 2011). These may be caused by multiple non-parallel neural tracts passing through a voxel (Vos et al., 2011; Jeurissen et al., 2013), or several tissue types present in a voxel (Pasternak et al., 2009; Metzler-Baddeley et al., 2012a).

Currently, the most common method in the analysis of DW-MRI data is diffusion tensor imaging (DTI; Basser et al., 1994a,b; Jones and Leemans, 2011; Tournier et al., 2011). The shortcoming of DTI is the inability to identify complex fiber configurations consisting of multiple fiber orientations (Alexander et al., 2001; Frank, 2001, 2002), present in 60–90% of WM voxels (Jeurissen et al., 2013). To overcome this, high angular resolution diffusion imaging (HARDI) methods (Tuch et al., 2002; Jansons and Alexander, 2003; Tournier et al., 2004, 2007; Tuch, 2004; Dell’Acqua et al., 2007; Descoteaux et al., 2007; Behrens et al., 2007) and methods based on diffusion spectrum imaging (DSI; Wedeen et al., 2005, 2008) have been developed. However, although able to identify complex fiber configurations, most of the HARDI methods do not account for PVEs caused by non-WM tissue, such as gray matter (GM) and CSF (Dell’Acqua et al., 2010; Metzler-Baddeley et al., 2012a).

The presence of non-WM PVEs is known in DW-MRI (Alexander et al., 2001; Pasternak et al., 2009; Dell’Acqua et al., 2010; Metzler-Baddeley et al., 2012a), but their effects in HARDI methods have not been widely studied. Diffusion in non-WM tissue is mostly isotropic within the resolution of DW-MRI (Dell’Acqua et al., 2010). Isotropic non-WM PVEs have been shown to affect DTI (Alexander et al., 2001; Pasternak et al., 2009) and tensor-derived measures (Metzler-Baddeley et al., 2012a). Pasternak and coworkers used constrained optimization of a bi-tensor model for “free water elimination” (FWE) in DTI (Pasternak et al., 2009), but they did not investigate GM PVEs. Metzler–Baddeley and coworkers used FWE to correct for CSF-contamination in tensor-derived measures in constrained spherical deconvolution (CSD) based tractography (Metzler-Baddeley et al., 2012a). Both fractional anisotropy (FA) and mean diffusivity (MD) were shown to increase in the presence of CSF-contamination (Pasternak et al., 2009). Moreover, diffusivity metrics were shown to be more sensitive to PVEs than anisotropy metrics (Metzler-Baddeley et al., 2012a). However, FWE-based approaches are not suitable for GM-contaminated regions.

In HARDI methods, very few studies account for the non-WM PVEs. The “ball and stick” model is the only method, which initially included an isotropic compartment and could be extended into multiple fiber orientations (Behrens et al., 2003, 2007; Jbabdi et al., 2012). In another study involving HARDI methods, isotropic PVEs were dampened by using adaptive regularization in the iterative Richardson–Lucy deconvolution algorithm (Dell’Acqua et al., 2010). Other methods that also account for isotropic compartments include diffusion basis spectrum imaging (Wang et al., 2011) and diffusion decomposition (Yeh et al., 2011; Yeh and Tseng, 2013). However, the non-WM PVEs are not taken into account and have not been studied earlier in CSD, one of the most popular, clinically feasible and readily available HARDI methods (Leemans et al., 2009; Tournier et al., 2012).

In this work, we perform simulations to assess non-WM PVEs in CSD (Tournier et al., 2004, 2007). This kind of comprehensive analysis has not been performed before, although the method is widely used and the consequences may be significant when studying the connectivity between GM regions. In addition, we analyze the proportion of voxels affected by isotropic PVEs, and present the fiber orientation distribution functions (fODFs) estimated with CSD in real data affected by non-WM PVEs.

 MATERIALS AND METHODS

We investigated the isotropic PVEs caused by non-WM tissue on fODFs estimated with CSD. DW signals were simulated with varying diffusion weightings, signal-to-noise ratios, fiber configurations, and tissue fractions. In addition, experiments with real data were performed.

ESTIMATION OF FIBER ORIENTATIONS WITH CONSTRAINED SPHERICAL DECONVOLUTION

In CSD, the full fODF is deconvolved from the DW signal using a kernel constructed from the single-fiber response function (RF), which can be estimated from the data (Tournier et al., 2004; Tax et al., 2014). During the deconvolution procedure, constraints are imposed to suppress the negative peaks in the fODF (Tournier et al., 2007, 2008). The number of distinct gradient directions limits the maximum order of the spherical harmonics (SH) decomposition, which can be used in the estimation in the fODF. However, the constraints used to suppress the negative peaks in the fODF can be exploited to estimate higher order solutions and thus, describe more complex fODFs. This is called super-resolved CSD (Tournier et al., 2007).

To find the peaks of the fODF estimated with CSD, a Newton optimization algorithm was used to extract the local maxima of the fODF directly based on the SH decompositions (Jeurissen et al., 2013). Optimization was initialized on a dense set of uniformly distributed spherical sample points. A threshold of 33% of the maximum amplitude of the fODF was used to discard small peaks. A maximum of six of the highest peaks were identified. The peaks were clustered around the peaks of the average fODF calculated over all simulation repetitions performed with the same parameter configuration. Peaks further away than half of the crossing angle (with an upper limit of 35°) from any of the peaks of the average fODF were not included in the clusters. A mean dyadic tensor was then used to derive the mean orientation for each of the identified fiber clusters (Basser and Pajevic, 2000; Jones, 2003). This orientation was then compared to the true orientations of the fiber bundles. Peaks in clusters that were less than half of the crossing angle (with an upper limit of 35°) from the true orientations were defined as true, and rest of the peaks, also if they were not assigned to a cluster, as false. From the true clusters, accuracy and precision (95th percentile confidence interval, CI) with respect to the orientation of the mean dyadic tensor were calculated.

 SIMULATION OF THE DW SIGNAL WITH PVES

Two crossing WM fiber configurations were simulated with equal weights. The orientation of the first fiber bundle was randomly selected, after which the orientation of the second fiber bundle was calculated in spherical coordinates with the defined crossing angle. The resulting angle was verified to be correct in each case.

Then, the DW signal was simulated separately for different tissue types, and the resulting signals were combined assuming no exchange between the compartments (Leemans et al., 2005). The number of gradient directions uniformly distributed on the unit hemisphere was 64 (Jones et al., 1999). To eliminate any bias caused by the gradient orientations, a different gradient set was used for each simulated DW signal. Signal from the specific WM fiber configurations was combined with isotropic CSF and GM compartments. In addition, air compartments were simulated to investigate only the effect of reduced signal of the WM compartment without any isotropic diffusion. Derived based on Basser and Jones (2002), the combined simulated DW signal S is:

S=(1−fisot)(ffibere−Trace(bDfiber1)+(1−ffiber)e−Trace(bDfiber2))+fisote−Trace(bDisot),

where fisot is the fraction of isotropic volume, ffiber is the fraction of the first fiber compartment with respect to the WM compartment, b is the b-matrix summarizing the attenuation in all three directions of the diffusion tensor (including information of the diffusion-weighting and the gradient orientations; Mattiello et al., 1997), and Dfiber1, Dfiber2 and Disot are the diffusion tensors of the two WM fibers and the isotropic compartment respectively (Basser et al., 1994b). The diffusion tensors were created with the following values. The MD for the simulation of different tissue types was 0.002 mm2/s for CSF, 0.0007 mm2/s for WM and GM (Dell’Acqua et al., 2010), and for air the signal was assumed to be zero. The FA was 0.8 for the WM signal and 0 for other tissue types. Rician noise was added to the combined DW signal. Finally, the DW signals were decomposed into an eighth-order series of SH (maximum possible order with the number of gradient orientations used).

 SIMULATION EXPERIMENTS

We performed simulation experiments to investigate the PVEs with different tissue compartments. Simulations and analyses of the simulation experiments were performed in Matlab (The MathWorks, Inc., Natick, MA, USA), by using dedicated software programmed by the authors.

The fraction of isotropic GM, CSF or air volume was varied from 0.00 to 0.95 with intervals of 0.05. We analyzed angles between the fiber populations in configurations ranging from 40° to 90° and with diffusion weightings (b-values) from 1000 to 3500 s/mm2. Signal-to-noise ratio (SNR) was calculated with respect to the non-diffusion weighted signal and simulated from 10 to 60, also generating a noiseless version of the DW signal. We performed 1000 repetitions with different noise realizations (resulting in Rician distributed data) for each parameter configuration. The fODFs were estimated from the simulated DW signals with CSD or super-resolved CSD using maximum orders of the SH from 4 to 14.

In addition to the isotropic volume fraction (VF) and PVE type (GM, CSF, or air) only one parameter at a time was investigated. The default values for the non-varying parameters were: b-value: 3000 s/mm2; angle between the crossing fiber configurations: 70°; SNR: 30. The default maximum order of the SH was 8 for CSD and 12 for super-resolved CSD.

 ACQUISITION AND ANALYSIS OF REAL DATA

High angular resolution DW data were acquired on a 3T MRI system with a 32-channel head coil. The subject gave written informed consent to participate in this study under a protocol approved by the local ethics committee. A single-shot echo-planar imaging (EPI) sequence was used with TR = 8100 ms, TE = 116 ms and 2.5 mm × 2.5 mm × 2.5 mm voxel size. The field of view (FOV) was 240 × 240 mm2 with a 96 × 96 acquisition matrix and the number of excitations (NEX) was 1. Fifty-four axial slices were imaged with 2.5 mm thickness and no gap. Diffusion sensitizing gradients with a b-value of 2800 s/mm2 were applied along 75 non-collinear directions. Ten images without diffusion weighting (b = 0 s/mm2) were acquired, of which one was acquired with reverse phase-encoding, for the purpose of EPI distortion correction. High-resolution anatomical T1-weighted images were acquired using a 3D magnetization-prepared rapid gradient-echo (MPRAGE) sequence (Mugler and Brookeman, 1990) with TR = 1900 ms, TE = 2.52 ms, TI = 900 ms and 1 mm × 1 mm × 1 mm voxel size (flip angle = 9° and NEX = 1). FOV was 250 mm × 250 mm × 176 mm with a 256 × 256 × 176 acquisition matrix.

The DW data were corrected for subject motion and eddy current induced distortions (Leemans and Jones, 2009; Andersson et al., 2012), and TOPUP was used to correct for EPI distortions (Andersson et al., 2003). The MRtrix package (J-D Tournier, Brain Research Institute, Melbourne, Australia, http://www.brain.org.au/software/; Tournier et al., 2012) was used for visualization of the real data. Tissue VFs for the DW data were estimated from the T1-weighted images, using a similar approach as presented by Smith et al. (2012).

The percentage of WM voxels affected by significant non-WM PVEs was estimated from real data. WM voxels were defined using a threshold of 25% WM tissue. The voxels with PVEs were estimated by using two threshold values: 25 and 10% non-WM volume.

 RESULTS

First, results of the simulation experiments are presented. Figure 1 shows the effects of isotropic non-WM VF in CSD (Figures 1A–D) and super-resolved CSD (Figures 1E,H). The bias and the 95% CI of the fiber orientations extracted with CSD are presented in Figures 1A,B. We also studied the effects on the number of correctly and falsely identified peaks (Figures 1C,D). The number of falsely identified peaks increased and the precision of the identified fiber orientations decreased, when the isotropic VF increased. The effects were stronger in GM than in CSF or air. However, the accuracy of the identified fiber orientations and the number of true peaks identified did not change until very high non-WM fractions. The similar performance with CSF and air PVEs using a high b-value indicates that the effect in CSF is mostly an SNR effect, which is clearly not the case in GM. Default values were used for the other parameters as specified in the methods section. The non-WM PVEs in super-resolved CSD, using up to 12th-order SH, started to affect the precision and the number of false peaks (Figures 1F,H) with lower fractions than when using up to eight-order SH. Accuracy remained high and was similar to the results when using SH up to eight-order. However, the ability to detect the two correct fiber orientations stayed higher with high isotropic fractions than while using SH up to eight order.

Effect of isotropic volume fraction with CSF, GM, and air to bias, 95% confidence interval (CI), and the number of correct and false peaks estimated with CSD using up to eighth-order SH (A–D) and super-resolved CSD using up to 12th-order SH (E–H; with diffusion weighting 3000 s/mm2, angle 70°, and SNR 30).

An illustration of the estimated fODFs based on only one noise realization per fraction is shown in Figure 2. The false peaks became more numerous and the correct peaks lost precision, when the isotropic VF increased. Next, the effects of varying maximum SH orders, diffusion weightings, SNRs, and angles between the two crossing fiber configurations were analyzed, while keeping the isotropic non-WM fraction constant at 0.5.

An illustration of fODFs estimated with CSD with GM partial volume increasing from 0 to 75% (with up to eighth-order SH, diffusion weighting 3000 s/mm2, angle 70°, and SNR 30). The blue lines correspond to the correct fiber orientations.

In Figure 3, the effects of maximum SH order on the non-WM PVEs are shown. The 95% CI and the number of false peaks increased when higher maximum SH orders were used. Bias was low with all orders except for the lowest maximum order 4 with GM PVEs. However, the correct peaks could be found properly with GM PVEs even with the lowest order, but not with CSF or air PVEs.

The effect of varying maximum SH orders with 50% non-WM partial volume (CSF, GM, and air) to bias (A), 95% confidence interval (CI) (B), and the number of correct (C) and false peaks (D) estimated with CSD (with diffusion weighting 3000 s/mm2, angle 70°, and SNR 30). For comparison, 100% WM measures are provided.

In Figure 4, the effects of varying diffusion weightings to the non-WM PVEs are shown. The 95% CI and the number of false peaks increased when diffusion weighting decreased. The number of false peaks was high and the precision of the correct peaks was low under GM PVEs compared to CSF, air, or 100% WM regardless of the diffusion weighting. The difference between CSF and air PVEs was visible only with very low diffusion weightings of 1000–1500 s/mm2.

The effect of varying diffusion weightings with 50% non-WM partial volume (CSF, GM, and air) to bias (A), 95% confidence interval (CI) (B), and the number of correct (C) and false peaks (D) estimated with CSD (with up to eighth-order SH, angle 70°, and SNR 30). For comparison, 100% WM measures are provided.

The effects of SNR on the non-WM PVEs are presented in Figure 5. Figures 5A–D show effects with 50% non-WM fractions and Figures 5E–H with 75% non-WM fractions. With 50% PVEs, increasing SNR improved the precision and reduced the number of false peaks identified. However, with 75% non-WM fractions, increasing SNR could not improve the situation with GM PVEs, and there were problems with precision also with high SNRs under CSF PVEs.

The effect of varying SNRs with 50% (A–D) and 75% (E–H) non-WM partial volume (CSF, GM, and air) to bias, 95% confidence interval (CI), and the number of correct and false peaks estimated with CSD (with up to eighth-order SH, angle 70°, and diffusion weighting 3000 s/mm2). For comparison, 100% WM measures are provided.

Figure 6 shows the effects of varying angle between the two crossing fiber configurations. With an angle of 40° between the two fiber configurations, the correct peaks could not be properly identified. However, with an angle of 50°, they could still be reliably detected without isotropic PVEs, but any type of non-WM volume caused a decrease in the fraction of the correct peaks identified (Figure 6C). With higher angles, the correct peaks were identified correctly and without more bias than in pure WM (Figure 6A). The precision of the identified fiber orientations and the number of false peaks identified improved when the angle between the fiber configurations increased (Figures 6B,D).

The effect of varying angle between the two crossing fiber configurations with 50% non-WM partial volume (CSF, GM and air) to bias (A), 95% confidence interval (CI) (B), and the number of correct (C) and false peaks (D) estimated with CSD (with up to eighth-order SH, diffusion weighting 3000 s/mm2, and SNR 30). For comparison, 100% WM measures are provided.

From real data, we estimated that 35.7% of WM voxels, defined to have at least 25% WM volume, had significant PVEs with non-WM tissue, also defined to be more than 25% VF. Lowering the non-WM tissue threshold to 10%, the proportion of WM voxels affected by PVEs increased to 46.8%. Of these voxels with non-WM PVEs, 96.0% were affected by PVEs with GM and 5.3% with CSF.

Figure 7 shows the fODFs estimated with CSD, using up to eight order SH, from real data overlaid on the WM tissue probability map of corona radiata extending towards cortical GM. The areas where WM interfaces with GM were affected both with CSD and super-resolved CSD. A large amount of voxels in the area had significant PVEs (gray-colored voxels), and perpendicular or spurious peaks appeared in the voxels with no apparent anatomical origin.

Illustration of the effects of isotropic non-WM PVEs with GM in superior corona radiata with real data. WM tissue probability map estimated from high-resolution anatomical MRI is visualized in the background, and fODFs estimated with CSD using up to eighth-order SH are overlaid.

Figure 8 shows the effect of CSF PVEs on the estimation of fODFs at the interface between the corpus callosum and CSF. Spurious orientations can be noticed, but they are much smaller in amplitude and the principal fiber orientation can still be clearly distinguished.

Illustration of the effects of isotropic non-WM PVEs with CSF in corpus callosum with real data. WM tissue probability map estimated from high-resolution anatomical MRI is visualized in the background, and fODFs estimated with CSD using up to eighth-order SH are overlaid.

 DISCUSSION

We studied the effects of isotropic non-WM partial volume on the fiber orientations estimated with CSD and super-resolved CSD by performing extensive simulations and real data experiments. CSD is a widely used method and knowledge about the implications of non-WM PVEs should be augmented.

Our results demonstrate that although CSD is efficient in the detection of PVEs caused by complex fiber configurations within a voxel, problems arise in the detection of the fODFs in the case of non-WM PVEs, which we estimated to be present in 35–50% of the WM voxels. As shown in Figure 1, the precision of the detected fiber orientations decreases and false peaks appear in the fODFs. This effect is most prominent for GM PVEs. The increase in bias with very high isotropic fractions may be at least partly caused by the inability to distinguish reliably both of the correct fiber orientations.

Part of these effects is due to the reduction of relative SNR in the WM tissue, which is caused by the decreased WM volume in the voxel, and not the isotropic diffusion properties of the non-WM tissue. PVEs with CSF tissue are mostly due to this effect, as shown by the similarity to air PVEs (Figure 1). Another part of the effects is caused by the isotropic diffusion, which invalidates the single fiber RF originally designed for pure WM. The more prominent PVEs in GM than in CSF or air are caused by this effect.

In addition, we showed that the PVEs increased when the maximum SH order increased. Therefore, the high maximum SH orders, although able to improve the angular resolution (Tournier et al., 2007), should be used with caution in the estimation of the fODFs under significant non-WM PVEs. Although the maximum angular frequency in the DW data is relatively low (Tournier et al., 2013), the fODFs contain higher angular frequencies, so a higher maximum SH order could still be useful in the estimation of the fODFs within pure WM regions. The use of lower diffusion weighting than generally used in HARDI sequences (i.e., less than 3000 s/mm2) increased the PVEs. Larger crossing angles could be detected with higher precision. With higher SNRs, moderate PVEs could be handled better, but high PVEs continued to decrease precision and increase the number of false peaks especially in GM.

Based on these results, we provide the following advice on how to operate CSD to maintain reasonable precision and number of false peaks under non-WM PVEs. Conditions with 95% CI lower than 20° and less than one identified false peak were considered reasonable. Thus, we suggest acquiring data with a high diffusion-weighting 2500–3000 s/mm2, and a reasonable SNR (~25–30). To extract the fiber orientations with CSD in regions with GM PVEs, we suggest using relatively low, from 6 to 8, maximum SH orders to minimize the loss in precision and the increase in the number of identified false fiber orientations. Nevertheless, the identified fiber orientations should be considered unreliable with higher than 60% GM and higher than 80% CSF VFs.

The isotropic PVEs, present in a significant proportion of WM voxels, lead to decreased precision and a high number of false peaks in the fODFs estimated with CSD, which in turn affects subsequent tractography algorithms, and may introduce false positives and hinder tract propagation into the cortex or near subcortical GM tissue. An algorithm already exists to discard tracts based on their anatomical feasibility and thus, only accept tracts that correctly propagate to the cortex (Smith et al., 2012). However, enabling the tracts to propagate properly into the cortex or adjacent to subcortical GM tissue would reduce the time needed for tracking and improve the precision of anatomically feasible tracts. Especially in connectomics, where reliable reconstruction of the fiber orientations profiles at the GM–WM interface is required to compute connectivity matrices, taking isotropic PVEs into account will be valuable.

Limitations of this study include the restriction to only one HARDI method, although it is one of the most commonly used ones (Metzler-Baddeley et al., 2012b; Emsell et al., 2013; Forde et al., 2013; Kristo et al., 2013; McGrath et al., 2013a,b; Reijmer et al., 2013a,b; Thompson et al., 2014). Previous studies indicate that the non-WM PVEs are present in DW-MRI in general (Alexander et al., 2001; Pasternak et al., 2009; Dell’Acqua et al., 2010; Metzler-Baddeley et al., 2012a). While some of the analysis methods already acknowledge or account for these PVEs (Behrens et al., 2003, 2007; Pasternak et al., 2009; Dell’Acqua et al., 2010; Wang et al., 2011; Yeh et al., 2011; Jbabdi et al., 2012; Yeh and Tseng, 2013), many of the currently used methods do not. For example, in CSD they have not yet been taken into account, and no detailed investigation about these effects had been performed previously. It is likely that also other methods which do not appropriately account for the non-WM PVEs will suffer from similar consequences. An additional limitation of this study is that there is no ground truth available in real data. Considering the clear effects demonstrated in the simulations, it is reasonable to assume that the spurious fiber orientations visible at the tissue interfaces are in fact false peaks also in real data. However, further experiments with real data are still necessary to completely understand the phenomenon and its effects in tractography. This would in turn help in the development of improvements for the fODF estimation with CSD, applicable also in real data, and thus allow improved tracking especially in the WM–GM interface.

In conclusion, we studied the effects of isotropic non-WM PVEs in CSD and found decreased precision and increased number of false peaks in the estimated fODFs. The effect was more pronounced with GM tissue. Considering the clear effects present in real and simulated data and the large proportion of WM voxels affected, it is important to take the non-WM PVEs into account in the extraction of fiber orientations with CSD. Therefore, we provide simple recommendations for the parameters used in the acquisition and the analysis, but acknowledge the need for more sophisticated methods to account for non-WM tissue in CSD.

 Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

This work was supported by the Fund for Scientific Research-Flanders (FWO), and by the Interuniversity Attraction Poles Program (P7/11) initiated by the Belgian Science Policy Office.

 REFERENCES Alexander A. L. Hasan K. M. Lazar M. Tsuruda J. S. Parker D. L. (2001). Analysis of partial volume effects in diffusion-tensor MRI. Magn. Reson. Med. 45 770–780.10.1002/mrm.1105 Andersson J. L. Skare S. Ashburner J. (2003). How to correct susceptibility distortions in spin-echo echo-planar images: application to diffusion tensor imaging. Neuroimage 20 870–888.10.1016/S1053-8119(03)00336-7 Andersson J. L. R. Xu J. Yacoub E. Auerbach E. Moeller S. Ugurbil K. (2012). “A comprehensive Gaussian process framework for correcting distortions and movements in diffusion images,” in Proceedings of the 20th Annual Meeting of ISMRM, Melbourne 2426. Basser P. J. Jones D. K. (2002). Diffusion-tensor MRI: theory, experimental design and data analysis: a technical review. NMR Biomed. 15 456–467.10.1002/nbm.783 Basser P. J. Mattiello J. LeBihan D. (1994a). MR diffusion tensor spectroscopy and imaging. Biophys. J. 66 259–267.10.1016/S0006-3495(94)80775-1 Basser P. J. Mattiello J. LeBihan D. (1994b). Estimation of the effective self-diffusion tensor from the NMR spin echo. J. Magn. Res. B 103 247–254.10.1006/jmrb.1994.1037 Basser P. J. Pajevic S. (2000). Statistical artifacts in diffusion tensor MRI (DT-MRI) caused by background noise. Magn. Reson. Med. 44 41–50.10.1002/1522-2594(200007)44:1<41::AID-MRM8>3.0.CO;2-O Basser P. J. Pajevic S. Pierpaoli C. Duda J. Aldroubi A. (2000). In vivo fiber tractography using DT-MRI data. Magn. Reson. Med. 44 625–632.10.1002/1522-2594(200010)44:4<625::AID-MRM17>3.0.CO;2-O Behrens T. E. J.Berg H. J. Jbabdi S. Rushworth M. F. S.Woolrich M. W. (2007). Probabilistic diffusion tractography with multiple fibre orientations: what can we gain? Neuroimage 34 144–155.10.1016/j.neuroimage.2006.09.018 Behrens T. E. J.Woolrich M. W. Jenkinson M. Johansen-Berg H. Nunes R. G. Clare S. (2003). Characterization and propagation of uncertainty in diffusion-weighted MR imaging. Magn. Reson. Med. 50 1077–1088.10.1002/mrm.10609 Conturo T. E. Lori N. F. Cull T. S. Akbudak E. Snyder A. Z. Shimony J. S. (1999). Tracking neuronal fiber pathways in the living human brain. Proc. Natl. Acad. Sci. U.S.A. 96 10422–10427.10.1073/pnas.96.18.10422 Dell’Acqua F. Rizzo G. Scifo P. Clarke R. A. Scotti G. Fazio F. (2007). A model-based deconvolution approach to solve fiber crossing in diffusion-weighted MR imaging. IEEE Trans. Biomed. Eng. 54 462–472.10.1109/TBME.2006.888830 Dell’Acqua F. Scifo P. Rizzo G. Catani M. Simmons A. Scotti G. (2010). A modified damped Richardson–Lucy algorithm to reduce isotropic background effects in spherical deconvolution. Neuroimage 49 1446–1458.10.1016/j.neuroimage.2009.09.033 Descoteaux M. Angelino E. Fitzgibbons S. Deriche R. (2007). Regularized, fast, and robust analytical Q-ball imaging. Magn. Reson. Med. 58 497–510.10.1002/mrm.21277 Emsell L. Leemans A. Langan C. Van Hecke W. Barker G. J. McCarthy P. (2013). Limbic and callosal white matter changes in euthymic bipolar I disorder: an advanced diffusion magnetic resonance imaging tractography study. Biol. Psychiatry 73 194–201.10.1016/j.biopsych.2012.09.023 Forde N. J. Ronan L. Suckling J. Scanlon C. Neary S. Holleran L. (2013). Structural neuroimaging correlates of allelic variation of the BDNF Val66met polymorphism. Neuroimage 90 280–289. 10.1016/j.neuroimage.2013.12.050 Frank L. R. (2001). Anisotropy in high angular resolution diffusion-weighted MRI. Magn. Reson. Med. 45 935–939.10.1002/mrm.1125 Frank L. R. (2002). Characterization of anisotropy in high angular resolution diffusion-weighted MRI. Magn. Reson. Med. 47 1083–1099.10.1002/mrm.10156 Jansons K. M. Alexander D. C. (2003). Persistent angular structure: new insights from diffusion magnetic resonance imaging data. Inverse Problems 19 103110.1088/0266-5611/19/5/303 Jbabdi S. Sotiropoulos S. N. Savio A. M. Graña M. Behrens T. E. (2012). Model-based analysis of multishell diffusion MR data for tractography: how to get over fitting problems. Magn. Reson. Med. 68 1846–1855.10.1002/mrm.24204 Jeurissen B. Leemans A. Jones D. K. Tournier J. D. Sijbers J. (2011). Probabilistic fiber tracking using the residual bootstrap with constrained spherical deconvolution. Hum. Brain Mapp. 32 461–479.10.1002/hbm.21032 Jeurissen B. Leemans A. Tournier J. D. Jones D. K. Sijbers J. (2013). Investigating the prevalence of complex fiber configurations in white matter tissue with diffusion magnetic resonance imaging. Hum. Brain Mapp. 34 2747–2766.10.1002/hbm.22099 Jones D. K. (2003). Determining and visualizing uncertainty in estimates of fiber orientation from diffusion tensor MRI. Magn. Reson. Med. 49 7–12.10.1002/mrm.10331 Jones D. K. (2008). Studying connections in the living human brain with diffusion MRI. Cortex 44 936–952.10.1016/j.cortex.2008.05.002 Jones D. K. (ed). (2010). Diffusion MRI: Theory, Methods, and Applications. Oxford: Oxford University Press.10.1093/med/9780195369779.001.0001 Jones D. K. Horsfield M. A. Simmons A. (1999). Optimal strategies for measuring diffusion in anisotropic systems by magnetic resonance imaging. Magn. Reson. Med. 42 515–525.10.1002/(SICI)1522-2594(199909)42:3<515::AID-MRM14>3.0.CO;2-Q Jones D. K. Leemans A. (2011). “Diffusion tensor imaging,” in Magnetic Resonance Neuroimaging eds Modo M. Bulte J. (New York: Humana Press) 127–144. Kristo G. Leemans A. Raemaekers M. Rutten G. J. Gelder B. Ramsey N. F. (2013). Reliability of two clinically relevant fiber pathways reconstructed with constrained spherical deconvolution. Magn. Reson. Med. 70 1544–1556.10.1002/mrm.24602 Leemans A. Jeurissen B. Sijbers J. Jones D. K. (2009). ExploreDTI: a graphical toolbox for processing, analyzing, and visualizing diffusion MR data. Proc. Int. Soc. Mag. Reson. Med. 17 3536. Leemans A. Jones D. K. (2009). The B-matrix must be rotated when correcting for subject motion in DTI data. Magn. Reson. Med. 61 1336–1349.10.1002/mrm.21890 Leemans A. Sijbers J. Verhoye M. Van der Linden AVan Dyck D. (2005). Mathematical framework for simulating diffusion tensor MR neural fiber bundles. Magn. Reson. Med. 53 944–953.10.1002/mrm.20418 Mattiello J. Basser P. JLe Bihan D. (1997). The b matrix in diffusion tensor echo-planar imaging. Magn. Reson. Med. 37 292–300.10.1002/mrm.1910370226 McGrath J. Johnson K. O’Hanlon E. Garavan H. Gallagher L. Leemans A. (2013a). White matter and visuospatial processing in autism: a constrained spherical deconvolution tractography study. Autism Res. 6 307–319.10.1002/aur.1290 McGrath J. Johnson K. O’Hanlon E. Garavan H. Leemans A. Gallagher L. (2013b). Abnormal functional connectivity during visuospatial processing is associated with disrupted organisation of white matter in autism. Front. Hum. Neurosci. 7:434. 10.3389/fnhum.2013.00434 Metzler-Baddeley C. O’Sullivan M. J. Bells S. Pasternak O. Jones D. K. (2012a). How and how not to correct for CSF-contamination in diffusion MRI. Neuroimage 59 1394–1403.10.1016/j.neuroimage.2011.08.043 Metzler-Baddeley C. Hunt S. Jones D. K. Leemans A. Aggleton J. PO’Sullivan M. J. (2012b). Temporal association tracts and the breakdown of episodic memory in mild cognitive impairment. Neurology 79 2233–2240.10.1212/WNL.0b013e31827689e8 Mori Svan Zijl P. (2002). Fiber tracking: principles and strategies: a technical review. NMR Biomed. 15 468–480.10.1002/nbm.781 Mugler J. P. Brookeman J. R. (1990). Three-dimensional magnetization-prepared rapid gradient-echo imaging (3D MP RAGE). Magn. Reson. Med. 15 152–157.10.1002/mrm.1910150117 Pasternak O. Sochen N. Gur Y. Intrator N. Assaf Y. (2009). Free water elimination and mapping from diffusion MRI. Magn. Reson. Med. 62 717–730.10.1002/mrm.22055 Reijmer Y. D. Leemans A. Brundel M. Kappelle L. J. Biessels G. J. (2013a). Disruption of the cerebral white matter network is related to slowing of information processing speed in patients with type 2 diabetes. Diabetes Metab. Res. Rev. 62 2112–2115.10.2337/db12-1644 Reijmer Y. D. Freeze W. M. Leemans A. Biessels G. J. (2013b). The effect of lacunar infarcts on white matter tract integrity. Stroke 10.1161/strokeaha.113.001321 Smith R. E. Tournier J. D. Calamante F. Connelly A. (2012). Anatomically-constrained tractography: improved diffusion MRI streamlines tractography through effective use of anatomical information. Neuroimage 62 1924–1938.10.1016/j.neuroimage.2012.06.005 Stejskal E. O. Tanner J. E. (1965). Spin diffusion measurements: spin echoes in the presence of a time-dependent field gradient. J. Chem. Phys. 42 28810.1063/1.1695690 Tax C. M. Jeurissen B. Vos S. B. Viergever M. A. Leemans A. (2014). Recursive calibration of the fiber response function for spherical deconvolution of diffusion MRI data. Neuroimage 86 67–80.10.1016/j.neuroimage.2013.07.067 Thompson D. K. Thai D. Kelly C. E. Leemans A. Tournier J. D. Kean M. J. (2014). Alterations in the optic radiations of very preterm children: perinatal predictors and relationships with visual outcomes. Neuroimage Clin. 4 145–153.10.1016/j.nicl.2013.11.007 Tournier J. Calamante F. Connelly A. (2007). Robust determination of the fibre orientation distribution in diffusion MRI: non-negativity constrained super-resolved spherical deconvolution. Neuroimage 35 1459–1472.10.1016/j.neuroimage.2007.02.016 Tournier J. Calamante F. Connelly A. (2012). MRtrix: diffusion tractography in crossing fiber regions. Int. J. Imaging Syst. Technol. 22 53–66.10.1002/ima.22005 Tournier J. Calamante F. Connelly A. (2013). Determination of the appropriate b value and number of gradient directions for high-angular-resolution diffusion-weighted imaging. NMR Biomed. 26 1775–1786.10.1002/nbm.3017 Tournier J. Calamante F. Gadian D. G. Connelly A. (2004). Direct estimation of the fiber orientation density function from diffusion-weighted MRI data using spherical deconvolution. Neuroimage 23 1176–1185.10.1016/j.neuroimage.2004.07.037 Tournier J. Yeh C. H. Calamante F. Cho K. H. Connelly A. Lin C. P. (2008). Resolving crossing fibres using constrained spherical deconvolution: validation using diffusion-weighted imaging phantom data. Neuroimage 42 617–625.10.1016/j.neuroimage.2008.05.002 Tournier J. D. Calamante F. Connelly A. (2010). “Improved probabilistic streamlines tractography by 2nd order integration over fibre orientation distributions,” in Proceedings of the International Society for Magnetic Resource Medicine Stockholm, Swedan. Abstract 1670. Tournier J. D. Mori S. Leemans A. (2011). Diffusion tensor imaging and beyond. Magn. Reson. Med. 65 1532–1556.10.1002/mrm.22924 Tuch D. S. (2004). Q-ball imaging. Magn. Reson. Med. 52 1358–1372.10.1002/mrm.20279 Tuch D. S. Reese T. G. Wiegell M. R. Makris N. Belliveau J. W. Wedeen V. J. (2002). High angular resolution diffusion imaging reveals intravoxel white matter fiber heterogeneity. Magn. Reson. Med. 48 577–582.10.1002/mrm.10268 Vos S. B. Jones D. K. Viergever M. A. Leemans A. (2011). Partial volume effect as a hidden covariate in DTI analyses. Neuroimage 55 1566–1576.10.1016/j.neuroimage.2011.01.048 Wang Y. Wang Q. Haldar J. P. Yeh F. C. Xie M. Sun P. (2011). Quantification of increased cellularity during inflammatory demyelination. Brain 134 3590–3601.10.1093/brain/awr307 Wedeen V. J. Hagmann P. Tseng W. Y. I.Reese T. G. Weisskoff R. M. (2005). Mapping complex tissue architecture with diffusion spectrum magnetic resonance imaging. Magn. Reson. Med. 54 1377–1386.10.1002/mrm.20642 Wedeen V. J. Wang R. P. Schmahmann J. D. Benner T. Tseng W. Y. I.Dai G. (2008). Diffusion spectrum magnetic resonance imaging (DSI) tractography of crossing fibers. Neuroimage 41 1267–1277.10.1016/j.neuroimage.2008.03.036 Yeh F. CTseng W. Y. I. (2013). Sparse solution of fiber orientation distribution function by diffusion decomposition. PLoS ONE 8:e75747. 10.1371/journal.pone.0075747 Yeh F. C. Wedeen V. JTseng W. Y. I. (2011). Estimation of fiber orientation and spin density distribution by diffusion deconvolution. Neuroimage 55 1054–1062.10.1016/j.neuroimage.2010.11.087