Human cerebral organoids recapitulate gene expression programs of fetal neocortex development

Abstract

Cerebral organoids-3D cultures of human cerebral tissue derived from pluripotent stem cells-have emerged as models of human cortical development. However, the extent to which in vitro organoid systems recapitulate neural progenitor cell proliferation and neuronal differentiation programs observed in vivo remains unclear. Here we use single-cell RNA sequencing (scRNA-seq) to dissect and compare cell composition and progenitor-to-neuron lineage relationships in human cerebral organoids and fetal neocortex. Covariation network analysis using the fetal neocortex data reveals known and previously unidentified interactions among genes central to neural progenitor proliferation and neuronal differentiation. In the organoid, we detect diverse progenitors and differentiated cell types of neuronal and mesenchymal lineages and identify cells that derived from regions resembling the fetal neocortex. We find that these organoid cortical cells use gene expression programs remarkably similar to those of the fetal tissue to organize into cerebral cortex-like regions. Our comparison of in vivo and in vitro cortical single-cell transcriptomes illuminates the genetic features underlying human cortical development that can be studied in organoid cultures.

Keywords: cerebral organoid; corticogenesis; neocortex; single-cell RNA-seq; stem cells.

Fig. 1.

Deconstructing cell composition in the fetal human neocortex. (A) scRNA-seq was performed on cells from two human neocortex specimens at 12–13 wpc. Schematic shows NPC types (APs, BPs) and neurons (N) enriched in zones within the human neocortex at midneurogenesis. AP, apical progenitor; BP, basal progenitor; CP, cortical plate; iSVZ, inner subventricular zone; oSVZ, outer subventricular zone; VZ, ventricular zone. (B) Heat maps show normalized correlation (Z-score) of single-cell transciptomes from human wpc 12 (light green) and wpc 13 (dark green) cerebral cortex with bulk RNA-seq data from laser-microdissected zones (left, 18) or FACS-purified cell types (right, 19) from the human neocortex at the same developmental stage. (C) Hierarchical clustering of scRNA-seq data reveals cell types in the human fetal cortex. Each row represents a single cell and each column a gene. Genes were discovered using PCA (SI Methods). The maximum correlation to bulk RNA-seq data from germinal zones is shown in the left sidebar. Cell-type assignment is shown on the right sidebar. Expression of genes used to classify APs, BPs, newborn neurons (N1), and maturing neurons (N2 and N3) are show to the right of the cell-type assignment bar. Top GO enrichments are shown above the heat map, with representative genes listed below.

Fig. S1.

PCA of fetal cerebral cortex. (A) PCA on all variable genes and all 226 cells was used to identify genes describing cortical cell populations. Each dot represents a cell that is color-coded in shades of green representing different experiments. (B) Hierarchical clustering and heat map visualization showing the expression of genes that have highest correlation and anticorrleation with PC1 and highest correlation with PC2. Note that the genes correlating with PC2 identified a single endothelial cell. Cells are shown in rows and genes in columns. (C) Hierarchical clustering and heat map visualization showing the discovery of five interneurons based on marker gene expression.

Fig. 2.

Reconstructing lineage relationships in the fetal neocortex. (A) Monocle reveals an AP–BP–neuron lineage that correlates with the zones of the developing neocortex. Cells (circles, colored based on cell type) are arranged in the 2D independent component space based on genes used to classify cells in Fig. 1C. The minimal spanning tree (gray lines) connects cells, with the black line indicating the longest path. (B) Monocle plots with single cells are colored based on the maximum correlation with bulk RNA-seq data from cortical zones (Far Left) or gene expression that distinguishes the lineage transitions (Middle Left, Middle Right, and Far Right). (C) Transcription factor (TF) correlation network during lineage progression. Shown are nodes (TFs) with more than three edges, with each edge reflecting a high correlation (>0.3) between connected TFs.

Fig. S2.

Gene expression signatures controlling fetal corticogenesis. (A) Cell lineage network based on pairwise correlations between cells and using the same gene set used in the Monocle analysis (Fig. 2 A and B). Cells are color-coded based on cell-type classification (Left) or maximum correlation to cortical zones (Right). Arrows and numbers mark topological features of the network that reflect cell proliferation and differentiation events inferred from gene expression signatures. (B) Ordering of scRNA-seq expression data according to the pseudotemporal position along the lineage revealed a continuum of gene expression changes from NPCs to neurons. Genes from the heat map depicted in Fig. 1C are shown. Each row represents a single cell and each column a gene. The cell-type assignment and maximum correlation to cortical zone are shown as sidebars. (C) Schematic illustrating lineage decisions in the fetal cortex, labeled with groups of genes characterizing lineage progression. For example, APs are defined by the expression of group A genes and self-renew through the expression of cell-cycle regulators shown in gene groups E and F.

Fig. 3.

Dissecting cerebral organoid cell composition using scRNA-seq. (A) scRNA-seq was performed on whole organoids dissociated at 33, 35, 37, 41, and 65 d after EB culture and four microdissected regions surrounding single ventricles from two organoids (day 53, r1, r2, ESC-derived; day 58, r3, r4, iPSC-derived). (B) The 30 d iPSC-derived organoid immunostained with proliferation marker MKI67 (magenta), neuronal marker DCX (doublecortin) (green), and DAPI (blue). (Scale bar, 100 μm.) Zoom to ventricle (asterisk) shows MKI67, DCX, and DAPI (Left), and NPC marker PAX6 (magenta), BP marker TBR2/EOMES (green), and mitosis marker phospo-vimentin (pVim; cyan) (Right). (Scale bar, 20 μm.) (C) Images of microdissected cortical regions (r1 and r2 from 53 d ESC-derived organoid; r3 and r4 from 58 d iPSC-derived organoid). Dotted lines show microdissection boundaries; asterisks mark ventricles. (D) PCA and unbiased clustering using t-SNE reveals cell populations within organoids. Shapes indicate experiments, and colors represent significant clusters. See Dissecting Cell Composition in Human Cerebral Organoids and SI Results for cluster descriptions. (E) Marker genes for each cluster. Cells are colored based on expression level. Cerebral cortex cells (c1, c2, c3, and c4) have high expression of FOXG1 and low expression of OTX2. (F) Violin plot shows FOXG1, NEUROD6, or OTX2 expression from each microdissected region compared with fetal cortex.

Fig. S3.

Supporting data for dissecting cerebral organoid cell composition using scRNA-seq. (A) Images from various stages of cerebral organoid development. From left to right, iPSC colonies (d0), EB (d6), neural induction (d10), 2 d after matrigel embedding (d12), 4 d after matrigel embedding (d14), and early organoids after 10 d in cerebral organoid media containing vitamin A (d24). For d12 and d14, Insets show stratified epithelium surrounding a ventricle. Boxes indicate areas shown at greater magnification, as indicated. (Scale bar, 100 μm.) (B and C) Cortical region from 35 d organoid immunostained for PAX6 and TBR2 (B) and DCX and DAPI (C). (D) AP dividing perpendicular to the apical (ventricular surface) immunostained for the ciliary marker ARL13B (magenta) and the centrosomal marker γ-Tubulin (TUBBG; green); nuclei were counterstained with DAPI (blue). (E) Violin plots show distribution of the transcript levels of the representative marker genes across all experiments. In situ hybridizations in E11.5 (Foxg1, Otx2, and Rspo2) and E13.5 (Dcn, Myt1l, and Neurod6) mouse brains show the patterns of gene expression along the rostro–caudal or apical–basal axes. Data are from the Allen Brain Institute Gene Expression Atlas (courtesy of the Developing Mouse Brain, Allen Brain Atlas). (F) More neurons were detected per organoid over the time course relative to progenitors. Neuron and progenitor count represents clusters 4 and 7, and clusters 1, 2, and 5 from B, respectively. Representative image of a 30- and 40-d iPSC-derived organoid showing the relative increase in DCX+ neurons during this time window. (Scale bar, 100 μm.) Inset shows progenitor and neuron stratification within an individual neurogenic zone. (Scale bar, 20 μm.) (G) Immunohistochemistry and transmission electron microscopy validation that mesenchymal cells secrete collagen and form a fibrous ECM surrounding the periphery of individual cortical regions. Black arrow points to a clathrin-coated pit, and the white triangles point to collagen fibrils. (EM image scale, 1 μm.) (H) Heat map showing expression of genes that distinguish cell type and regional identity of the different t-SNE clusters. The cell type or region that is marked by a given gene set is shown above the heat map. Cells are in rows and genes in columns (Dataset S2 and SI Methods). (I) Bar plots show the number of cells in each cluster over time for all whole organoid experiments. The color- code of each bar corresponds to the cluster number.

Fig. S4.

Characterizing organoid neuron heterogeneity. (A) Heat map visualization of variable genes describing neuronal heterogeneity in cerebral organoids. PCA on neuron clusters 4 (orange) and 7 (blue) from Fig. 3D was used to discover genes describing neuronal variation, and BackSPIN was used to cluster genes and cells to define different neuronal cell types. There were four major clusters of neurons (1a, 1b, 2a, and 2b), which are subpopulations of cluster 4 (1a and 1b) and 7 (2a and 2b), respectively. Note that four cells in cluster 2a coexpress GAD1, GAD2, SOX14, and ERBB4 and are likely interneurons derived from ventral telencephalon-like regions of the organoid. (B) Violin plots showing the distribution of representative marker gene expression for each neuron cluster. (C) Expression of representative markers from bulk RNAseq data along a time course in multiple regions of the developing human brain (sourced from the Allen Brain Atlas). (D) Gene expression averaged across all neuronal cells from each microdissected region shows that r1, r2, and r3 all express NFIA, NEUROD6, and the other cluster 4 markers, whereas r4 does not. These data support that r1, r2, and r3 contain neurons from organoid dorsal cortex-like regions.

Fig. S5.

Extended data showing similar gene expression profiles characterize lineage progression in organoid and fetal cerebral cortex. (A) Organoid cerebral cortex-like cells (c1, c2, c3, and c4 from Fig. 3) have differential correlation with bulk RNA-seq data from different laser-microdissected zones or FACS-purified cell types from fetal cerebral cortex (Fig. 1B). (B) Organoid cell lineage network based on pairwise correlations between cells. Cells are colored based on maximum correlation to cortical zones (Left) or cell type (Right). (C) Pseudotemporal cell ordering along the organoid lineage reveals gene expression changes from NPC to neuron. Genes with highest correlation and anticorrelation to PC1–3 are shown. Rows represent cells and columns genes. Maximum correlation to cortical zone and cell type is shown in the left sidebars. Top GO enrichments are shown above the heat map, with representative genes listed below. (D) Correlation network using the same TFs as in the fetal TF network reveals two highly connected subnetworks controlling AP proliferation/self-renewal and neuron differentiation. Shown are nodes (TFs) with more than two edges, with each edge reflecting a correlation between connected TFs that is greater than 0.3. (E) Bar plot shows the fraction of fetal and organoid cells that are assigned as AP, BP, or neuron. (F) The maximum and minimum correlation of each organoid cerebral cortex-like cell with any cell from the fetal neocortex plotted along the organoid pseudotime. The correlation is relatively even throughout the lineage. (G) Heat map visualization of fetal expression of the same genes discovered by PCA on organoid cells shown in C. (H) Heat map visualization of organoid expression of the same genes discovered by PCA on fetal cells, which was used in fetal cortex lineage analysis. For comparison, see Fig. 1C and Fig. S2B.

Fig. 4.

Similar gene expression profiles characterize lineage progression in organoid and fetal cerebral cortex. (A) Organoid AP–BP–neuron lineage. PCA on organoid dorsal cortex cells identified genes used for Monocle. The minimal spanning tree (gray lines) connects cells (circles, colored by cell type). Black line indicates the longest path (Dataset S3 and Fig. S5). (B) Monocle plot with cells colored by maximum correlation with bulk RNA-seq data from germinal zones (Far Left) or marker gene expression (Middle Left, Middle Right, and Far Right). (C) Scatter plot shows correlation (Pearson) between fetal and organoid average expression per cell type for marker TFs. (D) Heat map shows TF expression in organoid and fetal cells ordered by pseudotime. Top bars show cell type and maximum correlation with germinal zones. Each TF’s expression was averaged across cells of a given type (AP, BP, N1, N2, and N3), and the Pearson correlation between fetal and organoid cell types is shown to the right of the heat map. (E) Cell lineage network and dendrogram (Top Right) based on pairwise correlations between fetal (orange) and organoid (red) cells show that NPCs and neurons intermix.

Fig. S6.

Cerebral organoids recapitulate gene expression programs underlying cerebral cortex cell biology. (A) APs express genes involved in ECM production and sensing. The heat map shows expression for each ECM gene for organoid and fetal tissue cells in the order of the monocle lineage. For each gene, the correlation between organoid and fetal tissue lineages is plotted to the right of the organoid heat map. Note that the organoid dataset has proportionally more progenitors than the fetal tissue dataset, whereas the fetal tissue has more neurons. (B–D) A similar comparison as in A is shown for genes involved in aRG delamination (B), Delta/Notch signaling (C), and neurite outgrowth (D).

Fig. 5.

Genomic scans of disease, evolutionary, and chromatin signatures highlight genetic aspects of human corticogenesis that can be modeled in vitro. Shown is the covariation network using genes that have high correlation (>0.3) with TFs controlling the AP–BP–neuron lineage from Fig. 2C. Select TF nodes are highlighted to delineate the path. (i–iv) Panels show genes that have amino acid changes that are modHum (i, green), OMIM (ii, turquoise), hCondel (iii, blue), or haDHS (iv, orange). The percentage of cells that have a positive correlation (>0.4) between fetal and organoid cells is shown, with nodes colored based on the correlation coefficient (Dataset S4).

Fig. S7.

Gene expression differences between fetal and organoid cell types. (A and B) Volcano plots visualizing differential gene expression analysis between organoid and fetal APs (A) and organoid and fetal neurons (Ns; B). For each gene, the average difference between cells from each comparison is plotted against the power to discriminate between groups. Top-scoring genes for both metrics are labeled and genes in the top GO enrichment category “Response to organic substance” are highlighted in orange (Dataset S4). (C) Violin plots for some of the most differentially expressed genes show the distribution of gene expression across APs and neurons from fetal, ESC-derived organoid, and iPSC-derived organoid dorsal cortex. (D) Box plots show the distributions of the power of differentially expressed genes to discriminate between different cell populations. Shown in the two left-most box plots are comparisons between APs in organoid versus fetal tissue and neurons in organoid versus fetal tissue. As a reference, box plots are shown for comparisons between iPSC-derived and ESC-derived organoid APs and neurons, respectively. In addition, comparisons are shown for two randomly selected subgroups within fetal APs or fetal neurons, respectively.