The challenge of mapping the human connectome based on diffusion tractography

Tractography based on non-invasive diffusion imaging is central to the study of human brain connectivity. To date, the approach has not been systematically validated in ground truth studies. Based on a simulated human brain data set with ground truth tracts, we organized an open international tractography challenge, which resulted in 96 distinct submissions from 20 research groups. Here, we report the encouraging finding that most state-of-the-art algorithms produce tractograms containing 90% of the ground truth bundles (to at least some extent). However, the same tractograms contain many more invalid than valid bundles, and half of these invalid bundles occur systematically across research groups. Taken together, our results demonstrate and confirm fundamental ambiguities inherent in tract reconstruction based on orientation information alone, which need to be considered when interpreting tractography and connectivity results. Our approach provides a novel framework for estimating reliability of tractography and encourages innovation to address its current limitations.

shows the corresponding smoothed streamlines where the fiber endpoints are mostly distributed at the bundle end-points. This is called the "refined bundle".

Supplementary Figure 8. Restrictive definition of a valid connection by tractometer metrics.
Illustration of the corticospinal tract (CST) bundle volume (a) and endpoints distribution (b) in red overlay. As soon as a submitted streamline (c) has a single point exiting the bundle volume or not exactly terminating within the endpoints distribution, the streamline will be flagged as an invalid connection. Figure 9. Results using traditional scoring. In the Tractometer metrics, a valid connection is defined as a streamline that connects expected regions of interest without exiting the expected bundle mask. This definition is extremely demanding for current tractography algorithms (please note the altered axis scalings). When scoring based on these metrics, on average over all submissions, only 3% ± 2% of the submitted streamlines were valid (a). 13 ± 3 valid bundles were detected, accompanied by 376 ± 149 invalid bundles detected (b). The average volume overlap was 6% ± 10% (c). At best, a technique achieved 8% valid connections and 92% invalid connections.

Team 3 -Submissions 3_0 to 3_4 (5 submissions)
The diffusion images were processed in DSI Studio (http://dsistudio.labsolver.org). Eddy current correction with b-table rotation was conducted by registering each DWI to the b0 image. In submission 1 to 3. The DWI data were further upsampled to the T1 space for 3. In submission 1, 2, 4, the diffusion data were reconstructed using generalized q-sampling imaging 11 with default parameters assigned in DSI Studio. In submission 3 and 5, diffusion decomposition 12 was used to reconstruct the data. A deterministic fiber tracking algorithm was used to generate whole brain fiber tracking 13 . An angular threshold of 70 degrees was used. The step size and anisotropy threshold were determined automatically by the tool. A hierarchical clustering method 14 was applied. The 50 largest clusters were selected and merged, whereas the small clusters were removed. The final results were submitted in trk format.

Pipeline 1 corresponds to 3_3
Pipeline 2 corresponds to 3_0 Pipeline 3 corresponds to 3_4 Pipeline 4 corresponds to 3_2 Pipeline 5 corresponds to 3_1

Pipelines description
The fiber bundle analysis were composed for three essential segments: (1) eddy current distortion, fieldmap correction and motion correction, (2) tensor calculation, reconstruction the fiber and (3) fiber classification into 26 fiber bundles by EM-step, and merge into whole brain tractogram.
-Eddy Current Correction is correct for these distortions, artifacts (include shear, enhanced background, image intensity loss, and image blurring), and head motion. For those correction perform affine registration to a reference volume by FSL (FMRIB). EPI distortion and fieldmap correction using FSL's FUGUE processing.
-Start the tensor fitting and calculation of fractional anisotropy (FA) value. The fiber-tracking used the streamlines approach, optimized with nearest neighbor interpolation. The range of the fiber length are 20 to 500 (mm).
-Fiber classification into 26 fiber bundles by EM-step with ten features (include first point, middle point and last point of the coordinate (x, y, z) ; length of the fiber) are extracted into the EM (expectation-maximization) step. To correct the fiber pathway are automatic classification into 26 fiber bundles by EM algorithm and select valid bundles.

Q Ji, W E Reddick, J O Glass
Pipeline descriptions :  Eddy current and motion corrections.  FSL bedpost were run on the data.  17800 seeds voxels were chosen within white matter randomly. 10 fibers for each seed were generated using probabilistic fiber-tracking implemented in FSL 15 . Each streamline was recorded.  Spline processing and a filter were applied to each generated fiber.

Pipeline 2 : 5_0
 Eddy current and motion corrections  17800 seeds voxels were chosen within white matter randomly, a stream line based algorithm was used to track 10 fibers from each seed voxel.  Filter was applied on each fiber, the fibers below threshold were discarded. . Y Feng, C Gao, Y Wu, J Ma, R He, Q Li, C-F Westin

Motion correction
Correcting for subject motion & eddy current induced geometric distortions: Only bulk motion considered is the general object movement, the misalignment will cause blurring of the diffusion tensor image, as well as the related orientation of diffusion. The rigid registration with quaternion is applied to all gradient directions to align them with B0 image. The reorientation of the corresponding B-matrix is completed so that orientation information is correctly preserved 18 .

Distortion correction
Correcting for EPI/susceptibility distortions (shape and intensity): The phase direction distortion of EPI is corrected with the provided field map. The approach we implemented is based on the work of 19 . In this way, the geometry distortion correction of subject motion and EPI distortion can be performed in one interpolation step to minimize blurring effects. However, geometric distortion of EPI is accompanied by compression or stretching of voxels so that the voxel brightness is altered as well as voxel position. To facilitate intensity distortion correction, we complete the geometry distortion correction of subject motion and EPI distortion in two separate steps. When a local field gradient causes a great phase change across a voxel, the signal from that voxel is not displaced but lost all together due to signal dephasing. Afterwards, there are unrecoverable EPI signal losses.

3.Diffusion model beyond DTI.
1. Preprocessing of DWI data. The brain mask is obtained using Mrtrix by averaged DWI images with default threshold. Then, a tensor map and FA map are generated by DTI model. High anisotropy mask is extracted by FA map with FA threshold value 0.2. Finally, the WM mask is produced by the intersection of brain mask and high anisotropy mask. 2. Global consistency model. A spherical double-lobe basis function is used to form an overcomplete dictionary that guarantees the local sparsity of FOD. And a global consistency spatial model which incorporated a spatial priori information based on the Bayesian framework is developed using the coefficients of the basis function, which integrates the local fiber distributions in the neighbourhood, detailed in Wu et al. 20 . Due to the deconvolution-based process, the response function is estimated by the average of diffusion signal with FA threshold value larger than 0.7. The exponent enhancement number used in the experiment is 10.
3. Deterministic tracking. A HARDI based streamline tracking method which chooses the closed fiber peak orientation extracted from FOD of the last tracking step as the tracking direction is carried out using formula Vnext = f*Vcurr + (1-f)* Vprev with f=0.7. The angular threshold for contiguous step is 50 degree to avoid unreasonable turn. The step length is 1mm. The seed point position is random generated in per WM voxel. The range of the fiber length is between 50 mm and 500 mm, and fibers whose length out of range are exclude. In our result, 13871 fibers are generated totally and 10000 fiber trajectories are reserved for the Pipeline.
Pipeline 1-7_0  + Unspiking (in-house algorithm) + motion correction: We defined a graph where each of its n nodes represents a dMRI scan (excluding the B0 image). The edge weights are given by the angle between the b-vectors associated with their edge nodes. We computed the minimum spanning tree (MST) of the graph and rigidly registered all pairs of images corresponding to adjacent nodes in the MST by maximization of their mutual information (exactly n-1 registrations). The registration between any pair of scans can then be computed by following the path between their corresponding nodes in the tree and composing the transforms associated with the path's edges. We selected the center of the tree as reference and registered the rest of the images towards it. The average of the dw images has higher SNR and can then be registered to the B0. + Susceptibility-induced geometric distortion correction: We directly applied the deformation field defined by the off-resonance field provided by the organizers of the challenge. The bandwidth per voxel along the phase-encoding direction was also known to be 9.26 Hz per voxel. + denoising Non Local Spatial and Angular Matching (NLSAM) 21 with 5 angular neighbors + Spherical finite rate of innovation (soon to be in dipy) 22

Pipelines description
The different pipelines are outlined in Fig. 1 (numbered 10_0 to 10_19). A pipeline typically consists of the following steps, and we have tested different combinations of options for each step:  Rician noise correction: Noise removal by convolution of a linear minimum mean square error estimator using a Rician noise model with the diffusion weighted images (DWIs) using a 5 x 5 pixels Gaussian kernel 26,27 . The SNR was determined from the background noise and was found to be 20 for the b=0 image.
 Rigid registration: A rigid (6 degrees of freedom) coregistration technique based on mutual information was used to realign the DWIs to the first non-DWI 28,29 . It was wrongly assumed here that no B-matrix rotation was applied when generating the synthetic dataset, and that performing B-matrix rotation would therefore result in "overcorrection".
 Fieldmap correction: Image deformations due to field inhomogeneities were corrected based on the provided B0 field map. Signal accumulation and loss was not corrected 27 .
 Registration to T1: Non-rigid registration of the fractional anisotropy image to the T1 image 30 .
 CSD: Determination of the fiber orientation distribution function with constrained spherical deconvolution (lmax = 6) and recursive calibration of the response function (peak ratio threshold 0.01) 8,28,33 . In case outlier detection with REKINDLE was used in the previous step, rejected points were interpolated based on the robust diffusion tensor, so that an equal amount of directions per voxel could be used for CSD. We defined a graph where each of its n nodes represents a dwMRI scan (excluding the B0 image). The edge weights are given by the angle between the b-vectors associated with their edge nodes. We computed the minimum spanning tree (MST) of the graph and rigidly registered all pairs of images corresponding to adjacent nodes in the MST by maximization of their mutual information (exactly n-1 registrations). The registration between any pair of scans can then be computed by following the path between their corresponding nodes in the tree and composing the transforms associated with the path's edges. We selected the center of the tree as reference and registered the rest of the images towards it. The average of the dw images has higher SNR and can then be registered to the B0. + Susceptibility-induced geometric distortion correction: We directly applied the deformation field defined by the off-resonance field provided by the organizers of the challenge. The bandwidth per voxel along the phase-encoding direction was also known to be 9. Then we detected spiked slices as having voxels in Fourier space in abnormal range compared to other directions of the same slice. Raw data spiked slices were re-interpolated from the 6 neighboring directions in term of diffusion vector encoding. Only neighbors non-spiked same and adjacent slices were considered to avoid correcting with artifacted data. Interpolation used motion parameters to account for motion between directions acquired, and 6 interpolated slices were averaged.
Then, interpolation of motion was applied to despiked data first filtered with Non-Local Means.
FSL topup tool corrected for EPI inhomogeneity induced distortions using reversed phase encoded B0 38 .

Anatomical Imaging
T1 was analyzed with adapted freesurfer recon-all routine by skipping skullstriping step, replacing masking with AFNI 3dAutomask, and using template without skull for all operations to accomodate the simulated T1 scan provided.
Anatomical priors for Anatomically Constrained Tractography 39 as 5 Tissue Type (5TT) format for MRTrix was prepared from freesurfer surfaces and parcellation by computing partial volume maps for: -cortical and cerebellar gray matter, -subcortical gray matter,white matter, -csf, -abnormalities.
A mock gray matter region was set at the end of the brainstem to allow spinal streamlines to be found, and anterior (AC) and posterior commissure (PC) were added to white matter map.

-Tractography
We initiated tractography with a seeding dynamically using SIFT model 40 . The 2 nd order Integration over Fiber Orientation Distribution (iFOD2) algorithm 41 provided by MRtrix in conjunction with ACT framework was used to recover streamlines 39 .

-Submissions
Submission -13_0 Starting from the raw tractrogram we used SIFT to improve its correspondence with the FODs reconstruction. It helped as well to reduce the size of the file and we ended with 430 000 streamlines.

Submission -13_2
Starting from the raw tractogram, we filtered it using different rules: -Streamlines that cross corpus collusum (CC) cannot cross brainstem (BS).
-Streamlines that cross BS cannot pass through both hemispheres.
-Streamlines that cross CC should require some symmetry aspects: -The medial point of the streamline that cross the CC should be within the 2nd third of the fiber -Within a streamline, points that cross the CC mask have to be contiguous.
-Splitting streamline at the medial point that cross the CC we compute the minimum average directflip (mdf) (dipy) of one half mirrored through medial wall with the other half to estimate the symmetry of the original streamline. Every streamline which has a distance higher than 4 is removed.
Filtered streamlines were then selected through SIFT resulting in 450 000 streamlines.

Submission -13_1
Starting from submission 2 before SIFT algorithm we ask a specialist in neuro-anatomy to manually select some bundles of streamlines using endpoints and waypoints landmarks. We added to this twelve manually segmented bundles, streamlines crossing the CC, AC, PC as well as those going through respectively BS and cerebellum. Finally, we again used SIFT algorithm and ended with 310000 streamlines.

Submission -13_3
Last submission is a compilation of all streamlines that have been removed during filtering for submission 2. While we expect this submission to be poorly evaluated, it is of interest to check if filtering has not been too conservative and removed relevant filtered streamlines. This submissions contains 295 000 streamlines.

2) MRTRIX-Diffusion-weighted MRI white matter tractography
. a) The pre-processed Diffusion image, registered on the pre-processed T1, and its brain mask were converted into mrtrix format (.mif) (mrconvert command). Same applies for b-values and gradient directions, where z coordinate was inverted (-z), for respecting mrtrix conventions.
. b) Tensor components and fractional anisotropy mask were generated as described in mrtrix documentation.
. c) A mask of high anisotropy voxels was generated by eroding the FA mask with a threshold of 0.7.
. d) The response function SH coefficients were estimated from the DW signal in the single-fibre voxels, with a harmonic order of 4.
. e) The constrained spherical deconvolution (CSD) was performed with a maximum harmonic order set to 4.

3) CONVERSION OF TRACTOGRAMS TO ORIGINAL SPACE
. a) gen_unit_warp script was used for extracting nowarp images from Diffusion image.
. b) applywarp was used for applying inverse transformation produced during registration of Diffusion-T1-standard space (step 1.c).
. c) normalize_tracks script was applied on the tractograms (.tck) produced during step 2.f, for converting them in the original space.
. d) Tractograms normalized to original space, were checked by using the script provided by the challenge organizers (validate_tracts_space.py).

A R. Khan, W Hodges, S Alexander
Pipeline is an amalgamation of tools from FSL and the ITK-based BrightMatter neurosurgical planning software (Synaptive Medical). Tracts are in VTK format, with RAS reference.

DWI artifact removal (step1)
• Denoising (optional): NLMEANS algorithm from the dipy project, with automatic estimation of noise variance using a coarse mask of white-matter (WM).
• Susceptibility-derived distortions: we applied the tool FUGUE from FSL. The fieldmap is computed using an in-house phaseunwrapping algorithm over the supplied field-map.
• Head motion: we used a dwi-to-dwi rigid-registration algorithm using ANTs. Settings for the pipeline, and the methodology used to accordingly rotate the b-vectors is contributed to the nipype tool, and actual code found in: https://github.com/oesteban/nipype/tree/enh/ReorganizeWorkflo ws.

Segmentation (step2)
• Preliminary image registration: the image provided in the challenge dataset was linearly registered to the artifact free b0 image using FSL • White-matter mask registration: an in-house method for simultaneous segmentation and registration of dMRI images, is used to refine the fitting of three surfaces onto the dMRI data. The three surfaces are the boundaries of the ventricular system, the white-matter/gray-matter interface and the pial surface. Tissue fractions were then used to compute GM and WM labels in diffusion space.

Parcellation (step3)
• The MNI template was non-linearly registered to the undistorted B0 image using FSL's FNIRT tool. The transformation was then applied to the Harvard-Oxford cortical and subcortical atlases as well as the cerebellar label of the MNI structural atlas provided in the FSL toolbox. GM voxels were then labeled according to the closest values in these atlases, and assigning different labels for left and right hemispheres. Missing labels were removed, leading to a subject atlas with 112 regions.

Reconstruction (step4)
• The reconstruction was implemented using a formulation that solves the fibre configuration of all voxels of interest simultaneously and imposes spatial regularisation directly on the fibre space. This reconstruction method allows us to exploit information from the neighbouring voxels, translating the natural smoothness of the anatomical fibre tracts through the brain into a certain spatial coherence of the FOD (Fiber Orientation Distribution function) in neighbouring voxels 42 .

Tractography (step5)
For all methods below, tractography was restrained to the WM mask computed during the Preprocessing step.
. Gibbs : Gibbs global tractography was computed using MITK and the undistorted diffusion signal.
The following parameters were used : particles had a length of 3mm, width of 1mm and a weight of 0.0008, start and ending temperatures were set to 0.1 and 0.001 respectively, the balance between external and internal energy was set to 0, the minimum fiber length to 20mm, the curvature threshold to 60° and the random seeds to 100. Fibres were then extended by half a voxel on each end, and we extracted the fibres that started and ended inside the GM mask.
. Streamlines : streamline tractography was performed using the fibre peaks generated in the reconstruction step. Tractography was done using the FACT algorithm implemented in MRTRix in order to generate 100'000 fibres. Fibres were then extended by half a voxel in both ends, and we extracted the fibres that started and ended inside the GM.
. Probabilistic tractography : probabilistic tractography was computed on the Spherical Harmonics representation of the ODF from the Reconstruction step. We used the iFOD2 algorithm implemented in MRTrix to generate 100'000 fibres. We then extended the fibres by half a voxel on both ends, and extracted a set of fibres that started and ended inside the GM mas

Pipeline 2 -16_2
Steps 1 to 5 were followed as in pipeline 1, and the output of probabilistic tractography was clustered using Nibabel's QuickBundle tool, with a distance threshold of 4mm, before concatenating the tractograms into a single superset. The trackfile was named Superset3.trk

Pipeline 3 -16_4
Fibres in Superset1 were asigned a weight according to their contribution to the undistorted diffusion signal using COMMIT with default parameters 43 . The fibres were then replicated proportionaly to their assigned weight in a file called Superset1_COMMIT.trk Pipeline 4 -16_1 Steps 1 to 5 were followed as in Pipeline 1. Output tractograms from probabilistic and streamline tractography were clustered using QuickBundles, and the clustered tractograms concatenated with Gibb's output. Fibres in the resulting superset were asigned a weight according to their contribution to the undistorted diffusion signal using COMMIT with default parameters. The fibres were then replicated proportionaly to their assigned weight in a file called Superset2_COMMIT.trk

Pipeline 5 -16_3
Fibres in Superset3 were asigned a weight according to their contribution to the undistorted diffusion signal using COMMIT with default parameters. The fibres were then replicated proportionaly to their assigned weight in a file called Superset3_COMMIT.trk

DWI artifact removal (step1)
-Denoising (optional): NLMEANS algorithm from the dipy project, with automatic estimation of noise variance using a coarse mask of white-matter (WM).
-Susceptibility-derived distortions: we applied the tool FUGUE from FSL. The fieldmap is computed using an in-house phaseunwrapping algorithm over the supplied field-map.
-Head motion: we used a dwi-to-dwi rigid-registration algorithm using ANTs. Settings for the pipeline, and the methodology used to accordingly rotate the b-vectors is contributed to the nipype tool, and actual code found in: https://github.com/oesteban/nipype/tree/enh/ReorganizeWorkflo ws.

Segmentation (step2)
-Preliminary image registration: the image provided in the challenge dataset was linearly registered to the artifact free b0 image using FSL -White-matter mask registration: an in-house method for simultaneous segmentation and registration of dMRI images, is used to refine the fitting of three surfaces onto the dMRI data. The three surfaces are the boundaries of the ventricular system, the white-matter/gray-matter interface and the pial surface. Tissue fractions were then used to compute GM and WM labels in diffusion space.

Parcellation (step3)
-The MNI template was non-linearly registered to the undistorted B0 image using FSL's FNIRT tool. The transformation was then applied to the Harvard-Oxford cortical and subcortical atlases as well as the cerebellar label of the MNI structural atlas provided in the FSL toolbox. GM voxels were then labeled according to the closest values in these atlases, and assigning different labels for left and right hemispheres. Missing labels were removed, leading to a subject atlas with 112 regions.

Reconstruction (step4)
-The reconstruction was implemented using a formulation that solves the fibre configuration of all voxels of interest simultaneously and imposes spatial regularization directly on the fibre space. This reconstruction method allows us to exploit information from the neighboring voxels, translating the natural smoothness of the anatomical fibre tracts through the brain into a certain spatial coherence of the FOD (Fiber Orientation Distribution function) in neighboring voxels 42 .

Tractography (step5)
For all methods below, tractography was restrained to the WM mask computed during the Preprocessing step.
-Gibbs : Gibbs global tractography was computed using MITK and the undistorted diffusion signal.
The following parameters were used : particles had a length of 3mm, width of 1mm and a weight of 0.0008, start and ending temperatures were set to 0.1 and 0.001 respectively, the balance between external and internal energy was set to 0, the minimum fiber length to 20mm, the curvature threshold to 60° and the random seeds to 100. Fibres were then extended by half a voxel on each end, and we extracted the fibres that started and ended inside the GM mask.
-Streamlines : streamline tractography was performed using the fibre peaks generated in the reconstruction step. Tractography was done using the FACT algorithm implemented in MRTRix in order to generate 100'000 fibres. Fibres were then extended by half a voxel in both ends, and we extracted the fibres that started and ended inside the GM.
-Probabilistic tractography : probabilistic tractography was computed on the Spherical Harmonics representation of the ODF from the Reconstruction step. We used the iFOD2 algorithm implemented in MRTrix to generate 100'000 fibres. We then extended the fibres by half a voxel on both ends, and extracted a set of fibres that started and ended inside the GM mask.

Superset (step6) -17_2
-The three tractograms above were concatenated after clustering the probabilistic tractography and streamline outputs using Dipy's QuickBundle tool, with a distance threshold of 4mm. The concatenated tractogram was saved as Superset2.trk 1. Eddy Correction (FSL) 35 2. Rotation of b-vectors using the transformation from the FSL eddy correction output. 3. EPI correction (FUGUE) 35 o This was done using an in-house EPI correction code named 3 Dimensional Mutal Information (DMI) by Alex Leow 49 4. Hough Tractography 50,51 o Seeds: 29588 o White matter mask (threshold>0.2)  Grey and white matter segmentation and intensity bias correction of T1 image using unified segmentation approach 52 as implemented in SPM12. White matter probability map was then binarized using a threshold of 0.2 and down-sampled to 2mm to create the white matter mask. o Hough tractography algorithm: briefly, this is an global probabilistic fiber tracking approach based on the Hough transform. This allows sorting the 3D curves in the volumes by computing a score from diffusion images. The curves with higher scores represent the potential anatomical connections. o Result: ISMRM_2015_challenge_Pipeline1_track.trk 5. Fiber Length Threshold: 18mm -109mm o To remove small incorrect fibers from the periphery Pipeline 2: 19_2 1. Eddy Correction (FSL) 35 2. Rotation of b-vectors using the transformation from the eddy correction tool. 3. EPI correction (FUGUE) 35 4. Up-sampling to 1mm resolution using FSL's FLIRT. 5. Hough Tractography 50,51 o Seeds: 50000 o GM+WM binary mask  Grey (GM) and white matter (WM) segmentation and intensity bias correction of T1 image using unified segmentation approach 52 as implemented in SPM12. Grey and white matter probability maps were then binarized using a threshold of 0.2 and down-sampled to 2mm. The mask was generated as GM + WM. o Result: ISMRM_2015_challenge_Pipeline2_track.trk 5. Length Threshold: 18mm -109mm 6. ROI filtering to also remove longer incorrect fibers from the outside of the brain Pipeline 3: 19_0 1. Eddy Correction (FSL) 35 2. Rotation of b-vectors using the transformation from the eddy correction tool. 3. EPI Correction (FUGUE) 35 4. Registration to T1 and up-sampling to 1mm (flirt, normmi -12dof) 53 interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed using different combinations of preprocessing from FSL and ExploreDTI as described above to identify and the best reconstruction for each anatomical tract.
Tract_file2 : DTI_clean -20_1 -Denoising was applied to diffusion data using Overcomplete Local PCA 3 . -TOPUP/EDDY correction from FSL and RESTORE from ExploreDTI 28 preprocessing were included (motion distortion, eddy_current, EPI geometrical distortions, spike correction). -Data was processed using deterministic DTI tractography in ExploreDTI, with angle threshold of 45 degrees and FA stopping threshold of 0.15 -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed using different combinations of preprocessing from FSL and ExploreDTI as described above to identify and the best reconstruction for each anatomical tract. -An additional automatic filter was applied to remove spurious streamlines stopping and starting in either white matter regions, or starting and stopping in CSF.
Tract_file3 : DTI_cleanest -20_6 -Denoising was applied to diffusion data using Overcomplete Local PCA 3 . -TOPUP/EDDY correction from FSL and RESTORE from ExploreDTI 28 preprocessing were included (motion distortion, eddy_current, EPI geometrical distortions, spike correction). -Data was processed using deterministic DTI tractography in ExploreDTI, with angle threshold of 45 degrees and FA stopping threshold of 0.15 -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed using different combinations of preprocessing from FSL and ExploreDTI as described above to identify and the best reconstruction for each anatomical tract. -An additional automatic filter was applied to remove spurious streamlines stopping in deep white matter regions, or starting and stopping in CSF.

Tract_file1 : baseline -20_9
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 600 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 700 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 800 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations ranging from 600 -800 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed with and without denoising, using both conservative and liberal spherical deconvolution settings to identify and the best reconstruction for each anatomical tract.

Tract_file1 : baseline_clean -20_10
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 600 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines starting and stopping in either deep white matter regions or CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 700 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines starting and stopping in either deep white matter regions or CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 800 to optimise angle resolution and noise stability.
-Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines starting and stopping in either deep white matter regions or CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations ranging from 600 -800 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed with and without denoising, using both conservative and liberal spherical deconvolution settings to identify and the best reconstruction for each anatomical tract. -An additional automatic filter was applied to remove spurious streamlines starting and stopping in either deep white matter regions or CSF.
Tract_file1 : baseline_cleanest -20_8 -TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 600 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines stopping in deep white matter regions, or starting and stopping in CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI 28 preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 700 to optimise angle resolution and noise stability.
-Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines stopping in deep white matter regions, or starting and stopping in CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations = 800 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -An additional automatic filter was applied to remove spurious streamlines stopping in deep white matter regions, or starting and stopping in CSF.
-TOPUP/EDDY correction from FSL follewed by standard ExploreDTI preprocessing including (motion distortion, eddy_current and EPI geometrical distortions) -Data was processed using Spherical Deconvolution with a damped Richardson-Lucy algorithm was applied to the data as described in Dell'Acqua et al. 55 using fibre-response model alpha = 2.5 and number of iterations ranging from 600 -800 to optimise angle resolution and noise stability. -Tractography was performed as described in Dell'Acqua et al. 56 , with absolute FOD threshold of 0.05 and angle threshold of 45degrees -Semi-automatic pruning using manually selected regions of interest was applied to all datasets to filter out spurious streamlines and identify major white matter bundles. -The data was processed with and without denoising, using both conservative and liberal spherical deconvolution settings to identify and the best reconstruction for each anatomical tract. -An additional automatic filter was applied to remove spurious streamlines stopping in deep white matter regions, or starting and stopping in CSF.

Supplementary Note 2: Extended tractography experiments
Tractography on ground truth directions We ran additional tractography experiments directly on the ground truth field of fiber orientations. Those orientations were found using the ground-truth bundles, by extracting the main directions of bundles going through each voxel. The experiments were performed using multiple resoltuions of the ground truth vector field (2mm, 1.75mm, 1.5mm, 1.25mm, 1.0mm, 0.75mm and 0.5mm) and two independent implementations of deterministic streamline tractography methods (GT 1 and GT 2 ). GT 1 is an in-house development of the Sherbrooke Connectivity Imaging Lab. GT 2 is implemented and openly available in MITK (mitk.org). Both tractography algorithms used the orientations as input. The two methods used the following parameters for all datasets: GT1: -100k streamlines seeded in the white matter -Step size 0.2mm -Max angle between 2 steps: 45 degrees -Min/max lenght of streamlines: 10mm/300mm GT2: -100k streamlines seeded in the white matter -Step size 0.5 * voxel size -Max angle between 2 steps: 45 degrees -Min/max lenght of kept streamlines: 20mm to 200mm