Cell-Type-Specific Impact of Glucocorticoid Receptor Activation on the Developing Brain: A Cerebral Organoid Study

Human developmental trajectories are influenced by both genetics and the environment (1). An important environmental factor is activation of the glucocorticoid (GC) system during critical brain development periods, which is associated with cognitive, behavioral, and psychiatric health outcomes (1, 2). The glucocorticoid receptor (GR; NR3C1 gene) mediates key processes in fetal organ development, and its activation is essential for proper formation and maturation of many organs (3). Upon activation, GR translocates to the nucleus, where it functions as a transcription factor by binding to GC response elements (4) in many genes, as well as through protein-protein interactions with other transcription factors (5).

GCs are critical for development, and thus alterations profoundly affect fetal outcome. Synthetic GCs readily cross the placenta (3, 6) and are prescribed throughout pregnancy (7), for example, for maternal autoimmune conditions or fetal congenital adrenal hyperplasia (8). Higher doses are given at certain times, such as during postconception weeks 0–12 to prevent miscarriages (9). Most commonly, synthetic GCs are given to promote organ maturation as early as postconception week 22 in pregnancies at risk of premature delivery (10). Synthetic GC treatment has clear benefits for pregnancy and birth outcomes, but negative long-term outcomes are increasingly being reported (11). Accumulating evidence from human cohorts indicates that fetal synthetic GC exposure in mid to late pregnancy may result in adverse postnatal outcomes such as dysregulated hypothalamic-pituitary-adrenal axis activity and altered brain structure and development (12). While earlier studies (N=∼2,000) did not find major long-term effects of antenatal synthetic GCs on neurobehavioral outcomes (13), more recent findings from larger population-based cohorts reveal a significant increase in neurocognitive and behavioral risk later in life. For example, a Finnish study of children (N>600,000) followed up to age 8 (14) found increased prevalence of mental and behavioral disorders in exposed children. A Canadian study (N>500,000) found increases in neurocognitive disorders at age 5 in those exposed to antenatal synthetic GCs (15). Extended synthetic GC treatment for congenital adrenal hyperplasia has been linked to altered cognitive performance, including inattention, and increased fearfulness in childhood (8). Furthermore, elevated endogenous maternal GCs can affect brain development (2, 16–18), for example, with altered neonatal amygdala connectivity, differences in sensory processing and integration, and subsequent internalizing symptoms in girls (19). Together with the increasing prevalence of prenatal synthetic GC administration, this has a large societal impact and underlines the necessity of better understanding how prenatal GR activation outside the normal balance in the developing brain affects health outcomes. For this, model systems are necessary.

GR activation can be modeled robustly in vitro and in vivo and is thus amenable to system-wide investigation. Rodent studies have provided important insights into the role of GCs during development, including that synthetic GC exposure during gestation can affect neurodevelopment (20), and additionally can lead to depression-like phenotypes and other negative outcomes in adult offspring (21). Although important, these models have limitations. During brain development, important differences exist between species in the abundance, lineage complexity, and proliferative potential of certain neural progenitor cell types such as basal radial glia. Additionally, GR response elements and enhancers are not always conserved. Human induced pluripotent stem cells (hiPSCs) represent a promising tool for filling this gap. In addition to having human genetic background, hiPSC-derived cerebral organoids recapitulate the three-dimensional complexity and cell-to-cell interactions (22, 23) of fetal brain development during early to mid gestation (24–26), a period critical for neuron production, migration, connection, and differentiation (1).

Environmental impacts during these sensitive periods have rarely been studied. We aimed to test whether cerebral organoids are suitable for studying the impact of prenatal GC exposure on the developing brain. Here, GR activation with dexamethasone (Dex) in organoids led to cell-type-specific transcriptional responses mapping to genes that moderate risk for negative neurodevelopmental trajectories and psychiatric and behavioral traits, especially in neuronal cells.

Methods

Three hiPSC lines were used to generate cerebral organoids (27), starting with 9,000 cells. Line1 cells were reprogrammed from NuFF3-RQ male human newborn foreskin feeder fibroblasts (GSC-3404, GlobalStem) (28). Line2 cells were generated using a plasmid-based protocol for integration-free hiPSCs from peripheral blood mononuclear cells (29) from a female donor through the BeCOME study (30). Line3 hiPSCs were reprogrammed skin fibroblasts (HPS0076:409b2, RIKEN BRC cell bank, female) (31, 32). All donors gave informed consent.

Downstream analyses used the Scanpy package (33) unless stated otherwise, and MAST (34) was used for all single-cell differential expression analyses. To test for enrichment of selected gene sets in the differentially expressed genes (DEGs; false discovery rate–corrected p≤0.05), we used hypergeometric tests or permutations based on equal-sized gene sets with either the same mean expression distribution as the significant DEGs from each cell class or random sampling.

The study methods are described in more detail in the Supplementary Methods section of the online supplement.

Results

Cerebral Organoids Model the Early Developing Brain In Vitro and Express the Molecular Machinery for Glucocorticoid Response

Organoid transcriptomes were explored using bulk RNA sequencing (RNAseq) at seven successive developmental times from day 17 to day 158 (see Table S1 in the online supplement). We delineated a robust developmental trajectory at the whole transcriptome level, with age-defining genes clearly clustering in a progressive temporal manner (Figure 1A). Established cell type markers displayed predictable developmental patterns at the RNA (see Figure S1A,B and Table S2 in the online supplement) and protein level (see Figure S1C,D in the online supplement). We compared this data set with other RNAseq data sets including fetal and adult postmortem brain (35, 36) and hiPSC-derived in vitro differentiation models. We found that organoids older than 70 days clustered with early fetal brain samples (Figure 1B).

Acute Glucocorticoid Exposure in Cerebral Organoids Elicits Differential Expression of Established Glucocorticoid-Responsive Genes

We found the GR (NR3C1 gene) to be present at the protein level in a fraction of cells throughout the ventricle structure of organoids, from basal to apical, as well as outside the ventricle in stromal cells. Thus, GR is neither exclusive to any cell type nor ubiquitously expressed (Figure 1C). To determine the availability of the molecular machinery for GR activation in organoids, we tested expression of chaperones and co-chaperones of the GR protein complex (37) at the RNA level (Figure 1D; see also Table S1 in the online supplement). NR3C1 expression increased with organoid age, significantly between day 23 and day 40 (p<0.0001; fold change [FC]=1.29), plateauing after day 40 (Figure 1D). However, since cell type composition is also dynamic with age, we cannot conclude that there is an age effect on gene expression levels, as opposed to a change in the abundance of gene-expressing cells. Genes encoding the GR protein complex were abundantly and stably expressed across organoid development (Figure 1D).

To test GR activation dynamics, we used the synthetic agonist Dex. Organoids were grown until day 45, that is, after GR expression plateau, and tested with two acute paradigms (4 hours and 12 hours), with three different medically relevant Dex concentrations (10 nM, 100 nM, and 1000 nM) and vehicle control (dimethyl sulfoxide). Neither progenitor (SOX2, PAX6) nor neuronal cell markers (TUBB3, MAP2), nor the GR itself (NR3C1) were significantly Dex reactive at the RNA level (Figure 1E,F; see also Table S3 in the online supplement). This was contrary to transcripts known to be up-regulated by GCs: FKBP5 (38), SGK1 (39), TSC22D3 (37), and ZBTB16 (40). Most of the concentration-time combinations tested resulted in significant up-regulation of these genes (Figure 1E,F; see also Table S3). The strongest response was elicited by 12-hour stimulation (Figure 1F), with a 39-fold change for TSC22D3 and a 20-fold change for FKBP5 (Table S3). This combination (100 nM Dex/12 hours) was chosen for subsequent experiments.

Single-Cell RNA and Protein Characterization Shows Non-Cell-Type-Specific Expression and GR Activation in Organoids

We next aimed to understand cell-type-specific GC responses. We performed single-cell transcriptome sequencing in hiPSC Line1 at three time points with relevant GR expression: at day 30, day 60, or day 90 (N=4 each) in organoids exposed to Dex (100 nM/12 hours). After quality control procedures, we profiled 14,002 cells across time points and treatment conditions (day 30, N=5,035 cells; day 60, N=4,315; day 90, N=4,652). We applied graph-based cluster analysis using Scanpy on all Line1 organoid cells and identified 17 fine cell types (Figure 2A) and three coarse cell classes (Figure 2B) based on established marker genes (see the Supplementary Results section, Table S4, and Figures S2 and S3A,B in the online supplement) (22, 23, 41). When separating cells by organoid age, the same cell types or cell classes were identified and represented by the same markers, with the main difference being relative cell abundance indicative of advancing neuronal maturation (see Supplementary Results and Figure S3C,D in the online supplement).

Treated and untreated organoids yielded comparable cell numbers and distribution across organoids and cell types (see Figure S3E in the online supplement). NR3C1 (Figure 2C) had detectable expression in 15% of all cells. High Dex concentrations have been shown also to induce activation of the mineralocorticoid receptor (NR3C2 gene). Given that NR3C2 (Figure 2D) is detected in ≤2% of cells as opposed to 15%−34% of cells expressing NR3C1, consistent with other organoid (22, 42, 43) or fetal brain (42, 44) data, it is more likely that the majority of Dex-induced transcriptional responses are conferred by GR activation. This difference between the two glucocorticoid receptors was consistent when separating cells by organoid age (see Figure S4D–I in the online supplement). NR3C1 was distributed across cell types, although on average more abundantly in nonneural progenitors. Specifically, in the ES, S, and M clusters, >50% of cells expressed NR3C1 (Figure 2E). NR3C1 expression at specific organoid ages was not significantly different, either in all cells (see Figure S4J in the online supplement) or in individual clusters (i.e., mesenchymal stromal cells or radial glia; see Figure S4H,I in the online supplement).

For protein-level evidence of GR activation (37) in single cells, we applied 100 nM Dex/12-hour exposure in day-60 organoids and queried cellular localization using immunofluorescence and confocal microscopy. In the vehicle condition, GR was preferentially localized in the cytoplasm but also sometimes detectable in the nucleus (Figure 2F,G). Following Dex there was a significant (χ²=100.9, p<0.0001) shift in detectable GR toward the nucleus (Figure 2F,H), suggesting Dex-induced nuclear translocation, and thus GR activation. This was true regardless of cell type, with significant shift to nuclear GR following Dex happening in nonneural progenitors, neural progenitors, and neurons (see Figure S5 in the online supplement).

Cell-Type-Specific Responses to Acute Glucocorticoids Point to a Disruption of Neuronal Processes With Lasting Effects

Given the comparable cell types observed at day 30, day 60, and day 90 (see Figure S4 and Tables S5–S7 in the online supplement), we performed the main differential expression analyses in all Line1 cells combined (referred to as “Alldays”) to maximize statistical power. In the 17-cluster annotation, we found significant DEGs in nine clusters, ranging from one to 55 DEGs (q≤0.05; see Figure S6 and Table S8 in the online supplement). More detail on specific genes and dysregulation patterns is presented in the Supplementary Results section in the online supplement.

When analyzing differential expression within the three broader cell classes, we found 899 DEGs in nonneural progenitors, 407 DEGs in neural progenitors, and 322 DEGs in neurons (q≤0.05) (Figure 3A; see also Table S8 in the online supplement). The number of DEGs increased with increasing NR3C1 expression, and fold changes were higher in NR3C1-positive cells (Figure 3B; see also the Supplementary Results section in the online supplement). Forty-eight DEGs were significant in all three cell classes (Figure 3C), of which 33% were down-regulated, including neuronal genes FOXG1 and NEUROD6, and 67% were up-regulated, including Dex-responsive genes FKBP5 and TSC22D3 and the progenitor gene MGARP. Gene ontology analyses of these 48 shared genes revealed “forebrain development” and “cerebral cortex development” among the most significantly enriched (see Table S9 in the online supplement). When separating cells by organoid age, the fold change direction for all significant DEGs was always consistent with the direction in Alldays in at least two ages (q≤0.05; see Table S10 in the online supplement), and switched direction extremely rarely (≤2%) only in one age, supporting the decision to combine cells across organoid age. However, few genes were significant across all organoid ages (q≤0.05) reflective of power limitations due to cell abundance.

To test whether this acute stimulation has a sustained effect beyond the immediate GC response, we exposed organoids to 100 nM Dex for 12 hours and collected them up to 8 days after Dex removal. At the bulk RNA level, the effect of this acute GR activation lasted for >6 days (Figure 3D–G; see also Figure S7 and Table S3 in the online supplement). When exploring this effect on cell-type composition at the protein level, after an 8-day washout period, we found an increase in new CTIP2-positive neurons (p=0.0007; FC=2.5; deep layer V marker) and an increase in TBR1 that fell short of significance (p=0.07; FC=2.4; layer VI marker) in the ventricular zone. There was no significant difference in upper layer marker SATB2, as these neurons have a different birth and migration window. These findings are in line with lasting effects of the Dex-related up-regulation of progenitor markers such as PAX6 (see Figure S7) and evidence from mouse studies (45).

Dex-Regulated Transcripts in Neurons Are Enriched for Genes Implicated by GWAS for Behavioral Traits

To better gauge a potential impact of altered GR activation during neurodevelopment on phenotypes later in life, we mapped DEGs to gene lists implicated by over 700 genome-wide association studies (GWASs) in different quantitative traits and diseases recorded in the GWAS Catalog (46) (see Table S11 in the online supplement). We tested enrichment of DEGs from the three cell classes (q≤0.05) within genes significantly associated with GWAS traits using FUMA (47). To increase specificity, we reduced the test set to only those traits with at least five shared genes with our DEG lists. No significant enrichment emerged for nonneural progenitors. Neural progenitor DEGs were enriched for genes genome-wide associated with 38 traits, and neuron DEGs with 22 traits (enrichment q≤0.05) (Figure 4A,B; see also Table S11). GWAS phenotypes not related to brain function were labeled as “other.” Brain function–related traits were labeled “brain/behavior.” The latter were overrepresented in traits enriched for neuron DEGs, making up 59% of the traits with significant enrichment for this cell class (see Table S11). This distribution was significantly different from the full list of tested traits (N=59) (χ²=4.23, p=0.04) (Figure 4C). Traits classified as “brain/behavior” included “depression,” “neuroticism,” “chronotype,” and “adventurousness,” suggesting that GR activation in the developing brain can affect genes relevant for behavioral phenotypes and psychiatric disorders later in life. This was not significant for neural progenitors, where 89% of traits fell in the “other” category (Figure 4D; see also Table S11).

Since prenatal GC exposure has been associated with cognitive and behavioral outcomes (14, 15), we tested DEG enrichment among genes carrying common variants identified by the Cross-Disorder Group of the Psychiatric Genomics Consortium GWAS of eight mental illnesses (48). We found a significant enrichment only for neuron DEGs (q=0.00475 and odds ratio=3.91; see Table S12 in the online supplement). To test whether the enrichment results were reflective of cell-type-specific gene expression abundance, rather than strictly GR activation, we used a permutation-based approach. Gene sets comparable in size and mean expression distribution but not significantly differentially expressed (q≤0.05) were tested. A significant enrichment of these gene sets was observed only in 2/1,000 iterations (corrected empirical p=0.009) (Figure 4E; see also Table S12). This indicates that the enrichment of DEGs among genes in the Cross-Disorder GWAS is specific to the GR responsiveness of transcripts and not their expression level in neurons. This test was not significant in neural progenitors (Figure 4F; see also Table S12). An approach where permuted gene sets were randomly distributed yielded comparable results (see Table S12), further reinforcing the specificity of the DEG enrichment finding in neurons.

GR Activation in Neurons Uniquely Regulates Neurodevelopmental Disease-Related Genes

For neurodevelopmental phenotypes, we used gene lists from the DisGeNET database (49) associated with autism spectrum disorder (ASD), neurodevelopmental disorders, and height (as a non-brain-related control). We found a significant enrichment for ASD genes in DEGs from all three cell classes (Figure 4G–I; see also Table S13 in the online supplement). The most pronounced result emerged for neuron DEGs, where we found a significant enrichment with both ASD and neurodevelopmental disorders but not height (ASD: p=8.89 × 10⁻⁵, odds ratio=3.04; neurodevelopmental disorders: p=3.41 × 10⁻⁶, odds ratio=5.59; height: p=0.83, odds ratio=0.67) (Figure 4I; see also Table S13). Size- and expression-level-matched as well as random permuted gene set enrichments for ASD and neurodevelopmental disorders never reached a p value less than or equal to the nominal p value in the neurons (see Figure S8A,B in the online supplement), while this was not the case for height, making this finding specific to the Dex-regulated genes in this cell class (empirical q=0.009; see Table S13). We further validated these findings with genes carrying loss-of-function mutations associated with intellectual disability from the Developmental Brain Disorders Database (50) (see the Supplemental Results section, Table S14, and Figure S8C,D in the online supplement).

Validation in Day-90 Organoids of Two Additional Genetic Backgrounds Supports Disease Relevance of Neuronal Glucocorticoid Response

Since the most interesting findings emerged in neurons, and this cell class was most abundant at day 90, we generated validation data sets in same-age organoids derived from two additional hiPSC lines. When cells were classified as before (nonneural progenitors, neural progenitors, and neurons) (Figure 5A–C; see also Figure S9 in the online supplement), the proportions of neurons were different across Lines1–3 (Figure 5D), supporting previous reports of heterogeneity of organoids from different donor hiPSCs (24, 43). In neurons (Line1: N=1,831; Line2: N=1,772; Line3 N=2,612) (Figure 5D), we observed 641 DEGs in Line2 and 446 DEGs in Line3 (q≤0.05; see Table S15 in the online supplement), compared with 337 DEGs from the comparable analysis in Line1 at day 90. Of the significant DEGs, 32%, 26%, or 39%, respectively, for Line1, Line2, and Line3 were shared with at least one other line. Thirty-one genes were shared across all three lines, including progenitor staples such as PAX6, FABP7, and CRABP1, and neuronal genes such as C1orf61, GAD2, and NEUROD1. When focusing on the direction of differential expression, we found that 120 of the DEGs that were significant in Line1 (N=337) shared fold-change direction across Lines1–3, three times more than would be expected by chance (binomial distribution p<0.0001).

Despite differences in cell numbers and proportions (Figure 5D), a similar pattern emerged in the pathways and gene sets enriched. When testing for enrichment of DEGs from neurons (q≤0.05) within genes significantly associated with GWAS traits using FUMA, traits classified as “brain/behavior” were overrepresented in significantly enriched traits, representing 82%, 88%, and 74%, respectively, in Lines1–3 (Figure 5E–G; see also Table S16 in the online supplement). This distribution was significantly different from the full list of traits (Line1: χ²=30.97, p<0.0001 [Figure 5E]; Line2: χ²=78.42, p<0.0001 [Figure 5F]; Line3: χ²=35.41, p<0.0001 [Figure 5G]). Many of the shared enriched GWAS traits across cell lines were related to neuropsychiatry, such as “adventurousness,” “irritable mood,” and “feeling worry.” There were 12 shared enriched GWAS terms across at least two of Lines1–3 (see Table S17A in the online supplement), of which 11 (92%) belonged to the “brain/behavior” category, while “irritable mood” was shared across Lines1–3 (see Table S17B). The 12 shared enrichments were driven by 107 DEGs (Table S17C). TCF4, a transcription factor associated with multiple psychiatric disorders, was the strongest contributor, with 25 instances (trait-by-Line) where it drove the enriched GWAS trait (see Table S17C).

Line2 and Line3 also validated the finding of an enrichment in Cross-Disorder GWAS associations (see Table S18 in the online supplement). Finally, we also tested the relationship between DEGs and intellectual disability genes with loss-of-function mutations in the Developmental Brain Disorders Database (49), where all lines showed a similarly strong enrichment among neuron DEGs (Figure 5H; see also Table S18). Random permutations of nonsignificant genes did not reach a p value less than or equal to the nominal p value (see Table S18). Furthermore, these gene set enrichment findings were specific to neurons in all three lines (see Table S18), validating the finding of a unique profile of GC response in this mature cell class.

Discussion

In epidemiological and clinical studies, GR activation via synthetic GC administration during critical periods of human brain development is associated with negative health outcomes, such as endocrine dysfunction and brain structural and functional alterations, and emotional, behavioral, cognitive, or mental health problems (12, 51, 52). We explored whether hiPSC-derived cerebral organoids can improve our mechanistic understanding of in utero GC exposure. We found GR expression and activity across all cell types in this system. Single-cell RNA sequencing showed that transcriptional GC response was cell-type-specific, pointing to interference with neuronal differentiation and maturation processes. Transcriptional response in neurons was significantly enriched for genes associated with behavioral and neurodevelopmental traits, suggesting a possible convergence of environmental and genetic factors on the same cell-type-specific neurodevelopmental pathways. Thus, cerebral organoids open the possibility of interrogating human-specific genetic vulnerability to prenatal exposures.

The most profound GC response was in neurons, with decreased expression of neuronal-specific genes NEUROD6, FOXG1, TBR1, MYT1L, and NFIA but an increase of progenitor genes PAX6 and MGARP across several cell types and time. This confirms cell and animal study findings of GC involvement in neurogenesis by increasing proliferation while decreasing differentiation (39, 53, 54). In addition to acute effects after 12 hours, we showed Dex effects that lasted several days, both at the RNA level, on GR-responsive genes FKBP5 and TSC22D3 and progenitor genes PAX6 and MGARP, and at the tissue level, with an increase in new CTIP2-positive deep layer neurons. This is in line with mouse models treated with a single acute GC dose (45). Thus, this study further supports the evidence of dysregulated GR activation leading to cell-type-specific effects on neuronal maturation and differentiation.

Since genes associated with psychiatric risk are also important players during fetal development (55), we wanted to understand how these genes respond to GC exposure during neurodevelopment. In our data, GC-responsive transcripts were enriched among gene loci carrying genetic variants associated with neurodevelopmental, behavioral, and psychiatric phenotype associations. Uniquely for GR-activated transcripts in neurons, we found an enrichment of genes linked to behavioral phenotype GWAS findings and to neurodevelopmental disorders. Importantly, validation showed that despite cell-type and individual gene-level heterogeneity between organoids derived from three donor hiPSCs, the pathways of GC-responsive transcripts in neurons are reproducibly consistent. Many of the genes driving these convergent enrichments were transcription factors essential for neural cell specification (TCF4, AUTS2, PAX6, RELN, OTX2) and are implicated in neurodevelopment as well as in brain-related disorders (see Table S17B in the online supplement). A limitation of these transcriptome-wide analyses is that the mechanism by which activated GR induces these expression changes (DNA binding or protein-protein interactions with other transcription factors) was not investigated, and the timing of GC exposure allows for both direct and indirect transcriptional responses. Mechanistic investigation of direct versus indirect GR effects would be highly interesting.

Overall, our results highlight the interplay between environmental exposure and inherited genetic factors on neurodevelopment. Organoids allow exploration of the interaction between human-specific genetic variation and environmental exposures in a three-dimensional system. While we showed that organoids constitute a suitable model for in utero GC exposure, they have limits in modeling complex exposures, such as prenatal stress. Nonetheless, this study helps fill knowledge gaps on how prenatal environmental disturbances affect neurodevelopment, leading to negative outcomes. In the future, this system could help screen compounds to counteract GC effects on brain development.