Joint Multi-Fiber NODDI Parameter Estimation and Tractography Using the Unscented Information Filter

Tracing white matter fiber bundles is an integral part of analyzing brain connectivity. An accurate estimate of the underlying tissue parameters is also paramount in several neuroscience applications. In this work, we propose to use a joint fiber model estimation and tractography algorithm that uses the NODDI (neurite orientation dispersion diffusion imaging) model to estimate fiber orientation dispersion consistently and smoothly along the fiber tracts along with estimating the intracellular and extracellular volume fractions from the diffusion signal. While the NODDI model has been used in earlier works to estimate the microstructural parameters at each voxel independently, for the first time, we propose to integrate it into a tractography framework. We extend this framework to estimate the NODDI parameters for two crossing fibers, which is imperative to trace fiber bundles through crossings as well as to estimate the microstructural parameters for each fiber bundle separately. We propose to use the unscented information filter (UIF) to accurately estimate the model parameters and perform tractography. The proposed approach has significant computational performance improvements as well as numerical robustness over the unscented Kalman filter (UKF). Our method not only estimates the confidence in the estimated parameters via the covariance matrix, but also provides the Fisher-information matrix of the state variables (model parameters), which can be quite useful to measure model complexity. Results from in-vivo human brain data sets demonstrate the ability of our algorithm to trace through crossing fiber regions, while estimating orientation dispersion and other biophysical model parameters in a consistent manner along the tracts.

1. Introduction

Diffusion MRI (dMRI) is a non-invasive technique to study the microstructure of brain tissue. However, we need a mathematical model to interpret the diffusion weighted signal to study the microstructure of white matter fibers. Broadly, such models fall into two categories: parametric and nonparametric. The simplest parametric model is the diffusion tensor model, which describes a Gaussian estimate of the strength and diffusion orientation at each voxel (Basser et al., 1994; Alexander et al., 2002). While robust, this model can be inadequate in cases of mixed fiber presence or more complex orientations (Frank, 2002). To handle more complex diffusion patterns, parametric models have been introduced including mixtures of tensors and directional functions (Alexander et al., 2001; Anderson, 2005; Kreher et al., 2005; Parker and Alexander, 2005; Peled et al., 2006; Kaden et al., 2007; Zhang et al., 2007; Rathi et al., 2008). Several techniques attempt to reconstruct pathways based on these models. In this case, tractography is done by following the principal diffusion direction(s).

With multifiber models, care must be taken in fitting the model parameters to the recorded signal consistently, such that, the correlation in diffusion along the tract is accounted for while estimating the parameters (Malcolm et al., 2009a). As demonstrated in this paper, estimation based on information from previous estimates aids in this process.

Nonparametric techniques, unlike estimating a fixed number of fibers in parametric models, estimate an oriented distribution function (ODF) describing an arbitrary configuration of fibers. For this estimation, Tuch (2004) introduced Q-ball imaging to numerically compute the ODF via the Funk-Radon transform, and later, spherical harmonics were used to simplify the computation with an analytic form (Anderson, 2005; Hess et al., 2006; Ozarslan et al., 2006; Descoteaux et al., 2007). Another method to estimate the fiber ODF is to assume a model for the signal response of a single-fiber and use spherical deconvolution (Jansons and Alexander, 2003; Tournier et al., 2004; Jian and Vemuri, 2007; Kaden et al., 2007) to obtain a much sharper orientation profile. A good review of both parametric and nonparametric models and diffusion MRI in general can be found in Assemlal et al. (2011).

Recently, a model of diffusion was proposed by Zhang et al. (2012) called Neurite Orientation Dispersion Diffusion Imaging or NODDI, which accounted for the dispersion in orientation of the axonal fibers. In this work, the authors proposed an algorithm to estimate the parameters of the NODDI model assuming the existence of a single fiber pathway at each voxel. However, it is well known that most of the white matter tissue in the brain has crossing fibers (Tuch et al., 2002), which must be taken into account for proper analysis of the estimated microstructural parameters. Further, the correlation in diffusion (and consequently the estimated NODDI parameters) is not taken into account in the work of Zhang et al. (2012). Thus, estimating the orientation dispersion in the presence of multiple fiber crossing would be quite useful in analyzing the geometrical structure of white matter in healthy as well as diseased subjects.

1.1. Our Contributions

In this work, we propose to use the unscented information filter (UIF) based framework to perform joint NODDI parameter estimation and tractography. Earlier works have used the unscented Kalman filter (UKF) with a parametric multi-tensor model or a non-parametric spherical harmonic representation to do simultaneous model estimation and tractography (Malcolm et al., 2010; Baumgartner et al., 2012). The UIF has several advantages over the UKF as has been noted in Lee (2008), namely, it propagates the Fisher-information matrix as opposed to the state-covariance matrix leading to a significant reduction in computational load, while providing a more robust estimation of the estimated parameters. Further, the Fisher-information is calculated for each of the model parameters along the tract, providing the variance in the estimation of the model parameters.

Existing methods either use a region-of-interest (ROI) or an atlas within which the NODDI parameters are estimated and analyzed or follow the principal orientation to perform tractography. Thus, the fit is performed independently at each voxel disregarding the correlation of diffusion along the fiber path. In the UIF framework, we perform model estimation and tractography simultaneously to trace fiber tracts. The methodology assumes a “time varying" Gaussian distribution of the model parameters as the algorithm moves from one location to the next. This allows for consistent estimation of the NODDI parameters, while allowing for smooth estimation of the fiber tracts. In the current work, we demonstrate the performance of our method in the case of single fiber and 2-fiber NODDI model. To the best of our knowledge, this is a first extension of the NODDI model to the multi-fiber case within a tractography framework. Thus, the ability to trace multiple crossing fibers while robustly estimating the NODDI model parameters is one of the major contributions of this work. We demonstrate our technique on in-vivo human data set from the human connectome project (HCP) data set. We expect, that the proposed method will be quite useful to study the tissue microstructure in several disorders.

2. Approach

The main idea of our approach is to trace the local fiber orientations employing the three compartment NODDI model (intra-cellular, extra-cellular, and isotropic free-water) and using the estimation at previous positions to guide estimation at the current position. In a recursive fashion, the information filter estimates the model at the current position, moves a step in the most consistent direction, and then begins estimation again. In this case, the model parameters, which form the state of the filter, are assumed to have a Gaussian distribution, whose mean and covariance change at each step. Thus, the space of possible solutions (among the infinitely many possibilities) is greatly reduced leading to better accuracy in resolving individual orientations and yielding inherently smooth tracts despite the presence of noise and uncertainty. This is in complete contrast to existing methods, where the model parameters at each voxel are estimated from an arbitrary initialization or a large number of random initializations to finally obtain a solution that best fits the data. This could potentially lead to vastly different solution even in neighboring voxels since the space of possible solutions is large. In the proposed method, since each iteration begins with a near-optimal solution provided by the previous estimation, and since the step size is kept to be very small, the UIF quickly converges to the optimal solution and many local minima are naturally avoided.

2.1. Modeling Local Fiber Orientations

The three compartment NODDI model consists of an intra-cellular, an extra-cellular and an isotropic compartment (Zhang et al., 2012; Daducci et al., 2015). The normalized signal E can be written as:

$\begin{array}{ll}E=(1-V_{iso})(V_{ic}{E}_{ic}+(1-V_{ic})E_{ec})+V_{iso}{E}_{iso},& (1)\end{array}$

where E_ic and V_ic are the normalized signal and volume fractions of the intra-cellular compartment; E_ec is the normalized signal corresponding to the extra-cellular compartment; and E_iso and V_iso are the normalized signal and volume fractions of the isotropic compartment. The space bounded by the membrane of neurites is called the intra-cellular compartment. In this paper, we adopt a model where neurites are modeled as a set of cylinders of zero radius capturing the unrestricted nature of diffusion along the neurites and restricted diffusion perpendicular to them, i.e.,

$\begin{array}{l}{E}_{ic}(u)={\displaystyle {\int }_{{\mathbb{S}}^2}f}(n)e^{-b{d}_{\parallel }{(u·n)}^2}dn,\end{array}$

where b is the b-value, u is the gradient direction; f(n) is the orientation dispersion function around vector n; and d_∥ is intrinsic diffusivity along the cylinder. As in Zhang et al. (2012), the watson distribution is used to model the orientation dispersion function f : 𝕊² → ℝ, f(n) = M(1/2, 3/3, κ)⁻¹e^κ(m·n)2, where m is the mean orientation about which the dispersion is symmetric and κ is the concentration parameter, which determines the extent of orientation dispersion along the mean orientation; and M is the confluent hypergeometric function. The orientation dispersion index is given by $OD=\frac{2}{\pi }arctan\frac{1}{\kappa }$.

The space around the neurites is called the extra-cellular compartment. In this compartment, water diffusion is assumed to be hindered and is modeled as a Gaussian anisotropic tensor, with

$\log E_{ec}(u)=-b{u}^T\left({\displaystyle {\int }_{{\mathbb{S}}^2}f}(n)D(n)dn\right)u,$

where D(n) is a cylindrically symmetric diffusion tensor with principal orientation n, and d_∥ and d_⊥ are the coefficients of diffusion parallel and perpendicular to n. Both these parameters can expressed in terms of the intra-cellular fraction V_ic and the concentration parameter κ (see Zhang et al., 2012 for more details). The isotropic compartment is modeled with an isotropic Gaussian diffusion with diffusivity set to d_iso.

2.2. Estimating the NODDI Parameters

As described in Zhang et al. (2012), we fix the following parameters within the NODDI model: $d_{\parallel }=1.7\times 10^{-3}m{m}^2{s}^{-1}$ and $d_{iso}=3\times 10^{-3}m{m}^2{s}^{-1}$. Consequently, for a single fiber NODDI model, the free parameters are: V_ic, κ, m and V_iso. For the two-fiber NODDI model, we used the following formulation:

$\begin{array}{l}E=(1-V_{iso})(0.5\{V_{ic1}{E}_{ic1}+(1-V_{ic1})E_{ec1}\}\\ +0.5\{V_{ic2}{E}_{ic2}+(1-V_{ic2})E_{ec2}\})+V_{iso}{E}_{iso},& (2)\end{array}$

In this case, the free parameters to be estimated are: V_ic1, κ₁, m₁, V_ic2, κ₂, m₂, and V_iso. Given the measured signal at a particular voxel, we want to estimate the above model parameters that best explain the signal. We propose to achieve this using the unscented information filter (UIF), which is a recursive non-linear least squares estimator, giving the maximum likelihood estimate of the model parameters. Further, the Fisher information computed by the UIF provides a lower bound on the precision with which the model parameters can be estimated given the data. This statistic can be summarized into a single number as the estimated uncertainity at each point, by computing the matrix norm of the estimated covariance matrix. This uncertainity measure can be quite useful in removing unlikely fibers (with high uncertainity) from the tractography obtained. Note that, such information is typically not available using other model estimation methods.

As in streamline tractography, we treat the fiber as the trajectory of a particle which we trace out. At each step, we examine the measured signal at that position, use that measurement to update our model parameters within the filter, and propagate forward in the most consistent direction. To use a state-space formulation for estimating the model parameters, we need the following application-specific definition of four filter components:

1. The system state (x): which are the model parameters in our case,

2. The state transition function (s): how the model changes as we trace the fiber,

3. The observation function (h): how the signal appears given a particular model state,

4. The measurement (y): the actual signal obtained from the scanner.

For our state (x), we directly use the parameters of the NODDI model (in the case of 2-fibers):

$\begin{array}{ll}\text{x}=\left[V_{ic1}{\kappa }_1{\text{m}}_1{V}_{ic2}{\kappa }_2{\text{m}}_2{V}_{iso}\right].& (3)\end{array}$

For the state transition function s, we assume identity dynamics since the local fiber configuration does not undergo drastic change as it moves from one location to the next, when the step size is kept very small. The predicted signal is computed using: y = h[x] given by Equations (1) or (2) depending on the model used (1 fiber or 2-fiber NODDI) and our actual measurement is the actual signal interpolated directly on the diffusion weighted images at the current position.

It is important to note that we chose the unscented Information filter for its low computational complexity compared to the UKF. In the UKF formulation, one not only estimates the state, but also the state covariance matrix P, whose dimension is k × k, where k is the number of diffusion weighted gradients. At each step, the covariance matrix P is updated, which involves inversion of a large matrix that depends on k. This makes the filter computationally very expensive, especially in cases like the human connectome data, which has about 270 measurements. In contrast, in the UIF filter, only the Fisher-Information matrix is computed and propagated, which involves inversion of a matrix whose maximal dimension is substantially small, i.e., the length of the state vector. Thus, in the case of the two-fiber NODDI model, one only needs to invert the information matrix of size 11 × 11. This significantly increases the computational speed of the algorithm.

We use the unscented information filter with constraints, where the final solution is projected onto the physiological range of the parameters using a quadratic programming problem. For example, V_iso, V_ic are constrained to lie between [0, 1], while κ is always set to be positive. For a more thorough treatment of the UKF and UIF filters with constraints (see Lee, 2008; Malcolm et al., 2009b, 2010). In Appendix, we describe the mathematical equations for the prediction and measurement update steps of the UIF filter.

2.3. The tractography Algorithm

The UIF filter provides the maximum a-posteriori estimate of the model parameters, given the signal at each location. We embed this into a tractography framework, as described in detail in Malcolm et al. (2010). Briefly, we begin by initializing the UIF filter at each seed point by using the coarse-grid search method as described in Zhang et al. (2012). Subsequently, we run the UIF filter at each seed point, which provides an estimate of the model parameters, including the principal orientation(s) of the fiber bundle. A small step (using a fixed pre-determined step size) is taken along the direction of the principal fiber orientation that is most consistent as compared to the previous estimate. At sub-voxel locations, the signal is interpolated using an isotropic Gaussian kernel with its width (variance) given by the smallest voxel length in any direction. This interpolated signal is then used as the acquired measurement (y), within the UIF filter to estimate the model parameters as well as the Fisher information and the covariance matrix. In a loop, the fiber is then traced until a termination criteria is reached. In the present case, we used generalized fraction anisotropy (GFA) threshold of 0.08 and κ of 0.06 as the stopping criteria.

3. Experiments

We tested the proposed algorithm by tracing several fiber bundles from the in-vivo human data set obtained from HCP (Van Essen et al., 2013). All b-values of {1000, 2000, 3000}s/mm² were used to perform whole brain tractography using single fiber and two-fiber NODDI models.

On a 2.4 GHz processor with 16 cores, the UKF filter (with 2-fiber NODDI model) required about 72 hours of computational time for a single subject whole-brain tractography with 1 seed per voxel on the HCP data set with 270 gradient directions. With the same set of parameters, the UIF-based whole brain tractography required about 30 hours, an improvement by more than a factor of 2.

Whole brain tractography was performed using with the following parameter settings: expected rate of change of orientation q_m = 0.002, rate of change of κ (dispersion), q_κ = 0.015 and expected rate of change of intracellular and isotropic volume fraction was set to q_ic = q_iso = 0.0007. The parameter r_s is akin to a regularization parameter that should be set based on the expected noise level in the data. In the case of HCP data, we set it to 0.02, based on the noise level in the data. Thus, for higher noise in the data, r_s should be increased, which implicitly implies that the algorithm will trust the model more than the data. On the other hand, for high SNR data, r_s should be reduced to allow the tractography method to trust the data more than the expected model. For all the experiments, we set the step length parameter to 0.5 mm.

The following figures show the traced fiber bundles extracted using the white matter query language (WMQL) (Wassermann et al., 2013). In particular, Figure 1 shows the arcuate fasciculus traced using the 1-fiber NODDI model. Estimates of several diffusion measures of interest, such as, intracellular volume fraction, orientation dispersion, isotropic volume fraction, normalized mean squared error (NMSE) in fitting the data and uncertainity in the estimated parameters are shown along the tract with the standard single diffusion tensor based FA map shown in the background.

FIGURE 1

Figure 1. Visualization of intracelluar volume fraction, fiber dispersion, isotropic volume fraction, data fitting error and uncertainity in parameter estimation using 1-fiber noddi model in the arcuate fasciculus. The background slice is the single tensor FA map.

In Figure 2, the corticospinal tract is shown as traced using the 1-fiber NODDI model. In the coronal view, the background is a slice of fiber orientation obtained using the method in Daducci et al. (2015). Note the similarity in the obtained measurements using both the methods, i.e., high dispersion in the CSF and gray matter areas and low in the deep white matter regions.

FIGURE 2

Figure 2. Corticospinal tract showing orientation dispersion from two different views. In the sagittal view, the background is FA, while in the coronal view the 1-fiber NODDI dispersion map (Daducci et al., 2015) is shown in the background.

Figure 3 shows results for the 2-fiber NODDI model obtained using the UIF filter. The cortico-spinal tract (shown in color) and the superior-longitudinal fasiculus II (SLF-II) intersect in the centrum-semiovale region. As can be seen, the proposed algorithm is able to trace fibers through crossing regions, which is not possible using the 1-fiber NODDI model. Also, the lateral fibers of the corticospinal tract (CST) that go to the hand and face area are missing in the 1-fiber NODDI model, but can be nicely traced using the 2-fiber NODDI model.

FIGURE 3

Figure 3. Corticospinal tract (CST) and SLF-II traced using the 2-fiber NODDI model. Intersecting fibers are seen in the centrum-semiovale region.

Results shown in Figure 4 demonstrate how the proposed method can be used to trace the CST to the hand and face motor areas. The estimated model parameters, such as intracellular volume fraction (for both the fibers of the 2-fiber NODDI model), orientation dispersion and the isotropic volume fraction are also shown. Note once again, that we show the orientation dispersion and intracellular volume fraction for the second fiber as well, although the tracts were obtained by following the primary fibers. The results demonstrate the smooth and robust estimation of all the model parameters, including the ones for the second fiber.

FIGURE 4

Figure 4. Corticospinal tract traced using the 2-fiber NODDI model. In the background is a slice of orientation dispersion obtained from 1-fiber NODDI model of Daducci et al. (2015). Several NODDI specific measures are shown along the tracts along with the data fitting error, which is below 3% in most cases.

The UIF filter can also estimate the uncertainty in the estimation of the model parameters. This information can be used to removed unlikely fibers or false positives from the tractography, which offers a powerful way to automatically detect such fibers and remove them. One such example is shown in Figure 5 as pointed by the white arrows. Note that, to the best of our knowledge, only the UIF (and UKF) based tractography methods allow an inherent way to detect unlikely fibers.

FIGURE 5

Figure 5. Estimated uncertainty in the model parameters in the CST (left) and IOFF (right) fiber bundles. Fibers with high uncertainty (likely false positives), as marked by white arrows, can be easily removed by thresholding the estimated uncertainty.

Figure 6 shows a comparison of the traced inferior occipito-frontal fibers (IOFF) as traced by the 1-fiber (red) and 2-fiber (green) NODDI model. The 2-fiber NODDI manages to trace and connect a different part of gray matter region that is missed by the 1-fiber NODDI model. In fact the 1-fiber tracts only trace the medial portion of the lateral occipital cortex, whereas the 2-fiber NODDI tracts cover most of the lateral-occipital cortex as labeled by Freesurfer. Thus, the 2-fiber NODDI potentially provides a better estimate of fiber connectivity.

FIGURE 6

Figure 6. Comparision of results between 1-fiber NODDI model tractography(red) and 2-fiber NODDI model tractography(green) for the IOFF fibers.

4. Discussion and Limitations

In this work, we applied a new computationally efficient and numerically robust unscented information filtering framework for joint estimation of NODDI parameters and tractography. The proposed UIF filter provides a 2-fold improvement in processing speed, while computing the uncertainty in the estimated parameters, which is generally lacking in most existing methods. The method allows for tracing crossing fibers while accounting for the correlation in diffusion along the tracts. To the best of our knowledge, this is the first time, multi-fiber NODDI model has been used to estimate fiber dispersion index and intracellular-extracellular volume fractions for each of the crossing fibers separately, within a tractography framework. This can be quite useful in understanding the geometric properties of the white matter as it is traced along the tracts. The proposed algorithm is an open-source software and can be downloaded from https://github.com/pnlbwh/ukftractography, with the option of using the NODDI model as one of the possible choices.

We should however note that, the fiber dispersion dispersion (OD) index computed using the NODDI model is different than the one obtained using the methods in Savadjiev et al. (2010, 2012). In particular, the latter are computed from the fiber tracts or tensor fields at a macroscopic level, whereas the fiber dispersion obtained from NODDI is inherently microscopic. While there could be regions where they agree (some regions of white matter), yet in the gray matter the OD measures from NODDI is very different than the fiber dispersion computed using the method in Savadjiev et al. (2012).

Nevertheless, the proposed method has a few limitations. First, we assume equal volume fraction for each of the crossing fibers in our 2-fiber NODDI model, which may not be accurate. Second, the model fit to the data in the CSF areas is poor due to high noise in the data itself, as seen in Figure 7. However, the error in most white and gray matter areas is quite low, i.e., less than 2%. Another limitation of the current implementation of the proposed method is its inability to trace more than 2 fiber crossings. While an extension to trace 3 fiber model is straightforward, it can be done in areas which are a-priori known to have 3 fiber crossings. However, we believe that using a 3-fiber model for tracing all fibers would result in over-fitting of the data. Yet, we should note that the a majority of the white matter voxels have two crossings, and a very small region has more than 2 fiber crossings. Thus, the proposed method can be applied in most parts of the brain to trace fibers and estimate the dispersion index of each fiber separately, which is a significant improvement over the existing single-fiber based model of Zhang et al. (2012).

FIGURE 7

Figure 7. Arcuate Fasciculus traced using 1-fiber NODDI model. Background is the b = 0 image, where very bright regions indicate CSF areas. As seen, higher error in data fitting occurs only in the CSF areas, which is extremely noisy with free isotropic diffusion.

Another limitation, as is the case with most tractography algorithms, is the directional nature of the estimated tracts. For example, the tracts seeded in region A and reaching B, may not be obtained if seeding was done in region B. However, a typical way this particular problem is addressed is by seeding the whole brain and extracting tracts of interest from the whole brain tractography as has been done here and in most works using deterministic tractography.

Nevertheless, we believe that the present method allows to estimate parameters of the NODDI model along fiber tracts and allows to trace fibers through crossing regions. This could be useful in neuroscience studies to detect changes in white matter structure due to disease.

Author Contributions

YR: Conception, Design, Manuscript writing and data analysis. PR: Implementation, paper writing, data analysis.

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.

Acknowledgments

This project was funded by the National Institute of Health, R01MH097979 (PI: YR).

References

Alexander, A., Hasan, K., Tsuruda, J., and Parker, D. (2001). Analysis of partial volume effects in diffusion-tensor MRI. Magn. Reson. Med. 45, 770–780. doi: 10.1002/mrm.1105

PubMed Abstract | CrossRef Full Text | Google Scholar

Alexander, D., Barker, G., and Arridge, S. (2002). Detection and modeling of non-Gaussian apparent diffusion coefficient profiles in human brain data. Magn. Reson. Med. 48, 331–340. doi: 10.1002/mrm.10209

PubMed Abstract | CrossRef Full Text | Google Scholar

Anderson, A. (2005). Measurement of fiber orientation distributions using high angular resolution diffusion imaging. Magn. Reson. Med. 54, 1194–1206. doi: 10.1002/mrm.20667

PubMed Abstract | CrossRef Full Text | Google Scholar

Assemlal, H.-E., Tschumperlé, D., Brun, L., and Siddiqi, K. (2011). Recent advances in diffusion MRI modeling: angular and radial reconstruction. Med. Image Anal. 15, 369–396. doi: 10.1016/j.media.2011.02.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Basser, P., Mattiello, J., and LeBihan, D. (1994). MR diffusion tensor spectroscopy and imaging. Biophys. J. 66, 259–265. doi: 10.1016/s0006-3495(94)80775-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Baumgartner, C., Michailovich, O., Levitt, J., Pasternak, O., Bouix, S., Westin, C., et al. (2012). “A unified tractography framework for comparing diffusion models on clinical scans,” in Computational Diffusion MRI Workshop of MICCAI, Nice, 27–32.

Google Scholar

Daducci, A., Canales-Rodríguez, E. J., Zhang, H., Dyrby, T. B., Alexander, D. C., and Thiran, J.-P. (2015). Accelerated microstructure imaging via convex optimization (amico) from diffusion MRI data. NeuroImage 105, 32–44. doi: 10.1016/j.neuroimage.2014.10.026

PubMed Abstract | CrossRef Full Text | Google Scholar

Descoteaux, M., Angelino, E., Fitzgibbons, S., and Deriche, R. (2007). Regularized, fast, and robust analytical Q-ball imaging. Magn. Reson. Med. 58, 497–510. doi: 10.1002/mrm.21277

PubMed Abstract | CrossRef Full Text | Google Scholar

Frank, L. (2002). Characterization of anisotropy in high angular resolution diffusion-weighted MRI. Magn. Reson. Med. 47, 1083–1099. doi: 10.1002/mrm.10156

PubMed Abstract | CrossRef Full Text | Google Scholar

Hess, C., Mukherjee, P., Han, E., Xu, D., and Vigneron, D. (2006). Q-ball reconstruction of multimodal fiber orientations using the spherical harmonic basis. Magn. Reson. Med. 56, 104–117. doi: 10.1002/mrm.20931

PubMed Abstract | CrossRef Full Text | Google Scholar

Jansons, K., and Alexander, D. (2003). Persistent angular structure: new insights from diffusion MRI data. Inverse Prob. 19, 1031–1046. doi: 10.1088/0266-5611/19/5/303

CrossRef Full Text | Google Scholar

Jian, B., and Vemuri, B. (2007). A unified computational framework for deconvolution to reconstruct multiple fibers from diffusion weighted MRI. Trans. Med. Imaging 26, 1464–1471. doi: 10.1109/TMI.2007.907552

PubMed Abstract | CrossRef Full Text | Google Scholar

Kaden, E., Knøsche, T., and Anwander, A. (2007). Parametric spherical deconvolution: inferring anatomical connectivity using diffusion MR imaging. NeuroImage 37, 474–488. doi: 10.1016/j.neuroimage.2007.05.012

PubMed Abstract | CrossRef Full Text | Google Scholar

Kreher, B., Schneider, J., Mader, I., Martin, E., Hennig, J., and Il'yasov, K. (2005). Multitensor approach for analysis and tracking of complex fiber configurations. Magn. Reson. Med. 54, 1216–1225. doi: 10.1002/mrm.20670

PubMed Abstract | CrossRef Full Text | Google Scholar

Lee, D.-J. (2008). Nonlinear estimation and multiple sensor fusion using unscented information filtering. IEEE Sig. Process. Lett. 15, 861–864. doi: 10.1109/LSP.2008.2005447

CrossRef Full Text | Google Scholar

Malcolm, J., Shenton, M., and Rathi, Y. (2009a). “Neural tractography using an unscented Kalman filter,” in IPMI, 126–138.

PubMed Abstract | Google Scholar

Malcolm, J. G., Shenton, M. E., and Rathi, Y. (2009b). “Two-tensor tractography using a constrained filter,” in Medical Image Computing and Computer Assisted Intervention (MICCAI), 894–902.

PubMed Abstract | Google Scholar

Malcolm, J. G., Shenton, M. E., and Rathi, Y. (2010). Filtered multi-tensor tractography. IEEE Trans. Med. Imaging 29, 1664–1675. doi: 10.1109/TMI.2010.2048121

PubMed Abstract | CrossRef Full Text | Google Scholar

Ozarslan, E., Shepherd, T., Vemuri, B., Blackband, S., and Mareci, T. (2006). Resolution of complex tissue microarchitecture using the diffusion orientation transform (DOT). NeuroImage 31, 1086–1103. doi: 10.1016/j.neuroimage.2006.01.024

PubMed Abstract | CrossRef Full Text | Google Scholar

Parker, G., and Alexander, D. (2005). Probabilistic anatomical connectivity derived from the microscopic persistent angular structure of cerebral tissue. Philos. Trans. R. Soc. B 360, 893–902. doi: 10.1098/rstb.2005.1639

PubMed Abstract | CrossRef Full Text | Google Scholar

Peled, S., Friman, O., Jolesz, F., and Westin, C.-F. (2006). Geometrically constrained two-tensor model for crossing tracts in DWI. Magn. Reson. Med. 24, 1263–1270. doi: 10.1016/j.mri.2006.07.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Rathi, Y., Michailovich, O., Bouix, S., and Shenton, M. (2008). “Directional functions for orientation distribution estimation,” in International Symposium on Biomedical Imaging, 927–930. doi: 10.1109/isbi.2008.4541149

PubMed Abstract | CrossRef Full Text | Google Scholar

Savadjiev, P., Kindlmann, G. L., Bouix, S., Shenton, M. E., and Westin, C.-F. (2010). Local white matter geometry from diffusion tensor gradients. NeuroImage 49, 3175–3186. doi: 10.1016/j.neuroimage.2009.10.073

PubMed Abstract | CrossRef Full Text | Google Scholar

Savadjiev, P., Rathi, Y., Bouix, S., Verma, R., and Westin, C.-F. (2012). “Multi-scale characterization of white matter tract geometry,” in Medical Image Computing and Computer-Assisted Intervention–MICCAI 2012 (Springer), 34–41. doi: 10.1007/978-3-642-33454-2_5

PubMed Abstract | CrossRef Full Text

Simon, D., and Simon, D. (2006). Kalman filtering with inequality constraints for turbofan engine health estimation. IEE Proc. Control Theory Appl. 153, 371–378. doi: 10.1049/ip-cta:20050074

CrossRef Full Text | Google Scholar

Tournier, J.-D., Calamante, F., Gadian, D., and Connelly, A. (2004). Direct estimation of the fiber orientation density function from diffusion-weighted MRI data using spherical deconvolution. NeuroImage 23, 1176–1185. doi: 10.1016/j.neuroimage.2004.07.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Tuch, D. (2004). Q-ball imaging. Magn. Reson. Med. 52, 1358–1372. doi: 10.1002/mrm.20279

PubMed Abstract | CrossRef Full Text | Google Scholar

Tuch, D., Reese, T., Wiegell, M., Makris, N., Belliveau, J., and Wedeen, V. (2002). High angular resolution diffusion imaging reveals intravoxel white matter fiber heterogeneity. Magn. Reson. Med. 48, 577–582. doi: 10.1002/mrm.10268

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Essen, D. C., Smith, S. M., Barch, D. M., Behrens, T. E., Yacoub, E., Ugurbil, K., et al. (2013). The wu-minn human connectome project: an overview. Neuroimage 80, 62–79. doi: 10.1016/j.neuroimage.2013.05.041

PubMed Abstract | CrossRef Full Text | Google Scholar

Wan, E. A., and Van Der Merwe, R. (2000). “The unscented kalman filter for nonlinear estimation,” in Adaptive Systems for Signal Processing, Communications, and Control Symposium 2000. AS-SPCC. The IEEE 2000 (IEEE), 153–158. doi: 10.1109/asspcc.2000.882463

CrossRef Full Text | Google Scholar

Wassermann, D., Makris, N., Rathi, Y., Shenton, M., Kikinis, R., Kubicki, M., et al. (2013). “On describing human white matter anatomy: the white matter query language,” in Medical Image Computing and Computer-Assisted Intervention–MICCAI 2013 (Springer), 647–654. doi: 10.1007/978-3-642-40811-3_81

PubMed Abstract | CrossRef Full Text

Zhang, F., Goodlett, C., Hancock, E., and Gerig, G. (2007). “Probabilistic fiber tracking using particle filtering,” in MICCAI, 144–152. doi: 10.1007/978-3-540-75759-7_18

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, H., Schneider, T., Wheeler-Kingshott, C. A., and Alexander, D. C. (2012). Noddi: practical in vivo neurite orientation dispersion and density imaging of the human brain. Neuroimage 61, 1000–1016. doi: 10.1016/j.neuroimage.2012.03.072

PubMed Abstract | CrossRef Full Text | Google Scholar

Appendix

The Unscented Information Filter (UIF)

The state x_i of the UIF is given the unknown parameters of the single or multi-fiber NODDI model described in Equation (3). The initial state x₀ is obtained using a coarse grid search for the unknown parameters as described in Zhang et al. (2012). The initial covariance matrix $P_0\in {ℝ}^{n\times n}$ is assumed to be diagonal with all entries initially set to 0.01. Here n is the dimension of the state. There are a total of 2n+1 sigma points, with weights $w_0=\frac{k}{(n+k)}$, $w_i=\frac{1}{2(n+k)}$, and ζ = n+k. k is an adjustable scaling parameter set to 0.01 in all our experiments (Wan and Van Der Merwe, 2000).

1. Form weighted sigma points ${\text{X}}_t={\left\{w_i,{\text{x}}_i\right\}}_{i=0}^{2n}$ around current mean x_t and covariance P_t with scaling factor ζ

$\begin{array}{ll}{\text{X}}_0={\text{x}}_t{\text{X}}_i={\text{x}}_t+{\left[\sqrt{\zeta P_t}\right]}_i{\text{X}}_{i+n}={\text{x}}_t-{\left[\sqrt{\zeta P_t}\right]}_i& (\text{A}1)\end{array}$

2. Prediction Equations

$\begin{array}{rc}{\text{X}}_{t+1|t}^x=s\left[{\text{X}}_t\right]& (\text{A}2)\end{array}$

$\begin{array}{rc}{\overline{\text{x}}}_{t+1|t}={\displaystyle \sum _i}{w}_i{\text{X}}_{i,t+1|t}^x& (\text{A}3)\end{array}$

Predicted State covariance matrix

$\begin{array}{ll}{P}_{t+1|t}={\displaystyle \sum _i}{w}_i({\text{X}}_{i,t+1|t}^x-{\overline{\text{x}}}_{t+1|t}){({\text{X}}_{i,t+1|t}^x-{\overline{\text{x}}}_{t+1|t})}^T& (\text{A}4)\end{array}$

Information prediction equations

$\begin{array}{rc}{\text{Y}}_{t+1|t}={(P_{t+1|t})}^{-1}& (\text{A}5)\end{array}$

$\begin{array}{rc}{\overline{\text{y}}}_{t+1|t}={\text{Y}}_{t+1|t}{\overline{\text{x}}}_{t+1|t}& (\text{A}6)\end{array}$

3. Predict new sigma points for update, e.g.,

$\begin{array}{l}{\text{X}}_{0,t+1|t}={\overline{\text{x}}}_{t+1|t}{\text{X}}_{i,t+1|t}={\overline{\text{x}}}_{t+1|t}+{\left[\sqrt{\zeta P_{t+1|t}}\right]}_i\\ {\text{X}}_{i+n,t+1|t}={\overline{\text{x}}}_{t+1|t}-{\left[\sqrt{\zeta P_{t+1|t}}\right]}_i\hfill & (\text{A}7)\end{array}$

4. Measurement update equations for y_t are

$\begin{array}{rc}{\text{Z}}_{t+1|t}=h\left[{\text{X}}_{t+1|t}\right]& (\text{A}8)\end{array}$

$\begin{array}{rc}{\overline{\text{z}}}_{t+1|t}={\displaystyle \sum _i}{w}_i{\text{X}}_{i,t+1|t}& (\text{A}9)\end{array}$

$\begin{array}{rc}{P}_{xz}={\displaystyle \sum _i}{w}_i({\text{X}}_{i,t+1|t}-{\overline{\text{x}}}_{t+1|t}){({\text{Z}}_{i,t+1|t}-{\overline{\text{z}}}_{t+1|t})}^T& (\text{A}10)\end{array}$

$\begin{array}{rc}{H}_{t+1}={\text{Y}}_{t+1|t}{P}_{xz}& (\text{A}11)\end{array}$

$\begin{array}{rc}{i}_{t+1}=H_{t+1}^T{R}_{t+1}^{-1}{H}_{t+1}& (\text{A}12)\end{array}$

$\begin{array}{rc}{\text{I}}_{t+1}=H_{t+1}^T{R}_{t+1}^{-1}(z-{\text{z}}_{t+1|t}+{(P_{xz})}^T{\overline{\text{y}}}_{t+1|t})& (\text{A}13)\end{array}$

$\begin{array}{r}{P}_{t+1}={({\text{Y}}_{t+1|t}+{\text{I}}_{t+1})}^{-1}\\ {\text{X}}_{t+1}=P_{t+1}(i_{t+1}+{\overline{\text{y}}}_{t+1|t})& (\text{A}14)\end{array}$

Note that, the matrix I_t above, is the Fisher Information matrix and can be potentially used to determine the amount of information that the observed measurement z carries about each of the model parameters. In the above formulation, we have ignored the constraints on our model. This could result in non-realistic estimates, i.e., the volume fractions may not lie between 0 and 1. To enforce appropriate constraints, one can directly project any unconstrained state x onto the constrained subspace as shown in Simon and Simon (2006). In other words, we wish to find the state $\widehat{\text{x}}$ closest to the unconstrained state x which still obeys the constraints, $C\widehat{\text{x}}\le \text{d}$. For example, for the single-fiber NODDI model, the state is $\text{x}=\left[V_{ic}\kappa \text{m}{V}_{iso}\right]\in {ℝ}^{6\times 1}$. In this case, one can define the matrix C and d by:

$\begin{array}{ll}C=\left[\begin{array}{c}1000\\ -1000\\ 0-100\\ 0001\\ 000-1\end{array}\right],d=\left[\begin{array}{c}1\\ 0\\ 0\\ 1\\ 0\end{array}\right].& (\text{A}15)\end{array}$

Using P_t as a weighting matrix, this becomes a quadratic programming problem:

$\begin{array}{ll}{\displaystyle \underset{\widehat{\text{x}}}{min}}{(\text{x}-\widehat{\text{x}})}^T{P}_t^{-1}(\text{x}-\widehat{\text{x}})\text{\hspace{1em}subject to\hspace{1em}}C\widehat{\text{x}}\le \text{d}.& (\text{A}16)\end{array}$

This projection procedure is applied within the UIF algorithm to correct at every place where we move in the state-space: after spreading the sigma points X_t, after transforming the sigma points X_t+1|t, and after the final estimate x_t+1.