Multifiber pathway reconstruction using bundle constrained streamline

Computerized Medical Imaging and Graphics 46 (2015) 291–299

Contents lists available at ScienceDirect

Computerized Medical Imaging and Graphics journal homepage: www.elsevier.com/locate/compmedimag

Multifiber pathway reconstruction using bundle constrained streamline Chun-Yu Chu a,b,∗ , Jian-Ping Huang a,b , Chang-Yu Sun b , Yan-Li Zhang a , Wan-Yu Liu a,b , Yue-Min Zhu a,b a b

Metislab, Harbin Institute of Technology, Harbin, China CREATIS, CNRS UMR 5220, Inserm U1044, INSA Lyon, University of Lyon, Villeurbanne, France

a r t i c l e

i n f o

Article history: Received 12 January 2015 Received in revised form 19 June 2015 Accepted 28 July 2015 JEL classification: 70: 3-D reconstruction

a b s t r a c t Fiber tractography techniques in diffusion magnetic resonance imaging have become a primary tool for studying the fiber architecture of biological tissues both noninvasively and in vivo. Streamline tracking, as a simple and efficient tractography technique, is widely used to reconstruct fiber pathways. It is however very sensitive to noisy estimation of local fiber orientations. In this paper, we propose a bundle constrained streamline method to accurately reconstruct multifiber pathways. The method introduces neighboring fiber consistency constraint in the tracking process and reconstructs fiber pathways that have optimal tradeoff between consistency with local fiber orientation estimations and similarity with neighboring fiber segment orientations. Results on synthetic, physical phantom and real human brain DW images show that the proposed method allows regular fiber pathways to be reconstructed and outperforms existing techniques. © 2015 Elsevier Ltd. All rights reserved.

1. Introduction Diffusion magnetic resonance imaging (dMRI) tractography provides a primary tool for studying the fiber architecture of biological tissues both noninvasively and in vivo. The technique allows us to estimate neuronal tract pathways to determine both the structure and the connectivity of the brain. Diffusion tensor imaging (DTI) based streamline tracking was the first technique proposed and has been largely used to reconstruct fiber pathways in the human brain [1]. However, the DTI has been shown to fail in regions presenting intravoxel orientational heterogeneity such as fiber or fiber bundle crossing and branching. A number of more sophisticated models for estimating multiple fiber orientations (MFOs) within a voxel have been proposed in order to overcome the single orientation limitation of DTI. We cite, for instance, multitensor based models [2–4], high-order diffusion tensor [5], Q-Ball [6,7], diffusion spectrum imaging [8], spherical deconvolution [9–11], independent component analysis (ICA) from DTI [12], etc. More details can be found in [13,14].

∗ Corresponding author at: Metislab, Harbin Institute of Technology, 92 Xidazhi Street, Harbin, China. Tel.: +86 451 86403245. E-mail addresses: [email protected] (C.-Y. Chu), [email protected] (W.-Y. Liu). http://dx.doi.org/10.1016/j.compmedimag.2015.07.010 0895-6111/© 2015 Elsevier Ltd. All rights reserved.

To reconstruct multifiber pathways form the MFO system, several extensions of DTI-based streamline have been proposed for white matter tractography. A modified fiber assignment using the continuous tracking method was proposed to reconstruct tracts along the MFOs obtained from Q-Ball imaging [15]. Two-tensor based streamline tractography was successfully obtained in the regions of crossing fibers [16]. ICA based streamline was also developed to conduct human brain tractography from diffusion data [17]. In the above-mentioned streamline tracking techniques, fiber pathways are constructed based on local fiber orientation estimation in a step-by-step way. Therefore, the accuracy of the reconstruction is sensitive to the estimation of intravoxel fiber orientations and suffers from cumulative errors along the tracking path. Recently, global optimization based methods have been reported to reconstruct the most probable paths in the whole brain. These methods estimate fiber pathways by a cost function minimization or posterior probability maximization process. Some of the global methods consider the optimization of the path between two designated regions to produce optimal trajectories [18–21]. Others reconstruct full fiber tracts over the whole brain at once [22–25]. The global methods are less affected by the noisy estimation of local fiber orientations and the cumulative errors. However, they are generally computationally intensive and may not satisfy the requirement of practical uses.

292

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

Initialization Streamline Tracking Initialized Fiber Pathways Iteration BCS Based Fiber Pathway Reconstruction Fiber Pathway Smoothing

Fig. 2. The determination of fiber segment direction in BCS tracking.

Optimal Fiber Pathways

Fig. 1. Block diagram of the BCS method.

In this paper, we propose a bundle constrained streamline (BCS) method which introduces the orientation information of neighboring fibers in the process of the streamline tracking to overcome the sensitivity to the local uncertainty in the estimation of intravoxel fiber orientations. Based on the assumption of local fiber consistency (i.e. the neighboring fiber pathways are similar in local orientation), the method reconstructs fiber pathways that have optimal tradeoff between consistency with local fiber orientation estimations and similarity with neighboring fiber segment orientations. 2. Methodology The proposed fiber pathway reconstruction method is depicted in Fig. 1. It consists essentially of two main elements: fiber pathway reconstruction according to BCS tracking rules and fiber pathway smoothing. 2.1. Streamline tracking The streamline tracking method reconstructs fiber pathways in a step-by-step way. The process starts with a predefined point in the region containing fibers, called the seed point, and iteratively moves the point along local diffusion orientations until a stopping criterion is fulfilled. In general, the stopping criteria include the exiting from the region of interest, intersection at voxels having low anisotropy values and large local curvature of the trajectory. The iterative tracking process can be descripted by the updating formula:

tracking, we determine the next fiber segment position pi+1 according to the local fiber orientation estimation as well as the mean fiber orientation of all or part of the fiber sets , where each fiber passes through the neighborhood of seed point p0 with radius RD . The determination of the next fiber segment position pi+1 is illustrated in Fig. 2. Assuming that multiple fiber bundles probably exist in , we divide  into one or more fiber sets denoted by 1 , 2 , . . ., NC , where NC is the number of sets. In the present study, to divide , we use a simple agglomerative hierarchical clustering, where a distance threshold HTC between two subnodes of hierarchical clustering tree is used as the criterion for forming clusters. The distance metric between two fibers is defined as the arithmetic mean of the distances between the corresponding points on two fibers. In Fig. 2, the curved surface Sk is a part of the sphere at pi with radius st (step size), and any point q ∈ Sk satisfies the angle between q–pi and di+1 (i.e. (pi − pi−1 )/st ) is smaller than an angle threshold  s . We determine the next fiber segment position pi+1 as the solution of the optimization problem: pi+1 = arg min [ED (q) + EF (q)] ,

(2)

q∈Sk

where  is the control parameter which provides a tradeoff between “diffusion potential energy” ED (q) and “flow potential energy” EF (q). The “diffusion potential energy” reflects the dissimilarity between the local fiber orientation estimate dD (main diffusion direction) and the fiber segment direction di = (q − pi )/st at pi , and is defined as



ED (q) = 1 − dTD di

2

.

(3)

where st is the step size, pi+1 the next position determined by the previous position pi , and di the local diffusion orientation at ith step. Due to the influence of noise in dMRI data on local diffusion orientations, the streamline method often generates inaccurate fiber reconstruction.

Note that multiple fiber orientations at pi are possible. In this case, dD is the nearest orientation to di . Especially, if there is no fiber orientation at pi (invalid voxel or with very low fractional anisotropy value), we set ED (q) = ∞. The “flow potential energy” measures the dissimilarity between the fiber segment direction di and the mean flow directions dFj of fiber sets j (j = 1, 2, . . ., NC ) at pi , and is defined as

2.2. BCS tracking rules

EF (q) = 1 − dTi dF

Like the streamline, the BCS tracking also reconstructs fiber pathways in a step-by-step way. However, in each step of BCS

where dF is one of dFj which is the nearest to di . Especially, if there is no valid dFj , we set EF (q) = ∞.

pi+1 = pi + st di

(1)



2

,

(4)

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

Fig. 3. Computation of the mean orientation of a fiber bundle.

The computation of the mean flow direction dFj is shown in Fig. 3. The plane Sp passes through pi and perpendicular to di–1 ; the kth fiber in j passes through Sp at point pjk , and pik is the fiber segment direction of the kth fiber at pjk . We compute dFj as

Nj dFj =

    pjk − pi  , arccos dTjk di−1  d jk   ,   fw pjk − pi  , arccos dTjk di−1 

f k=1 w

Nj

k=1

(5)

293

Algorithm 1 BCS tracking Step size st , parameters RD , HTC , ,  s , Rd , smoothing radius rs , Input: number of iterations Niter , seed points Pseeds , and MFO field. Output: Fiber pathway set . 1: ←Streamline tracking For it = 1 to Niter do 2: For ∀F ∈ ˚ do 3: 4: Remove F from  5: P0 ←The seed point of F 6: ←Find neighboring fibers of p0 with radius Rd 7: {1 , 2 , . . ., NC } ←Agglomerative hierarchical clustering i←0 8: 9: Repeat 10: Sk ← A part of the sphere at pi with radius st (Fig. 2) 11: Uniformly select 100 points Q = {q1 , q2 , . . ., q100 } on Sk 12: For ∀q ∈ Q do 13: Compute energy ED (q) using (2) 14: Compute energy EF (q) using (3) 15: End For 16: pi+1 ← arg min [ED (q) + EF (q)] q∈Q

17: 18: 19: 20: 21: 22: 23:

i←i+1 Until ED (q) + EF (q) = / ∞, ∀q ∈ Q F ← {p0 , p1 , . . ., pi } Smoothing F using (6) Add F to  End For End For

where Nj designates the number of fibers in j and fw is a weighting function defined by

⎧ ⎨ 1 − x , 0 ≤ x ≤ R , y ≤ s d Rd , fw (x, y) = ⎩ 0, otherwise

2.5. Determination of parameters (6)

where x is the distance between two points, y is the angle between two fiber segments, and Rd and  s are the parameters allowing omitting the fibers that obviously deviate from pi . In practice, to reduce the computational complexity of (2), we uniformly select 100 points on Sk , and choose the point having the minimum energy as the solution of (1). The stop criterion of BCS tracking is ED (q) + EF (q) = ∞, ∀q ∈ Sk .

2.3. Fiber pathways smoothing We smooth fiber pathways using simple neighborhood averaging method. To a fiber pathway F are linked a series of points p1 , p2 , . . ., ps , . . ., pn , where ps denotes the seed point of F. Each point is smoothed with odd span rs as

pi =

⎧ ⎪ ⎨



i+(rs −1)/2

⎪ ⎩ j=i−(rs −1)/2 ps ,

pj rs

, i= / s

.

(7)

i=s

The position of seed point is preserved for the subsequent analysis.

2.4. Procedure of BCS tracking The BCS tracking procedure is described as follows:

The proposed method involves several parameters, so we briefly discuss the impact of the parameters in this section. We suggest that some parameters are chosen according to the smallest spacing between two neighboring voxels sm . Empirically, we give our default choice of these parameters used in the present work. The step size st should be reasonably small. However, if st is too small, it will lead to increased computational workload. As a default, we choose st = 0.5sm . The radius RD is used to determine the fiber set  around a seed point. The parameter RD should be reasonably large such that  contains a sufficient number of fibers. An excessively large RD will not lead to bad results but will increase numerical complexity. Our default choice is RD = 2sm . The parameter HTC is a distance threshold used to divide . A too small HTC may lead to a fiber set  that will be divided into too many clusters, while a too large HTC may lead to one cluster containing several fiber bundles diverging in different directions. In both the extreme cases, the reconstruction error will become relatively higher. In our experiments, we set HTC = 30 mm. The parameter  is an important parameter which ensures the tradeoff between diffusion potential and flow potential. The choice of this parameter should take into account the noise level of data. In the case of high level noise, a large  would be advised because of the low reliability of the local fiber orientation estimations. As a default, we choose  = 6. The parameter  s defining the fiber segment direction is limited to a reasonable extent. Our default choice is  s = 30 degrees. The parameter Rd combined with  s determines which fibers in j are to be used to compute the mean flow direction. As a default, we choose Rd = sm degrees. The radius rs controls the level of fiber smoothing. Our default choice is rs = 3. For the number of iterations Niter , a choice of 3 to 5 is generally sufficient. A too large Niter will lead to unnecessary numerical complexity. Our default choice is Niter = 5.

294

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

Fig. 4. The fiber pathways of the simulated data: (a) ground truth, (b) reconstructed using the streamline method, (c) reconstructed using the BCS method. (d) Zoom of the red box region in (b). The red, black, and yellow ellipses in (d) respectively indicate the premature termination, divagation, and inflection of the fiber pathways. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

3. Experiments and results The proposed method is evaluated on synthetic, physical phantom and real human brain data. In all our experiments, the local fiber orientation is estimated using the constrained compressed sensing (CCS) method [3]. As a default, we set the fractional anisotropy threshold as 0.1 and angle threshold as 30 degrees for the streamline method.

The seed points were marked with the small red circles and the intravoxel fiber structures were shown as the background. It can be seen that, at the SNR of 17, the streamline method is very sensitive to the local fiber orientation estimation, which results in some fibers showing premature termination (marked by the red ellipse in Fig. 4d), divagation (marked by the black ellipse in Fig. 4d) or inflection (marked by the yellow ellipses in Fig. 4d). However, the fibers from the BCS method are well consistent with the ground truth, without premature termination, divagation, or inflection at all.

3.1. Synthetic data The diffusion-weighted (DW) data of size 16 × 16 were synthesized from the multitensor model with b-value of 1000 s/mm2 and 30 gradient directions. The diffusion tensor eigenvalues were set to (1.2, 0.3, 0.3) × 10−3 mm2 /s. The fiber architecture of the simulated data (noise free) is shown in Fig. 4a, which consists of two fiber bundles form a crossing configuration. Rician noise is added to generate the noisy data with a signal-to-noise ratio (SNR) of 17. In the present paper, the SNR of the simulation data was defined as the ratio of nondiffusion-weighted signal intensity to the noise standard deviation. In the simulation experiment, we set the parameters rs = 2 and Niter = 3 for the BCS method. Fig. 4b and c shows the fiber pathways reconstructed using the streamline and BCS methods, respectively.

3.2. Physical phantom data We used the diffusion phantom data published for the Fiber Cup at MICCAI 2009 [26], specifically the data with b-value of 1500 s/mm2 , 64 gradient directions, imaging matrix of 64 × 64 × 3, and isotropic spatial resolution of 3 mm. The phantom data were designed to compare the performance of fiber tractography algorithms. It simulates a coronal section of the human brain, containing several fiber crossing and kissing configurations with different curvatures. In the Fiber Cup challenge, a total of 16 seed points and the corresponding 16 fiber pathways (see Fig. 7a) were supplied as the ground truth to compare the performance of different methods.

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

295

Table 1 The empirical parameters, range of optimal parameters and that used in the present work. Parameter

RD (mm)

Rd (mm)

 s (degree)

HTC (mm)

Niter



Empirical parameters Range of optimal value Used in our work

12 ≥6 6

3 ≥3 3

30 30–40 30

24 21–30 27

5 ≥3 5

10 6–8 8

Three metrics were used for the contest: spatial metric, tangent metric, and curve metric [27]. Some seed points are selected manually and marked as the small red circles in Fig. 6. In particular, the selected seed points contain the 16 seed points initially designed for the contest. A mask has been defined to exclude the background region for

better visualization and reducing computation workload. In all the experiments on the phantom data, we set the fractional anisotropy threshold to 0.01 for the streamline method. The parameters of the BCS method were optimized in terms of the three metrics used for the contest. The optimization process starts with a set of empirical parameters (the first row in Table 1), and varies one

Fig. 5. The impact of main parameters on the 16 fiber pathways reconstruction, in terms of the spatial metric, tangent metric and curve metric.

296

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

Fig. 6. Fiber tracking results of the phantom data: (a) streamline method, (b) BCS method.

Fig. 7. The fiber pathways of the phantom data from the 16 seed points: (a) the ground truth, (b) reconstructed using the winner’s method, (c) using BCS method.

Fig. 8. The fiber pathways of the real data: (a) reconstructed using streamline method, (b) using BCS method.

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

297

Fig. 9. Three regions of interest.

parameter at each time. Fig. 5 shows the impact of main parameters on the results. In some cases, we observe the curve metric inconsistent with the spatial metric and tangent metric. For example, in Fig. 5d the spatial metric and tangent metric indicate better performance when HTC = 21 to 30 mm, while the curve metric indicates the contrary results. Generally, the spatial metric reflects more adequately the similarities between pathways [28]. Table 1 summarizes the range of optimal parameters and that used in the present work. Fig. 6 shows the fiber tracking results obtained using the streamline and BCS methods with the selected seed points. By visually comparing the fiber tracking results, the BCS method is obviously superior to the streamline method. Fig. 7b and c shows the results of the winner’s method [29] reported in the Fiber Cup and our BCS method, respectively. It can be seen that the fiber pathways, reconstructed by both methods, grossly follow the ground truth. However, close inspections reveal that the winner’s method tends to generate fiber pathways that confuse the spatial position of adjacent fibers with different seed points, such as Fiber 3 and Fiber 4. Furthermore, the Fiber 12 was incorrectly reconstructed using the winner’s method, whereas the shapes of all the fiber pathways reconstructed using the BCS method were similar to the ground truth. Table 2 gives the quantitative results of the BCS method and the winner’s method in terms of the three metrics. The best values have been highlighted by a bold font. The scores indicate that the BCS method has better performance.

3.3. Real human brain data The real human brain data come from the open data collection at http://hdl.handle.net/1926/1687. The datasets were acquired on a 3 T GE system using an echo planar imaging (EPI) DTI Tensor sequence. Each DW volume consists of 86 contiguous axial slices of size 144 × 144. The acquisition parameters were: TR = 17,000 ms, TE = 78 ms, FOV = 24 cm, slice thickness = 1.7 mm, diffusion sensitivity factor b = 900 s/mm2 , and diffusion gradient directions = 51. In the real data experiment, we randomly select about 20 thousands seed points in the whole brain. Fig. 8 shows the fiber pathways reconstructed using the streamline and BCS methods, respectively. The fiber pathways are colored by their local orientation. The reconstructions contain approximately 13 thousands fiber pathways longer than 20 mm. The mean fiber-length of the two results is 57 mm and 76 mm, respectively. To analyze the details of the reconstructions, we define three regions of interest (Fig. 9) and extracted fiber bundles of the corpus callosum, the corticospinal tract, and the cingulate cortex, respectively. The MITK Diffusion toolkit [30] was used for the fiber extraction and visualization work. Fig. 10 shows the results of fiber extraction. For the corpus callosum, the streamline method has less fibers going to the apical cortex and the most lateral projections are missing, while our BCS method shows more widespread projections and covers larger cortical areas in both hemispheres. For the corticospinal tract and cingulate cortex reconstruction, our BCS method exhibits better local consistence and many more fibers than the streamline. Furthermore, the results reconstructed using

Table 2 Quantitative results of the BCS method on the Fiber Cup phantom compared to the contest winner’s quantitative results. The best values have been highlighted by a bold font. Fiber

Reisert et al. [23] Spatial metric (mm)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 Mean Mean* *

2.24 2.37 4.98 2.18 1.98 4.25 5.62 2.11 2.61 5.78 3.36 17.02 4.66 2.56 2.16 5.81 4.36 3.51

BCS method Tangent metric (degree)

Curve metric (mm−1 )

Spatial metric (mm)

Tangent metric (degree)

Curve metric (mm−1 )

9.33 12.38 6.35 5.37 6.20 8.10 11.44 8.26 6.73 12.10 4.85 46.74 12.75 14.74 4.25 7.54 11.07 8.69

0.025 0.040 0.021 0.014 0.020 0.026 0.023 0.043 0.017 0.026 0.013 0.033 0.070 0.069 0.011 0.021 0.030 0.029

2.35 2.49 1.36 2.20 1.45 1.17 2.83 4.75 3.62 3.41 2.86 1.98 1.21 4.71 3.26 3.86 2.72 2.77

7.13 6.29 5.84 5.69 7.12 6.09 9.94 13.01 7.17 11.26 7.83 7.77 7.38 15.30 9.48 6.73 8.38 8.42

0.019 0.023 0.020 0.025 0.025 0.033 0.041 0.023 0.022 0.025 0.023 0.025 0.036 0.047 0.032 0.022 0.028 0.028

The mean value was computed without using the Fiber 12.

298

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299

Fig. 10. The results of fiber extraction: (a) for the reconstruction of streamline, (b) BCS method. From top to bottom: fiber bundles of the corpus callosum, the corticospinal tract, and the cingulate cortex, respectively.

our BCS method were more consistent with known neuroanatomy in comparison with tracts from the streamline. 4. Discussion and conclusion We have proposed a BCS tracking method for reconstructing multifiber pathways based on multiple fiber orientations in a voxel. The method presents the particularity of achieving optimal tradeoff between consistency with local fiber orientation estimations and similarity with neighboring fiber segment orientation. The results on synthetic, physical phantom and real human brain data consistently demonstrated that the proposed BCS tracking method enables regular fiber pathways to be reconstructed and presents improved performance over the commonly used streamline and contest winner’s method. In comparison with the streamline method, the proposed BCS method does not exhibit premature termination, divagation, or inflection. The reason is that the streamline method tracks fiber pathways strictly based on the inaccurate MFO estimation caused by noise in the data or improper model fitting. In contrast, the proposed BCS method introduces the neighboring fiber consistency constraint to overcome the problems of premature termination and divagation and smooth the fiber pathways to diminish the inflection. The computational complexity of the BCS method depends on the fiber counts of initial reconstruction. Generally, more fibers mean more complexity of searching neighboring fibers and calculating the mean fiber orientation. In this case, we suggest searching neighboring fibers in a random subset of initial reconstruction as a compromise proposal. The BCS method is based on the assumption that fiber bundles are locally smooth. Therefore, at voxels where the local smoothness

assumption is obviously violated, such as trauma or pathology region, the neighboring fiber orientation information may lead to inaccurate reconstruction. In this case, the results obtained should be interpreted with precaution since fiber pathways could be altered. Conflict of interest statement None of the authors declare any conflicts of interest. Acknowledgments This work was supported by the National Natural Science Foundation of China (61271092), International S&T Cooperation Project of China, Ministry of Science and Technology of China (2007DFB30320), the Program PHC-Cai Yuanpei 2012, the Applied Technology Research and Development Program of Heilongjiang Province, Office of Science and Technology in Heilongjiang Province (GC13A311), the applied technology research and development program for Harbin Creative Talents (2014RFXXJ102), the Fundamental Research Funds for the Central Universities (HIT. NSRIF. 2015021), the Scientific Research Startup Fund for HIT Future Faculty, and the French ANR under ANR-13-MONU-0009-01.

References [1] Mori S, van Zijl PCM. Fiber tracking: principles and strategies—a technical review. NMR Biomed 2002;15(7–8):468–80. [2] Tuch DS, Reese TG, Wiegell MR, Makris N, Belliveau JW, Wedeen VJ. High angular resolution diffusion imaging reveals intravoxel white matter fiber heterogeneity. Magn Reson Med 2002;48(Oct (4)):577–82.

C.-Y. Chu et al. / Computerized Medical Imaging and Graphics 46 (2015) 291–299 [3] Landman B, Bogovic J, Wan H, Elshahaby FEZ, Bazin P-L, Prince JL. Resolution of crossing fibers with constrained compressed sensing using diffusion tensor MRI. NeuroImage 2012;59(Oct (3)):2175–86. [4] Ramirez-Manzanares A, Rivera M, Vemuri BC, Carney P, Mareci T. Diffusion basis functions decomposition for estimating white matter intravoxel fiber geometry. IEEE Trans Med Imaging 2007;26(8):1091–102. [5] Ozarslan E, Mareci TH. Generalized diffusion tensor imaging and analytical relationships between diffusion tensor imaging and high angular resolution diffusion imaging. Magn Reson Med 2003;50(Nov (5)):955–65. [6] Tuch DS. Q-ball imaging. Magn Reson Med 2004;52(Dec (6)):1358–72. [7] Descoteaux M, Angelino E, Fitzgibbons S, Deriche R. Regularized, fast, and robust analytical Q-ball imaging. Magn Reson Med 2007;58(Sep (3)):497–510. [8] Wedeen VJ, Hagmann P, Tseng W-YI, Reese TG, Weisskoff RM. Mapping complex tissue architecture with diffusion spectrum magnetic resonance imaging. Magn Reson Med 2005;54(Dec (6)):1377–86. [9] Tournier J-D, Calamante F, Connelly A. Robust determination of the fibre orientation distribution in diffusion MRI: non-negativity constrained super-resolved spherical deconvolution. NeuroImage 2007;35(May (4)):1459–72. [10] Tournier J-D, Calamante F, Gadian DG, Connelly A. Direct estimation of the fiber orientation density function from diffusion-weighted MRI data using spherical deconvolution. NeuroImage 2004;23(Nov (3)):1176–85. [11] Kaden E, Knösche TR, Anwander A. Parametric spherical deconvolution: inferring anatomical connectivity using diffusion MR imaging. NeuroImage 2007;37((Aug) 2):474–88. [12] Singh M. Estimation of multiple fiber orientations from diffusion tensor MRI using independent component analysis. IEEE Trans Nucl Sci 2005;52(Feb (1)):266–73. [13] Hasan KM, Walimuni IS, Abid H, Hahn KR. A review of diffusion tensor magnetic resonance imaging computational methods and software tools. Comput Biol Med 2011;41(12):1062–72. [14] Assemlal H-E, Tschumperlé D, Brun L, Siddiqi K. Recent advances in diffusion MRI modeling: angular and radial reconstruction. Med Image Anal 2011;15(Aug (4)):369–96. [15] Chao Y-P, Chen J-H, Cho K-H, Yeh C-H, Chou K-H, Lin C-P. A multiple streamline approach to high angular resolution diffusion tractography. Med Eng Phys 2008;30(Oct (8)):989–96. [16] a Qazi A, Radmanesh A, O’Donnell L, Kindlmann G, Peled S, Whalen S, et al. Resolving crossings in the corticospinal tract by two-tensor streamline

[17]

[18] [19]

[20]

[21]

[22]

[23] [24] [25]

[26]

[27]

[28] [29] [30]

299

tractography: method and clinical assessment using fMRI. NeuroImage 2009;47(Aug (Suppl. 2)):T98–106. Singh M, Wong C-W. Independent component analysis-based multifiber streamline tractography of the human brain. Magn Reson Med 2010;64(Dec (6)):1676–84. Jbabdi S, Woolrich MW, Andersson JLR, Behrens TEJ. A Bayesian framework for global tractography. NeuroImage 2007;37(Aug (1)):116–29. Wu X, Xu Q, Xu L, Zhou J, Anderson AW, Ding Z. Genetic white matter fiber tractography with global optimization. J Neurosci Methods 2009;184(Nov (2)):375–9. Wu X, Xie M, Zhou J, Anderson AW, Gore JC, Ding Z. Globally optimized fiber tracking and hierarchical clustering—a unified framework. Magn Reson Imaging 2012;30(May (4)):485–95. Frindel C, Robini M, Schaerer J, Croisille P, Zhu Y-M. A graph-based approach for automatic cardiac tractography. Magn Reson Med 2010;64(Oct (4)): 1215–29. Kreher BW, Mader I, Kiselev VG. Gibbs tracking: a novel approach for the reconstruction of neuronal pathways. Magn Reson Med 2008;60(Oct (4)): 953–63. Reisert M, Mader I, Anastasopoulos C, Weigel M, Schnell S, Kiselev V. Global fiber reconstruction becomes practical. Neuroimage 2011;54:955–62. Mangin J, Fillard P, Cointepas Y, Le Bihan D. Toward global tractography. NeuroImage 2013;80:290–6. Fillard P, Poupon C, Mangin J-F. A novel global tractography algorithm based on an adaptive spin glass model. International conference on medical image computing and computer-assisted intervention 2009;12(Pt 1):927–34. Poupon C, Rieul B, Kezele I, Perrin M, Poupon F, Mangin J-F. New diffusion phantoms dedicated to the study and validation of high-angular-resolution diffusion imaging (HARDI) models. Magn Reson Med 2008;60(Dec (6)):1276–83. Fillard P, Descoteaux M, Goh A, Gouttard S, Jeurissen B, Malcolm J, et al. Quantitative evaluation of 10 tractography algorithms on a realistic diffusion MR phantom. NeuroImage 2011;56(May (1)):220–34. Aranda R, Rivera M, Ramirez-Manzanares A. A flocking based method for brain tractography. Med Image Anal 2014;18(Apr (3)):515–30. Reisert M, Mader I, Anastasopoulos C, Weigel M, Schnell S, Kiselev V. Global fiber reconstruction becomes practical. NeuroImage 2011;54(Jan (2)):955–62. Fritzsche KH, Neher PF, Reicht I, van Bruggen T, Goch C, Reisert M, et al. MITK diffusion imaging. Methods Inf Med 2012;51(5):441–8.