Microstructural and functional gradients are increasingly dissociated in transmodal cortices

Abstract

While the role of cortical microstructure in organising neural function is well established, it remains unclear how structural constraints can give rise to more flexible elements of cognition. While nonhuman primate research has demonstrated a close structure–function correspondence, the relationship between microstructure and function remains poorly understood in humans, in part because of the reliance on post mortem analyses, which cannot be directly related to functional data. To overcome this barrier, we developed a novel approach to model the similarity of microstructural profiles sampled in the direction of cortical columns. Our approach was initially formulated based on an ultra-high–resolution 3D histological reconstruction of an entire human brain and then translated to myelin-sensitive magnetic resonance imaging (MRI) data in a large cohort of healthy adults. This novel method identified a system-level gradient of microstructural differentiation traversing from primary sensory to limbic regions that followed shifts in laminar differentiation and cytoarchitectural complexity. Importantly, while microstructural and functional gradients described a similar hierarchy, they became increasingly dissociated in transmodal default mode and fronto–parietal networks. Meta-analytic decoding of these topographic dissociations highlighted involvement in higher-level aspects of cognition, such as cognitive control and social cognition. Our findings demonstrate a relative decoupling of macroscale functional from microstructural gradients in transmodal regions, which likely contributes to the flexible role these regions play in human cognition.

About the Authors:

Casey Paquola

Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Validation, Visualization, Writing – original draft, Writing – review & editing

* E-mail: [email protected] (BCB); [email protected] (CP)

Affiliation: McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada

ORCID logo http://orcid.org/0000-0002-0190-4103

Reinder Vos De Wael

Roles Data curation, Methodology, Writing – review & editing

Affiliation: McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada

Konrad Wagstyl

Roles Methodology, Software, Writing – review & editing

Affiliation: McGill Centre for Integrative Neuroscience, McGill University, Montreal, Canada

ORCID logo http://orcid.org/0000-0003-3439-5808

Richard A. I. Bethlehem

Roles Formal analysis, Visualization, Writing – review & editing

Affiliation: Autism Research Centre, Department of Psychiatry, University of Cambridge, Cambridge, United Kingdom

ORCID logo http://orcid.org/0000-0002-0714-0685

Seok-Jun Hong

Roles Data curation, Methodology, Writing – review & editing

Affiliation: McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada

ORCID logo http://orcid.org/0000-0002-1847-578X

Jakob Seidlitz

Roles Resources, Validation, Writing – review & editing

Affiliations Developmental Neurogenomics Unit, National Institute of Mental Health, Bethesda, Maryland, United States of America, Brain Mapping Unit, University of Cambridge, Department of Psychiatry, Cambridge, United Kingdom

ORCID logo http://orcid.org/0000-0002-8164-7476

Edward T. Bullmore

Roles Resources, Validation, Writing – review & editing

Affiliation: Brain Mapping Unit, University of Cambridge, Department of Psychiatry, Cambridge, United Kingdom

ORCID logo http://orcid.org/0000-0002-8955-8283

Alan C. Evans

Roles Methodology, Resources, Writing – review & editing

Affiliations McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada, McGill Centre for Integrative Neuroscience, McGill University, Montreal, Canada

Bratislav Misic

Roles Methodology, Software, Writing – review & editing

Affiliation: McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada

ORCID logo http://orcid.org/0000-0003-0307-2862

Daniel S. Margulies

Roles Methodology, Visualization, Writing – review & editing

Affiliation: Frontlab, Institut du Cerveau et de la Moelle épinière, Paris, France

ORCID logo http://orcid.org/0000-0002-8880-9204

Jonathan Smallwood

Roles Investigation, Supervision, Writing – review & editing

Affiliation: York Neuroimaging Center, University of York, York, United Kingdom

Boris C. Bernhardt

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

* E-mail: [email protected] (BCB); [email protected] (CP)

Affiliation: McConnell Brain Imaging Centre, Montreal Neurological Institute and Hospital, McGill University, Montreal, Canada

ORCID logo http://orcid.org/0000-0001-9256-6041

Introduction

A core principle of neuroscience is that brain structure governs ongoing function. For decades, nonhuman primate work has confirmed the intrinsic relationship between microstructure and macrolevel function [1,2]. While the role of structure in cortical processing is well characterised in the sensory motor domain, it is less clear how this constraint gives rise to more flexible elements of cognition. In an attempt to describe how cortical areas are able to take on more abstract functional roles, Mesulam (1998) postulated a hierarchical axis of large-scale cortical organisation and connectivity, referred to as the ‘sensory–fugal gradient’ [3]. This axis describes gradual transitions at the whole-cortex level, running from primary sensory and motor regions involved in externally focussed computations toward transmodal cortices, in which neural responses are not segregated by modality and operations are increasingly decoupled from perception [4,5] (see S1 Table for nomenclature). Unlike sensory cortices, transmodal cortices have a less hierarchical organisation with dense interconnectivity and top-down projections that tend to jump synaptic levels and that allow spatially distributed areas to respond flexibly to different types of information [5–10].

The current study systematically examined the interplay between microstructural constraints and functional connectivity in humans and its contribution to high-level cognition and behaviour. Studying cortical microstructure in humans, neuroanatomists have traditionally used cell-staining techniques to map spatial variations in both cyto- and myeloarchitectural features of post mortem brains [11–15]. Extending upon this work, the ‘structural model’ proposes a tight coupling of microstructural similarity with increased probability of interareal connectivity, which was, however, primarily formulated on non-human primate data [1,16, see 17 for a recent review]. Although post mortem methods are the gold standard for describing microstructure per se, they cannot be mapped directly to function in vivo, making it hard to directly quantify how these metrics relate to neural function. With the advent of high-field magnetic resonance imaging (MRI), it has become possible to probe microstructural properties of different cortical regions in the living human brain. In particular, myelin-sensitive imaging contrasts can differentiate regions with distinct myeloarchitectural profiles at an individual level [18–22]. In parallel, resting-state functional MRI analysis can identify highly reproducible intrinsic networks formed by cortical areas [23–28]. These studies have highlighted that transmodal cortex is largely composed of two spatially distributed yet functionally cohesive networks—the fronto–parietal network thought to respond to the demands of the moment in a flexible way [29–31] and the default mode network that depends on abstract information from self-generated memory and thought processes [32–34]. How flexible functional connectivity profiles are underpinned by microstructural differentiation remains to be seen. We expected a common sensory–fugal gradient to define microstructural and functional organisation, but we also hypothesised that local differences would signify weaker hierarchical constraints in areas responsible for higher-order cognitive processes.

Core to our analysis was the formulation of a systematic approach that modelled cortico–cortical networks from similarity of microstructure profiles sampled across thousands of points and in the direction of cortical columns. The model was first developed on an ultra-high–resolution post mortem 3D histological reconstruction of an entire human brain [35], and we show robust evidence for a principal spatial axis of gradual cytoarchitectural variation running from primary sensory to transmodal areas, recapitulating established spatial trends in laminar differentiation and cytoarchitectural complexity [3,36]. We then translated our approach to myelin-sensitive MRI in a large cohort of healthy adults, showing consistency in vivo and across individuals [37]. In addition to showing correspondence between histological and MRI-derived topographies, microstructure-based gradients only partially converged with macroscale functional topographies derived from task-free functional connectome analysis obtained in the same subjects [28]. In fact, while primary sensory regions served as a common anchor of microstructural and functional gradients, the microstructural axis depicted a progression toward limbic cortices, while its functional counterpart traversed toward default mode and fronto–parietal networks. Critically, meta-analytic decoding revealed that dissociation of functional from microstructural gradients was related to patterns of higher-order thought, such as a cognitive control or social cognition. Together, our analyses on the convergence and divergence of spatial trends in microstructure and function support the hypothesis that reductions in hierarchical constraints in transmodal cortex is a central mechanism underlying flexible cognitive functions.

Results

Formulation of the histology-based microstructure profile covariance analysis

We modelled cortico–cortical microstructural similarity across a 100-μm resolution Merker-stained 3D histological reconstruction of an entire post mortem human brain (BigBrain; https://bigbrain.loris.ca/main.php [35]) (Fig 1A). Staining intensity profiles, representing neuronal density and soma size by cortical depth, were generated along 160,000 surface points (henceforth, vertices) for each hemisphere (Fig 1B). Profile residuals, obtained after correcting intensity profile data for the y coordinate to account for measurable shifts in intensity in anterior-to-posterior direction (S1 Fig), were averaged within 1,012 equally sized, spatially contiguous nodes [38]. Pairwise correlations of nodal intensity profiles, covaried for average intensity profile, were thresholded at 0, and positive edges were log transformed to produce a microstructure profile covariance (MPCHIST) matrix; in other words, MPCHIST captures cytoarchitectural similarity between cortical areas (see S2 Fig for distribution of values).

[Figure omitted. See PDF.]

Fig 1. Histology-based MPCHIST analysis.

(A) Pial (purple) and WM (yellow) surfaces displayed against a sagittal slice of the BigBrain (left) and with the midsurface (blue) in magnified view (right). (B) Mean and SD in residual intensity at each node are displayed on the cortex (left). Cortex-wide intensity profiles were calculated by systematic intensity sampling across intracortical surfaces (rows) and nodes (columns). (C) The MPCHIST matrix depicts node-wise partial correlations in intensity profiles, controlling for the average intensity profile. Exemplary patterns of microstructural similarity from S1, ACC, V1, and the temporal pole. Seed nodes are shown in white. Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). ACC, anterior cingulate cortex; HIST, histology-based; MPC, microstructure profile covariance; S1, primary somatosensory; V1, primary visual; WM, white matter.

https://doi.org/10.1371/journal.pbio.3000284.g001

The pipeline was optimised with respect to the number of intracortical surfaces based on matrix stability (see Methods and S3 Fig). While microstructural similarity had a small but significant statistical relationship with spatial proximity (adjusted R2 = 0.02, P < 0.001), similar findings were obtained after correcting for geodesic distance.

The principal gradient of microstructural similarity reflects sensory–fugal neurostructural variation

Diffusion map embedding, a nonlinear dimensionality reduction algorithm [39] recently applied to identify an intrinsic functional segregation of cortical regions based on resting-state functional MRI [28], was applied to the histology-based MPCHIST matrix (Fig 2A). The relative positioning of nodes in this embedding space informs on (dis)similarity of their covariance patterns. The first principal gradient (G1HIST), accounting for 14.5% of variance, was anchored on one end by primary sensory and motor areas and on the other end by transmodal association and limbic cortices (Fig 2B; see S4 Fig for the second gradient and S5 Fig for results on inflated cortical surfaces). G1HIST depicted the most distinguishable transition in the shape of microstructure profiles (Fig 2B, right). Regions of the prefrontal cortex (green) expressed an intracortical profile that was closest to the cortex-wide average. Extending outward from the centre of G1HIST, sensory and motor regions (blue-purple) exhibited heightened cellular density around the midsurface, whereas paralimbic cortex (red) displayed specifically enhanced density near the cortical borders. For further validation of the biological basis of G1HIST, we mapped independent atlases of laminar differentiation [40] and cytoarchitectural class [13,41] onto the BigBrain midsurface (Fig 2C).

[Figure omitted. See PDF.]

Fig 2. The G1HIST of the histology-based MPCHIST.

(A) Identification: the MPCHIST matrix was transformed into an affinity matrix, which captures similarities in MPCHIST between nodes; this affinity matrix was subjected to diffusion map embedding, a nonlinear compression algorithm that sorts nodes based on MPCHIST similarity. (B) Variance explained by embedding components (left). The first component, G1HIST, describes a gradual transition from primary sensory and motor (blue) to transmodal and limbic areas (red), corresponding to changes in intensity profiles, illustrated with the mean residual intensity profile across 10 discrete bins of the gradient (right). (C) Spatial associations between G1HIST and levels of laminar differentiation (left; [40]) and cytoarchitectural taxonomy (right; [13,41]), ordered by median. Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). G1, first principal gradient; HIST, histology-based; MPC, microstructure profile covariance.

https://doi.org/10.1371/journal.pbio.3000284.g002

Multiple linear regression analyses showed that levels of laminar differentiation and cytoarchitectural taxonomy each accounted for 17% of variance in G1HIST (S2–S3 Tables). Strongest predictors were idiotypic (β = 0.06, P < 0.001) and paralimbic (β = 0.05, P < 0.001) laminar differentiation levels and limbic (β = 0.07, P < 0.001) as well as motor (β = 0.06, P < 0.001) classes in the cytoarchitectural model, demonstrating the cytoarchitectural distinctiveness of regions at the extremes of G1HIST.

Development of in vivo MPC analysis

The MPC approach was adapted to in vivo data in single individuals using T1w/T2w MRI, a ratio indexing cortical microstructure and myelin [18] shown to recapitulate sensory–fugal transitions [42]. Multimodal surface-matched T1w/T2w images, pial surfaces, and white matter surfaces were obtained from the minimally preprocessed S900 release of the Human Connectome Project (HCP) [37,43]. We selected a total of 219 unrelated subjects and grouped these randomly into a Discovery sample (n = 110, 66 females, mean ± SD age = 28.8 ± 3.8 years) and a Replication sample (n = 109, 62 females, mean ± SD age = 28.5 ± 3.7 years). For each individual, we systematically generated intracortical surfaces using the same equivolumetric transformation algorithm as for the histological data set [44,45] and aggregated whole-cortex intensity profiles across 64,984 linked vertices that were subsequently parcellated into 1,012 contiguous nodes [38]. We computed pairwise partial correlations between nodal intensity profiles (controlling for the average intensity profile), kept only positive correlations, and log transformed the result to produce a cortex-wide microstructure profile covariance matrix (MPCMRI; see S6 Fig for distribution of values). A group average MPCMRI matrix was calculated across all participants in the Discovery sample. While microstructural similarity estimated in vivo was stronger between proximal nodes, the variance accounted for by geodesic distance was low (adjusted R2 = 0.04, β = −2.57, P < 0.001).

Diffusion map embedding revealed a principal gradient of microstructural differentiation accounting for 13.7% of variance (G1MRI, Fig 3B; see S7 Fig for the second gradient). In line with its histological counterpart, G1MRI was anchored on one end by primary sensory areas and on the other end by limbic regions. Cortex-wide analysis demonstrated a high correlation of G1HIST and G1MRI at a group level (r = 0.63, P < 0.001) as well as at an individual level (0.41 < r < 0.62, mean r = 0.52 ± 0.03), driven by the close spatial correspondence of gradient extremes (Fig 3C, left) and was unrelated to interindividual variance in cortical curvature (r = −0.01, P = 0.867). G1MRI depicts increasing mean myelin content as well as a gradual transition in the relative myelin content around the midsurface (Fig 3B, right). Microstructural profiles of the prefrontal cortex (orange) again resembled the cortex-wide average; however, comparing node ranks across both modalities revealed a shift in prefrontal regions toward the transmodal anchor in G1MRI (Fig 3C, right). This effect appeared to be driven by a downward shift of lateral occipital–parietal areas toward the sensory anchor, owing to heavy myelination relative to their cytoarchitectural complexity [46]. Laminar differentiation and cytoarchitectural taxonomy accounted for 44% and 37% of variance in G1MRI, respectively (S4 and S5 Tables). As in the histological analysis, the paralimbic (β = 0.10, P < 0.001) and idiotypic (β = 0.07, P < 0.001) levels were the strongest predictors within the laminar differentiation model, while motor (β = 0.15, P < 0.001), limbic (β = 0.07, P < 0.001), and primary sensory (β = 0.10, P < 0.001) classes were strong predictors within the cytoarchitectural model.

[Figure omitted. See PDF.]

Fig 3. In vivo MPCMRI.

(A) Left hemisphere pial, mid, and WM surfaces superimposed on a T1w/T2w image (left); whole-cortex intensity profiles were calculated by systematic sampling across surfaces (rows) and vertices and then averaged with each node (columns). Mean at each node (centre right); MPCMRI matrix depicts node-wise partial correlations in intensity profiles, covaried for mean whole-cortex intensity profile (right). (B) Normalised angle matrix sorted by the principal gradient (left); variance explained by the diffusion-embedding components (left centre) and the principal gradient (right centre); mean residual intensity profiles within 10 discrete bins of the gradient (right). (C) Similarity of histological and in vivo gradients (G1HIST, G1MRI) shown in a density plot (left) and node-wise rank differences shown on the cortical surfaces (right). (D) Associations of G1MRI to levels of laminar differentiation [40] and cytoarchitectural class [13,41] ordered by median. In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). GI, first principal gradient; HCP, Human Connectome Project; HIST, histology-based; MPC, microstructure profile covariance; MRI, magnetic resonance imaging; WM, white matter.

https://doi.org/10.1371/journal.pbio.3000284.g003

Correspondence of microstructural similarity with macroscale functional organisation

To examine the role of microstructural similarity in macroscale functional organisation, we generated a group average resting-state functional connectome across the Discovery subsample and derived gradients with diffusion map embedding. As shown previously [28], the principal functional gradient (G1FUNC) extends from primary sensory and motor networks, through dorsal attention and salience networks, to finally culminate in the transmodal core composed of fronto–parietal and default mode networks (Fig 4A). At the group level, cortex-wide analyses demonstrated moderate-to-high correlations of G1MRI with G1FUNC (r = 0.52, P < 0.001; Fig 4B) and G1HIST with G1FUNC (r = 0.31, P < 0.001), illustrating a common topography of MPC and functional connectivity. Additionally, the topographies of microstructure and function were more closely related than node-to-node correspondence of microstructural similarity with functional connectivity (in vivo: r = 0.10, P < 0.001; histology-based: r = 0.11, P < 0.001), supporting the utility of connectivity-informed dimensionality reduction techniques to reveal common principles of sensory–fugal cortical organisation [47].

[Figure omitted. See PDF.]

Fig 4. Cross-modal correspondence of the MPCMRI and intrinsic functional gradients.

(A) Transformation from individual functional connectomes (left) to a group average normalised angle matrix (centre) to diffusion-embedding components (right). (B) The group-level G1MRI (left), group-level G1FUNC (centre), and density plots depicting the correlation between the gradients (right). (C) Consistency across three example subjects and (D) interindividual variability of the gradients and cross-modal correspondence in the Replication data set. In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). FUNC, functional; HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.g004

The Replication data set underwent identical processing procedures as the Discovery data set. At a group level, G1MRI derived from the Replication and Discovery cohorts were highly correlated (r = 0.98, P < 0.001). High correlations between G1MRI and G1FUNC were also evident at the individual subject level, following alignment of individual functional and microstructural gradients to templates built from the group Discovery data set (see Methods for details). For every participant, we observed a significant correlation between G1MRI and G1FUNC (0.37 < r < 0.61, P < 0.001; Fig 4C and 4D).

In addition to studying their commonalities, we assessed unique topographic features of the modality-specific gradients. Cortex-wide nodal rank comparisons highlighted an upward shift in the position of the prefrontal cortex and precuneus in functional gradient space, relative to the microstructure-derived gradient, and a converse downward shift of the posterior inferior temporal and midcingulate cortices (Fig 5A; see S8 Fig for replication at vertex level). These shifts were independent of regional differences in the interindividual variance of curvature (r = 0.01, P = 0.741). Notably, the greatest rank differences were evident in transmodal cortices, suggesting a specific dissociation of function from microstructure in these higher-order regions. Inspecting the distribution of nodes in intrinsic functional communities [26] along each gradient (Fig 5B), we noted that while the sensory–fugal gradient was overall preserved across all modalities, different sensory and transmodal networks occupied extremes in each modality. G1MRI extended from somatomotor to limbic networks, whereas G1FUNC extended from visual to transmodal default mode networks. Comparing the average node rank of each functional community between modalities and across individuals (S6 Table), we noted that the default mode and fronto–parietal networks shifted to the apex of G1FUNC, reflecting segregation of higher-order communities during rest despite their similar myeloarchitecture.

[Figure omitted. See PDF.]

Fig 5. Divergent representations of the cortical hierarchy derived from microstructure and function.

(A) Differences in nodal ranks between G1MRI (blue) and G1FUNC (red). (B) Radar plot depicting the difference in mean node ranks of functional communities [26] between G1MRI (blue) and G1FUNC (red), with 95% confidence intervals calculated across individuals presented with dotted lines; stacked bar plots depicting the proportion of each bin accounted for by intrinsic functional communities. (C) Meta-analysis maps for diverse cognitive terms were obtained from Neurosynth [28]. We calculated node-wise z-statistics, capturing node-term associations, and calculated the centre of gravity of each term along G1FUNC and G1MRI. The density plots depict the mean difference in the centre of gravity of meta-analysis maps in G1FUNC and G1MRI space across subjects. In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). FUNC, functional; GI, first principal gradient; HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.g005

Our final analysis examined whether the functional topography diverged from microstructure specifically in cortical areas involved in abstract, perceptually decoupled functions. We conducted a meta-analysis using the Neurosynth database and estimated the centre of gravity across a set of diverse cognitive terms [28] along G1FUNC relative to G1MRI (S5C Fig). Top terms exhibiting the strongest upward shift from G1MRI to G1FUNC involved multidomain, integrative functions, such as ‘working memory’, ‘verbal’, ‘cognitive control’, and ‘social control’, whereas the strongest downward shifts involved higher-order visual processing, such as ‘face’ and ‘visual perception’.

Robustness of the MPC and resultant gradients

A series of replication analyses assessed the robustness of the MPC framework.

1. Histology: G1HIST was highly consistent across variations in thresholding, parcellation, intracortical surface number, and BigBrain voxel resolution (S3 Fig). Furthermore, varying the α parameter for diffusion map–embedding algorithm from 0 and 1 (in increments of 0.1) resulted in virtually identical gradients (all r > 0.99). We also compared diffusion map derived gradients to more conventional graph theoretical characterisations of the MPCHIST matrix. Specifically, we applied Louvain modularity detection [48], incorporating a procedure that maximized consensus across fine- and coarse-grained decompositions (by varying the tuning parameter γ from 0.5–1.5) across 100 repetitions. We resolved three modules (0.49 < Q < 0.71; S9 Fig), which occupied distinct positions along the first two gradients. Diffusion map embedding thus represents a continuous and complementary approach to order microstructure profile covariance, which emphasises gradual transitions within and between discrete cortical areas that are fundamental to the hierarchical organisation of layer-specific cortical projections [49].

2. In vivo: As for their histological counterpart, G1MRI solutions were robust against variations in processing parameters, including matrix thresholding, parcellation scheme, and surface number (S10 Fig). Again, variation of the α parameter resulted in virtually identical gradients (all r > 0.99). Application of Louvain community detection identified only two modules in the microstructure profile covariance matrix, which concisely halved G1MRI (0.35 < Q < 0.41, S11 Fig). We replicated G1MRI in an independent data set of healthy late adolescents/young adults [50] that underwent magnetisation transfer (MT) imaging (spatial correlation: r = 0.78, P < 0.001; S12 Fig) and in a cohort of 17 healthy adults scanned at our imaging centre in whom quantitative T1 relaxation data was available (spatial correlation: r = 0.81, P < 0.001; S13 Fig), demonstrating robustness of the approach to acquisition site, acquisition type, and surface construction.

Discussion

Cortical areas identified by classical neuroanatomical studies represent discrete regions embedded within gradual transitions in cyto- and myeloarchitecture [11–13,51]. Cortical gradients are nearly ubiquitous across microstructural and functional domains of the mammalian neocortex, for which mounting evidence supports a common and overarching ‘sensory–fugal’ organisation [49,52]. Capitalising on a ultra-high–resolution 3D reconstruction of an entire post mortem human brain [35], we sampled microstructure profiles and utilised unsupervised techniques to identify smooth transitions in cytoarchitectural composition. We translated our approach from histology to myelin-sensitive in vivo MRI and recovered similar microstructural gradients, which were consistent across individuals. Collectively, our findings support a common sensory–fugal axis of interareal cortical differentiation across histology and in vivo microstructure. Importantly, we also observed a divergence between microstructural and functional gradients toward transmodal default mode and fronto–parietal networks. While these functional networks were distributed across the microstructural gradient, they were localised to the top of the functional gradient. Based on meta-analytical decoding, this divergence was found to relate to aspects of human cognition such as cognitive control and social cognition, supporting views that these regions may be relatively untethered from hierarchical constraints that allow them to take on more flexible cognitive roles.

To bridge different scales of human brain organisation [53], we developed a cortico–cortical network model based on microstructural similarity between areas. Microstructure profile generation, based on state-of-the-art equivolumetric surface construction techniques [44], provided a simplified recapitulation of cellular changes across the putative laminar structure of the cortex. Importantly, our covariance framework generated cytoarchitecturally grounded networks, which were sensitive to laminar thickness as well as cell size and density. The current study revealed gradual cytoarchitectural transitions across the folded neocortex with primary sensory areas on one end of the gradient and cingulate and parahippocampal cortices on the other end. These architectonically unique regions anchor the observed macroscopic sensory–fugal axis of microstructural differentiation, reflecting a continuum through distinct levels of laminar differentiation and cytoarchitectural complexity from the koniocortical primary sensory areas via unimodal association areas toward dysgranular paralimbic cortices [40,54]. The findings strongly support prior evidence on the successive steps of change in cortical architecture and connectivity in reptiles [55], monotremes [56], marsupials [57,58], cetaceans [59], prosimians [60], squirrel monkeys [61], rhesus monkeys [4,62,63], and humans [36,49,64]. Despite differences in cytoarchitectural nomenclature across prior research (see also S1 Table), our findings also mirror the watercolour illustrations of cortical gradients created by Bailey and von Bonin (1951), which highlight the distinctiveness of these gradient anchors.

We directly translated our approach to myelin-sensitive in vivo MRI [18–22], recovering highly consistent topographic maps. Despite there not being a universal 1:1 mapping between cyto- and myeloarchitecture, they are closely related and have complimentary aspects of the same cortical wiring system [14,65,66]. The size, density, and position of pyramidial neurons determine the total intracortical myelin content as well as myelin banding. Histological studies have shown how this microscale correspondence confers similar macroscale topographies, whereby laminar differentiation and mean intracortical myelin gradually change along the sensory–fugal gradient [54,67–69]. Our findings extend upon this work by showing that the global topography of myeloarchitectural similarity, which involves both mean myelin and myelin banding [66], is strongly related to macroscale cytoarchitectural differentiation.

Microstructural and functional gradients were consistently anchored by primary sensory regions on the lower end, recapitulating their consistent location at the bottom of functional and structural hierarchies in primates [3,7,28,36,70,71]. In primary regions, where neural processes are strongly constrained by sensorimotor contingencies, information processing is driven by extrinsic inputs and interactions with the outside world [72,73] as well as intrinsic signalling molecules involved in axonal guidance, cell adhesion, regional circuit formation, and thus cortical patterning [74]. The impact of these forces decreases with synaptic distance, producing a sensory–fugal axis of neurostructural and functional differentiation along the cortical mantle. In fact, given the default mode network’s location as the most distant functional system relative to the location of primary sensory sulci [28], it may be relatively untethered by the influence of intrinsic signalling molecules and extrinsic activity [75] and thus expresses less unique microstructural profiles, as observed here and in classical neuroanatomical studies [11–13].

By showing a close correlation between individualised microstructural gradients and corresponding functional connectome gradients [28], our findings expand the structural model of connectivity to living humans. Originally formulated in nonhuman primates [1,2] and later demonstrated in rodents [76] and cats [77,78], the structural model predicts increased projections between areas with more similar internal architecture. This principle underpins concomitant variations in cytoarchitecture and patterns of connectivity along the sensory–fugal axis and, despite providing a parsimonious account of cortical connectivity, had not been previously explored in living humans or using functional connections. A recent study carried out diffusion MRI network construction and showed an association between structural connectivity strength in selected regions to microstructural similarity information obtained from BigBrain [79]. While this study demonstrated the feasibility of integrating BigBrain data with in vivo connectomics, our study had a different scope. Firstly, we formulated a novel step-by-step procedure to derive networks of microstructural similarity, both based on BigBrain as well as microstructural MRI, bringing a new network modality to connectomics. Of note, MPC networks have whole-cortex coverage, and their construction accounts for curvature using state-of-the-art equivolumetric surface transformation techniques. As we have shown, the MPC procedure can be consistently translated from histology to high-resolution in vivo imaging, thus allowing a direct comparison of microstructural similarity networks with those based on functional connectivity and thus probes structure–function associations in the same subjects.

In contrast to a close functional coupling of microstructurally similar areas at the low end of the cortical gradient, the functional connectivity between higher-order regions may not be as strongly constrained by hierarchical principles. For example, prefrontal and parietal association regions bilaterally project feedforward and feedback connections [80], which can be densely interdigitated [81]. The diverse connectivity profiles in transmodal cortices are likely related to heightened synaptic plasticity [69] that enables more flexible reconfigurations of functional relationships. In contrast, synaptic plasticity is lower in sensory cortices in which quick, accurate response to external stimuli necessitate a constrained hierarchical organisation [3]. Another plausible mechanism for reduced influence of microstructural homophily in higher-order networks may be the increased relevance of long-range structural connections [82]. Such wider-ranging structural connections may accommodate the more distributed spatial layout and diverse functional roles of transmodal cortices. Conversely, primary sensory and motor regions exhibit more locally clustered short-range connectivity profiles related to microstructural similarity [83,84], likely in accordance with their more specialized and stable functional roles.

Our approach was robust with respect to variations in algorithmic and analytical choices. While the singular nature of the BigBrain data set prohibited replication of the histological pipeline, we demonstrated consistency of neuroimaging-derived microstructural gradients across three cohorts with unique myelin-sensitive contrasts [18,20,22,37,50]. Importantly, the in vivo approach can serve as a lower resolution yet biologically meaningful extension of the histological work. In fact, myelin-based gradients exhibited a comparable association with laminar differentiation and cytoarchitectural complexity as the histology-based gradient, as shown to be the case by earlier post mortem work [66–69]. Albeit replication of our findings at a vertex-wise level, the millimetre resolution of in vivo imaging may not always capture detailed spatial features visible on histology. With the emergence of algorithms that detect cortical laminae based on histological data [85] and increasing availability of ultra-high–field MRI scanners (at field strengths of 7T and higher), we expect that the histological and in vivo MPC approach may be further refined in future work, providing an even more direct bridge between the micro- and macroscales of human brain organisation. As it offers a hierarchy-dependent reference frame to examine the interplay of microstructural and system-level network mechanisms in single individuals, the proposed framework may also be advantageous in neurodevelopmental research and complement existing covariance network–mapping approaches that have shown promise in studying typical and atypical development [86–91]. As such, it represents a powerful and freely available (https://github.com/MICA-MNI/micaopen/tree/master/MPC) tool to investigate coordinated changes in cortical microstructure paralleling the emergence and maturation of large-scale functional networks during development [92,93] and may conversely provide insight into the structural underpinnings of atypical network configurations in complex neurodevelopmental disorders [94].

We close by considering the significance of the observed progressive dissociation between structural and functional topographies in cortical areas involved in multidomain, integrative processing. Components of the fronto–parietal network have been shown to guide behaviour in an adaptive manner, changing in line with external demands [29] and possibly substantiated by its ability to dynamically reconfigure its functional connectivity [95]. Although the functional role of the default mode network as the putative apex of the cortical functional hierarchy in primates remains subject to debate [96], it is evidently involved in a broad class of memory-driven operations, involving self-referential and simulative thought processes and some degree of abstraction [32–34]. The interesting possibility that reduced hierarchical constraints enable functional diversity and flexibility was also supported by ad hoc meta-analysis, suggesting that microstructural and functional gradients decouple in regions contributing to processes such as working memory, social cognition, and cognitive control. Such a hypothesis provides a potential mechanistic account for why some of the more creative acts of the human mind emerge through the interaction of the two most dominant yet structurally diverse functional systems [97].

Methods

Histology-based MPC

Histological data acquisition and preprocessing.

An ultra-high–resolution Merker-stained 3D volumetric histological reconstruction of a post mortem human brain from a 65-year-old male was obtained from the open-access BigBrain repository on February 2, 2018 (https://bigbrain.loris.ca/main.php [35]). The post mortem brain was paraffin-embedded, coronally sliced into 7,400 20-μm sections, silver-stained for cell bodies [98], and digitised. Manual inspection for artefacts (i.e., rips, tears, shears, and stain crystallisation) was followed by automatic repair procedures, involving nonlinear alignment to a post mortem MRI, intensity normalisation, and block averaging [99]. 3D reconstruction was implemented with a successive coarse-to-fine hierarchical procedure [100]. We downloaded the 3D volume at four resolutions, with 100-, 200-, 300-, and 400-μm isovoxel size. We primarily analysed 100-μm data and used 200-, 300-, and 400-μm data to assess consistency of findings across spatial scales. Computations were performed on inverted images, on which staining intensity reflects cellular density and soma size. Geometric meshes approximating the outer and inner cortical interface (i.e., the GM/CSF boundary and the GM/WM boundary) with 163,842 matched vertices per hemisphere were also available [101].

Histology-based MPC analysis.

1) Surface sampling. We systematically constructed 10–100 equivolumetric surfaces in steps of 1 between the outer and inner cortical surfaces [45]. The equivolumetric model compensates for cortical folding by varying the Euclidean distance ρ between pairs of intracortical surfaces throughout the cortex to preserve the fractional volume between surfaces [44]. ρ was calculated as follows for each surface:(1)in which α represents fraction of the total volume of the segment accounted for by the surface, while Aout and Ain represent the surface area of the outer and inner cortical surfaces, respectively. Next, vertex-wise microstructure profiles were estimated by sampling intensities along linked vertices from the outer to the inner surface across the whole cortex. In line with previous work [85], layer 1 was approximated as the top 10% of surfaces and removed from the analysis due to little inter-regional variability. Note, however, that findings were nevertheless virtually identical when keeping the top 10% of surfaces. To reduce the impact of partial volume effects, the deepest surface was also removed. Surface-based linear models, implemented via SurfStat for Matlab (http://mica-mni.github.io/surfstat) [102], were used to account for an anterior–posterior increase in intensity values across the BigBrain due to coronal slicing and reconstruction [35], whereby standardised residuals from a simple linear model of surface-wide intensity values predicted by the midsurface y coordinate were used in further analyses.

2) MPC matrix construction. Cortical vertices were parcellated into 1,012 spatially contiguous cortical ‘nodes’ of approximately 1.5 cm2 surface area, excluding outlier vertices with median intensities more than three scaled median absolute deviations away from the node median intensity. The parcellation scheme preserves the boundaries of the Desikan Killany atlas [38] and was transformed from conte69 surface to the BigBrain midsurface via nearest neighbour interpolation. Nodal intensity profiles underwent pairwise Pearson product–moment correlations, controlling for the average whole-cortex intensity profile. MPCHIST for a given pair of nodes i and j was thus(2)in which rij is the Pearson product–moment correlation coefficient of the BigBrain intensity profiles at nodes i and j, ric the correlation coefficient of the intensity profile at node i with the average intensity profile across the entire cortex, and rjc the Pearson correlation of the intensity profile at node j with the average intensity profile across the whole brain. The MPC matrix was thresholded above zero, and remaining MPC values were log-transformed to produce a symmetric 1,012 ☓ 1,012 MPCHIST matrix. The in-house developed code for MPC construction is available online (https://github.com/MICA-MNI/micaopen/tree/master/MPC).

3) Parameter estimation. The optimal surface number was determined based on the stability of the MPC matrix. This procedure involved (repeatedly and randomly) dividing the vertex intensity profiles within each node into two groups and constructing two MPC matrices, then calculating the Euclidean distance between them. The procedure was repeated 1,000 times. Although the MPC matrix instability was robust to variations in surface number, the 18-surface solution exhibited a noticeable local minimum MPC instability in the studied range (10–100 surfaces) and was used in subsequent analyses (S2 Fig). Notably, the MPC gradient was similar using two finer grained solutions (i.e., 54 and 91 surfaces), in which local minima were observed as well. More details on the origins of the stability statistic in clustering algorithms may be found elsewhere [103].

4) Relation to spatial proximity. To determine whether MPCHIST was not purely driven by spatial proximity, we correlated MPCHIST strength with geodesic distance for all node pairs. The latter was calculated using the Fast Marching Toolbox between all pairs of vertices, then averaged by node (https://github.com/gpeyre/matlab-toolboxes/tree/master/).

Histology-based MPC gradient mapping.

In line with previous studies [28,104], the MPCHIST matrix was proportionally thresholded at 90% per row and converted into a normalised angle matrix. Diffusion map embedding [39], a nonlinear manifold learning technique, identified principal gradient components, explaining MPCHIST variance in descending order (each of 1 × 1,012). In brief, the algorithm estimates a low-dimensional embedding from a high-dimensional affinity matrix. In this space, cortical nodes that are strongly interconnected by either many suprathreshold edges or few very strong edges are closer together, whereas nodes with little or no intercovariance are farther apart. The name of this approach, which belongs to the family of graph Laplacians, derives from the equivalence of the Euclidean distance between points in the embedded space and the diffusion distance between probability distributions centred at those points. Compared to other nonlinear manifold learning techniques, the algorithm is relatively robust to noise and computationally inexpensive [105,106]. Notably, it is controlled by a single parameter α, which controls the influence of the density of sampling points on the manifold (α = 0, maximal influence; α = 1, no influence). In this and previous studies [28,104], we followed recommendations and set α = 0.5, a choice that retains the global relations between data points in the embedded space and has been suggested to be relatively robust to noise in the covariance matrix. Gradients were mapped onto BigBrain midsurface visualised using SurfStat (http://mica-mni.github.io/surfstat) [102], and we assessed the amount of MPCHIST variance explained. To show how the principal gradient in MPCHIST (G1HIST) relates to systematic variations in microstructure, we calculated and plotted the mean microstructure profiles within ten equally sized discrete bins of G1HIST.

Relation of G1HIST to laminar differentiation and cytoarchitectural taxonomy.

We evaluated correspondence of G1HIST to atlas information on laminar differentiation and cytoarchitectural class. To this end, each cortical node was assigned to one of four levels of laminar differentiation (i.e., idiotypic, unimodal, heteromodal, or paralimbic) derived from a seminal model of Mesulam, which was built on the integration of neuroanatomical, electrophysiological, and behavioural studies in human and nonhuman primates [40] and one of the seven Von-Economo/Koskinas cytoarchitectural classes (i.e., primary sensory, secondary sensory, motor, association 1, association 2, limbic, or insular) [13,41]. In the case of laminar differentiation maps, assignment was done manually; in the case of cytoarchitectural classes, we mapped previously published Von Economo/Koskinas classes [91] to the BigBrain midsurface with nearest neighbour interpolation and assigned nodes to the cytoarchitectural class most often represented by the underlying vertices. Finally, we estimated the contribution of level of laminar differentiation (D, a categorical variable) and cytoarchitectural class (C, a categorical variable) to the principal gradient G1HIST of the MPCHIST within two separate multiple regression models:(3)(4)

We evaluated model fit via adjusted R2 statistics and unique variances explained by each predictor (β).

In vivo MPC

MRI data acquisition and preprocessing.

We studied data from 219 unrelated healthy adults from the minimally preprocessed S900 release of the HCP [43]. The Discovery data set included 110 individuals (66 females, mean ± SD age = 28.8 ± 3.8 years) and the Replication data set 109 (62 females, mean ± SD age = 28.5 ± 3.7 years). MRI data were acquired on the HCP’s custom 3T Siemens Skyra equipped with a 32-channel head coil. Two T1w images with identical parameters were acquired using a 3D-MPRAGE sequence (0.7 mm isotropic voxels, matrix = 320 × 320, 256 sagittal slices; TR = 2,400 ms, TE = 2.14 ms, TI = 1,000 ms, flip angle = 8°; iPAT = 2). Two T2w images were acquired using a 3D T2-SPACE sequence with identical geometry (TR = 3,200 ms, TE = 565 ms, variable flip angle; iPAT = 2). Four rs-fMRI scans were acquired using multiband accelerated 2D-BOLD echo-planar imaging (2 mm isotropic voxels, matrix = 104 × 90, 72 sagittal slices; TR = 720 ms, TE = 33 ms, flip angle = 52°; mb factor = 8; 1,200 volumes/scan). Participants were instructed to keep their eyes open, look at fixation cross, and not fall asleep. While T1w and T2w scans were acquired on the same day, rs-fMRI scans were split over two days (two scans/day).

Structural and resting-state functional MRI data underwent HCP’s minimal preprocessing [43,107,108]. For structural MRI, images underwent gradient nonlinearity correction. When repeated scans were available, these were coregistered and averaged. Following brain extraction and readout distortion correction, T1w and T2w images were coregistered using rigid body transformations. Subsequently, nonuniformity correction using T1w and T2w images was applied [109]. Preprocessed images were nonlinearly registered to MNI152 space, and cortical surfaces were extracted using FreeSurfer 5.3.0-HCP [110–112], with minor modifications to incorporate both T1w and T2w [18]. Cortical surfaces in individual participants were aligned using MSMAll [113,114] to the hemisphere-matched conte69 template [115]. T1w images were divided by aligned T2w images to produce a single volumetric T1w/T2w image per subject [18]. Notably, this contrast nullifies inhomogeneities related to receiver coils and increases sensitivity to intracortical myelin.

For rs-fMRI, the timeseries were corrected for gradient nonlinearity and head motion. The R-L/L-R blipped scan pairs were used to correct for geometric distortions. Distortion-corrected images were warped to T1w space using a combination of rigid body and boundary-based registrations [116]. These transformations were concatenated with the transformation from native T1w to MNI152 to warp functional images to MNI152. Further processing removed the bias field (as calculated for the structural image), extracted the brain, and normalised whole-brain intensity. A high-pass filter (>2,000s FWHM) corrected the timeseries for scanner drifts, and additional noise was removed using ICA-FIX [117]. Tissue-specific signal regression was not performed [118,119]. We finally transformed these rs-fMRI to native space and sampled time-series at each vertex of the MSMAll-registered [113,114] midthickness cortical surfaces.

In vivo MPC analysis.

We estimated MPC in the in vivo data set (MPCMRI) in the same manner as MPCHIST with the only adjustment that intensity profiles were not corrected for y coordinates; instead, the contrast reversal for T1w and T2w data was used to correct for inhomogeneity as part of the HCP minimal processing pipeline. We generated equivolumetric surfaces between the outer and inner cortical surfaces (see Eq 1) and systematically sampled T1w/T2w values along 64,984 linked vertices from the outer to the inner surface across the whole cortex. In turn, MPCMRI can be denoted as an extension of Eq 2, in which MPCMRI(i,j) for a given pair of nodes i and j is defined by(5)in which s is a participant and n is the number of participants. We systematically evaluated matrix stability with four to 30 intracortical surfaces and selected 14 surfaces as the most stable solution (S10 Fig).

In vivo MPCMRI gradient: Relation to the sensory–fugal gradient and the histological MPC gradient.

As for the histological data, diffusion map embedding derived a principal gradient (G1MRI) from the group average MPCMRI matrix. Correspondence of the in vivo G1MRI to the histological G1HIST gradient was estimated via Spearman rank correlation between spatially matched nodes. To localise differences, we calculated the difference in rank of each node. To dispel potential confounds due to regional differences in the variance in curvature across participants, we calculated average curvature of each node for each subject, then estimated the correlation between node-wise standard deviation in curvature and the G1HIST and G1MRI difference map. As before, we assessed the contribution of level of laminar differentiation and cytoarchitectural class to the G1MRI via multiple regression.

Correspondence of microstructure and functional connectivity

Individual functional connectomes were generated by averaging preprocessed timeseries within nodes, correlating nodal timeseries and converting them to z scores. For each individual, the four available resting state scans were averaged at the matrix level. Then, a group average functional connectome was calculated across the Discovery cohort. Correlation coefficients were calculated between the group average functional connectome and the MPCHIST and MPCMRI matrices. The group average functional connectome was proportionally thresholded at 90% per row, transformed into a cosine similarity matrix, transformed into a normalised angle matrix, then diffusion map embedding was applied, producing G1FUNC. We calculated the correspondence of the G1FUNC with G1HIST and G1MRI with Spearman rank correlations. The differences in the gradients were localised by comparing node ranks across the whole cortex and within functional communities. Differences were also calculated using vertex-wise gradients to ensure the effects were not confounded by averaging disparate microstructural or functional profiles within nodes. Again, we also assessed associations to across-subject variance in curvature and the G1MRI and G1FUNC difference map. Seven functional communities were mapped onto the conte69 surfaces from a previous parcellation [26] with nearest neighbour interpolation from fsaverage5, then nodes were assigned to the functional community most often represented by the underlying vertices. To aid interpretation of the modality-specific gradients, G1MRI and G1FUNC were discretised into 50 equally sized bins, and we calculated the proportion of each bin accounted for by each functional community, then performed seven paired t tests contrasting average node rank of a functional community in G1MRI and G1FUNC across individuals.

We assessed which cognitive faculties are related to dissociations between G1MRI and G1FUNC using meta-analytic maps from Neurosynth [120]. Neurosynth (http://www.neurosynth.org) combines automated text mining and meta-analytical techniques to produce probabilistic mappings between cognitive terms and spatial brain patterns. We downloaded and surface projected meta-analytic z-statistic maps of 24 terms covering a wide range of cognitive functions, which correspond to the topic names defined by Margulies and colleagues (2016). We discretised G1MRI and G1FUNC into five-percentile bins and, for each term, calculated the mean z-statistic within each bin. From this, we deduced the centre of gravity of each term within gradient space and calculated the differences between the G1MRI centre of gravity and G1FUNC centre of gravity for each term using subject-specific gradients.

Robustness of the MPC approach

Individualised MPCMRI gradients and relation to the individual-specific functional hierarchy.

Inter-individual consistency of G1MRI was assessed in the Replication data set. The MPCMRI pipeline was deployed at an individual level, thus resolving individualised MPCMRI matrices and gradients. Additionally, the diffusion map embedding was employed on functional connectomes to derive individual functional gradients [28]. To ensure the spatial correspondence of individual gradients, the individual gradients from the Replication data set underwent Procrustes linear alignment to the Discovery data set group average embedding. Individual cross-modal coupling was calculated as the Spearman rank correlation between G1MRI and G1FUNC.

Robustness of G1HIST.

We assessed the robustness of G1HIST by altering pipeline parameters, repeating MPCHIST generation and diffusion map embedding, then calculating the Pearson correlation of the modified G1HIST with the original gradient. In particular, we evaluated variable matrix thresholds (i.e., 70%–95%, in steps of 1%), alternative surface number in which MPCHIST matrix instability reached a local minima (i.e., 54- and 91-surface solutions), the voxel resolution of the BigBrain volume (100–400 um), and spatial scale (i.e., vertex- versus parcel-wise construction). For the latter, we correlated nodal gradient values with the median vertex within each parcel.

Robustness of G1MRI.

We repeated the robustness procedures reported for the histological gradients in the in vivo data set, including variation of thresholding level, parcellation usage, and surface number. Here, the pipeline was repeated with 23 surfaces, pertaining to local minima in the in vivo MPCMRI matrix instability.

Independent replications of G1MRI.

We replicated the in vivo gradient in two independent data sets based on two additional myelin sensitive magnetic resonance imaging contrasts.

1. Quantitative T1. We implemented the MPC approach on 17 healthy adults (5 females, mean ± SD age = 28.1 ± 6.1, 2 left-handed) for whom quantitative T1 relaxation time–mapping (qT1) images were available. All participants gave informed consent, and the study was approved by the local research ethics board of the Montreal Neurological Institute and Hospital. MRI data was acquired on a 3T Siemens Magnetom Prisma-Fit with a 64-channel head coil. A submillimetric T1-weighted image was acquired using a 3D-MPRAGE sequence (0.8 mm isotropic voxels, 320 x 320 matrix, 24 sagittal slices, TR = 2,300 ms, TE = 3.14 ms, TI = 900 ms, flip angle = 9°, iPAT = 2), and qT1 data was acquired using a 3D-MP2RAGE sequence (0.8 mm isotropic voxels, 240 sagittal slices, TR = 5,000 ms, TE = 2.9 ms, TI 1 = 940 ms, T1 2 = 2,830 ms, flip angle 1 = 4°, flip angle 2 = 5°, iPAT = 3, bandwidth = 270 Hz/px, echo spacing = 7.2ms, partial Fourier = 6/8). The combination of two inversion images in qT1 mapping minimises sensitivity to B1 inhomogeneities [121] and provides high intrasubject and intersubject reliability [122].

Cortical surfaces were extracted from the T1-weighted scans using FreeSurfer 6.0 [110–112], and 14 equivolumetric intracortical surfaces were generated [45]. qT1 was registered to Freesurfer native space using a boundary-based registration [116], and FreeSurfer native space was registered to standard conte69 space using Caret5 landmark-based registration [115]. We used the former to sample qT1 intensity values along the intracortical surfaces and the latter to resample the evaluated surfaces to a common space with 64,984 matched vertices. As in the main approach, we averaged vertex-wise intensity profiles within 1,012 nodes [38], computed pairwise partial correlations between nodal intensity profiles (controlling for the average intensity profile), kept only positive correlations, and log transformed the result to produce a MPCMRI-QT1 matrix. Finally, we generated a group-average MPCMRI-QT1 matrix and applied diffusion map embedding. The similarity of G1MRI- QT1 to the original G1MRI was measured with a node-wise Spearman rank correlation.

1. 2. MT. In an open data set of 297 healthy young adults (149 female, mean ± SD age = 19.1 ± 2.9; [50]), we studied intracortical depth profiles of MT. MT is a validated measure of myelination [22] and was available in the form of average intensity values along eight equidistant intracortical surfaces (10%–90% in 10% intervals) within 308 cortical regions (for further details, see [92]). We similarly applied the MPC framework to the MT profiles, then averaged the MPCMRI-MT matrix across the group and applied diffusion map embedding. We measured the Spearman rank correlation between G1MRI-MT and the original G1MRI, which was recalculated with 308 matched cortical areas.

Supporting information

S1 Fig. Horizontal view of the volumetric reconstruction of the BigBrain with pial, mid, and white matter surfaces projected on the left hemisphere.

Notably, we corrected for the linear relationship between intensity values and midsurface y coordinate (r = −0.68, P < 0.001), which existed due to coronal slicing and reconstruction of the BigBrain. Histological data is openly available as part of the BigBrain initiative (https://bigbrain.loris.ca/main.php).

https://doi.org/10.1371/journal.pbio.3000284.s001

(TIF)

S2 Fig. Distribution of values in MPCHIST matrix.

(Left) Frequency of r values calculated by Pearson product moment correlation coefficient of the nodal intensity profiles, controlling for the average intensity profile. (Right) Frequency of positive z values following log transformation of r values. Histological data is openly available as part of the BigBrain initiative (https://bigbrain.loris.ca/main.php). HIST, histology-based; MPC, microstructure profile covariance.

https://doi.org/10.1371/journal.pbio.3000284.s002

(TIF)

S3 Fig. Robustness of G1HIST to parameter variation.

(A) MPCHIST matrix instability using between 10 and 100 intracortical surfaces. G1HIST was consistent regardless of the number of intracortical surfaces used, as shown by the strong spatial correlation of the 18-, 54- and 91-surface solutions (all r > 0.97, all P < 0.001). (B) Correlation matrix depicting the high correspondence of G1HIST solutions with 70–95% row-wise matrix thresholding (0.95 < r < 1, all P < 0.001). (C) Estimation of G1HIST from 20488 vertices resulted in a consistent G1HIST to the 1,012 parcel construction pipeline (r = 0.87, P < 0.001). (D) The 200-μm, 300-μm, and 400-μm resolution BigBrain data sets were characterised as lower resolution replications, and G1HIST was found to be highly correlated across these resolutions (all r > 0.81, all P < 0.001). Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). G1, first principal gradient; HIST, histology-based; MPC, microstructure profile covariance.

https://doi.org/10.1371/journal.pbio.3000284.s003

(TIF)

S4 Fig. First two principal components of MPCHIST.

(A) The second principal component, accounting for 12.7% of variance in MPCHIST components, is projected on the BigBrain midsurface. (B) Scatterplot depicting the first two embedding gradients, with corresponding probability density functions. The second gradient divides the lower-order areas of the first gradient, insomuch that somatomotor and primary visual areas (red) are separated from ventral prefrontal areas and secondary visual areas (blue). Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). HIST, histology-based; MPC, microstructure profile covariance.

https://doi.org/10.1371/journal.pbio.3000284.s004

(TIF)

S5 Fig. Principal gradients and node-rank differences projected on inflated cortical surfaces.

Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). HCP, Human Connectome Project.

https://doi.org/10.1371/journal.pbio.3000284.s005

(TIF)

S6 Fig. Distribution of values in MPCMRI matrix.

(Left) Frequency of r values calculated by Pearson product–moment correlation coefficient of the nodal intensity profiles, controlling for the average intensity profile. (Right) Frequency of positive z values following log transformation of r values. In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s006

(TIF)

S7 Fig. First two principal components of MPCMRI.

(A) The second principal component, accounting for 11.7% of variance in MPCMRI components, projected on the conte69 midsurface. (B) Scatterplot depicting the first two embedding gradients, with corresponding probability density functions. The second gradient divides the higher-order areas of the first gradient, insomuch that the cingulate, orbitofrontal cortex and the inferior temporal gyrus (red) are separated from the prefrontal cortex, precuneus, temporo–parietal junction, and superior temporal gyrus (blue). In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s007

(TIF)

S8 Fig. Differences in vertex ranks between G1MRI (blue) and G1FUNC (red).

In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). FUNC, functional; G1, first principal gradient; HCP, Human Connectome Project; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s008

(TIF)

S9 Fig. Community structure of the MPCHIST.

(A) MPCHIST matrix sorted by community membership. (B) Modular decomposition of MPCHIST projected on the BigBrain midsurface. (C) Mean ± SD of residual intensity profile for each module, after correction for the midsurface y coordinate and demeaning. (D) Boxplot depicts the unique positions of modules along the first two principal gradients. Histological data is openly available as part of the BigBrain repository (https://bigbrain.loris.ca/main.php). HIST, histology-based; MPC, microstructure profile covariance.

https://doi.org/10.1371/journal.pbio.3000284.s009

(TIF)

S10 Fig. Robustness of G1MRI to parameter variation.

(A) MPCMRI matrix instability using between 4 and 30 intracortical surfaces. G1MRI was consistent regardless of the number of intracortical surfaces used, as shown by the strong spatial correlation of the 14- and 23-surface solutions (r = 0.98, P < 0.001). (B) Correlation matrix depicting the high correspondence of G1MRI solutions with 70%–95% row-wise matrix thresholding (0.91 < r <1, all P < 0.001). (C) Estimation of G1MRI from 20,464 vertices resulted in a consistent G1MRI to the 1,012-parcel construction pipeline (r = 0.98, P < 0.001). In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). G1, first principal gradient; HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s010

(TIF)

S11 Fig. Community structure of the MPCMRI.

(A) MPCMRI matrix sorted by community membership. (B) Modular decomposition of MPCMRI projected on the conte69 midsurface. (C) Mean ± SD of residual intensity profile for each module, after demeaning. (D) Boxplot depicts the unique positions of modules on the first, but not the second, principal component. In vivo imaging data is openly available as part of the HCP S900 release (https://www.humanconnectome.org/study/hcp-young-adult/document/900-subjects-data-release). HCP, Human Connectome Project; MPC, microstructure profile covariance; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s011

(TIF)

S12 Fig. Independent replication of G1MRI, using MT data.

(A) Mean MT intensity across subjects projected onto the cortical surface. (B) Gradient ordered normalised angle matrix. (C) Variance explained by embedding components. (D) G1MRI-MT projected onto conte69 midsurface. (E) Scatterplot depicting the strong correlation between the G1MRI (reconstructed with the 308 parcellation scheme) and G1MRI-MT (r = 0.79, P < 0.001). MT metadata was acquired from the github repository (https://github.com/KirstieJane/NSPN_WhitakerVertes_PNAS2016). G1, first principal gradient; MRI, magnetic resonance imaging; MT, magnetisation transfer.

https://doi.org/10.1371/journal.pbio.3000284.s012

(TIF)

S13 Fig. Independent replication of G1MRI, using qT1 images.

(A) Mean qT1 intensity across subjects projected onto the cortical surface. (B) Gradient ordered normalised angle matrix. (C) Variance explained by embedding components. (D) G1MRI-QT1 projected onto conte69 midsurface. (E) Scatterplot depicting the strong correlation between the G1MRI and G1MRI-QT1 (r = 0.81, P < 0.001). qT1 metadata was acquired in-house and is available on the Github (https://github.com/MICA-MNI/micaopen/tree/master/MPC). G1, first principal gradient; MRI, magnetic resonance imaging; qT1, quantitative T1 relaxometry.

https://doi.org/10.1371/journal.pbio.3000284.s013

(TIF)

S1 Table. Comparison of nomenclature and parcellation of cytoarchitectural groupings.

https://doi.org/10.1371/journal.pbio.3000284.s014

(PDF)

S2 Table. Individual predictors in a multiple regression model of G1HIST by levels of laminar differentiation.

G1, first principal gradient; HIST, histology-based.

https://doi.org/10.1371/journal.pbio.3000284.s015

(PDF)

S3 Table. Individual predictors in a multiple regression model of G1HIST by classes of cytoarchitecture.

G1, first principal gradient; HIST, histology-based.

https://doi.org/10.1371/journal.pbio.3000284.s016

(PDF)

S4 Table. Individual predictors in a multiple regression model of G1MRI by levels of laminar differentiation.

G1, first principal gradient; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s017

(PDF)

S5 Table. Individual predictors in a multiple regression model of G1MRI by classes of cytoarchitecture.

G1, first principal gradient; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s018

(PDF)

S6 Table. Statistical outcome of paired t tests between G1MRI and G1FUNC within each functional community, taken across individuals.

FUNC, functional; G1, first principal gradient; MRI, magnetic resonance imaging.

https://doi.org/10.1371/journal.pbio.3000284.s019

(PDF)

Acknowledgments

The authors would also like to express their gratitude to the teams at the Foschungszentrum Julich and the Montreal Neurological Institute, who made the BigBrain data set available. Furthermore, we thank Dr. Nicola Palamero-Gallagher for helpful and inspiring discussions.

Citation: Paquola C, Vos De Wael R, Wagstyl K, Bethlehem RAI, Hong S-J, Seidlitz J, et al. (2019) Microstructural and functional gradients are increasingly dissociated in transmodal cortices. PLoS Biol 17(5): e3000284. https://doi.org/10.1371/journal.pbio.3000284

References

1. Barbas H. Pattern in the laminar origin of corticocortical connections. J Comp Neurol. 1986 Oct 15;252(3):415–22. pmid:3793985

2. Barbas H, Rempel-Clower N. Cortical structure predicts the pattern of corticocortical connections. Cereb Cortex. 1997;7(7):635–46. pmid:9373019

3. Mesulam M-M. From sensation to cognition. Brain. 1998;121:1013–52. pmid:9648540

4. Jones EG, Powell TPS. An Anatomical Study of Converging Sensory Pathways within the Cerebral Cortex of the Monkey. Brain. 1970;93(4):793–820. pmid:4992433

5. Mesulam M. The evolving landscape of human cortical connectivity: Facts and inferences. Neuroimage. 2012 Oct 1;62(4):2182–9. pmid:22209814

6. Goldman-Rakic PS. Topography of Cognition: Parallel Distributed Networks in Primate Association Cortex. Annu Rev Neurosci. 1988;11(1):137–56.

7. Felleman DJ, Van Essen DC. Distributed hierarchical processing in the primate cerebral cortex. Cereb Cortex. 1991;1(1):1–47. pmid:1822724

8. Selemon LD, Goldman-Rakic PS. Common Cortical and Subcortical Targets of the Dorsolateral Prefrontal and Posterior Parietal Cortices in the Rhesus Monkey: Evidence for a Distributed Neural Network Subserving Spatially Guided Behavior. Vol. 8, The Journal of Neuroscience. 1988.

9. Cavada C, Goldman-Rakic PS. Posterior parietal cortex in rhesus monkey: II. Evidence for segregated corticocortical networks linking sensory and limbic areas with the frontal lobe. J Comp Neurol. 1989 Sep 22;287(4):422–45. pmid:2477406

10. Cavada C, Goldman-Rakic PS. Posterior parietal cortex in rhesus monkey: I. Parcellation of areas based on distinctive limbic and sensory corticocortical connections. J Comp Neurol. 1989 Sep 22;287(4):393–421. pmid:2477405

11. Brodmann K. Vergleichende Lokalisationslehre der Grosshirnrinde in ihren Prinzipien dargestellt auf Grund des Zellenbaues. Leipzig: Barth JA; 1909.

12. Vogt C, Vogt O. Allgemeinere Ergebnisse unserer Hirnforschung. J für Psychol und Neurol. 1919;25(Suppl. 1):273–462.

13. Von Economo C, Koskinas G. Die Cytoarchitektonik der Hirnrinde des erwachsenen Menschen. Berlin: Springer; 1925.

14. Palomero-Gallagher N, Zilles K. Cortical layers: Cyto-, myelo-, receptor- and synaptic architecture in human cortical areas. NeuroImage. 2017 Aug 12;

15. Flechsig PE. Anatomie des menschlichen Gehirns und Rückenmarks auf myelogenetischer Grundlage—Paul Emil Flechsig—Google Books. Leipzig: Georg Thieme; 1920.

16. Beul SF, Barbas H, Hilgetag CC. A Predictive Structural Model of the Primate Connectome. Sci Rep. 2017;7.

17. García-Cabezas MÁ, Zikopoulos B, Barbas H. The Structural Model: a theory linking connections, plasticity, pathology, development and evolution of the cerebral cortex. Brain Struct Funct. 2019 Feb 9;1–24.

18. Glasser MF, Van Essen DC. Mapping human cortical areas in vivo based on myelin content as revealed by T1- and T2-weighted MRI. J Neurosci. 2011 Aug 10;31(32):11597–616. pmid:21832190

19. Bock NA, Kocharyan A, Liu J V., Silva AC. Visualizing the entire cortical myelination pattern in marmosets with magnetic resonance imaging. J Neurosci Methods. 2009 Dec 15;185(1):15–22. pmid:19737577

20. Stüber C, Morawski M, Schäfer A, Labadie C, Wähnert M, Leuze C, et al. Myelin and iron concentration in the human brain: A quantitative study of MRI contrast. Neuroimage. 2014 Jun 1;93(P1):95–106.

21. Geyer S, Weiss M, Reimann K, Lohmann G, Turner R. Microstructural Parcellation of the Human Cerebral Cortex–From Brodmann’s Post-Mortem Map to in vivo Mapping with High-Field Magnetic Resonance Imaging. Front Hum Neurosci. 2011;5(19):1–7.

22. Schmierer K, Tozer DJ, Scaravilli F, Altmann DR, Barker GJ, Tofts PS, et al. Quantitative magnetization transfer imaging in postmortem multiple sclerosis brain. J Magn Reson Imaging. 2007 Jul;26(1):41–51. pmid:17659567

23. Biswal B, Zerrin Yetkin F, Haughton VM, Hyde JS. Functional connectivity in the motor cortex of resting human brain using echo-planar MRI. Magn Reson Med. 1995;34(4):537–41. pmid:8524021

24. Raichle ME. The brain’s default mode network. Annu Rev Neurosci. 2015/05/06. 2015;38:433–47. pmid:25938726

25. Smith SM, Fox PT, Miller KL, Glahn DC, Fox PM, Mackay CE, et al. Correspondence of the brain’s functional architecture during activation and rest. Proc Natl Acad Sci U S A. 2009;106(31):13040–5. pmid:19620724

26. Yeo BTT, Krienen FM, Sepulcre J, Thomas Yeo BT, Krienen FM, Sepulcre J, et al. The organization of the human cerebral cortex estimated by intrinsic functional connectivity. J Neurophysiol. 2011 Sep;106(3):1125–65. pmid:21653723

27. Cauda F, D’Agata F, Sacco K, Duca S, Geminiani G, Vercelli A. Functional connectivity of the insula in the resting brain. Neuroimage. 2011 Mar 1;55(1):8–23. pmid:21111053

28. Margulies DS, Ghosh SS, Goulas A, Falkiewicz M, Huntenburg JM, Langs G, et al. Situating the default-mode network along a principal gradient of macroscale cortical organization. Proc Natl Acad Sci. 2016;113(44):12574–9. pmid:27791099

29. Duncan J. The multiple-demand (MD) system of the primate brain: mental programs for intelligent behaviour. Trends Cogn Sci. 2010 Apr;14(4):172–9. pmid:20171926

30. Dosenbach NUF, Fair DA, Miezin FM, Cohen AL, Wenger KK, Dosenbach RAT, et al. Distinct brain networks for adaptive and stable task control in humans. Proc Natl Acad Sci. 2007 Jun 26;104(26):11073–8. pmid:17576922

31. Corbetta M. Frontoparietal cortical networks for directing attention and the eye to visual locations: identical, independent, or overlapping neural systems? Proc Natl Acad Sci U S A. 1998 Feb 3;95(3):831–8. pmid:9448248

32. Vatansever D, Menon DK, Stamatakis EA. Default mode contributions to automated information processing. Proc Natl Acad Sci U S A. 2017 Nov 28;114(48):12821–6. pmid:29078345

33. Braga RM, Buckner RL. Parallel Interdigitated Distributed Networks within the Individual Estimated by Intrinsic Functional Connectivity. Neuron. 2017 Jul 19;95(2):457–471.e5.

34. Spreng RN, Grady CL. Patterns of Brain Activity Supporting Autobiographical Memory, Prospection, and Theory of Mind, and Their Relationship to the Default Mode Network. J Cogn Neurosci. 2010 Jun;22(6):1112–23. pmid:19580387

35. Amunts K, Lepage C, Borgeat L, Mohlberg H, Dickscheid T, Rousseau M-E, et al. BigBrain: An Ultrahigh-Resolution 3D Human Brain Model. Science (80-). 2013;340(6139):1472–5.

36. Sanides F. The Cyto-myeloarchitecture of the Human Frontal Lobe and its Relation to Phylogenetic Differentiation of the Cerebral Cortex. J Hirnforsch. 1964;6(5):269–82.

37. Van Essen DC, Smith SM, Barch DM, Behrens TEJ, Yacoub E, Ugurbil K. The WU-Minn Human Connectome Project: An overview. Neuroimage. 2013;80:62–79. pmid:23684880

38. Hong S-J, Bernhardt BC, Gill RS, Bernasconi N, Bernasconi A. The spectrum of structural and functional network alterations in malformations of cortical development. Brain. 2017 Aug 1;140(8):2133–43. pmid:28899007

39. Coifman RR, Lafon S, Lee AB, Maggioni M, Nadler B, Warner F, et al. Geometric diffusions as a tool for harmonic analysis and structure definition of data: Multiscale methods. Proc Natl Acad Sci. 2005;102(21):7432–7. pmid:15899969

40. Mesulam M-M. Behavioral neuroanatomy: Largescale networks, association cortex, frontal syndromes, the limbic system, and hemispheric specialization. In: Principles of Behavioral and Cognitive Neurology. 2000. p. 1–120.

41. Vértes PE, Rittman T, Whitaker KJ, Romero-Garcia R, Váša F, Kitzbichler MG, et al. Gene transcription profiles associated with inter-modular hubs and connection distance in human functional magnetic resonance imaging networks. Philos Trans R Soc B Biol Sci. 2016;371(1705).

42. Burt JB, Demirtaş M, Eckner WJ, Navejar NM, Ji JL, Martin WJ, et al. Hierarchy of transcriptomic specialization across human cortex captured by structural neuroimaging topography. Nat Neurosci. 2018 Sep 6;21(9):1251–9. pmid:30082915

43. Glasser MF, Sotiropoulos SN, Wilson JA, Coalson TS, Fischl B, Andersson JL, et al. The minimal preprocessing pipelines for the Human Connectome Project. Neuroimage. 2013 Oct 15;80:105–24. pmid:23668970

44. Waehnert MD, Dinse J, Weiss M, Streicher MN, Waehnert P, Geyer S, et al. Anatomically motivated modeling of cortical laminae. Neuroimage. 2014 Jun 1;93:210–20. pmid:23603284

45. Wagstyl K, Paquola C, Bethlehem R, Huth A. kwagstyl/surface_tools: Initial release of equivolumetric surfaces. 2018.

46. Lewis JW, Van Essen DC. Mapping of architectonic subdivisions in the macaque monkey, with emphasis on parieto-occipital cortex. J Comp Neurol. 2000 Dec 4;428(1):79–111. pmid:11058226

47. Mars RB, Passingham RE, Jbabdi S. Connectivity Fingerprints: From Areal Descriptions to Abstract Spaces. Trends Cogn Sci. 2018 Sep 18;22(11):1026–37. pmid:30241910

48. Lancichinetti A, Fortunato S. Consensus clustering in complex networks. Sci Rep. 2012;2.

49. Goulas A, Zilles K, Hilgetag CC. Cortical Gradients and Laminar Projections in Mammals. Trends Neurosci. 2018 Jul 3;

50. Kiddle B, Inkster B, Prabhu G, Moutoussis M, Whitaker KJ, Bullmore ET, et al. Cohort Profile: The NSPN 2400 Cohort: a developmental sample supporting the Wellcome Trust NeuroScience in Psychiatry Network. Int J Epidemiol. 2018 Feb 1;47(1):18–19g. pmid:29177462

51. Bailey P, Bonin G von. The isocortex of man. Urbana: University of Illinois Press; 1951.

52. Fulcher BD, Murray JD, Zerbi V, Wang X-J. Multimodal gradients across mouse cortex. Proc Natl Acad Sci. 2019 Mar 5;116(10):4689–95.

53. van den Heuvel MP, Yeo BTT. A Spotlight on Bridging Microscale and Macroscale Human Brain Architecture. Vol. 93, Neuron. 2017. p. 1248–51. pmid:28334602

54. Sanides F. Functional Architecture of Motor and Sensory Cortices in Primates in the Light of a New Concept of Neocortex Evolution *. In: Noback CR, Montagna W, editors. The Primate Brain. New York; 1970. p. 136–208.

55. Dart RA. The dual structure of the neopallium: Its history and significance. J Anat. 1934 Oct;69(P1):3.

56. Abbie AA. The excitable cortex in the Monotremata. Sydney: Australian Journal of Experimental Biology and Medical Science; 1938.

57. Abbie AA. Cortical lamination in a polyprotodont marsupial, Perameles nasuta. J Comp Neurol. 1942 Jun 1;76(3):509–36.

58. Abbie AA. The excitable cortex in Perameles, Sarcophilus, Dasyurus, Trichosurus and Wallabia (Macropus). J Comp Neurol J Comp Neurol. 1940;72(3):469–87.

59. Morgane PJ, Glezer II, Jacobs MS. Comparative and Evolutionary Anatomy of the Visual Cortex of the Dolphin. In Springer, Boston, MA; 1990. p. 215–62.

60. Sanides F, Krishnamurti A. Cytoarchitectonic subdivisions of sensorimotor and prefrontal regions and of bordering insular and limbic fields in slow loris (Nycticebus coucang coucang). J Hirnforsch. 1967;9(3):225–52. pmid:6077605

61. Sanides F. The architecture of the cortical taste nerve areas in squirrel monkey (Saimiri sciureus) and their relationships to insular, sensorimotor and prefrontal regions. BRES Brain Res. 1968;8(1):97–124.

62. Mesulam M ‐Marsel Mufson EJ. Insula of the old world monkey. Architectonics in the insulo‐orbito‐temporal component of the paralimbic brain. J Comp Neurol. 1982 Nov 20;212(1):1–22. pmid:7174905

63. Pandya DN, Sanides F. Architectonic parcellation of the temporal operculum in rhesus monkey and its projection pattern. Zeitschrift Anat und Entwicklungsgeschichte. 1973;139(2):127–61.

64. Sanides F. Die Architektonik des menschlichen Stirnhirns zugleich eine Darstellung der Prinzipien seiner Gestaltung als Spiegel der stammgeschichtlichen Differenzierung der Grosshirnrinde. Berlin: Springer; 1962.

65. Braitenberg V, Schüz A. Cyto- and Myeloarchitectonics: Two Aspects of the same Reality. In: Cortex: Statistics and Geometry of Neuronal Connectivity. Berlin, Heidelberg: Springer Berlin Heidelberg; 1998. p. 151–8.

66. Hellwig B. How the myelin picture of the human cerebral cortex can be computed from cytoarchitectural data. A bridge between von Economo and Vogt. J Hirnforsch. 1993;34(3):387–402. pmid:8270790

67. Braitenberg V. A note on myeloarchitectonics. J Comp Neurol. 1962 Apr 1;118(2):141–56.

68. Barbas H, Pandya DN. Architecture and intrinsic connections of the prefrontal cortex in the rhesus monkey. J Comp Neurol. 1989 Aug 15;286(3):353–75. pmid:2768563

69. García-Cabezas M, Joyce MKP, John YJ, Zikopoulos B, Barbas H. Mirror trends of plasticity and stability indicators in primate prefrontal cortex. Eur J Neurosci. 2017 Oct;46(8):2392–405. pmid:28921934

70. Fulcher BD, Murray JD, Zerbi V, Wang X-J. Multimodal gradients across mouse cortex.

71. Huntenburg JM, Bazin PL, Goulas A, Tardif CL, Villringer A, Margulies DS. A Systematic Relationship Between Functional Connectivity and Intracortical Myelin in the Human Cerebral Cortex. Cereb Cortex. 2017;27(2):981–97. pmid:28184415

72. Rakic P. Specification of cerebral cortical areas. Science (80-). 1988;241(4862):170–6.

73. O’Leary DDM, Chou S-J, Sahara S. Area Patterning of the Mammalian Cortex. Neuron. 2007 Oct 25;56(2):252–69. pmid:17964244

74. Fukuchi-Shimogori T, Grove EA. Neocortex patterning by the secreted signaling molecute FGF8. Science (80-). 2001;294(5544):1071–4.

75. Buckner RL, Krienen FM. The evolution of distributed association networks in the human brain. Vol. 17, Trends in Cognitive Sciences. 2013. p. 648–65. pmid:24210963

76. Goulas A, Uylings HBM, Hilgetag CC. Principles of ipsilateral and contralateral cortico-cortical connectivity in the mouse. Brain Struct Funct. 2017 Apr 6;222(3):1281–95. pmid:27497948

77. Hilgetag CC, Grant S. Cytoarchitectural differences are a key determinant of laminar projection origins in the visual cortex. Neuroimage. 2010 Jul 1;51(3):1006–17. pmid:20211270

78. Beul SF, Grant S, Hilgetag CC. A predictive model of the cat cortical connectome based on cytoarchitecture and distance. Brain Struct Funct. 2015 Nov 26;220(6):3167–84. pmid:25062666

79. Wei Y, Scholtens LH, Turk E, van den Heuvel MP. Multiscale examination of cytoarchitectonic similarity and human brain connectivity. Netw Neurosci. 2018 May 29;1–34.

80. Medalla M, Barbas H. Diversity of laminar connections linking periarcuate and lateral intraparietal areas depends on cortical structure. Eur J Neurosci. 2006 Jan;23(1):161–79. pmid:16420426

81. Goldman-Rakic PS, Schwartz ML. Interdigitation of contralateral and ipsilateral columnar projections to frontal association cortex in primates. Science. 1982 May 14;216(4547):755–7. pmid:6177037

82. Hagmann P, Cammoun L, Gigandet X, Meuli R, Honey CJ, Van Wedeen J, et al. Mapping the structural core of human cerebral cortex. Friston KJ, editor. PLoS Biol. 2008 Jul 1;6(7):1479–93.

83. Oligschläger S, Huntenburg JM, Golchert J, Lauckner ME, Bonnen T, Margulies DS. Gradients of connectivity distance are anchored in primary cortex. Brain Struct Funct. 2017 Jul;222(5):2173–82. pmid:27807628

84. Sepulcre J, Liu H, Talukdar T, Martincorena I, Thomas Yeo BT, Buckner RL. The organization of local and distant functional connectivity in the human brain. Sporns O, editor. PLoS Comput Biol. 2010 Jun 10;6(6):1–15.

85. Wagstyl K, Lepage C, Bludau S, Zilles K, Fletcher PC, Amunts K, et al. Mapping Cortical Laminar Structure in the 3D BigBrain. Cereb Cortex. 2018 Jul 1;28(7):2551–62. pmid:29901791

86. Alexander-Bloch A, Raznahan A, Bullmore E, Giedd J. The Convergence of Maturational Change and Structural Covariance in Human Cortical Networks. J Neurosci. 2013;33(7):2889–99. pmid:23407947

87. Zielinski BA, Gennatas ED, Zhou J, Seeley WW. Network-level structural covariance in the developing brain. Proc Natl Acad Sci U S A. 2010;107(42):18191–6. pmid:20921389

88. Bonilha L, Tabesh A, Dabbs K, Hsu DA, Stafstrom CE, Hermann BP, et al. Neurodevelopmental alterations of large-scale structural networks in children with new-onset epilepsy. Hum Brain Mapp. 2014 Aug;35(8):3661–72. pmid:24453089

89. Palaniyappan L, Park B, Balain V, Dangi R, Liddle P. Abnormalities in structural covariance of cortical gyrification in schizophrenia. Brain Struct Funct. 2015 Jul 26;220(4):2059–71. pmid:24771247

90. Bethlehem RAI, Romero-Garcia R, Mak E, Bullmore ET, Baron-Cohen S. Structural Covariance Networks in Children with Autism or ADHD. Cereb Cortex. 2017 Aug 1;27(8):4267–76. pmid:28633299

91. Seidlitz J, Váša F, Shinn M, Romero-Garcia R, Whitaker KJ, Vértes PE, et al. Morphometric Similarity Networks Detect Microscale Cortical Organization and Predict Inter-Individual Cognitive Variation. Neuron. 2018;97(1):231–247.e7.

92. Whitaker KJ, Vértes PE, Romero-Garcia R, Váša F, Moutoussis M, Prabhu G, et al. Adolescence is associated with genomically patterned consolidation of the hubs of the human brain connectome. Proc Natl Acad Sci. 2016 Aug 9;113(32):9105–10. pmid:27457931

93. Grydeland H, Walhovd KB, Tamnes CK, Westlye LT, Fjell AM. Intracortical Myelin Links with Performance Variability across the Human Lifespan: Results from T1- and T2-Weighted MRI Myelin Mapping and Diffusion Tensor Imaging. J Neurosci. 2013;33(47):18618–30. pmid:24259583

94. DiMartino A, Fair DA, Kelly C, Satterthwaite TD, Castellanos FX, Thomason ME, et al. Unraveling the miswired connectome: A developmental perspective. Vol. 83, Neuron. 2014. p. 1335–53. pmid:25233316

95. Cole MW, Reynolds JR, Power JD, Repovs G, Anticevic A, Braver TS. Multi-task connectivity reveals flexible hubs for adaptive task control. Nat Neurosci. 2013 Sep 28;16(9):1348–55. pmid:23892552

96. Buckner RL, Margulies DS. Macroscale Cortical Organization and a Default-Like Transmodal Apex Network in the Marmoset Monkey. bioRxiv. 2018 Sep 12;415141.

97. Beaty RE, Kenett YN, Christensen AP, Rosenberg MD, Benedek M, Chen Q, et al. Robust prediction of individual creative ability from brain functional connectivity. Proc Natl Acad Sci U S A. 2018 Jan 30;115(5):1087–92. pmid:29339474

98. Merker B. Silver staining of cell bodies by means of physical development. J Neurosci Methods. 1983;9(3):235–41. pmid:6198563

99. Lepage C, Mohlberg H, Pietrzyk U, Amunts K, Zilles K, Evans A. Automatic repair of acquisition defects in reconstruction of histology slices of the human brain. In: 16th Annual Meeting of the Organization for Human Brain Mapping (OHBM). Barcelona; 2010.

100. Mohlberg H, Tweddell B, Lippert T, Amunts K. Workflows for Ultra-High Resolution 3D Models of the Human Brain on Massively Parallel Supercomputers. In Springer, Cham; 2016. p. 15–27.

101. Lewis LB, Lepage C, Marc F, Zilles K, Amunts K, Evans AC. BigBrain: Initial Tissue Classification and Surface Extraction. In: Organisation for Human Brain Mapping. Hamburg; 2014.

102. Worsley K, Taylor J, Carbonell F, Chung M, Duerden E, Bernhardt B, et al. SurfStat: A Matlab toolbox for the statistical analysis of univariate andmultivariate surface and volumetric data using linear mixed effects modelsand random field theory. In: Human Brain Mapping. 2009.

103. Lange T, Roth V, Braun ML, Buhmann JM. Stability-Based Validation of Clustering Solutions. Neural Comput. 2004;16(6):1299–323. pmid:15130251

104. Vos de Wael R, Larivière S, Caldairou B, Hong S-J, Margulies DS, Jefferies E, et al. Anatomical and microstructural determinants of hippocampal subfield functional connectome embedding. Proc Natl Acad Sci U S A. 2018 Sep 24;115(40):10154–9. pmid:30249658

105. Von Luxburg U. A tutorial on spectral clustering. Stat Comput. 2007 Dec 22;17(4):395–416.

106. Tenenbaum JB, De Silva V, Langford JC. A global geometric framework for nonlinear dimensionality reduction. Science (80-). 2000 Dec 22;290(5500):2319–23.

107. Glasser MF, Coalson TS, Robinson EC, Hacker CD, Harwell J, Yacoub E, et al. A multi-modal parcellation of human cerebral cortex. Nature. 2016;536(7615):171–8. pmid:27437579

108. Smith SM, Beckmann CF, Andersson J, Auerbach EJ, Bijsterbosch J, Douaud G, et al. Resting-state fMRI in the Human Connectome Project. Neuroimage. 2013;80:144–68. pmid:23702415

109. Sled JGG, Zijdenbos APP, Evans ACC. A nonparametric method for automatic correction of intensity nonuniformity in MRI data. IEEE Trans Med Imaging. 1998 Feb;17(1):87–97. pmid:9617910

110. Fischl B, Sereno MI, Dale AM. Cortical surface-based analysis. II: Inflation, flattening, and a surface-based coordinate system. Neuroimage. 1999/02/05. 1999;9(2):195–207. pmid:9931269

111. Dale AM, Fischl B, Sereno MI. Cortical surface-based analysis. I. Segmentation and surface reconstruction. Neuroimage. 1999;9(2):179–94. pmid:9931268

112. Fischl B, Sereno MI, Tootell RB, Dale AM. High-resolution intersubject averaging and a coordinate system for the cortical surface. Hum Brain Mapp. 2000/01/05. 1999;8(4):272–84. pmid:10619420

113. Robinson EC, Jbabdi S, Glasser MF, Andersson J, Burgess GC, Harms MP, et al. MSM: A new flexible framework for multimodal surface matching. Neuroimage. 2014;100:414–26. pmid:24939340

114. Robinson EC, Garcia K, Glasser MF, Chen Z, Coalson TS, Makropoulos A, et al. Multimodal surface matching with higher-order smoothness constraints. Neuroimage. 2018;167:453–65. pmid:29100940

115. Van Essen DC, Glasser MF, Dierker DL, Harwell J, Coalson T. Parcellations and hemispheric asymmetries of human cerebral cortex analyzed on surface-based atlases. Cereb Cortex. 2012;22(10):2241–62. pmid:22047963

116. Greve DN, Fischl B. Accurate and robust brain image alignment using boundary-based registration. Neuroimage. 2009;48(1):63–72. pmid:19573611

117. Salimi-Khorshidi G, Douaud G, Beckmann CF, Glasser MF, Griffanti L, Smith SM. Automatic denoising of functional MRI data: Combining independent component analysis and hierarchical fusion of classifiers. Neuroimage. 2014;90:449–68. pmid:24389422

118. Murphy K, Fox MD. Towards a consensus regarding global signal regression for resting state functional connectivity MRI. Neuroimage. 2017 Jul 1;154:169–73. pmid:27888059

119. Vos de Wael R, Hyder F, Thompson GJ. Effects of Tissue-Specific Functional Magnetic Resonance Imaging Signal Regression on Resting-State Functional Connectivity. Brain Connect. 2017;7(8):482–90. pmid:28825320

120. Yarkoni T, Poldrack RA, Nichols TE, Van Essen DC, Wager TD. Large-scale automated synthesis of human functional neuroimaging data. Nat Methods. 2011 Aug 26;8(8):665–70. pmid:21706013

121. Marques JP, Kober T, Krueger G, van der Zwaag W, Van de Moortele P-FF, Gruetter R. MP2RAGE, a self bias-field corrected sequence for improved segmentation and T1-mapping at high field. Neuroimage. 2010 Jan 15;49(2):1271–81. pmid:19819338

122. Haast RAM, Ivanov D, Formisano E, Uludaǧ K. Reproducibility and Reliability of Quantitative and Weighted T1 and T2* Mapping for Myelin-Based Cortical Parcellation at 7 Tesla. Front Neuroanat. 2016;10:112. pmid:27917112