Gyral and sulcal connectivity in the human cerebral cortex

Abstract

Introduction

The remarkable trajectory on the sophistication of basic neuroanatomical studies has provided a comprehensive assessment of cortico-cortical connectivity over the last 150 years (e.g. Dejerine and Dejerine-Klumpke 1895; Schmahmann and Pandya 2009) and allowed understanding connectivity at different levels of scale in nonhuman primates, most notably, macaques (e.g. Rockland and Pandya 1979; Barbas 1986; Schmahmann et al. 2007; de Schotten et al. 2012; Harriger et al. 2012; Morecraft et al. 2012; Goulas et al. 2014; Morecraft et al. 2015; Girard et al. 2020) and humans (e.g. Catani et al. 2003; Budde and Annese 2013; Vergani et al. 2014; Galinsky et al. 2016; Jitsuishi et al. 2017; Bullock et al. 2019). Coverging evidence from antero- and retrograde tracing along viral tracing experiments in macaques led to the formation of the so-called structural model of connectivity (Barbas and Rempel-Clower 1997), which is comprised of 3 components: (i) the description of feedforward and feedback projections between cortical areas (e.g. Jones and Powell 1970; Barbas 1986; Cavada and Goldman-Rakic 1989a, 1989b; Markov et al. 2014), (ii) a hierarchical processing from primary over uni- and multi-modal to limbic cortices (Mesulam 1998; Felleman and van Essen 1991; Chanes and Barrett 2016; Garcia-Cabezas et al. 2019); (iii) a communication model based on predictive coding (Rao and Ballard 1999; Friston 2010) embedded in an adaptive Bayesian context (Tucker and Luu 2021).

However, for a long-time, besides the description of tracts by Dejerine and Dejerine-Klumpke (1895), a comprehensive description of global, whole-brain tract architecture in humans has been lagging behind. The development of diffusion-based magnetic resonance imaging (dMRI), however, has added an unprecedented view on brain’s tract architecture during the last 20 years (e.g. Hagmann et al. 2007; Gong et al. 2009; Cheng et al. 2012; Wedeen et al. 2012; Sotiropoulos and Zalesky 2017). In this context, the data on tract topology obtained from previous studies in macaques have been essential as dMRI still cannot demonstrate formally connectivity. This limitation not withstanding, the rapid evolution of image acquisition and data-analytic methods (review: Yeh et al. 2021; Andersson et. al. 2003; Robinson et al. 2010; Jeurissen et al. 2014; Smith et al. 2016; Guevara et al. 2017; Buchanan et al. 2020; Khalilian et al. 2021) now readily allow studying the topology of tracts throughout life and under the condition of brain diseases (Bassett and Bullmore 2009), albeit at a much lower spatial resolution than histological assessment.

Because of the inherent limitations of imaging studies, several have tried to validate them by comparing fiber-tracking results between tracer-studies and with observations derived from diffusion-based methods acquired in macaques and humans. Schmahmann et al. (2007) identified the long association tracts using diffusion spectrum imaging in macaques. Quantitative properties of tracer-based network connections were analyzed by Harringer et al. (2012). Quantitative comparisons between tracer- and diffusion-based tractography in macaques (Donahue et al. 2016; Girard et al. 2020) demonstrated a considerable agreement between methods (i.e. area under the curve of the receiver operating characteristic up to 0.78). Cross-species imaging studies of humans and macaques pointed out differences between fiber tracts (de Schotten et al. 2012), network properties (Goulas et al. 2014), and gyral and sulcal connectivity (Zhang et al. 2020).

Due to the limited spatial resolution, most imaging-based studies have focused on assessing the large, relatively invariant long-range fiber bundles (>60 mm length) containing several tracts (Hagmann et al. 2007; Gong et al. 2009; Zhu et al. 2012; Bullock et al. 2019; Zhang et al. 2020; Yeh et al. 2021), which, however, comprise only 10% of all fibers (Braitenberg 1974; Braitenberg and Schüz 1998). Association fibers (u-fibers, 20–60 mm length) that connect neighboring gyri and 30% of all fibers (Vergani et al. 2014; Galinsky et al. 2016; Bullock et al. 2019) are now amenable for in vivo detection (review: Guevara et al. 2020; Zhang et al. 2014; Guevara et al. 2017). The majority of cortical connections are short (<20 mm length), reside mostly within cortical layers I, IV, and V (Vogt and Vogt 1919), and can only be studied with high-field MRI systems at this time (Leuze et al. 2014; Rowley et al. 2015).

To analyze the organization of the complex cortical networks, most typically, authors used graphs as an abstraction (e.g. Iturria-Medina et al. 2008; He and Evans 2010). In this approach, cortical regions correspond to nodes and edges to the tracts between nodes that carry a measure of connectivity between regions, which can be binary (connected or not) or weighted (e.g. reflect connection strength). An extensive number of metrics derived from graph theory were suggested to describe and analyze various properties of the network, either at the level of the whole-brain or within single nodes (review: Rubinov and Sporns 2010). There is converging evidence that connectivity is not uniform across regions in terms of connection strength and range (van den Heuvel and Sporns 2011b; Zhang et al. 2020). Local hubs have rather few, local connections, while global hubs have more, including distant connections, for which they have been called the “rich club” of the network (van den Heuvel and Sporns 2011a). Importantly, using a database of tract-tracing studies in macaques, Harriger et al. (2012) established a similar connectivity pattern as found in humans, which was confirmed by Goulas et al. (2014) in a comparative study.

The analytic approach sketched above albeit very informative lacks anatomical precision in terms of the location and structural properties of the cortical regions. Typically, cortical regions are defined by adopting any of the common atlases (e.g. (Tzourio-Mazoyer et al. 2002; Desikan et al. 2006) based on anatomical and/or functional constraints. After the quantitative analysis of the graph network properties, results are interpreted in terms of their nodal region labels as a whole. Interestingly, there are only 2 studies that aimed at differentiating connectivity between gyral and sulcal areas, which were conducted by a similar group of authors. Nie et al. (2012) pointed out that gyral regions are more strongly connected than sulci. Zhang et al. (2020) found that regions where 2 or more gyri meet (so called “3-hinge” regions) were more strongly connected than nearby gyral regions and pointed out their potential role in cortical patterning from simulations (Zhang et al. 2018).

However, the well-known considerable individual variability of cortical patterning is met with further variability of cortical connectivity (e.g. Cheng et al. 2012; Zhu et al. 2012; Chavoshnejad et al. 2021), rendering a limited understanding of the governing principles behind cortical connectivity. In our previous work, we presented a data-driven approach to determine anatomical cortical patterns with high resolution based on developmental cortical pattern formation (Kruggel 2018; Kruggel and Solodkin 2019; Kruggel and Solodkin 2020). Indeed, neuroanatomical studies of fetal development established locations and time points at which sulci appear on the cortical surface (Cunningham 1892; Chi et al. 1977; Nishikuni and Ribas 2012). Recent imaging-based studies provided converging evidence that these sulcal roots (Regis et al. 1995) are relatively invariant landmarks across subjects (Im et al. 2010) and are genetically determined (le Guen et al. 2018; Kruggel and Solodkin 2020). The earlier a sulcal root appears during development, the deeper and larger the sulcal basin around it (Kruggel 2018; Kruggel and Solodkin 2019). We found that basins cluster into 7 communities, which are further divided into 13 in some subjects (Kruggel 2018). The individual variability of the cortical patterning is governed by typically 2 sulcal roots per community that are found at relatively invariant locations across subjects (Kruggel and Solodkin 2019), called centers of low variability (CLV). We note that structurally defined communities (Kruggel 2018) at least partially correspond with results from clustering cortical modules based on connectivity (Hagmann et al. 2008). One important advantage of our topological classification is that it allows for a symbolic matching of regions between all subjects without the need for registration.

In this study, we aimed at elucidating how cortical connectivity complements the current understanding of cortical pattern formation. Instead of employing a common atlas, we used local macro-structural properties (surface curvature and geodesic depth) to define individual sulcal and gyral regions-of-interest. Due to their defining properties as landmarks, we further distinguished central sulcal regions at invariant (CLV) locations from peripheral ones. During this study, we found a considerable similarity between the structurally defined communities and the connectivity defined modules, although they are not identical. Because community borders coincide with gyral crowns, we decided to distinguish border from inner gyral regions. Thus, we analyzed properties of cortico-cortical connections in these 4 region types based on individual segmentations. We were surprised to find clearly distinguishing features between regions that were highly similar across subjects, which allowed addressing distinctly defined roles to region types in the cortical network. In addition, features helped in identifying differences in network properties between hemispheres and sexes.

Material and methods

Subjects and imaging data

From the 1113 subjects in the “1200 Subjects Release” of the Human Connectome Project (Human Connectome Project 2017) we selected 1061 subjects, for which anatomical and diffusion-weighted data sets were available. Unprocessed T₁-, T₂-, and diffusion-weighted (DW) magnetic resonance (MR) images were used in this study. For detailed acquisition information, refer to the release document (Human Connectome Project 2017).

Segmentation of cortical features

Procedures for the generation of white matter/grey matter (WM/GM) brain surfaces and basin segmentations were identical to those described previously (Kruggel and Solodkin 2019). As a result, we obtained triangulated meshes of the WM/GM interface, for each subject and both hemispheres, that contained vertex-wise basin labels. Individual hemispheric meshes can be mapped to a unit sphere (Kruggel 2008) and linearly registered across the cohort (Kruggel 2018). Basin labels were re-sampled on a common icosahedral mesh (ico7, 163842 vertices, corresponding to a spatial resolution of about 0.56 mm|$^2$|⁠) that matched to resolution of the original data. Thus, from each vertex in this space, the best matching surface location in a specific subject and hemisphere can be referenced. We collected an image of (i) the basin label, (ii) the geodesic depth, (iii) the local curvature expressed as the shape index, and (iv) the cortical myelin ratio across the |$n_v = 163842$| vertices and both hemispheres of the 1061 subjects in this cohort.

Processing of DW imaging data

The 3 B-shell series of paired L/R blip DW MR scans were corrected for subject motion and susceptibility distortions (Andersson et al. 2003). Next, voxel-wise estimates of the orientation distribution function of water motion were computed using the constrained spherical deconvolution (CSD) method (Jeurissen et al. 2014). Finally, tracking of fiber connections was performed via a probabilistic method (Smith et al. 2016).

Seeds for the tracking process were generated by the following procedure: The voxel-wise fractional anisotropy was estimated from the CSD data above, masked by the averaged |$B_0$|-images of the preprocessed DWI data. The high-resolution T₂-weighted anatomical scan was registered with the |$B_0$|-images in DWI space. Both scans were combined and segmented using a Gaussian mixture model into 4 classes (roughly corresponding to compartments cerebrospinal fluid (CSF), WM, GM, and other). The WM compartment was extracted as a binary object in voxel space. The hemispheric surface mesh with vertex-wise basin labels was adapted to this WM object using linear registration followed by nonlinear deformation. Finally, basin labels were transferred from the surface mesh to the closest voxel on the surface of the WM object. A careful assessment of the matching quality is important, in order not to bias the generation of seed voxels for tracking in sulcal vs. gyral areas. Thus, we obtained a set of basin-labeled voxels on the WM/GM interface in DWI space per hemisphere. We used on average 85 000 seeds per hemisphere and 100 repetitions per seed for tracking. Tracks that did not end at a labeled surface voxel were discarded. This process yielded 2–5 million tracks per hemisphere. We filtered out some spurious fibers based on a length/distance criterion, i.e. long fibers that connected 2 cortical areas nearby.

As a result, we obtained a set of fiber tracks, for each subject and hemisphere, connecting WM/GM surface locations with known basin labels.

Clustering of fiber connections into communities

Using fiber tracks and basin labeling, we computed hemisphere-wise fiber connectivity matrices. Consider subject |$a$| with |$n_a$| distinct basin labels and allocate a matrix |$C(n_a,n_a)$|⁠. Increment element |$C(i,j)$| if a fiber starts at basin label |$i$| and ends at basin label |$j$|⁠. Normalize |$C$| by the number of connections. Thus, entries in this matrix represent the connection probability between regions.

Fig. 1

Left: Vertex-wise variability of fiber connectivity. Regions in magenta show low inter-subject variability, separated by rims of higher variability (green–red). Middle: Segmentation into 6 fiber communities. Right: Mapping of fiber communities onto an “inflated” left hemisphere.

Open in new tabDownload slide

Segmentation of gyral and sulcal regions

During this study, we found distinctive differences in the cortical connectivity between gyral and sulcal regions. So we decided to introduce a finer grained segmentation of an hemispheric surface than the basin segmentation described above.

An individual hemispheric surface was segmented into sulcal and gyral areas using the scale-invariant shape index, which is positive in concave and negative in convex regions. Peaks of gyral crowns were determined as local maxima of the shape index. Sulcal roots were determined as local minima of the shape index with a minimum geodesic depth of 5 mm (Fig. 2, top left). Results were filtered such that peaks had a distance of at least 8 arc degrees. Starting from these peaks and roots, a curvature-guided watershed region growing procedure was used to segment the hemispheric surface into the gyral and sulcal patches.

Fig. 2

Region segmentation in an example subject. Top left: Borders of fiber communities (black) overlaid onto geodesic depth (magenta: gyral crowns; green: sulci). Top right: Borders of fiber communities (black) overlaid onto CLV (colored). Below left: Region segmentation into (i) border gyral regions (shades of red); ii) inner gyral regions (shades of yellow); (iii) core sulcal regions (shades of blue); (iv) peripheral sulcal regions (shades of green). Below right: Segmented regions mapped onto the corresponding WM/GM surface.

Open in new tabDownload slide

By construction, borders between fiber communities coincide with gyral ridges. Borders were extracted and dilated on the spherical mesh to bands of 8 vertices width. Gyral patches that coincide with borders were denoted as border gyral regions, the others as inner gyral regions. In our previous work (Kruggel and Solodkin 2019), we found that the deepest sulcal roots show much less variability in space across subjects, called CLV (Fig. 2, top right). Thus, we distinguished sulcal patches that coincide with CLVs as core sulcal regions, the others as peripheral sulcal regions (Fig. 2, below left).

Consider the frontal community (Fig. 1, right, in red). An example mapping of segmented regions on an individual hemisphere (Fig. 2, below right) shows community borders (in red) along the superior rim of the hemisphere and the precentral gyrus, inner gyral regions (in yellow), peripheral sulcal regions in green, and core sulcal regions (in blue) in the inferior section of the precentral sulcus and the posterior part of the superior frontal sulcus. Region statistics were compiled in Table 1.

Measures of network connectivity

Subject-wise connectivity matrices were recomputed according to the regions obtained above. Connectivity matrices are equivalent to a weighted bidirectional network graph, in which nodes corresponds to cortical regions, and edge weights to the connection probability between regions. Regions were characterized by 13 nodal network metrics. The first 10 were weighted nodal metrics as compiled by Rubinov and Sporns (2010):

• Degree corresponds to the number of regions connected to region |$i$|⁠:

  $$\begin{align}& d_i = \sum_j^n C(i,j)> 0. \end{align}$$

  (2)
• Strength measures the cumulative probability of connections to region |$i$|⁠:

  $$\begin{align}& s_i = \sum_j^n C(i,j). \end{align}$$

  (3)
• Betweenness centrality corresponds to the fraction of all shortest paths that connect through region |$i$|⁠:

  $$\begin{align}& b_i = \frac{1}{(n-1)(n-2)} \sum_j^n \frac{\# p_i(j,k)}{\# p(j,k)}, \end{align}$$

  (4)

  where |$\# p(j,k)$| is the number of shortest paths between regions |$j,k$| and |$\# p_i(j,k)$| those that pass through region i. Regions with higher betweenness centrality have more control of the network.
• The efficiency of a region’s connectivity is computed as:

  $$\begin{align}& e_i = \frac{1}{d_i(d_i-1)} \sum_j^n \left ( \frac{C(i,j) \, C(i,k)}{l(p_i(j,k))} \right )^{1/3}, \end{align}$$

  (5)

  where |$d_i$| corresponds to the degree of region |$i$| and |$l(p_i(j,k))$| to the length of the shortest path between regions |$j,k$| that contains only neighbors of region |$i$|⁠. Efficiency represents the resilience at region |$i$|⁠: how well connectivity is provided by its neighbors when it is removed.
• Clustering coefficient represents the embeddedness of a region in the network:

  $$\begin{align}& c_i = \frac{2}{d_i(d_i-1)\, C_{max}} \sum_{j,k} \left ( C(i,j) \, C(j,k) \, C(k,i) \right )^{1/3}, \end{align}$$

  (6)

  where |$d_i$| corresponds to the degree of region |$i$|⁠, |$C_{max}$| to the maximum connection probability in the network, and |$C(i,j), C(j,k), C(k,i)$| to the connectivity of regions with triangle connections in the one ring at region |$i$|⁠.
• Small-worldness measures the locality of a region’s connectivity:

  $$\begin{align}& w_i = \frac{c_i \, l_r}{l_i \, c_r} \quad\mbox{with}\quad l_i = \frac{\sum l(p_i(j,k))}{\# p_i(j,k)}, \end{align}$$

  (7)

  where |$c_i$| corresponds to the clustering coefficient of region |$i$| and |$l_i$| to the average length of shortest path connecting through region |$i$|⁠, in relation to measures of a random graph |$r$|⁠. A smallworld topology has a high clustering and a short path length.
• Strength z-score corresponds to the normalized strength of a region within its community:

  $$\begin{align}& z_i = \frac{(s_i - \overline{s_j})}{\sigma(s_j)} \quad\mbox{with}\quad s_j = \sum_{k \in M_i} C(j,k), \end{align}$$

  (8)

  where |$s_i$| corresponds to the connection strength of region |$i$| in community |$M_i$|⁠, |$\overline{s_j}$| and |$\sigma (s_j)$| to the mean and standard deviation of the within community strength.
• Participation coefficient measures the relative connection strength of region |$i$| between communities:

  $$\begin{align}& y_i = 1 - \frac{\sum_m s_m^2}{s_i^2} \quad\mbox{with}\quad s_m = \sum_{k \in M_m} C(i,k), \end{align}$$

  (9)

  where |$s_m$| corresponds to the connection strength of region |$i$| to those in community |$m$| and |$s_i$| to the strength of region |$i$|⁠.

• S-level decomposition: In the graph that represents a connectivity matrix, the smallest connectivity |$s_{min}$| between any 2 nodes was determined. All edges in the graph with a connectivity less or equal than |$s_{min}$| were removed. A level counter for all remaining nodes was increased. This process was repeated until the graph was empty. Thus, the s-level of a region represents its ranked strength in the network.

• S-core decomposition: The graph was pruned as described above, but the highest value of |$s_{min}$| was kept for each node. Thus, the s-core of a region corresponds to its weakest connectivity in the network.

We added the following metrics:

• Density represents the quotient of the connection strength and surface area of a region.

• Reach measures a region’s connectivity, weighted by the length of shortest paths:

  $$\begin{align}& r_i = \frac{\sum_j^n C(i,j) l(p(i,j))}{\sum_j^n C(i,j)}. \end{align}$$

  (10)

• Myelin ratio is the quotient of the |$T_1$|- by |$T_2$|-weighted intensity, averaged across the cortical segment of a region. It is understood that this ratio is proportional to the fraction of cortical myelin (Glasser and van Essen 2011; Norbom et al. 2020).

• The cortical thickness is the distance between the CSF/GM and GM/WM interfaces, averaged across the cortical segment of a region.

Note that the last 2 metrics measure were not derived from DW images and do not represent connectivity per se.

We acknowledge some redundancy in this list of connectivity measures. Several metrics assess similar properties, as indicated by high correlation coefficients between them. However, we found that each measure provided a valuable contribution to the overall characterization of a region’s connectivity. We decided to (i) prefer probabilistic over deterministic tracking; (ii) use continuous over binary network metrics; (iii) avoid using arbitrary thresholds on network metrics. These decisions match recommendations of a recent publication by an expert panel (Yeh et al. 2021).

Because ranges of measures differed across metrics, and most metrics did not follow a Gaussian distribution, we resorted to a nonparametric normalization. Region-wise results of a specific metric obtained in a hemisphere of a specific subject were ranked and z-transformed. Region-wise z-scores can be (i) averaged per region type, leading to 4 z-scores per metric and hemisphere; (ii) mapped to the reference surface, and averaged across all subjects in this cohort, leading to a spatial map of a metric. We considered this approach as viable because it was based on sufficiently large sample of |$\approx $|920 regions. Transformed metrics were now Gaussian distributed with similar ranges across metrics, which makes using operations such as principal component analysis viable.

Results

In this section, we described (i) general properties of connectivity for the cortical region types; (ii) provided a deeper analysis of the cortical network and its local properties; and (iii) analyzed inter-hemispheric and sex-related differences in brain connectivity.

Connectivity of regions

The connectivity between region types was expressed as the fraction of total streamlines between them (Table 2). Border gyral regions mostly connected to core sulcal and other border regions. An example of border-to-border gyral connectivity is the aslant tract that bypasses under the frontal community from the medial community to the insula. Similarly, inner gyral regions connected to core sulcal and other inner regions. Inner-to-inner gyral connections likely correspond to local within community subcortical u-fibers. Besides the reciprocal connections to gyral regions, core sulcal regions predominately connected to other core sulcal regions. Because core regions were sparse (12% of all, with 2–3 per community), this finding spurred the hypothesis that they form nonlocal (between-community) connections. Finally, peripheral sulcal regions connected almost equally to all other regions, indicating a role as local communication integrators. Results were highly similar for both hemispheres. Functional segregation, facilitated by high structural clustering coefficient (Rubinov and Sporns 2010), is widely acknowledged to be a fundamental characteristic of the cortex.

Next, we assessed the distribution of connectivity vs. fiber path length, as given by the shortest path between 2 regions in the anatomical network. Each path segment corresponds to roughly 10 mm length. In a first approximation, the distribution over all region types followed the “power law” described before (Ercsey-Ravasz et al. 2013), which was actually composed of a mixture of log-normal distributions. However, when separating for specific region types, distributions showed considerable differences. For a finer analysis, we compared the amount of connections between specific region types at different path lengths using t-tests. Results were compiled in Fig. 3 as t-values vs. path length. Values are positive if the first sample had stronger connections at a specific path length.

Fig. 3

Strengths of fiber connections vs. path length in the left (black) and right (red) hemisphere. Left: t-value of comparison between (i) gyral and (ii) sulcal fibers. Right: t-value of comparison between (i) border vs. (ii) inner gyral fibers. t-values are positive if compartment (i) and negative if compartment (ii) had stronger connections. Values above and below the dotted lines were statistically significant (p < 0.05).

Open in new tabDownload slide

From the comparison of gyral vs. sulcal regions (left panel), we found more gyral connections at a mid-range distance (path length 5–10), but more sulcal connections at shorter and longer distances. This peak of the gyral distribution likely corresponded to the local and regional u-fibers that originated and ended in gyral regions, whereas most long (between-community) fibers originated in sulcal regions. Further differentiation between border and inner gyral regions (right panel) showed that connectivity of border regions was stronger for path lengths greater than 6. The connectivity was always stronger in core than peripheral sulcal regions (not shown). In summary, inner gyral and peripheral sulcal regions had predominately short, local connectivity, border gyral regions mostly mid-range (within- and between-community) connectivity, whereas core sulcal regions dominated in nonlocal connectivity.

Network metrics in regions

Network metrics of the individual fiber connectivity were computed in the 4 region types, as described in Section 2.6. Region-wise measures were ranked, z-transformed, and averaged by region type and hemisphere (Table 3). Surprisingly, metrics were strikingly different between region types. We ranked regions by their metrics (t-tests, < indicates |$p < 1e-4$| after Bonferroni correction).

Core sulcal regions had the highest scores in 10 of 12 metrics. They showed a high fiber density with a strong local connectivity (degree, strength, efficiency, clustering, z-score, and small-worldness). However, they also had considerable remote connectivity (s-core, s-level, reach). Together, they form the “rich club” of local hubs. In contrast, peripheral sulcal regions scored lowest in 10 of 12 metrics, with higher scores in the participation coefficient, characterizing them as local information integrators.

Border gyral regions scored highest in terms of betweenness centrality and participation coefficient, but second lowest in most other network metrics, which characterized these regions as between-community connectors, in line with their anatomical location. Finally, inner gyral regions scored second highest in 8 of 12 metrics, consistently below core sulcal regions. Their low participation coefficient, but high small-worldness and efficiency characterized them as local processing units.

In addition to these 12 metrics based on DW data, we used the cortical myelin ratio and cortical thickness as complimentary metrics that were based on the available |$T_1$|- and |$T_2$|-weighted images. We found strong differences between region categories: The myelin ratio and thickness decreased from border and inner gyral to core and peripheral sulcal regions.

As discussed in Section 2.6, some metrics are highly correlated. To better understand the relevance of region properties, we used a principal component analysis to project results into a 3D space, retaining 86.3% of the variance on the left, and 82.1% on the right side. Component plots (Fig. 4) showed distinct clusters for the 4 region types. Note that the clusters spread decreased with the proportion of regions per type. The most compact cluster (peripheral sulcal regions, in blue) was based on 44% of all regions, whereas the most extensive cluster (core sulcal regions, in green) was based on only 10% of all regions. This may indicate that rather measures are “noisy” than region properties differ between individuals.

Fig. 4

Principal component analysis of the 14 metrics in the left (top) and right hemisphere (below). Colors denote region types: border gyral (black), inner gyral (red), core sulcal (green), peripheral sulcal regions (blue).

Open in new tabDownload slide

Loadings on principal components (Table 4) demonstrated a predominant influence of strength- and density-related metrics on component 1 (density, degree, strength, betweenness centrality, z-score, s-level, and s-core decomposition); of locality-related metrics on component 2 (efficiency, clustering coefficient, small-worldness, participation, fiber reach); and anatomical metrics (myelin ratio, cortical thickness) on component 3. Thus, component 1 distinguished between border gyral and core sulcal vs. inner gyral and peripheral sulcal regions in terms of connection strength. Component 2 distinguished between inner gyral and core sulcal vs. border gyral and peripheral sulcal regions in terms of connection locality. Gyral regions were distinguished from sulcal regions in terms of thickness and myelin content by component 3.

Finally, we clustered PCA-transformed metrics using a Gaussian mixture model (Scrucca et al. 2017). The best model had 4 classes on both sides (as indicated by the Bayesian information criterion). Comparing the clusters against the region type in a confusion matrix demonstrated an excellent match (normalized mutual information left: 0.957; right: 0.959). Note that the same level of discrimination can be achieved with only 5 metrics (strength, clustering coefficient, reach, myelin ratio, and cortical thickness), indicating the redundancy of the metrics. In conclusion, regional metrics can be used to predict the anatomically defined region type.

Hemispheric asymmetries in region connectivity

Strongly significant inter-hemispheric asymmetries in region connectivity were found for most metrics (Table 5).

Interestingly, most asymmetries were leftward for inner gyral and peripheral sulcal regions, and rightward for border gyral and core sulcal regions. In the component analysis (Fig. 4), core sulcal (green) and border gyral regions (black) were more distinct from inner gyral regions (red) on the right side. On the left side, core sulci had a weaker connectivity (PC 1: L-R = -0.3687), but a much higher locality (PC 2: L-R = +0.5097), whereas peripheral sulci showed little differences. Border gyri were more strongly connected (PC 1: L-R = +0.209), whereas inner gyri had a weaker locality (PC 2: L-R = -0.3036). Results of network metrics (Table 3) demonstrated that core sulci were strongly connected and had much higher scores in the s-level and s-core decomposition on the right side, indicating stronger remote than local connectivity.

Sex-related differences in region connectivity

Sex-related differences in network and anatomical metrics were compiled in Table 6.

Interestingly, differences in network metrics for border gyral and core sulcal regions were generally not significant. In all significant differences for inner gyral regions, females scored lower than males, whereas for peripheral sulcal regions, females scored higher than males. Thus, in the context of the ranking of network metrics (Table 3) differences between inner gyral and peripheral sulcal regions were smaller for females than males. The component analysis (Table 4) pointed rather towards sex-related differences in local connection strength (PC 1) than in the network architecture (PC 2).

Discussion

This study found 4 cortical region types with highly distinct anatomical and network features. Statistical properties strongly differed between region types, but were highly similar across subjects (Fig. 4). Note that region types were defined on hemisphere-wide averages, e.g. each hemisphere contributed only data point per region type to this figure. It is conceivable that these region types provide a meaningful abstraction of the individual neuroanatomy, as will be discussed in the following. First, we provide a more comprehensive understanding of the cortico-cortical connectivity in terms of their anatomical basis and developmental aspects, leading to an anatomically driven model of cortical connectivity.

Role of region types in the brain network

In summary, we interpreted the role of the 4 region types as:

• With a high betweenness and participation coefficient but low small-worldness, gyral regions at community borders were characterized as connector nodes between fiber communities, in agreement with their anatomical definition. These regions were strongly connected to core sulcal regions (and other border regions) and had the highest proportion of mid-range fibers. In particular, border gyral regions had the highest myelin ratio and cortical thickness.

• Core sulcal regions scored highest in 10 of 12 network metrics, featuring a strong local connectivity to inner and border gyri of the same community. These regions had the highest fraction of connections to other core regions, the highest fraction of nonlocal fibers, and high s-core and s-level coefficients, which all indicated a strong remote connectivity to core sulci in other communities. With 2 or 3 such regions per fiber community, they assume a central role as communication hubs within a community and beyond. These regions are located at the deepest sulcal roots, which develop first during gestation, are genetically determined, and show a low inter-subject variability. Together, they form the “rich club” and backbone of the communication network.

• In contrast, sulcal regions in the periphery scored lowest in 10 of 12 network metrics, but with higher ranks for the participation coefficient. They had mostly short fibers that were equally connected to the other region types, characterizing them as local information integrators. Peripheral sulcal regions had the lowest myelin ratio and cortical thickness.

• Gyral regions inside a fiber community ranked second highest in 8 of 12 metrics, but consistently below core sulcal regions. They had a high fraction of predominately short connections to border gyral and core sulcal regions, and the lowest participation coefficient, characterizing them as regions with predominantly local connectivity.

The communication roles of the 4 region types were compiled in a scheme of a cortical sub-module (Fig. 5). Note that this scheme was simplified in its anatomical composition, because the number and spatial configuration of regions depended on the location of the submodule and differed between subjects (compared with Fig. 2, bottom row).

Fig. 5

Simplified schema of 2 (remote) communities, consisting of border gyral (B, red), core sulcal (C, blue), inner gyral (I, yellow), and peripheral sulcal regions (P, green). Line widths approximately denote the connection strength between regions.

Open in new tabDownload slide

The bi-directional arrows in the model should not be understood as reciprocal connections as documented in tracer experiments (e.g. Jones and Powell 1970; Barbas 1986; Cavada and Goldman-Rakic 1989a, 1989b; Markov et al. 2014). Diffusion-based methods cannot distinguish between source and target of fiber connections. We differentiated results in Table 2 between source and target in order to demonstrate that tracking from either direction yielded similar results. This equivalence rendered potential bias due to seed voxel selection and tracking procedure as unlikely.

Roles of sulci and gyri in anatomical connectivity

There is mounting evidence that sulci and gyri serve different roles in brain functioning (review: Jiang et al. 2021; Amiez et al. 2019; Ge et al. 2017; Liu et al. 2019). Surprisingly, the distinction of gyri and sulci in terms of their connectivity has received little attention so far, both in tracer- and diffusion-based studies, although differences in anatomical properties (e.g. cortical thickness, myelin content) are well documented. We only found 2 (not independent) diffusion-based studies (Nie et al. 2012; Zhang et al. 2020) that explicitly differentiated between sulci and gyri in terms of connectivity. Both reported remarkably lower connection strength and density in sulci compared with that in gyri. In a first approximation, we confirmed this qualitative difference, but not to the extent as reported by Nie et al. (2012), likely due to differences in their analytic approach. Convergent evidence in the neuroanatomical literature documented a significant but not substantial lower fiber density in sulci vs. gyri (e.g. Smart and McSherry 1986; Braitenberg and Schüz 1998; Ercsey-Ravasz et al. 2013; Galinsky et al. 2016; Jitsuishi et al. 2017; Mortazavi et al. 2017).

Previous studies did not distinguish between sulcal region types. We found that core sulcal areas take a central role in the brain communication network. Because core regions comprised only 20% of sulcal areas, their distinguishing role was likely missed. Their special role is neurobiologically plausible. Sulcal roots appear early in development (Nishikuni and Ribas 2012; Kruggel and Solodkin 2019), are genetically determined (le Guen et al. 2018; Kruggel and Solodkin 2020), and are the most invariant across subjects (Im et al. 2010; Kruggel and Solodkin 2019). In order to achieve regularity and consistency during development (White et al. 2010), it is useful to draw the formation of long-range connections—that are “expensive” to form (Bullock and Sporns 2012; Ercsey-Ravasz et al. 2013)—to early defined landmark locations. This pre-allocation of long fibers does not impede the formation of other (e.g. inter-gyral) connections, due to the layered organization of fiber bundles. The shortest fibers are the most superficial, whereas longer fibers travel at deeper levels of the white matter. Indeed, the early formation of long inter-sulcal connections may facilitate and guide the routing of other fibers between (gyral) structures that develop later.

We also distinguished between gyral regions on borders and inside communities. The finding that border gyral regions have markedly different network properties was, at a second glance, less surprising. Due to their distinguished location, these regions, we suggest, serve a role as “staging areas” where information is relayed between neighboring communities. Often, border regions directly connect to core sulcal areas, e.g. portions of the uncinate fasciculus connect a core sulcal area in the fronto-orbital cortices with border regions between the 2 temporal communities close to the temporal pole. We also found the highest myelin ratio and thickest cortices in border areas, which we suggest is as an independent corroboration of their proposed distinctive role.

As an alternative to our border regions, we assessed “3-hinge” regions (Deng et al. 2014; Zhang et al. 2018; Zhang et al. 2020). Their locations were readily available to our workflow, because 3 basins meet at such locations. We could qualitatively confirm their conclusion that these regions were more strongly connected than other gyral regions nearby. However, the distinction between our border and inner gyral regions was statistically much stronger than theirs. We point out that “3-hinge” locations are often located at the “end units” of large sulci that coincide with community borders. Thus, we assumed that the statistically weaker differences between “3-hinge” and other gyral regions was a consequence of diluting statistical power of the distinct properties of our border and inner gyral regions.

Inter-hemispheric and sex-related differences

It is much expected that inter-hemispheric and sex-related differences of region-wise properties were statistically weaker than those between region types. In the context of the point clouds in Fig. 4, differences between points within the same cloud were assessed here. In particular, we found that core sulcal and border gyral regions were more distinct from inner gyral regions on the right side, which may lead to the conclusion that regions on the right hemisphere show a higher degree of anatomical specialization. Because we compared averaged properties per hemisphere, it must be emphasized that a localized assessment may reveal stronger differences in specific brain areas. In this context, Iturria-Medina et al. (2011) found that the right hemisphere had a higher efficiency and stronger connectivity, which was in line with our results.

Sex-related differences were more minute. Interestingly, network metrics for border gyral and core sulcal regions were similar. Whereas inner gyral and peripheral sulcal regions were more distinct in males, females tended to have a weaker overall connectivity, but a stronger nonlocal connectivity, as indicated by a lower small-worldness, confirming results of a previous study (Tunc et al. 2016).

An anatomically driven model of cortical connectivity

Recent studies about the allocation of cortical sectors (Puelles et al. 2019) and mechanisms that predate gyrification (Rash et al. 2019) have added significant insight into the neurobiology of early cortical development. Our results demonstrated differences in the properties of gyri and sulci and the spatial heterogeneity of cortical connectivity that presumably develop after gestation week 20 until the first years of life. At this time, early models do not explain or predict the neuroanatomically well-documented temporal sequence of gyrification events during the second trimester (Nishikuni and Ribas 2012; Kruggel and Solodkin 2019) and its influence on the formation of cortical connections.

Our results are compatible with the structural model of connectivity (Barbas and Rempel-Clower 1997; Garcia-Cabezas et al. 2019) because this model does not explicitly differentiate between gyri and sulci or subtypes thereof. It is possible that their microscopically defined “cortical type” is related to our macroscopically defined “region type”. If whole cortical maps of the cortical type become available (Garcia-Cabezas et al. 2020), a quantitative comparison between both models appears feasible. In addition, our model scheme (Fig. fig-scheme) does not conflict with the hierarchical aspect of directed information flow of the structural model, because it concerns the organization of local connectivity and does not cover directional aspects of information flow.

At the highest level of structural organization we found regions cluster into spatially compact communities (or modules) based on an abundance of predominantly local connections, in line with several previous reports (e.g. Rubinov and Sporns 2010; Bullock and Sporns 2011b; van den Heuvel et al. 2012). Connectivity-based communities largely aligned with anatomically defined communities that we described earlier (Kruggel 2018), especially along anatomical constraints (e.g. the rims of a hemisphere).

Within a community, a hierarchical organization was proposed into (i) locally connected “small-world” regions; (ii) “feeder” regions, and (iii) “rich-club” regions that access the backbone communication network of the brain (Bullock and Sporns 2012; van den Heuvel and Sporns 2012). This hierarchically modular organization was found to facilitate limited sustained network activity (Kaiser and Hilgetag 2010) and may support the critical range of human physiology (Kitzbichler et al. 2009). The modular architecture devised from our results (Fig. 5) fits nicely due to their high s-level and s-core decomposition scores and their low small-worldness. The fact that these regions are most invariant across subjects and appear early in development allows them to provide a consistent base or scaffolding to the cortex’s modular structure. In addition, their distributed and densely interconnected network may play a central role in facilitating efficient global functional integration and information flow in the cortex (Bullock and Sporns 2012; Samu et al. 2014). We addressed border gyral regions with a role as communicators with adjacent communities, although some of their connections can be long (e.g. the aslant tract or sections of the arcuate fasciculus). Peripheral sulcal and inner gyral regions take roles as feeders and local connectors.

We included about 20 core sulcal locations (i.e. 2–4 per community) that comprised the most invariant across subjects. Each core and its adjacent regions corresponded to 30–70 cm|$^2$| of cortical area. The path length between peripheral/inner regions and the closest core sulcal was never greater than 2 (i.e. <50 mm). As indicated in Fig. 5, together they form a sub-modular “small-world” of regional communication.

Note that the 4 region types with distinct structural and connection properties may serve as a meaningful abstraction of the individual cortical variability. Our initial studies of pattern types (“motifs”) around a specific region indicated that there were only a few and that their frequency matched with the scheme of Fig. 5. When framing the individual cortical surface into sub-modules, it appears that the inter-individual variability is considerably reduced. This is facilitated by construction, because each sub-module was centered around a core sulcal area at strongly consistent locations across subjects. It is conceivable that much of the inter-individual variability is confined within a sub-module and is related to individual differences in the number and spatial configuration of the 4 region types defined in this work. Thus, the problem of finding corresponding cortical locations between subjects may be reduced to finding a region-to-region mapping between homologue sub-modules, which can be solved by computationally efficient methods.

Finally, we note that our work opens up a straightforward approach for a meaningful assessment and classification of fiber tracts, in terms of the connected region types and their locality (i.e. within a sub-module, module, or across modules).

Acknowledgment

Authors would like to thank the Research Cyberinfrastructure Center at UC, Irvine for their technological support.

Conflict of interest statement: None declared.

References