FGF/MAPK signaling sets the switching threshold of a bistable circuit controlling cell fate decisions in ES cells

Intracellular transcriptional regulators and extracellular signaling pathways together regulate the allocation of cell fates during development, but how their molecular activities are integrated to establish the correct proportions of cells with particular fates is not known. Here we study this question in the context of the decision between the epiblast (Epi) and the primitive endoderm (PrE) fate that occurs in the mammalian preimplantation embryo. Using an embryonic stem (ES) cell model, we discover two successive functions of FGF/MAPK signaling in this decision. First, the pathway needs to be inhibited to make the PrE-like gene expression program accessible for activation by GATA transcription factors in ES cells. In a second step, MAPK signaling levels determine the threshold concentration of GATA factors required for PrE-like differentiation, and thereby control the proportion of cells differentiating along this lineage. Our findings can be explained by a simple mutual repression circuit modulated by FGF/MAPK signaling. This may be a general network architecture to integrate the activity of signal transduction pathways and transcriptional regulators, and serve to balance proportions of cell fates in several contexts.


Introduction
To ensure the faithful development of multicellular organisms, cell fate decisions in populations of undifferentiated cells have to be tightly balanced. It is now well established that transcriptional networks and extracellular signals together control these decisions, but how their interactions determine the proportions of cells differentiating along particular lineages is often not known. This question is of particular importance in one of the first cell fate decisions in the mammalian preimplantation embryo, where a small number of inner cell mass (ICM) cells have to reliably populate both the epiblast (Epi) lineage that will give rise to the embryo proper, as well as the primitive endoderm (PrE) lineage, which differentiates into tissues that function in patterning and nutrient supply of the embryo (Rossant and Tam, 2009). The factors underlying this fate decision have been studied in embryos and embryonic stem (ES) cells (Cho et al., 2012;Fujikura et al., 2002;Shimosato et al., 2007), clonal derivatives of ICM cells that are biased towards the Epi fate but harbor a latent PrE differentiation potential (Beddington and Robertson, 1989;Morgani et al., 2013). Mutant analysis indicates that a transcriptional network centered on the transcription factor Nanog marks and defines the Epi fate (Chambers et al., 2007;Frankenberg et al., 2011), while a network centered on GATA transcription factors underlies the PrE fate (Bessonnard et al., 2014;Schrode et al., 2014). Nanog and GATA6 are co-expressed in early ICM cells, but their expression patterns become mutually exclusive as cells commit to specific fates (Plusa et al., 2008), suggesting mutually repressive interactions between the two programs (Frankenberg et al., 2011;Schrode et al., 2014). Studies using genetic mutants and pharmacological inhibitors have furthermore shown that the FGF/MAPK signaling pathway promotes the PrE fate at the expense of the Epi fate in embryos (Kang et al., 2012;Nichols et al., 2009;Yamanaka et al., 2010), and that it is required for PrE-like differentiation of ES cells (Cho et al., 2012). How the activities of the transcriptional networks are integrated with the activity of the FGF/MAPK signaling pathway, and how these inputs together control the proportion of cells differentiating along either lineage has not been systematically investigated.
Recently, a mathematical model for the decision between the Epi and the PrE fate has been proposed (Bessonnard et al., 2014), in which Nanog and GATA repress each other, and reinforce their own expression through direct positive feedback. This defines a dynamic system with three stable states in which cells either express GATA6 or Nanog alone, or co-express the two markers. In this model, FGF/MAPK signaling both promotes GATA6 expression and inhibits Nanog expression, and differences in FGF/MAPK signaling between cells have been proposed to underlie fate choice from the co-expression state (Bessonnard et al., 2014). While this model is consistent with static phenotypes of wildtype embryos and genetic mutants, the gene expression dynamics proposed by this model have not directly been tested. It is also not clear whether all proposed links are required to explain the behavior of the genetic circuit underlying this cell fate decision, and which one of the two inputs into the system -signaling and transcription factor activity -bears most on the fate decision.
Addressing these open questions requires quantitative modulation of the inputs into the genetic circuit regulating fate choice, and following its dynamics in single cells in real time. Here we achieve this by transiently expressing fluorescently tagged GATA Fate Decisions in ES Cells 4 factors in ES cells carrying live reporters for the Epi and the PrE fate. This allows us to recreate a state of co-expression of Epi and PrE determinants akin to the state of ICM cells in the embryo, and to follow the resolution of this state in real time. We find that cells rapidly exit the co-expression state towards one of two mutually exclusive states, i.e. the system is bistable. PrE-like differentiation occurs in cells exposed to GATA factor levels above a threshold, and the function FGF/MAPK signaling is to set this threshold dose. This provides a mechanism through which both transcription factor activity and signaling can tune the proportions of cells with specific fates. Recapitulating the dynamic behavior of the circuit in silico only requires mutual repression between the transcriptional networks underlying the Epi and the PrE fates without any positive feedback loops, and a single repressive input of MAPK signaling on the Epi-specific program. This data-based model for the Epiversus-PrE fate decision, much simpler than previously proposed models, will serve as a basis to guide further experimental and theoretical exploration of this crucial fate decision of mammalian embryogenesis. Furthermore, our finding that FGF/MAPK signaling can balance the proportions of alternative fates in cell populations by setting the response threshold of a regulatory network to a transcription factor input is a novel principle for this signaling pathway which might be relevant in developing tissues beyond the ICM.

An ES-cell model system to investigate PrE-like fate choice in culture
To model in culture the transition from Nanog;GATA6 co-expression to mutually exclusive expression of Epi and PrE markers that characterizes the Epi-versus-PrE fate decision (Plusa et al., 2008), we used a doxycycline-inducible system to transiently express GATA6-FLAG in ES cells (Beard et al., 2006;Mulvey et al., 2015;Wamaitha et al., 2015) (Fig. 1A). Individual cells co-expressed inducible GATA6-FLAG and endogenous Nanog protein after a 6 h doxycycline pulse (Fig.   1B). 24 h after doxycycline removal, the cells had degraded the exogenous GATA6-FLAG, but a subset now stained positive for the endogenous PrE marker GATA4 ( Fig. 1C). Virtually all GATA4-positive cells were negative for Nanog staining, suggesting that following GATA6/Nanog co-expression, ES cells transition to one of two mutually exclusive states, marked by the expression of Epi-and PrE markers, respectively. This is similar to the behavior of ICM cells, and suggests that a previously reported stable state of co-expression of Nanog and endogenous GATA factors (Bessonnard et al., 2014) is not accessible in our system. Consistent with previous studies (Fujikura et al., 2002;Mulvey et al., 2015;Shimosato et al., 2007)

ES cell culture conditions affect the expression of endogenous PrE markers following a GATA4-mCherry pulse
For an induced transcription factor to trigger a specific differentiation program, this program needs to be molecularly accessible. We therefore set out to determine culture conditions for which transient GATA4-mCherry expression led to efficient expression of endogenous GATA6. While in the presence of feeders and 15% serum a 6 h pulse of GATA4-mCherry expression resulted in approx. 10% GATA6-positive cells one day later, this proportion dropped to approx. 1.5% GATA6-positive cells when cells were cultured without feeders in 10% serum (supplementary material Fig. S3A, B), even though GATA4-mCherry was efficiently induced in both conditions (supplementary material Fig. S3C), and cells were positive for the pluripotency marker Nanog (supplementary material Fig. S3A). Next, we pre-cultured cells in 2i + leukemia inhibitory factor (LIF), a condition reported to promote extraembryonic differentiation potential (Morgani et al., 2013), before simultaneous addition of doxycycline and transfer into serum-containing medium. This increased the proportion of GATA6-positive cells induced by a 6 h doxycycline pulse from 11.3 ± 1.8 % (mean ± standard deviation (SD)) for 1 day of pre-culture in 2i + LIF to 51.7 ± 9.8 % for 7 days of pre-culture ( Fig. 2A -C). Because the duration of the preculture in 2i + LIF also affected the levels of GATA4-mCherry expression induced by doxycycline (Fig. S4A), we determined the ratio between the fraction of GATA6positive cells one day after a 6 h doxycycline pulse and the fraction of GATA4-mCherry-positive cells immediately after the pulse as a measure for the efficiency of PrE-like differentiation. This ratio plateaus at approximately 55% after 3 days of preculture in 2i + LIF (Fig. 2C).
To assess the influence of the components of 2i + LIF, we removed each of them from the complete 2i + LIF medium or added them individually to serum-containing medium during 3 days of pre-culture. All conditions led to an increase in the

Fig. 3: Transient GATA4-mCherry expression induces stable expression of PrE marker genes. A Immunostaining for Venus and GATA6 protein in Gata6:H2B-Venus reporter cells 24 h after a 6 h pulse of doxycycline-induced GATA4-mCherry expression. Scale bar, 50 µm. B Flow cytometry of Gata6-reporter cells treated as in A. C Flow cytometry detecting Gata6:H2B-Venus expression at indicated time-points after a 6 h pulse of GATA4-mCherry induction. D Percentages of Venus high cells at different times after a 6 h doxycycline pulse. Data points represent mean ± SD from three independent experiments. E Flow cytometry of Venus high cells sorted 18 h after the end of a 6 h doxycycline pulse (dark blue), cultured for 48 h and analyzed for Venus expression. F H2B-Venus fluorescence intensity values of individual Gata6-reporter cells tracked in time-lapse movies during and following a 6 h doxycycline pulse (shaded area). Colorcode is informed by hierarchical clustering based on H2B-Venus expression levels.
Small panels on the right show traces for each cluster separately. See also supplementary material Fig. S6 and Movie 1.

A threshold level of GATA4-mCherry controls PrE-like differentiation
We then asked whether the flipping of the bistable switch that we had identified To correlate GATA4-mCherry input levels more precisely with subsequent fate choice in single cells, we performed time-lapse imaging of GATA4-mCherryinducible cells during and after a doxycycline pulse, followed by immunostaining for Nanog and GATA6 (supplementary material Movie 2, Fig. S9A, B). We found that most cells with GATA6-positive progeny had experienced higher GATA4-mCherry expression levels than cells with Nanog-positive progeny (purple and green datapoints in Fig. 4E, F, and in supplementary material Fig. S9C). We used receiver operating characteristics (ROC) analysis (Fawcett, 2006; for details see supplementary material) to assess how well the two classes of cells could be separated by a threshold value of GATA4-mCherry expression. Plotting the ratio of correctly and incorrectly separated events over the total number of cells (true positive ratio TPR, and false positive ratio, FPR) for varying threshold values gives a characteristic curve for a single time-point ( Fig. 4G); the larger the area under this curve (AUC), the better the differentiation outcomes can be separated or predicted from GATA4-mCherry expression levels. The AUC increased quickly upon doxycycline addition and reached a plateau between 0.8 and 0.9 after approx. 3 h (Fig. 4H). Similar results were obtained when we used cumulative instead of instantaneous GATA4-mCherry expression measurements (not shown), suggesting that non-systematic measurement errors are not a major limit of predictive power. The optimal prediction threshold that maximizes the difference between TPR and FPR tracked the expression dynamics of the GATA4-mCherry

FGF/MAPK signaling modulates the proportion of cells with PrE-like differentiation
In the mouse embryo, both GATA factors and FGF/MAPK signaling are required to establish PrE identity. Having shown above that inhibiting MAPK signaling prior to doxycycline-induced GATA expression increases the proportion of cells with PrE-like differentiation, we next wanted to test how MAPK signaling affected the decision to embark on PrE-like differentiation during and after the GATA pulse. MAPK activity required for PrE-like differentiation was almost completely saturated in serum-free medium, possibly through autocrine FGF signaling (supplementary material Figs S11, S12), prompting us to tune the levels of Erk phosphorylation following removal of the pre-culture medium with subsaturating doses of PD03 (Fig. 5A). Partial inhibition of MAPK signaling during the 6 h doxycycline pulse and the 24 h chase period reduced the fraction of GATA6 positive cells, but not the expression levels of GATA6 in individual differentiated cells (Fig. 5B, C), with a quasi-linear relationship between the level of Erk phosphorylation and the percentage of differentiating cells (Fig. 5D).
We obtained similar results using the FGF receptor inhibitor PD173074 (supplementary material Fig. S12), indicating that most of the MAPK activity relevant for PrE-like differentiation of ES cells is triggered by FGF ligands, consistent with literature reports (Kunath et al., 2007). FGF/MAPK signaling levels therefore control the fraction of cells that embark on the PrE-like differentiation path.
To investigate how partial Mek inhibition affected the GATA4-mCherry threshold required for PrE-like differentiation, we performed time-lapse imaging and cell tracking for maximal and reduced MAPK signaling in parallel (Fig. 5E, F), using a PD03 concentration that led to a significant reduction of the number of differentiating cells without inducing cell death (supplementary material Movie 3). ROC analysis gave similar AUC values for both conditions, indicating that differentiation can be predicted based on GATA4-mCherry expression levels with similar confidence at different signaling levels (Fig. 5G). However, the optimal prediction threshold was consistently increased upon partial Mek inhibition (Fig. 5F, H). We conclude that MAPK signaling levels set the GATA4-mCherry threshold dose required to trigger differentiation.
We noticed that the distribution of GATA4-mCherry expression levels in differentiating and non-differentiating cells changed upon partial inhibition of signaling (Fig. 5E, F). In addition to setting the transcription factor threshold dose, partial Mek inhibition therefore appears to modulate heterogeneities in the population that affect PrE-like differentiation.

A simple mutual repression circuit recapitulates the experimentally observed gene expression dynamics
To gain insights into the formal nature of the interactions between signaling and transcriptional regulators we then sought to identify the minimal circuit model of the components of the decision machinery that would recapitulate our data. The irreversible, switch-like behavior of our system indicates the presence of positive feedback in the underlying genetic network. Because Nanog directly represses Gata6 (Singh et al., 2007), and GATA expression led to rapid repression of Nanog in our system, we chose a network of two mutually repressive nodes, Gata and Nanog, as a minimal system with net positive feedback to formalize a bistable genetic switch (Cherry and Adler, 2000;Plahte et al., 1995;Snoussi, 1998;Thomas, 1981) (Fig. 6A, see supplementary material for a detailed description of the model). This system is described by two coupled ordinary differential equations that account for the dynamics of Nanog (N) and endogenous GATA (G) as markers for the Epi and PrE programs in individual cells, respectively: A third equation models the externally supplied pulse of GATA (G X ) that drives the endogenous circuit: To reflect the experimentally observed heterogeneous expression of exogenous GATA factors, we varied the maximum transcription rate D of exogenous GATA between cells (for details see supplementary material) . This was the only source of cell-to-cell variability in our model. As initial conditions, we endowed cells with high levels of Nanog and no GATA to reflect pre-culturing in the presence of PD03.
To assess the dynamics of the endogenous circuit described by this model we plotted the nullclines dG/dt=0 and dN/dt=0 for the specific set of parameters used. Two of the three equilibrium states defined by the intersections of the nullclines are stable and correspond to the fully differentiated GATA-positive and the undifferentiated Nanogpositive state, respectively (Fig. 6B). A boundary in the phase space (dashed line in endogenous Gata6 gene (Fig. 6C, compare to Fig. 3G), and exogenous GATA4-mCherry (Fig. 6D, compare to Fig. 4E), suggesting this simple mutual repression circuit is sufficient to capture essential dynamics of the experimental system.
To further test the model, we compared the dynamics of Nanog expression in silico and in vivo. In model simulations, Nanog expression levels first decreased rapidly in all cells from the initial conditions chosen to represent the effects of the pre-culture regime towards lower steady state levels, before differences in Nanog expression levels in differentiating and non-differentiating cells became apparent (Fig. 6E).
These expression dynamics were observed experimentally in cells expressing a translational Nanog fusion reporter (Filipczyk et al., 2013) following transient expression of GATA4-mCherry (supplementary material Movie 4, Fig. 6F). This agreement between Nanog expression dynamics in silico and in vivo further supports the idea that a simple mutual repression circuit is sufficient to capture the dynamics of the system.

Inhibition of the Epi program by MAPK signaling controls the proportion of cells with PrE-like differentiation
To pinpoint the main mechanism by which FGF/MAPK signaling controls the fraction of cells with PrE-like differentiation, we considered two simple extensions of the model, one in which FGF/MAPK signaling promotes expression of the PrE program (Fig. 7A), and an alternative model in which signaling inhibits the Epi program (Fig. 7B). In both cases a reduction of signaling led to a simulated increase in the number of cells in the Nanog-positive peak and a decrease of cells in the GATA-positive peak (Fig. 7A, B, middle). However, expression levels in the respective positive peaks changed in distinct ways depending on the type of signaling input (Fig. 7A, B,  Finally, we sought to develop a visual representation of the system's dynamics in different signaling regimes by estimating the path-integral quasi-potential surfaces (Bhattacharya et al., 2011) of the system for two different signaling levels (Fig. 7E,  F). This representation highlights two basins of attraction corresponding to the Nanog-positive state and the GATA positive state, respectively. A reduction in signaling bends the ridge that separates the basins towards the GATA-positive state, making its basin of attraction narrower and shallower relative to that of the Nanogpositive state (Fig. 7E, F).
We conclude that a simple mutual repression circuit is sufficient to capture the dynamic hallmarks of the Epi-versus-PrE fate decision, and that through repression of the Epi program, FGF/MAPK signaling sets the relative sizes of the basins of attraction corresponding to the two fates defined by this circuit, allowing signaling to regulate the proportions of cells adopting either fate. respectively. We detect a well-defined threshold level of exogenous GATA factor expression required to flip this switch and induce differentiation, and find that MAPK signaling sets this threshold dose. This decision is therefore a strongly regulated process that is largely determined by few well-defined transcriptional and signaling inputs.

The accessibility of the PrE program depends on ES cell culture conditions
We find that cells cultured in the presence of serum are refractory to PrE-like differentiation upon induced GATA expression, but responsiveness to doxycyclineinduced GATA factors can be restored by extended exposure to GSK3 or Mek inhibitors, e.g. in 2i medium. One interpretation of this finding is that ES cells cultured in serum are strongly biased towards embryonic fates, and as a consequence have blocked the PrE-like differentiation program. In line with this idea, ES cells grown in the presence of serum display higher levels of repressive chromatin marks on a subset of promoters, including the Gata6 promoter, than cells grown in 2i + LIF (Marks et al., 2012). Furthermore, the transcriptome of ICM cells resembles more closely that of ES cells cultured in 2i medium than that of ES cells cultured in serum (Boroviak et al., 2014). This suggests that pre-culture in 2i brings ES cells to a molecular state mirroring that of early ICM cells, from which, upon induced GATA expression, the decision between the Epi and the PrE fate can be taken similar to the situation in the embryo.

Extraembryonic fate choice is determined by the output of a simple mutual repression circuit
Our finding that precise measurements of GATA4-mCherry expression levels allow  (Albeck et al., 2013;Aoki et al., 2013), and we expect they will be functionally relevant for PrE-like differentiation of ES cells.

Integration of signaling into the mutual repression circuit serves to balance proportions of cell fates in developing tissues
The mathematical model of a mutual repression circuit has previously been applied to describe the dynamics of the switch between the lysogenic and lytic phases of the lifecycle of bacteriophage lambda (Ptashne, 2004), and a genetically engineered toggle switch circuit in E. coli (Gardner et al., 2000). Our work is one of the first experimentally supported examples demonstrating that this network can be used to formalize the decision between two fates during mammalian development. Extending the model with a signaling input allows for dynamic control of the sizes of the basins of attraction corresponding to the different states of the bistable system. The mammalian preimplantation embryo may harness this property to balance the proportion of Epi and PrE cells. The initial expression of transcriptional regulators driving lineage choice is stochastic, possibly due to the mechanisms that control gene expression in the early embryo (Dietrich and Hiiragi, 2007;Ohnishi et al., 2014), and the resulting heterogeneous distributions of transcription factor concentrations in ICM cells bias cells towards specific fates (Xenopoulos et al., 2015). It has been shown that lineage commitment occurs non-synchronously in the cells of the ICM, and that the first cells to commit are fated towards the epiblast (Grabarek et al., 2012). Because Epi cells produce Fgf4 (Frankenberg et al., 2011;Nichols et al., 1998), Fgf4 levels will reflect the number of Epi-committed cells and act on the as-yet uncommitted cells. By modulating the bistable switch operating in these cells, this process may ultimately place the appropriate number of uncommitted cells in the basin of attraction corresponding to the PrE fate. FGF/MAPK signaling might thus acts as feedback mechanism to balance proportions of two distinct cell fates in populations (Lander et al., 2009). It will be interesting to see whether this new principle applies to differentiation decisions beyond those in the preimplantation embryo.
All cell lines used in this study were based on the KH2 ES cell line (Beard et al., 2006). Engineering of ES cells is described in more detail in the supplementary material. Transgene expression was induced by adding 500 ng/ml doxycycline to the culture medium.

Flow cytometry
Cells for flow cytometry were trypsinized and either analyzed immediately or fixed for 15 minutes in 3% PFA/PBS. Intracellular antigens were stained in suspension using the same primary and secondary antibodies as used for immunostaining. mCherry fluorescence was measured on a BD Fortessa Flow cytometer, all other flow cytometric analysis was performed using a Beckman Coulter CyAn ADP analyzer.
Cell sorting was done on a Beckman Coulter MoFlo. To estimate peak positions, histograms were smoothened, followed by detection of local maxima with customwritten Python scripts.
Detection was performed using fluorescently labeled secondary antibodies and scanning in a LI-COR Odyssey system. Intensities of bands were quantified in ImageStudio (LI-COR).

Time-lapse imaging and cell tracking
Time-lapse imaging was performed in DMEM based medium without phenol red, supplemented as detailed above. We used a Zeiss Axiovert M200 microscope equipped with a SOLA LED light source, an Andor iXON Ultra 888 EMCCD camera and a heated stage with CO 2 supply. Hardware was controlled by MicroManager software (Edelstein et al., 2001). Time-lapse movies were acquired using a 40x longworking distance lens. See Supplemental Experimental Procedures for details on image analysis.

Mathematical modeling
Numerical simulations of the model were implemented in Python language. Parameter values used in the simulations are given in supplementary material Table S1. For details on the model see supplementary material.

Author contributions
CS, PR and AMA conceived the study; CS performed experiments; CS, PR and JPM analyzed the data; PR developed the mathematical model; CS wrote the manuscript with input from PR and AMA; all authors approved the manuscript.