Variable Outcomes in Neural Differentiation of Human PSCs Arise from Intrinsic Differences in Developmental Signaling Pathways

Abstract

Directed differentiation of human pluripotent stem cells varies in specificity and efficiency. Stochastic, genetic, intracellular, and environmental factors affect maintenance of pluripotency and differentiation into early embryonic lineages. However, factors affecting variation in in vitro differentiation to defined cell types are not well understood. To address this, we focused on a well-established differentiation process to cerebral cortex neural progenitor cells and their neuronal progeny from human pluripotent stem cells. Analysis of 162 differentiation outcomes of 61 stem cell lines derived from 37 individuals showed that most variation occurs along gene expression axes reflecting dorsoventral and rostrocaudal spatial expression during in vivo brain development. Line-independent and line-dependent variations occur, with the latter driven largely by differences in endogenous Wnt signaling activity. Tuning Wnt signaling during a specific phase early in the differentiation process reduces variability, demonstrating that cell-line/genome-specific differentiation outcome biases can be corrected by controlling extracellular signaling.

Keywords: Hh signalling; Wnt signalling; cortical differentiation; dual SMAD inhibition; hPSC variation; human puripotent stem cell; neural induction; patterning of the cortex; regional identity specification; transcriptional profiling.

Graphical abstract

Figure 1

Gene Expression Profiling in 84 Directed Differentiations Highlights Broad Transcriptional Similarity and Specific Differences in Expression of Regional Brain Genes (A) Protocol used to differentiate cortical cultures from PSCs. The early and late stages analyzed are highlighted. (B) Hierarchical clustering of gene expression from 84 early-stage differentiations profiled with Codeset1. Clusters are named early cluster 1 (EC1)–EC5. Highly expressed cortical development genes are indicated with white arrowheads. Variation was observed in expression of transcripts specific to the telencephalon (FOXG1), the ventral telencephalon (LHX8, LHX6, NKX2-1, DLX1, and DLX5), the hindbrain (HOXA2 and HOXB2), and the dorsal telencephalon (cortex) (EMX1, EMX2, and EOMES), indicated with black arrowheads. (C) Replicating the patterns observed in (B), genes associated with specific brain regions are highly variable across differentiations in a second independent dataset of 65 early-stage differentiations profiled with Codeset2. See also Figure S1.

Figure 2

Analysis of Variation in Gene Expression at Early Stage Reveals Differences in Spatial Identity Reflecting In Vivo Dorsoventral and Rostrocaudal Axes (A) Principal-component analysis (PCA) of early-stage gene expression data (84 samples, 174 genes); samples are plotted along the two components explaining most gene expression variation and colored by cluster. Caudal outliers are circled. EC, early cluster. (B) Gene contributions to variation in the early-stage dataset plotted using the Z-scored loadings for PC1 and PC2 as coordinates. Highest contributors (absolute Z scores >2) are labeled. (C) The distribution of highest loading genes in (B) is correlated with the gene expression patterns along the dorsoventral and rostrocaudal axes in vivo. (D) Expression of selected high-loading genes along PC1–PC2. (E) Mapping to 5 regions of the E11.5 Allen Developing Mouse Brain Atlas based on correlated expression of variable genes. RSP, rostral secondary prosencephalon. (F) Expression heatmap of selected high-loading genes in individual differentiations. See also Figure S2.

Figure 3

Developmental Patterning Axes Are the Main Drivers of Variation in Late-Stage Differentiations (A) PCA of late-stage gene expression data (44 samples, 171 genes) with samples colored by cluster. LC, late cluster. (B) Gene contributions to variation in the late-stage dataset plotted using the Z-scored loadings for PC1 and PC2 as coordinates. Highest contributors (absolute Z scores >2) are labeled. (C) The distribution of highest loading genes in (B) is correlated with the gene expression patterns along the dorsoventral and rostrocaudal axes in vivo and resembles the first two principal components in the early-stage data. (D) Expression of selected high-loading genes along PC1–PC2. (E) Mapping to 5 regions of the E15.5 Allen Developing Mouse Brain Atlas based on correlated expression of variable genes. RSP, rostral secondary prosencephalon. (F) Expression heatmap of selected high-loading genes in individual differentiations and interpretation summary showing variable contribution of different regional identities to composition of differentiations in three late-stage clusters. See also Figure S3.

Figure 4

Early-Stage Gene Expression Is Predictive of Late-Stage Expression and Differentiation Outcomes (A) Heatmap of Pearson correlation coefficients between expression of late-stage cluster markers and early-stage gene expression in 44 time pairs. Only early-stage genes with an absolute correlation coefficient greater than 0.75 for at least one of the late-stage markers are included. (B) Correspondence between early-stage and late-stage cluster assignmentS mapped for 44 differentiations when two, three, or four clusters are considered. (C) Late-stage cluster outcome was used to establish acceptable thresholds of expression of predictive genes at the early stage: low DLX5 and FOXG1 expression were used to classify differentiations as partially caudalized. For the remaining rostral differentiations, high NKX2-1 expression and low PAX6 expression were used to classify differentiations as highly ventralized, while samples with high PAX6 and high NKX2-1 expression were classified as partially ventralized. Correspondence of empirical classification at early-stage with late-stage cluster classification is mapped on the right. See also Figures S4 and S5.

Figure 5

Cell-Line-Specific Variation in Spatial Identity of Neural Differentiations (A) Classification of differentiations from early-stage data as dorsalized (green), partially ventralized (orange), highly ventralized (red), and partially caudalized (blue) using thresholds determined in Figure 4C. Additionally, samples with extremely low DLX5 and FOXG1 expression (belonging to EC5) were classified as highly caudalized (purple). (B) Differentiation outcome frequency plotted by PSC line for lines with at least two differentiations. Names are formatted as individual.clone; numbers indicate separate differentiations per line. Asterisks indicate lines for which the Wilson 95% confidence interval of the difference in dorsal differentiation frequency compared to overall frequency does not include zero. (C) Differentiation outcome frequency plotted by genotypes for which at least three cell lines were induced. MAPT Ex10+16, frontotemporal-dementia-causing exon 10 splicing mutation in gene encoding tau protein; TS21, trisomy 21 (Down syndrome); AD, Alzheimer’s disease. Frequencies were normalized to account for the variable number of differentiations per line. See also Figure S6.

Figure 6

Hedgehog and Wnt Signaling Display Different Temporal Dynamics in PSC Lines with Inherent Tendencies to More Cortical or Ventral Differentiation Outcomes (A) Early-stage differentiations classified as highly or partially ventralized express higher levels of SHH compared to dorsalized and caudalized differentiations, as well as levels of Hedgehog signaling readouts PTCH1, GLI1, and GAS1 consistent with increased pathway activity. Dorsalized differentiations express higher levels of Wnt signaling readouts AXIN2 and TNFRSF19 compared to ventralized differentiations but lower AXIN2 compared to partially caudalized differentiations (pairwise Welch’s t test, false discovery rate [FDR]-corrected p values: ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001; ns, not significant). Only comparisons to the dorsalized class are shown. Differentiations per group: dorsalized, 86; highly ventralized, 23; partially ventralized, 23; and partially caudalized, 14. (B) Early-stage differentiations from separate PSC lines vary in average expression of Hedgehog and Wnt signaling pathway activation, consistent with different tendencies in regional patterning (one-way ANOVA: ^∗p < 0.05; ^∗∗p < 0.01). Error bars represent standard error; n = 2–6 differentiations per line. (C) Gene expression time course during differentiation for selected genes associated with forebrain regions and Hedgehog or Wnt signaling. Differentiations from ventral-prone lines (partially ventral iPSC22.1 and highly ventral GMESC01.1, iPSC01.1, iPSC06.1, and iPSC14.1, n = 8–11) were compared to differentiations of a dorsal-prone line (iPSC21.1, n = 1–2). Profiles represent average gene expression, and error bars represent standard deviation. Significance shown for dorsal versus ventral comparison at 17 and ∼35 dpi (Welch’s t test, FDR-corrected p values: ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001).

Figure 7

Inherent Tendency to Ventralization Is Largely Rescued by a Brief, Specific Phase of Wnt Signaling Activation (A) Diagram of Hedgehog pathway components targeted by small molecules in (B) and (C). (B) Treatment with purmorphamine (1 μM) between 7 and 17 dpi results in a more ventralized gene expression profile at 33 dpi (line iPSC22.1, n = 2). (C) Treatment with Hedgehog inhibitor cyclopamine (1 μM) between 7 and 17 dpi has no observable effect on dorso-ventral gene expression in a highly ventralized line at ∼35 dpi (iPSC14.1, n = 3). (D) Summary diagram of Wnt pathway components targeted by small molecules in (E)–(G). (E) Treatment with Wnt inhibitor IWP2 (2 μM) between 0 and 12 dpi results in a more ventralized gene expression profile at 33 dpi (line iPSC22.1, n = 2). (F) Treatment with Wnt activator CHIR99021 (1 μM) between 13 and 17 dpi significantly increases cortex-associated gene expression and decreases MGE-associated expression at ∼35 dpi in differentiations of 4 ventral-prone lines compared to vehicle treatment (GMESC01.1, iPSC01.1, iPSC06.1, and iPSC22.1) (one-sample Student’s t test, mu = 0, n = 7, FDR-corrected p values: ^∗p < 0.05; ^∗∗p < 0.01; ^∗∗∗p < 0.001). Top panels show combined trend; bottom panels show breakdown by PSC line. (G) Clustering of ∼35-dpi differentiations from ventral-prone lines treated between 7 and 17 dpi with either vehicle or 1 μM CHIR99021. Treatment stimulating Wnt/β-catenin signaling results in shift in classificationof ventralized differentiations to more dorsalized clusters. (H) Model for outcome of differentiation of distinct cell lines. All error bars represent standard error. See also Figure S7.