Kurt G Schilling
email: kurt.g.schilling.1@vumc.org

Chantal M W Tax

Francois Rheault

Bennett Landman

Adam Anderson

Maxime Descoteaux

Laurent Petit

Prevalence of white matter pathways coming into a single diffusion MRI voxel orientation: the bottleneck issue in tractography

Keywords: tractography, white matter, fiber pathways, bottleneck, tractometry, crossing fibers

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Two paths diverged from a single orientation, And streamlines could not travel both Else it be a false positive, long it stood And looked down one as far it could

To which cortex should it approach? Diffusion MRI fiber tractography is currently the only tool to map the long-range structural brain connectivity in vivo. However, there are a number of limitations and ambiguities that affect the ability of tractography to accurately map the connections of the brain [START_REF] Schilling | Challenges in diffusion MRI tractography -Lessons learned from international benchmark competitions[END_REF]. At the voxel level, significant attention has been given to the "crossing fiber problem" [START_REF] Alexander | Mathematics of crossing fibers[END_REF][START_REF] Tournier | The Biophysics of Crossing Fibers[END_REF][START_REF] Tuch | Diffusion MRI of complex neural architecture[END_REF]. This problem typically refers to the situation when two or more differently oriented fiber bundles are located in the same imaging voxel, which causes a partial volume effect that can lead to ambiguous or incorrect estimates of fiber orientation and subsequent failure of tractography [START_REF] Wheeler-Kingshott | About "axial" and "radial" diffusivities. Magnetic resonance in medicine[END_REF]. Crossing fibers have been shown to occur in a majority of the voxels in the brain [START_REF] Jeurissen | Investigating the prevalence of complex fiber configurations in white matter tissue with diffusion magnetic resonance imaging[END_REF], and for the last decade the crossing fiber problem has been cited as the major limitation that diffusion tractography faces, with a vast number of algorithms [START_REF] Daducci | Quantitative comparison of reconstruction methods for intra-voxel fiber recovery from diffusion MRI[END_REF] and papers referring to this problem [START_REF] Schilling | Challenges in diffusion MRI tractography -Lessons learned from international benchmark competitions[END_REF][START_REF] Alexander | Mathematics of crossing fibers[END_REF][START_REF] Tournier | The Biophysics of Crossing Fibers[END_REF][START_REF] Tuch | Diffusion MRI of complex neural architecture[END_REF][START_REF] Schilling | Histological validation of diffusion MRI fiber orientation distributions and dispersion[END_REF][START_REF] Schilling | Can increased spatial resolution solve the crossing fiber problem for diffusion MRI? NMR in biomedicine[END_REF][START_REF] Knösche | Validation of tractography: Comparison with manganese tracing[END_REF][START_REF] Mani | Acceleration of high angular and spatial resolution diffusion imaging using compressed sensing with multichannel spiral data[END_REF][START_REF] Zhang | Axon diameter mapping in crossing fibers with diffusion MRI[END_REF][START_REF] Tournier | Direct estimation of the fiber orientation density function from diffusion-weighted MRI data using spherical deconvolution[END_REF][START_REF] Rheault | Common misconceptions, hidden biases and modern challenges of dMRI tractography[END_REF][START_REF] Schilling | Limits to anatomical accuracy of diffusion tractography using modern approaches[END_REF][START_REF] Jeurissen | Diffusion MRI fiber tractography of the brain[END_REF][START_REF] Poulin | Tractography and machine learning: Current state and open challenges[END_REF][START_REF] Tax | Localizing and characterizing single fiber populations throughout the brain[END_REF][START_REF] Descoteaux | Deterministic and probabilistic tractography based on complex fibre orientation distributions[END_REF][START_REF] Donahue | Using Diffusion Tractography to Predict Cortical Connection Strength and Distance: A Quantitative Comparison with Tracers in the Monkey[END_REF][START_REF] Behrens | Probabilistic diffusion tractography with multiple fibre orientations: What can we gain?[END_REF][START_REF] Behrens | Characterization and propagation of uncertainty in diffusion-weighted MR imaging[END_REF][START_REF] Dyrby | Validation strategies for the interpretation of microstructure imaging using diffusion MRI[END_REF][START_REF] Dyrby | Validation of in vitro probabilistic tractography[END_REF]. The identification and characterization of this problem led to a fundamental paradigm shift in diffusion processing, moving the field beyond classical diffusion tensor imaging, and has led to the development of a number of algorithms capable of resolving crossing fibers [2-4, 13, 25, 26].

A more recently described limitation of fiber tractography is the "bottleneck problem" [START_REF] Maier-Hein | The challenge of mapping the human connectome based on diffusion tractography[END_REF]. In contrast to the crossing fiber problem, bottlenecks occur at a global level when multiple fiber populations converge towards a narrow region, temporarily aligning and sharing the same orientation and trajectory, before re-emerging from the bottleneck region [START_REF] Schilling | Challenges in diffusion MRI tractography -Lessons learned from international benchmark competitions[END_REF][START_REF] Rheault | Common misconceptions, hidden biases and modern challenges of dMRI tractography[END_REF]. Current tractography algorithms cannot adequately choose the correct pathway upon reemerging, which leads to generation of a potentially large number of false positive pathways [START_REF] Maier-Hein | The challenge of mapping the human connectome based on diffusion tractography[END_REF], and limits the ability to use tractography as a tool to explore potential connections and fiber pathways of the brain. While significant efforts have gone into solving the crossing fiber problem, the bottleneck problem has received far less attention. A thorough characterization and investigation of bottleneck locations and prevalence may highlight the extent of this problem, and much like the crossing fiber problem, cause a paradigm shift in tractography in order to solve this issue.

In this work, we utilize well-known and well-characterized white matter fiber bundles extracted using automated tools, to quantify how often they overlap within the same imaging voxels, and also how often the overlap occurs within the same voxel while also sharing the same dominant orientation. These locations represent known bottleneck regions for tractography, and indicate areas in the brain where a number of white matter pathways converge, and where tractography may lead to incorrect or erroneous estimates of brain connectivity.

Nomenclature

Here, we aim to clarify nomenclature that will be used in this study to describe our methodology and results. First, a bundle, or fiber bundle, is a group of streamlines that is created from a diffusion MRI dataset and is intended to represent a specific white matter pathway of the brain. Bundles contain streamlines with start and end points generally belonging to the same brain territories, respectively. In this study, we create bundles using two common white matter atlases and tractography dissection techniques that are informed by prior anatomical knowledge and contain pathways for which there is extensive evidence of their existence. Thus, in this study we are analyzing only anatomically plausible bundles.

Next, a voxel represents orientational information in threedimensional space. In MRI, the size of voxels is on the order of millimeters, whereas axons of the brain have diameters on the scale of micrometers, and a single voxel can contain hundreds of thousands of axons. The directionality of axons within a voxel can be summarized by the fiber orientation distribution (FOD). The FOD is a continuous function over a sphere and can be visualized as a 3D directional histogram where peaks, or local maxima, are assumed to point parallel to the direction of axons. The FOD can be segmented into discrete elements based on peaks, or lobes, that are considered to be representative of a specific orientation of a set of axons within each voxel. These fiber elements are referred to as fixels [START_REF] Raffelt | Investigating white matter fibre density and morphology using fixel-based analysis[END_REF], and are parameterized by the mean Because the UF and IFOF diverge in Voxel C, this is an example of a multi-orientation voxel and multi-bundle voxel, with one fixel classified as a single-bundle fixel and the other a multibundle fixel. Voxel D highlights the fanning of the IFOF, which results in a multi-orientation voxel and single-bundle voxel, and both fixels are single-bundle fixels. Finally, both the IFOF and UF pass through voxel E following the same orientation, thus voxel E is a single-orientation voxel but multi-bundle voxel, and also a multi-bundle fixel. This fixel, and thus also this voxel, represents a bottleneck for tractography. orientation of fibers within the lobe of the FOD. Thus, there can be multiple fixels in a single voxel, with the advantage that we can now assign a specific property or index to each fixel, or orientation, within a voxel.

For this study, it is necessary to clarify or define six classifications, written out in Table 1 and displayed as a cartoon in Figure 1. If a FOD has only a single peak, or local maxima, we classify it as a single-orientation (fiber) voxel, and if it has greater than one peak we call this a multiorientation voxel [START_REF] Jeurissen | Investigating the prevalence of complex fiber configurations in white matter tissue with diffusion magnetic resonance imaging[END_REF]. In contrast to simply counting peaks, we also count the number of bundles passing through the same voxel, then characterizing as a single-bundle voxel and multi-bundle voxel. Similarly, we count the number of bundles associated with each fixel and, characterize the fixel as a single-bundle fixel, or multi-bundle fixel. As described above, a fixel is usually used to describe a single fiber bundle element. However, we hypothesize that a single fixel may be associated with several bundles, thus creating bottlenecks for tractography. A bottleneck 'region' then is a cluster of coherent multi-bundle fixels.

Methods

Data

We utilized data from 25 healthy subjects in the Human Connectome Project (HCP) S1200 release [START_REF] Van Essen | The WU-Minn Human Connectome Project: an overview[END_REF]. The HCP protocol (custom 3T Siemens Skyra) included T1-weighted images acquired using a 3D MPRAGE sequence (TE = 2.1 ms, TR = 2400 ms, flip angle = 8 deg, FOV = 224×224mm, acquisition, voxel size = 0.7mm isotropic). Diffusion images were acquired using a single-shot EPI sequence, and consisted of three b-values (b = 1000, 2000, and 3000 s/mm 2 ), with 90 directions (and 6 b=0 s/mm 2 ) per shell (TE = 89.5 ms, TR = 5520 ms, slice thickness = 1.25 mm, flip angle = 78 degrees, FOV = 210*180, voxel size = 1.25mm isotropic). Data pre-processing included correction for susceptibility distortions, subject motion, and eddy current correction.

Processing

All processing was performed with the MRTrix3 software package [START_REF] Tournier | MRtrix3: A fast, flexible and open software framework for medical image processing and visualisation[END_REF]. Multi-shell, multi-tissue constrained spherical deconvolution (dwi2fod) was used to estimate the white matter FOD for each subject using a group averaged white matter response function computed using the "dhollander" algorithm [START_REF] Dhollander | Unsupervised 3-tissue response function estimation from single-shell or multi-shell diffusion MR data without a coregistered T1 image[END_REF]. Next, we generated a study-specific unbiased FOD template (population_template) for population-based analysis.

From both the population FOD template and subject-specific FODs, we extracted peaks (sh2peaks) and fixels (fod2fixel), with a threshold of 0.1 times the maximum peak amplitude to remove spurious peaks. We then counted the number of peaks (i.e., lobes of the FOD) and defined those voxels which have only a single peak to be a single-orientation voxel, and those with greater than one peak to be multi-orientation voxels.

Bundle segmentation

We utilized two common, automated, pipelines for white matter bundle extraction, Recobundles [START_REF] Garyfallidis | Recognition of white matter bundles using local and global streamline-based registration and clustering[END_REF] and TractSeg [START_REF] Wasserthal | Combined tract segmentation and orientation mapping for bundle-specific tractography[END_REF][START_REF] Wasserthal | TractSeg -Fast and accurate white matter tract segmentation[END_REF]. These are both informed by prior anatomical knowledge in order to generate bundles representative of wellcharacterized, and well-validated, white matter pathways of the brain. All analysis is performed separately for both techniques to show that results generalize across slight deviations in the number and definitions of white matter pathways, and the techniques used to extract them.

Recobundles is based on an atlas of 78 bundles [START_REF] Yeh | Populationaveraged atlas of the macroscale human structural connectome and its network topology[END_REF], although 12 bundles are cranial nerves outside of the cerebrum and brainstem. Whole brain tractography was performed using anatomically constrained probabilistic tractography with the iFOD2 propagation algorithm (tckgen), to generate 25 million streamlines, which were filtered based on the diffusion signal using the SIFT algorithm [START_REF] Smith | SIFT: Spherical-deconvolution informed filtering of tractograms[END_REF] (tcksift) to 2 million streamlines. Bundle recognition was performed following streamline linear registration to the HCP842 [START_REF] Yeh | Populationaveraged atlas of the macroscale human structural connectome and its network topology[END_REF] bundle template (dipy_slr) and bundle recognition using default parameters of the RecoBundles algorithm [START_REF] Garyfallidis | Recognition of white matter bundles using local and global streamline-based registration and clustering[END_REF].

TractSeg [START_REF] Wasserthal | TractSeg -Fast and accurate white matter tract segmentation[END_REF] is a tool based on convolutional neural networks that is trained to create tract orientation maps and segmentations of end regions, which can be used to perform probabilistic bundle-specific tractography [START_REF] Wasserthal | Combined tract segmentation and orientation mapping for bundle-specific tractography[END_REF]. We implemented the processing pipeline provided at (https://github.com/MIC-DKFZ/TractSeg) with the MRTrixderived FODs as input, in order to generate 72 bundles per subject.

Assigning bundles to voxels and fixels Figure 2 visualizes the procedure used to assign bundles to voxels and bundles to fixels. For each bundle (Figure 2,A; N=78 RecoBundles; N=72 TractSeg), a tract density map was created (tckmap) by counting the number of streamlines coursing through each voxel (Figure 2,B). Density maps were thresholded at 5% of the maximum streamline density per bundle in order to create binary segmentations indicating the voxel-wise profile of each bundle. Next, a fixel density map (Figure 2,C) was created (tck2fixel) by counting the number of streamlines that are most closely aligned to each fixel. Again, fixel-density maps were thresholded at 5% of the maximum density in order to create binary segmentations indicating the fixel-wise profile of each bundle. Now, the number of known bundles per voxel, as well as number of known bundles per fixel, can be counted. A single fixel that is associated with multiple bundles represent a potential bottleneck. Voxel-based and fixel-based analysis is performed both at the level of the individual in subject-space, and also at the population-level in template-space.

Figure 2. Assigning bundles to voxels and fixels. Each segmented white matter bundle (A) was assigned to each voxel (B) based on thresholding of the density map, as well as assigned to each fixel (C) within a voxel based on matching orientations of streamlines. This allowed us to query both the number of known bundles per voxel, and number of known bundles per fixel. White matter bundles were derived from TractSeg (N=72 bundles) and Recobundles (N=78 bundles).

Results

Bottleneck Prevalence

Investigating fixels throughout the white matter, we first ask what is the prevalence of multi-bundle fixels? Figure 3 shows maps of the number of bundles assigned to each fixel, visualized in coronal, sagittal, and axial views, for both TractSeg bundles (top) and Recobundles (bottom). Most noticeably, most fixels in the white matter, in all orientations, are associated with multiple bundles. In fact, many regions have groups of oriented fixels associated with 7+ unique bundles. Figure 4 quantifies the number of bundles assigned to each individual fixel. The results confirm the qualitative observation that most fixels contain multiple bundles converging in a given orientation, with greater than 50% of fixels in Recobundles and greater than 70% of fixels in TractSeg containing greater than a single fiber population. In general, TractSeg bundles show a higher prevalence of multi-bundle fixels (i.e., bottlenecks) than Recobundles. In summary, a majority of fixels in the brain that contain known fiber bundles act as bottleneck regions for tractography.

Bottleneck Locations

Next, we ask where do the most important bottleneck regions occur? These are regions in which groups of fixels with similar orientation exhibit a convergence of the largest number of pathways. Highlighted bottleneck regions, e.g., 7 or more bundles, are visualized in Figures 5678. In all Figures, the fixels are color-coded by the number of fibers traversing through each orientation, and exemplar bundles converging in each region are displayed.

The first displayed region with the highest number of bottlenecks (Figure 5) is located in the deep white matter of the occipital lobe oriented in the anterior-posterior direction, covering a wide dorsal-ventral expanse and including the stratum sagittale. A large number of individual WM pathways, each with unique starting and/or ending connections, converge through this region with the same orientation, including the posterior part of the inferior fronto-occipital fasciculus (IFOF), inferior longitudinal fasciculus (ILF), middle longitudinal fasciculus (MdLF), parieto-thalamic and occipito-thalamic (optic radiations, OR) connections, parieto-striatal and parieto-occipital pontine tract (POPT), and several subdivisions or segmentations of the striato-cortical connections, and splenium of the corpus callosum. Thus, while most pathways terminate throughout the occipital lobe, these fibers can project to sub-cortical nuclei, temporal or frontal lobes, or to the contralateral hemisphere as commissural fibers.

Figure 6 further highlights the bottleneck problem in this region, and illustrates how different tracts may share a similar location AND orientation, yet have unique start and end points. We filtered the previously described bundles using a single region of interest (a 2x2x2 voxel cube i.e., a 2.5mm isotropic region), and show just the streamlines from each bundle that traverse this area (top row). While the full extent of each pathway does not traverse this region, a large and coherent subset of each bundle does, all oriented in the anterior-posterior direction. To simplify the illustration, a representative streamline is shown for each filtered bundle (bottom row), exemplifying the bottleneck problem: there is a large combinatorial number of possible pathways that traverse through this voxel following this single welldefined orientation.

The second main bottleneck region is the convergence of superior-inferior oriented fibers converging and traversing throughout the brainstem, from the mid-brain to the medulla (Figure 7). This includes a number of ascending and descending fibers projection pathways, corticopontine fibers arising from the cortex, and cerebellar tracts. Again, fibers traversing through this region all share the same dominant orientation, yet end throughout the extent of the cortex, subcortex, spinal cord, and cerebellum.

The third region with highest number of bottlenecks occurs in the superior-inferior oriented fibers of the posterior limb of the internal capsule (Figure 8). This region contains pathways such as the corticospinal tract, frontal and parietal pontine fasciculi, striato-postcentral and striato-precentral bundles, superior thalamic radiations toward the parietal, precentral, and postcentral cortices. While most of these fibers project from the mid-brain and nuclei, projections cover the expanse of the parietal and frontal cortices, with many projecting onward dorsally towards the superior frontal gyrus.

Single-bundle and multi-bundle voxels

Rather than assessing fixels, we can also ask what is the prevalence of single-bundle and multi-bundle voxels, and where do known bundles overlap? Thus, this overlap can contain both bundles oriented in the same direction, and different directions within a voxel. Figure 9 quantifies the number of bundles per voxel for both Recobundles and

. Illustration of the bottleneck in the anterior-posterior oriented white matter of the occipital lobe. Fiber bundles from Figure 5 (same color scheme) were filtered to select only streamlines (top) traversing a small 2x2x2 voxel region of interest (arrow). A single representative streamline from each sub-bundle is also shown (bottom). This example emphasizes that streamlines belonging to many fiber bundles may traverse through the same small region, in the same orientation.

TractSeg bundles, and visualizes bundle overlap in the template space. It is clear that a majority of voxels in the white matter contain multiple, overlapping bundles. In fact, a number of voxels in the brain again contain as many as 7 or more unique overlapping bundles, with overlap of these bundles frequently occurring in regions that parallel the bottleneck regions -the centrum semiovale and posterior corona radiata, and also the posterior limb and retrolenticular limb of the internal capsule. Thus, most white matter voxels (>65-75%) contain overlap of multiple bundles.

Single-orientation and multi-orientation voxels

Finally, we ask what is the prevalence of single-orientation and multi-orientation voxels (i.e., the prevalence of the "crossing fiber problem"), and where do single-and multiorientation voxels occur? Figure 10 quantifies and visualizes the prevalence of multi-orientation voxels in both an individual subject and across the population. In agreement with previous literature, our results show that a majority of voxels in the white matter (>60%) contain multiorientation voxels. These voxels are prevalent throughout the entire white matter, with more complex (for example >2 peaks) crossings in the centrum semiovale and cerebellum. Voxels with only a single peak are prevalent in the corpus callosum and internal capsule, as well as near the crowns of various gyri.

One particularly interesting finding is that many of the identified bottlenecks above are often associated with singleorientation voxels -the internal capsule, mid-brain, and less frequently in the deep white matter of the occipital lobe.

Discussion

In this study, we investigated the prevalence and locations of bottleneck regions in the brain, which present obstacles in our ability to build anatomically correct maps of the human brain using diffusion tractography. We find that not only do a majority of voxels contain multiple bundles, but also that a majority of individual FOD peaks within a voxel (i.e., fixels) are associated with multiple bundles. As much as 50-70% of fixels in the brain contain multiple fiber bundles traversing through them. This is based on the use of "well-known" bundles, representing only a lower-bound that can only go up as we understand and map all existing bundles in the brain. The convergence of bundles into a nearly parallel funnel, and subsequent convergence, may lead to a combinatorial number of possible pathways that tractography algorithms may choose to take, leading to the generation of possibly anatomically non-existent pathways. While this may not present a problem for bundlespecific tractography with the use of manually placed priors [START_REF] Chamberland | Active delineation of Meyer's loop using oriented priors through MAGNEtic tractography (MAGNET)[END_REF], bundle templates [START_REF] Hansen | Pandora: 4-D white matter bundle population-based atlases derived from diffusion MRI fiber tractography[END_REF][START_REF] Rheault | Bundle-specific tractography with incorporated anatomical and orientational priors[END_REF], or machine learning [35, 41,

Figure 7. Bottleneck region in the superior-inferior oriented white matter of the brain-stem (arrow) contains a large number of white matter bundles with unique starting and ending connections. Hot-cold colormap ranges from 0 to 7 bundles. Pathways (derived from Recobundles) from top to bottom, left to right: CST, FPT, LL, MLL, CTT, OPT, TPT, SCP (for full names see acronyms at end of document).

42], mapping the entirety of the human connectome ultimately strives to map all of the unique bundles of the brain. The identified bottleneck regions clearly present obstacles to this mapping when the true connections are not known a priori. From a whole-brain connectome, streamlines which pass through these fixels without explicit knowledge of their existence should be suspected to be false positive connections generated by the process. Over and above problems created in connectomics studies, this also highlights the problems of proposing the existence of a new or unique pathway from diffusion MRI alone, even if the pathway is reproducible across scans and subjects. We posit

Figure 8. Bottleneck region in the superior-inferior oriented white matter of the internal capsule (arrow) contains a large number of white matter bundles with unique starting and ending connections. Hot-cold colormap ranges from 0 to 7 bundles. Pathways (derived from TractSeg) from left to right, top to bottom: T_PREM, T_PAR, STR, ST_PREF, ST_POST, POPT, FPT, CST (for full names see acronyms at end of document).

Figure 9. Many voxels in the white matter contain multiple known bundles. Prevalence of voxels with 1 to 7+ bundles is quantified for Recobundles (A) and visualized in template space (B), and also quantified for TractSeg bundles (C) and visualized overlaid in template space (D).

Note that many voxels contain 0 bundles (i.e, are not associated with known bundles in our atlas) and are thus not quantified. that orthogonal information in the form of blunt dissection, tracers in animals, and alternative contrasts is necessary for inferences which demand highly specific tractography.

This bottleneck problem has further implications on quantitative tractography, or tractometry, where the alongtract profile of measures along multiple tracts allows a comprehensive characterization of white matter. Most methods map voxel-wise values along the tract (such as FA) [START_REF] Yeatman | Tract profiles of white matter properties: automating fiber-tract quantification[END_REF], which we know are affected by crossing fibers, while recent studies map fixel-wise values (such as apparent fiber density), along the tract [START_REF] Chamberland | Detecting microstructural deviations in individuals with deep diffusion MRI tractometry[END_REF][START_REF] Chamberland | Dimensionality reduction of diffusion MRI measures for improved tractometry of the human brain[END_REF]. Yet, we show here that these measures are still not yet truly specific to one bundle. Finally, studies that use global methods of filtering or quantification [START_REF] Smith | SIFT: Spherical-deconvolution informed filtering of tractograms[END_REF][START_REF] Daducci | COMMIT: Convex optimization modeling for microstructure informed tractography[END_REF] can map an average streamlinespecific estimate of diffusion or relaxometry along the tract [START_REF] Barakovic | Resolving bundle-specific intra-axonal T2 values within a voxel using diffusion-relaxation tract-based estimation[END_REF], however, these will still be affected by partial volume at different points along the streamline profile.

We additionally confirm findings from previous studies which indicate that most voxels have multiple peaks in their FOD. Previous studies have estimated anywhere from 60-90% of the brain contains "crossing fibers" [START_REF] Descoteaux | Deterministic and probabilistic tractography based on complex fibre orientation distributions[END_REF][START_REF] Behrens | Probabilistic diffusion tractography with multiple fibre orientations: What can we gain?[END_REF], estimates which vary with signal to noise ratio, diffusion sensitivity, and image resolution [START_REF] Jeurissen | Investigating the prevalence of complex fiber configurations in white matter tissue with diffusion magnetic resonance imaging[END_REF][START_REF] Schilling | Can increased spatial resolution solve the crossing fiber problem for diffusion MRI? NMR in biomedicine[END_REF]. Our estimated fraction is at the low-end of literature values, but we utilize multi-shell deconvolution as opposed to previous work, which may result in more regularized reconstruction and less spurious peaks. Regardless, this suggests that most voxels have more than one bundle traversing its location. Over and above this, we find that most peaks within a voxel have more than one bundle traversing in its direction. Thus, even solving the crossing fiber problem does not solve the bottleneck problem that may cause a large number of ambiguous, false positive pathways. Here, we propose that characterizing and describing the prevalence of this problem should lead to the development of methods which may alleviate or overcome these obstacles, much like that done for the crossing fiber problem. It is likely that advances in anatomical knowledge in combination with innovation in streamline generation and innovation in streamline filtering, will be needed to mitigate this problem.

Bottleneck locations

We have described the highest bottleneck regions in this study. Specifically, we highlighted the deep white matter of the occipital lobe, the brainstem, and the internal capsule. These regions included a number of association, projection, and commissural fibers, all with unique trajectories and fundamentally different structural connections. While these were the most visually apparent 'hot spots', it is clear that a majority of fixels in the brain are associated with multiple white matter fiber pathways. From a tractography perspective, these regions may cause ambiguous connectivity estimates, yet, anatomically, these locations may have significant functional relevance, representing the intersection or merging of the many anatomo-functional highways of the brain.

Importantly, these bottlenecks are almost certainly an underestimation of the true prevalence and extent of this problem. We choose segmentation techniques which reconstruct only known anatomical pathways of the brain (72 and 78 bundles, respectively) for which there is broad agreement on their existence. Several other segmentation techniques exist which suggest the existence of a much greater number of unique pathways in the brain, however we have chosen to perform a conservative estimation, in addition to a conservative thresholding of density and/or streamlines to highlight the prevalence of this problem. Thus, these numbers are derived from a conservative estimate of the complexity of the brain, and represent the bare minimum of the number of convergent bundles. However, the absolute quantification itself may be biased towards the pathways from the chosen techniques and the atlases these are based upon. Additionally, while the bundles defined by these techniques have proven reliability, they are not themselves immune to the problems posed by bottlenecks. It is possible that improvements in tracking will refine the definitions of the bundles, as well as identify new pathways; while this may change estimates of the number of bundles in a given fixel, it won't eliminate the basic feature of the coincidence of bundles in common pathway segments.

Adding to known atlases

While the primary aim of this study was to characterize bottlenecks, in which >1 bundle passes through a fixel, we made the interesting ancillary observation (Supplementary Figure 1) that many fixels within the white matter were not associated with any bundles defined in our utilized atlases. We note that these 'zero-bundle' fixels were not included in the analysis because there is a lack of information regarding these orientations. While many of these could be spurious peaks or isolated voxels in white matter, larger expanses of coherently oriented FODs could be regions which may represent underexplored white matter pathways which can eventually be added to our repertoire and collection of bundles. While some regions are easily explained, for example due to a lack of TractSeg bundles in the cerebellum, other regions occur both along and across several gyral blades, as well as the relatively underexplored system of Ufibers and local association fibers. Further reasons for this could be the thresholding of densities and streamlines and parameter configurations for the bundle segmentation techniques. Additional atlases or bundle segmentation procedures may include pathways through many of these regions, which would solve the missing-bundle problem, but would likely also increase the prevalence of bottlenecks. This is a potential limitation of the current study -the choice of atlases (bundle segmentation procedures). We have purposefully chosen to only include pathways for which there is broad agreement on their existence, and techniques which have been validated and well-utilized in the field. There would also be the potential to define pathways and bottleneck regions using whole-brain connectivity, or atlases which are derived from clusters or reproducibility of large datasets [START_REF] Siless | AnatomiCuts: Hierarchical clustering of tractography streamlines based on anatomical similarity[END_REF][START_REF] Guevara | Automatic fiber bundle segmentation in massive tractography datasets using a multi-subject bundle atlas[END_REF][START_REF] Ros | Atlas-guided cluster analysis of large tractography datasets[END_REF], however, these are potentially confounded by the bottleneck problem itself (i.e., contain false positive, yet reproducible, pathways), and an extensive comparison of algorithms and segmentation methods is beyond the scope of the current study.

"Single Fiber Populations" and microstructure

In addition to tractography, the diffusion signal that results from a single fiber population (i.e., the fiber response function) has applications towards tissue microstructural modeling. The response function can be used to estimate the FOD and tissue microstructural properties including diffusivities, compartment sizes, and orientation dispersions. Typically, this response function is derived from studying regions of low complexity, specifically, regions that are considered single fiber populations and contain only a single peak in the FOD [START_REF] Tournier | Robust determination of the fibre orientation distribution in diffusion MRI: non-negativity constrained super-resolved spherical deconvolution[END_REF][START_REF] Dhollander | Unsupervised 3-tissue response function estimation from single-shell or multi-shell diffusion MR data without a coregistered T1 image[END_REF][START_REF] Tax | Recursive calibration of the fiber response function for spherical deconvolution of diffusion MRI data[END_REF]. Here, however, we can see that even if a voxel or region contains a dominant orientation and high anisotropy, a majority of these regions are composed of multiple, distinct fiber pathways, that may have varying densities, sizes, and distributions of axons. For example, if we were to use the traditional definition of a 'single fiber population' (i.e., our single-orientation voxel) we would find that only 35% and 27% contain just a single fiber bundle passing through them (Supplementary Figure 2). Thus, these so-called single fiber regions are very often multi-bundle regions. Truly single-fiber and single-bundle regions are rare, even with our limited selection of known bundles used in this study, and in our case, occur in the cingulum and specific gyral blades (Supplementary Figure 3). While some works show that biological differences between fiber populations are negligible in the response function formulism [START_REF] Christiaens | On the need for bundlespecific microstructure kernels in diffusion MRI[END_REF], there is some evidence that the response functions do vary across pathways [START_REF] Tax | Localizing and characterizing single fiber populations throughout the brain[END_REF][START_REF] Schilling | Histologically derived fiber response functions for diffusion MRI vary across white matter fibers-An ex vivo validation study in the squirrel monkey brain[END_REF][START_REF] Howard | Joint modelling of diffusion MRI and microscopy[END_REF], which may lead to variation in estimates of FODs and subsequent microstructure.

Nomenclature

Here, we have also introduced slightly different nomenclature than past literature. As described above a single-orientation voxel has traditionally been called a single-fiber voxel, whereas a multi-orientation voxel has been called crossing-fiber voxel. Clearly, a voxel with only one orientation is not limited to only containing the presence of fibers from a single bundle, hence the new description. Additionally, a fixel traditionally refers to a "specific fiber bundle within a specific voxel" [START_REF] Raffelt | Investigating white matter fibre density and morphology using fixel-based analysis[END_REF][START_REF] Honnedevasthana Arun | Characterisation of white matter asymmetries in the healthy human brain using diffusion MRI fixel-based analysis[END_REF][START_REF] Genc | Development of white matter fibre density and morphology over childhood: A longitudinal fixel-based analysis[END_REF][START_REF] Raffelt | Connectivity-based fixel enhancement: Whole-brain statistical analysis of diffusion MRI measures in the presence of crossing fibres[END_REF]], yet we have shown again that a fixel (which is truly a segmented lobe, or orientation, from the FOD), is also not limited to a single specific fiber bundle, and in fact, likely contains axons from several bundles. We are not proposing to change the use and discussion of these elements throughout the field, but rather chose clarifying nomenclature to remove ambiguity in this specific study.

Conclusion

In this work, we investigated the prevalence of bottleneck regions, or where multiple white matter pathways of the brain converge and subsequently diverge. Our results indicate that most white matter contains multiple overlapping bundles, and individual orientations within a voxel are associated with multiple bundles. These findings have profound implications for tractography analysis which aims to map unknown connections across the brain, and strengthen the awareness of limitations or challenges facing these image processing techniques. institutional Research Chair in Neuroinformatics and NSERC Discovery.

Figure 1 .

 1 Figure 1. Nomenclature. Two bundles, the UF and IFOF, are used to highlight classifications of voxels (A-D), and fixels within the voxels. Voxels in A and B are examples of single-orientation voxel and single-bundle voxels and also single-bundle fixel. Because the UF and IFOF diverge in Voxel C, this is an example of a multi-orientation voxel and multi-bundle voxel, with one fixel classified as a single-bundle fixel and the other a multibundle fixel. Voxel D highlights the fanning of the IFOF, which results in a multi-orientation voxel and single-bundle voxel, and both fixels aresingle-bundle fixels. Finally, both the IFOF and UF pass through voxel E following the same orientation, thus voxel E is a single-orientation voxel but multi-bundle voxel, and also a multi-bundle fixel. This fixel, and thus also this voxel, represents a bottleneck for tractography.

Figure 3 .

 3 Figure 3. There is a high prevalence of bundles assigned to each fixel in the brain. Results are shown as vectors, colored by the number of associated bundles, and averaged across the population. TractSeg results are shown on top, Recobundles on bottom. Fixels with more than one bundle traversing through them represent bottleneck regions for tractography.

Figure 5 .

 5 Figure 5. Bottleneck region in the anterior-posterior oriented white matter of the occipital lobe (arrows) contains a large number of white matter bundles with unique starting and ending connections.Colormap ranges from 0 to 7+ bundles. Pathways (derived from TractSeg) from left to right, top to bottom: CC7, ST_par, OR, ST_OCC, POPT, MdLF, T_PAR, IFO (for full names see acronyms at end of document).

Figure 4 .

 4 Figure 4. Number of bundles assigned to fixels in the brain. Most fixels had greater than one bundle traversing through their designated orientation. Note that fixels which were assigned to 0 bundles are not shown. Y-axis is shown as a fraction of fixels. Error bars represent variation across the studied population.

Figure 6

 6 Figure 6. Illustration of the bottleneck in the anterior-posterior oriented white matter of the occipital lobe. Fiber bundles from Figure 5 (same color scheme) were filtered to select only streamlines (top) traversing a small 2x2x2 voxel region of interest (arrow). A single representative streamline from each sub-bundle is also shown (bottom). This example emphasizes that streamlines belonging to many fiber bundles may traverse through the same small region, in the same orientation.

Figure 10 .

 10 Figure 10. There is a high prevalence of multi-orientation voxels throughout the brain. Prevalence of voxels with one peak, two peaks, and three peaks is quantified for a single subject (A) and visualized overlaid on an anatomical image (B) and also averaged across the population (C) and visualized overlaid across the population template (D).Note that the number of peaks is discrete on a single subject but continuous when averaged across all subjects.

Acknowledgments

This work was supported by the National Institutes of Health under award numbers R01EB017230, and T32EB001628, and in part by ViSE/VICTR VR3029 and the National Center for Research Resources, Grant UL1 RR024975-01. CMWT was supported by a Sir Henry Wellcome Fellowship (215944/Z/19/Z) and a Veni grant (17331) from the Dutch Research Council (NWO). MD was supported by the

Data Availability

All resulting voxel-wise and fixel-wise bundle overlaps, for both TractSeg and Recobundles, is available upon request.

Appendix

The bundles resulting from each bundle-segmentation pipeline are given as a list below.

Recobundles


Francesco Cappuzzello
email: cappuzzello@lns.infn.it

Luis Acosta

Clementina Agodi

Carmen Altana

Paulina Amador-Valenzuela

Ismail Boztosun

Giuseppe A Brischetto

Salvatore Calabrese

Daniela Calvo

Vittoria Capirossi

Diana Carbone

Manuela Cavallaro

Efrain R Chávez Lomeli

Irene Ciraldo

Mauro Cutuli

Franck Delaunay

Haris Djapo

Canel Eke

Paolo Finocchiaro

Suna Firat

Maria Fisichella

Marcilei A Guazzelli Da Silveira

Aylin Hacisalihoglu

Felice Iazzi

Gaetano Lanzalone

Laura La Fauci

Roberto Linares

Nilberto H Medina

Mauricio Moralles

José R B Oliveira

Athena Pakou

Luciano Pandola

Horia Petrascu

Federico Pinna

Giuseppe Russo

Roberto B B Santos

Onoufrios Sgouros

Selçuk O Solakci

George Souliotis

Vasileios Soukeras

Alessandro Spatafora

Domenico Torresi

Salvatore Tudisco

Aydin Yildirin

Vinicius A B Zagatto

Recent experimental activity on heavy-ion induced reactions within the NUMEN project

The possibility to use a special class of heavy-ion induced direct reactions, such as double charge exchange reactions, is discussed in view of their application to extract information that may be helpful to determinate the nuclear matrix elements entering in the expression of neutrinoless double beta decay half-life. The strategies adopted in the experimental campaigns performed at INFN -Laboratori Nazionali del Sud are briefly described, emphasizing the advantages of the multi-channel approach to nuclear reaction data analysis.

Introduction

Charge exchange (CE) reactions are nuclear processes in which the individual balance of neutrons and protons inside the projectile and the target is modified, while keeping the overall number of nucleons unchanged. In the isospin formalism the CE promote isovector excitations in both the projectile and target. CE are thus selective tools to emphasize the role of isospin in both nuclear structure and reaction studies. Relevant aspects of single charge exchange (SCE) induced by light projectile have been discussed in several articles, among which excellent reviews are in ref. [START_REF] Alford | Nucleon charge-exchange reactions at intermediate energy[END_REF] [2] [START_REF] Taddeucci | The (p, n) reaction as a probe of beta decay strength[END_REF]. For heavy-ion induced CE we refer to the recent review [START_REF] Lenske | Heavy ion charge exchange reactions as probes for nuclear β-decay[END_REF]. An especially appealing aspect is the possibility to use CE to investigate the spinisospin-dependent response of nuclei. Recently, particular interest has been given to double charge exchange (DCE) processes [5] [6] [START_REF] Cappuzzello | The NUMEN Project: Shedding Light on Neutrinoless Double Beta Decay by Heavy-Ion Nuclear Reactions[END_REF] as they can give access to the nuclear response to second order isospin operators, which also promote double beta decays. SCE and DCE studies have a genuine nuclear physics interest, since they reflect a key aspect of nuclear dynamics. However, beyond that, there is an even broader impact of such investigations because of their relevance for nuclear beta decay, thus connecting strong and weak interactions. In this context, the NUMEN and NURE projects at INFN-LNS laboratory in Catania [START_REF] Cappuzzello | The NUMEN project: NUclear Matrix Elements[END_REF] [9] aim at exploring heavyion induced single and double charge exchange (DCE) reactions in a systematic approach. The main goal is to extract information on nuclear matrix elements for SCE and DCE processes and to identify the possible connections with double beta decay. Detailed information about these projects are found in recent articles, which we refer to [10] [11] [12] for both scientific and technological aspects. In the present paper, the global strategy and methodologies of the above mentioned programmes, the challenges of the experimental activities together with the status of the already performed experimental runs are briefly outlined.

The motivation

Neutrinoless double beta decay (0νββ) is an hypothetical process in which a "parent" nucleus is spontaneously transformed into a "daughter" isobar by the conversion of two neutron/protons into two protons/neutrons. The process is accompanied by the emission of two electrons (or positrons), which guarantee also electric charge conservation. If detected, 0νββ would provide evidence for the Majorana nature of the neutrino (i.e. the equivalence of neutrino and antineutrino) and give access to the neutrino mass [13] [14] [15] [16]. Under the hypothesis that the neutrino, electrons and nuclear degrees of freedom can be separated in 0νββ, the decay rate explicitly depends on the nuclear transition probability, accounted for in the so-called Nuclear Matrix Element (NME). The 0νββ NMEs are very sensitive to nuclear physics involved in decay, since they probe the nuclear wave functions of the initial and final states. In particular, the nuclear response to high-momentum second-order isospin operators is probed in 0νββ, similarly to what happen in DCE reactions [START_REF] Cappuzzello | Nuclear Response to Second-Order Isospin Probes in Connection to Double Beta Decay[END_REF]. However, there are no experimental methods to directly measure 0νββ NMEs and accurate theoretical calculations are challenging, leading to discrepancy factors larger than two, as presently reported in literature [13]. Thus experimental inputs from SCE and DCE studies are useful to help evaluate the 0νββ NMEs and to constrain the calculations. In this context, an experimental campaign has started at the INFN-Laboratori Nazionali del Sud in Catania using the MAGNEX large acceptance magnetic spectrometer [START_REF] Cappuzzello | The MAGNEX spectrometer: Results and perspectives[END_REF]. It focuses on the study of heavy-ion induced DCE reactions and competing reaction channels involving the nuclei of interest for 0νββ decay. The possible connections of DCE measured crosssections with double beta decay NMEs and, consequently, with the still unsolved neutrino puzzles are among the goals of the present studies.

Heavy-ion induced double charge exchange reactions

Challenging experimental issues must be addressed to measure heavy-ion induced DCE cross sections. Indeed one needs to detect heavy ions with good isotopic separation and energy resolution in a wide angular range, including zero-degree, in order to distinguish transitions to individual states and explore a wide momentum transfer range. In addition, the rather tiny DCE cross sections (down to few nb) demand for high experimental sensitivity, which strongly depends on the rejection capability against unwanted spurious events. These latter may be generated by competing reaction processes, some of which are more likely to occur than DCE. For this reason highly selective particle identification is a necessary prerequisite for these experiments. In addition, spurious data associated to a wrong determination of detection parameters could casually be misinterpreted as valid DCE data. Characterizing the relevant background sources can be subtle and require stringent conditions for the experimental set-up as well as a careful study of the detection response [START_REF] Calabrese | Analysis of the background on cross section measurements with the MAGNEX spectrometer: The (20Ne,20O) Double Charge Exchange case[END_REF]. Recently, accurate cross sections for heavy-ion DCE reactions have been measured at the INFN-LNS laboratory, thanks to the use of the MAGNEX spectrometer. The accurate description of heavy-ion induced DCE is a demanding task also for nuclear theory. From the reaction theory side initial and final state interactions of projectile and target systems need to be under control and the reaction amplitude must be described properly. In this context direct reaction theory can provide the suitable toolbox. Nuclear structure theory is asked to provide accurate input for modelling both projectile and target [20]. Recently, new theoretical developments for the reaction mechanism involved in heavy ion induced reactions have been achieved alongside. Presently this activity is in progress within the NUMEN project [21] [22] [23] [START_REF] Lenske | Nuclear Matrix Elements for Heavy Ion Sequential Double Charge Exchange Reactions[END_REF].

The experiments

The NUMEN experiments are performed at INFN-LNS using the K800 Superconducting Cyclotron to accelerate beams and the MAGNEX large acceptance magnetic spectrometer for the detection of the ejectiles. MAGNEX is a large acceptance magnetic device consisting of a large aperture vertically focusing quadrupole and a horizontally bending dipole magnet. A hybrid Focal Plane detector (FPD) featuring a gaseous tracker and a wall of silicon pad sensors is mounted at the focal plane to detect the reaction ejectiles [START_REF] Torresi | An upgraded focal plane detector for the MAGNEX spectrometer[END_REF] [START_REF] Cavallaro | The low-pressure focal plane detector of the MAGNEX[END_REF]. MAGNEX was designed to investigate processes characterized by very low yields, guaranteeing the identification of heavy ions with quite high mass (ΔA/A∼1/160), angle (Δθ∼0.2°) and energy resolutions (ΔE/E∼1/1000), within a large solid angle (Ω∼50 msr) and momentum range (-14% < Δp/p < +10%). Highresolution measurements for quasi-elastic processes, characterized by differential cross-sections as low as tens of nb/sr, were already performed with this setup [27] [28]. A crucial feature is the implementation of a technique of trajectory reconstruction, based on differential algebraic techniques and the accurate mapping of the spectrometer magnetic fields [29] [30] [31] [32], which solves the equation of motion of each detected particle to 10 th order [33]. This is a unique characteristic of MAGNEX, which guarantees the above mentioned performances and its relevance in the worldwide scenario of heavy-ion physics.

The experimental activity with accelerated beams proposed and presently in progress consists of two main classes of experiments, aiming at the exploration of the two directions of isospin lowering τ -τ -and rising τ + τ + , characteristic of β -β -and β + β + decays, respectively. In particular, the β + β + direction in the target excitation modes is investigated using a 18 O 8+ beam and measuring the ( 18 O, 18 Ne) DCE transitions, together with other reaction channels involving same beam and target. Similarly, the β -β -direction is explored via the ( 20 Ne, 20 O) reaction, using a 20 Ne 10+ beam and detecting the reaction products from the DCE channel and from the other open reaction channels. In each case, the experiments are organized in order to explore together with the DCE reactions of interest, all the relevant quasi-elastic processes occurring in the projectile-target collision. In this way, we can get supplementary information from the detailed study of the reaction network, profiting of the strong selectivity of quasi-elastic processes. The nucleus-nucleus interaction can be extracted from elastic scattering data and applied to the evaluation of the distorted waves in each direct reaction channel. Inelastic scattering data give selective to the role of the deformation in the coupled channel reaction network equations. Single nucleon transfer reactions shed light on the role of single particle configurations in the involved nuclear many body states. Two nucleon transfer reactions probe pairing correlations, while SCE give access to 1particle-1hole excitations. All of these information are provided in a consistent experimental framework and analysed in a multi-channel approach, based on quantum scattering theory of coupled channels. A unique opportunity is thus searched by NUMEN to build an over-constrained set of nuclear structure information, promoting the development of a new approach to 0νββ NMEs calculations. Investigations of the two classes of experiments have been started in the recent years with the present facility. These studies highlighted the strengths and the limiting aspects of the adopted technique, providing first valuable information on relevant nuclear structure aspects for a number of double beta decay emitter candidates [START_REF] Cappuzzello | The NUMEN project: NUclear Matrix Elements[END_REF]. This activity has also given precious information to guide the major upgrade of the facility toward much higher beam intensity, which is presently going on.

Experiments with 18 O beam

For the experiments of this class, the reaction channels of our interest are listed below:

• ( ) one-neutron stripping reaction One of the main challenges of such experiments is the measurement at very forward angles, including zerodegree. This is performed by placing the spectrometer with its optical axis at +3° with respect to the beam axis. Thanks to its large angular acceptance, the range -2° < θlab < 9° is explored. The MAGNEX quadrupole and dipole magnetic fields are set in order that the incident beam, after passing through the magnets, reaches a region besides the FPD. For this class of experiments, in fact, the incident ion beam ( 18 O 8+ ) has higher magnetic rigidity (Bρ) than the ejectiles of interest (namely fully stripped 18 Ne, 18 F, 20 Ne, 19 F, 16 O, 17 O ions). The beam stops in a Faraday cup, placed in the high-Bρ region besides the FPD, which measures the collected charge for each run. 116 Sn, 76 Se and 48 Ti are the targets of interest for 0νββ, which have already been explored by NUMEN, via ( 18 O, 18 Ne) reaction at 15 and 22 AMeV. The purpose was to study the 116 Sn → 116 Cd, 76 Se → 76 Ge and 48 Ti → 48 Ca transitions and the competing channels as shown in Fig. 2. A dedicated run was also performed on 12 C target, in order to compare reaction models with an easier light system. The reduction and analysis of the collected data from these campaigns is presently in progress. First results have been already published in ref. [34] [35] [36] [37].

Experiments with 20 Ne beam (β -β -direction)

In the class of experiments triggered by 20 Ne 10+ beams, the reaction channels we are interested are the following:

• ( 20 Ne, 20 21 Ne) one-neutron pickup reaction. For these experiments, the incident beam ( 20 Ne 10+ ) has a lower magnetic rigidity compared to to the ejectiles of interest. Thus, for a fixed magnetic field setting, the beam is more bent than the ejectiles of interest. The spectrometer optical axis is typically placed at -3°, thus the covered angular range is -8° < θlab < +3°. The quadrupole and dipole magnetic fields are set in order that the 20 Ne 10+ beam reaches the low-Bρ region besides the FPD. The 20 Ne 9+ and 20 Ne 8+ components of the beam emerging for the target could be a source of background in these experiments. They are significantly reduced by the use of a supplementary stripper foil mounted after the target itself [38]. The systems already experimentally explored using the ( 20 Ne, 20 O) reaction at 15 AMeV are the 116 Cd target (to study the 116 Cd → 116 Sn transition), the 130 Te (for the 130 Te → 130 Xe) and the 76 Ge (for the 76 Ge → 76 Se). The data reduction and analysis is in progress and the first results from some reaction channels have been recently published [39]

Fig. 2

 2 Fig. 2 Scheme of the reaction channels already studied within the NURE and NUMEN experiments.

ConclusionsA

  systematic study of heavy-ion induced DCE reactions and of the other reaction channels characterized by the same projectile and target is in progress at INFN-LNS. The main goal is to investigate the nuclear response to DCE reactions for all the isotopes candidate for 0νββ decay towards the determination of 0νββ NMEs.

Acknowledgments

This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No 714625).