Tempo and mode of gene expression evolution in the brain across primates
Abstract
Primate evolution has led to a remarkable diversity of behavioral specializations and pronounced brain size variation among species (Barton, 2012; DeCasien and Higham, 2019; Powell et al., 2017). Gene expression provides a promising opportunity for studying the molecular basis of brain evolution, but it has been explored in very few primate species to date (e.g. Khaitovich et al., 2005; Khrameeva et al., 2020; Ma et al., 2022; Somel et al., 2009). To understand the landscape of gene expression evolution across the primate lineage, we generated and analyzed RNA-seq data from four brain regions in an unprecedented eighteen species. Here, we show a remarkable level of variation in gene expression among hominid species, including humans and chimpanzees, despite their relatively recent divergence time from other primates. We found that individual genes display a wide range of expression dynamics across evolutionary time reflective of the diverse selection pressures acting on genes within primate brain tissue. Using our samples that represent a 190-fold difference in primate brain size, we identified genes with variation in expression most correlated with brain size. Our study extensively broadens the phylogenetic context of what is known about the molecular evolution of the brain across primates and identifies novel candidate genes for the study of genetic regulation of brain evolution.
Editor's evaluation
This is an important study that represents a significant contribution to our understanding of how gene expression in the primate brain has evolved across the extant primate phylogeny. It provides solid evidence for potential links between gene expression variation and brain size, although these are somewhat limited by the focus only on adult brains, since many key changes likely occur during development. Nevertheless, both the taxonomically broad data set and the analysis are likely to be of broad interest to the evolutionary biology, anthropology, and comparative neuroscience communities.
https://doi.org/10.7554/eLife.70276.sa0Introduction
Primates are distinguished from other mammals by their large brains relative to body size (Boddy et al., 2012; Martin, 1981; Smaers et al., 2021). Among the diversity of primate species, there is remarkable variation in behavioral specializations, including differences in social structure, spatial, dietary and visual ecology, and locomotion (Powell et al., 2017; Barton, 2012; DeCasien and Higham, 2019). Despite the impressive array of cognitive attributes displayed by primates, molecular and cellular studies investigating various aspects of brain evolution tend to sample from a small number of species to address questions of how humans are unique. In large part, the emphasis on human brain evolution is warranted. Humans are unmatched in possessing exceptionally large brains and unparalleled cognitive abilities, such as language (Konopka and Roberts, 2016; Rilling, 2014). While valuable, the limited number of species included in prior research lacks a comprehensive perspective of the phylogenetic context in which the human brain evolved within the diversity of primates.
Although researchers have used a variety of approaches to assess whether the human brain is unique (Stout and Hecht, 2017) and how it might have evolved (Hrvoj-Mihic et al., 2013; Sousa et al., 2017a), the full potential for using gene expression to evaluate patterns of brain evolution in primates has not yet been met. Upon observing the remarkable similarity between human and chimpanzee protein sequences, King and Wilson, 1975 proposed that the basis of the physical and behavioral phenotypic differences between these two species must be found in changes within gene regulatory regions that drive expression. Previous studies have explored how changes in regulatory regions can influence gene expression but have often sampled various organs from species across broad spans of evolutionary time, such as mammals or vertebrates (Brawand et al., 2011; Breschi et al., 2016). In studies focusing on gene expression in primate brain tissues, research has mostly focused on the neocortex and cerebellum in human, chimpanzee, and rhesus macaque (Babbitt et al., 2010; Blekhman et al., 2010; Khaitovich et al., 2005; Khaitovich et al., 2004; Khrameeva et al., 2020; Konopka et al., 2012; Ma et al., 2022; Somel et al., 2009; Sousa et al., 2017b). However, new insights can be gained by sampling at greater neuroanatomical resolution from a broader array of primates. Examining gene expression of the brain from a more comprehensive landscape empowers novel inquiry in primate brain evolution, including questions pertaining to the sources of variation that drive expression differences, rates of expression change across the primate phylogenetic tree, and genes that correlate with brain size across primates.
In the current study, we sampled prefrontal cortex (PFC), primary visual cortex (V1), hippocampus (HIP), and lateral cerebellum (CBL) from 18 primate species, the broadest diversity of primates sampled in any study of gene expression in the brain to date, including species from several rarely-studied lineages. Our dataset represents 70–90 million years (Perelman et al., 2011) of primate evolution, providing a more thorough understanding of the evolution of gene expression across primates and allowing for an unprecedented view of how gene expression in the brain has changed over time across all major clades of primate phylogeny.
Results
Most variation can be explained at the species level, not by brain region
To understand how gene expression in the brain has evolved across the primate lineage, we generated and analyzed RNA-seq data from 18 primate species (including five hominoids, four cercopithecoids, four platyrrhines, and five strepsirrhines, with 1–3 biological replicates) across four brain regions, including PFC, V1, HIP, and CBL (Figure 1, Supplementary file 1). The transcriptomes and gene models were assembled de novo (Haas et al., 2013) (see Materials and methods). We quantified the expression of 15,017 orthologs within hominoids, and 3432 on-to-one orthologs across all 18 species (Supplementary files 2 and 3). Variation in interspecific mammalian gene expression has been shown to be less pronounced than that observed across samples from different organs, reflecting the diversity of underlying organ physiology (Brawand et al., 2011; Sudmant et al., 2015). Furthermore, it has been reported that the rate of gene expression divergence evolved more slowly within the cerebral cortex and cerebellum compared to other organ systems from developmentally distinct germ layers (Brawand et al., 2011; Khaitovich et al., 2005). To explore the variability of gene expression from distinct regions of the brain across our broad sampling of primates, we constructed a pairwise distance matrix of the 500 most variable protein-coding genes based on the standard deviation of expression across samples (Methods). This subset of genes was enriched with glycoproteins, signal peptides, and plasma membrane proteins, with roles in immune function, molecular trafficking, and cell signaling. Using this distance matrix, we performed a principal coordinates analysis (PCoA) on data from all brain regions. Because our samples represent disparate regions of the same organ, we expected less variation to be attributed to brain regions than primate species or taxa, reflecting the similarity in physiology of brain tissues. Unsurprisingly, the variation from our complex gene expression dataset is represented across multiple axes of the PCoA (Supplementary file 4).
We plotted the first three axes and created polygons around data derived from samples sharing a common taxa (Figure 2a–c) or region (Figure 2d–f). As predicted, taxon assignment explains a large amount of variation to the dataset, with clear trends emerging independent of brain region. We find the greatest divergence in expression patterns among hominoid and strepsirrhine species, while there is more similarity observed among cercopithecoids and platyrrhines. The hominoids, displaying the greatest level of diversity of any primate phylogenetic group, demonstrate variation that is particularly apparent along Axis 1 and largely driven by human and chimpanzee expression patterns (Figure 2A, B). Strepsirrhines also exhibit a large amount of variation, especially apparent along Axis 2, which can mostly be attributed to the three species of lorises. When Axis 1 and 2 of the PCoA are plotted on the same bivariate plot (Figure 2a), the hominoids display more variation than the strepsirrhines by about 24% (Supplementary file 5). However, a large portion of the variation in the strepsirrhines is attributed to evolutionary divergence over about 63 million years (since the last common ancestor of lemurs and lorises), whereas the variation within hominoids has largely accrued over only 9 million years (since humans and chimpanzees shared an ancestor with gorillas). The hominoid and strepsirrhine samples represent similar variation in terms of sex and life stage, suggesting that these factors do not account for the variability seen in these taxa. Therefore, a remarkable finding of this analysis is how much variation is represented by hominoids, despite the fact that this lineage represents a much shorter evolutionary divergence time.
The cerebellum differs most significantly from the other three sampled brain regions
We observed trends in gene expression by brain region that are predominantly seen along Axis 3 (Figure 2e–f). Here, we observe that the variation attributed to CBL is distinguished from that of PFC, V1, and HIP, which are very similar in their distributions. The fact that CBL differs in its pattern of gene expression from other brain regions (Hawrylycz et al., 2015; Hawrylycz et al., 2012; Itõ, 2012) is not surprising given that it is the only sampled region that develops from a different part of the embryonic neural tube (namely, the hindbrain vs. forebrain) and exhibits a neuronal packing density (predominantly glutamatergic granule cells) that far exceeds these other brain regions (Azevedo et al., 2009; Herculano-Houzel, 2011; Itō, 2012). Enrichments for cerebellar gene expression reveal expression changes to categories such as ‘Cell Surface Receptor Signaling Pathway,’ ‘Cell Projection Organization,’ and ‘Wnt Signaling Pathway.’ This is true for humans in relation to chimpanzees as well as other species with deeper evolutionary relationships. (Supplementary files 6 and 7). Notably, genes involved in cell migration and cell surface receptor signaling Emera et al., 2016 have been shown to mediate the cell-cell interactions necessary for axon guidance (Koropouli and Kolodkin, 2014). Although all the brain regions surveyed show some enrichment for categories related to cell signaling, the cerebellum shows a unique increase in signaling activity-related terms (both in the number of enrichment categories and the degree of change associated with these categories). This is indicative of region-specific differential expression.
Gene expression profiles accurately replicate phylogenetic relationships
To determine whether gene expression profiles can reconstruct known phylogenetic relationships among primates, we built expression phenograms (Methods) for all brain regions combined (Figure 2—figure supplement 2) and each region separately (Figure 2—figure supplement 3). Samples that were derived from individuals of the same species tended to be grouped together, regardless of brain region, revealing that inter-individual differences are minor compared to other sources of variation. Gene expression profiles also replicated the phylogenetic relationships of closely related species (e.g. humans and chimpanzees; pig-tailed and rhesus macaques) when all regions were considered, but these relationships became less phylogenetically structured in the phenograms constructed using expression data from individual brain regions. All neighbor-joining phenograms accurately represent cercopithecoids and strepsirrhines as monophyletic groups; however, expression data produces paraphyletic groups of hominoid and platyrrhine species. This result potentially reflects the fact that taxa with longer periods of independent evolution (i.e. strepsirrhines) are more likely to show divergent patterns of gene expression than more closely related groups. Meanwhile, a more dense sampling of cercopithecoids (three individuals per species), permits a fairly accurate reconstruction of this taxon.
Differential expression in the context of the Ornstein-Uhlenbeck model
Previous studies have used a variety of different approaches to model gene expression changes over time (Brawand et al., 2011; Perry et al., 2012). Here, we used a recently described Ornstein–Uhlenbeck (OU) model to analyze neutral and conserved processes as determined by changing gene expression levels (Chen et al., 2019). OU processes have been proposed to model gene expression evolution as they model both drift and stabilizing selection (Rohlfs et al., 2014). Previous studies have shown how models that incorporate stabilizing selection are more accurately able to predict gene expression evolution in mammals than models that account only for neutral drift (Chen et al., 2019). Across all expressed one-to-one orthologs represented in the sampled primates, we found that ~15–20% of genes show differential expression across all species-to-species comparisons (q-value <0.05). As expected, the relative amount of differential expression increases over evolutionary time, both between species and clades (Figure 3 and Figure 3—figure supplement 1 ). However, both the species and clade-wise comparisons show larger numbers of differentially expressed (DE) genes in the comparisons when strepsirrhines are included in the contrast, a result of the more than 60 million years of independent evolution of this taxon. In the prefrontal cortex, we see slightly smaller DE in humans as compared to siamang and baboon; however, in total number, this DE is not appreciably different (DE Human-Siamang n=251, DE Human-Baboon n=255) to other comparisons (DE Human-Rhesus Macaque n=315). Upset plots allow us to further analyze these unique patterns of differential expression (Figure 3—figure supplement 2). In the cerebellum, the lemur, baboon, and human samples are particularly unique and show relatively higher levels of DE compared to other species. However, in the hippocampus, visual cortex, and prefrontal cortex, the chimpanzees appear to show higher DE than the human samples. This suggests again that the cerebellum is a particularly unique brain structure and that chimpanzee gene expression is significantly different in the other three brain regions studied, warranting further analysis.
Human-specific enrichment for metabolic processes, neural development, and gene regulation
When comparing gene expression in the PFC of humans relative to other primates, human PFC shows an enrichment of metabolic processes, including ‘regulation of cellular metabolic process’ and ‘regulation of macromolecule metabolic process’ (Supplementary files 6 and 7). Comparing human and chimpanzee PFC reveals that categories that support neural growth and development (e.g. ‘neuron projection morphogenesis,’ ‘cell morphogenesis involved in neuron differentiation’), gene regulation, and metabolic processes are enriched as differentially expressed in human PFC relative to chimpanzees.
In addition to examining human and chimpanzee data in isolation, we also analyzed human-specific changes in expression within the context of other outgroup species. Using the Ornstein-Uhlenbeck process as a model of continuous trait evolution across our 18-species primate phylogeny, we again observe similar categories of enriched processes for the human prefrontal cortex in comparison to that of chimpanzees. These terms include ‘Nervous System Development,’ ‘Neurogenesis,’ ‘Glial Cell Differentiation,’ ‘Neuron Projection Morphogenesis,’ ‘Regulation of Gene Expression,’ and ‘Metabolism.’ To further validate our results, we also looked at the human PFC in comparison to other primate species. In analyzing differential expression between the human and siamang PFC, we note that similar trends for enrichment are also found, such as ‘Neural Growth and Development’ and ‘Metabolic Processes,’ and ‘Gene Regulation.’ Of interest, under the category of ‘Positive Regulation of Transcription by RNA polymerase II’ we find several genes that appear to be upregulated in humans compared to siamang: APP (amyloid precursor protein, related to plaque formation in Alzheimer’s disease) as well as PRKN (found to be causal in Parkinson’s disease) (Funayama et al., 2023). This supports the idea that these enrichments are human-specific, have relevance to important human neurodegenerative disease states, and are not a reflection of changes occurring within the chimpanzee lineage.
Broader species comparisons show similar trends across evolutionary time
Beyond human and chimpanzee comparisons, we also note many interesting broader temporal trends observed from the EVEE-based differential expression analysis. When examining Hominini (humans and chimpanzees) compared to other Great Apes, terms related to ‘Regulation of Metabolic Process,’ ‘Nervous System Development,’ and ‘Biosynthetic processes’ are all enriched within the hippocampus. Meanwhile, ‘Negative Regulation of Synaptic Transmission’ was enriched in the other ape species. In comparing the PFC of the Hominoid clade to that of the Strepsirrhine clade, we found that ‘Neuron Development and Differentiation’ was enriched. Overall, among various species and clade comparisons, there is a general trend of decreasing specificity in enrichment categories over increasing evolutionary time. Using PFC data, we found that the relationship between humans and chimpanzees, some of the closest relatives in our dataset, shows terms related to ‘Synapse Assembly,’ ‘Regulation of Glial Cell Differentiation,’ ‘Regulation of Astrocyte Differentiation,’ ‘Axonogenesis,’ and ‘Neuron Projection Morphogenesis.’ Looking at a more distantly related species pair, the human and rhesus macaque comparison shows enrichment for terms related to ‘Cell Growth,’ ‘Cell Development,’ ‘Biological Regulation,’ ‘Neuron Projection Development,’ ‘Regulation of Neurogenesis,’ and ‘Positive Regulation of RNA Biosynthetic Processes.’ The most distantly related species, humans and lemurs, have overall the largest number of differentially expressed genes, and with that, the broadest categories of enrichment. These include categories such as ‘Regulation of Developmental Processes,’ ‘Regulation of Nervous System Development,’ ‘Cell Development,’ and ‘Multicellular Organism Development’.
Our PCoA analyses showed that gene expression in brain regions sampled from the lorises (i.e. the slender loris, slow loris, and pygmy slow loris) diverged from other strepsirrhines, and other primates more generally (Figure 2—figure supplement 1). When strepsirrhines are compared to other primates in differential expression analyses, transcription factors, and other genes involved in gene transcription and translation and multiple biosynthetic pathways involved in cellular metabolism are among the categories of DE genes (Supplementary files 6 and 7).
Evolutionary rates of expression change across clades and brain regions
Using the OU model, we found that individual genes exhibit wide variation in expression dynamics across the primate lineage (Figure 4a). Enrichments for genes showing low variation or stabilizing selection (q=0.05) reveal categories related to transport and cellular localization (GO Biological Processes, Supplementary files 6 and 7). In contrast, genes that are less constrained or neutrally evolving (q>0.05) have a number of processes related to neuron morphogenesis, plasticity, and cell death. Yet, unlike sequence evolution, gene expression is not linear across evolutionary time but a saturation point in pairwise comparisons of gene expression is reached due to stabilizing selection pressures. Here, we find that pairwise expression differences between humans and the other species increasingly diverge with evolutionary distance in all brain regions sampled (Figure 4b); however, these pairwise comparisons do not seem to saturate with evolutionary time across the primate comparisons (Chen et al., 2019). We note that the saturation of pairwise expression differences from humans may be found at a phylogenetic node ancestral to primates (Figure 3—figure supplement 1).
Expression of a majority of genes evolves under stabilizing selection
We utilized EVEE-tools, developed in the context of the Ornstein-Uhlenbeck process of continuous trait evolution, to classify genes as primarily under the effects of stabilizing selection vs. a model of neutral drift (Chen et al., 2019). In this analysis, we found that across tissues, on average, 64% of genes fit better under stabilizing selection (64% CBL, 59% HIP, 72% PFC, and 60% V1; FDR-corrected q-value calculated via the BH procedure to correct for multiple hypothesis testing; FDR threshold of 5% to determine significance). In the context of the OU model of continuous trait evolution, we found that the overall phylogenetic signal in brain expression divergence was slightly smaller than observed with edgeR over the entire combined dataset (EVEE: 4–6.7% across all single species comparisons using a logFC of +/-2; edgeR: 15–20%). However, the relative amount of differential expression does increase gradually with evolutionary distance (as expected based on results in edgeR). The only exception to this is the human and chimpanzee comparison, which shows considerable variation in differential expression across the four brain regions. Representing such a short period of evolutionary divergence, this increase in DE suggests direction selection within that lineage.
Evidence of potential directional selection in human and chimpanzee data
To better understand the outlier gene expression in this dataset, and overall to gain insight into genes that may be subject to directional selection pressures, we again used the EVEE-tools OU-based model to score our dataset for outlier expression. This requires the determination of the evolutionary mean and variance for each gene across our entire expression dataset, from which we then can compare to individual species expression data. Using this method, we found a small subset of genes to be defined as having an expression that deviates from the optimal OU-distribution (Z-score >2 or < –2, p-value <0.05) (Figure 4Figure 4—figure supplement 1). It is important to note that a significant FDR (<5%) was not reached in this dataset, among all species comparisons. This is expected based on previous applications of this package, in which a mammalian dataset was also unable to reach an FDR below 18% (Chen et al., 2019). This is likely a reflection of the limitations of our phylogeny, and suggests that future projects should aim to sample even more broadly.
In our outlier expression analysis, we found that patterns of outlier gene expression occur in a species and tissue-specific manner (Supplementary file 10). For example, in the PFC, the chimpanzee and marmoset samples appear to have the highest number of outlier genes deviating from average expression patterns (compared to human, siamang, baboon, rhesus macaque, and marmoset samples). We found this to be true regardless of how the dataset was normalized in order to determine the average expression for each gene (via defining a reference species). This is particularly interesting and, in combination with differential expression analysis, highlights the chimpanzee PFC as a particularly divergent structure. In contrast to the PFC, outlier analysis in the CBL reveals that the human and lemur samples have the highest number of outlier genes.
Upon gene set category enrichment analysis, we find that many of these genes that are deemed ‘outliers’ are related to functions in development, transcription, nervous system development, neurogenesis, and metabolism. For example, the chimpanzee PFC shows a significant upregulation in genes involved in energy storage and transfer, such as ETFA and NDUFS4, which are both involved in electron transport for ATP generation, as well as ANKH, implicated in phosphate transport (Szeri et al., 2020; Henriques et al., 2021; Shil et al., 2021) We also see that the chimpanzee PFC shows unique expression patterns in genes related to synaptic activity and neurotransmitter release, including significant downregulation of GABRA4 and SYN1 (Fassio et al., 2011; Fan et al., 2020 Fan et al., 2020; Fassio et al., 2011). In contrast, the human CBL and PFC both display a significant upregulation of genes related to Amyloid protein production (APP), a major component of many neurodegenerative diseases with functions in synaptic signaling (O’Brien and Wong, 2011). Unique to the human CBL we also see enrichments in genes related to neurogenesis and synaptic activity, including SDK1, FZD5, and CDH10 (Redies et al., 2012; Slater et al., 2013; Bagot et al., 2016) This is suggestive of directional selection pressures occurring in a tissue-specific manner and encourages future investigation of these outlier genes.
Implications of using humans as a reference species
We continue to use humans as a reference species in these analyses as, compared with other primates, humans have exceptionally large brain sizes and unique cognitive abilities. However, we do recognize that there are some implications for having humans as a reference, especially given our data that would suggest human gene expression as being largely different from the rest of our primate dataset. To address this, we repeated analyses using EVEE-tools by including two additional reference species: siamang and rhesus macaque. The percentage of genes that fit better under the model of stabilizing selection (in comparison to neutral drift) is not statistically different from those observed when using humans as a reference (on average, 58–64% across all four major brain regions). We additionally looked at pairwise species comparisons to determine if the general trends of directional selection and differential gene expression were comparable to the human-reference data, and again confirmed that the effects of the reference species used for normalization here are negligible in terms of analyzing differential expression across our primate tree.
Correlation of gene expression to brain size and their change over evolutionary time
Because an increase in absolute brain size is one of the most striking characteristics of humans, we asked what subset of genes is correlated with this important biological trait across the primate tree. Our dataset provides a unique opportunity to evaluate this question, since the average brain size across our study species varies by ~190 fold. Using multivariate analysis (Methods), we defined gene lists most strongly correlated with brain size in each brain region (Supplementary file 8). Results indicate that the same genes rank among those with the strongest positive correlation in PFC, V1, and HIP. CBL also shares some of these same genes but includes more variation among the genes most strongly correlated to overall brain size than the other three brain regions (Figure 5, Supplementary file 9), potentially reflecting the more recent expansion of the CBL relative to the rest of the neocortex (Miller et al., 2019; Smaers et al., 2018).
Discussion
We performed RNA-seq on four brain regions from 18 primate species, representing the broadest sampling for any gene expression study in primate brain tissue to date. Through more representative sampling of primate species, we found substantial variation in gene expression levels within the hominoid and strepsirrhine lineages, with the diversity among hominoids particularly impressive due to the recent divergence of this taxon. Using the OU model, we found that a substantial proportion of genes showed differential expression across species. The relative amount of differential expression increased over evolutionary time, both between species and clades. Additionally, when comparing gene expression across broader species and clade comparisons, we observed trends related to brain development, nervous system regulation, and cellular metabolism.
Our findings point to human-specific enrichment for metabolic processes, neural development, and gene regulation. The considerable diversity of gene expression in human and chimpanzee brain tissue has profound implications for understanding the distinct evolutionary processes that have acted upon the brains of the ancestral species of these two lineages. We observed a wide variety of expression dynamics of individual genes in the pairwise comparisons of humans to other primate species. Enrichments for genes under stabilizing selection indicated processes related to transport and cellular localization, while less constrained genes were associated with neuron morphogenesis, plasticity, and cell death. Importantly, gene expression evolution did follow a linear pattern, but did not reach a saturation point due to stabilizing selection pressures as seen in other studies (Chen et al., 2019). Less constrained, neutrally evolving patterns appeared to be the most prevalent pattern in each brain region studied, preventing a saturation point of stabilizing pressures to be reached within Primates with increased phylogenetic distance from humans. Lastly, we identified genes that are correlated with brain size across all major primate taxa, providing candidates for further inquiry.
Our deeper analysis of gene expression has revealed evolutionary patterns that were inaccessible with a more limited sampling of primate brain tissue. We anticipate that the candidate genes and data provided by this study will serve as a resource for many other lines of inquiry into human and non-human primate brain evolution.
Materials and methods
Biological sample collection and RNA extraction
Request a detailed protocolThe sample includes brain tissue from human and nonhuman primates. All samples were obtained from adult individuals free from known neurological disease. If available, the right hemisphere was preferentially sampled. Human brain samples were obtained from the National Institute for Child Health and Human Development Brain and Tissue Bank for Developmental Disorders at the University of Maryland (Baltimore, MD). Chimpanzee brain tissue was obtained from the National Chimpanzee Brain Resource (https://www.chimpanzeebrain.org, supported by NIH grant NS092988). All other sources of brain tissue are listed in Supplementary file 1.
From each individual, we sampled four regions of the brain, including PFC, V1, HIP, and CBL. PFC was sampled from the frontal pole, corresponding to Brodmann’s area 10 in humans. In other primates, the PFC region sampled more broadly encompassed prefrontal cortical areas but was limited to the most anteriorly projecting part of the frontal pole. All V1 samples were dissected around the calcarine sulcus to include primary visual cortex (Brodmann’s area 17). Samples from both PFC and V1 contained all cortical layers and a small amount of underlying white matter (<10%). The HIP was sampled from the medial aspect of the temporal lobe and included all hippocampal subfields. The CBL was sampled from the most laterally projecting region of lateral hemisphere in all primates. In humans, the CBL region corresponded to Crus I or Crus II. CBL samples contained all layers of cerebellar cortex and a small amount of underlying white matter (<10%). Each sample was briefly homogenized using a Tissuelyzer (Qiagen), and the total RNA was isolated using an RNAeasy kit (Qiagen) with a DnaseI treatment.
Library preparation and sequencing
Request a detailed protocolSingle-end RNA-seq libraries were made using the NEBNext mRNA Library Prep Reagent Set for Illumina. Libraries were prepared in batches of 4–8 samples of randomly sampled species and brain regions. Library sizes were checked on the Bioanalyzer (Agilent). RNA-seq libraries were multiplexed on the NextSeq500 (Illumina) in the Genomics Resource Laboratory at the University of Massachusetts Amherst, also randomly distributed across NextSeq500 runs. All fastq files have been submitted to the SRA: https://www.ncbi.nlm.nih.gov/bioproject/?term=PRJNA639850.
Mapping and transcriptome analysis
Request a detailed protocolSequencing reads were assembled into species-specific transcriptomes (containing the reads from all four brain regions) using Trinity (Grabherr et al., 2011). With Trinity assembly, fragments of the original RNA reads for each species are compiled and clustered into groups based on sequence similarity, which eventually are extended to reconstruct full-length transcripts (Grabherr et al., 2011). These transcripts were then blasted against the human blast nt database, with the alignment thresholds for the top hits from different clades listed in Supplementary file 2, similar to the approach in Perry et al., 2012. Individual libraries were then mapped to the species-specific transcriptome using bowtie (Langmead et al., 2009) in RSEM, and count tables were generated using RSEM (Li and Dewey, 2011). Orthology assignments were additionally checked using the Ensemble one-to-one orthology alignments as a guide for the subset of species with a publicly available genome. Transcriptome quality was assessed using BUSCO to determine assembly and annotation completeness (Seppey et al., 2019), Supplementary file 3. We recognize that some of the species transcriptomes show relatively lower BUSCO completeness scores (namely the Slender Loris at 35.8% complete). We hypothesize that this is likely a reflection of the limited tissue sampling in this dataset of only brain tissue. Previous studies have shown that transcriptome assemblies from single tissue regions on average have lower completeness scores than assemblies composed of reads from a variety of tissue types (Simão et al., 2015). This is likely a reflection of tissue-specific gene expression. As evidence of this, we further analyzed reads deemed as ‘missing’ by BUSCO and found that many of these showed little to no expression across the human brain (Uhlén et al., 2015) and Human Protein Atlas proteinatlas.org. We do not make any major conclusions about the Loris species in this manuscript and thus do not believe these BUSCO scores significantly affect the conclusions made in this manuscript. We only consider Loris data in concert with Lemurs, which by comparison have much more complete transcriptomes.
Distance-based data analyses (PCoA and phenograms)
We performed principal coordinates analyses (PCoA) based on a pairwise distance matrix of all 137 samples. The distance matrix was comprised of the top 500 most variably expressed protein-coding genes. Pairwise distances were calculated by leading log2 fold change, providing a symmetrical representation of the expression ratio centered around 0 (i.e. log2(2)=1 while log2(0.5) = –1) (Robinson et al., 2010). Creating the distance matrix and plotting the PCoA were performed using the plotMDS.DGElist function (based in limma) in the edgeR package in R. Although variation is represented across more than 20 axes (Supplementary file 4), the first three axes were plotted to compare patterns across primate taxa and brain region sampled (Figure 2, Figure 2—figure supplements 1–3). Polygons overlap the data points representing taxa or brain region. The area of each polygon was computed using functions in the sp package of R. The chull() function was used to define the points around the perimeter of each polygon, and the Polygon() function calculated the area of each. Relative areas of each polygon are listed in Supplementary file 5. Figure 3 displays the same data as the PCoA but uses an array of colors allowing the data from each individual species to be visualized.
The same log2 fold change distance matrix was then used to create phenograms representing the similarity of gene expression profiles among samples. The minimum distance neighbor-joining function in the ‘ape’ package of R created a tree based on the method proposed by Saitou and Nei, 1987. The boot.phylo function estimated the reliability of given nodes of the tree by resampling over 1000 iterations. Although our objective in this analysis was to investigate patterns of evolution across the primate order and the brain regions, our sample included multiple individuals from the same species. By treating these samples separately, our analyses represent both within- and between-species variability in gene expression over time.
Analyzing differential expression
Request a detailed protocolCounts were filtered and normalized using edgeR (Robinson et al., 2010), with any multispecies comparisons using the GLM functionality (McCarthy et al., 2012). Gene Ontology enrichments were performed using the DAVID gene ontology tool 6.8 (Huang et al., 2009a; Huang et al., 2009b) and g:Profiler (Reimand et al., 2016). Supplementary files 6 and 7 shows results from ordered g:profiler enrichments (g:GOSt) performed on DE genes where q<0.05 (note, this is not ranked on polarity of expression, just absolute change) with all genes expressed in this study used as background.
Differential Expression was also analyzed using the package EVEE-tools to incorporate the Ornstein-Uhlenbeck model (Chen et al., 2019). Evolutionary means and variance values were calculated for each gene across the entire phylogeny for individual brain-region datasets as well as the entire RNA-seq bulk dataset in the context of the Ornstein-Uhlenbeck model. Differential expression was determined using a multivariate OU model. Regimes were defined by species/clades of interest (i.e. in human vs. chimpanzee comparisons, the human samples represent one regime in one model to compare to all other samples, represented by a second regime. In comparison, the chimpanzee samples would be represented as an additional model and a separate regime.) p-values were calculated to represent the fit of each OU model (i.e. Human or Chimp-specific expression patterns) against a Brownian motion model. These were then corrected for multiple hypothesis testing using the Benjamini-Hochberg FDR procedure. We used an FDR threshold of 0.05 to define significance. Additionally, Akaike and Bayesian Information Criterion (AIC and BIC) scores were calculated for each gene for each model, and only genes with AIC and BIC scores significant against the null were considered in further analysis. Directionality of differential expression was further determined by comparing estimated mean expression levels for each regime in each model.
We also looked at differential expression beyond pairwise comparisons, again using EVEE-tools, to see how each individual species differed across the entire dataset. With this, one species of interest was treated as a single regime while all 17 other species were grouped as a second regime. The same criteria as above were used to determine significant differential expression, and this data was utilized to construct UpSet plots using the UpsetR package (Conway et al., 2017). We included only relevant species in these UpSet plots (representing each major primate clade) to simplify the graphs. Each graph is also separated by brain region, similar to previous analyses (Figure 3—figure supplement 2).
We also analyzed sources of variation in our data to determine how significant of an effect species differences have compared to other factors, including primate families and individual variation. For this, we conducted a correlation analysis utilizing the non-parametric Mann-Whitney U (MWU) test on the differentially expressed genes between humans and chimpanzees to other primate species. We focused this analysis on only those species with three individuals per brain region (olive baboon, rhesus macaque, and lemurs). Using the set of human-chimpanzee DEGs, we calculated the Spearman rank-based correlation coefficient between each species to either human or chimpanzee expression. We determined whether or not these correlation coefficients were significantly different across species and across brain regions using the Mann-Whitney U test. We determined that there was no significant correlation between any of the three species to human or chimpanzee expression in any of the four brain regions (Figure 3—figure supplement 1; Yapar et al., 2021). This suggests that the expression profiles of these more distantly related primates are equally similar to human and chimpanzee expression patterns. We also conducted additional Analysis of Covariates (ANCOVA) to confirm that other factors, namely age and sex, were not significant sources of variation in our gene expression analyses. These analyses, along with our PCoA plots, show that taxon identity and brain region are the two most significant determinants for DE. Additionally, differences in samples from the same individual are defined by differences at the level of brain region. ANCOVA analyses showed minor residual effects that we deemed as random and were not further analyzed for the purposes of this manuscript.
Outlier expression was determined using the scoreGenes.R script from the EVEE-Tools script suite (Chen et al., 2019). We compared the total dataset-normalized mean expression of a single gene to that of a single select species, in the context of the overall mean expression and variance across the dataset. We specifically analyzed outlier expression in a brain region-specific manner, looking at the datasets subset by individual brain region. To normalize the expression of the entire dataset, a single species was selected to be used as a reference in TMM normalization. For all outlier analyses except for humans, the human samples were used as a reference. For the human outlier analysis, the rhesus samples were used as a reference for normalization. Importantly, we tested the use of different reference species for dataset normalization and did not find a significant difference in the number of outlier genes and the enrichment categories associated. Z-scores were calculated for genes whose expression patterns fit an OU model (in comparison to the null model of Brownian Motion) and whose evolutionary mean is above 5 CPM for the entire dataset (Chen et al., 2019). In order to determine significant outliers, an FDR threshold of 0.05 was again employed, however at this level limited significance was found. Supplementary file 10 shows the results of this analysis for the CBL and PFC brain regions.
Phylogenetic and evolutionary distance analysis
Request a detailed protocolCategorical enrichments for the contrasts between species and clades are in Supplementary file 6. The phylogenetic tree for primates was downloaded from the UCSC Genome Browser (30 primate species) (Kent et al., 2002). Distances between species were extracted using the Environment for Tree Exploration Toolkit (Huerta-Cepas et al., 2010). The residuals and mean squared expression differences of all orthologs across 18 species were found using the package EVEE (Chen et al., 2019), and in all contrasts, humans were used as the reference species. We then analyzed the subset of genes showing either broadly defined conserved (low) or neutral (higher) variation across species (low <q = 0.05, high q>0.05), with categorical enrichments for these two groups in Supplementary file 7.
Comparisons of the heterogenous tissues used and single-cell gene expression data
Request a detailed protocolIn any study that derives results from homogenized tissue samples, the composition heterogeneity of the samples may drive differences in gene expression (Montgomery and Mank, 2016). To address this issue, we compared the expression of our tissue samples to recent studies that have performed single-cell RNA-seq on neurons and astrocytes. RNA-seq data from primary neurons and astrocytes were obtained from NCBI’s Gene Expression Omnibus (GEO) and processed in the same manner as the tissue samples for all human samples. These included four hippocampal astrocytes, four cortical astrocytes, and one cortical neuron from Zhang et al., 2016 (GEO accession number GSE73721) and three pyramidal neuron samples isolated from an unspecified brain region by the ENCODE project (ENCODE Project Consortium, 2012; Davis et al., 2018) (accession numbers GSM2071331, GSM2071332, and GSM2071418). Only genes with counts greater than zero in all samples and (CPM)>1 in all 23 samples were included in this analysis (n=7111). A PCoA was made from a distance matrix of the top 500 most variably expressed genes by the pairwise biological coefficient of variation (method = “bcv”) across samples (Robinson et al., 2010). Creating the distance matrix and plotting the PCoA were performed using the plotMDS.DGElistfunction in the edgeR package in R. The PCoA of our human samples in comparison to primary neurons and astrocytes suggests that our heterogeneous tissue samples are not biased to contain more neurons or astrocytes as compared to each other (i.e. one tissue is not biased within this small sample set) (Figure 5—figure supplement 1), and is consistent with other neural cell and brain tissue comparisons (Khrameeva et al., 2020).
Additionally, prior knowledge of variation across primates in cell type composition of the brain is informative in interpreting bulk RNA-seq data. It is well known that neuron densities tend to decrease as brain size increases (Sherwood et al., 2020). This suggests that larger brains accommodate a smaller number of neurons per unit volume compared to smaller brains. However, it is interesting to note that other structural elements, such as astrocytes (Munger et al., 2022), microglia (Dos Santos et al., 2020), and synapses (Sherwood et al., 2020), exhibit a relatively invariant density per unit volume across species. Despite changes in brain size, these essential components involved in neural communication and support maintain a consistent presence, emphasizing their crucial role in brain function regardless of the species' brain size. Thus, based on this previous research, in cross-species comparisons of bulk tissue, we can tentatively interpret differences in gene expression to reflect generally similar proportions of major cell types, except for neurons, which are expected to decrease in proportion with larger brain sizes.
Brain size analyses
Request a detailed protocolAverage species endocranial volumes (ECVs) were obtained primarily from Kamilar and Cooper, 2013 and Isler et al., 2008 from the mean of male and female volumes. ECVs were used since reliable brain size data does not exist for all species samples. The data for human ECV (also averaged from male and female data points) was previously published by Coqueugniot and Hublin, 2012. Isler et al., 2008 reported a minimal error when ECV was transformed to brain mass using the correction factor of 1.036 g/ml (Stephan, 1960), and we used this conversion to obtain brain mass estimates from ECVs for all species (Supplementary file 8).
Within each brain region, we performed a PCA on the species average gene expression of the 500 most variable genes by standard deviation, using the prcomp function in R. For each regional PCA, 14 PCs were required to account for about 90% of the variance in gene expression. We performed multiple regression analyses to determine which of the PCs could predict brain size. Using all 14 PCs accounted for at least 95% of the variation in each brain region. Akaike information criterion was applied. However, it was noted that for each brain region, PC2 was the most predictive of brain size by low (regional adjusted R2 values for PC2 against brain size were: PFC, 0.42; V1, 0.56; HIP, 0.50; CBL, 0.36). The 500 genes and their loadings on PC2 are listed in Supplementary file 9. Across the four sampled regions, we find general uniformity in the extent to which individual genes affect brain size, but the cerebellum displays the most unique signature of the regions sampled (Figure 5).
Data availability
Sequencing data have been deposited in the Short Read Archive: BioProject PRJNA639850.
-
NCBI BioProjectID PRJNA639850. Tempo and mode of gene expression evolution in the brain across primate tree.
-
NCBI Gene Expression OmnibusID GSE73721. RNA-Seq of human astrocytes.
-
NCBI Gene Expression OmnibusID GSE78331. single cell RNA-seq from pyramidal cell (ENCLB928LID).
References
-
Equal numbers of neuronal and nonneuronal cells make the human brain an isometrically scaled-up primate brainThe Journal of Comparative Neurology 513:532–541.https://doi.org/10.1002/cne.21974
-
Both noncoding and protein-coding RNAs contribute to gene expression evolution in the primate brainGenome Biology and Evolution 2:67–79.https://doi.org/10.1093/gbe/evq002
-
Embodied cognitive evolution and the cerebellumPhilosophical Transactions of the Royal Society of London. Series B, Biological Sciences 367:2097–2107.https://doi.org/10.1098/rstb.2012.0112
-
Sex-specific and lineage-specific alternative splicing in primatesGenome Research 20:180–189.https://doi.org/10.1101/gr.099226.109
-
Age-related changes of digital endocranial volume during human ontogeny: results from an osteological reference collectionAmerican Journal of Physical Anthropology 147:312–318.https://doi.org/10.1002/ajpa.21655
-
The encyclopedia of DNA elements (ENCODE): data portal updateNucleic Acids Research 46:D794–D801.https://doi.org/10.1093/nar/gkx1081
-
Primate mosaic brain evolution reflects selection on sensory and cognitive specializationNature Ecology & Evolution 3:1483–1493.https://doi.org/10.1038/s41559-019-0969-0
-
SYN1 loss-of-function mutations in autism and partial epilepsy cause impaired synaptic functionHuman Molecular Genetics 20:2297–2307.https://doi.org/10.1093/hmg/ddr122
-
Molecular genetics of parkinson’s disease: contributions and global trendsJournal of Human Genetics 68:125–130.https://doi.org/10.1038/s10038-022-01058-5
-
Full-length transcriptome assembly from RNA-Seq data without a reference genomeNature Biotechnology 29:644–652.https://doi.org/10.1038/nbt.1883
-
Canonical genetic signatures of the adult human brainNature Neuroscience 18:1832–1844.https://doi.org/10.1038/nn.4171
-
Evolution, development, and plasticity of the human brain: from molecules to bonesFrontiers in Human Neuroscience 7:707.https://doi.org/10.3389/fnhum.2013.00707
-
ETE: a python environment for tree explorationBMC Bioinformatics 11:24.https://doi.org/10.1186/1471-2105-11-24
-
Endocranial volumes of primate species: scaling analyses using a comprehensive and reliable data setJournal of Human Evolution 55:967–978.https://doi.org/10.1016/j.jhevol.2008.08.004
-
Phylogenetic signal in primate behaviour, ecology and life historyPhilosophical Transactions of the Royal Society of London. Series B, Biological Sciences 368:20120341.https://doi.org/10.1098/rstb.2012.0341
-
Regional patterns of gene expression in human and chimpanzee brainsGenome Research 14:1462–1473.https://doi.org/10.1101/gr.2538704
-
Semaphorins and the dynamic regulation of synapse assembly, refinement, and functionCurrent Opinion in Neurobiology 27:1–7.https://doi.org/10.1016/j.conb.2014.02.005
-
Differential expression analysis of multifactor RNA-Seq experiments with respect to biological variationNucleic Acids Research 40:4288–4297.https://doi.org/10.1093/nar/gks042
-
Comparative analysis of astrocytes in the prefrontal cortex of primates: insights into the evolution of human brain energeticsThe Journal of Comparative Neurology 530:3106–3125.https://doi.org/10.1002/cne.25387
-
Amyloid precursor protein processing and alzheimer’s diseaseAnnual Review of Neuroscience 34:185–204.https://doi.org/10.1146/annurev-neuro-061010-113613
-
A molecular phylogeny of living primatesPLOS Genetics 7:e1001342.https://doi.org/10.1371/journal.pgen.1001342
-
Re-evaluating the link between brain size and behavioural ecology in primatesProceedings. Biological Sciences 284:20171765.https://doi.org/10.1098/rspb.2017.1765
-
Cadherins and neuropsychiatric disordersBrain Research 1470:130–144.https://doi.org/10.1016/j.brainres.2012.06.020
-
g:Profiler-a web server for functional interpretation of gene lists (2016 update)Nucleic Acids Research 44:W83–W89.https://doi.org/10.1093/nar/gkw199
-
Comparative primate neurobiology and the evolution of brain language systemsCurrent Opinion in Neurobiology 28:10–14.https://doi.org/10.1016/j.conb.2014.04.002
-
Modeling gene expression evolution with an extended ornstein-uhlenbeck process accounting for within-species variationMolecular Biology and Evolution 31:201–211.https://doi.org/10.1093/molbev/mst190
-
The neighbor-joining method: a new method for reconstructing phylogenetic treesMolecular Biology and Evolution 4:406–425.https://doi.org/10.1093/oxfordjournals.molbev.a040454
-
Busco: Assessing genome assembly and annotation completenessMethods in Molecular Biology 1962:227–245.https://doi.org/10.1007/978-1-4939-9173-0_14
-
Ndufs4 ablation decreases synaptophysin expression in hippocampusScientific Reports 11:10969.https://doi.org/10.1038/s41598-021-90127-4
-
The evolution of mammalian brain sizeScience Advances 7:eabe2101.https://doi.org/10.1126/sciadv.abe2101
-
Evolution of the human nervous system functionStructure, and Development. Cell 170:226–247.https://doi.org/10.1016/j.cell.2017.06.036
-
Methodische studien uber den quantitativen vergleicharchitektonischer struktureinheiten des gehirnsZeitschrift Fur Wissenschartliche Zoologie 164:143–172.
Article and author information
Author details
Funding
National Science Foundation (BCS-1750377)
- Courtney C Babbitt
National Institutes of Health (T32 GM135096)
- Katherine Rickelton
James S. McDonnell Foundation
https://doi.org/10.37717/220020293- Chet C Sherwood
National Institutes of Health (NS-092988)
- Chet C Sherwood
National Science Foundation (SMA-1542848)
- Chet C Sherwood
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
We would also like to acknowledge our funding from NSF BCS-1750377, the Wenner Gren Foundation, James S McDonnell Foundation (220020293), NSF INSPIRE (SMA-1542848), NIH (NS-092988), and a fellowship to KR from NIH T32 GM135096. https://doi.org/10.37717/220020293.
Copyright
© 2024, Rickelton et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics
-
- 290
- downloads
-
- 2
- citations
Views, downloads and citations are aggregated across all versions of this paper published by eLife.
Download links
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Further reading
-
- Chromosomes and Gene Expression
- Developmental Biology
About 70% of human cleavage stage embryos show chromosomal mosaicism, falling to 20% in blastocysts. Chromosomally mosaic human blastocysts can implant and lead to healthy new-borns with normal karyotypes. Studies in mouse embryos and human gastruloids showed that aneuploid cells are eliminated from the epiblast by p53-mediated apoptosis while being tolerated in the trophectoderm. These observations suggest a selective loss of aneuploid cells from human embryos, but the underlying mechanisms are not yet fully understood. Here, we investigated the cellular consequences of aneuploidy in a total of 125 human blastocysts. RNA-sequencing of trophectoderm cells showed activated p53 pathway and apoptosis proportionate to the level of chromosomal imbalance. Immunostaining corroborated that aneuploidy triggers proteotoxic stress, autophagy, p53-signaling, and apoptosis independent from DNA damage. Total cell numbers were lower in aneuploid embryos, due to a decline both in trophectoderm and in epiblast/primitive endoderm cell numbers. While lower cell numbers in trophectoderm may be attributed to apoptosis, aneuploidy impaired the second lineage segregation, particularly primitive endoderm formation. This might be reinforced by retention of NANOG. Our findings might explain why fully aneuploid embryos fail to further develop and we hypothesize that the same mechanisms lead to the removal of aneuploid cells from mosaic embryos.
-
- Chromosomes and Gene Expression
- Developmental Biology
Transcription often occurs in bursts as gene promoters switch stochastically between active and inactive states. Enhancers can dictate transcriptional activity in animal development through the modulation of burst frequency, duration, or amplitude. Previous studies observed that different enhancers can achieve a wide range of transcriptional outputs through the same strategies of bursting control. For example, in Berrocal et al., 2020, we showed that despite responding to different transcription factors, all even-skipped enhancers increase transcription by upregulating burst frequency and amplitude while burst duration remains largely constant. These shared bursting strategies suggest that a unified molecular mechanism constraints how enhancers modulate transcriptional output. Alternatively, different enhancers could have converged on the same bursting control strategy because of natural selection favoring one of these particular strategies. To distinguish between these two scenarios, we compared transcriptional bursting between endogenous and ectopic gene expression patterns. Because enhancers act under different regulatory inputs in ectopic patterns, dissimilar bursting control strategies between endogenous and ectopic patterns would suggest that enhancers adapted their bursting strategies to their trans-regulatory environment. Here, we generated ectopic even-skipped transcription patterns in fruit fly embryos and discovered that bursting strategies remain consistent in endogenous and ectopic even-skipped expression. These results provide evidence for a unified molecular mechanism shaping even-skipped bursting strategies and serve as a starting point to uncover the realm of strategies employed by other enhancers.