1 Introduction

Diffusion magnetic resonance imaging (dMRI) is the only non-invasive tool to obtain information about the neural architecture within the white matter of the brain in vivo using the tractography process [1, 2]. It provides detailed information about white matter (WM) in many clinical applications like neurosurgical planning [3], post-surgery evaluations [4], and many aspects of assessment and study of neurological diseases such as multiple sclerosis [5, 6], Alzheimer’s disease [7], and schizophrenia [8].

Currently, fiber tracking is most commonly implemented using the principal diffusion directions of the diffusion tensor. All imaging systems, even clinical systems, now include a fiber tracking module available for clinicians and in common use among neuroscientists. However, an important limitation of the diffusion tensor model is the Gaussian diffusion assumption, which implies that there can only be a single-fiber population per voxel [9,10,11], posing significant problems for the diffusion tensor imaging (DTI) fiber tractography and the interpretation of DTI integrity metrics. The resolution of current typical diffusion-weighted (DW) data is much larger than the diameter of neural axons. Thus, there are multiple fiber tract orientations within a given voxel. Recent studies have demonstrated that more than 90% of white matter voxels contain more than one fiber population [12]. Thus, the diffusion tensor model is inadequate in the majority of white matter voxels as a result of partial volume effects.

A great variety of solutions have been proposed to tackle the problem of intravoxel fiber orientation estimation. The specificity of reconstructed tracts and diffusion features can be improved by achieving a high spatial resolution [13] or a high angular resolution [14]. However, there are practical limitations in increasing the angular resolution and the spatial resolution of the acquired data directly, such as a reduced signal-to-noise ratio (SNR) and a prolonged scanning time which limits the use of DW-MRI in clinical practice.

High angular resolution diffusion imaging (HARDI) [15] has become one of the most used methods to better characterize the complexity of water motion, thereby overcoming the limitations of the DTI technique. This technique is based on sampling the signal on a single or many spheres of the q-space. Among the HARDI techniques, there are the QBI (Q-Ball Imaging) method [16, 17] and the spherical deconvolution one [18, 19]. A common objective among these methods is the estimation of white matter fiber orientations within each voxel and the ability to study the brain connectivity. However, although able to identify complex fiber configurations, HARDI techniques require a very large number of DW images to estimate the diffusion signal when fiber tracts cross [20]. Increasing the angular resolution in HARDI technique is a challenging task because of the trade-off between the angular resolution and the acquisition time.

Using a reduced set of data combined with a post-processing method is considered as an alternative. This solution has been used quite successfully in other fields of medical research involving the cardiac medical image analysis and the electroencephalographic (EEG) domain. Wang et al. [45] used a novel meshless deformable model (MDM) to analyze the 3D cardiac motion and strain. The volumetric MDM represents the object with points cloud. The deformation is performed only by control points that are obtained by the intersections of three stacks of orthogonal tagging planes in the myocardium. The experimentation showed that the MDM based on only control points well recover the myocardium motion. Zhang et al. [44] applied the meshfree representation using element-free Galerkin method. This framework uses a limited number of nodes without explicit connectivity to present the entire analysis domain. It was validated in different cardiac motion applications. As a result, this framework provides a way to avoid the complicated meshing procedures while preserving the accuracy. J\(\ddot {\mathrm {a}}\)ger et al. [48] used the radial basis function (RBF) interpolation to derive 64-channel EEG data from only 32-channel recording. This approach permits a reconstruction of missing data which were lost due to technical problems caused by a broken electrode or physiological artifacts. The difference between the measured data and the predicted data is used to determine the errors measured. The experimentation demonstrates the advantages of the RBF interpolation in EEG reconstruction compared to the commonly used interpolation methods.

In HARDI domain, some works have proposed to use the compressive sensing (CS) technique by exploiting a sparse representation of the data to shorten the acquisition time. Michailovich et al. [21] used this technique to represent the HARDI signals in the basis of spherical ridgelets by using a relatively small number of representation coefficients. Sparse regularization has been also applied to accelerate the HARDI technique [23]. Tristn-Vega et al. [22] proposed to represent the probabilistic orientation distribution function (ODF) in the frame of spherical wavelets (SW). Many researchers focus on proposing novel sampling schemes to reduce the number of measurements used in the HARDI technique [24, 25]. Cho et al. [26] propose to use the hemispherical encoding scheme with the cross-term correction to reduce the scan time.

In this work, we propose to explore the angular interpolation approach in the gradient direction domain to construct the diffusion signal from a reduced number of gradient directions. These missing gradient directions will be obtained according to their neighborhood from a reduced set of data on the sphere of the q-space. The interpolation of diffusion MRI before fiber reconstruction has been used, but mainly for registering DW images [27]. It is, also, applied to improve the diffusion MRI spatial resolution and, therefore, to obtain a higher image resolution before fiber reconstruction using a diffusion tensor model [28,29,30]. The problem of interpolation on a spherical grid has already been tackled in physics, meteorology, and climatology applications [31]. A little bit of works used the angular interpolation in the field of DW-MRI with the aim of reducing the acquisition time. The focus of the current study was to propose a complete system of tractography and to demonstrate the advantages of using angular interpolation in tractography process. This work is based on the preliminary results presented in [32]. The evaluation of our method was performed using the data and results from the ISMRM tractography challenge 2015, including 96 tractography methods for a quantitative comparison. A standard method commonly used in the literature to perform fiber-tracking will be chosen and the evaluation will be done following the Tractometer approach. We will also use in vivo human brain data collected by clinical imaging protocols to demonstrate the potential of our fiber detection technique on real data, where we recover multiple fiber crossings on known regions of interest in the human brain. Finally, we will preform qualitative comparisons between our approach and the widely employed in clinical studies (DTI model). Therefore, the main contributions of this paper can be summarized as follows:

  1. 1.

    Developing a full fiber tractography pipeline based on an angular interpolation in gradient direction domain,

  2. 2.

    Quantitative evaluation of the performance of our pipeline on the phantom of the ISMRM 2015 Tractography Challenge,

  3. 3.

    Quantitative comparison of all methods on common diffusion MRI metrics,

  4. 4.

    Quantitative and qualitative examination of the proposed pipeline in vivo data.

Thus, this paper will be structured as follows: we present in the second section the steps of the followed process to obtain the brain connectivity: from the HARDI acquisition technique to the tractography method. Moreover, we present the ground truth dataset specifically designed for the contest followed by the evaluation methodology. We will, also, describe the real data that were acquired from a healthy adult volunteer. In the third section, we present the experimental results. In the first part, we compare our approach to the results of the ISMRM 2015 Tractography Challenge. In the second part, we present the experiments in vivo brain DW-MRI data by specifying all the steps leading to the tractography. We discuss the results in the fourth section.

2 Methods

In this section, we present, first of all, the steps of the followed process to obtain the brain connectivity. We describe, also, the used database and the evaluation criteria. We, finally, present real data acquired from a healthy adult volunteer.

2.1 Tractography processing pipeline

We present in this part the proposed tractography pipeline. It incorporates four steps. The first step consists on acquiring a reduced set of DW images. Then, the DW images are estimated in the missing directions, using the proposed interpolation method. In the third step, we estimate the fiber orientation distribution (FOD) using the spherical deconvolution (SD) methods. It provides an accurate depiction of fiber crossing. Finally, we use these orientations across the whole brain to reconstruct the fiber bundles using tractography methods.

In this work, the main contribution concerns the second step. We developed an angular interpolation algorithm that allows us to estimate DW image intensity at a given location in any direction from the acquired DW images. Figure 1 shows a workflow of the overall proposed pipeline.

Fig. 1
figure 1

The proposed tractography pipeline

Full size image

2.1.1 Step 1: Data acquisition

The first step consists on acquiring a serie of diffusion-weighted (DW) signals by applying N gradients of different strengths and directions. It should be noted that these N directions are uniformly distributed on the sphere of the q-space. The application of a single pulsed gradient produces one DW image that corresponds to a single point in the q-space, as shown in Fig. 2.

Fig. 2
figure 2

Explanation of the diffusion-weighted image acquisition procedure in the HARDI technique

Full size image

In the proposed framework, only a subset of the total number of gradient directions is acquired, which reduces the scan time significantly while making the technique robust to estimate the diffusion signal. In previous work [32, 35], we showed through numerical experiments, that as few as 34 directions of an original set of 64 diffusion gradients are sufficient for recovering intravoxel fiber structure. We used a general mixed integer linear programming (MILP) framework proposed by Cheng et al. [33] to design “uniform” sampling schemes. Then, we use a non-local mean filter with Rician noise correction to denoise the DW data before HARDI reconstruction. This method was shown to have the desired effect of correction for the noise bias without blurring-out figure crossing information [34].

In a diffusion image, each voxel has a signal which results from the motion of a large number of water molecules. Hence, in a voxel, the diffusion signal Sa{i, j} is obtained by the collection of the magnitudes in all diffusion-weighted images of the brain. It is defined in our case using 34 directions by the following equation:

$$ S_{a}\{i,j\} = \{E_{q_{1}}(i,j),E_{q_{2}}(i,j),...E_{q_{34}}(i,j)\} $$
(1)

such as:

$$ E_{q_{d}}(i,j)=I_{q_{d}}(i,j)/I_{q_{0}}(i,j) $$
(2)

where 1 ≤ d ≤ 34, (i, j) corresponds to the voxel location in a DW image, \(I_{q_{d}}(i,j)\) is the measured signal obtained from the acquired DW image in the direction qd, and \(I_{q_{0}}(i,j)\) is the unweighted one obtained without any applied diffusion gradients.

2.1.2 Step 2: Interpolation

In the second step, we explore the interpolation in the diffusion gradient direction domain as a method to estimate the HARDI signal from a reduced set of gradient directions. We describe an overview of the proposed approach steps as illustrated in Fig. 3.

Fig. 3
figure 3

A flow chart of the DW images estimation using the proposed approach

Full size image

Let denote by {Ii}i= 1..34 the set of acquired images and by {Iei}i= 1..30 the set of missing images to be estimated. For an image Ii or Iei, its corresponding direction on the sphere of the q-space will be denoted respectively by P(Ii) and P(Iei). Having acquired 34 DW images, we estimate in the second step the DW images in the missing directions. We propose to compute the Delaunay triangulation diagram of the 34 points subset. A missing direction P(Ie) corresponding to an image Ie will have three neighboring points (P(I1), P(I2), P(I3)) that correspond to the points of the nearest triangle from the Delaunay triangulation. Then, the missing DW images are estimated. We compute the DW value of each voxel as a linear combination of its determinated neighboring points DW values.

$$ Ie(i,j)=\alpha_{1}I_{1}(i,j)+\alpha_{2}I_{2}(i,j)+\alpha_{3}I_{3}(i,j) $$
(3)

For 1 ≤ t ≤ 3, we formulate the value of αt as follows to characterize the participation of the point P(It) in the estimation of the missing DW value of the point P(Ie).

$$ \alpha_{t}= \frac{\text{Area}(P(Ie),P(I_{m}),P(I_{n}))}{\text{Area}(P(I_{1}),P(I_{2}),P(I_{3}))} $$
(4)

where (m, n) ∈{{1, 2, 3} × {1, 2, 3}}, t ≠ m, t ≠ n and n ≠ m.

Area(P(I1), P(I2), P(I3)) corresponds to the area of the triangle composed by the three points (P(I1), P(I2), P(I3)). After estimating the DW signal in each voxel using the proposed approach, we reconstruct the DW image. This procedure is repeated to estimate all the DW images in the missing directions.

In HARDI technique, we used a full set of directions (N = 64) uniformly distributed on the sphere to avoid issues with imperfections in the uniformity of the DW gradient directions [43]. Thus, the number of missing directions is equal to 30 directions (64 − 34). The diffusion signal in each voxel will be:

$$ S_{e}\{i,j\} = \{Ee_{q_{1}}(i,j),Ee_{q_{2}}(i,j),...Ee_{q_{30}}(i,j)\} $$
(5)

such as:

$$ Ee_{q_{d}}(i,j)=Ie_{q_{d}}(i,j)/I_{q_{0}}(i,j) $$
(6)

where 1 ≤ d ≤ 30, \(Ie_{q_{d}}(i,j)\) is the measured signal obtained from the estimated DW image in the direction qd, (i, j) and \(I_{q_{o}}(i,j)\) are the same terms used in Eq. 2.

Hence, the final diffusion signal will be the combination of the acquired signal Sa and the one estimated using the proposed approach Se. Our approach uses the twin concepts of HARDI technique and linear interpolation to reduce acquisition time. Therefore, within our approach, we seek an optimal balance between the number of samples and reconstruction accuracy. Importantly, we apply this method on DW images before the diffusion modeling step, removing the limitation of applicability to a specific model.

After recovering the HARDI signal, we will estimate the spherical function that characterizes the intravoxel fiber orientation.

2.1.3 Step 3: FOD estimation

We estimate the fiber orientation distribution (FOD) using the spherical deconvolution (SD) methods. It provides an accurate depiction of fiber crossing. Precise local reconstruction of the principal diffusion directions is critical to obtain accurate fiber tracts for visualization. The spherical deconvolution methods use the HARDI data to estimate the orientation of multiple intravoxel fiber population within the regions of complex white matter architecture. In fact, a good local estimation of the configuration in each voxel permits a good fiber tracking.

2.1.4 Step 4: Tractography

Finally, we use these orientations across the whole brain to reconstruct the fiber bundles using tractography methods. In order to better deal with fiber crossings in tractography, HARDI-based techniques are needed. Fiber tracking algorithms can be categorized as deterministic [37] and probabilistic [38]. The first tractography method uses the local fiber orientations to delineate the whole trajectory. Starting from a seeding point, the fiber trajectory is computed in a step-wise fashion based on a local measurement orientation until termination criteria are met. Stopping criteria are, usually, a brain mask, low fractional anisotropy (FA) values, or maximum curvatures. As opposed to deterministic fiber tracking, probabilistic fiber tracking was proposed to explore the connection probability. The output of such algorithms is a probability map of connection between the given seed voxel and other voxels of the brain.

2.2 HARDI databases

Our first objective was to quantitatively evaluate our proposed pipeline against a known ground truth. Thus, we evaluate the performance of our methods on the 2015 ISMRM challenge dataset using the Tractometer connectivity metrics. We, first, show how to simulate the DW data. Secondly, we define the global metrics included in Tractometer system. Afterwards, we will focus on real data to demonstrate the performance of the proposed approach. We will describe the acquisition protocol of DW real images.

2.2.1 Realistic brain-like phantom data

- Data generation :

In this part, we will first briefly review the construction of a realistic diffusion ground truth dataset. The evaluation of analysis methods remains challenging due to the lack of gold standards and validation frameworks. To do so, we used a phantom with a known connectivity. We employed a version of the ISMRM 2015 Tractography Challenge data simulated with a b-value of 3000s/mm2 and 64 non-zero gradient directions. The other simulation parameters are the same as for the original phantom. We generate the DW data of the realistic diffusion MRI phantom using the simulation software FiberFox [39].

This phantom is developed for the evaluation of structural connectivity pipelines. It consists on a set of 25 major fiber bundles. Figure 4 illustrates these 25 bundles that serve as ground truth models. This phantom contains realistic fiber configurations: crossing, kissing, and bending. The ground truth is available and enables to correctly compare tractography algorithms. For details about the ISMRM tractography challenge 2015 and the phantom, we can refer to Maier-Hein et al. [40] and the corresponding challenge homepage.Footnote 1

Fig. 4
figure 4

Illustration of the used phantom: a Data generation process and b its 25 fiber bundles

Full size image
- Evaluation criteria :

As one of the aims of this work is to quantitatively assess the quality of the tractograms obtained on the phantom dataset, we score the results using global connectivity metrics included in the Tractometer system [41]. Table 1 summarizes the metrics, with abbreviations and definitions.

Table 1 Metrics and definitions [41]
Full size table

2.2.2 Human brain data

We used in vivo data, that were acquired from a healthy adult volunteer using a 3-T Siemens system. We describe below the acquisition protocol, which is also summarized in Table 2. We consider diffusion images acquired on single spherical shell in q-space.

Table 2 Summary of the used parameters for data acquisition
Full size table

HARDI sequence is characterized primarily by two parameters: the number of DW directions and the b-value at which the DW images are acquired. A recent study [43] suggests that this number should be at least equal to 45 to capture all the angular features of the DW signal at b = 3000 s/mm2. Thus, the DW images were acquired with 64 uniformly distributed gradient directions at b-values of 3000 s/mm2. Acquisition time in diffusion MRI increases with the number of diffusion-weighted images that need to be acquired. The total scan time of the HARDI acquisition was approximately 15 min in duration. The acquisition time represents a major problem of the HARDI technique. In fact, its duration is long and thus, this technique is clinically infeasible.

3 Results

We focus, in this section, on the experimentation procedure. The first part of the experimentation is conducted on simulated DW images of a brain-like geometry. The goal of this first part of study is to perform a quantitative evaluation based on the Tractometer connectivity metrics. We, first, compare the reconstruction using the full set of DW images and using the proposed interpolation approach. Then, we compare the performance of the proposed algorithm against several alternative solutions. The second part of the experimentation is performed on real data. We develop two tests to evaluate the performance of our pipeline. Firstly, we, quantitatively, compare the DW images estimated in each direction with the real acquired data. We, then, qualitatively evaluate our method by fiber tracking on in vivo data.

3.1 Results on simulated phantom data

3.1.1 34 versus 64 directions

We perform, in this part, a quantitative evaluation of the proposed approach for the reconstruction of the brain connectivity. We, therefore, quantitatively compare between the reconstruction with the proposed approach and the one obtained with the original set of directions (64 directions) using global connectivity metrics. We will establish probabilistic tractography using constrained spherical deconvolution with MRtrix software [36] to evaluate the potential of the proposed approach in connectivity. The used method in MRtrix involves the estimation of the fiber orientation distribution (FOD) using the CSD method [18, 19], which requires acquired data using the high angular resolution DW imaging (HARDI) strategy [42]. The FOD estimated using CSD provides an accurate estimation of the amount of white matter fibers aligned with any given direction. For this reason, FOD-based tractography algorithms were used. It allows the tracking through regions of crossing fibers.

Based on the known ground truth bundles, we calculated the global metrics.

In agreement with the results obtained in our previous work [32], no significant difference was found between the two approaches. Clearly, Table 3 shows that the reconstruction of fibers with the proposed approach is similar to the one obtained with the full directions. The method was able to reconstruct 23/25 valid bundles in the ISMRM tractography challenge phantom, while reconstructing 53/2450 invalid bundles. We obtained a good ratio VB/IB. The proposed approach provides a good compromise between accuracy and computational cost.

Table 3 Comparison between the results obtained using classical probabilistic CSD tractography and the proposed pipeline
Full size table

The obtained results show that our novel method enables reducing the number of measurements commonly used in DW-MRI acquisition. This method has proved also its potential in brain connectivity.

3.1.2 Comparison with state-of-the-art tractography pipelines

This experiment allows us to compare our approach with all 96 original submissions of the tractography challenge, comprising a large variety of tractography pipelines with different pre-processing, local reconstruction, tractography, and post-processing algorithms. The quantitative results of our algorithm are summarized in the Table 4 and the final results are accessible on the Tractometer website.Footnote 2

Table 4 Quantitative results of the proposed algorithm compared to the mean scores over all 96 benchmark submissions
Full size table

The results demonstrate the vast potential of interpolation for fiber tractography. Compared to the mean scores over all 96 benchmark submissions, the proposed approach performed better in terms of valid bundles and invalid bundles. Table 4 shows that on average 21 out of 25 VBs were identified by the participating teams. Most teams (61.46%) find the valid bundles (VB= 23) but also find many invalid ones [40]. The percentage of invalid bundles (53/2450) is rather low compared to the majority of benchmark algorithms. Submitted tractograms contained an average of 88 ± 59 IBs. This fact demonstrates the inability of current state-of-the-art tractography algorithms to control for false positives. Our approach outperformed more than half of the challenge submissions (67 Submissions). Thus, we enhance Trade-off between valid and invalid which leads to better connectivity results. These experiments showed that our method is especially efficient to resolve regions where fiber bundles are crossing, and still well handle other fiber bundle configurations. Figure 5 reports the quality of the reconstructions in terms of valid an invalid bundles. A good pipeline tractography should find a high number of valid bundles and a low number of invalid bundles. As can be appreciated from Fig. 5, the proposed pipeline clearly outperforms a large selection of methods of state-of-the-art according to the valid and invalid bundles criteria.

Fig. 5
figure 5

Scores of the tractograms obtained using the proposed approach (PA)(red color) in comparison to the original challenge submissions (gray) in terms of Valid/Invalid Bundles

Full size image

The large number of invalid bundles that were identified by most state-of-the-art algorithms compared to our obtained score (53 IB) demonstrates the advantage of our pipeline.

Additionally, we quantitatively compare the performance of the proposed approach in terms of valid/invalid connections (VC/IC). Figure 6 shows the obtained results with the proposed method compared to the 96 submissions. Here, a good method must ensure a low number of IC and a large number of VC. The obtained results in term VC/IC is very near to the average of all submissions. It is also important to note that there are many methods that outperform the preposed approach according to this criteria. Table 5 shows the number of the original challenge teams that are outperformed by the presented approach with respect to the mean scores over all 96 submissions.

Fig. 6
figure 6

Scores of the tractograms obtained using the presented method (red color) in comparison to the original challenge submissions (gray) in terms of valid/invalid connections

Full size image
Table 5 Quantitative results of the proposed algorithm compared to the mean scores over all 96 benchmark submissions
Full size table

No submission contained all 25 VBs, but 10 submissions (10.4%) recovered 24 VBs, and 59 submissions detected 23 VBs. In the case of the IB score, 67 submissions could be outperformed. The results indicate that interpolation methods can be successfully applied to DWI datasets for mining anatomical details that are normally seen only at higher resolution, which will aid in tractography and microstructural mapping of tissue compartments.

3.2 Results on human brain data

3.2.1 Test 1: DW image estimation

In order to provide a meaningful basis for quantitative comparison, in the absence of ground truth connectivity information, we use the full (64 directions) high angular resolution diffusion-weighted imaging (HARDI) signal acquired by the MRI system as a ground truth. We estimated the DW signal in each voxel by the interpolation at 30 directions with knowledge of the 34 directions. Figure 7 displays an example of an axial slice of the generated DW image by the MRI system and the estimated DW image using the proposed approach, in the same gradient direction.

Fig. 7
figure 7

Exemplary illustration of the estimation in real data: a Acquired DW image of axial slice. b Estimated DW image of axial slice with the proposed approach

Full size image

Visual evaluation shows results that are in line with the first experiment. We, furthermore, computed the mean error in each voxel calculated for each sample j ∈{1..30} such as:

$$ \text{error}=\frac{|S^{\text{acquired}}_{j}-S^{\text{estimated}}_{j}|}{S^{\text{acquired}}_{j}} $$
(7)

We consider only the voxels inside the brain. Figure 8 shows the obtained error for all slices, displayed as a box plot.

Fig. 8
figure 8

Box plot representing the error of the different slices on real data acquisition

Full size image

The median value is indicated by the location of the line that divides the box into two parts. The medians of all slices are approximately at the same level. We obtained a very low mean error for all slices. As can be appreciated, the reconstruction of the diffusion signal from undersampled data using the proposed method yields accurate results.

3.2.2 Test 2: Probabilistic fiber tracking

Here, we, qualitatively, illustrate that the gain in angular resolution remarkably improves the fiber tracts. Whole-brain probabilistic tractography was finally performed using MRtrix software and we carry out the same procedure for the human brain data. The objective is the validation of the tractography algorithms based on the interpolation in gradient direction domain with a focus on their ability to correctly reconstruct crossing regions, between two and three different fiber populations, in the known regions of interest in the human brain.

The tractography reconstruction is realized by both the original signal (Fig. 9a) and the one constructed by the proposed approach (Fig. 9b). We compared the reconstructions obtained from undersampled data (34 directions) to those with fully-sampled data (64 directions), considering this latter as ground truth. In the absence of ground truth connectivity information, this coronal slice was used to illustrate the qualitative analysis of the fiber reconstruction results.

Fig. 9
figure 9

Qualitative comparison between fiber tracts reconstructed using a the original dMRI data and b the proposed method

Full size image

We can denote from Fig. 9 that the performance of the proposed method is comparable to the gold standard using 64 directions. The obtained results show the accuracy of the proposed approach to perform the brain connectivity. In fact, the trajectory and the orientation of the white matter fibers is preserved compared with the reconstruction using the full 64 directions.

The results were also compared to current standard tensor tractography algorithms [46] using 34 DW images. We examined a well-known region in the white matter containing crossing fibers (the region highlighted by the blue box in Fig. 10) consisting of intersecting callosal fibers (red), corticospinal fibers (blue), and association fibers (green) as revealed by Wedeen et al. [47] using the reference standard technique (diffusion spectrum magnetic resonance imaging (DSI) tractography).

Fig. 10
figure 10

Qualitative comparison between fiber tracts reconstructed using a the proposed method and b the diffusion tensor model

Full size image

Figure 10b illustrates the inability of the DTI model to resolve the intersection between fibers. The proposed approach clearly outperforms the standard diffusion tensor model. We illustrate the additional information that can be gained from the proposed approach: note that in this case, the tensor-based algorithm does not identify the full extent of the target white matter pathways. Our approach permits tracking through regions containing crossing fibers, thus addressing the well-known limitations of the diffusion tensor.

Our results show that fiber crossing models are able to reveal connections between more brain areas than the simple tensor model. The proposed approach can reliably resolve crossing white matter tracts. In contrast, the diffusion tensor model often fails to accurately represent the microstructure at this junction and only depicts the orientation of the CC.

Thus, linear interpolation combined with high angular resolution diffusion imaging (HARDI) technique achieves significant improvement. It displayed more compact fiber bundles compared with the DTI-based results in all cases.

We are able to reconstruct fibers with a clinically acceptable number of measurements.

4 Discussion

We have proposed an integral framework to perform brain connectivity. It incorporates an interpolation step on DW images before the diffusion modeling to accurately describe crossing fiber bundles from a reduced number of measurements. The motivation for this choice of approach was to create a solution comparable to a clinical acquisition in terms of number of directions and acquisition time.

Importantly, we highlight that our proposed approach can be used with any diffusion model and any tractography algorithm, removing the limitation of applicability to a specific model such as diffusion tensor.

We have demonstrated its improved accuracy through experiments on phantom based on a simulated human brain data set with ground truth tracts and real HARDI data. We used the metrics included in the Tractometer system to evaluate tractograms reconstructed from interpolated data: Valid Bundles (VB), Invalid Bundles (IB), Valid Connections (VC), Invalid Connection (IC), and No Connections (NC).

The capability of the presented approach to estimate complex sub-voxel fiber geometries was successfully demonstrated in the first experiment conducted on simulated DWI of a brain-like geometry. We confirm that as few as 34 directions are sufficient for the reconstruction of the white matter fibers compared with the full 64 directions. The obtained results show that the reconstruction of fibers with the proposed approach are similar to the one obtained with the full directions.

We quantify the performance of our methods on the 2015 ISMRM challenge dataset and evaluate using the Tractometer connectivity metrics. In doing so, we can compare ourselves to the 96 original challenge submissions. Our proposed pipeline achieves competitive results to other challenge submissions. Overall, we have shown that our pipeline provides a good ranking in terms of valid/invalid bundles.

The proposed pipeline successfully recovered the most valid bundles (23 VBs, recall that there are 25 VBs in the phantom). In general, the presented approach yields one of the best performances in terms of valid bundles, while keeping a low invalid bundles.

One interesting aspect that catches the eye is the fact that the presented approach outperforms a large number of submissions in term of valid/invalid ratio, which is a known issue of current fiber tractography approaches. However, we should also mention that in term of VC/IC many methods outperform the proposed approach but it is very near to the average of all submissions.

On real human data, we present a quantitative comparison between the DW images acquired by the MRI system and those estimated using the proposed approach. We, furthermore, computed the mean error in each voxel. We obtained a very low mean error for all slices. We, then, performed a visual evaluation between tracts reconstructed from the original (64 directions) and interpolated dMRI data (34 directions). The second experimental part compares the proposed method with tractography obtained from DTI in vivo human data. DTI cannot accurately describe the microstructure in complex white matter voxels that contain more than one fiber population, due to intersecting tracts. Qualitatively, we show that fiber tracts reconstructed with the proposed approach reveal exquisite details beyond that are achievable with the diffusion tensor model using the same number of DW images. This approach improves the ability to resolve different diffusion directions within the same voxel that result from crossing axonal bundles. The proposed approach potentially reduces the scan time of HARDI technique from 15 min to approximately 7 min while maintaining its precision and accuracy, thus driving the acquisition cost of HARDI closer to DTI.

One limitation of our study is that the proposed pipeline was applied on a one healthy volunteer. We have to apply the technology developed in this work on several subjects in different brain regions and on clinical applications such as traumatic brain injury and multiple sclerosis to confirm the validity and clinical utility of our approach.

5 Conclusions

In conclusion, angular interpolation has a particular role in diffusion data with substantially reduced number of diffusion encoding directions, which allows a reduced scanning time. The experimentations are conducted on two steps. In the first one, tests are performed on realistic brain-like phantom with a known ground truth. We thoroughly evaluated the performance of the presented framework in comparison to 96 state-of-the-art tractography pipelines. The experimentations are realized in the second step on in vivo human brain data collected by clinical imaging protocols.

The obtained results show the accuracy of the proposed approach to perform the brain connectivity. We show that our method improves the angular resolution, as well as fiber crossing discrimination. Thus, our method is a new approach that combines the advantages of the estimation of multiple intravoxel fiber populations of HARDI with the clinical feasibility of routinely used DTI image data acquisition in clinical practice.

While the presented results are promising, there are still some challenges to address. Further work is necessary to quantify and improve the performance of the presented approach. The ability to recover reliable and accurate intravoxel fiber distributions within the human brain is promising and opens new perspectives for doing extensive validation on several subjects in different brain regions and we plan to apply the technology developed in this work to clinical applications such as traumatic brain injury and multiple sclerosis to confirm the validity and clinical utility of our approach. Further investigations are needed to test more complex interpolation techniques. For clinically orientated future investigations, the most important steps will be the integration of the novel approach in the MRI system besides conventional DTI-based tractography.