Fiber length profiling: A novel approach to structural brain organization

&NA; There has been a recent increased interest in the structural connectivity of the cortex. However, an important feature of connectivity remains relatively unexplored; tract length. In this article, we develop an approach to characterize fiber length distributions across the human cerebral cortex. We used data from 76 participants of the Adult WU‐Minn Human Connectome Project using probabilistic tractography. We found that connections of different lengths are not evenly distributed across the cortex. They form patterns where certain areas have a high density of fibers of a specific length while other areas have very low density. To assess the relevance of these new maps in relation to established characteristics, we compared them to structural indices such as cortical myelin content and cortical thickness. Additionally, we assessed their relation to resting state network organization. We noted that areas with very short fibers have relatively more myelin and lower cortical thickness while the pattern is inverted for longer fibers. Furthermore, the cortical fiber length distributions produce specific correlation patterns with functional resting state networks. Specifically, we find evidence that as resting state networks increase in complexity, their length profiles change. The functionally more complex networks correlate with maps of varying lengths while primary networks have more restricted correlations. We posit that these maps are a novel way of differentiating between ‘local modules’ that have restricted connections to ‘neighboring’ areas and ‘functional integrators’ that have more far reaching connectivity.


Introduction
In this paper, we introduce a less studied, but important, feature of the brain; maps of cortical tract length (CTL), which we define as the length of tractography streamlines initiated from a specific point on the graywhite matter interface.
The first motivation for this work arose from current approaches to diffusion weighted MRI data analysis. Diffusion magnetic resonance imaging (dMRI) is an imaging modality that has most prominently been used to extract biologically significant information about brain microstructure by measuring the diffusion of water in the brain (Mori and Tournier, 2013;Assaf and Basser, 2005;Zhang et al., 2012;Caspers and Axer, 2017).
Examples of these metrics range from the simple fractional anisotropy (FA) and mean diffusivity (MD) (Basser and Pierpaoli, 1996) to the more complex biophysical models such as NODDI (Zhang et al., 2012). Beyond microstructure, dMRI data can be used to model fiber directions in every voxel and then use these orientation distribution functions to model possible pathways that axons might follow to connect to different brain regions (white matter tracts); this method is called tractography (Jeurissen et al., 2017). To date, the features of interest derived from tractography have been either the pathways the streamlines traverse (virtual dissections) (Catani et al., 2002;Catani and Thiebaut de Schotten, 2008), target regions of the streamlines (Sporns, 2009) or the similarity between streamlines (parcellation or ordering studies) (Johansen-Berg et al., 2004;Bajada et al., 2017b;Eickhoff et al., 2015;Cerliani et al., 2012). In a similar way to how FA and MD are used to infer information about tissue microstructure from the diffusion signal (Jones, 2009), this study uses tractography streamline lengths as a measure of a biologically useful measurement to infer which areas of cortex are involved in local or distant connections.
To this end, we used probabilistic tractography to compute streamline length distributions across each vertex of the cerebral cortex. This procedure then generated a cortical map per fiber length. The result was a unique contrast that highlights areas that have a relatively high amount of streamlines of a particular length versus areas with low amounts of streamlines propagating from them.
The second motivation for this study was to explore the neuroscientific implications of different regions being dominated by short, medium and long range fibers. An inspiration for this work comes from the field of network science. As described by (Milgram, 1967), a network can be thought of as "small-world" if every node in the network is connected to every other by a relatively short path. While sometimes disputed (Hilgetag and Goulas, 2016), during the past decade, many experiments have demonstrated small-world properties of the human brain (Bassett and Bullmore, 2016;Sporns and Honey, 2006). A hypothetical network consisting solely of short-range fibers would only be able to efficiently communicate with their anatomical neighbors. However, we know from decades of neuroanatomical work that the brain also contains long and medium fibers (described as far back as in the 19th Century) (Bajada et al., 2017a;Bajada et al., 2015;Catani and Thiebaut de Schotten, 2008;Catani and Mesulam, 2008;Menjot de Champfleur et al., 2013;Schmahmann et al., 2007;Makris et al., 2009;Dejerine and Dejerine-Klumpke, 1895).
It has been posited that the role of these longer fibers is to decrease the average path length in the brain. One can conceptualize such a network as one that contains functionally segregated modules, integrated with one another through nodes that are connected by long intermodular fibers. There is evidence that certain regions process information much like local modules , while others have been shown to be hubs that integrate information from across other brain areas to process higher cognitive functions (Margulies et al., 2016;Lambon Ralph et al., 2017).
There is another possible hypothesis for the function of these longrange fibers (Betzel and Bassett, 2018) argue that these long range fibers have a relatively minor role in the reduction of the average path length (and hence the generation of a small-world network). Indeed, if the role of these long range fibers was simply a reduction in path length, their distribution (start and termination) could be random (Betzel and Bassett, 2018) propose that the function of these fibers are to increase the diversity of inputs and outputs to specific brain regions. While it is not the aim of this article to enter theoretical discussions as to whether long range fibres contribute primarily to the "small-worldness" of the brain or to the diversification of inputs, our results clearly show that long range fibres have a clear pattern and show a predilection for specific cortical areas.
This study incorporates concepts derived from classical discussions about functional segregation and integration (Deco et al., 2015) alongside ideas from network science (Sporns, 2013), and approaches them from the perspective of brain mapping . We used our newly developed contrast to identify cortical regions with relatively more short range fibers that receive and transmit information to and from their anatomical neighbors (local modules). We also identified areas that contain longer range fibers which may allow for small-worldness, or the diversification of information, of the brain's network. These areas have the potential to communicate with disparate brain areas (integrators). In these ways, understanding the patterns of tract length across the cortex provides a novel characteristic of brain organization.
Our study is the first to demonstrate the full range (distribution) of length maps. A recent study assessed vertex-wise mean streamline lengths and a single measure of short-range and long-range fibers (Padula et al., 2017), it can be clearly seen ( Fig. 1) that streamline length per voxel has a complex distribution which is not fully captured by a single summary statistic of streamline length per voxel. Furthermore, average white-matter streamline lengths have been used previously as structural markers in the Average Path Length Maps (APMs) (Pannek et al., 2011), but once again, these reduce the information to a single statistic, the mean.

Methods
This project used diffusion weighted MRI data acquired by the WU-Minn Human Connectome Project (S900 release) (Van Essen et al. 2012, 2013Sotiropoulos et al., 2013), that were pre-processed using the minimal processing pipeline (Glasser et al., 2013;Jenkinson et al. 2002Jenkinson et al. , 2012Fischl, 2012). Further processing, detailed below, was then undertaken in order to produce the cortical tract length maps.

Participant selection
The HCP study population includes 1200 healthy adult subjects. As the cohort includes a significant proportion of siblings which included (monozygotic and dizygotic) twins, the final sample had to be reduced by excluding twins and non-twin siblings as the aim of this study was to provide normative data for CTL. Furthermore, in order to ensure as homogeneous a population as possible, participants that fit the following criteria were also excluded from the study: endocrine abnormalities including hyper and hypothyroidism, a handedness score below zero (people with a left handed tendency), color vision abnormalities, illegal drug use, history of psychiatric problems, hypertensive individuals, alcohol detected by a breathalyser, data with quality control (QC) issues A or B. All data was obtained from the S900 release and were processed using version 3.1 of the HCP pipelines. The final cohort included the first 76 individuals (age range 22-35; 52 female) of the remaining list of subjects after exclusions. In addition to the original HCP study and data collection being approved by the Washington University IRB, the current study with retrospective usage of the acquired data, including the restricted data, was approved by the local ethics committee of the Fig. 1. Histograms of the first 90 grayordinates (coordinates of gray matter vertices) across length showing the predominance of short streamlines. The colormap is normalized at every length bin; this reflects the fact that each CTL map represents a length bin across grayordinates and each map is presented with a normalized colormap (see Figs. 4 and 5). Furthermore, the maximum and minimum colormap values are truncated to the 2nd and 98th percentile. We use the same colormap throughout all the images in this paper unless stated otherwise. See BALSA scene https://balsa.wustl.edu/5KGM. University of Duesseldorf (Germany).

Image acquisition and preprocessing
The data, T1w 3D MPRAGE (0.7 mm), T2w, and diffusion weighted images, were acquired on a customized 3T scanner by the HCP group using advanced MR acquisition techniques as described in the HCP overview article  with further details found in the following articles Feinberg et al., 2010;Setsompop et al., 2012;Xu et al., 2012;Sotiropoulos et al., 2013). Images were defaced and anonymised using the algorithm in (Milchenko and Marcus, 2013). DICOM to nifti conversion was carried out using the MRICron utility dcm2nii (http://www.nitrc.org/projects/mricron).

Local diffusion modeling and tractography
Further processing was carried out using the MRTrix3 (RC 3) software package . The majority of the processing was carried out on the JURECA general-purpose supercomputer at the Juelich Supercomputing Center (Jülich Supercomputing Centre, 2016).
The local model of the diffusion signal was computed with the multishell-multi-tissue constrained spherical deconvolution algorithm (spherical harmonics order 8) (Jeurissen et al., 2014;Tournier et al., 2007). Gray matter, white matter and cerebrospinal fluid tissue probability maps were created from the T1 weighted images using the FSL FAST (Zhang et al., 2001) and deep gray matter structures were estimated using FSL FIRST (Patenaude et al., 2011) in order to perform anatomically constrained (ACT) tractography . This approach filters out any streamlines that terminate in regions that are anatomically implausible such as the middle of the white matter. Note that another approach to performing anatomically constrained tracking has also been described in the literature (Girard et al., 2014). A primary point of difference between these two methods is in the subcortical gray matter terminations of the anatomically constrained tracts. While the  approach, used in this study, provides a binary subcortical mask that allows streamlines to enter but not exit, the (Girard et al., 2014) approach does not require a binary subcortical mask.
Cortical, vertex-wise, unidirectional, probabilistic tractography using MRtrix3 was carried out (ACT; IFOD2; max angle 45 ; max distance 250 mm; downsample factor 0.25; step size 0.625; backtrack) Tournier et al., 2010). A maximum of 100K streamline propagation attempts were made from spherical regions of interest (ROIs) with a 1 mm radius centered at each vertex location of the subject specific white matter, 32k MSMAll HCP surfaces in the subjects real physical volume space . Streamline propagation was stopped either if the maximum amount of streamlines (100K) were propagated or if 5000 streamlines met the inclusion criteria determined by the ACT algorithm . This procedure generated one track file per vertex.
The streamlines were then normalized to a group template using tcknormalise. Since the MNI space has been found to overestimate the mean subject volume by differing amounts in the x, y, and z directions (Abdollahi et al., 2014), the subjects T1w images were normalized to a template that preserves the average volume of the subjects. The code used for dedrifting the MNI warps (that are provided by the HCP group) can be found on github at: https://github.com/claudebajada/mni_ dedrifting/.
The normalization compensates for head size and shape. While structural measures such as myelin and cortical thickness are generally done in the subjects' own, anatomical, space, since we are averaging histogram bins of a global metric (tract length), we wanted to ensure that, as far as possible, tracts of a particular length correspond in anatomical meaning, hence, the normalization approach was used. The dedrifting procedure ensures that the lengths obtained are correct, on average, and are not systematically affected by the MNI drift.

Length mapping
The lengths of every streamline in the track files were explicitly calculated by measuring the Euclidean distance between each streamline vertex using tckstats from the MRtrix3 package along the entire streamline in the dedrifted MNI space.
Histograms of streamline length were then created per cortical vertex and stacked along the grayordinate axis. Each row (along the length axis) in the stacked histogram represents a histogram of streamline lengths for a single vertex (or grayordinate; a "coordinate" or vertex in the gray matter as described by the HCP), while each column (along the grayordinate axis) represents a map of streamlines of a particular length.
Mean and median group CTL maps (excluding data points where tractography failed) were created across individuals and saved into ".spec" and ".scene" files for easy viewing in the connectome workbench viewer. A reliability map was also generated to determine the proportion of successful streamlines per vertex across all participants. This indicates the level of data that is available to support the value at any particular vertex (Fig. 2).

Null maps and the estimation of reliable tract lengths
Since tractography can be influenced by both microstructure as well as macroanatomy (Reveley et al., 2015), we assessed at which lengths the CTL map patterns change from ones that are only influenced by macro-anatomy and hence would help us determine which tract lengths are likely to be affected by macroanatomy. In order to explore this, we performed tractography where the local voxel model was replaced by a null distribution (NullDist2 -MRtrix3) (Morris et al., 2008). All other parameters were held the same as in the CSD tractography. This effectively means that the streamlines undergo a random walk constrained only by the anatomical properties of the brain and the basic tractography algorithm settings. These models were used to approximate the length at which the CSD tractography length distribution differed from the null tractography distributions. For this, each null CTL map was correlated with its corresponding CSD tractography map. The length at which a sharp decrease in correlation was taken as the point at which streamlines were not likely to be affected by cortical macrostructure (see Fig. 3).

Myelin, cortical thickness and resting state networks
To understand the role of CTL maps in the structural and functional relationships in the brain, we associated CTL with measures of structure (myelin and cortical thickness) and function (resting state networks). The data for magnetic resonance derived myelin (Glasser and Van Essen, 2011), cortical thickness (Glasser et al., 2013) and resting state networks (12 independent components -from a total of 15 with three excluded due to one being a subcortical network and two others not clearly definable) were acquired from the HCP S1200 group average data release -see relevant papers for detailed acquisition protocols (Van Essen et al. 2012, 2013Glasser et al., 2013;Jenkinson et al. 2002Jenkinson et al. , 2012Fischl, 2012).
For the myelin and cortical thickness relationships, individual subject length maps were correlated (Pearson Correlation) to their respective myelin and cortical thickness map. The results were Fisher-z transformed and averaged to produce a group profile along with a 95% confidence interval of the mean. This produced a mean correlation across individuals.
The same procedure was carried out for the resting state networks, though as resting state networks were defined at group level for reliable extraction of the networks, the individual CTL maps were correlated with these group resting state network independent component maps.

Visualization
The CTL maps are shown both on averaged mid-thickness cortical surfaces (Fig. 4) and on template flat maps (Fig. 5). While the midthickness surfaces are two dimensional surfaces embedded in three dimensions (real world coordinates), the flat maps have the advantage of being truly two dimensional (i.e. the z dimension has zero value). We thus populated the z-dimension of the flat maps with the CTL value (see Fig. 5 and BALSA scene: https://balsa.wustl.edu/Bv55).
The flat map version of the data is not intended to elucidate the anatomy, as this can be adequately done in on the midthickness maps (Fig. 4). The aim of the flat map visualization is to highlight the normalizations of colormaps that occur at each length bin. Indeed, due to the abundance of short range fibers as compared to long range fibers, in order to discern any contrast at higher lengths, each map needs to be min-max normalized (standard from the connectome workbench -colormap additionally truncated at 2nd and 98th percentile to exclude any potential extreme maxima and minima). This gives a misleading impression, for example, that the anterior temporal lobe had more long range fibers than it had short-range fibers. By populating the z-dimension of the flat length maps and using the min-max normalize colormap, we are able to demonstrate both the within-length contrast and the between-length comparison (for an even more intuitive insight, see the BALSA file: https://balsa.wustl.edu/file/show/nP0j). Since the patterns between the mean and median CTL maps are qualitatively similar, a decision was made to use the mean maps for visualization in this article. However, the median maps are available as BALSA database.
All figures found in this article along with the data used to create them can be downloaded from the BALSA database at: https://balsa. wustl.edu/study/show/1K3l

Results
Cortical Tract Length (CTL) maps shed light upon which cortical brain regions have a relatively high amount of structural connectivity across the spectrum of lengths (2 mm -250 mm).
CTL maps were generated along with null model tractography length maps. A reliability map was also generated to indicate the vertices where less than 5000 streamlines were generated (see Fig. 2). The null maps were compared to the probabilistic maps to disentangle the effects of brain structure on tract length. We saw a correlation with maps of connections <40 mm and a slight anti-correlation with maps of connections >40 mm (see Fig. 3).
This likely reflects the fact that the streamlines from the deep sulci (in the probabilistic tractography) do not penetrate the deep white matter (Reveley et al., 2015) while in the null model, the streamlines propagated from the deep sulci are the longest.
The overall pattern indicates that short connections occur ubiquitously across the cortex but have peaks in primary areas. The maps exhibited 'network-like' properties where major parts of the cortex had a relatively small contribution of streamlines of certain lengths while other regions showed high densities. The mid-range fibers showed distributions across secondary areas while the longer connections distributed primarily across tertiary association areas. Fig. 4 shows examples of the CTL maps across the short, medium and long fiber ranges.

Short range fibers
The shorter fibers are the most numerous and also most well distributed across the cortex (see Fig. 4b). The regions that demonstrate  The red lines show the correlation between the empirical length maps and the null model length maps with the thickness of the line representing the 95% confidence interval of the mean Fisher's z-score across all participants. Note that the length at which there is a sudden decrease in correlation with the null model is approximately 40 mm. See BALSA scene https://balsa.wustl.edu/gxjx. relatively higher concentrations in primary areas. The shortest (<40 mm) fibers on the mean image show peaks in the primary visual, somatomotor, auditory and parahippocampal cortices. At 40 mm, fronto-temporal regions along with some medial and lateral occipital areas become prominent (Fig. 4a).
The shortest fibers were found to be greatly affected by the gyral patterns in the cortex. Since up to 40 mm in length, the null tractography length maps correlated with the empirical tractography results (see Fig. 3). Due to this, we do not consider streamlines of below 40 mm to be reliably disentangled from macroanatomical information and we thus are   5. Exemplars of short, medium and long range fiber distributions in both hemispheres plotted on flat maps with the z dimension representing the streamline count. One notes that the short range fibers are all above the baseline (black flatmap visible under the short range image), the medium range fibers include regions (in red) that rise well above baseline, in the long range fibers, even the regions with the highest fiber count are barely above baseline (zero fibers). The same colormaps and color ranges were used as in the above figure. See BALSA dataset (https://balsa.wustl.edu/study/show/1K3l) and BALSA scene https://balsa.wustl.edu/Bv55. Note that these images are not designed to elucidate anatomy, this can be done with ease in Fig. 4. This figure provides a way to visualize the actual CTL values across different lengths. advise caution in interpreting the shortest fibers.

Medium range fibers
The medium range maps (Fig. 4c) tended to exhibit a "network-like" appearance, where large swathes of cortex have relatively low fiber density while specific areas show high fiber densities. Some known functional areas such as ones in the language network (with the inferior frontal gyrus, inferior parietal lobule and perisylvian temporal cortex), and limbic network (cingulum and orbitofrontal regions) formed part of the mid range fiber maps.

Long range fibers
Long range maps (Fig. 4d) tended to be restricted to the so-called tertiary association cortices. These fibers have very low absolute numbers (see Fig. 5d) and are relatively more dense in cortical regions where tracts may cross the centrum semiovale.

Resting state network length profiles
The resting state network profiles split into three main groups representing, primary, secondary and higher-order cortical networks. As one progresses across the rows in Fig. 6, the RSNs become increasingly more complex in their function. This is associated with a widening repertoire of lengths that the RSN recruits. While the primary networks recruit a narrow band of short range fibers, complex networks such as the frontoparietal networks recruit a large spectrum of lengths from 40 mm up to about 160 mm.

Structural length profiles
When correlating myelin and cortical thickness to the CTL maps, the profiles tended to present mirror opposite findings. While myelin tended to be correlated to the shorter length maps and anticorrelated to the longer ones, the opposite was generally true for cortical thickness.
Furthermore, the myelin and cortical thickness length profiles showed similar patterns to specific resting state network length profiles. For example, at approximately 100 mm, where myelin is lowest, the cortical thickness peaks and the somatosensory length profile has a pronounced dip, while the default mode and cortical thickness profiles have a small peak (see Fig. 7). There is also a peak in the myelin length profiles in the very short range fibers (0-40 mm), corresponding with a peak in the primary visual network.

Discussion
This study provided the first encompassing solution to map white matter tract lengths across the human cortex using diffusion MR tractography. This approach maps features of tractography streamlines to provide new contrast. Previous techniques such as tract density imaging (TDI) , Tract Orientation Density Imaging (TODI)  and Average Pathlength Maps (APMs) (Pannek et al., 2011) are all examples of such mapping techniques. Our study, however, is the first to investigate the full extent of streamline lengths across the entire cerebral cortex.
Our findings show that streamlines of different lengths are not evenly distributed across the cortical surface; while shorter ranged streamlines were more abundant than their long ranged counterparts, i.e. the fibers form patterns, across all lengths, that, despite us not measuring any direct connectivity, remind one of known structural and functional networks.
With regards to resting state networks, primary networks, such as the visual network, correlated with a relatively narrow band of length maps, while as cognitive functions become more complex, the networks tended to correlate with a wider range of maps (Fig. 6). These results are consistent with findings from the primate literature (Markov et al., 2012).
We further noted that cortical areas that emanate short range fibers (<40 mm) have relatively high myelin content and lower cortical thickness while this pattern inverts at longer lengths. We note, that Fig. 6. Resting state networks length profiles arranged, roughly moving from primary networks to secondary and tertiary (third row) networks. The line thickness indicates the 95% confidence interval of the mean. See BALSA scene href{https://balsa.wustl.edu/l5Vr. certain components of the length profiles of the resting state networks are also mirrored in peaks in the myelin and cortical thickness length profiles.
Finally, from the technical side, we showed that cortical streamlines below 40 mm (in the dedrifted MNI space) were potentially influenced by macro-anatomical features.

Structure and function
This work proposed a novel way to characterize cortical organization. The length maps circumscribed regions that have been described within the functional and clinical literature since the 19th century (Broca, 1861;Kelly et al., 2010;Noonan et al., 2013). The pattern of the CTLs across the cerebral cortex shows that fiber concentrations are often distributed across the cortex in patterns that strikingly resemble known functional networks. The fiber densities peak in specific anatomical regions in the literature such as Broca's area (~80 mm) (Fig. 4c), the motor cortex (~50 mm), the precuneus and the temporal pole (~160 mm). These localizations of fiber densities suggested that fibers of different lengths may contribute to different functional networks.
In order to test whether our length maps corresponded well to known functional networks, we explored the similarity of length maps with different functional resting state networks and found that lower and higher-order networks exhibited very different length profiles. As the network increased in functional complexity, the length profiles also become more complex. While lower-order networks showed one or two dominating peak fiber lengths, higher-order networks showed a contribution from all lengths (Fig. 6).
This could be interpreted in such a way that while specific length maps correspond well to relatively low-level functional tasks, complex integrative functional networks require regions that connect across multiple distances and hence correlate with multiple maps. These patterns may arise due to the fact that primary networks mainly act as local modules that do not share much information to anatomically distinct cortical areas. On the other hand, complex functional networks presumably need to integrate information from across the entire cortex in order to produce complex cognitive experiences (Margulies et al., 2016).
We then investigated the association between fiber length and structural measures such as myelin content and cortical thickness. Indeed, Cortical myelination has long been one of the features that have been used to parcellate the cerebral cortex, most famously by Oskar and Cecile Vogt. The classical developmental work by Flechsig provided evidence that primary cortices myelinate earlier than secondary and tertiary ones (Flechsig, 1901). Furthermore, these early myelinating regions tend to also be the most densely myelinated ones .
With regards to cortical thickness, while there is less of a tradition of using cortical thickness as a feature for parcellating the cortex, recent work by  have included it as a feature in their multimodal parcellation. Indeed, we propose that fibre length may vary as a function of cortical thickness since a thicker cortex may allow more fibres to emanate from it, thereby increasing the chance of there being long range fibres.
It follows, that changes in such features such as myelin and cortical thickness likely correspond to some changes in other features of the underlying neuron populations within areas of cortex, we were interested in whether the length of the neuronal populations emanating from specific regions related to the cortical thickness and myelination. Indeed, our data suggest that while shorter fibers are positively correlated with the myelin maps, longer fibers are positively correlated with the cortical thickness maps.
When viewed in light of the findings from the resting state networks, our results may indicate that primary networks involve a narrow band of short range fibers with higher myelination and thinner cortex while higher cognitive networks involve a wider band of longer fiber populations which tend to emanate from thicker cortex with lower myelination.
Our results generally mirror a study examining local and distant cortical connections using functional MRI (Sepulcre et al., 2010). Their results are consistent with our findings that a relative increase in short range connectivity underlie primary networks while association areas require more distant connectivity. We posit that the association areas require long range fibers in order to obtain information from diverse regions (Betzel and Bassett, 2018).
These results have also been shown to mirror a second pattern found by (Mueller et al., 2013), in which primary areas (of low fiber length) tend to be less variable across individuals in terms of functional connectivity measures. A later study also found that this variability extends to structural connectivity measures . This variability was correlated with amount of evolutionary cortical expansion (Mueller et al., 2013), tentatively suggesting that areas of predominantly shorter fibers may be more evolutionary conserved than areas of longer ranged fibers.

Caveats regarding the technical interpretation of the results
Diffusion MRI is a versatile in-vivo technique that allows for noninvasive investigation of properties of white matter. CTL mapping uncovers biologically relevant information using this technique. However, some technical aspects to the above study should be kept in mind when using the CTL maps.
First and most importantly, when interpreting the CTL maps, one is tempted to discuss the localizations in the maps as the termination of large association, projection or commissural tracts such as the uncinate fasciculus and the corpus callosum. This, however, is not necessarily the case (and is not designed to be the case). The maps actually show the relative density of streamlines with a given length. Hence, one region may have streamlines that contribute to two or more tracts. Furthermore, a single tract, especially if it fans (such as the uncinate fasciculus or inferior fronto-occipital fasciculus), may not all be localized at a single length. In fact, specific named fasciculi such as the ones described above may have their own characteristic length distributions. This point has been discussed as far back as the 19th century by the works of Joseph and Augusta Dejerine where they allude to the fact that fasciculi are constituted from fibers of differing length (Dejerine and Dejerine-Klumpke, 1895;Bajada et al., 2017a). Furthermore, since each map of streamlines of a particular length is a histogram bin across all cortical vertices, a value at any particular length must be interpreted relative to all other values. Finally, some tracts, for example the corticospinal tract, are truncated at the brainstem. This would lead to some areas that have a "brain tract length" rather than a "true tract length" (due to the truncation at the medulla oblongata).

Length effects and spurious streamlines
The effect of streamline length on tractography experiments is well known and widely reported in the literature (Jeurissen et al., 2017). In tractography, error propagates, hence the longer the streamline, the less reliable is its trajectory. Indeed, a recent article (Maier-Hein et al., 2017) showed that long streamlines can 'hitchhike' from one real fiber bundle to another creating an artificial bundle. While this limitation affects any technique that uses tractography such as our CTL maps, we are potentially more resilient to this error than other studies because we are not interested in the exact trajectory of the streamlines.
While our proposed approach is likely to be relatively robust to this limitation, since we are not interested in precise terminations of fibers, only their lengths, the effect may still be apparent in some regions of the cortex. For example, streamlines that are propagated from the parietal and temporal lobes tend to spread into the centrum semiovale. This is a region of high uncertainty which could cause streamlines to diverge from an anatomically plausible path and increase the apparent streamline length. As a consequence of this, at the extreme of length, the CTL maps, particularly where the absolute streamline count is very low (see Fig. 5c), may highlight a methodological limitation rather than true anatomy.
There have been techniques described recently (such as SIFT, SIFT2 and COMMIT) that attempt to filter out spurious streamlines based on microstructure as measured by the diffusion local models (Daducci et al., 2015;Smith et al., 2013). These techniques are gaining in popularity, however, we have decided against implementing them in the current study, primarily in order to establish the approach with a simple, tried and tested approach to tractography. Fiber filtering algorithms can be used in later studies in order to assess the impact of such techniques on length mapping.

Macrostructural effects, microstructural effects, and the gyral bias
Tractography is affected by both the underlying data, the physical boundaries that form part of the stopping criteria of tractography and possible biological features inherent in brain anatomy. Indeed, these independent effects may contribute to a biasing of the tractography data.  describe three factors that contribute to the gyral bias, the first relates to the fact that local model fiber orientations seldom point towards the gyral banks. Indeed, they show that the principal orientations in the gyral crown is perpendicular to the surface while along the bank and the sulcal fundus, they are almost parallel to the surface.
Secondly, tractography algorithms prefer fibers that have a smaller angular deviation across steps, this will bias streamlines to continue on straights up through a gyrus rather than terminating in a bank. This is most likely influenced by the superficial white matter fiber system (Ufibers) (Reveley et al., 2015).
Finally (Van Essen et al., 2014), show that a patch of white matter surface at the gyral crown directly corresponds to a much larger patch directly above it, while the inverse is true in the sulcal fundus (Donahue et al., 2016) examined the effect of cortical folding in relation to the gyral bias. Along the same lines of , they predicted that, since gyral crowns have increased gray matter per unit area while sulcal troughs have decreased gray matter as compared to a flat surface, the gyral bias of streamline densities should be predicted by these gyrification patterns. Indeed, the authors found that although regional differences were found in the unidirectional cortex to cortex tractography, these did not correlate with the calculated folding predicted bias. The tractography in this study, however, cannot be directly compared to that by (Donahue et al., 2016) primarily since, our current study generated up to a maximum of 100k streamlines but only kept the first 5k successful streamlines, while (Donahue et al., 2016) seeded 5k streamlines irrespective of how many streamlines were successful.
In the current study we attempted to quantify the length at which the macrostructure stops influencing the tractography. We do this by correlating CSD length maps to length maps generated with null data. The null maps can be thought of as tractography limited only by the cortical folding patterns rather than by the underlying data. We found that at 40 mm in the standard, dedrifted, space, the length patterns obtained by the null tractography and CSD tractography stop being correlated. Indeed, at higher lengths, there is a slight anticorrelation. We use this result to urge caution in the interpretation of lengths below 40 mm. Despite, this result, however, we still describe the shortest maps (<40 mm) and leave it up to the reader to judge the confidence that they have in the results.
While some approaches have been tested in order to mitigate the effect of the gyral bias (St-Onge et al., 2018), the work is still very recent and would requires further investigation before a standardised, accepted approach is formulated.

Future perspective
We believe that CTL maps can be used to increase both our basic knowledge about the structural and functional organization of the brain, and as a clinical research tool. One such basic question that may be answered is which areas (ones with short, medium or long fibers) produce greater deficit when lesioned.
Furthermore, features of the CTL maps can also be correlated with other characteristics such as cognitive performance, genetics and environmental factors.