Research Article

Spatiotemporal transcriptomic divergence across human and macaque brain development

See allHide authors and affiliations

Science  14 Dec 2018:
Vol. 362, Issue 6420, eaat8077
DOI: 10.1126/science.aat8077

Structured Abstract

INTRODUCTION

Improved understanding of how the developing human nervous system differs from that of closely related nonhuman primates is fundamental for teasing out human-specific aspects of behavior, cognition, and disorders.

RATIONALE

The shared and unique functional properties of the human nervous system are rooted in the complex transcriptional programs governing the development of distinct cell types, neural circuits, and regions. However, the precise molecular mechanisms underlying shared and unique features of the developing human nervous system have been only minimally characterized.

RESULTS

We generated complementary tissue-level and single-cell transcriptomic datasets from up to 16 brain regions covering prenatal and postnatal development in humans and rhesus macaques (Macaca mulatta), a closely related species and the most commonly studied nonhuman primate. We created and applied TranscriptomeAge and TempShift algorithms to age-match developing specimens between the species and to more rigorously identify temporal differences in gene expression within and across the species. By analyzing regional and temporal patterns of gene expression in both the developing human and macaque brain, and comparing these patterns to a complementary dataset that included transcriptomic information from the adult chimpanzee, we identified shared and divergent transcriptomic features of human brain development. Furthermore, integration with single-cell and single-nucleus transcriptomic data covering prenatal and adult periods of both species revealed that the developmental divergence between humans and macaques can be traced to distinct cell types enriched in different developmental times and brain regions, including the prefrontal cortex, a region of the brain associated with distinctly human aspects of cognition and behavior.

We found two phases of prominent species differences: embryonic to late midfetal development and adolescence/young adulthood. This evolutionary cup-shaped or hourglass-like pattern, with high divergence in prenatal development and adolescence/young adulthood and lower divergence in early postnatal development, resembles the developmental cup-shaped pattern described in the accompanying study by Li et al. Even though the developmental (ontogenetic) and evolutionary (phylogenetic) patterns have similar profiles, the overlap of genes driving these two patterns is not substantial, indicating the existence of different molecular mechanisms and constraints for regional specification and species divergence.

Notably, we also identified numerous genes and gene coexpression modules exhibiting human-distinct patterns in either temporal (heterochronic) or spatial (heterotopic) gene expression, as well as genes with human-distinct developmental expression, linked to autism spectrum disorder, schizophrenia, and other neurological or psychiatric diseases. This finding potentially suggests mechanistic underpinnings of these disorders.

CONCLUSION

Our study provides insights into the evolution of gene expression in the developing human brain and may shed some light on potentially human-specific underpinnings of certain neuropsychiatric disorders.

Concerted ontogenetic and phylogenetic transcriptomic divergence in human and macaque brain.

Left: Human and macaque brain regions spanning both prenatal and postnatal development were age-matched using TranscriptomeAge. Right: Phylogenetic transcriptomic divergence between humans and macaques resembles the developmental (ontogenetic) cup-shaped pattern of each species, with high divergence in prenatal development and adolescence/young adulthood and lower divergence during the early postnatal period (from perinatal to adolescence). Single-cell transcriptomics revealed shared and divergent transcriptomic features of distinct cell types.

Abstract

Human nervous system development is an intricate and protracted process that requires precise spatiotemporal transcriptional regulation. We generated tissue-level and single-cell transcriptomic data from up to 16 brain regions covering prenatal and postnatal rhesus macaque development. Integrative analysis with complementary human data revealed that global intraspecies (ontogenetic) and interspecies (phylogenetic) regional transcriptomic differences exhibit concerted cup-shaped patterns, with a late fetal-to-infancy (perinatal) convergence. Prenatal neocortical transcriptomic patterns revealed transient topographic gradients, whereas postnatal patterns largely reflected functional hierarchy. Genes exhibiting heterotopic and heterochronic divergence included those transiently enriched in the prenatal prefrontal cortex or linked to autism spectrum disorder and schizophrenia. Our findings shed light on transcriptomic programs underlying the evolution of human brain development and the pathogenesis of neuropsychiatric disorders.

The development of the human nervous system is an intricate process that unfolds over a prolonged time course, ranging from years to decades, depending on the region (16). Precise spatial and temporal regulation of gene expression is crucial for all aspects of human nervous system development, evolution, and function (613). Consequently, alterations in this process have been linked to psychiatric and neurological disorders, some of which may exhibit primate- or human-specific manifestations (11, 1418). However, our ability to explain many aspects of human nervous system development and disorders at a mechanistic level has been limited by our evolutionary distance from genetically tractable model organisms, such as the mouse (15, 16, 1922), and by a lack of contextual and functional interpretations of polymorphisms and disease-associated variations in the human and nonhuman primate (NHP) genomes (11, 17, 21, 23). Moreover, neither the extent of molecular changes underlying human-specific differences nor the specific developmental programs affected by these changes have been thoroughly studied.

The rhesus macaque (Macaca mulatta) is the most widely studied NHP in neuroscience and medicine (2426). The macaque nervous system parallels the human nervous system with its complex cellular architecture and extended development, and thereby offers a unique opportunity to study features of neurodevelopment that are shared and divergent between the two closely related primates. Furthermore, studies of post mortem NHP tissues provide a unique opportunity to validate results obtained using post mortem human tissue, especially those from critical developmental periods that can be confounded by ante mortem and post mortem factors and tissue quality. Finally, substantial advances in transgenic and genome-editing technologies now allow the possibility of creating more precise genetic models for human disorders in macaques (2426). This will facilitate the interrogation of the effects of specific gene mutations in a model that is closer to the human brain than any other experimental animal.

Comparative transcriptomic profiling offers unbiased insight into conserved and clade- or species-specific molecular programs underlying cellular and functional development of the human nervous system (2731). However, a systematic characterization of the spatial and temporal transcriptomic landscapes of the macaque brain at the region-specific and single-cell levels, as well as the identification of shared and divergent features between humans and macaques, are lacking. Data and analyses such as we present here should provide both retrospective and prospective benefits to the fields of neuroscience, evolutionary biology, genomics, and medicine.

Study design, data generation, and integrated analysis

RNA sequencing (RNA-seq) data were obtained from bulk tissue (366 samples from 26 prenatal and postnatal brains) or single cells/nuclei (113,274 cells or nuclei from two fetal and three adult brains) from post mortem rhesus macaque specimens. Both tissue and single cell/nucleus datasets were subjected to multiple quality control measures (figs. S1 to S6 and tables S1 and S2) (32). Tissue-level samples covered the entire span of both prenatal and postnatal neurodevelopment (Fig. 1, A and B, and table S1) and included 11 areas of the cerebral neocortex (NCX), hippocampus (HIP), amygdala (AMY), striatum (STR), mediodorsal nucleus of thalamus (MD), and cerebellar cortex (CBC). Subject ages ranged from 60 post-conception days (PCD) to 11 postnatal years (PY) and were matched by age and brain region to 36 human brains from an accompanying study (33) and five adult chimpanzee brains from a previous study (34) (Fig. 1A). To investigate the contribution of different factors to the global transcriptome dynamics, we applied unsupervised clustering and principal components analysis, which revealed that age, species, and regions contributed more to the global transcriptomic differences than did other tested variables (figs. S3 and S4).

Fig. 1 Conserved and divergent transcriptomic features of human and macaque neurodevelopmental processes.

(A) Plot depicting the real age (x axis) and the age predicted by TranscriptomeAge (y axis) of human, chimpanzee, and macaque. Macaque (164 PCD) and human (266 PCD) births are shown as green and red dashed lines, respectively. (B) Schematic showing human developmental periods as described in Kang et al. (29) and the matched macaque developmental and chimpanzee adult datasets. Each line corresponds to one macaque or one chimpanzee specimen and the corresponding predicted age when compared to human neurodevelopment. PCD, post-conception day; PY, postnatal year. The asterisk indicates the extension of the early fetal period, in which early fetal macaques (60 PCD) cluster with midfetal humans. (C) The weight (W) of five transcriptomic signatures in the developing human (solid line) and macaque (dashed line) NCX and the respective association with neurodevelopmental processes. In signature 1 (neurogenesis), the arrow indicates the point at which the signature reaches the minimum in humans (red) and macaques (green). The asterisk indicates the same as in (B). In transcriptomic signatures 2, 3, 4, and 5, arrows indicate the point at which the signatures reach the maximum in humans (red) and macaques (green). Note that for transcriptomic signatures 2 and 3 (neuronal differentiation and astrogliogenesis), there is a synchrony between humans and macaques, whereas for transcriptomic signatures 4 and 5 (synaptogenesis and myelination), there is heterochrony between the species, with acceleration in human synaptogenesis and delay in human myelination. Prefrontal cortical areas are plotted in red, primary motor cortex in orange, parietal areas in green, temporal areas in blue, and primary visual cortex in gray. MFC, medial prefrontal cortex; OFC, orbital prefrontal cortex; DFC, dorsolateral prefrontal cortex; VFC, ventrolateral prefrontal cortex; M1C, primary motor cortex; S1C, primary somatosensory cortex; IPC, inferior posterior parietal cortex; A1C, primary auditory cortex; STC, superior temporal cortex; ITC, inferior temporal cortex; V1C, primary visual cortex. (D) Cell type enrichment is shown for each signature. P values adjusted by Benjamini-Hochberg procedure are plotted (with ranges indicated by size of dots); significance is labeled by color (red, true; gray, false). H, human; M, macaque; eNEP/RGC, embryonic neuroepithelial progenitor/radial glial cell; eIPC, embryonic intermediate progenitor cell; eNasN, embryonic nascent neuron; ExN, excitatory neuron; InN, interneuron; Astro, astrocyte; OPC, oligodendrocyte progenitor cell; Oligo, oligodendrocyte; Endo, endothelial cell; VSMC, vascular smooth muscle cell.

To explore cell type origins of tissue-level interspecies differences, we conducted single-cell RNA-seq (scRNA-seq) on 86,341 cells from six matching regions of two 110-PCD fetal macaque brains [i.e., the dorsolateral prefrontal neocortex (DFC, also called DLPFC), HIP, AMY, STR, MD, and CBC] and single-nucleus RNA-seq (snRNA-seq) of 26,933 nuclei from three adult macaque DFCs (8, 11, and 11 PY; tables S2 and S3) (32). These data were complemented by 17,093 snRNA-seq samples from adult humans [see (33)] as well as two scRNA-seq datasets from embryonic and fetal human NCX (33, 35). In the six fetal macaque brain regions, we identified 129 transcriptomically distinct clusters of cell types (i.e., 19 in DFC, 20 in HIP, 25 in AMY, 22 in STR, 20 in MD, and 23 in CBC) (figs. S7 to S12 and tables S3 and S4). In the adult human DFC (fig. S13) and adult macaque DFC (fig. S14), we identified 29 and 21 transcriptomically distinct cell types, respectively (tables S3, S5, and S6). Alignment of our macaque fetal data with the adult single-nucleus data revealed hierarchical relationships and similarities between major cell classes, reflecting their ontogenetic origins and functional properties (fig. S15). Cell clusters were categorized by their gene expression patterns and assigned identities commensurate with their predicted cell type and, in the case of human adult neocortical excitatory neurons, their putative laminar identity. Although the majority of cell clusters were composed of cells derived from all brains, we found a few clusters in subcortical regions (AMY, 2 of 25 clusters; CBC, 1 of 23 clusters; STR, 1 of 22 clusters) that included cells from a single donor brain. This might be due to variations in dissection, age (even though both fetal macaques were 110 PCD, a 3- to 4-day variation remains), individual differences, and other technical bias. We used the single-cell datasets in this and the accompanying study (33) to deconvolve tissue-level RNA-seq data, identify temporal changes in cell type–specific signatures, analyze differences in cell types and their transcriptomic profiles, and conduct cell type enrichment analyses.

Similarities and differences in the spatiotemporal dynamics of the human and macaque brain transcriptomes

Unsupervised hierarchical clustering and principal components analysis of bulk tissue revealed common principles of transcriptomic regional architecture across development in macaques and humans (figs. S3 and S4). Among macaque regions, these analyses showed distinct and developmentally regulated clustering of the NCX (combination of 11 areas), HIP, and AMY, with CBC exhibiting the most distinctive transcriptional profile—an observation shared with our complementary study in humans (27, 29, 30, 33, 36). A hierarchical clustering of both fetal and postnatal NCX areal samples revealed their grouping by topographical proximity and functional overlap, similar to those relationships that we observed in the human brain (fig. S3). Thus, these results show that the transcriptomic architecture of the macaque brain is regionally and temporally specified and reflects conserved global patterns of ontogenetic and functional differences that are also found in humans.

To explore species similarities and differences in the spatiotemporal dynamics of the brain transcriptome, we used the XSAnno computational framework (37) to minimize biases in comparative data analyses arising from the disparate quality of gene annotation for the two species. We created common annotation sets of 27,932 and 26,514 orthologous protein-coding and noncoding mRNA genes for human-macaque and human-chimpanzee-macaque comparisons, respectively (fig. S2) (32). Next, we developed TranscriptomeAge, an algorithm to unbiasedly predict the equivalent ages of human and macaque samples on the basis of temporal transcriptomic changes (32). We chose to optimize this model for age-matching the aforementioned 11 neocortical areas, which are highly similar in terms of their transcriptomes, cellular composition, and developmental trajectories when compared to other brain regions [see (33)]. TranscriptomeAge confirmed transcriptomic similarities in both species coinciding with major prenatal and postnatal developmental phases, including fetal development, infancy, childhood, and adulthood (Fig. 1, A and B, and figs. S16 to S18). However, we identified two human developmental periods where alignment suggested that they are transcriptomically distinct from macaques and/or are especially protracted. First, 60-PCD macaque specimens [which correspond to the human early fetal period (29) according to the Translating Time model (38)] were most closely aligned with midfetal human samples (102 to 115 PCD, i.e., 14.5 to 16.5 post-conception weeks). This suggests that, transcriptomically, human brain development is protracted even at early fetal periods. Second, we found that 2-, 3.5-, 4-, 5-, and 7-PY macaque specimens, of which at least the youngest should chronologically match to human childhood (39), did not align with any of our human specimens from early or late childhood [1 to 12 PY, or periods 9 and 10 according to (29)] but did align with adolescent and adult humans (Fig. 1, A and B). Consistent with previous morphophysiological and behavioral studies (5), these results indicate that macaques lack global transcriptomic signatures of late childhood and/or that humans have a prolonged childhood relative to macaques (Fig. 1, A and B).

Species differences in the timing of concerted neurodevelopmental processes

We hypothesized that the observed developmental differences between humans and macaques might be grounded on transcriptomic changes in concerted biological processes in developmental timing (i.e., heterochrony). By decomposing the gene expression matrix of human neocortical samples, we identified five transcriptomic signatures underlying neocortical development (32). Using top cell type–specific genes derived from our prenatal single-cell and adult single-nucleus data, we analyzed cell type enrichment of each of the five signatures, and ascribed them to neurogenesis, neuronal differentiation, astrogliogenesis, synaptogenesis, and oligodendrocyte differentiation and myelination (Fig. 1, C and D, and fig. S19). To determine whether the transcriptomic signatures we identified were correctly assigned, we compared their developmental patterns to the timing of major human neurodevelopmental processes, expression trajectories of key genes previously implicated in those processes, and trajectories of cell type proportions identified by the deconvolution of tissue-level data (figs. S19 and S20). We found that the developmental trajectories of genes associated with neuronal differentiation, synaptogenesis, and myelination, as well as the cell type proportions of fetal human or macaque excitatory neurons, astrocytes, and oligodendrocytes, matched those of the corresponding transcriptomic signatures (fig. S20). Moreover, the identities we assigned to these transcriptomic signatures were confirmed by comparison of transcriptomic signatures to independently generated nontranscriptomic data predicting the start and end of human neocortical neurogenesis (for neurogenesis) (40) and to data measuring the number of doublecortin (DCX)–immunopositive nascent neurons in the human hippocampus throughout development and adulthood (for neuronal migration and initial differentiation) (41), developmental variation in synaptic density in the human cortex (for synaptogenesis) (42), and myelinated fiber length density (for myelination) (43) (fig. S19).

Next, we analyzed how the shape of the five transcriptomic trajectories was conserved across the 11 neocortical areas within each species and between species. Analysis of their trajectories within each species revealed that the shape of a given trajectory is similar across neocortical areas (Fig. 1C and fig. S17). However, the transcriptomic trajectories associated with oligodendrocyte differentiation and myelination exhibited a prominent temporal shift (asynchrony) across neocortical areas in both species (fig. S17). Between species, myelination and, to a lesser extent, synaptogenesis exhibited species differences in the shapes of these trajectories; the myelination transcriptomic signature progressively increased in the human NCX beginning from late fetal development through adulthood without reaching an obvious plateau until 40 PY, but in the macaque NCX the myelination signature reached a plateau around the first postnatal year (Fig. 1C). This corresponds to early childhood in human neurodevelopment [window 6 or period 10 according to (33) and (29), respectively] and is consistent with histological studies and reflective of previously reported hierarchical maturation of neocortical areas (4347). Similarly, we corroborated synchronous or concurrent transcriptomic patterns of neocortical synaptogenesis by analyzing previously collected data on synaptic density in multiple areas of the macaque NCX (48) (fig. S19). However, we observed that the synaptogenesis transcriptomic trajectory peaked earlier in humans than in macaques, at the transition between late infancy and early childhood (Fig. 1C). In addition, expression trajectories of genes induced by neuronal activity—a process critical for synaptogenesis—also showed drastic increases during late fetal development and infancy, and, like the synaptogenesis trajectory, displayed a concurrent or synchronous shape across neocortical areas [see (33)]. Interestingly, the developmental transcriptomic profile of DCX (a marker of nascent, migrating neurons) showed that macaques maintain higher expression in the hippocampus throughout postnatal development and adulthood; this suggests that postnatal neurogenesis is more prominent in the macaque hippocampus than in the human hippocampus, as recently shown (fig. S19) (49). Thus, both species exhibited distinct transcriptomic signatures of neoteny, such as prolonged myelination in humans and prolonged postnatal hippocampal neurogenesis in macaques. Together, these data suggest that the temporal staging of major neurodevelopmental processes, in particular with myelination beginning in primary areas before association neocortical areas, is a conserved feature of primate development, although the temporal progression of certain processes is heterochronic.

Concerted ontogenetic and phylogenetic transcriptomic divergence

After matching the global transcriptome by age between the two species, we analyzed regional differences in gene expression (heterotopy) within each species. By adopting Gaussian-process models to accommodate the spatiotemporal correlations of gene expression (32), we found that the developmental cup-shaped or hourglass-like pattern of transcriptomic interregional differences we observed in humans (33) is also present in macaque neocortices and other brain regions (Fig. 2, A and B, and fig. S21), with greater differential expression between regions observed during early and midfetal ages preceding this period and subsequent young adulthood. Notably, two brain regions—CBC and STR—exhibited greater differences, relative to other brain regions, beginning immediately after birth, rather than beginning during childhood or adolescence (fig. S21). This suggests that the development of the primate forebrain may be constrained by unique developmental or evolutionary influences, which led us to investigate the gene expression patterns, developmental processes, and cell types underlying this transcriptomic phenomenon.

Fig. 2 Ontogenetic interregional transcriptomic differences display a cup-shaped pattern in humans and macaques.

(A and B) The interregional difference was measured as the average distance of each neocortical area to all other areas in the human (A) and macaque (B) neocortices across development. The upper-quartile interregional difference among all genes is plotted; the color scale indicates magnitude. The gray planes represent the transition from prenatal to early postnatal development (late fetal transition) and from adolescence to adulthood. (C) The number of coexpression modules that display gradient-like expression (anterior to posterior, posterior to anterior, medial to lateral, temporal lobe–enriched) and enrichment in primary areas or enrichment in association areas in each developmental phase. Left, human modules; right, macaque modules. (D) Donut plots depicting the modules from (C) that exhibited species-distinct interregional differences. The expression pattern of each species-distinct module is shown for humans (top) and macaques (bottom). Color scales indicate expression level of the genes in each module. Prenatal modules show a human-distinct anterior-to-posterior expression gradient (left); macaque-distinct early postnatal modules show enrichment in primary or association areas (center); and a macaque-distinct adult module is enriched in association areas, especially in MFC (right). HS, human (Homo sapiens) module; MM, macaque (Macaca mulatta) module.

To do so, we considered three phases of brain development mirroring major transitions in the cup-shaped pattern: prenatal development, early postnatal development, and adulthood. Between these three phases are two transitional periods: a steep late fetal transition (33) and a more moderate transition between childhood/adolescence and adulthood. We performed weighted gene coexpression network analysis (WGCNA) independently for each phase and species, resulting in Homo sapiens (HS) and macaque (Macaca mulatta, MM) modules (32) (table S7), with analyses conducted on 11 neocortical areas; this allowed us to identify discrete spatiotemporal expression patterns that otherwise might be co-mingled as a result of the highly disparate nature of CBC and other non-neocortical regions. Within the prenatal phase, we found 12 modules consisting of genes exhibiting spatial expression gradients along the anterior-posterior (8 modules) and medial-lateral (1 module) axes of the NCX and broadly reflecting prospective neocortical areal topography (Fig. 2C). For example, prenatal modules HS85 and HS87 exhibited prefrontal/frontal-enriched graded expression in the human brain, tapering to lowest expression in the temporal and occipital lobes (Fig. 2D). Furthermore, prenatal modules, such as HS15 and MM57, had their highest expression restricted to the temporal lobe (table S8 and figs. S22 and S23) during prenatal development.

In contrast to the prenatal phase, modules identified from early postnatal development (i.e., infancy, childhood, and adolescence) in either species did not exhibit anterior-to-posterior or medial-to-lateral expression gradients. Rather, the greater regional synchrony characterizing gene expression in this phase yielded differences organized not around topography but between primary and association areas of the NCX (Fig. 2C, figs. S24 and S25, and table S9). This suggests that the gradient-like transcriptomic patterns arising during prenatal development are superseded by myelination and neuronal activity–related processes postnatally, which may differentiate the separation between primary and association areas. Early postnatal modules such as MM42, MM24, and MM23, among others, exhibited greater expression in primary areas such as the primary motor cortex (M1C), primary auditory cortex (A1C), and primary visual cortex (V1C) than in association areas such as DFC and ventrolateral prefrontal cortex (VFC) (Fig. 2D).

The transition to young adulthood was marked by another decrease in interregional differences, but this reduction was not as pronounced as in the late fetal transition, nor were interregional patterns of gene expression markedly different in the adult. Thus, gene expression differences between primary and association areas continued to drive regional variation in both adult humans and macaques (Fig. 2, C and D, figs. S26 and S27, and table S10). Gene Ontology (GO) enrichment analysis using the top variant genes in each period, with all genes expressed in each period as background, indicated differential enrichment of biological processes associated with different cell populations across areas and time. As observed in the accompanying human study (33) and commensurate with the developmental trajectories of the observed transcriptomic signatures, the functional terms enriched prenatally were generally related to neurogenesis and neuronal differentiation, whereas early postnatal and adult functional terms were enriched for processes related to synaptogenesis and myelination (fig. S28).

We next sought to determine whether the regional-specific expression patterns of coexpression modules detected in human brains correlated with their expression patterns in macaque brains, and vice versa (32). We found that two human prenatal modules contained genes exhibiting a pronounced anterior-to-posterior gradient in the human NCX, HS85 and HS87, but these genes did not exhibit enriched expression in the macaque prefrontal cortex (Fig. 2D and table S8). Among genes in these modules were RGMA and SLIT3, two genes encoding axon guidance molecules (50), and BRINP2 and CXXC5, which encode proteins involved in retinoic acid signaling (51), potentially implicating this signaling pathway—critical for early brain development and neuronal differentiation (51)—in the patterning of the human prefrontal cortex. We also observed that several modules in macaque postnatal development that did not correlate well with human modules (MM23, MM24, MM26, and MM42) were enriched for genes that are expressed in oligodendrocytes (Fig. 2D, fig. S24, and table S9) and were up-regulated in all primary areas of macaque NCX relative to association areas. Conversely, genes in these modules were up-regulated in humans only in M1C and A1C, but not in primary somatosensory cortex (S1C) or V1C (fig. S24 and table S9). Integration with our multi-regional database of the adult chimpanzee transcriptome (34) indicates that the macaque gene expression pattern, rather than the human gene expression pattern, may be unique among these species (fig. S29). Many of the species-specific patterns of diversification between primary and association areas that we observed during early postnatal development were preserved in adult modules of both species (fig. S26), with some notable exceptions. For example, the adult macaque module MM25 exhibited up-regulation in association areas in both species, but prominent up-regulation in the medial prefrontal cortex (MFC) and down-regulation in V1C were observed only in macaques (Fig. 2D, fig. S26, and table S10).

These findings reaffirm a conserved framework in primate neocortical development and function (21), including a topographic basis for transcriptomic differences during prenatal development and functional relationships postnatally. Our analyses also suggest that interregional and interspecies differences in oligodendrocyte development and myelination, particularly during early postnatal development, mediate key aspects of transcriptomic variation both within and among species.

Heterotopic changes in human and macaque brain transcriptomes

We next investigated the transcriptomic divergence between humans and macaques for each brain region across development. We found that the developmental phases exhibiting high levels of interregional differences within each species (i.e., prenatal development and young adulthood) also displayed greater transcriptomic divergence between the two species, revealing a concerted phylogenetic (evolutionary) cup-shaped pattern (Fig. 3A). This phylogenetic cup-shaped pattern divided neurodevelopment into the same three phases as the regional ontogenetic (developmental) cup shape (Fig. 3A). However, unlike the ontogenetic (developmental) cup-shaped pattern, where CBC, MD, and STR disproportionally exhibited more intraspecies differences than NCX, HIP, and AMY, all regions appeared to exhibit a relatively similar amount of interspecies differences (Fig. 3A). Interestingly, interspecies differences among neocortical areas were distinct enough to provide clear clustering of topographically and functionally related prefrontal areas [i.e., MFC, orbital prefrontal cortex (OFC), DFC, and VFC], particularly during prenatal development, or topographically distributed nonvisual primary areas (i.e., M1C, S1C, and A1C) in adulthood. Prospective areas of the prefrontal cortex, which underlie some of the most distinctly human aspects of cognition, were more phylogenetically distinct than other neocortical areas during early prenatal development (Fig. 3A and fig. S30). Together, these findings suggest that the evolutionary and developmental constraints acting on the brain transcriptome, in particular the NCX, may share some overlapping features.

Fig. 3 Transcriptomic divergence between humans and macaques throughout neurodevelopment reveals a phylogenetic cup-shaped pattern.

(A) Interspecies divergence, measured as absolute difference in gene expression, between humans and macaques in each brain region throughout development (coded as in Fig. 2A). The upper-quartile divergence among all genes is plotted. The gray planes represent the transition from prenatal to early postnatal development (late fetal transition, left) and from adolescence to adulthood (right). (B) Venn diagrams displaying the number of differentially expressed genes (DEX, top) or genes with differential exon usage (DEU, bottom) between humans and macaques in at least one brain region during prenatal development, early postnatal development, and adulthood. (C) Bubble matrix with examples of genes showing global or regional interspecies differential expression. Brain regions displaying significant differential expression between humans and macaques are shown with black circumference. Red circles show up-regulation in humans; blue, up-regulation in macaques. Circle size indicates absolute log2 fold change. (D) Percentage of overlap between genes showing the highest interspecies divergence in each region (driving the evolutionary cup-shaped pattern) and genes with the largest pairwise distance between brain regions in prenatal, early postnatal, and adult human and macaque brains (driving the developmental cup-shaped pattern). The result is plotted using a variable number of the highest-ranked genes based on interregional difference and interspecies divergence. Data are means ± SD across regions.

To gain insight into the transcriptomic programs driving phylogenetic divergence across neocortical areas, we conducted a functional annotation of the top 100 genes driving the observed variation along the first principal component (PC1). We found that interspecies divergence in the prenatal prefrontal cortex could be explained by an enrichment of genes related to cell proliferation [false discovery rate (FDR) < 10−5]. This indicated that the observed interspecies divergence in the prefrontal cortex was likely due to a different proportion of progenitor cells in the early fetal human prefrontal tissue samples (fig. S30). In contrast, during postnatal development, PC1 separated prefrontal areas and the inferior temporal cortex (ITC) from the other neocortical areas. This pattern was mainly driven by genes associated with myelination-associated categories (FDR < 0.05; fig. S30) and genes associated with synaptic transmission (FDR < 0.05; fig. S29). Although speculative, these observations potentially link the expansion of the human prefrontal cortex, the wealth of human-specific connectivity made possible by that extension, and the altered patterns of myelination we observe between humans and macaques.

Confirming the observed regional diversification in each species, postnatal development displayed the lowest number of differentially expressed genes between species; most of these (89.3%) were also differentially expressed in adulthood, the phase where we observed the greatest number of interspecies differentially expressed genes (Fig. 3B and table S11). Genes differentially expressed between humans and macaques exhibited distinct patterns of spatiotemporal divergence (Fig. 3C) and showed diverse functional enrichment (table S12). Although 229 genes (2.6%) displayed up- or down-regulation in all the sampled brain regions throughout development and adulthood, others were specifically up- or down-regulated in a subset of brain regions and/or during a particular developmental phase.

To test whether genes with differential expression between humans and macaques showed distinct conservation profiles, we compared values of dN/dS (the ratio of nonsynonymous to synonymous substitution rates) for the whole set of genes differentially expressed in any of the 16 brain regions in at least one of the three developmental phases (32). We found that the differentially expressed genes between humans and macaques also show significantly higher dN/dS values associated with higher evolutionary rates than the remaining protein-coding genes (Wilcoxon-Mann-Whitney P = 2.2 × 10−8, n = 4429 genes). This result was also observed when we focused on the genes differentially expressed in prenatal development (P = 3.7 × 10−11, n = 2380 genes), early postnatal development (P = 4.5 × 10−24, n = 1765 genes), or adulthood (P = 1.0 × 10−6, n = 3837 genes) separately. Moreover, these higher dN/dS values for differentially expressed genes remained highly significant in all the brain regions and developmental phases analyzed, highlighting the consistent association between interspecies transcriptional variation and gene evolution.

Integration with our complementary dataset generated on adult chimpanzee brains (34) revealed that 531 (10.6%), 507 (12.9%), and 1079 (13.9%) genes differentially expressed between species in prenatal development, early postnatal development, and adulthood, respectively, showed human-specific expression in the same brain region in the adult brain. Several genes among those exhibiting species- or human-specific patterns of gene expression were developmentally and regionally regulated. PKD2L1, a gene that encodes an ion channel (52), exhibited human-specific up-regulation only postnatally (Fig. 3C). Conversely, TWIST1, a gene encoding a transcriptional factor implicated in Saethre-Chotzen syndrome (53), showed human-specific down-regulation only postnatally (Fig. 3C). In contrast, MET, a gene linked to autism spectrum disorders (54), showed human-specific up-regulation in the prefrontal cortex and STR postnatally (Fig. 3C). PTH2R, a gene encoding the parathyroid hormone 2 receptor, exhibited macaque-distinct up-regulation in the prenatal NCX but human-distinct up-regulation in the adult NCX, and immunohistochemistry showed that PTH2R is enriched in excitatory neurons (fig. S31). These results show that at least some of the tissue-level interspecies differences we observed are due to changes at the level of specific cell types. Furthermore, even though the ontogenetic and phylogenetic patterns have similar profiles, the overlap of genes driving these two patterns is not substantial (Fig. 3D), indicating the existence of different molecular mechanisms and constraints for regional specification and species divergence.

To gain a more complete understanding of the interspecies transcriptomic differences, we performed an analysis of interspecies differential exon usage as a conservative way of exploring the impact of putative differential alternative splicing. We detected largely similar numbers of genes containing differentially used exons between species in all developmental phases (32) (table S13), with 1924 genes showing interspecies differential exon usage in at least one brain region during the prenatal phase, 1952 during the early postnatal phase, and 1728 during adulthood (Fig. 3B and fig. S32). In our set of differentially used exonic elements, non–protein-coding regions were overrepresented (P < 2.2 × 10−16, χ2 independence test), with 4705 of the 5372 differentially used exonic elements in noncoding regions. This enrichment was especially strong for non–untranslated region (UTR) exonic elements belonging to non–protein-coding transcripts from protein-coding genes and 5′ UTR regions (P < 2.2 × 10−16), but was also significant for 3′ UTR regions (P = 1.81 × 10−11) and non-UTR exonic elements from non–protein-coding genes (P = 0.02364); these results suggest that posttranscriptional regulation may contribute to species differences at the exon level.

Phylogenetic divergence in transcriptional heterotopic regulation

Because transcription factors can regulate the expression of multiple genes, the differential expression we observed between species in different brain regions might be mediated in part by differential expression of a relatively small number of transcription factors. To assess this possibility, we searched for transcription factor binding sites (TFBSs) that were enriched in the annotated promoters of interspecies differentially expressed genes for each brain region and developmental stage in our analysis (32). We found that the binding sites for 86 transcription factors were enriched among interspecies differentially expressed genes; 7 of these 86 transcription factors were differentially expressed between humans and macaques (table S14). RUNX2 was differentially expressed between humans and macaques in the prenatal HIP, PAX7 in the early postnatal AMY, STAT6 in the prenatal NCX, STAT4 in the early postnatal and adult NCX, SNAI2 in the adult CBC, and EWSR1 and NEUROD1 in the adult NCX. Although these enriched motifs were found in only a relatively small proportion of the promoters of the interspecies differentially expressed genes (table S15), expression changes of almost 30% of the differentially expressed genes in the NCX can be explained solely by the transcription factors STAT4, EWSR1, and NEUROD1, which have been previously implicated in neuronal development (55) and brain disorders (56, 57). This suggests that species differences in the expression levels of influential transcription factors could be phenotypically relevant.

To substantiate the possibility that these transcription factors might regulate interspecies differences in gene expression, we next conducted an independent analysis that integrated epigenomic data. We used previously published data on macaque-human differential regulatory elements (active promoters and enhancers) in several regions of adult brains (58). Using region-matched (i.e., NCX, STR, MD, and CBC) aspects of this dataset, we performed TFBS enrichments for the regions defined as up-regulated in humans as well as those down-regulated in humans relative to macaques (32) (tables S16 to S18). As before, we then compared TFBSs enriched among regulatory elements differentially detected in humans and macaques with the transcription factors differentially expressed in a given area or region between species. We observed a higher number of differentially expressed transcription factors associated with binding sites selective for epigenetic loci down-regulated in humans (17, 6, 6, and 1 for NCX, CBC, MD, and STR, respectively) than for loci up-regulated in humans (3, 1, and 1 for NCX, CBC, and MD, respectively). Moreover, 86% of promoters associated with interspecies differentially expressed genes in the NCX contained TFBSs for transcription factors that were differentially expressed between species in the NCX. The same was true for 33% of all differentially expressed genes retrieved from the CBC, 29% for the differentially expressed genes in the MD, and 8.5% of the differentially expressed genes in the STR.

Analysis of epigenomic data (58) in matched brain regions and developmental stages showed that all TFBSs enriched in differentially expressed genes were also found to be enriched in differential regulatory elements. The good agreement between the two independent datasets supports the regulatory relevance of these differentially expressed TFBSs in driving the expression changes of other differentially expressed genes.

Diversity and cell type specificity of species differences

To explore whether cell type–specific transcriptomic changes account for the interspecies divergence observed at the tissue level, we tested the enrichment of human up-regulated genes in human single cells and human down-regulated genes in macaque single cells. Furthermore, we used prenatal scRNA-seq data for prenatal differentially expressed genes and adult snRNA-seq data for the early postnatal and adulthood periods (Fig. 4, A and B, and fig. S33). In all prenatal neocortical areas, human up-regulated genes were enriched in neural progenitors, indicating that the human NCX may possess more neural progenitors at matched time points relative to macaque counterparts, although we cannot completely exclude the possibility that a lack of macaque samples matching human early fetal samples (Fig. 1, A and B) might contribute to this observation, despite the efforts we made to minimize the effects of sampling bias between species by fitting a Gaussian-process model. In contrast, macaque up-regulated genes were enriched in multiple subtypes of excitatory and inhibitory neurons in all neocortical areas (Fig. 4A). Interestingly, a specific subtype of excitatory neurons (i.e., ExN2) was enriched for the macaque up-regulated genes only in prefrontal areas. In the postnatal and adult NCX, human up-regulated genes were enriched in a single population of likely upper-layer excitatory neurons (ExN2b), which was not described in a recent snRNA-seq study of the adult human NCX (59). Conversely, postnatally up-regulated macaque genes were enriched in multiple subtypes of excitatory neurons (Fig. 4B). Interspecies differentially expressed genes in non-neocortical brain regions of the prenatal brain were also enriched in specific cell types (fig. S33). For example, genes displaying interspecies differential expression in HIP and CBC were enriched in a population of oligodendrocyte progenitor cells (OPCs) and external granular layer transition to granule neuron (EGL-TransGraN) cells, respectively. Furthermore, genes showing interspecies differential expression in HIP, AMY, STR, and CBC were enriched in a population of microglia (fig. S34).

Fig. 4 Cell type specificity of species differences.

(A) Cell type enrichment for differentially expressed genes up- or down-regulated in human neocortical areas. Enrichment of genes up-regulated in humans or macaques was tested using single cells from prenatal human NCX (33) or macaque DFC, respectively. The plot shows –log10 P values adjusted by Benjamini-Hochberg procedure averaged across all neocortical areas (NCX), prefrontal areas (PFC), and non-prefrontal areas (nonPFC). Significance (average −log10 P > 2) is labeled by color (red, true; gray, false). (B) Same as (A) for early postnatal and adult periods. (C) Cell type enrichment of selected genes showing human-distinct up- or down-regulation in adult brain regions or neocortical areas (34). Preferential expression measure is plotted to show the cell type specificity.

By integrating our single-cell datasets with a tissue-level transcriptomic dataset of adult human, chimpanzee, and macaque brains (34), we identified the cell type enrichment of several genes showing human-specific up- or down-regulation in NCX or all brain regions relative to chimpanzees and macaques. For example, CD38 was found to be down-regulated in all human brain regions and enriched in astrocytes (Fig. 4C). This gene encodes a glycoprotein that is important in the regulation of intracellular calcium, and its deletion leads to impaired development of astrocytes and oligodendrocytes in mice (60). CLUL1, a gene reported to be specifically expressed in cone photoreceptor cells (61), showed human-specific up-regulation in all brain regions and was enriched in oligodendrocytes and astrocytes. TWIST1 exhibited human-specific down-regulation in all neocortical areas postnatally and was enriched in upper-layer excitatory neurons (Fig. 4C). Conversely, PKD2L1 is up-regulated in NCX postnatally and was enriched in putative deep-layer excitatory neurons (Fig. 4C). MET exhibited human-specific up-regulation in the prefrontal cortex and STR postnatally and was enriched in upper-layer excitatory neurons (Fig. 4C).

Shared and divergent transcriptomic features of homologous cell types

To test whether the observed differential expression between humans and macaques was due to differences in cell type composition or due to transcriptomic differences between homologous cell types, we performed a comparative analysis between human and macaque cell types of prenatal and adult dorsolateral prefrontal cortices. The correlation between human and macaque cell types showed that all human cell types had a close homolog in macaques, and vice versa (Fig. 5, A and B). Nonetheless, we identified genes showing interspecies differential expression in homologous cell types (Fig. 5C). To avoid biases inherent to high variation in scRNA-seq or snRNA-seq, we filtered out genes that did not display differential expression between species at the tissue level and only included genes that exhibited enrichment in cell types where they showed interspecies differential expression [preferential expression measure > 0.3 (32)].

Fig. 5 Shared and divergent transcriptomic features of homologous cell types between humans and macaques.

(A) Dendrogram and heat map showing diversity and correlation of prenatal cell types within and between the two species. The human single cells were from (33). (B) Dendrogram and heat map showing diversity and correlation of adult cell types within and between the two species. (C) Cell type specificity of interspecies differentially expressed genes based on the single cell/nucleus information. Blue, human down-regulated genes; red, human up-regulated genes.

We identified 14 differentially expressed genes in prenatal development and 41 differentially expressed genes in adulthood (Fig. 5C). For example, TRIM54, which encodes a protein implicated in axonal growth (62), was down-regulated in human prenatal neocortical excitatory neurons (Fig. 5C). VW2CL, which encodes a protein associated with α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid (AMPA)–type glutamate receptors (63), was down-regulated in prenatal human neocortical interneurons. SLC17A8 (aka VGLUT3), which encodes vesicular glutamate transporter 3, is up-regulated in human postnatal somatostatin-positive interneurons (InN8). Overall, we found that human DFC cell types showed high correlation with macaque DFC cell types and that only a small set of genes displayed differential expression between these homologous cell types (Fig. 5C). Thus, the interspecies differences identified at the tissue level are likely to result from variations in cellular diversity, abundance, and, to a lesser extent, transcriptional divergence between cell types.

Heterochronic changes in human and macaque brain transcriptomes

The observed heterotopic differences may result, in part, from changes in the timing of gene expression, or heterochrony. To identify such heterochronic differences, we created a Gaussian process–based model [TempShift (32)] and applied this model independently to human and macaque gene expression datasets. To maintain consistency with earlier analyses, we focused our analysis on 11 neocortical areas, which had similar transcriptomic signatures relative to other brain regions [see (33)]. We identified genes with interregional temporal differences within neocortical areas of each species and aggregated them into 36 regional clusters (RCs; fig. S35 and table S19). For both human and macaque brains, analysis of all heterochronic genes revealed greater interareal differences during prenatal periods than at early postnatal or adult ages (fig. S36). In addition, although we observed differences in interareal heterochrony between the early postnatal phase and the adult phase in humans, we did not observe these differences in macaques (fig. S36). This suggests that interregional synchrony in macaques precedes that in age-matched humans, possibly reflecting the protracted development of the human brain during childhood and the earlier plateauing of myelination-associated processes in macaque postnatal development (Fig. 1C and fig. S19).

Analysis of the regional clusters revealed further insights into shared and species-distinct aspects of neurodevelopment. For example, we identified five regional clusters (RC4, 21, 26, 29, and 34) enriched for genes expressed selectively by neural progenitors that exhibited temporal differences between human neocortical areas (fig. S35). Each of these clusters exhibited a gradient whereby a decrease in expression in central regions of the prenatal NCX preceded a decrease at the anterior and posterior poles, suggesting increased progenitor populations or a prolonged neurogenic period in the prefrontal cortex as well as superior temporal cortex (STC), ITC, and V1C. However, although we observed similar temporal gradients in macaques for RC4, 26, and 29, neither RC21 nor RC34—the modules exhibiting the sharpest delay in the posterior NCX—exhibited a similar central-to-polar gradient in macaques (Fig. 6A). Conversely, RC10 and RC12 exhibited an inverse gradient in humans, with decreased expression in the prefrontal NCX, STC, ITC, and V1C preceding a decrease in the central cortex. These modules, which are enriched in astrocytes, did not exhibit a similar gradient in macaques (Fig. 6A and fig. S35). This indicates that even though the transcriptomic signature associated with astrogliogenesis showed a global synchronicity between species (Fig. 1C and fig. S19), a smaller group of genes enriched in astrocytes displayed heterochrony between species.

Fig. 6 Heterochronic expression of regional and interspecies gene clusters.

(A) Clusters of genes exhibiting species-distinct regional heterochronic expression patterns in human and macaque brains at various prenatal periods and adulthood. The timing of expression of genes in the cluster is represented by a color scale (blue, earlier expression; red, later expression). Prenatal heterochronic regional clusters RC21 and RC34 show earlier expression in human prenatal frontoparietal perisylvian neocortical areas (M1C, S1C, and IPC) and enrichment in neural progenitors. RC10 is composed of genes with earlier expression in the human prenatal prefrontal cortex and enrichment in astrocytes. These observed regional expression patterns are not present in the macaque prenatal NCX. Adult heterochronic cluster RC25 shows earlier expression in primary areas of the macaque cortex and enrichment for genes associated with oligodendrocytes. (B) A network of 139 interspecies heterochronic genes (blue) is enriched for targets of putative upstream transcriptional regulators that include those encoded by eight genes of the same network (red) and TWIST1 (green), a transcription factor with interspecies heterotopic expression (fig. S34). Arrows indicate direction of regulation. (C) Top five canonical pathways enriched among interspecies heterochronic genes in at least one neocortical area. The dashed red line corresponds to P = 0.01. (D) Cluster EC14 shows interspecies heterochronic expression, exhibits a delayed expression specifically in the human prenatal prefrontal cortex, and is enriched for genes selectively expressed by intermediate progenitor cells (IPC).

Despite the global enrichment of heterochronic genes in prenatal development (fig. S36), we also identified clusters exhibiting higher interregional differences in postnatal development and adulthood. One example is RC25, a cluster enriched for oligodendrocyte markers that exhibited a pattern of early expression in primary motor and somatosensory areas in the macaque NCX but not the human NCX (Fig. 6A). This finding corroborates myelination-related regional asynchrony (because primary areas myelinate earlier) as well as interspecies heterochrony in oligodendrocyte maturation and myelination-associated processes. Reflective of the cup-shaped pattern of regional variation in global development, the regional clusters also suggest the asynchronous maturation of prenatal areas, a gradual synchronization during early postnatal development in both species, and additional postnatal and adult differences driven in part by myelination.

We next applied TempShift to identify genes exhibiting interspecies heterochronic divergence. Among 11 neocortical areas, we identified approximately 3.9% of coding and noncoding mRNA genes (1100 of 27,932 analyzed orthologous genes) exhibiting interspecies heterochronic expression in at least one neocortical area. We then used Ingenuity Pathway Analysis (Qiagen) to assess upstream transcriptional regulation of heterochronic genes. We found that the differential expression of 139 interspecies heterochronic genes could be explained by as few as eight co-regulated heterochronic transcriptional regulators (Fig. 6B) (32), plus one transcription factor with heterotopic expression (down-regulated in the postnatal human NCX) between species, TWIST1 (fig. S37). A majority (90 of 139) of these putative target genes of the nine transcriptional regulators exhibited accelerated expression in the human NCX. As mentioned above, humans exhibit an accelerated heterochronic pattern for the synaptogenesis transcriptomic signature; the presence of FOS, a neuronal activity–regulated gene, as one of the hubs of this transcriptional network indicates that this accelerated synaptogenesis likely drives the accelerated expression of several genes in the human NCX. Furthermore, an ontological analysis of the genes with heterochronic expression revealed an enrichment for functional categories such as “axonal guidance signaling,” “glutamate receptor signaling,” and “CREB signaling in neurons” (Fig. 6C), which suggests that heterochronic processes include molecular pathways related to axon guidance and synaptic activity.

We next identified 15 evolutionary clusters (ECs) on the basis of the 1100 heterochronic genes displaying interspecies neocortical heterochronic expression patterns (table S20). Among the evolutionary clusters, EC14 exhibited a delayed expression in the human dorsolateral prefrontal cortex and was enriched for intermediate progenitor cell (IPC) markers (Fig. 6D and fig. S38), in agreement with the progenitor cell population differences we observed previously in the prefrontal cortex, indicating that this neocortical prefrontal area likely has a protracted neurogenesis relative to macaques. Similarly, the species-distinct maturation gradients of neural progenitors, astrocytes, and oligodendrocytes also support observations we made concerning interspecies heterotopy. These results were supported by selective validation of the expression profiles of heterochronic genes; using droplet digital polymerase chain reaction, we selected five genes with different developmental profiles across regions and species (figs. S39 to S43), which enabled us to confirm the expression profiles of these genes as well as to ensure that our observations were not the result of biases introduced by TranscriptomeAge.

Species difference in spatiotemporal expression of disease genes

Next, we investigated whether genes associated with risk for neuropsychiatric disorders exhibited differences in their spatiotemporal expression between humans and macaques. We focused our analysis on genes linked to autism spectrum disorders (ASD) and other neurodevelopmental disorders (NDD), attention deficit hyperactivity disorder (ADHD), schizophrenia (SCZ), bipolar disorder (BD), major depressive disorder (MDD), Alzheimer’s disease (AD), and Parkinson’s disease (PD) in previous genetic studies or through our integrative analysis from the accompanying study (33) (table S21). We next sought to determine whether the expression of genes associated with these neuropsychiatric disorders were enriched in any particular developmental phase. Consistent with previous studies associating the midfetal time frame with specific high-confidence ASD (hcASD) genes (64), we found that a larger group of hcASD genes were more highly expressed in the prenatal brains than in the early postnatal and adult brains in both species (fig. S44). In contrast, AD-related genes were more highly expressed in the early postnatal and adult brains than in the prenatal brains in both species (fig. S44). Other groups of disease-related genes did not show any obvious global difference across development. We identified genes with heterochronic or heterotopic expression between the two species that are associated with ASD (6 and 0, respectively), non-hcASD NDD (56 and 14, respectively), and SCZ (45 and 14, respectively) (Fig. 7). This finding potentially suggests the involvement of species-specific aspects in the etiology of ASD, NDD, and SCZ. Unsupervised hierarchical clustering of SCZ-associated genes with heterotopic expression yielded five obvious spatiotemporal clusters, three of which exhibited species differences exclusively during prenatal development (fig. S44). NDD-associated genes with heterotopic expression did not yield any obvious spatiotemporal clusters. Of the prenatal clusters, cluster 1 showed enrichment in the prefrontal cortex, cluster 3 in the temporal cortex, and cluster 2 in both the frontal and temporal cortices, in humans; in macaques, cluster 4 displayed an enrichment in the postnatal and adult frontal cortex, and cluster 5 exhibited a similar enrichment in the adult prefrontal cortex (Fig. 7D).

Fig. 7 Heterotopic and/or heterochronic expression of disease-associated genes between humans and macaques.

(A) Bar plot depicting the number of genes associated with autism spectrum disorder (ASD; hc, high confidence), neurodevelopmental disorders (NDD), attention deficit hyperactivity disorder (ADHD), schizophrenia (SCZ), bipolar disorder (BD), major depressive disorder (MDD), Alzheimer’s disease (AD), and Parkinson’s disease (PD) that display heterochronic divergence between humans and macaques. (B) Bubble matrix showing the heterochronic expression of ASD- and SCZ-associated genes. Blue represents earlier expression in humans; red represents earlier expression in macaques. (C) Bar plot depicting the number of genes associated with neuropsychiatric disorders that exhibit heterotopic divergence between humans and macaques. The 14 SCZ-associated genes that displayed heterotopy are grouped into five clusters on the basis of their spatiotemporal expression profiles (fig. S41). (D) Donut plots exhibiting the centered expression of the five SCZ-associated heterotopic clusters in prenatal development, early postnatal development, and adulthood. Clusters that are not significantly divergent between species in each period are gray and do not have a black border. Red indicates high expression; blue indicates low expression.

Further analysis revealed that the ASD-associated genes SHANK2 and SHANK3, which encode synaptic scaffolding proteins at the postsynaptic density of excitatory glutamatergic synapses, exhibited earlier expression in the macaque NCX and other brain regions relative to humans (Fig. 7B). Commensurate with a role for these proteins in neural circuit development, and in agreement with analyses suggesting the involvement of neocortical projection neurons in the etiology of ASD, these two genes also became progressively more expressed across prenatal ages in both humans and macaques (fig. S45). SCZ-associated genes displaying interspecies heterochrony included GRIA1, a glutamate ionotropic receptor AMPA-type subunit that has different expression trajectories in MFC and OFC relative to other neocortical areas, and that is expressed earlier in human VFC, M1C, S1C, IPC, and STC (Fig. 7B and fig. S45).

These evolutionary changes in the spatiotemporal expression of certain disease-associated genes might therefore imply transcriptional underpinnings for potential human-specific aspects of neuropsychiatric disorders. For example, the presence of human-distinct heterochrony in synapse-related proteins associated with ASD, coupled with the lack of obvious heterotopic expression in hcASD genes, may suggest that conserved neurodevelopmental programs common to primate species are uniquely shifted temporally in some areas in the human brain, potentially implicating key developmental periods, places, and cell types involved in disease etiology. Similarly, the heterochronic and heterotopic changes we associated with SCZ—in particular, those affecting the prenatal prefrontal and temporal cortices—may be involved in human-specific aspect of disease etiology.

Given the importance of UTRs and other noncoding regions in the regulation of gene expression as well as disease, we next explored differences in exon usage between species in genes associated with neuropsychiatric disorders. We observed that 413 genes with differentially expressed exonic elements were linked to the studied diseases. Moreover, we detected 35 disease genes showing differentially used exonic elements with predicted binding sites (65) for microRNAs (miRNAs) independently associated with central nervous system diseases (66) (table S22). Several of these genes (e.g., GRIN2B, BCL11B, and NKPD1) were potentially targeted by a large number of disease-associated miRNAs (fig. S46), and gene-miRNA interactions have already been experimentally validated for 11 of the 35 genes we identified, according to miRTarBase (67) (table S23). For example, we detected differential exon usage of BCL11B, a gene involved in the development of medium spiny neurons (68), between humans and macaques in the adult STR (fig. S46). However, although BCL11B shows lower expression in the human STR than in the macaque STR, the exonic element containing the 3′UTR of BCL11B was itself not differentially expressed. This observation suggests that overexpression in macaques is associated with an alternative isoform containing a shorter 3′UTR region. This shorter 3′UTR lacks predicted binding sites for various miRNAs, including members of the brain-specific miR-219 family, which have been experimentally shown to interact with BCL11B mRNA (69). Together, these findings indicate that certain genes associated with neuropsychiatric disorders exhibit changes in the timing of their expression, location, and splicing pattern between human and NHP brains, and thus may lead to species differences in disease pathogenesis.

Discussion

In this study, we present a comprehensive spatiotemporal transcriptomic brain dataset of the macaque brain. Our integrative and comparative analysis involving complementary humans and adult chimpanzees (33, 34) revealed similarities and differences in the spatiotemporal transcriptomic architecture of the brain and the progression of major neurodevelopmental processes between the two species. For example, we have identified shared and divergent transcriptomic features among homologous brain regions and cell types. We found transcriptomic evidence suggesting that human childhood is especially protracted relative to that of macaques. It has long been recognized that the development of the human brain is prolonged relative to that of other NHPs, and that this slower rate of maturation expands the period of neural plasticity and capacity for learning activities, memory, and complex sensory perception, all processes necessary for higher-order cognition (14, 14, 28). We also found that, relative to macaques, the early periods of human fetal neurodevelopment are transcriptomically distinct and protracted. A similar observation of early neurodevelopmental protraction was recently observed in vitro, in neural progenitors derived from pluripotent cells of human and NHPs (70). However, we also identified cases of neoteny in macaques, such as the protracted postnatal expression of DCX in the hippocampus, likely reflecting differences in neurogenesis between the two species, as recently shown (49).

We found that global patterns of spatiotemporal transcriptomic dynamics were conserved between humans and macaques, and that they display a highly convergent cup-like shape. The sharpest decrease in interregional differences occurs during late fetal ages and before birth; this is likely a consequence of reorganizational processes at this developmental period rather than extrinsic influences due to birth and subsequent events (i.e., respiratory activity or other developmentally novel stimuli). Interestingly, after this transitional period, diversification of neocortical areas appears to be driven mainly by differences between primary and association areas. In addition to these largely conserved broad developmental patterns of interregional differences, we identified numerous genes and gene modules with human-distinct heterochronic or heterotopic expression. These patterns involved brain regions such as the developing prefrontal areas, which are central to the evolution of distinctly human aspects of cognition and behavior (1921). Surprisingly, we also found that developmental phases exhibiting high levels of interregional differences (i.e., early to midfetal periods and young adulthood) were also less conserved between the two species. The coincident convergence of the ontogenetic and phylogenetic cups during the late fetal period and infancy is strikingly distinct from the previously reported phylogenetic transcriptomic hourglass-like pattern that occurs during the embryonic organogenetic period (71, 72).

Genes with divergent spatiotemporal expression patterns included those previously linked to ASD, SCZ, and NDD. These species differences in the expression of disease-associated genes linked to synapse formation, neuronal development, and function, as well as regional and species differences in synaptogenesis and myelination, might have implications for the overall development of neural circuitry and consequently human cognition and behavior. These observations are possibly relevant for recent NHP models of neuropsychiatric disease, such as the SHANK3-deficient macaque model (73), which might therefore not be capable of fully capturing human-distinct aspects of SHANK3 regulation during neurodevelopment.

Our study reveals insights into the evolution of gene expression in the developing human brain. Future work on the development patterns and the functional validation of the genes we report to have heterotopic and/or expression patterns between humans and macaques will likely shed some light on potentially human-specific underpinnings of certain neuropsychiatric disorders.

Materials and methods

Sixteen regions of the macaque brain spanning from early prenatal to adulthood were dissected using the same standardized protocol used for human specimens and described in the accompanying study by Li et al. [(33); see also (32)]. The macaque brain regions and developmental time points matched human brain regions and time points analyzed in (33). The sampled homologous brain regions were identified using anatomical landmarks provided in the macaque brain atlas (74). An overview of dissected brain regions is provided in fig. S1. The Translating Time model (38) was used to identify equivalent time points between macaque and human prenatal development. The list of macaque brains used in this study and relevant metadata are provided in tables S1 and S2. Macaque studies were carried out in accordance with a protocol approved by Yale University’s Committee on Animal Research and NIH guidelines.

We performed tissue-level RNA extraction and sequencing of all 16 regions, scRNA-seq of dorsolateral prefrontal cortex (DFC), hippocampus (HIP), amygdala (AMY), striatum (STR), mediodorsal nucleus of the thalamus (MD), and cerebellar cortex (CBC) of midfetal macaques, and snRNA-seq of DFC of adult macaques. Single cell/nucleus sample processing was done with 10X Genomics and sequencing was done with Illumina platforms.

For tissue-level analysis, we generated annotations of human-macaque orthologs using the XSAnno pipeline, and matched the developmental age of human and macaque samples based on their respective transcriptome using our algorithm TranscriptomeAge. We also developed TempShift, a method based on a Gaussian-process model, to reveal the interregional differences, interspecies divergence, and genes with heterotopic and heterochronic expression. We also queried differentially expressed genes for enrichment in transcription factor binding sites using findMotifs.pl, and analyzed interspecies differential exon usage using the R package DEXSeq.

The single cell/nucleus data were first analyzed by cellranger for decoding, alignment, quality filtering, and UMI counting. After that, data were further analyzed with Seurat according to its guidelines, and cell types were clustered for classification with SpecScore.R. To perform direct comparisons between human and macaque at the single-cell level, we focused on the homologous genes between these species and aligned monkey and human cells together to further analyze interspecies divergence of homologous cell types (fig. S47). We used MetageneBicorPlot function to examine the correlation of neuronal and glial cell subtypes, and we employed correlation analysis to detect the correspondence of excitatory neuron and interneuron subtypes. Finally, we did functional enrichment of disease-associated genes in both tissue-level and single-cell datasets.

Supplementary Materials

www.sciencemag.org/content/362/6420/eaat8077/suppl/DC1

Materials and Methods

Figs. S1 to S47

Tables S1 to S27

References (7598)

References and Notes

  1. See supplementary materials.
Acknowledgments: We thank M. Horn, G. Sedmak, M. Pletikos, D. Singh, G. Terwilliger, and S. Wilson for assistance with tissue acquisition and processing. We also thank A. Duque for using equipment from MacBrainResource (NIH/NIMH R01MH113257). Funding: Data were generated as part of the PsychENCODE Consortium, supported by NIH grants U01MH103392, U01MH103365, U01MH103346, U01MH103340, U01MH103339, R21MH109956, R21MH105881, R21MH105853, R21MH103877, R21MH102791, R01MH111721, R01MH110928, R01MH110927, R01MH110926, R01MH110921, R01MH110920, R01MH110905, R01MH109715, R01MH109677, R01MH105898, R01MH105898, R01MH094714, P50MH106934 awarded to S. Akbarian (Icahn School of Medicine at Mount Sinai), G. Crawford (Duke University), S. Dracheva (Icahn School of Medicine at Mount Sinai), P. Farnham (University of Southern California), M. Gerstein (Yale University), D. Geschwind (University of California, Los Angeles), F. Goes (Johns Hopkins University), T. M. Hyde (Lieber Institute for Brain Development), A. Jaffe (Lieber Institute for Brain Development), J. A. Knowles (University of Southern California), C. Liu (SUNY Upstate Medical University), D. Pinto (Icahn School of Medicine at Mount Sinai), P. Roussos (Icahn School of Medicine at Mount Sinai), S. Sanders (University of California, San Francisco), P. Sklar (Icahn School of Medicine at Mount Sinai), M. State (University of California, San Francisco), P. Sullivan (University of North Carolina), F. Vaccarino (Yale University), D. Weinberger (Lieber Institute for Brain Development), S. Weissman (Yale University), K. White (University of Chicago), J. Willsey (University of California, San Francisco), P. Zandi (Johns Hopkins University), and N.S. Also supported by BFU2017-86471-P (MINECO/FEDER, UE), U01 MH106874 grant, Howard Hughes International Early Career, 3P30AG021342-16S2 (H.Z.); Obra Social “La Caixa” and Secretaria d’Universitats i Recerca and CERCA Programme del Departament d’Economia i Coneixement de la Generalitat de Catalunya (GRC 2017 SGR 880) (T.M.-B.); a Formació de Personal Investigador fellowship from Generalitat de Catalunya (FI_B00122) (P.E.-C.); La Caixa Foundation (L.F.-P.); a Juan de la Cierva fellowship (FJCI-2016-29558) from MICINN (D.J.); and NIH grants MH109904 and MH106874, the Kavli Foundation, and the James S. McDonnell Foundation. Author contributions: A.M.M.S. and N.S. designed the study and procured and dissected all samples; A.M.M.S. and Y.I.K. performed all tissue-level experiments and validations; T.G. and M.S. performed single-cell experiments; M.S., D.J.M., and M.Y. performed single-nucleus experiments; Y.Z. developed TranscriptomeAge and TempShift under supervision of H.Z.; Y.Z., M.L., and G.S. analyzed the data; P.E.-C., D.J., L.F.-P., and T.M.-B. performed transcription factor enrichment, differential exon usage analyses, and evolutionary conservation analyses; and Y.Z., A.M.M.S., F.O.G., and N.S. wrote the manuscript with input of the other authors. Competing interests: Authors have no conflict of interest. Data and materials availability: Data are available at NCBI BioProjects (accession number PRJNA448973) and via Synapse in psychencode.org. All algorithms, packages, and scripts are available at evolution.psychencode.org. Supplement contains additional data.
View Abstract

Subjects

Navigate This Article