ABSTRACT
Naïve epiblast cells in the embryo and pluripotent stem cells in vitro undergo developmental progression to a formative state competent for lineage specification. During this transition, transcription factors and chromatin are rewired to encode new functional features. Here, we examine the role of mitogen-activated protein kinase (ERK1/2) signalling in pluripotent state transition. We show that a primary consequence of ERK activation in mouse embryonic stem cells is elimination of Nanog, which precipitates breakdown of the naïve state gene regulatory network. Variability in pERK dynamics results in heterogeneous loss of Nanog and metachronous state transition. Knockdown of Nanog allows exit without ERK activation. However, transition to formative pluripotency does not proceed and cells collapse to an indeterminate identity. This outcome is due to failure to maintain expression of the central pluripotency factor Oct4. Thus, during formative transition ERK signalling both dismantles the naïve state and preserves pluripotency. These results illustrate how a single signalling pathway can both initiate and secure transition between cell states.
INTRODUCTION
Understanding how cells traverse between molecularly and functionally distinct identities is a central question in developmental and stem cell biology. Pluripotent mouse embryonic stem cells (ESCs) reside in a naïve state. To undergo germline and somatic differentiation, they must first transition to a formative state that is competent to respond to lineage induction (Hayashi et al., 2011; Kinoshita et al., 2021; Mulas et al., 2017). In the naïve state, ESCs undergo continuous symmetrical self-renewal while retaining pre-implantation epiblast identity and the ability to colonise blastocyst chimaeras (Bradley et al., 1984). Stable propagation of naïve ESCs is achieved using ‘2i’ culture conditions, comprising small molecule inhibition of MEK upstream of ERK1 (also known as MAPK3) and ERK2 (MAPK1) (ERK1/2: henceforth ERK) and partial inhibition of glycogen synthase kinase 3 (GSK3) (Ying et al., 2008). Leukaemia inhibitory factor (LIF) further increases self-renewal efficiency. Withdrawal of 2i and culturing cells in only N2B27 basal media opens the path to differentiation via transition to formative pluripotency (Smith, 2017). Formative pluripotent cells have distinct molecular, cell biological and developmental properties compared with naïve cells. They have different transcription factor dependencies, gene expression signatures, chromatin modifications, DNA methylation, metabolism and cell surface mechanics (De Belly et al., 2021; Hayashi et al., 2011; Kalkan et al., 2017a; Kinoshita et al., 2021; Mulas et al., 2017; Smith, 2017; Wang et al., 2021; Yu et al., 2021). Most importantly, although naïve cells cannot respond productively to lineage-inducing signals, formative cells have gained the competence to undergo efficient germline and germ layer lineage specification (Hayashi et al., 2011; Kinoshita et al., 2021; Mulas et al., 2017). Similar findings have been demonstrated for human naïve pluripotent stem cells (Rostovskaya et al., 2019).
Multiple studies have established that the ERK pathway acting downstream of autocrine FGF4 is a key driver of the pluripotency transition. Perturbations in components of the FGF/ERK pathway prevent or significantly impede differentiation (De Belly et al., 2021; Betschinger et al., 2013; Burdon et al., 1999; Cheng et al., 1998; Findlay et al., 2013; Kunath et al., 2007; Leeb et al., 2014; Molotkov et al., 2017; Sangokoya and Blelloch, 2020; Stavridis et al., 2007; Yang et al., 2012). Inhibition of FGF/ERK thus empowers the robust derivation and propagation of ESCs (Batlle-Morera et al., 2008; Ying et al., 2008). Conversely, increasing ERK activity by modulating negative feedback regulators, such as the ribosomal S6 kinase (RSK) family, leads to a more rapid and synchronised progression from naïve to formative cell states (Nett et al., 2018; Yang et al., 2012). A recent study using optogenetic stimulation of FGFR1 has further confirmed that ERK activity propels naïve state exit (Arekatla et al., 2023). How this is achieved remains unclear, however.
Here, we combine chemical modulation of ERK activity with transcription factor perturbations to examine the requirements for both naïve state exit and progression to formative pluripotency.
RESULTS AND DISCUSSION
Kinetics of naïve state exit depends on activation of ERK
To monitor the effect of ERK activity on exit from the naïve ESC state, we employed two assay systems: (1) colony formation on replating in 2i/LIF (Dunn et al., 2014); (2) expression of the Rex1:GFPd2 (RGd2) fluorescent reporter, specific for the naïve state (Kalkan et al., 2017a). We first examined cells withdrawn from 2i for 30 h in N2B27 basal media with or without the MEK inhibitor [MEK(i)] PD0325901 to block ERK activation. In agreement with previous results (Nett et al., 2018), we saw that, although cells in N2B27 alone lost capacity for ESC colony formation, this capacity was maintained in MEK(i) cultures (Fig. S1A). Cells in MEK(i) also retained expression of the naïve RGd2 reporter (Fig. S1B). To examine the effect of increased ERK signalling, we blocked negative regulation by RSK proteins using the pan-RSK inhibitor [RSK(i)] BI-D1870 (Nett et al., 2018). Following treatment with RSK(i) for 30 h, cells formed fewer ESC colonies (Fig. S1A) and showed more extensive downregulation of the RGd2 reporter than cells kept in N2B27 only (Fig. S1B). Thus, we confirmed that ERK signalling activity regulates the rate of exit from the naïve ESC state (Nett et al., 2018).
Mouse ESCs are sustained by cooperative and partially redundant activity of a suite of transcription factors that act in conjunction with the constitutive pluripotency factors Oct4 (encoded by Pou5f1) and Sox2 (Dunn et al., 2014; Niwa, 2018). Among these, Klf4 and Nanog have been reported to be directly phosphorylated and potentially destabilised by pERK (Dhaliwal et al., 2018; Kim et al., 2014). We used siRNA to assess whether depletion of either factor, or of Esrrb, which lies downstream of Nanog (Festuccia et al., 2012), would allow naïve state exit by 30 h without ERK signal (Fig. 1A). siKlf4 showed no effect on RGd2 expression in the presence of MEK(i). In contrast, treatment with siNanog or siEsrrb resulted in RGd2 downregulation (Fig. 1B; Fig. S1C) and furthermore reduced colony formation on replating (Fig. 1C). Therefore, the ability of MEK(i) to delay naïve state exit requires Nanog or Esrrb. We also examined the consequences of knockdown on naive state exit when ERK signalling is active (Fig. 1D). Both siNanog and siEsrrb resulted in faster downregulation of RGd2 compared with scramble control siRNA or siKlf4 (Fig. 1E). These data indicate that clearance of Nanog and Esrrb is rate-limiting for exit from the naïve ESC-state.
Identifying transcription factors downstream of ERK. (A) Strategy to identify transcription factors downstream of ERK. (B) Percentage of RGd2-positive cells determined by flow cytometry, 30 h after 2i withdrawal. Unpaired, two-tailed Student's t-test with respect of N2B27/MEK(i)+siNeg (N=6). (C) Representative images of colony formation in 2i/LIF conditions of 30 h cells. (D) Strategy to identify rate-limiting steps during exit from naïve pluripotency. (E) Percentage of RGd2-positive cells determined by flow cytometry following siRNA treatment in N2B27 conditions. Unpaired, two-tailed Student's t-test with respect of N2B27/MEK(i)+siNeg (N=6). (F-H) Box plots showing percentage of Nanog- and Esrrb-positive cells over time in N2B27/DMSO under endogenous ERK activation (F), in N2B27/MEK(i) (G) and in N2B27/RSK(i) (H) conditions (N=2, six images/condition). Two-sided Wilcoxon test. Box plots show median values (middle bars) and first to third interquartile ranges (boxes); whiskers indicate 1.5× the interquartile ranges; dots indicate outliers. (I) Percentage of Nanog-positive cells at 20 h under different concentrations of MEK(i) and RSK(i) (N=2, six images/condition). Unpaired, two-tailed Student's t-test with respect to DMSO. (J) Expression of naïve genes by RT-qPCR at 30 h cells under indicated conditions. Expression was normalised to Actb, then to the maximum for each gene (N=3). Data are mean±s.d. **P<0.005, ***P<0.0005. ns, not significant.
Identifying transcription factors downstream of ERK. (A) Strategy to identify transcription factors downstream of ERK. (B) Percentage of RGd2-positive cells determined by flow cytometry, 30 h after 2i withdrawal. Unpaired, two-tailed Student's t-test with respect of N2B27/MEK(i)+siNeg (N=6). (C) Representative images of colony formation in 2i/LIF conditions of 30 h cells. (D) Strategy to identify rate-limiting steps during exit from naïve pluripotency. (E) Percentage of RGd2-positive cells determined by flow cytometry following siRNA treatment in N2B27 conditions. Unpaired, two-tailed Student's t-test with respect of N2B27/MEK(i)+siNeg (N=6). (F-H) Box plots showing percentage of Nanog- and Esrrb-positive cells over time in N2B27/DMSO under endogenous ERK activation (F), in N2B27/MEK(i) (G) and in N2B27/RSK(i) (H) conditions (N=2, six images/condition). Two-sided Wilcoxon test. Box plots show median values (middle bars) and first to third interquartile ranges (boxes); whiskers indicate 1.5× the interquartile ranges; dots indicate outliers. (I) Percentage of Nanog-positive cells at 20 h under different concentrations of MEK(i) and RSK(i) (N=2, six images/condition). Unpaired, two-tailed Student's t-test with respect to DMSO. (J) Expression of naïve genes by RT-qPCR at 30 h cells under indicated conditions. Expression was normalised to Actb, then to the maximum for each gene (N=3). Data are mean±s.d. **P<0.005, ***P<0.0005. ns, not significant.
We examined the effect of ERK perturbations on Nanog and Esrrb protein by immunofluorescence staining. In N2B27 medium alone, downregulation of both proteins occurred between 12 h and 18 h (Fig. 1F; Fig. S1D). MEK(i) delayed loss of Nanog protein for at least 22 h (Fig. S1E) but had a lesser effect on loss of Esrrb, which proceeded in ∼50% of cells from 18 h (Fig. 1G). Therefore, Esrrb downregulation can occur in part independently of ERK or Nanog. Enhancing ERK activity with RSK(i) treatment resulted in earlier downregulation of Nanog protein with little change in Esrrb (Fig. 1H). To determine whether Nanog responsiveness is dose-dependent, we titrated MEK(i) and RSK(i) and measured the percentage of Nanog-positive cells at 20 h. We found that higher levels of ERK activity resulted in a lower number of Nanog-positive cells (Fig. 1I). Treatment with both MEK(i) and RSK(i) phenocopied MEK(i)-only, as expected if the effect of RSK(i) on Nanog is through the ERK pathway (Fig. 1I).
These results suggest that ERK signalling leads directly to depletion of Nanog protein. We also analysed mRNA expression of a panel of naive ESC transcription factors (Fig. 1J). These markers are barely expressed after 30 h in N2B27. In MEK(i) we saw reduction in mRNA for all except Esrrb, which increased slightly (Fig. 1J). When we combined MEK(i) with siNanog, Esrrb mRNA was lost, in line with the observed exit from the naïve state. Transcriptional regulation of Esrrb is multifaceted, including input from Oct4, Sox2 and Tcf3 (Festuccia et al., 2018a; Martello et al., 2012; Whyte et al., 2013). However, Esrrb has been shown to be a direct transcriptional target of Nanog (Festuccia et al., 2012). Persisting Nanog protein in MEK(i) conditions is therefore likely to be the main driver of maintained Esrrb transcription and delay in naïve state exit.
Our findings reveal that during pluripotency transition the proximal effect of acute ERK activation on the naïve transcription factor network is to reduce Nanog protein, which in turn diminishes Esrrb transcription. This sequence of events is predicted to collapse the naïve gene regulatory network and propel exit from the naïve state (Dunn et al., 2014; Festuccia et al., 2018b).
ERK activation dynamics correlates with loss of Nanog protein in individual cells
The temporal dynamics of ERK activation and of Nanog downregulation are heterogenous in ESC populations (Arekatla et al., 2023; Reimann et al., 2023). We investigated whether the events are correlated in single cells. To capture the dynamics of ERK activity over the timescale of cell state transitions, we monitored expression of a direct transcriptional target, Spry4. We employed ultrasensitive bioluminescence live imaging, which has low phototoxicity, and used luciferase reporters with fast folding and degradation kinetics (Mandic et al., 2017). Using ESCs carrying a Spry4-Fluc transcriptional reporter (Phillips et al., 2019) and a Nanog::Nluc protein fusion knock-in, we imaged at 10 min intervals. Withdrawal of 2i/LIF resulted in dynamic expression of the Spry4-Fluc reporter, consistent with ERK activation profiles observed in ESCs and other cells (Fig. 2A,B; Figs S2A,B and S6) (Albeck et al., 2013; Aoki et al., 2013; Nett et al., 2018). Of interest, we observed only one Spry4-Fluc peak in most cell division cycles.
Live imaging of ERK activity and Nanog protein dynamics. (A) Representative time series images of Spry4-Fluc, Nanog::Nluc signal after withdrawal of 2i culture conditions. Red, orange and yellow circles show different cells across selected frames. (B) Representative trace of Nanog::Nluc and Spry4-Fluc over ∼38 h of differentiation. (C) Cross-correlation function plot between Spry4 promoter activity and Nanog protein for the combined data (all traces joined end-to-end) in N2B27/DMSO, alongside controls in which we measured the cross-correlation between Nanog and a randomised Spry4 signal (rand. signal) or between Nanog and Spry4 traces time shifted by a randomised time (rand. time). (D) Comparison of maximum cross-correlation between Spry4 promoter and Nanog protein in individual traces against the ‘randomised signal’ control in N2B27/DMSO conditions. Black shows mean±s.d. Wilcoxon rank sum test (or unpaired two-samples Wilcoxon test) (n=44, two independent imaging experiments). (E) Comparison of peak lag time at maximum cross-correlation for all the cells that show a significant maximum cross-correlation above the ‘randomised signal’ control. Data values were compared with the lag time observed between Nanog and Spry4 traces time shifted by a random time (‘randomised time’). Black shows mean±s.d. Statistics show comparison of the distributions using Kolmogorov–Smirnov (n=28, two independent imaging experiments). (F) Autocorrelation function plot between Spry4 promoter activity and Nanog protein for the combined data (all traces joined end-to-end) in N2B27/RSK(i), alongside controls, as in C. (G) Comparison of maximum cross-correlation between Spry4 promoter and Nanog protein in individual traces against the randomised control in N2B27/RSK(i) conditions. Black shows mean±s.d. Wilcoxon rank sum test (or unpaired two-samples Wilcoxon test) (n=38, two independent imaging experiments). (H) Comparison of peak lag time at maximum cross-correlation for all the cells that show a significant maximum cross-correlation above the ‘randomised signal’ control. Data values were compared with the lag time observed between Nanog and Spry4 traces time shifted by a random time (‘randomised time’). Black shows mean±s.d. Statistics show comparison of the distributions using Kolmogorov–Smirnov (n=35, two independent imaging experiments).
Live imaging of ERK activity and Nanog protein dynamics. (A) Representative time series images of Spry4-Fluc, Nanog::Nluc signal after withdrawal of 2i culture conditions. Red, orange and yellow circles show different cells across selected frames. (B) Representative trace of Nanog::Nluc and Spry4-Fluc over ∼38 h of differentiation. (C) Cross-correlation function plot between Spry4 promoter activity and Nanog protein for the combined data (all traces joined end-to-end) in N2B27/DMSO, alongside controls in which we measured the cross-correlation between Nanog and a randomised Spry4 signal (rand. signal) or between Nanog and Spry4 traces time shifted by a randomised time (rand. time). (D) Comparison of maximum cross-correlation between Spry4 promoter and Nanog protein in individual traces against the ‘randomised signal’ control in N2B27/DMSO conditions. Black shows mean±s.d. Wilcoxon rank sum test (or unpaired two-samples Wilcoxon test) (n=44, two independent imaging experiments). (E) Comparison of peak lag time at maximum cross-correlation for all the cells that show a significant maximum cross-correlation above the ‘randomised signal’ control. Data values were compared with the lag time observed between Nanog and Spry4 traces time shifted by a random time (‘randomised time’). Black shows mean±s.d. Statistics show comparison of the distributions using Kolmogorov–Smirnov (n=28, two independent imaging experiments). (F) Autocorrelation function plot between Spry4 promoter activity and Nanog protein for the combined data (all traces joined end-to-end) in N2B27/RSK(i), alongside controls, as in C. (G) Comparison of maximum cross-correlation between Spry4 promoter and Nanog protein in individual traces against the randomised control in N2B27/RSK(i) conditions. Black shows mean±s.d. Wilcoxon rank sum test (or unpaired two-samples Wilcoxon test) (n=38, two independent imaging experiments). (H) Comparison of peak lag time at maximum cross-correlation for all the cells that show a significant maximum cross-correlation above the ‘randomised signal’ control. Data values were compared with the lag time observed between Nanog and Spry4 traces time shifted by a random time (‘randomised time’). Black shows mean±s.d. Statistics show comparison of the distributions using Kolmogorov–Smirnov (n=35, two independent imaging experiments).
We examined the relationship between Spry4 activation and Nanog protein downregulation. After smoothing to remove noise, we used a simple set of ordinary differential equations to calculate the Spry4 promoter activity for each Spry4-Fluc trace (Fig. S2C; see Materials and Methods for details). We created continuous traces by adding the measurements made from each cell end-to-end (Fig. S2D). We then measured the cross-correlation between activity of the Spry4 promoter and Nanog protein level. As controls, we randomised the Spry4 signal in two ways: first, we randomised the Spry4 signals to measure correlations between Spry4 promoter activity and Nanog downregulation that could be attributable to noise; second, we assigned random time shifts to the Spry4 traces recorded (schematics shown in Fig. S2E). Cross-correlation for the real data is higher than for either of the controls, meaning that the measured Spry4 and Nanog signals are correlated above noise levels and there is a consistent time delay between the two signals (Fig. 2C). We repeated the analysis for individual traces and observed the same trend (Fig. 2D,E; Fig. S2F,G). The average lag time is ∼80 min, indicating that activation of the Spry4 promoter precedes Nanog downregulation. We repeated the analysis for RSK(i) treated cells and observed a stronger correlation at the level of the combined dataset (Fig. 2F) as well as in individual cells (Fig. 2G; Fig. S2H). Interestingly, RSK(i) treatment, which leads to a more sustained peak of pERK activity (Fig. S2B), decreased the average delay (lag) between Spry4 promoter activation and Nanog downregulation to 20 min (Fig. 2H; Fig. S2I). The fact that the lag is short, and not evident in all cells, suggests that Nanog downregulation might not require transcriptional activation.
Overall, these analyses indicate that the level or dynamics of pERK that activate Spry4 transcription in individual cells also initiate downregulation of Nanog protein. We surmise that the observed metachronous exit from the naïve state is determined by cell-to-cell variability in pERK activation.
Ongoing ERK activation is required for entry into formative pluripotency
On release from 2i, ERK signalling downregulates the Nanog/Esrrb axis and allows cells to exit the naïve state (Fig. 3A). When MEK(i) is maintained to block ERK activation, cells do not upregulate formative genes by 30 h (Fig. 3Aii,B). However, depletion of Nanog allows us to enforce naïve state exit in the absence of ERK activity. We investigated whether exit in these conditions is sufficient for transition to the formative state (Fig. 3Aiii). We found that, upon Nanog knockdown and collapse of the naïve transcription factor network, cells maintained in MEK(i) remained viable but failed to upregulate formative genes, notably including Otx2, a key effector of formative pluripotency (Fig. 3B) (Acampora et al., 2012; Buecker et al., 2014). This indicates an ongoing requirement for ERK signalling to execute the transition and establish formative pluripotency.
ERK is required for entry into formative state. (A) Schematics summarising the effect of differentiation in control (i), MEK(i) (ii) and MEK(i)+siNanog (iii) conditions. (B) Expression of formative genes by RT-qPCR under different conditions at 30 h. Expression was normalised to Actb, then to the maximum for each gene. Heatmap of means (N=3). (C) Experimental set up and results to rescue expression of formative genes following withdrawal of MEK(i). Expression normalised to Actb. Data are mean±s.d. (N=4). Unpaired, two-tailed Student's t-test: ***P<0.0001. n.s., not significant. (D) Representative immunostaining for Sox1 expression on day 4 of neural differentiation of siNanog cells treated with DMSO or MEK(i). Scale bar: 50 µm.
ERK is required for entry into formative state. (A) Schematics summarising the effect of differentiation in control (i), MEK(i) (ii) and MEK(i)+siNanog (iii) conditions. (B) Expression of formative genes by RT-qPCR under different conditions at 30 h. Expression was normalised to Actb, then to the maximum for each gene. Heatmap of means (N=3). (C) Experimental set up and results to rescue expression of formative genes following withdrawal of MEK(i). Expression normalised to Actb. Data are mean±s.d. (N=4). Unpaired, two-tailed Student's t-test: ***P<0.0001. n.s., not significant. (D) Representative immunostaining for Sox1 expression on day 4 of neural differentiation of siNanog cells treated with DMSO or MEK(i). Scale bar: 50 µm.
We investigated whether MEK(i)+siNanog cells may be stalled in transition and poised for formative progression upon ERK activation. We therefore removed MEK(i) after 30 h to release ERK signalling. However, we saw no increase in expression of Otx2 or other formative marker genes (Fig. 3C). We also assessed whether N2B27/MEK(i)+siNanog cells may be diverted into an extra-embryonic fate. However, we did not detect significant expression of trophoblast markers Cdx2, Eomes or Hand1 (Fig. S3A). Similarly, we could not detect appreciable expression of endoderm markers Sox17, Gata6 or Gata4, ruling out both hypoblast and definitive endoderm differentiation (Fig. S3B).
We considered whether MEK(i)+siNanog cells might experience precocious differentiation, truncating or bypassing the formative stage. To test this possibility, siNanog cells were treated with MEK(i) before withdrawing the inhibitor and measuring expression of Sox1, a neural specific transcription factor, on day 4. Culture of siNanog cells with no exposure to MEK(i) resulted in efficient upregulation of Sox1 on day 4 (Fig. 3D). This confirms that Nanog is not required for cells to gain competence for neural lineage specification. However, treatment of siNanog cells with MEK(i) for 30 h substantially reduced upregulation of Sox1 (Fig. 3D). Furthermore, by removing MEK(i) at different time points, we observed that treatment for longer than 18 h was sufficient to significantly impair neural differentiation (Fig. S3C). Therefore, naïve state exit in the presence of MEK(i) does not facilitate neural differentiation.
In summary, MEK(i)+siNanog cells downregulate ESC markers and exit the naïve state, but they neither upregulate formative genes nor divert to extra-embryonic identity. These findings confirm that ERK signalling is necessary after exit for progression to the formative state and indicate that, without ERK, input cells become stranded in an indeterminate state.
ERK inhibition does not dysregulate the global chromatin landscape upon naïve state exit
ERK signalling has been shown to have effects on chromatin and transcription (Tee et al., 2014). Therefore, we asked whether MEK(i)+siNanog cells were disabled by a global disruption of the chromatin landscape. Using Cut&Tag (Kaya-Okur et al., 2019) we examined the distribution of H3K27me3, H3K4me3 and H3K27ac histone modifications. During normal transition, most H3K27me3 promoter peaks are maintained at 30 h (1472/1796 peaks), and only 7-10% differ (127/1796 lost or 197/1796 acquired) (Fig. S4A). In MEK(i)+siNanog the vast majority of H3K27me3 promoter peaks are similarly maintained (1401/1472 peaks). Therefore, the failure in formative transition is unlikely to be due to global dysregulation of the repressive chromatin landscape. However, nearly 50% of the differential sites were alternatively modulated, with 85/127 sites failing to lose H3K27me3 and 117/197 sites failing to gain H3K27me3 (Fig. S4B,C). Amongst promoter regions that do not lose H3K27me3, many are associated with ERK responsive genes such as Dusp6 or with formative genes like Pou3f1 (Fig. S4F). A small number of promoters (114), including Myc, show aberrant gain of H3K27me3 only in MEK(i)+siNanog conditions (Fig. S4C). H3K27me3 redistribution at enhancers was similar to that observed at promoters, with most peaks unchanged in either condition and, of those that do change, only a small proportion behave differently in MEK(i)+siNanog (Fig. S4D).
H3K27ac deposition is more dynamic during transition. Among detected peaks, more than half are lost during exit (5016/9674 peaks) with only a quarter (2512/9674) maintained, and ∼20% new loci (2146/9674; Fig. S4E). This pattern was not dramatically different in MEK(i)+siNanog-treated cells, although ∼30% of ‘new’ sites (680/2146) did not gain H3K27ac (Fig. S4E), including the ERK target and formative marker Lef1 (Fig. S4F).
Overall, these data do not indicate global dysregulation of histone modifications due to ERK pathway inhibition. Instead, the data show that many ERK response genes and formative genes show marks of silencing [high H3K27me3 and lower H3K4me3 and H3K27ac (Fig. S4F)]. A set of genes also fails to gain H3K27me3 and lose H3K4me3 and H3K27ac (Fig. S4G), but we could not detect any clear pattern in the genes surrounding these differentially bound sites. These differences are consistent with, and potentially consequent to, reduced transcription of those genes. However, higher H3K27me3 deposition at promoters and enhancers of many formative genes may contribute to the inability to restore expression by removing MEK(i) after 30 h (Fig. 3C).
ERK activity is required to maintain Oct4 in the naïve to formative pluripotency transition
The transcription factor Oct4 (encoded by Pou5f1) is the mainstay of pluripotency, expressed through all stages in vivo and in vitro (Schöler et al., 1990; Wu and Schöler, 2014). In the preceding analyses of histone modifications, we noticed specific changes at the Oct4 promoter, which prompted further investigation. In MEK(i)+siNanog, we detected loss of H3K4me3 and, to a lesser extent, of H3K27ac at the Pou5f1 locus at 30 h (Fig. 4A). Oct4 expression is normally maintained throughout the naïve to formative transition. However, in line with disappearance of H3K4me3, we observed that MEK(i)+siNanog cells did not express Oct4 mRNA (Fig. 4B) and also lost Oct4 protein (Fig. 4C,D). Furthermore, release from MEK(i) did not reactivate Oct4 expression (Fig. 4B). We similarly saw that MEK(i)+siEsrrb cells lost Oct4 protein (Fig. S5A). To determine when cells lost Oct4, we determined the percentage of immunopositive cells over time. In N2B27, Oct4 protein levels reduced between ∼30-48 h after 2i withdrawal, as cells began to upregulate Sox1 and enter the neural lineage (Fig. 4C) (Mulas et al., 2017). In contrast, in MEK(i)+siNanog cells Oct4-negative cells were detected as early as 18 h after 2i-withdrawal and expression was extinguished throughout the population by 30 h (Fig. 4C,D). In cultures treated with MEK(i) only, Oct4 protein was downregulated in a fraction of cells from 30-36 h (Fig. 4C,E; Fig. S5B).
Failure to maintain Oct4 expression contributes to defective entry into the formative state. (A) Distribution of H3K27me3, H3K4me3 and H3K27ac along the Pou5f1 locus at 30 h. (B) mRNA levels determined by RT-qPCR of Pou5f1 (Oct4) in control, N2B27/MEK(i)+siNanog cells and after withdrawal of MEK(i) (N=4). (C) Percentage of Oct4-positive cells as determined by immunofluorescence over time in N2B27/DMSO+siNegative, N2B27/MEK(i)+siNanog and N2B27/MEK(i)+siNegative conditions. Pairwise comparison with respect of N2B27/DMSO+siNegative. (D) Representative images showing immunofluorescence for Oct4 protein in N2B27/DMSO+siNegative and N2B27/MEK(i)+siNanog conditions at 24 h. (E) Representative image showing immunofluorescence for Oct4 protein in N2B27/DMSO and N2B27/MEK(i) conditions at 36 h. (F) Rescue of expression of formative genes by constitutive expression of Oct4. Experimental set up (top left) and results. Graphs show gene expression determined by RT-qPCR, normalised to Actb (N=4). Data are mean±s.d. Unpaired, two-tailed Student's t-test: *P<0.05, **P<0.001, ***P<0.0001. n.s., not significant. Scale bars: 100 μm.
Failure to maintain Oct4 expression contributes to defective entry into the formative state. (A) Distribution of H3K27me3, H3K4me3 and H3K27ac along the Pou5f1 locus at 30 h. (B) mRNA levels determined by RT-qPCR of Pou5f1 (Oct4) in control, N2B27/MEK(i)+siNanog cells and after withdrawal of MEK(i) (N=4). (C) Percentage of Oct4-positive cells as determined by immunofluorescence over time in N2B27/DMSO+siNegative, N2B27/MEK(i)+siNanog and N2B27/MEK(i)+siNegative conditions. Pairwise comparison with respect of N2B27/DMSO+siNegative. (D) Representative images showing immunofluorescence for Oct4 protein in N2B27/DMSO+siNegative and N2B27/MEK(i)+siNanog conditions at 24 h. (E) Representative image showing immunofluorescence for Oct4 protein in N2B27/DMSO and N2B27/MEK(i) conditions at 36 h. (F) Rescue of expression of formative genes by constitutive expression of Oct4. Experimental set up (top left) and results. Graphs show gene expression determined by RT-qPCR, normalised to Actb (N=4). Data are mean±s.d. Unpaired, two-tailed Student's t-test: *P<0.05, **P<0.001, ***P<0.0001. n.s., not significant. Scale bars: 100 μm.
We investigated whether maintained expression of Oct4 was sufficient to restore entry into the formative state in MEK(i)+siNanog conditions. For this, we used ZHBTc4 ESCs, in which both Oct4 alleles are inactivated and Oct4 is produced from a regulatable transgene (Niwa et al., 2000). Expression of the transgene is repressed by addition of doxycycline (DOX) (Fig. 4F). Without DOX (Oct4 ON), ZHBTc4 cells efficiently upregulated formative genes after withdrawal from 2i (Fig. 4F). If DOX is added to silence Oct4 expression, ZHBTc4 cells become trophectoderm-like (Niwa et al., 2000). In this condition formative genes are not activated, consistent with a requirement for Oct4 for formative transition. We then treated ZHBTc4 cells with MEK(i)+siNanog. When DOX was added to silence Oct4, cells failed to upregulate formative genes. In contrast, without DOX (Oct4 ON) MEK(i)+siNanog ZHBTc4 cells displayed upregulation of key formative transcription factors Otx2 and Pou3f1. Fgf5 was still not expressed, however; likely because it is a direct target of ERK (Kalkan et al., 2017b).
We conclude that ERK1/2 activity becomes essential for continued Oct4 expression as cells exit the naïve state and that, without Oct4, cells are unable to upregulate formative transcription factors and proceed with transition. Collectively these results demonstrate that ERK1/2 signalling drives ESC transition both by destabilising the naïve state transcription factor network and by sustaining the platform pluripotency factor Oct4. The detailed mechanisms of these two distinct effects remain to be determined. It is known that pERK can phosphorylate Nanog (Brumbaugh et al., 2014; Kim et al., 2014) but whether this is sufficient to acutely destabilise the protein has yet to be determined. pERK has also been shown to phosphorylate METTL3 (Sun et al., 2020) and stimulate m6A methylation of mRNAs including Nanog (Batista et al., 2014), which could alter translation efficiency as well as RNA stability. How pERK contributes to maintaining Oct4 gene expression requires further investigation, but we speculate that this may be linked to the enhancer switching that is known to occur at the Pou5f1 locus during formative transition (Choi et al., 2016; Tesar et al., 2007; Yeom et al., 1996).
In summary, our findings show how variable dynamics of ERK activation in ESCs lead to metachronous dissolution of naïve pluripotency and that ongoing ERK signalling preserves expression of the anchor factor(s) Oct4 to secure cell state transition.
MATERIALS AND METHODS
Culture of mouse ESCs
Mouse embryonic stem cells were routinely cultured in 2i {1 μM PD0325901 [MEK(i)], 3 μM CHIR99021} or 2i/LIF (10 ng/ml LIF Qkine) following previously published protocols (Mulas et al., 2019). Cells were dissociated using Accutase (Millipore, SCR005) and pipetted to obtain a single cell suspension. The cell suspension was diluted in 5× wash buffer [DMEM/F12, 0.03% bovine serum albumin (BSA) Fraction V (Thermo Fisher Scientific)], spun down and resuspended in fresh 2i or 2i/LIF medium before cell counting. Cells were plated at 1.5×104 cells/cm2. For maintenance, cells were plated on dishes coated with 0.1% gelatine for at least 15 min at 37°C (Sigma-Aldrich, G1890), whereas all exit and differentiation experiments were performed in plates coated with 10 μg/ml laminin in PBS for at least 30 min at 37°C (Sigma Millipore, CC095). For all experiments, batch tested, home-made N2B27 was used (N2B27.BV; Mulas et al., 2019). Unless indicated, the following inhibitors were used: MEK(i) PD0325901, 1 μM; RSK(i) BI-D1870 (Axon Medchem, Axon 1528), 3 μM. For experiments with the Zhbtc2 line, 1 μg/ml doxycycline (Sigma-Aldrich) was added to the media when indicated.
siRNA treatment and exit from naïve state
All siRNA experiments were performed in 24-well plates in technical duplicates and biological triplicates, following the protocols indicated in Mulas et al. (2019). For each well, we plated 30,000 cells, 20 nM of total siRNA and 0.5μl Lipofectamine RNAiMax (Thermo Fisher Scientific) in 500 μl 2i. After overnight incubation (<18 h), cells were either collected for RNA extraction and cDNA synthesis to determine the efficiency of knockdown or further differentiated. For differentiation experiments, cells were gently washed with 1 ml PBS before adding fresh media (e.g. N2B27±inhibitors). For each gene we used two siRNAs (10 nM each). Samples were analysed at the indicated time points. The siRNAs used are listed in Table S1 and have been previously validated in Dunn et al. (2014).
RT-qPCR
RNA was extracted using QIAGEN RNeasy or ReliaPrep RNA Miniprep System (Promega). cDNA synthesis was performed using SuperScript III (Invitrogen) following the manufacturer's protocols. Quantitative PCR (qPCR) was performed using TaqMan (Thermo Fisher Scientific), Sybr Green (Thermo Fisher Scientific) or UPL (Roche) technology. A list of primers and probes is provided in Table S1.
Flow cytometry
Cells were dissociated using Accutase until colonies detached, before adding an equal volume of cold PBS/2% foetal bovine serum and dissociated into a single cell suspension using pipettes. Samples were stored on ice until ready to be analysed. Data from at least 8000 cells were collected on a LSR Fortessa Analyzer (BD Biosciences). Forward and side scatter were used to identify cells, while pulse width was used to gate for single cells. To quantify the effect of perturbations on the percentage of Rex1::GFPd2-positive cells, for each biological replicate we set the gate at ∼50% 30 h N2B27/DMSO control. All other samples were analysed with respect of that gate (Fig. S1C).
Immunostaining and quantification
Samples were fixed with cold 4% paraformaldehyde (Santa Cruz Biotechnology) at room temperature for 10 min and washed twice with PBS. Samples were stored in PBS at 4°C until all time points could be stained in parallel. Samples were permeabilised and blocked in block buffer (0.3% donkey serum, 0.15% Triton X-100 in PBS) for at least 2 h at room temperature. Primary antibodies were diluted in block buffer and incubated overnight at 4°C (Table S1). Samples were washed three times with PBST (0.15% Triton X-100 in PBS) for 15 min each time. Alexa-conjugated secondary antibodies raised in donkey (Thermo Fisher Scientific, H+L Highly Cross-Adsorbed Secondary Antibodies) were diluted 1:1000 (2 µg/ml) in block buffer alongside a nuclear stain (DAPI or Hoechst 33342) and incubated for 2 h at room temperature. Samples were washed three times in PBST for 15 min each before storing samples overnight in PBS at 4°C.
Images of fixed samples were acquired using a Leica DMI4000 at 10× or 20× magnification. For each biological repeat of a time series, all images were acquired in a single imaging session. For each well, we acquired three to eight images in a grid pattern. Image processing and segmentation was performed in CellProfiler (https://github.com/CMulas/IF_quantification). Nuclear stain was used to segment nuclei before determining the intensity on the remaining channels. The intensity measures were then analysed using R code. We used a receiver operating characteristic (ROC) curve to determine a binary threshold between positive and negative control cells. For Nanog and Esrrb stainings, day 3 differentiated cells stained with both primary and secondary antibodies were used as negative staining controls (confirmed with siNanog and siEsrrb treatment), whereas 2i cultures were used as positive controls. For Oct4 stainings, secondary-only stainings were used as negative controls, and 2i cultures were used as positive controls. For each image, the percentage of positive cells for each channel was calculated. The statistical method is indicated in the figure legend. Code and example data can be found here: https://github.com/CMulas/IF_quantification.
Live bioluminescent imaging
Sample preparation
Calibration cells (PGK-Nluc-Fluc; Mandic et al., 2017) and cells carrying Spry4-Fluc transcriptional reporter (Phillips et al., 2019) and Nanog::Nluc fusion were routinely cultured in 2i/LIF as described above. The Spry4-Fluc construct contains a splice acceptor site followed by an IRES and an Bsd/F2A/NLSluc cassette and has a half-life of 1.56 h (Phillips et al., 2019). The Nanog:Nluc targeting construct was generating using the previously validated targeting construct for Sox2 (Strebinger et al., 2019) in which 5′ and 3′ homology arms flank a Nluc-loxP-P2A-Puro-sfGFP-loxP cassette. Integration of Nanog::Nluc was initially verified by GFP and Nluc expression, and finally by PCR following excision of the loxP cassette. Using cycloheximide treatment we determined that the reporter had a half-life of 3.02 h. Imaging experiments were performed in FluoroDishes (WPI, FD35-100) coated with CellAdhere (now discontinued) for 30 min in the incubator before washing gently twice with Dulbecco's phosphate buffered saline (DPBS) with calcium and magnesium (Thermo Fisher Scientific/Gibco, 14040133). Then 100,000 cells/cm2 reporter cells and 10,000 calibration cells were plated overnight in 2i medium, before washing gently with DPBS with calcium and magnesium and changing media to N2B27 supplemented with DMSO or RSK(i), 0.5 mM Luciferin (to visualise Fluc, NanoLight Technology, 306A) and 1:2000 RealTime-Glo MT Cell Viability Assay Substrate (to visualise Nluc, first diluted 1:1 in DMSO, Promega G9711). The plate was allowed to equilibrate for ∼30 min before transferring to the microscope and setting up. Therefore, the image acquisition started 45 min post-media change.
Image acquisition
Bioluminescence live imaging was performed on an Olympus LV200 with an EMCCD camera at 60× (oil-immersion) under environmental control (37°C, 7% CO2) using binning 1×1 (512×512) and photon imaging mode 5×. We acquired two images per channel (Fluc: 600 nm LP filter; Nluc: 460/36 nm filter), one with a shorter integration time to measure the bioluminescence signal of calibration cells, one at longer integration time to measure bioluminescent signal of the reporters. The acquisition times for Fluc were 20 ms and 310 ms, and for Nluc were 20 ms and 240 ms.
Image processing and quantification
Image quantification was performed semi-automatically. The background was subtracted for each image to remove optical aberrations. Next, we generated a maximum intensity mask merging both Fluc and Nluc channels to identify objects. We used a custom Fiji plugin to identify regions of interest (ROIs, 4 px diameter circles) in the maximum intensity mask before measuring the intensity in the Nluc and Fluc channel. Cell division was annotated manually. We also measured the background around each ROI for each frame (empty field 10 px diameter circle next to the cell ROI).
Signal normalisation
As Nluc signal intensity increases over time because of the progressive increase in available substrate in cells (Mandic et al., 2017), for each experiment we included a 1:10 dilution of calibration cells. Calibration cells constitutively express an Nluc-Fluc fusion construct (Mandic et al., 2017). For each slide, we calculated the normalisation ratio by dividing the Fluc over the Nluc signal for the calibration cells before applying a Loess regression. Next, we applied a rolling average to smooth the background signal. The Fluc reporter signal was normalised by subtracting the raw signal to the smooth background, and denoised by fitting a smoothing spline regression. The Nluc signal was normalised by subtracting the raw Nluc signal to the background and dividing by the normalisation ratio and denoised using a low pass filter.
Reconstructing promoter activity
We used a previously established approach to infer promoter activity from the de-noised Spry4 reporter readout (Kannan et al., 2018). Briefly, this method relies on a series of ordinary differential equations to infer the promoter activity given a resulting output (in this case bioluminescence signal) and a set of fixed parameters. For simplicity, we converted all time units to slices (=10 min) The key parameters are the Fluc degradation rate [d, from Phillips et al. (2019)=0.07 slices−1], the translation rate [β, estimated from the average translation rate in ESCs, Ingolia et al. (2011)=5.59 slices−1] and the mRNA degradation rate (d.r, average=0.075 slices−1). We assumed that the magnitude of the maturation rate of bioluminescent proteins is negligible.
Autocorrelation between Spry4 promoter activity and Nanog downregulation
First, we filtered the dataset to exclude time points when Nanog is already downregulated. This removed bias between fast and slow differentiating cells. To obtain a combined dataset, we joined the traces of each cell end-to-end. Next, we generated two control datasets. For the randomised control dataset, the Fluc promoter signal was completely randomised for each cell. For the time shifted control dataset, a random frame shift was added to the Fluc promoter signal. Next, we calculated the autocorrelation function (ACF) between the smoothened Nanog::Nluc signal and the smoothened Fluc signal as well as the randomised and time-shifted Fluc controls. To determine whether there was a significant correlation, we compared the maximum ACF for each trace of the data against the randomised control (noise) and performed a Wilcoxon test. To determine if there was a consistent and significant lag (time relationship) between the Fluc promoter signal and Nanog downregulation, we compared the lag of the curves which show significant autocorrelation (above noise) against the time shifted controls (where the time relationship between the two signals should be random). We compared the distributions using a Kolmogorov–Smirnov test.
CUT&Tag
Sample preparation
Each batch consisted of 12 samples processed in parallel. In each batch, the three samples [2i, N2B27/DMSO+siControl and N2B27/MEK(i)+siNanog] were processed in parallel for a given antibody, as well as two negative controls and one positive control (typically H3K27me3). Each sample consisted of a single well of a 12-well plate. Cells were differentiated as described above and processed following the benchtop CUT&Tag protocol v2 (Kaya-Okur et al., 2019; dx.doi.org/10.17504/protocols.io.z6hf9b6) with minor modifications: all Eppendorfs were pre-coated with PBS +0.1% BSA for 30 min and all buffers were kept on ice before use. Nextera primers were used to generate libraries following the benchtop CUT&Tag recommended protocol. Antibodies and primers used are listed in Table S1. Libraries were pooled and sequenced using 50 bp pair-end Nextera (Illumina).
Alignment and normalisation
Samples were processed in Galaxy (Afgan et al., 2022) workflows available in https://github.com/CMulas/CUT-Tag_analysis). FASTA files were aligned using Bowtie2 (Langmead and Salzberg, 2012) against mouse (mm10) and Escherichia coli K12 (spike in control). BAM files were then sorted, removing unaligned reads and PCR duplicates. For each antibody, all samples were normalised to 1× using bamCoverage (deepTools; Ramírez et al., 2016).
Differential binding analysis
For each sample, peaks were called using SEACR (Meers et al., 2019). To determine differential binding for a given antibody, we first generated a file of all possible binding regions across samples. To generate this list, all individual lists of peaks were concatenated, sorted and merged. Next, we centred the peaks and trimmed all regions to 1000 bp. Finally, we generated a count matrix of the total signal over the 1000 bp regions for each sample. Differential analysis was performed using EdgeR (Robinson et al., 2010). For H3K27me3, we divided regions into promoters and enhancers based on overlap with the transcriptional start site (TSS). Peaks were annotated using PAVIS (Huang et al., 2013). The code and list of differentially bound sites can be found in https://github.com/CMulas/2024-Mulas_NanogERK.
Acknowledgements
We thank Peter Humphries and Darren Clements for assistance with microscopy, and Brian Hendrich, Nichola Reynolds, Masaki Kinoshita, Harry Leitch, Jenny Nichols and Lawrence Bates for discussions. Elena Corujo-Simon provided TS cells for positive control RT-qPCR. Pre-conjugated pA-Tn5 was a kind gift from the Henikoff lab. We also thank the Genomics and Tissue Culture facility at the Wellcome–MRC Cambridge Stem Cell Institute for support.
Footnotes
Author contributions
Conceptualization: C.M.; Methodology: C.M., D.M.S.; Software: C.M.; Validation: M.S.; Formal analysis: C.M.; Investigation: C.M., M.S., S.I.S., C.H.; Resources: D.M.S.; Data curation: C.M.; Writing - original draft: C.M., A.S., K.J.C.; Writing - review & editing: C.M., S.I.S., D.M.S.; Visualization: C.M.; Supervision: C.M., A.S., K.J.C.; Project administration: C.M., A.S., K.J.C.; Funding acquisition: A.S., K.J.C.
Funding
C.M. was funded by the King's Prize Fellowship at the Randall Centre for Cell and Molecular Biology, King's College London and a travel grant by the British Society of Developmental Biology/The Company of Biologists. A.S. is a Medical Research Council Professor (G1100526/2). This work was also funded by a Leverhulme Trust grant (RPG-2016-418 to A.S. and K.J.C.) and by a European Research Council grant (‘CellFateTech’ 772798 to K.J.C.). The Wellcome–MRC Cambridge Stem Cell Institute receives core support from Wellcome and the Medical Research Council. Deposited in PMC for immediate release.
Data availability
Cut&Tag sequencing data was deposited in ArrayExpress and can be found under accession number E-MTAB-14031. The Cut&Tag analysis notebooks can be found in https://github.com/CMulas/2024-Mulas_NanogERK. The pipelines to quantify immunofluorescence can be found in https://github.com/CMulas/IF_quantification.
Peer review history
The peer review history is available online at https://journals.biologists.com/dev/lookup/doi/10.1242/dev.203106.reviewer-comments.pdf
References
Competing interests
The authors declare no competing or financial interests.