Transcriptional profiles of supragranular-enriched genes associate with corticocortical network architecture in the human brain

Corticocortical Connectivity Networks Have Distinct Molecular Profiles of Human Supragranular Enriched Genes.

The primary gene set [Human Supragranular Enriched (HSE) set, n = 19 genes] (Table 1) was selected based on in situ hybridization images available from the Allen Institute Human Transcriptional Atlas, as well as laminar characterizations in Zeng et al. (18). Genes in this set included those whose laminar expression in the cerebral cortex is either (i) unique to humans relative to mouse, and expressed in upper layers of the cortex, or (ii) predominantly enriched in layer V or layers V/VI in mouse, but shift to predominately layer III or layers II/III in humans.

We used previously published fcMRI parcellations to explore how transcriptional profiles vary within and across human corticocortical networks. fcMRI has proven to be a powerful, if indirect, method for estimating network topography in the human brain (26, 27). Connectivity estimates using fcMRI are broadly consistent with known anatomical systems as estimated from nonhuman primate anatomical tracing, with some exceptions (28–33). A parcellation of 17 fcMRI networks based on resting-state data analyzed in Yeo et al. (31) was split into sets of components according to their locations and extent on the cortical surface, resulting in a set of 114 cortical regions composed of roughly symmetric territories in the left and right hemispheres.

Our analyses revealed that distinct molecular profiles of the HSE set define different classes of cerebral cortex. These dissociations are consistent with groupings of large-scale corticocortical networks. Fig. 1A shows the correlation matrix for the HSE genes from microarray data available from six human postmortem brains released by the Allen Institute (human.brain-map.org/). The correlation matrix reflects the similarity of transcriptional profiles of brain regions distributed across the cortex. For each postmortem brain, atlas coordinates of cerebral samples [based on the Montreal Neurological Institute (MNI) atlas] were used to assign each to one of the 114 17-network components (see Methods and Fig. S1 for details on network assignments and naming conventions). Matrix rows and columns correspond to the 114 regions. Fig. 1B shows the resting-state fcMRI correlation matrix for the same 114 components (27) for comparison. The transcriptional profiles in Fig. 1A fall into clear groupings: regions within auditory, somato/motor, and visual networks tend to have similar molecular profiles, which are distinct from regions within paralimbic and association networks. Primary motor and somatosensory samples falling on the precentral and postcentral gyrus, respectively, had highly similar transcriptional expression of the genes examined here, consistent with previous reports (34), and are grouped together for the remainder of the paper. Regions within paralimbic and certain association networks tend to also have similar molecular profiles to each other.

As an illustration of transcriptional similarity across spatially distributed locations, brain samples falling in the superior temporal gyrus at or near the primary auditory cortex, and belt/parabelt auditory association regions have similar transcriptional profiles to nearby somato/motor network regions, as well as to more distant occipital regions (Fig. 1C), relative to neighboring regions in the temporoparietal and middle temporal gyrus that belong to higher-order association networks. As another example of similar transcriptional profiles spanning long distances, regions in the default network across the prefrontal, parietal, temporal, and cingulate cortices have higher similarity to each other than to association regions that fall outside of that network (Fig. S2). These examples highlight that transcriptional similarity occurs not only between regions of the same cortical class (e.g., primary sensory/sensorimotor regions), but also between regions that connect together to form distributed anatomical systems.

Similarity of transcriptional profiles was assessed by testing how many pairs of networks (17 × 16/2 = 136) were significantly positively or negatively correlated. For the HSE set, 103 of 136 network pairs have significant gene expression profile correlations [either positively or negatively correlated; q < 0.05, false-discovery rate (FDR) corrected] (Fig. S3A).

Comparison with Other Gene Sets.

Multiple comparison gene sets were analyzed. Some comparison sets share relevant properties in common with the HSE set: (i) genes that have conserved expression in supragranular layers in both mouse and human (Conserved Supragranular set, n = 14 genes) (18), (ii) genes previously shown to be cortically enriched relative to subcortical structures in humans relative to macaque monkeys and chimpanzees (Human Cortically Enriched set, n = 20 genes) (2), and (iii) genes that were previously associated with anatomic connectivity in rodents and could be matched to homologous genes in the microarray dataset (Rodent Connectivity set, n = 381 genes) (24). The Conserved Supragranular set contains genes that are enriched in layers II and III in both mouse and human (18). Genes in the Human Cortically Enriched set are selectively up-regulated in human cerebral cortex relative to chimpanzee and macaque cerebral cortex (2), and therefore also have human-specific enrichment in the cerebral cortex. The Rodent Connectivity set contains a large number of genes previously shown to predict anatomical connectivity in rodents (24). A fourth comparison set consisted of genes whose expression was shown to be distinctly associated with four human cortical networks (Human/Mouse Connectivity set, n = 136) (25). The expression pattern of this gene set also follows structural connectivity in homologous regions in the mouse (25).

Permutation testing was used to determine whether correlations of transcriptional profiles are stronger for the HSE set than for these alternatives. Seventy-six of the 153 edges in the matrix were significantly stronger for the HSE set relative to the Conserved Supragranular set, 77 of 153 for the HSE vs. Human Cortically Enriched set, and 103 of 153 for the Rodent Connectivity set (q < 0.05, FDR-corrected) (Fig. S3B). For the reverse contrasts, 15 of 153, 16 of 153, and 1 of 153 pairs were significant, respectively (q < 0.05, FDR-corrected). Twenty of 153 network pairs were significantly stronger in the HSE vs. Human/Mouse Connectivity set (q < 0.05, FDR-corrected; 0 of 153 were stronger in the reverse contrast). The HSE and Human/Mouse Connectivity set share four genes in common. Because this overlap comprises 21% of the HSE set, it is not surprising that similar correlation structure exists between these two sets. The Human/Mouse Connectivity set additionally shares two genes in common with the Conserved Supragranular Enriched set and four genes in common with the Rodent Connectivity set.

An open question is whether the within- and across-network transcriptional coexpression observed with the HSE set is present in other sets comprising a wider ontology. A collection of 1,454 gene sets obtained from the Broad Institute’s MSigDB was used to test this possibility.

We used the widely used network-based statistic (NBS) 35) to compare the HSE set against all 1,458 alternatives (Methods). NBS reflects the extent to which connected components (sets of nodes connected by suprathreshold links) are present in a graph. NBS is particularly suited for determining whether the graph structure of two groups (here, two gene sets) differ in a coherent manner, which is not reflected in univariate testing of each link independently. Across a range of tested correlation thresholds (0.1–0.4), the HSE maximal connected component size was significantly larger than that of alternatives (all P < 0.05), with the exception of the Human/Mouse Connectivity set, in which the maximal component size of the HSE set was larger but the difference did not reach statistical significance (P = 0.07–0.12 across thresholds). Overall, the HSE set exhibits a significantly stronger network coexpression pattern than all of the alternative gene sets tested (with P < 0.05 for 99.9%, or 1,457 of 1,458 sets).

Although significantly smaller than HSE, 2 of the 1,454 MSigDB sets had relatively large connected components (Fig. S4C). Eleven other MSigDB sets appeared to have large connected components, but a single gene common to all 13 MSigDB sets, CARTPT (cocaine and amphetamine related transcript, prepropeptide), likely drives this result (Discussion and Fig. S4B).

Gene-Expression Covariance Patterns Are Consistent Across Individuals.

Fig. 2 shows representative individual subject correlation matrices using the HSE set. Brain samples are arranged in the same order as in the group-averaged data in Fig. 1. Each individual has a different matrix size reflecting the number of brain samples available from each network. The resulting matrices consistently show a “blocked” pattern, whereby strong positive correlations exist between the visual and somato/motor networks, as well as to the auditory cortex. Strong negative correlations exist between these sensory/motor regions and paralimbic and certain association networks, particularly the ventral attention and default network regions.

Supragranular Molecular Profiles Are More Similar Within than Across Networks.

Gene-expression profiles are similar both within and between cortical networks for the HSE set. For example, regions within visual networks have similar molecular profiles to each other and also to regions in the somato/motor networks and auditory cortex. Regions within paralimbic networks have similar expression profiles to each other, as well as to regions in default and ventral attention networks. Within-network correlations are consistently higher than out-of-network correlations (Fig. 3).

Association networks have the lowest average difference in within- and out-of-network correlation strengths. This reflects relatively greater heterogeneity in transcriptional profiles in the regions assigned to these networks (Figs. 1A and 3). The 17-network fcMRI parcellation identifies multiple subnetworks within the control network (31) that each has distributed components across prefrontal, parietal, and temporal association cortices (e.g., orange and maroon networks in Fig. 1C). Regions within one subnetwork (Fig. 1C, orange) tend to have more similar transcriptional expression to sensory/motor regions and the dorsal attention network, whereas another control subnetwork (Fig. 1C, maroon) tends to have more similar expression to paralimbic, ventral attention/salience, and default networks (Fig. 1A, see also Fig. 7).

Gene expression is generally expected to be highly consistent across the cerebral cortex relative to other regions of the brain (4, 36, 37). Fig. 4 summarizes mean expression values for nondetrended and detrended (achieved by subtracting the mean expression across all cortical brain samples for each probe) data for the HSE set in one representative individual (H0351.1016). When nondetrended data are used, gene expression is highly similar across networks, consistent with expectations (4, 34). The detrended means for each gene (Fig. 4, Right) illustrate the differences in relative expression of these genes across networks. For example, the relative expression levels for these genes in the paralimbic network are an almost complete reversal of their expression in the visual network. Although the exact ordering of the 19 genes differs across individuals, the effect of detrending was consistent: KCNC3, COL6A1, SYT2, SCN4B, and VAMP1 had the highest relative expression in the visual network across individuals, whereas SCN3B, COL24A1, PRSS12, C1QL2, and TPBG had the lowest relative expression in the visual network. Across individuals this pattern was also consistently reversed for the paralimbic and association networks.

Effect of Spatial Gradients on Molecular Profile Similarity.

Graded expression along a principal axis is common in the adult cerebral cortex and often takes the form of rostrocaudal molecular gradients (8). To assess whether physical distance accounted for transcriptional similarity, the MNI atlas coordinates of brain samples were sorted along the rostrocaudal axis for two representative individuals (Fig. 5). In each case, gene sets showed a “checkerboard” pattern, whereby positive correlations appeared both on and off the diagonal. The extreme off-diagonal positive correlations indicate similarity in gene expression between certain brain samples located far apart along the rostrocaudal axis. Negative correlations close to the diagonal indicate that regions closer together in terms of rostrocaudal location have dissimilar molecular profiles. Transcriptional variation of the HSE genes is not captured by a simple rostrocaudal gradient.

More generally, transcriptional similarity is expected to be negatively correlated with spatial distance (34). Given that corticocortical connectivity is densest at short distances (38), it could be that spatial proximity is a primary driver of transcriptional similarity within and across networks for the HSE genes. Fig. 6 and Figs. S5–S7 plot the correlation strength in gene expression between pairs of brain samples as a function of the Euclidean distance separating them for each individual, for HSE and four alternative gene sets. These results indicate that the influence of physical distance on transcriptional profile correlation is negative, consistent with expectations (r = −0.30 to −0.51 for the HSE set). Beyond the expected negative spatial trend, color-coding the sample pairs reveals the additional influence of network identity (Fig. 6 and Figs. S5–S7). At equivalent spatial distances, cross-network region pairs of visual, somato/motor, or auditory cortical regions have correlations that fall above the trend line, whereas region pairs between the sensory and association or paralimbic cortex tend to cluster below the trend line.

Cross-network correlations were directly compared by binning pairs of samples across five equally spaced Euclidean distance bins. At equivalent distances, transcriptional similarity was consistently higher among visual, somato/motor, and auditory cortical pairs than between those regions and certain association network regions (default, ventral attention) or paralimbic regions (e.g., red versus magenta points in Fig. S5) (two-sample t tests, P’s < 0.001 across individuals, family-wise error rate-corrected). Exceptions were the closest and most extreme distance bins, where there were insufficient numbers of pairs in each category for most individuals. Comparison of sensory-association or sensory-paralimbic pairs to association-paralimbic pairs at equivalent distance bins was also significant (magenta vs. blue points in Fig. S5) (two-sample t tests, P’s < 0.001 across individuals, family-wise error rate-corrected).

This separation of sample pair correlations by cortical type or network assignment was strongest in the HSE set (Fig. 6 and Fig. S5), with diminished separation in other gene sets (Figs. S6 and S7). These results indicate that spatial proximity does not uniformly affect the tendencies of brain samples to have similar or dissimilar transcriptional profiles: at equivalent spatial distances, paralimbic and certain associational networks are more likely to have similar expression profiles to each other than to auditory, visual, or somato/motor profiles.

Degree of Connectivity Predicts Transcriptional Similarity.

Unimodal and hierarchically organized sensory pathways tend to form the densest connections to adjacent areas (16, 38–40). In contrast, heteromodal association areas tend to support longer-range connections and to participate in distributed systems (15, 41, 42).

Neuroimaging approaches have similarly demonstrated that the relative balance of local versus distant connectivity is heterogenous across the cerebral mantle in humans (13). Unimodal sensory cortical regions display preferentially local connectivity, whereas higher order heteromodal association regions distributed across the cortex have preferentially distant connectivity. The analysis in Sepulcre et al. (13) relied on fMRI-based functional connectivity measures with limited spatial resolution and is an indirect measure of anatomical connectivity. Notably, task interactions can shift the balance of local to distant connectivity, indicating dynamic influences on these estimates. However, the broad connectivity principles observed in that study converge with basic observations from nonhuman primate tracing (14, 15).

To test the assumption that relative differences in expression of these HSE genes are associated with differences in corticocortical connectivity, Fig. 7 examines whether relative differences in the degree of local and long-range functional connectivity provide further insight into how the gene profiles group together.

Fig. 7A shows the relative preference for distant versus local connectivity plotted across the human cortical surface from Sepulcre et al. (13) (Supporting Information). The surface map shows the direct subtraction between the degree of positive correlations outside of and within local neighborhoods of voxels across the brain. Warmer colors indicate regions that have preferentially long-distance coupling, whereas cool colors indicate regions with preferentially local coupling. Fig. 7B plots the group-averaged HSE set correlation matrix for 17 networks. Values along the rows indicate the average distant–local score for each network. Networks with similar connectional profiles also tend to have similar gene-expression profiles. Fig. 7C shows the correlation between group-averaged HSE gene profiles for all cortical components (Fig. 1A) and the absolute difference in their distant–local degree connectivity score. The negative correlation (ρ = −0.38; P < 0.001) indicates that transcriptional profile similarity tends to decrease with more divergent distant–local degree connectivity between pairs of brain regions.

Transcriptional Similarity Within and Across Networks.

Transcriptional profiles were more similar within-network than out-of-network boundaries (Figs. 1C and 3 and Fig. S2). In many cases, cross-network transcriptional similarity follows expectations from anatomic connectivity derived from nonhuman primate tract tracing (e.g., Fig. S2). Network groupings, such as those produced from fcMRI estimates, may be best conceived as a means for identifying regions of the cortex with similar connectional properties. The transcriptional similarity between these regions, then, could reflect those shared properties regardless of whether they are directly coupled. For example, note that functional connectivity estimates do not group the default network with the ventral attention/salience network (Fig. 1B), yet their transcriptional profiles for the HSE genes are similar (Figs. 1A and 6). This result is a reminder that transcriptional similarity may primarily be indexing broad classes of cortex and not necessarily markers of individual networks. In addition to their connectivity differences, these regions also tend to show protracted development, have distinct metabolic profiles, and are more variable between individuals (11).

Further subdivisions within the broad classes discussed here (sensory/motor, paralimbic, association) could reveal additional structure of how transcriptional profiles vary across the cerebral cortex. For example, transcriptional profiles of the HSE genes are not uniform across all of the association cortex: in certain cases neighboring regions of the prefrontal cortex have different transcriptional profiles, as do neighboring regions in the parietal and temporal association cortices. The present results show that grouping by network is useful for revealing where these transitions occur, and whether they tend to occur in concert across distributed, interconnected regions.

The profiles of regions in the dorsal attention network take transitional forms between the sensory/motor profiles and association cortex profiles (Fig. 3). The putative homolog of the dorsal attention network in nonhuman primates forms a canonical pathway linking the primary visual cortex to extrastriate visual areas, the middle temporal area complex (MT+), posterior parietal cortex, and the frontal eye fields (43, 44). Recent anatomical tracing in mice suggest that the visual processing pathway may be similarly organized into distributed, dorsal, and ventral processing streams in that species (45). This finding would suggest that it is a conserved pathway that may have been present in the rodent-primate common ancestor. The dorsal attention network is a distributed network involved in sensory-motor integration and may represent an early, evolved prototype network that has been expanded and elaborated in the distinct cortical networks that are observed in primates.

The present results show that regions in the dorsal attention network have transcriptional profiles that align most closely with gene-expression profiles of sensory/motor regions, but are also similar to expression profiles of some association regions. In particular, dorsal attention network regions have similar transcriptional expression of the HSE genes to regions that participate in systems involved in executive control (control network) (Figs. 1A and 6) (46, 47). Consistent with this pattern, human functional connectivity analyses show that the dorsal attention network couples to both extrastriate visual regions, as well as to executive control regions distributed across prefrontal, parietal, and occipito-temporal cortices (Fig. 1) (48).

In other instances, transcriptional agreements between regions belonging to different networks signify an important departure from the pattern observed from fcMRI estimates. The most striking departure is the similarity in expression for the HSE genes between regions belonging to visual networks and those belonging to the somato/motor networks and auditory cortex (Figs. 1, 2, 4, and 7). Here, gene expression is similar despite the lack of strong coupling between these networks indicated by human fcMRI estimates (Fig. 1B), lack of direct strong anatomic connectivity in nonhuman primates (49, 50), and dissimilarity in cytoarchitecture. Transcriptional similarity of supragranular genes in sensory and motor zones may reflect common tendencies of these zones to preferentially support local connectivity (Fig. 7). Thus, a reasonable speculation is that these transcriptional profiles are marking broadly distinct types of cortical architecture, which tend to be shared within a network, but are also shared across networks when the properties of the distinct networks are similar, such as may be the case for sensory and motor zones that are distributed across the cortex. Sensory and motor zones may be similar in the sense that they have preferential local corticocortical connectivity, or in the sense that they lack the distinct alternative properties of association networks. Either form of similarity—shared properties or shared absence of properties—could lead to their transcriptional profiles being similar relative to other cortical zones.

Transcriptional Variation Beyond Spatial Gradients.

Gradients in the expression of key transcription factors are strong drivers of cortical arealization and patterning during development (51–53). In adults, gene expression is fairly uniform across cortical regions (8, 18, 34). However, in adult macaque monkeys, a number of genes with maximal variation in expression across cortical areas are expressed in rostrocaudal gradients. This finding may reflect persistent traces of gradients in developmental origin and timing during cortical neurogenesis (8).

Although spatial proximity is likely to be a primary factor constraining gene expression across the cerebral cortex (54, 55), the current analysis presents a different perspective: relative expression patterns of a small number of genes that are uniquely expressed in human upper cortical layers compared with mouse are correlated in distributed patterns across the cerebral cortex. Their spatial distributions across the cortex are closely linked to cortical subtypes and their associated connectivity profiles. Although the dissociations were strongest for the HSE genes, other gene sets often showed weaker variants of the same phenomenon: in general, sensory and motor regions have more similar gene profiles to themselves and each other, and paralimbic and association regions have more similar expression profiles to themselves and each other.

Functional Ontology of Human Supragranular Enriched Genes.

Some genes within the HSE set have known functional roles in the formation or maintenance of corticocortical connectivity. For example, neurofilament protein-heavy (NEFH) predominantly localizes to corticocortical projecting neurons (56). Laminar density of this protein is higher in V1 and in extrastriate visual areas relative to prefrontal areas, particularly in layer III (57). NEFH is also present in layer III long-range corticocortical projection neurons in the association cortex, but less so in paralimbic regions (56). Consistent with these protein-expression patterns, transcriptional expression of NEFH was consistently higher in the visual network and lower in the paralimbic network (Fig. 4). γ-Synuclein (SNCG) regulates lipid metabolism in the brain (58) and plays a role in synaptic plasticity and formation (59, 60) and neurofilament network integrity (61). Other genes have a less-obvious role in connectivity. µ-Crystallin (CRYM), a cytosolic thyroid hormone-binding protein, is one such unexpected candidate; its precise physiological role in the brain is not yet known. Interpretations of our results must therefore take into account the incomplete record of genetic functional ontologies.

Across individuals, the genes most consistently up-regulated in the visual and somato/motor networks include SYT2, SCN4B, COL6A1, and VAMP1. COL6A1, which encodes an extracellular matrix protein of the collagen family, regulates myelination in the peripheral nervous system (62). It is also one of a set of genes that are up-regulated in humans relative to chimpanzees (6), but the reason for its laminar distribution in the human cerebral cortex is presently unknown. The differential VAMP1 expression in sensory and somato/motor regions, versus paralimbic and association cortices, is corroborated by laminar differences evident in the in situ hybridization patterns: VAMP1 is ubiquitously expressed across layers in primary visual cortex (BA 17), but is restricted and enriched in layer III pyramidal neurons in higher-order areas, including the temporal and dorsolateral prefrontal cortex (human.brain-map.org).

Genes with highest expression in the association and paralimbic cortices relative to sensory and motor regions tended to be consistent across individuals. In particular, COL24A1, SCN3B, PRSS12, and C1QL2 consistently rank lowest in relative gene expression in the visual network, and highest in the paralimbic and default network regions. The function of these genes within the context of supragranular laminar circuitry has not been well characterized in model systems, because by definition they are not prominently expressed in upper layers in mouse. However, the C1QL family is thought to play a critical role in synapse formation and maintenance of synapses between climbing fibers from the inferior olive and Purkinje cells in the cerebellum (63), as well as in postnatal synapse elimination between retinal ganglion cells and the lateral geniculate nucleus (64). It is also expressed in dentate gyrus granule cells (65). The inclusion of several genes implicated in sodium, calcium, and potassium ion channels suggests that the laminar remodeling in the human cerebral cortex compared with the mouse may have been accompanied by changes to the balance of excitability within local and distributed cortical circuits.

Alternative Gene Sets.

Overall, the HSE gene set has a significantly stronger association to corticocortical network connectivity properties than 99.9% of the alternative gene sets (1,457 of 1,458) tested. Beyond direct comparisons to the HSE set, the MSigDB collection presented an opportunity to discover other gene sets with potentially meaningful correlation structure in the context of distributed connectivity networks. Two MSigDB sets, Neurotransmitter Secretion and Regulated Secretory Pathway, retained high correlation structure even when CARTPT—an outlier in other MSigDB sets—was omitted (Fig. S4). The Neurotransmitter Secretion set is a subset of 13 of 15 genes in the Regulated Secretory Pathway set. There is no overlap between these sets and the HSE gene set; we consider these genes interesting candidates for further cross-species cerebral cortical characterization.

Intriguingly, CARTPT appears to be selectively expressed by human upper layer cortical pyramidal neurons, and is expressed in the prefrontal, temporal association, as well as visual cortex. In mouse, Cartpt is also primarily expressed in layer III in the cerebral cortex, and appears particularly enriched in primary sensory areas in that species. Cartpt expression has been best characterized in the hypothalamus. Its function there appears to be related to leptin regulation and shaping the negative response to appetite (66). The function of this gene in the cerebral cortex is not well understood.

In addition to the MSigDB collection, we considered four additional curated gene sets. The enrichment of HSE genes in upper layers likely involves multiple cell types beyond pyramidal projection neurons. That transcriptional coexpression of conserved upper markers (i.e., the Conserved Supergranular set) also, but to a lesser extent, dissociates network types is expected. Transcriptional covariation of genes forming the Human/Mouse Connectivity set was recently shown in an elegant study (25) to dissociate four human connectivity networks on the basis of a metric of higher within- than across-network association. Four of 19 of the HSE genes (21%) were also contained in this set. Several genes in the Rodent Connectivity and Conserved Supragranular sets are also shared with the Human/Mouse Connectivity set, which was generally enriched for genes implicated in ion channel function and in synapses.

The enrichment of these genes in human corresponds with structural connectivity in homologous regions in mouse (25). Thus, the Human/Mouse Connectivity set may preferentially index genes associated with connectivity in cortical regions that are conserved across species. Consistent with this, the four genes (SYT2, NEFH, SV2C, SCN4B) that overlap with the HSE set tend to be up-regulated in the visual and somato/motor networks and down-regulated in paralimbic and association networks. The HSE set captures the main properties of the Human/Mouse Connectivity set, but additionally reveals substantial cross-network transcriptional similarity in humans that reflects network properties across sensory/motor, paralimbic, and association cortical classes.

Microarray Datasets.

Publicly available gene-expression datasets for human postmortem brains (n = 6) were obtained from the Allen Institute for Brain Science. Each individual’s dataset includes normalized microarray expression data for 58,692 probes measured from dissected samples from each postmortem brain. Microarray data were filtered to retain samples from the cerebral cortex (n = 4 left hemisphere only, n = 2 both left and right hemispheres). Before dissection, each postmortem brain was MRI-scanned, and T1-weighted volumes were anatomically segmented using FreeSurfer’s recon-all pipeline (freesurfer.net). Individual subject anatomical volumes were registered to the MNI coordinate space. For more details on microarray and structural imaging datasets, refer to the Allen Human Brain Atlas technical white papers (help.brain-map.org/display/humanbrain/Documentation).

Gene Probe Sets.

Gene sets were obtained from published lists. A primary source consisted of 1,454 annotated gene sets that belong to the Molecular Signatures Database maintained by the Broad Institute. MSigDB gene sets are named by their Gene Ontology term and consist of genes annotated by that term. The sets fall into three broad categories: Biological Process (a recognized series of events or molecular functions), Cellular Component (describing locations at the level of subcellular or macromolecular complexes), and Molecular Function (the functions of a gene product). Four additional alternative gene sets were analyzed: (i) Conserved Supragranular set, n = 14 genes (18), (ii) Human Cortically Enriched set, n = 20 genes (2), (iii) Rodent Connectivity set, n = 381 genes (24), and (iv) Human/Mouse Connectivity set, n = 136 genes (25). Gene information for MSigDB sets can be found online (software.broadinstitute.org/gsea/msigdb/index.jsp). The remaining four comparison sets are available in Dataset S1.

Functional Connectivity Datasets.

Consensus maps of 7 and 17 functional connectivity networks from resting-state data from 1,000 healthy young adults in the Brain Genomics Superstruct Project (71) were used in the present study. Participants provided written informed consent in accordance with guidelines set by institutional review boards of Harvard University or Partners Healthcare. The Yeo et al. (31) maps are publicly available for download (https://surfer.nmr.mgh.harvard.edu/fswiki/CorticalParcellation_Yeo2011), and are highly similar to alternative estimates (72). The method for clustering has been described elsewhere (31). Briefly, for each participant, the Pearson's product moment correlation was computed between each surface vertex (n = 18,715) and 1,175 uniformly distributed cortical regions of interest (ROIs). The “connectivity profile” of each surface vertex was its functional coupling to these ROIs. Each participant's 18,715 × 1,175 matrix of correlations was binarized to retain the top 10% of correlations before summing across subjects to obtain an overall group estimate P. Therefore, the ith row and jth column of the matrix P was the number of subjects whose correlations between the ith vertex and jth ROI are within the top 10% of correlations (within each individual subject). In other words, each matrix component took on integer values from 0 to 1,000. The connectivity profiles were clustered using a mixture of von Mises–Fisher distributions. For more details, see refs. 31 and 73. Because our previous analyses (31) identified solutions with 7 and 17 network clusters to be particularly stable, these were adopted for the present study.

Network Assignments.

A volumetric 17-network parcellation in MNI space from Yeo et al. (31) was split into sets of components according to their locations and extent on the cortical surface. Specifically, a set of 114 cortical regions composed of roughly symmetric territories in the left and right hemispheres were defined by selecting vertices on a spherical inflation of each consensus map with respect to network boundaries, sulcal patterns, and confidence maps. “Confidence maps” refer to estimates of instability of network assignments and suggest either uncertainty about proper network assignment in a given region, or anticipate further subdivisions of networks (31). A similar approach to obtain sets of subregions within cortical networks has been described elsewhere (74–76).

The cerebral cortical brain samples were assigned to network components by associating the MNI coordinate of each sample with the corresponding network label at that coordinate. Brain samples that fell on the borders between network components—or whose neighboring voxels did not have the same network assignment—were excluded. Most of the excluded regions fell in subcortical structures (thalamus, basal ganglia) that were given cortical slab labels in the Allen Institute ontology (see Fig. S1 for locations of included and excluded samples). Filtering subcortical structures was necessary as differences between cortical and subcortical transcriptional expression dominate those within the cerebral cortex (34).

For visualization, the components were divided into seven groups (default, control, paralimbic, salience/ventral attention, dorsal attention, somato/motor, visual), which broadly correspond to major networks discussed in the literature. The paralimbic network includes paralimbic cortical regions, such as the medial surface of the temporal lobe, orbitofrontal cortex, and subgenual cingulate cortex. Note that the somato/motor network clusters both pre- and postcentral gyrus, and also likely includes the loci of primary auditory core and belt areas in the temporal lobe (31). Although strong anatomic connections between motor and somatosensory regions are expected (77), we could not resolve from the fcMRI estimates whether the auditory cortex is truly functionally coupled to these zones using fcMRI approaches (31). However, the transcriptional profile data are able to make this distinction, so transcription data are separated by whether they likely were obtained from motor (precentral gyrus) or primary somatosensory areas (postcentral gyrus). Similarly, the close proximity of somatosensory regions, particularly S2, to the auditory cortex and spatial limitations inherent to fMRI-based approaches induces clustering of auditory regions to somato/motor regions because of partial voluming of the fcMRI data (see figure 14 of ref. 31). Transcriptional profiles within auditory regions, which do not suffer from this spatial blurring, are assessed separately from the somato/motor network in the present results.

Molecular Profile Correlation Analysis.

When multiple probes in the microarray data existed for a given gene, expression values were averaged across probes. For each gene set, and for each postmortem brain, pairwise correlations of expression values of the probes were computed for all possible pairs of network components.

Because gene profiles are in general highly similar across the cerebral cortex (4), expression values were detrended by subtracting the mean expression value across samples before correlation analyses. Because the presence of outliers in a gene set can potentially inflate or deflate the correlation coefficient, rank-signed Spearman’s ρ is used to generate correlation coefficients for each gene set (see Fig. S4A for a comparison between Pearson’s product-moment correlation and Spearman’s ρ). Rows and columns of the resulting cross-correlation matrices were arranged so that components that belong to the same functional connectivity network are grouped together. For brains in which data from both hemispheres were available (H0351.2001 and H0351.2002), left hemisphere and right hemisphere components were interleaved within each network. For visualization, correlations were thresholded at ρ = ± 0.1.

To obtain a group summary of molecular profile agreements within and across networks, the correlations for all regions belonging to the same network (excluding self-correlations at the diagonal) were averaged. These were averaged across individuals. Sixteen of the 114 regions did not contain brain samples from any of the individuals and were excluded from analysis, leaving 98 regions.

To determine whether correlations of gene expression between pairs of the 17 networks are significantly different from zero, average expression of each gene was computed for each network by averaging across all brain samples that fell within a given network of the 17-network parcellation (Fig. S3A). The labels of the HSE genes corresponding to the averaged expression values were permuted independently for each network (retaining the same permutation order across subjects for a given permutation). The resulting correlation matrix was computed to construct the null distribution (10,000 permutations). This procedure breaks the correlation between samples while retaining all other data structure. Correlations surviving an FDR-corrected value of q < 0.05 are reported. See Supporting Information for additional details.

To compare gene sets, a null distribution was obtained by randomly exchanging the 19 HSE genes with those of the alternative set and computing the statistic of interest using the permuted gene sets (Fig. S3B). For example, for the Rodent Connectivity set, the 19 HSE genes were pooled with the 381 Rodent Connectivity genes. The pool was then randomly and repeatedly (10,000 permutations) split into two sets of 19 and 381 genes, and the absolute difference between the correlation matrices of the two gene sets was computed. The resulting “difference matrix” was compared with the difference matrix computed from the two original gene sets. Only correlations surviving an FDR-corrected value of q < 0.05 are interpreted. See Supporting Information for additional details.

In addition to comparing individual elements in the correlation matrices of two gene sets, we also used a global metric called network-based statistic (35) to examine whether there is significant correlation structure in a connectivity matrix of interest. Briefly, NBS is a nonparametric statistical test used to isolate components (sets of vertices connected by suprathreshold edges) of an n × n (here 98 × 98, where each element is a network component from the 17 network parcellation) undirected connectivity matrix that differs between two populations (in the present case, between two gene sets). NBS is a procedure to control family-wise error rate in the weak sense. Specifically, the difference matrix between the HSE and a given alternative gene set is thresholded, and the largest connected component is computed. This process is repeated for the permuted difference matrices created by the pooling and random splitting procedure described above, and component sizes between the original difference matrix and the permuted difference matrix are compared.