1. Computational and Systems Biology
  2. Developmental Biology
Download icon

Rapid changes in morphogen concentration control self-organized patterning in human embryonic stem cells

  • Cited 0
  • Views 870
  • Annotations
Cite this article as: eLife 2019;8:e40526 doi: 10.7554/eLife.40526

Abstract

During embryonic development, diffusible signaling molecules called morphogens are thought to determine cell fates in a concentration-dependent way. Yet, in mammalian embryos, concentrations change rapidly compared to the time for making cell fate decisions. Here, we use human embryonic stem cells (hESCs) to address how changing morphogen levels influence differentiation, focusing on how BMP4 and Nodal signaling govern the cell-fate decisions associated with gastrulation. We show that BMP4 response is concentration dependent, but that expression of many Nodal targets depends on rate of concentration change. Moreover, in a self-organized stem cell model for human gastrulation, expression of these genes follows rapid changes in endogenous Nodal signaling. Our study shows a striking contrast between the specific ways ligand dynamics are interpreted by two closely related signaling pathways, highlighting both the subtlety and importance of morphogen dynamics for understanding mammalian embryogenesis and designing optimized protocols for directed stem cell differentiation.

Editorial note: This article has been through an editorial process in which the authors decide how to respond to the issues raised during peer review. The Reviewing Editor's assessment is that all the issues have been addressed (see decision letter).

https://doi.org/10.7554/eLife.40526.001

Introduction

Mammalian development depends crucially on diffusible signaling molecules called morphogens, that are thought to determine cell fates in a concentration-dependent manner (Green et al., 1992; Wilson et al., 1997; Wolpert, 1969), and protocols for directed stem cell differentiation are based on this picture (Chambers et al., 2009; Loh et al., 2016; McLean et al., 2007; Mendjan et al., 2014). However, in the vertebrate embryo, expression patterns of these morphogens change rapidly, simultaneous with large-scale cell movements, and therefore individual cells experience substantial changes in morphogen levels during differentiation (Arnold and Robertson, 2009; Balaskas et al., 2012; Dessaud et al., 2007; Dubrulle et al., 2015; Kinder et al., 2001; van Boxtel et al., 2015). Duration of ligand exposure must logically be a relevant parameter, as was confirmed in a number of contexts, including Activin/Nodal signaling in early Zebrafish development and Sonic Hedgehog (Shh) signaling in the mouse neural tube (Dessaud et al., 2007; Sako et al., 2016). However, duration is only one of a large number of features of a dynamic signal. It is unknown whether the precise time course of ligand exposure plays a role in cell fate decisions, and if so, whether different pathways interpret signaling histories differently. In analogy with human speech, which enables sophisticated communication by relying on temporal modulation of a single mode of signaling (sound), it is possible that complex information is encoded in developmental signals by temporal modulation to enable a range of different responses to a single pathway. In addition to ligand concentration and duration (‘integral’), cells may also be sensitive to ligand rates of change (‘derivative’), and it has been suggested that adaptive signaling pathways allow cells to perform this derivative computation (Sorre et al., 2014; Yi et al., 2000).

These concepts have begun to be explored in mammalian cell culture systems. For the Nκfb pathway, an intricate relation between ligand dynamics, signaling response and target gene activation was found (Hoffmann et al., 2002; Kellogg et al., 2017; Selimkhanov et al., 2014). Further demonstrating the importance of changing ligand concentrations, a class of ERK target genes was recently shown to be activated more efficiently by pulses than sustained signaling in 3T3 cells (Wilson et al., 2017), and we have shown that the response to TGFβ in C2C12 mouse myoblasts reflects ligand rate of increase (Sorre et al., 2014; Warmflash et al., 2012). However, these ideas have not been applied to mammalian development or differentiation of pluripotent stem cells, and their relevance to developmental patterning remains unexplored.

Gastrulation is the first differentiation event of the embryo proper, when the germ layers are formed and the body axes are established. Nodal and BMP4 are morphogens crucial for gastrulation in vertebrates (Winnier et al., 1995). In the early mammalian embryo, BMP4 is required for both initiating gastrulation and specifying the dorsal-ventral axis, while Nodal maintains the pluripotent epiblast and is subsequently required for mesoderm and endoderm differentiation (Conlon et al., 1994). Each pathway has distinct receptor complexes that phosphorylate specific signal transducers, known as receptor-Smads, which then complex with the shared cofactor Smad4 and translocate to the nucleus to activate target genes (Figure 1a) (Zinski et al., 2018). Both Nodal and BMP4 have been claimed to act in a concentration-dependent manner based on classic Xenopus experiments (Green et al., 1992; Gurdon et al., 1994; Wilson et al., 1997). However, those experiments allow alternative interpretations (Heemskerk and Warmflash, 2016), and the role of BMP4 and Nodal ligand dynamics has not been investigated.

Micropatterned human embryonic stem cells (hESCs) self-organize into reproducible spatial domains corresponding to each of the germ layers and were recently established as a method to recapitulate human gastrulation in vitro (Warmflash et al., 2014). This system can be easily manipulated and observed. Further, in contrast to the micropatterned colonies, when hESCs are grown more sparsely, their response to exogenous signals is uniform and not dependent on secondary signals, allowing for dissection of the dynamics of response (Nemashkalo et al., 2017). Thus, the response of cells to dynamic signals can be systematically investigated in sparse culture, and this information can be used to unravel the complexity of self-organized pattern formation in micropatterned colonies. In this study, we take this approach and use hESCs to evaluate the role of changing concentrations of BMP4 and Activin/Nodal in cell-fate decisions associated with gastrulation. Unexpectedly, we find an important role for rapid concentration changes in Nodal pathway response, while the BMP pathway responds to concentration and duration of ligand exposure more directly.

Results

SMAD4 signaling response of hESCs to BMP4 is sustained but to Activin is adaptive

To investigate how sudden increases in BMP4 and Nodal levels are interpreted by hESCs, we performed live imaging of hESCs with GFP:SMAD4 in the endogenous locus (Nemashkalo et al., 2017), and quantified signaling strength as the ratio of nuclear and cytoplasmic SMAD4 intensity (Figure 1—figure supplement 1a). We found that a sudden increase in BMP4 leads to sustained SMAD4 signaling (Figure 1b,d), consistent with our previous work (Nemashkalo et al., 2017). In contrast, the response to addition of Nodal and its substitute Activin is strongly adaptive, that is transient, and returns to a signaling baseline of around 20% of the response peak (Figure 1c,d, Figure 1—figure supplement 1i, j), similar to the previously observed response to TGFβ in C2C12 cells (Sorre et al., 2014; Warmflash et al., 2012). The response to recombinant Nodal was weak at all doses, likely reflecting the low quality of recombinant Nodal. To elicit a response in Nodal required 2 μg/ml of recombinant protein, whereas the response to Activin was saturated by 5 ng/ml (Figure 1f, Figure 1—figure supplement 1i, j). However, when the responses to Nodal and Activin were normalized to their respective maxima, their dynamics were identical (Figure 1—figure supplement 1j). Given our focus on the dynamics of signaling, the extremely similar dynamics in response to Nodal and Activin, and the impracticality of performing experiments with Nodal given the concentrations required, we used Activin in all further experiments. The peak and baseline signaling levels, but not the timescale of adaptation, depended on the Activin dose (Figure 1f). For BMP4, there was a sharp response so that even low doses gave a nearly full initial response, however, the dose affected the duration of signaling in a way consistent with ligand depletion at low doses (Figure 1e, Appendix 1). Immunofluorescence staining for receptor-Smads revealed that SMAD1/5/8 activation mirrors the SMAD4 response to BMP4, while in response to Activin, SMAD2/3 nuclear localization adapts less than SMAD4 to about 60% of peak response (Figure 1—figure supplement 1e,g,h). In each case, the response was due to the added ligand and not to induced secondary signaling through the other pathway, as addition of the BMP inhibitor Noggin had no effect on the dynamics in response to Activin, while addition of the Activin inhibitor SB431542 slightly increased the response to BMP (Figure 1g,h). Thus, endogenous BMP has no effect on the Nodal response, while Nodal signaling has a repressive effect on BMP. We note that some target genes may respond to SMAD2/3 without SMAD4 (Levy and Hill, 2005), and that since the dynamics of SMAD2/3 and SMAD4 are different, care is required in interpreting the dynamics of individual genes. Hereafter, when we refer to Nodal signaling as reflected by the SMAD4 reporter, we mean ‘SMAD4-dependent Nodal signaling’.

Figure 1 with 3 supplements see all
SMAD4 signaling response of hESCs to BMP4 is sustained while that to Activin is adaptive.

(a) BMP and Nodal pathways share the signal transducer Smad4. (b, c) hESCs expressing GFP:SMAD4 at 0, 1 and 6 hr after treatment with BMP4 (b) or Activin (c). Scalebar 30 μm (d) GFP:SMAD4 average nuclear:cytoplasmic intensity ratio after treatment with BMP4 (blue) or Activin (red). Error bars represent standard error. Ncells ~ 700, distributions shown in (Figure 1—figure supplement 1b–c). (e) SMAD4 response to different doses of BMP4 shows decline at low doses with a dose-dependent time scale, suggesting ligand depletion. Doses in graph legend are in ng/ml. (f) SMAD4 signaling response to different doses of Activin shows fixed time scale of adaptation. (g) Quantification of GFP:SMAD4 nuclear to cytoplasmic ratio in response to either Activin alone or together with the BMP inhibitor Noggin (h) Quantification of GFP:SMAD4 nuclear to cytoplasmic ratio in response to either BMP alone or together with the Activin/Nodal inhibitor SB431542.

https://doi.org/10.7554/eLife.40526.002

Adaptive signaling is caused by negative feedback controlling sequestration

Although TGFβ signaling through SMAD4 has been shown to adapt in C2C12 cells (Sorre et al., 2014; Warmflash et al., 2014), the molecular mechanisms remain unclear. A previous attempt to uncover relevant genes through a genome-wide siRNA screen uncovered a large number of potential regulators of signaling, but none of the knockdowns completely blocked adaptation, suggesting redundancy in adaptation mechanisms (Deglincerti et al., 2015). As an alternative approach, we used FRAP at different times after Activin treatment to better understand how cells modulate nuclear SMAD4 levels over the course of the signaling response (Figure 2a–c). In particular, we performed measurements before ligand treatment, after 1 hr of treatment when signaling peaks, and after 6 hr when signaling has adapted (Figure 2c). In all cases, we observed recovery on a timescale of minutes following nuclear photobleaching, which confirmed that SMAD4 continuously shuttles between the nucleus and cytoplasm even in the absence of ligand stimulation (Figure 2d). We also observed a reduced recovery rate after 1 hr of Activin treatment, consistent with previous work on TGFβ signaling which suggested that nuclear accumulation of SMAD4 is caused by reduced nuclear export (Nicolás et al., 2004; Schmierer and Hill, 2005). Importantly, however, we found that the recovery rate is not restored during adaptation (Figure 2d), showing that cells do not revert to a pre-signaling state. This excludes upstream mechanisms of adaptation as might be caused by secretion of extracellular feedback inhibitors or depletion of receptors or R-Smads (Vizán et al., 2013). Moreover, mathematical modeling showed that if adaptation is caused by depletion of an upstream component such as receptors, then the magnitude of adaptation is governed by the ratio of the degradation rates of active and inactive receptors (Appendix 1). These degradation rates also control the timescales for adaptation and recovery so that strong adaptation necessitates a large difference in timescales. However, we do not see these different time scales in our pulse experiments below (see Figure 5), which show that the refractory period after ligand exposure is similar to the time scale of adaptation. In contrast, models with negative feedback causing inhibition of downstream signaling are capable of explaining all features of our data including the observed time scales (Appendix 1).

Adaptive Activin response is not a return to the pre-stimulus state and is explained by a model that includes sequestered SMAD4 populations.

(a) Photobleaching and recovery of nuclear Smad4 at 2 hr after Activin treatment. Scalebar 5 μm. (b) Exponential fit to recovery of nuclear fluorescence after bleaching yields amplitude A and recovery rate k. Intensity drop at t = 0 shows bleaching event. (c) Photobleaching was performed on untreated cells, at the peak response to Activin, and after adaptation. (d) Boxplot of distribution of recovery rates. Recovery rates at peak signaling and after adaptation are significantly smaller than for untreated cells (t-test p < 10−6), the difference in recovery rate between peak signaling and adapted state is not significant. N > 12 cells for each FRAP condition. (e) Boxplot of distribution of recovery amplitudes. (f) Nuclear to cytoplasmic intensity ratio after bleaching (R’) is systematically smaller than nuclear to cytoplasmic intensity ratio before bleaching (R). (g) Cartoon of mathematical model for Smad4 localization, with sequestered populations of Smad4 that are confined to either nucleus or cytoplasm, and free Smad4 shuttling with import/export rates kin/kout. (h) Cartoon demonstrating that this model explains results in (f). (i) Changes in Smad4 sequestration in nucleus (blue) and cytoplasm (red) determined through model fitting. Error bars in i and j represent error propagation of standard errors in measured parameters over N > 12 cells as described in Appendix 1. (j) Nuclear export rates, which given the fixed nuclear import rate of the model directly reflect the measured exchange rates k.

https://doi.org/10.7554/eLife.40526.007

To understand the mechanism of this inhibition, we first considered fitting our FRAP data to a model in which the entire population of SMAD4 is free to shuttle between the nucleus and cytoplasm. In this case, bleaching reduces the total pool of observable GFP-SMAD4 molecules but the nuclear to cytoplasmic intensity ratio depends only on the kinetic constants, and therefore would recover to the same value following bleaching. However, we observed that this ratio systematically decreases after nuclear bleaching (Figure 2f). This suggests a more general model that includes sequestered SMAD4, which may move within but not between the nuclear and cytoplasmic compartments. Because production and degradation are slow compared to shuttling, recovery from bleaching comes only from redistribution of unbleached molecules. The nuclear sequestered population is not exported to the cytoplasm, while the cytoplasmic sequestered population is unavailable to enter the nucleus and replenish the population of unbleached molecules there, leading to a reduced nuclear to cytoplasmic ratio after equilibration of the free dark and fluorescent molecules through exchange (Figure 2g,h, Appendix 1). The parameters obtained from fitting this model to our data suggest that initial accumulation of SMAD4 in the nucleus reflects both a lower export rate and a reduction in cytoplasmic sequestration (Figure 2i,j). Adaptation, however, does not result from altering shuttling kinetics, but instead reflects reduced nuclear sequestration and increased cytoplasmic sequestration. Thus, taken together, our FRAP data and mathematical modeling suggest a mechanism of adaptation that relies on negative feedback and acts by modulating sequestration rather than nuclear exchange rates.

Transcriptional dynamics of differentiation targets follows SMAD4 dynamics

Next, we evaluated transcriptional dynamics downstream of BMP4 and Activin using qPCR, which showed that BMP targets are stably induced (Figure 3a), while differentiation targets of Nodal show adaptive transcription on a timescale consistent with SMAD4 signaling (Figure 3b). Moreover, shared targets of the pathways were found to be transcribed adaptively in response to Activin and stably in response to BMP4 (Figure 3c, Figure 3—figure supplement 1a). In contrast, the transcription of NODAL, WNT3, and their inhibitors LEFTY1 and CER1 were sustained upon Activin treatment (Figure 3d). Molecularly, the two classes of transcriptional dynamics in response to Activin may reflect differential requirements for SMAD4 signaling levels with lower levels required to maintain the targets with sustained dynamics so that these are continuously transcribed due to the baseline signaling following adaptation. Alternatively, transcription of these genes may require only SMAD2/3 activation, which is more sustained than that of SMAD4 (Figure 1—figure supplement 1e,g,h). The differences in expression of these sets of targets are not due to differences in mRNA stability as mRNAs for stably expressed genes were found to decline rapidly upon pathway inhibition with SB431542 indicating a need for ongoing signaling to maintain expression (Figure 3—figure supplement 1g).

Figure 3 with 1 supplement see all
Transcription of BMP targets and Nodal differentiation targets reflects SMAD4 dynamics, while other Nodal targets show sustained transcription.

(a, b) qPCR measurements of transcriptional response to BMP4 treatment (a) and of differentiation targets to Activin (b) y-axes show relative CT values. (c) Transcription of the shared Activin/BMP4 target HAND1 after BMP4 (blue) or Activin (red) treatment. (d) Non-adaptive response to Activin of ligands and inhibitors involved in initiating the primitive streak. (e) Transcriptional response to Activin under pluripotency maintaining conditions (red) and mesendoderm differentiation conditions (blue) of Activin target LHX1 (e) and joint Activin/Wnt target EOMES (f). Error bars represent standard deviations over three replicates. Logarithms are base 2.

https://doi.org/10.7554/eLife.40526.009

The sustained transcription of Nodal and Wnt pathway ligands and inhibitors may be required to activate the positive feedbacks between the Nodal and Wnt pathways, which are known to be involved in establishing the primitive streak, the region of the mammalian embryo where mesoderm and endoderm form (Ben-Haim et al., 2006). This suggests a picture where stable transcription of the ligands and inhibitors allows for the establishment of signaling patterns in the embryo, while cells receiving these signals to differentiate interpret them according to their dynamics. Several other genes not related to mesendoderm differentiation were also found to be stably induced by Activin (Figure 3—figure supplement 1b–d).

The measurements above were performed by adding Activin to mTeSR1 media which contains high levels of FGF. Activin/Nodal signaling plays a dual role in the early embryo and is involved in both maintaining pluripotency, and differentiation to primitive streak fates. It maintains epiblast pluripotency in combination with FGF (James et al., 2005; Vallier et al., 2005), and while differentiation genes are transiently induced by Activin treatment in the presence of only FGF, cells do not robustly differentiate. In contrast, Activin induces primitive streak fates when Wnt is present. We asked whether target genes also respond adaptively to Activin/Nodal signaling during this differentiation. Initial transcriptional response was found to be qualitatively similar with or without Wnt activation, although shared targets of Wnt and Activin such as the mesendodermal markers EOMESODERMIN (EOMES) and BRACHYURY (BRA) showed a stronger response and were stably reactivated following adaptation on longer time scales (Figure 3e,f, Figure 3—figure supplement 1e,f).

BMP4 response reflects concentration, Activin response reflects rate of concentration increase

Sustained response to BMP4 suggests sensing of ligand concentration. In contrast, adaptive response to Activin with dose-dependent amplitude suggests sensing of ligand rate of increase. We directly tested whether cells are sensitive to the rate of increase of Activin, but not BMP4, by slowly raising ligand concentrations (concentration ramp) and comparing the response with the response to a single step to the same final dose (Figure 4a,d). If cells are primarily sensitive to ligand doses, the step and ramp should eventually approach the same final activity, while if cells are sensitive to the rate of ligand increase, the response to the ramp should be reduced. As expected, BMP4 signaling responses to the ramp and step approached each other as ramp concentration increased, while Activin signaling was dramatically reduced in the case of the ramp (Figure 4b,e). Moreover, transcriptional dynamics of the shared target HAND1 matched the signaling pattern and showed dramatically reduced transcription in response to the Activin but not the BMP ramp (Figure 4c,f). As noted in the discussion of Figure 1e, the BMP response is switch-like with increasing dosage, and this was reflected by the ramp approaching the levels of the step within the first few increases. Similar results were obtained for other adaptive Activin targets, in contrast to non-adaptive Activin targets which, as expected, also showed sustained transcription in response to the ramp (Figure 4—figure supplement 1). Finally, we varied the rate of the ramp of Activin and found that intermediate rates of increase also yielded intermediate signaling responses (Figure 4g,h).

Figure 4 with 1 supplement see all
BMP4 response reflects concentration, but Activin response reflects rate of concentration increase.

(a, d) Ligand concentration over time for slow ramp (blue) or sudden step (red) in BMP4 (a) or Activin (d). (b, e) SMAD4 signaling response to BMP4 (b) or Activin (e). (c, f) Transcriptional response of HAND1 to concentration ramp versus step for BMP4 (c) or Activin (f). Error bars in qRT-PCR data (c, f) represent standard deviations over three replicates. (g, h) Activin response to different ramp rates and step sizes. In (b, e, h) small hourly wiggles are artifacts of performing media changes and do not reflect actual signaling responses.

https://doi.org/10.7554/eLife.40526.012

Repeated rapid increases in Activin/Nodal enhance differentiation to primitive streak fate

Morphogens control cell fate, and the dependence of transcription of mesodermal and endodermal genes on Activin rate of increase suggests rapid Activin increase may boost differentiation to these fates. We hypothesized that exposing cells to repeated rapid increases by pulsing the level of Activin could enhance this effect. To rigorously test this hypothesis, we grew cells in differentiation conditions and compared pulses that switch between high and low doses of Activin with a sustained high Activin dose while performing media changes at the same times (‘dummy pulses’) (Figure 5a). Each pulse of Activin elicited a strong response in the translocation of SMAD4 to the cell nucleus, while no such responses were seen in response to the dummy pulses (Figure 5b). Differentiation to mesendooderm as marked by Bra expression was also enhanced under pulsed conditions (Figure 5c,d), despite the reduced integrated ligand exposure. Continuous exposure to lower doses of Activin and treating cells with Activin for the same total time as the three pulses but in a single pulse showed that the effect is specific to pulsed Activin and not a consequence of simply reducing the integrated Activin exposure (Figure 5e–g, Figure 5—figure supplement 1).

Figure 5 with 1 supplement see all
Repeated rapid increases in Activin/Nodal enhance differentiation to primitive streak fate.

(a) Schematic of pulse experiment, graph shows ligand concentration, controls receive media changes at the same time as the pulsed well. (b) SMAD4 signaling profile in response to Activin pulses (red), high Activin (yellow), and no Activin (blue). (c) Immunofluorescence staining for BRA after high constant Activin or pulsed Activin. Scalebar 100 μm. (d) Distribution of BRA expression per cell (Ncells per condition ~6 × 103) determined from immunofluorescent images (c). (e, f, g) Dose response series showing BRA expression monotonically increases with Activin dose, and therefore the effect of pulses is not due to reduced average Activin exposure. (e) Immunofluorescence staining for BRA after 34 hr differentiation with different doses of Activin. (f) Distributions of BRA intensity per cell in the images containing (d). (g) Cumulative distributions of BRA intensity.

https://doi.org/10.7554/eLife.40526.015

Rapid changes in endogenous Nodal signaling occur during self-organized patterning

To test whether rapid concentration increases are relevant to endogenous Nodal during embryonic patterning, we turned to micropatterned colonies of hESCs treated with BMP4. These colonies differentiate in a spatial pattern with reproducible rings of extraembryonic cells and all three germ layers, and represent a model for the patterning associated with gastrulation in the human embryo (Heemskerk and Warmflash, 2016; Warmflash et al., 2014). Nodal signaling is required for the formation of mesoderm and endoderm within these colonies, and both small molecule inhibition and genetic knockout of Nodal drastically reduce mesendoderm differentiation (Chhabra et al., 2018; Warmflash et al., 2014). For the dynamics described here to be relevant, Nodal signaling should evolve rapidly compared to the timescale for adaptation, and we should observe rapidly changing signaling patterns with the GFP:SMAD4 reporter. In contrast, if cells read a stable gradient of Nodal protein during patterning, as posited by classic models, we would expect correspondingly stable patterns in signaling and the adaptive dynamics would not be relevant.

To distinguish these hypotheses, we used live imaging to observe SMAD4 signaling in micropatterned colonies during the 42 hr in which these patterns form (Figure 6a-c). The cells initially respond uniformly to the BMP4 treatment. The response then is restricted to the colony edge by approximately 12 hr (Figure 6a,b), and this pattern is maintained until approximately 25 hr. Beginning at 25 hr, a wave of increased nuclear SMAD4 spreads inward again from the edge (Figure 6a,c). SMAD4 convolves the BMP and Nodal responses, and we hypothesized that the stable response at the edge of the colony represents BMP signaling while the inward moving wave results from Nodal signaling. To test this, we looked at the pathway specific SMAD1 and SMAD2/3 (Figure 6d-g). These confirmed that the BMP-specific active SMAD1 signaling remains restricted to the edge (Figure 6f), while the Nodal transducer SMAD2/3 activity spreads rapidly inward from the colony edge between 24 and 36 hr (Figure 6e). The SMAD2/3 and SMAD4 signal transducers reveal a rapidly evolving Nodal distribution with a wavefront that moves through the colony at approximately one-cell diameter (10 μm) per hour. The velocity of the Nodal wavefront suggests that individual cells see Nodal levels increase rapidly in 1 hr, consistent with the time scale of ligand increase required for strong response to exogenous ligand in the previous experiments. Importantly, this wave of signaling activates Bra expression in its wake (Figure 6g,h). This indicates that rapid increases in endogenous Nodal signaling are associated with mesoderm differentiation in this model of human gastrulation.

Figure 6 with 1 supplement see all
Rapid changes in endogenous Nodal signaling occur during self-organized patterning.

(a) Average radial profile of SMAD4 signaling over time (kymograph) in micropatterned colonies after BMP4 treatment (N = 4 colonies). (b) SMAD4 signaling in single colony at 40 hr. (c) Radial SMAD4 signaling profiles at discrete times from 20 hr to 40 hr. (d) Immunofluorescence staining 40 hr after BMP4 treatment for pSMAD1, SMAD2/3 and BRA. (e, f, g) Normalized radial profiles of SMAD2/3 (e) pSMAD1 (f) and BRA (g) averaged over N > 5 colonies per time. (h) Half-maximum versus time for SMAD2/3 (blue), SMAD4 (purple) and pSMAD1 (red) and BRA expression domain defined by a threshold of at least 20% of maximal expression (yellow). In all panels, error bars represent standard deviations taken over different colonies.

https://doi.org/10.7554/eLife.40526.018

Discussion

Our work shows that morphogens in the mammalian embryo do not act in a purely concentration-dependent matter and has revealed an important role for Activin/Nodal rate of change in specifying cell fate as a consequence of adaptive signaling. Adaptive signaling could serve to restrict the response to a narrow competence window or to separate the multiple roles of a single morphogen by using distinct dynamics to selectively activate different target genes, effectively expanding the information content of the morphogen gradient (Selimkhanov et al., 2014). It is currently unclear whether adaptive signaling is an intrinsic feature of the Activin/Nodal pathway or is context dependent. A recent study found that Activin target genes are induced more stably when cells are also treated with Wnt, but that the dynamics of Smad2 translocation were not affected (Yoney et al., 2018). For Wnt, costimulation with Activin or BMP switches the dynamics of Wnt signaling from transient to sustained, and it will be interesting to determine whether Activin/Nodal signaling is sustained in some contexts (Massey et al., 2019).

Activin and Nodal are generally thought of as text-book morphogens acting in a concentration-dependent manner (Gilbert, 2014). This is based on experiments in which an artificial gradient is created by placing an Activin-soaked bead in an intact animal cap or in which dissociated Xenopus animal cap cells are exposed to a range of Activin concentrations (Green et al., 1992; Gurdon et al., 1994). It is important to note that these experiments are not inconsistent with our findings. In both cases, the highest concentration corresponds to the highest rate of concentration increase, and our results also show the peak of the transient signaling response induced by a step increase is dose-dependent. This fact could underlie the differences in cell fates observed in dissociated animal cap cells.

Our results are in contrast to those for another well-characterized system, the patterning of the neural tube by Shh, which has also been shown to be adaptive and regulates cell fate decisions through a described gene regulatory network (Balaskas et al., 2012). In that case, adaptation depends on upstream inhibition and results in a relatively constant integrated signaling output in which time and duration can be interchanged (Dessaud et al., 2007). This interchangeability between concentration and duration can be achieved by upstream feedback such as degradation of activated receptors. Intuitively, increasing ligand concentration increases the response but also leads to more rapid degradation of receptors limiting the response in time.

Several lines of reasoning argue that adaptation in Activin/Nodal signaling follows a different mechanism, with consequences for cell fate. The signaling behavior of our system is qualitatively different, as there is no tradeoff between concentration and duration of strong signaling. Longer exposure to ligand does not increase the duration of signaling because the system has already adapted, and reactivating the pathway can only be achieved by pulsing rather than extending the duration of ligand exposure. Interestingly, there is a low level of baseline signaling following adaptation which is maintained by continued ligand exposure. Given the distinct classes of target genes we found responding to each aspect of the signal (baseline and adaptive pulse), this may demonstrate the principle that dynamic signals relay more information through a single pathway, with one set of target genes responding to signal duration, and another to rapid signal increase.

In addition to the qualitative behavior being inconsistent with a system that integrates signal, like the neural tube, our quantitative data and mathematical modeling rule out the mechanism that would most naturally implement such a behavior, namely adaptation through degradation. First, our FRAP measurements show that the adapted state is kinetically distinct from the pre-signaling state, while the degradation model would predict a return to the same state. Second, we show analytically in Appendix 1 that a model based on upstream feedback cannot account for our dynamic measurements because in order to adapt through receptor degradation, the timescale for recovery must be slower than that for adaptation, which is not what we find in our pulsing experiments. It will be interesting to elucidate the genetic network that interprets adaptive Nodal signaling and to compare it with the one that interprets Shh.

Our conclusions regarding signaling dynamics rely on nuclear localization of Smads as a proxy for signaling activity. The validity of this measure is supported by strong correlation between SMAD2 nuclear localization and pSMAD2 levels, as well as the fact that transcriptional activity of a number of direct targets mirrors nuclear SMAD4 both in this paper and in previous work (Warmflash et al., 2012). Although some target genes do not follow this time course, it is important to note that target gene dynamics represent those of signaling filtered through a specific promoter (Dubrulle et al., 2015). For example, a promoter with high affinity for the signal transducer will saturate at low levels and will not reflect fluctuations in signaling that remain above these levels. Although it was shown that blocking nuclear export of Smad4 does not interfere with function of the pathway (Biondi et al., 2007), this does not invalidate our conclusions, as blocking nuclear export of Smad4 destroys the strong correlation between pathway activity and nuclear localization by artificially keeping Smad4 in the nucleus following the termination of transcription.

Although relatively unexplored in the context of development, adaptive signaling is a common feature of biological systems. Well-studied examples in bacteria include the chemotactic response where adaptive signaling allows cells to sense the rate of change and move to areas of higher nutrients (Block et al., 1983), and the σb response where it allows cells to sense the rate of increase of environmental stress (Young et al., 2013). In mammalian cells, stimulation with TNFα leads to a transient pulse of NFκB signaling, followed by a low, sustained baseline, and this has been suggested to allow specificity in target gene activation with different targets responding to different features of the dynamics (Hoffmann et al., 2002). Similarly, DNA damage has been shown to give rise to stereotyped pulses of p53 signaling due to negative feedback (Lahav et al., 2004).

It will be important to understand how our results fit with current protocols for directed stem cell differentiation which typically perform a single media change per day and therefore do not take advantage of potential enhancements resulting from optimizing ligand dynamics. Routinely, ranges of ligand concentration are explored to optimize differentiation protocols. Our results suggest a different approach that starts with characterizing the signaling dynamics during each cell fate decision and then optimizing the relevant dynamics of ligand presentation. As such, we believe the results presented here will serve as a basis for a dynamic understanding of embryonic patterning and stem cell differentiation.

Materials and methods

Key resources table
Reagent type
(species) or resource
DesignationSource or referenceIdentifiersAdditional information
Cell line (Homo sapiens)ESI017ESIBIORRID:CVCL_B85
Cell line (H. sapiens)RUES2Ali Brivanlou (Rockefeller)RRID:CVCL_B810
AntibodyGoat polyclonal anti-BrachyuryR and D SystemsRRID:AB_2200235(1:400)
AntibodyRabbit monoclonal anti-Sox2Cell Signaling TechnologiesRRID:AB_1904142(1:200)
AntibodyMouse monoclonal anti-NanogBD BiosciencesRRID:AB_1645598(1:400)
AntibodyMouse monoclonal anti-Smad2/3BD BiosciencesRRID:AB_398162(1:100)
AntibodyRabbit monoclonal anti-pSmad1Cell Signaling TechnologiesRRID:AB_2493181(1:100)
Software, algorithmImage processing and data analysis codeThis studyhttps://github.com/idse/stemcells/ commit a5ee164

Cell lines

The cell lines used were ESI017 and RUES2 GFP:Smad4 RFP:H2B. ESI017 cells were obtained directly from ESIBIO while RUES2 were a gift of Ali Brivanlou (Rockefeller University). The identity of these cells as pluripotent cells was confirmed via triple staining for pluripotency markers OCT4, SOX2, and NANOG. All cells were routinely tested for mycoplasma contamination and found negative.

Cell culture and differentiation protocols

The cell lines and their maintenance are described in Nemashkalo et al. (2017). For all experiments except micropatterning, cells were seeded at a low density of 6 × 104 cells per cm2 and grown with rock-inhibitor Y27672 (10 uM; StemCell Technologies). This ensured uniform response to exogenous ligand and minimized the effect of secondary endogenous signaling. Unless otherwise indicated in the figures, experiments in Figures 13 were done in mTeSR1 medium (StemCell Technologies) (also referred to as pluripotency conditions), and cells were treated with 50 ng/ml Activin A (R and D Systems) or 50 ng/ml BMP4 (R and D Systems). Noggin was used at 500 ng/ml. SB431542 at 10 μM. Differentiation conditions in Figure 3e,f are defined as Essential six medium (Gibco) +3 μM CHIR99021 (StemCell technologies).

Imaging and image analysis

Live imaging was done on an Olympus/Andor spinning disk confocal microscope with a 40×, NA 1.25 silicon oil objective. Immunofluorescence data for Figure 5 were collected using a 20×, NA 0.75 objective on an Olympus IX83 inverted epifluorescence microscope. Fixed micropatterned colonies for Figure 6 were imaged using an 30x NA 1.05 silicon oil objective on an Olympus FV1200 laser scanning confocal microscope. All image analysis used Ilastik (Sommer et al., 2011) for initial segmentation and custom written MATLAB code, available at https://github.com/idse/stemcells for further analysis (Heemskerk, 2019; copy archived at https://github.com/elifesciences-publications/stemcells). Figure 1—figure supplement 1a shows how Smad4 nuclear to cytoplasmic ratio, our measure for signaling intensity, is defined from the segmentation by subtracting the mean value in the background mask from nuclear and cytoplasmic masks for each cell. Smad2/3 signaling was similarly measured by nuclear to cytoplasmic ratio. Since pSmad1 and Bra are purely nuclear, these were measured by nuclear intensity, normalized by DAPI to correct for intensity variation due to optics.

FRAP

FRAP experiments were done on an Olympus FV1200 laser scanning confocal microscope using a 100x NA 1.49 oil objective. To deal with cell movement recovery, curves were obtained by reading out a polygon defined by the interpolation between the initial bleach window and a manually defined polygonal nuclear mask in the final frame. Supplemental Information describes the mathematical model for nuclear localization of Smad4 which is shown schematically in Figure 2g,h, and the inference of model parameters shown in Figure 2i,j from the data.

qPCR

For qPCR experiments, ESI017 cells were grown in 24-well plates. RNA was extracted using Ambion RNAqueous-Micro Total RNA Isolation Kit and cDNA synthesis was performed with Invitrogen SuperScript Vilo cDNA Synthesis Kit according to the manufacturer’s instructions. Measurements were performed with SYBR green and the primers in the Table 1. ATP5O was used for normalization in all experiments.

Table 1
qPCR primers used in this study.
https://doi.org/10.7554/eLife.40526.023
ATP5OACTCGGGTTTGACCTACAGCAAAATGAACGGACAGAACCG
BRATGCTTCCCTGAGACCCAGTTGATCACTTCTTTCCTTTGCATCAAG
CERACAGTGCCCTTCAGCCAGACTACAACTACTTTTTCACAGCCTTCGT
GATA3TTCCTCCTCCAGAGTGTGGTAAAATGAACGGACAGAACCG
HAND1GTGCGTCCTTTAATCCTCTTCGTGAGAGCAAGCGGAAAAG
ID2GCAGCACCTCATCGACTACAAATTCAGAAGCCTGCAAGGA
ID4CCCTCCCTCTCTAGTGCTCCGTGAACAAGCAGGGCGAC
LEFTY1ACCTCAGGGACTATGGAGCTCAGGAGAAATGGCCAATTGAAGGCCAGG
LHX1TCCCCAATGGTCCCTTCTCCGTAGTACTCGCTCTGGTAATCTCC
MIXL1CCGAGTCCAGGATCCAGGTACTCTGACGCCGAGACTTGG
NANOGCCGGTCAAGAAACAGAAGACCAGACCATTGCTATTCTTCGGCCAGTTG
NODALATGCCAGATCCTCTTGTTGGAGACATCATCCGCAGCCTAC
NOGCATGAAGCCTGGGTCGTAGTTCGAACACCCAGACCCTATC
OCT4GGGCTCTCCCATGCATTCAAACCACCTTCCCTCCAACCAGTTGC
SOX2CCATGCAGGTTGACACCGTTGTCGGCAGACTGATTCAAATAATACAG
TBR2/EOMESCACATTGTAGTGGGCAGTGGCGCCACCAAACTGAGATGAT
WNT3CTCGCTGGCTACCCAATTTGAGCCCAGAGATGTGTACTGC

Pulses

For pulse and dose response experiments in Figure 5, differentiation was done in Essential six medium +1 uM CHIR99021 +20 ng/ml bFGF (Life Technologies)+10 uM Y27672 with 30 ng/ml Activin added as indicated in the figures. Time between pulses was always 5 hr to allow the pathway to relax. Duration of individual pulses was chosen for experimental convenience and lengths of pulses shown in Figure 5 are 6, 10, and 8 hr. Controls were subjected to media changes at the same time. During media changes, cell were washed three times with PBS.

Micropatterning

For micropatterned colonies, we followed the protocol in Deglincerti et al. (2016) using the chemically defined medium mTeSR1. For fixed micropatterns, we used the CYTOO Arena EMB chip and analyzed the 800 um colonies, while for live imaging we used a CYTOO 96-well plate RW DS-S-A, which has 700 um colonies. For the analysis in Figure 6, radial profiles of SMAD1 and SMAD2/3 were normalized to have the lowest signaling level be zero and the highest be one in at each time, as minimal and maximal levels were similar at each time and only their spatial distribution varied. SMAD4 was normalized such that its maximum over all positions interior in the colony to the most interior half-maximum of pSMAD1 intensity was equal to one. This position was chosen for normalization as it reflects the peak of Nodal-dependent SMAD4 signaling. Finally, because BRA levels change substantially, BRA at all times was normalized to the maximum of the latest time (48 hr).

Immunostaining

Fixing and immunostaining of cells was done as described in Nemashkalo et al. (2017). Table 2 lists the antibodies that were used.

Table 2
Antibodies used for immunofluorescence in this study.
https://doi.org/10.7554/eLife.40526.024
ProteinSpeciesDilutionCatalog no.Vendor
BRAGoat1:400AF2085R and D Systems
SOX2Rabbit1:2005024SCell Signaling Technology
NANOGMouse1:400560482BD Biosciences
SMAD2/3Mouse1:100610842BD Biosciences
pSMAD1Rabbit1:10013820Cell Signaling Technology

Replicates, sample sizes, and error bars

All experiments were performed at least twice. Data shown are from representative experiments, except for FRAP data, where data from multiple experiments are pooled because of the difficulty of performing FRAP on a sufficient number of cells in a given experiment. Error bars on single-cell data (SMAD4 signaling dynamics, pSMAD1, SMAD2/3, BRA, and HAND1 immunostainings) are standard error over cells. Full distributions are also provided to show cell-to-cell variability. Error bars on qRT-PCR data are over technical replicates due to the difficulty of quantitatively comparing biological replicates, likely owing to differences in culture densities, timing between seeding and ligand stimulation, and other culture variables, however, multiple biological replicates were performed in all cases. Error bars in micropatterning experiments are standard deviation over colonies. For single-cell imaging experiments, at least 600 cells were measured for each condition, for FRAP experiments at least 12 cells were measured for each condition, and for micropatterning experiments at least five colonies were analyzed for each condition.

Supplemental information

Supplemental information includes four figures, three movies, and an Appendix on mathematical modeling.

Appendix 1

‘Rapid changes in morphogen concentration control self-organized patterning in human embryonic stem cells’

FRAP measurements and inference of model parameters

Kinetic model

In this section we formulate a kinetic model for Smad4 localization and work out its consequences, the following sections will show that the assumptions on which our model is based are consistent with the data.

The simplest model for Smad4 involves a single population of freely shuttling Smad4 whose total Tf changes slowly compared with the time scale of shuttling and whose localization is determined entirely by nuclear import and export (Schmierer and Hill, 2005). The nuclear free Smad4 Nf and cytoplasmic free Smad4 Cf then obey

(1) N˙f=kinCfkoutNf,Tf=Nf+CfC˙f=N˙f,

where kin and kout denote the nuclear import and export rate, respectively. Combining these equations we obtain

(2) N˙f=kinTf-kNf,kkin+kout,

from which we see that the system approaches equilibrium exponentially with a time scale 1/k. In equilibrium, the free populations are given by

(3) Nfeq=Tfkinkin+kout=κTf,Cfeq=Tfkoutkin+kout=(1-κ)Tf,κkinkin+kout

We therefore see that in this model, the nuclear to cytoplasmic ratio is given by

(4) NC=κ1-κ=kinkout,

which is independent of the amount of Smad4 and therefore unaffected by bleaching. In other words, the prediction is that after recovery from photobleaching the nuclear to cytoplasmic intensity will have the same value as before bleaching. This is not what we find (Figure 2f).

We therefore consider a more general model in which we assume that nuclear and cytoplasmic Smad4 consist of both free and sequestered populations

(5) N=Nf+Ns,C=Cf+Cs,

and where the free part cycles between nucleus and cytoplasm. We define sequestered as unable to shuttle between nucleus and cytoplasm, but possibly still free to diffuse within the compartment. It is important that this is different from the definition in Schmierer and Hill (2005) where sequestration is taken to imply reduced mobility within each compartment. We are agnostic about the molecular mechanism of sequestration. Moreover, slow shuttling of the sequestered populations and slow interconversion between free and sequestered (i.e. changes in Tf) do not affect the conclusions as long as the time scales are long compared to the time scale of shuttling.

What we measure is fluorescence intensity which corresponds to Smad4 density, related to total Smad4 by the volume N=nVn,C=cVc, so the measured nuclear to cytoplasmic ratio is given by

(6) rnc=αR=αNs+κTfCs+(1-κ)Tf,αVcVn

Measurements show that total Smad4 turns over very slowly (half life greater than 12 hr), so we can still consider it fixed on the time scale of our experiments

(7) T=Tf+Ns+Cs.

Moreover, we are only interested in the relative size of the Smad4 populations, so we can set T=1, and end up with the equilibrium nuclear to cytoplasmic ratio

(8) r=αNs+κ(1NsCs)Cs+(1κ)(1NsCs).

If nuclear exchange is much faster than the rate of change of any of the parameters in the expression above, we will always find the system in equilibrium and we can understand changes in nuclear to cytoplasmic ratio through the above expression. Since we do not expect α to change rapidly in response to Activin or BMP4, we conclude that an increase in r can be achieved through increase in Ns, decrease in Cs, increase in kin, or decrease in kout. As we will show in the next section, this model has a enough freedom to fit our data and makes some basic predictions that we verified.

Fluorescence recovery after photobleaching

To determine the kinetics of Smad4 shuttling, we used FRAP, which turns part of the Smad4 dark, allowing us to observe the equilibration dynamics of the complement. Nuclear bleaching initially sets the observable part of N=0, while cytoplasmic bleach sets C=0. Denoting quantities after bleaching with a prime, in the case of nuclear bleach we have Ns=0, while the new total free part Tf is the original cytoplasmic free part Tf=Cf=(1-κ)Tf. Substituting this into the expression in the previous section yields the new equilibrium quantities and from (6) yields the prediction that the nuclear to cytoplasmic ratio goes down after nuclear bleach if there is sequestered Smad4

(9) R=NC=AN1-ACR=κ(1-κ)TfCs+(1-κ)2TfR

This inequality becomes an equality only when both Cs=0 and Ns=0. The approach to equilibrium is found by solving (2). Normalizing by the pre-bleach level, the nuclear fluorescence recovery is

(10) N(t)N=κ(1κ)TfN(1ekt)AN(1ekt),

while the cytoplasmic recovery after nuclear bleach is given by

(11) C(t)C=1C(Cs+(1κ)2Tf+κ(1κ)Tfekt)BC+ACekt.

From fitting these functions to the measured FRAP curves we obtain AN,BC,AC and k. Because of the normalization two of these values are not independent: AC+BC=1, so we are left with three measured values: k, AN and AC, while we have four remaining unknown parameters: Ns,Cs,kin,kout. The nuclear to cytoplasmic ratio r provides an additional measurement, but also introduces the additional unknown α.

Although one might have thought bleaching the cytoplasm would provide additional information, it is a simple exercise to check that it does not. Specifically using cb to label quantities after cytoplasmic bleaching, and defining

(12) Ccb(t)C=ACcb(1-e-kt),Ncb(t)N=BNcb+ANcbe-kt,

we find

(13) ACcb=AC,ANcb=AN.

The conclusion is that we have one more parameter than we have measurements, leaving the model underdetermined. We resolve this by assuming kin does not depend on the signaling state of the cell, consistent with the literature (Schmierer and Hill, 2005), and hold it fixed as we solve for the other variables from the data using the equations

(14) AN=κ(1κ)(1NsCs)Ns+κ(1NsCs),K=kin+kout,AC=κ(1κ)(1NsCs)Cs+(1κ)(1NsCs),r=αNs+κ(1NsCs)Cs+(1κ)(1NsCs)

Although for any kin there is a solution to these equations, sensible solutions have 0Ns,Cs1 and kout>0 and fortunately this requirement strongly constrains the value of kin. With the corrections and fitting procedures discussed in the following sections we find that for sensible solutions we must have 6.6×10-4s-1kin1.1×10-3s-1, and that within this range the qualitative behavior of the inferred parameters does not change. We therefore fixed kin in the middle of this range at kin=8.9×10-4s-1.

Practicalities

In practice it turns out to be difficult to bleach the nucleus perfectly and also to prevent the cytoplasm from bleaching somewhat. Due to fast diffusion within the cytoplasm and media, a uniform drop in cytoplasmic intensity can be observed at the time of nuclear photobleaching. Therefore, to represent the fraction of each compartment that is bleached, we define the nuclear and cytoplasmic bleach factors

(15) bn=In-IbgIn-Ibg,bc=Ic-IbgIc-Ibg

Here I and I denote the intensity right before and after photobleaching with subscript indicating the compartment. The background intensity is taken to be the mean intensity in areas without cells. We denote with a tilde, quantities that do not take this correction into account, and define their relationship with the original variables. Now the free fraction after bleaching is more complex

(16) T~f=((1-κ)bc+κbn)Tf.

Moreover, the nucleus starts with intensity bnTf so the recovery amplitude is proportional to the difference between the new equilibrium and this initial value. The recovery curves become

(17) N~(t)N=A~NN(1ekt)+B~NN=κ(T~fbnTf)N(1ekt)+bn,

and

(18) C~(t)C=A~NCekt+B~NC=(1κ)(bcTfTf)Cekt+(1κ)Tf+bcCsC.

The fit parameters A,B are not constrained to sum to one, with B now containing additional information about the degree of bleaching for both compartments

(19) B~NN=bn,A~NC+B~NC=bc.

The equations relating the measured amplitudes to the inferred parameters are now multiplied by a factor bc-bn,

(20) ANNA~NN=(bc-bn)ANN,ANCA~NC=(bc-bn)ANC,

while the baseline for the cytoplasm transforms as

(21) BNCB~NC=(bc-bn)BNC+bn.

There appears to be no expression for nuclear to cytoplasmic ratio after recovery R~ in terms of only R, but we can get the corrected R from plugging the corrected amplitudes (A, rather than the measured A~) into (9).

Parameter inference and sensitivity

To infer NS,Cs,kout we used a weighted least square fit, which minimized

(22) E(kin,kout,Cs,Ns)=iEi2,Ei=pi(p)-p¯iSpi,pi(p)=pi(kin,kout,Cs,Ns)

using the built-in MATLAB function lsqnonlin. Here p denotes the measurements in (14), for example p1=AN, while p denotes the inferred parameters. p¯ denotes the mean of p, and Sp the standard error in the mean of p. The values of the inferred parameters can be found in Figure 1 of the main manuscript. The main feature is that adaptation is caused by changes in Ns,Cs rather than kin,kout.

Error bars in the inferred parameters were determined from the error in the input parameters by error propagation. Specifically, if we write (22) in matrix notation, with δp=p-p¯ and S the diagonal matrix of standard errors, then

(23) E=δpTS-1δp=δpTJTS-1Jδp,Jij=pipj,

where J is the Jacobian, so the errors in the inferred parameters are given by the square root of the diagonal elements of

(24) S=(JTS-1J)-1.

To understand the sensitivity of the measured quantities to changes in the inferred parameters, we can define a sensitivity matrix

(25) Nij=Jijpjpi

which gives the relative change in the measured parameter given for a given relative change in the inferred parameter. Ordering the measurements (column) as (AN,AC,k,R) and the parameters as kout,Cs,Ns, for the parameter values in Figure 1 the sensitivity matrix is given by

(26) Nijuntreated=(0.220.250.460.810.650.0410.86000.580.40.42),Nijpeak=(0.0150.110.350.570.350.0280.76000.730.330.32)
(27) Nijadapted=(0.190.070.0930.821.20.0180.74000.720.750.073)

This effectively decomposes the error in the inferred parameters and for example tells us that AC mostly constrains kout and Cs, but is not sensitive to Ns.

Mathematical modeling of adaptation mechanisms

In this section we extend the purely kinematic model of the previous section to a dynamic model; that is, we explore the predictions of different models for the mechanism of adaptation. In the section titled 'Receptor degradation cannot explain Activin/Nodal signaling dynamics', we show that models where adaptive signaling is caused by depletion of upstream pathway components are not only inconsistent with our FRAP measurements, as noted in the main text, but are also inconsistent with the approximately equal timescales of adaptation and recovery that we observed in our pulse experiments. In contrast, we show in the section titled 'Feedback and feed-forward models are consistent with observed Activin/Nodal signaling dynamics', that either negative feedback or incoherent feed-forward models can account for all observed features of the signaling response. Finally, we show in the section titled 'Ligand degradation is sufficient to explain BMP dose response', that the BMP4 signaling response can be accounted for by ligand degradation.

Minimal model

We consider models for signaling in which a signal transducer can be localized to either the nucleus or the cytoplasm. We assume that the production and degradation of the signal transducer are slow compared to the signaling dynamics so that the total amount of signal transducer is fixed. We adopt the notation from the simple model above ('Kinetic model'). Denoting the total amount of transducer by Tf and the amount in the nucleus by Nf, we can write down a simple equation for the nuclear fraction as:

(28) Nf˙=kin(I)(Tf-Nf)-kout(I)Nf

where kin and kout are the rates of nuclear import and export, respectively, which can both depend on the concentration of an upstream component, I, that can be considered the ligand, the activated receptor, or the receptor-associated Smads. As discussed in the previous section, we choose kin to be independent of I and assume that kout is reduced by I with half-saturation constant K,

(29) kin(I)=kin,kout(I)=kout(0)K+I.

At steady state

(30) N¯f(I)=kinTfk,

where k is the sum of the in and out rates k=kin+kout(I). If the input goes from I0 to I, then the dynamics of Nf are given by

(31) Nf(t)=N¯f(I0)+(1-e-kt)(N¯f(I)-N¯f(I0)).

That is, Nf exponentially approaches its new steady state level with time scale 1/k. As is to be expected, in this simple model, there is no adaptation.

Receptor degradation cannot explain Activin/Nodal signaling dynamics

A simple mechanism of adaptation would be if the upstream component denoted by I was degraded upon activation. To be concrete, we assume that I represents activated receptors, Ii represents inactivated receptors, and that receptors are activated by a ligand with concentration L, where we have subsumed the rate constant for activation into L for simplicity. Then

(32) I˙=LIi-(γ+δ)I,I˙i=β-(δi+L)Ii+γI,

where β represents production of receptors, γ represents return of receptors to the inactive state, and δ,δi respectively represent the degradation rates of activated and inactivated receptors. At steady state, we have

(33) I¯=βLδ(L+δi)+γδi,Ii¯=β(γ+δ)δ(L+δi)+γδi.

Note that when δi=0, I¯=β/δ becomes independent of L so that the model is perfectly adaptive in this case. However, in this case I¯i=β(γ+δ)/(δL) and so diverges as L goes to 0. At non-zero δi, there is a ligand-dependent baseline that eventually saturates at the level β/δ. We first solve the homogeneous equation

(34) 𝒙˙=M𝒙,𝒙=[IIi],M=[-(δ+γ)Lγ-(L+δi)],

where I=I-I¯ and Ii=Ii-I¯i. The solutions to the homogeneous equation can be written as

(35) 𝒙(t)=A𝒗+eλ+t-B𝒗-eλ-t,

where the eigenvalues are given by the solutions to the characteristic polynomial for M

(36) λ±=12[-(L+γ+δ+δi)±(L+γ+δ+δi)2-4[δ(L+δi)+γδi],

and the (non-normalized) eigenvectors are

(37) 𝒗±=[1(λ±+γ+δ)/L].

By taking the gradient of the term inside the square root in the eigenvalue equation (Equation. 36) in the space of the four parameters Lδ, δi, and γ, it is straightforward to show the minimum of this term is 0 and therefore it is always positive. The system is therefore over-damped and will always decay to the equilibrium value.

To proceed further, we consider the case γ=0 which corresponds to receptor binding being irreversible. This is a reasonable assumption as off-rates for ligand-receptor binding are generally small (Sako et al., 2010), and further this is only way to achieve a timescale of adaptation that is independent of the ligand levels, consistent with our data. In this case, steady state levels reduce to

(38) I¯=βδLL+δi,Ii¯=βL+δi,

and the eigenvalues to

(39) λ+=-δ;λ-=-(L+δi).

The eigenvalue λ- will set the timescale for the signaling response following ligand exposure, and λ+ will give the timescale to adapt, which depends only on the degradation constant for activated receptors and is therefore not dependent on the ligand concentration.

We now solve for the dynamics of the non-homogeneous system. We assume that initially L=0, and the system is in equilibrium at this value. That is, I(0)=I¯(L=0)=0 and Ii(0)=I¯i(L=0)=β/δi. These initial conditions yield the equations

(40) 0=I¯+A-B;βδi=I¯i+L+δi-δLB.

Solving these for A and B yields

(41) A=βL(δ-δi)δδi(L+δi-δ),B=βL2δi(L+δi)(L+δi-δ).

In order for the system to adapt, the maximum value that I reaches must be significantly higher than the baseline. Denoting the time at which this maximum is achieved by tmax, we have

(42) I˙(tmax)=Aλ+eλ+tmax-Bλ-eλ-tmax=0,

which yields

(43) tmax=(δ-(L+δi))-1lnδ-δiL.

This implies that as δiδ, tmax. This is logical because if the degradation rates of activated and inactive receptors are the same, there will be no adaptation and the signaling will reach its maximum at t=. Adaptation therefore requires δi<δ.

Defining the maximum relative to baseline, D, as a measure of adaptation, and plugging in the values for A, B, and λ±, we obtain after some algebra

(44) DI(tmax)-I¯I¯

Given that 0<δi<δ and tmax>0, we find that

(45) D<δδi.

This implies that to explain the observed D4 (Figure 1d), we must have δi<δ/4. Now consider if the ligand is removed so that L returns to 0. The system will regain the ability to respond to ligand following adaptation once the number of inactive receptors recovers. From Equation. (39), the rate of recovery in this case is λ-=δi, and recall from the same equation that the rate of adaptation is δ. With δiδ this means that the recovery from adaptation (refractory period) will take significantly longer than the time to adapt. This is not consistent with our data on Activin pulses in Figure 5a,b which shows that a full response can obtained with 5 hr between pulses, similar to the time required to adapt initially. This model is thus excluded as inconsistent with our data. As discussed in the main text and previous section it is also not consistent with our FRAP measurements which show that the adapted state is kinetically distinct from the pre-stimulation state. A model with degradation of upstream components would predict similar kinetics in these two cases.

Feedback and feed-forward models are consistent with observed activin/nodal signaling dynamics

Feedback and feed-forward models provide alternative mechanisms of adaptation. In these models signaling activates an additional component which then causes inactivation of Smad4 and re-localization to the cytoplasm. This component could either be a Smad4 target (feedback model) or induced in parallel such as if it required only Smad2 or was induced by non-canonical signaling (feed-forward model). To be general, we extend the model from above ('Minimal model') to

(46) x˙=aLK+L+bNf-dx Nf˙=kin-koutK+LNf-cxNf

The variable Nf again represents the nuclear Smad4. The variable x represents a component that inhibits Smad4 with the strength of this inhibition governed by the parameter c. This component can either be induced directly by the ligand L (feed-forward) or by Smad4 (feedback), and the strengths of these inductions are controlled by the parameters a and b, respectively. x decays with rate d.

This model is not analytically tractable due to the non-linear term, so we examined whether we could reproduce the data presented in this paper for a variety of different combinations of a and b. In particular, we considered whether we could accurately reproduce the dose response and the comparison between ramps and steps. Appendix 1—figure 1 shows the dose response of Smad4 (Nf) to changing L. Although both the feed-forward and feedback models are capable of reproducing the features of the experimentally determined dose response, the feedback model does so more robustly. While many different values of the feedback strength (b) give the qualitatively correct behavior, when the feed-forward induction (a) is too strong, the curves tend to cross in the dose response due to a lowering of the steady-state value of Nf with increasing ligand. Both models are also capable of reproducing the comparison of the step versus the ramp (Appendix 1—figure 2) with the same caveat that high values of a produced baselines that were reduced compared to the starting levels.

Appendix 1—figure 1
Simulated dose responses for the indicated values of a and b.

The x-axis indicates time and the y-axis nuclear Smad4, both are in arbitrary units. The values of L used were 0, 0.5, 1, 2, 5, and 10. The values of the other parameters are kin=1, kout=20, c=1, d=0.1.

https://doi.org/10.7554/eLife.40526.027
Appendix 1—figure 2
Simulated step and ramp responses for the indicated values of a and b.

Blue lines indicate a step increase in ligand to L=10 while red lines indicate a ramp where the same level was achieved in 10 steps which were spaced by 0.5 time units each. The x-axis indicates time and the y-axis nuclear Smad4, both are in arbitrary units Other parameters are the same as in Appendix 1—figure 1.

https://doi.org/10.7554/eLife.40526.028
Ligand degradation is sufficient to explain BMP dose response

Although the receptor degradation model is not consistent with Activin response, if irreversible receptor binding is taken to imply joint degradation of ligand/receptor complex in the presence of a finite amount of ligand, it can account for the observed dependence of Smad4 signaling on BMP4 dose. Without assuming irreversible binding, we can supplement Equations (32) with an equation accounting for the free ligand

(47) L˙(t)=-L(t)Ii(t)+γI(t)

Integrating the model shows that duration of signaling depends on the initial ligand level, and larger supersaturating doses show sustained signaling for longer times (Appendix 1—figure 3), as experimentally observed for BMP4 (Figure 1e).

Appendix 1—figure 3
Simulated dose-response for the ligand depletion model, showing that this accounts for the BMP4 dose-response.

The parameters used are γ=β=δ=1 and δi=0.1. The initial conditions were I=0, Ii=1 and L=L0, and the value of L0 was varied as indicated in the legend.

https://doi.org/10.7554/eLife.40526.029

References

  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
    Adaptation kinetics in bacterial chemotaxis
    1. SM Block
    2. JE Segall
    3. HC Berg
    (1983)
    Journal of Bacteriology 154:312–323.
  6. 6
  7. 7
  8. 8
    A primary requirement for nodal in the formation and maintenance of the primitive streak in the mouse
    1. FL Conlon
    2. KM Lyons
    3. N Takaesu
    4. KS Barth
    5. A Kispert
    6. B Herrmann
    7. EJ Robertson
    (1994)
    Development 120:1919–1928.
  9. 9
  10. 10
  11. 11
  12. 12
  13. 13
    Developmental Biology
    1. SF Gilbert
    (2014)
    Sinauer Associates Incorporated.
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
    The organizer of the mouse gastrula is composed of a dynamic population of progenitor cells for the axial mesoderm
    1. SJ Kinder
    2. TE Tsang
    3. M Wakamiya
    4. H Sasaki
    5. RR Behringer
    6. A Nagy
    7. PP Tam
    (2001)
    Development 128:3623–3634.
  22. 22
  23. 23
  24. 24
  25. 25
    Synergy with tgfβ ligands switches WNT pathway dynamics from transient to sustained during human pluripotent cell differentiation
    1. J Massey
    2. Y Liu
    3. O Alvarenga
    4. T Saez
    5. M Schmerer
    6. A Warmflash
    (2019)
    PNAS, 10.1101/406306.
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
    Ilastik: interactive learning and segmentation toolkit
    1. C Sommer
    2. C Straehle
    3. U Kothe
    4. FA Hamprecht
    (2011)
    8th IEEE International Symposium on Biomedical Imaging (ISBI 2011). pp. 230–233.
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
    Concentration-dependent patterning of the Xenopus ectoderm by BMP4 and its signal transducer Smad1
    1. PA Wilson
    2. G Lagna
    3. A Suzuki
    4. A Hemmati-Brivanlou
    (1997)
    Development 124:3177–3184.
  42. 42
  43. 43
  44. 44
  45. 45
  46. 46
  47. 47
  48. 48

Decision letter

  1. Roel Nusse
    Reviewing Editor; Stanford University, United States
  2. Naama Barkai
    Senior Editor; Weizmann Institute of Science, Israel
  3. Kyle M Loh
    Reviewer; Stanford University, United States
  4. Roy Wollman
    Reviewer; University of California, San Diego, United States

In the interests of transparency, eLife includes the editorial decision letter, peer reviews, and accompanying author responses.

[Editorial note: This article has been through an editorial process in which the authors decide how to respond to the issues raised during peer review. The Reviewing Editor's assessment is that all the issues have been addressed.]

Acceptance notification:

The manuscript addresses the important question whether hESCs respond to changing concentrations of growth factors, in particular how BMP4 and Nodal control the various cell-fate decisions that hESCs take. It is shown that the response to BMP4 depends on the absolute concentration of the factor. In contrast, the response to the BMP-related Nodal signal depends on the rate of change in concentration. These findings will guide other investigators to consider the dynamics of signaling pathways and their concentrations in experiments manipulating stem cells differentiation.

Decision after peer review:

Thank you for submitting your article "Rapid changes in morphogen concentration control self-organized patterning in human embryonic stem cells" for consideration by eLife. Your article has been reviewed by three peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Naama Barkai as the Senior Editor. The following individuals involved in review of your submission have agreed to reveal their identity: Kyle M Loh (Reviewer #1) and Roy Wollman (Reviewer #3). Reviewer #2 remains anonymous.

By measuring the dynamics of an endogenously tagged Smad4 fluorescent reporter, the authors provide evidence for adaptive responses of theTGF-β pathway to Nodal/Activin morphogens, They demonstrate that hESCs at sparse densities can sense the rate of change of ligand concentration. This work presents a quantitative approach to an important developmental pathway and should be of broad interest.

Major concerns:

The reviewers make a number of suggestions for further improvement of the manuscript. The comments are mostly restricted to the interpretation of the experiments and would require few if any new experimental data. We suggest taking these comments into account in resubmitting this work to eLife.

Separate reviews (please respond to each point):

Reviewer #1:

The temporal kinetics with which signaling pathways operate have generally been overlooked. Previous work by Warmflash and colleagues in cultured C2C12 cells suggested that TGFβ signaling is transient (i.e., cells respond to TGFβ ligands only for several hours) and adaptive (i.e., cells do not sense the absolute concentration of TGFβ ligands but rather they sense the rate of change, with large and sudden changes in TGFβ concentrations leading to the strongest responses) (Warmflash et al., 2012 and Sorre et al., 2014). However, whether C2C12 cells are a physiological model is questionable. In recent work with human embryonic stem cells (hESCs) it seems that a closely related signaling pathway, BMP, shows quite different dynamics: BMP4 treatment of hESCs leads to sustained (not transient) BMP signaling for many hours (Nemashkalo et al., 2017).

In the present submission, Heemskerk et al. perform a side-by-side comparison of the temporal dynamics of TGFβ and BMP signaling in hESCs. Monitoring nuclear translocation of SMAD4 by live imaging, they confirm that BMP4 treatment leads to sustained SMAD4 nuclear translocation in hESCs (thus reproducing their prior results reported in Nemashkalo et al., 2017) whereas ACTIVIN treatment leads to transient SMAD4 nuclear localization (Figure 1B-D). Indeed, many BMP target genes show sustained upregulation in response to BMP4, whereas TGFβ target genes show a mixed pattern of either transient or sustained upregulation in response to ACTIVIN (Figure 2). Using sophisticated microfluidic methods to control the rate of ligand delivery, they suggest that TGFβ (but not BMP) signaling relies on the change of ligand concentration (i.e., its time-derivative but necessarily the absolute concentration; Figure 3) and further that there is a traveling wave of TGFβ activity in BMP4-treated hESC colonies of defined sizes (Figure 4). By virtue of the transient and adaptive nature of TGFβ signaling, the authors conclude by suggesting that temporal pulses of ACTIVIN treatment might be best to efficiently upregulate BRACHYURY (Figure 4D).

Overall, the suggestion that TGFβ and BMP signaling kinetics (on the scale of hours) could be quite different is intriguing. The Discussion is excellently written, and appropriately emphasizes that a mechanistic understanding of the kinetics with which these pathways signal could greatly enrich both stem cell and developmental biology.

Major critiques:

1) SMAD4. The authors measure the activity of the TGFβ as well as BMP pathways by principally assessing a single metric: SMAD4 nuclear localization. However it is controversial whether SMAD4 localization as a single metric can accurately describe the total levels of TGFβ and BMP pathway activity. There are some reports that only a subset of TGFβ target genes require SMAD4 (e.g., Levy and Hill, 2005). If there is SMAD4-independent TGFβ signaling, it is therefore plausible that certain TGFβ target genes could still be transcribed (i.e., the TGFβ pathway is activated) in the window of time when SMAD4 is not in the nucleus (i.e., in the "adaptation" phase that the authors identify). If so, SMAD4 nuclear localization might only partially readout TGFβ signaling, which could therefore complicate the interpretation of the authors' data. Indeed, the authors show that upon ACTIVIN exposure for 5 hours (a timepoint by which SMAD4 has largely vacated the nucleus; Figure 1D), a number of TGFβ target genes (LEFTY1, NODAL, CER1 and WNT3) still continue to be transcribed (Figure 2D), indicating the possibility that SMAD4 nuclear localization does not completely reflect the degree/extent of TGFβ signaling. Could the authors perform luciferase reporter assays of ACTIVIN-treated hESCs to examine whether TGFβ transcriptional responses are truly transient? It would be preferable for the authors to perform live imaging of endogenously-expressed, labeled SMAD2/3 or SMAD1/5/8 proteins to examine their subcellular localization upon ACTIVIN or BMP treatment, respectively, but this is unreasonable to ask the authors to conduct in a short span of time. Barring those experiments, the authors might want to confine their conclusions to refer to "SMAD4-dependent TGFβ and BMP signaling" as opposed to "TGFβ and BMP signaling". For instance, would it be more accurate to specifically say that SMAD4-dependent TGFβ signaling (as opposed to TGFβ signaling altogether) is transient and adaptive?

2) SMAD4 signal-to-noise ratio. Their live imaging data suggests that upon ligand stimulation, the increase in nuclear SMAD4 is modest (2-fold or less; Figure 1D). This is rather surprising, because one might anticipate a massive change in SMAD4 levels is needed to drive large transcriptional changes. The modest SMAD4 signal-to-noise ratio also raises questions about the authors' subsequent and quite detailed analysis of temporal dynamics, adaptation and so on. Why is the signal-to-noise ratio so modest? One possibility is that there are high nuclear SMAD4 levels in the steady-state (i.e., there is high SMAD4 background), perhaps because the authors' hESCs are grown in mTeSR1 (Materials and methods), which contains 2 ng/mL of TGFβ1 (Ludwig et al., 2006; Nature Cell Biology), and would be expected to lead to some level of basal nuclear SMAD4. What if the authors transiently pulsed their hESCs with a TGFβ inhibitor before performing the various ACTIVIN and BMP stimulations, would this reduce their background and thus enhance the signal-to-noise ratio?

Minor critiques:

1) Figure 1H-J should be more clearly labeled and explained.

2) "Error bars" erroneously merged into 1 word (Figure Legends).

3) Data in the Nemashkalo et al., 2017, paper (by the same authors) seems to already suggest that BMP4 leads to sustained SMAD4 signaling, therefore it seems appropriate for Heemskerk et al. to cite this earlier work, as their Figure 1 seems to be confirming-and extending beyond-this earlier paper.

4) There seem to be some strange artifacts in certain subpanels. For instance why does the orange line in Figure 3B (reflecting nuclear SMAD4 localization after sudden exposure to high BMP) jitter at regular intervals? It makes sense that the blue lines in Figure 3B,E would jitter (as they might correlate with each step of the ligand concentration increase) but not the orange line in Figure 3B.

5) Similarly, in Figure 4B (red line: pulses), there seem to be discontinuous spikes in SMAD4 whenever the media is changed.

6) The authors conclude that BMP4 induces sustained, whereas ACTIVIN induces transient, SMAD4 signaling, yet in Figure 4 they suggest that BMP4 might induce endogenous TGFβ signaling. One alternate interpretation of the authors' data is that maybe BMP4 itself induces a transient SMAD4 response, but that by secondarily inducing TGFβ signaling, this latter wave of TGFβ might lead to sustained SMAD4.

7) Another oddity is in Figure 3C. If hESCs sense the absolute concentration of BMP4 ligand, one would expect that at low levels of BMP, they would respond weakly and at higher levels of BMP, they would respond strongly. Yet in Figure 3C, at the 2-3 hr timepoint, the "ramp" cells (blue line, exposed to low BMP4 at that point in time) express HAND1 mRNA almost as highly as the "step" cells (orange line, exposed to high BMP4 at that point in time). The only other interpretation I can think of is that maybe the low BMP4 levels in the blue line (at 2-3 hrs) are super-threshold already.

Reviewer #2:

Heemskerk et al. provide evidence for adaptive responses of the TGF-β pathway to Nodal/Activin morphogens. The authors elegantly measured dynamics of an endogenously tagged Smad4 fluorescent reporter that transduces both BMP and Nodal/Activin input signals, and demonstrate that at sparse densities hESCs can sense the rate of change of ligand concentration. Using mathematical models and experimental measurements of endogenously tagged Smad4, they rule out some mechanisms of pathway adaptation, and propose that Nodal/Activin adaptation occurs through reduced nuclear sequestration and/or increased cytoplasmic sequestration. Under sparse growth conditions, they also find that endogenous target genes can selectively respond to pulsatile and sustained modes of Smad4 activation resulting from Activin and BMP inputs, respectively. The authors suggest that sensing the rate of change of ligand concentration, demonstrated in sparse cultures, could be used by the mammalian embryo during mesodermal patterning, and use an in vitro model of human gastrulation to test this hypothesis. In this system, they correlate a moving front of Smad4 nuclear localization with the timing and spatial position of Brachyury expression, a transcription factor involved in mesodermal differentiation.

Overall, this study brings a quantitative perspective to a central developmental pathway, and should help to further strengthen a growing sense that dynamics are critically important in the behavior of many developmental signaling pathways. It suggests that controlling ligand dynamics could facilitate directed differentiation protocols. It should therefore be of broad interest, and we support publication.

Below, we describe several issues that could be addressed to strengthen the work in a revised manuscript:

Major comments:

Figures 1 and Figure 1—figure supplement 1F together show that BMP4 triggers a more stable nuclear localization of BMP-responsive Smad1, while activin results in adaptive nuclear localization of Nodal-responsive Smad2/3. However, the cross-pathway dynamics are unclear. Does BMP4 affect Smad2/3 and does Activin/Nodal affect Smad1, and if so, with what dynamics? This control is important to support the authors claim that the difference in dynamics results from distinct pathway-specific transducers.

Figure 2, B and D, present gene sets that show qualitatively different dynamic response to Activin. The authors suggest the interesting possibility that these sets have different sensitivities to Smad4. However, an alternative explanation could be a difference in RNA half-lives, with the monotonic responses in Figure 2D resulting from stable mRNAs that accumulate over time, effectively integrating the response. Is it possible to discriminate between these explanations (e.g. by examining nascent or intronic RNA), or at least discuss these two possibilities?

Figure 3 presents an elegant comparison between ramps (staircases) and steps of BMP and Activin. The results are striking. However, two issues complicate interpretation. First, the ramp and step differ in total integrated signaling over time, making it difficult, from this single experiment, to disentangle the roles of ligand concentration and dynamics. The simplest way to address this might be to also include a step of half the amplitude, or to do a ramp at double the slope. Second, based on Figure 1 and Figure 1—figure supplement 1B-C, which shows dose-dependent responses to BMP concentration, we would expect the ramp response in Figure 3B to depend on ligand dose and therefore to grow slowly over time. Why does it step up so rapidly?

In Figure 3, rate responsive behavior is demonstrated by comparing a step to a ramp. Does the Nodal pathway act as a quantitative sensor of the rate of increase of Activin, or does it exhibit some kind of rate threshold? This could be answered by showing how response depends quantitatively on the ramp slope. A mathematical model of the adaptation process could help to provide insight into the expected behavior here.

Rate-responsiveness in other systems. Work in bacteria (Young et al., 2013) showed similar rate-responsive behavior in the sigB stress response system, where two kinds of stresses activate sigB in different ways, with one input showing rate-responsive behavior, assayed by temporal "ramps" of different slopes. In that system the mechanism of rate responsiveness has been determined. It would be interesting to discuss the similarity (and differences) in signal processing behaviors between these very different biological systems.

In Figure 4, the authors suggest that during development the moving Smad4 signaling front is responsible for the activation of Brachyury expression. Although strong evidence for adaptative and rate sensing responses to Nodal/ActivinA ligands is provided for sparse hESCs cultures, the micropattern experiment cannot distinguish between a rate sensing mechanism and a pure concentration sensing mechanism. Even though Brachyury responds with higher expression to Activin pulses in sparse hESC cultures, these two experimental conditions differ in both the time scales and the number of pulses cells could perceive (i.e. in sparse culture cells are exposed to multiple Activin pulses separated by ~4.5 hours, whereas in micropattern the data is consistent with only one wave of pSmad2/3, one wave of pSmad1 and one front of Smad4 moving at 1 cell diameter (10µm) per hour, and traveling at least 150 µm to fully cover the spatial domain of brachyury expression. Since no data on dynamics of endogenous Nodal expression is provided, it remains unknown if, over time, ligand concentration is locally increasing (or oscillating) on a timescale relevant to Brachyury induction. As result, in the micropattern, in which responses to BMP and endogenous Nodal co-occur, it is unclear if cells are sensing absolute ligand concentration and/or a rapid increase in ligand concentration. Experimentally determining Nodal dynamics in micropatterns remains difficult, and perturbing its dynamics would require significant cell line engineering. Thus, an alternative avenue to support the claim that in the gastrulating human embryo expression of fate determinants results from rate sensing of ligand increments could be to computationally compare whether patterning triggered by adaptation to Nodal achieved by rate sensing rate qualitatively or quantitatively differs from patterning the mesoderm through sensing absolute ligand concentrations. What new or different emergent properties result from sensing rates as opposed to just concentrations?

Relation to other dynamic encoding systems: This study provides an important advance in understanding the role of dynamics in core mammalian signaling systems. But the relationship to other work on dynamic signal encoding in signaling pathways, as reviewed for example in Purvis and Lahav (2013), Behar and Hoffmann (2010), etc. is minimal. Specifically, it would be useful to know which aspects of the Nodal system in hES resemble or differ from phenomena observed in other pathways or other cell types. For example, different target programs are activated by sustained vs pulsatile activation of p53 and Notch, among other pathways. Also, a recent report (Frick et al., 2017 PNAS) studied the response of C2C12 cells to TGFb ligands, which belong to the Activin/Nodal family and are also transduced by Smads2/3, and showed that the pathway senses fold changes. How does this fold change detection relate to the adaptive rate-responsive behavior observed here? Discussion of these issues would help to connect this study to the broader literature and unify results from multiple studies.

Minor Comments:

In Figure 2, please include the base of the log on the y axis: is it 2, e, or something else?

In Figure 1E there is an orange line near the y-axis. It remains unclear if this line represents the fact that after photobleaching Smad4 fluorescence does not recover to the original amplitude.

Reviewer #3:

The paper by Idse et al. includes fantastic work by the authors to explore fundamental properties of mammalian signaling in an interesting developmental context using an innovative experimental system.

I am very happy to participate in this "peer-review trial in which the authors decide how to respond to the peer-review comments." and to learn that "the article will almost certainly be published".

I do not have any additional experiment to suggest, the authors did a good job providing the adequate controls and I have little doubt in the correctness of the results they present.

Perhaps the biggest suggestion I can make with regards to the work is that it is a bit too dense and hard to grasp. Overall the work follows three parts: 1) An analysis of mechanism for adaptation using innovative combination of FRAP and parameter identification of ODE model. 2) analysis of gene expression dynamics response to dynamic SMAD4, and 3) implication of the differences in response dynamics on overall self-organization patterns in the in-vitro gastrulation model. I will refer to these as parts 1-3 below for brevity.

I am not sure why the authors restricted themselves to only four figures, although I can make an educated guess;). To make the paper easier to follow I strongly recommend to take advantage of eLife lack of restriction on length and number of figures and use the space to elaborate much more about the work.

In part 1 the author only use two paragraphs on main text to explain a lot of stuff. To understand the adaptation and the FRAP experiments one has to think about two different math models that use a parameter fitting procedure to identify underlying mechanisms. The material is mostly in the supp text, but supp text should not replace main text rather provide more details for completeness. I recommend to move some of the material from supp mat to the main text, explain the FRAP procedure, show the alternative fits for different parameters etc. This part can easily take 2-3 figures.

Similarly, part II seemed rushed. The authors show how some genes follow SMAD4 dynamics and some do not. But the sparse text doesn't allow them to go into details, what is different about these genes? Is it RNA halflife? Even without additional experiments, few simple models would be very helpful to understand Figure 2A-D. I am not sure I got the point of panels 2EF, so either explain that in more details or remove those. The different gene expression response to different stimuli regimes (Figure 3) is very nice with its demonstration of the role of SMAD4 dynamics.

Finally, the last part is very short and descriptive, The author measure different SMADs over space and time, show that there are different responses and that these change gene expression. I was left a bit confused as to why these differences arise? What aspect of self-organization is in play to cause these patterns. The fact that all these results are described in a single paragraphs of 322 words might be part of the problem.

Minor Comments:

1) The use of SMAD4 signaling in figure titles should be replaced with y-axis label saying SMAD4_{nuc}

2) The cartoon in Figure 2I is confusing, why don't all the no sequestered dots recover?

https://doi.org/10.7554/eLife.40526.031

Author response

Reviewer #1:

[…] Overall, the suggestion that TGFβ and BMP signaling kinetics (on the scale of hours) could be quite different is intriguing. The Discussion is excellently written, and appropriately emphasizes that a mechanistic understanding of the kinetics with which these pathways signal could greatly enrich both stem cell and developmental biology.

Major critiques:

1) SMAD4. The authors measure the activity of the TGFβ as well as BMP pathways by principally assessing a single metric: SMAD4 nuclear localization. However it is controversial whether SMAD4 localization as a single metric can accurately describe the total levels of TGFβ and BMP pathway activity. There are some reports that only a subset of TGFβ target genes require SMAD4 (e.g., Levy and Hill, 2005). If there is SMAD4-independent TGFβ signaling, it is therefore plausible that certain TGFβ target genes could still be transcribed (i.e., the TGFβ pathway is activated) in the window of time when SMAD4 is not in the nucleus (i.e., in the "adaptation" phase that the authors identify). If so, SMAD4 nuclear localization might only partially readout TGFβ signaling, which could therefore complicate the interpretation of the authors' data. Indeed, the authors show that upon ACTIVIN exposure for 5 hours (a timepoint by which SMAD4 has largely vacated the nucleus; Figure 1D), a number of TGFβ target genes (LEFTY1, NODAL, CER1 and WNT3) still continue to be transcribed (Figure 2D), indicating the possibility that SMAD4 nuclear localization does not completely reflect the degree/extent of TGFβ signaling. Could the authors perform luciferase reporter assays of ACTIVIN-treated hESCs to examine whether TGFβ transcriptional responses are truly transient?

We thank the reviewer for this suggestion. We have attempted to perform luciferase assays using two different TGFb responsive promoters described in the literature (3TP and CAGA), however, neither showed a strong response to Activin treatment in hESCs, indicating that there are likely cell-type specific cofactors necessary for these promoters. In agreement with this, we found by qPCR that the PAI1 gene, from which the 3TP promoter was derived, is also not induced by Activin treatment in hESCs. We agree with the reviewer that this is an important issue. We believe there are two possible, not mutually exclusive, interpretations of the gene expression data. (1) Most transcription is Smad4 dependent but low levels of Smad4 are sufficient for some genes. These genes are maintained by the baseline in Smad4 following adaptation and yield sustained transcription. (2) Only a subset of genes are Smad4 dependent, while others can be induced by Smad2 alone. Since Smad2 is more sustained, these genes are sustained. Of course, one of the explanations may be true for some genes and the other for other genes. We attempted to distinguish between these possibilities by using CRISPR/Cas9 to knock out Smad4 but we found that we were unable to maintain pluripotent Smad4-KO cells. We believe that the dataset presented accurately captures the fact that there are multiple different kinetic behaviors of downstream target genes. We have expanded our discussion to better reflect the different possible explanations.

It would be preferable for the authors to perform live imaging of endogenously-expressed, labeled SMAD2/3 or SMAD1/5/8 proteins to examine their subcellular localization upon ACTIVIN or BMP treatment, respectively, but this is unreasonable to ask the authors to conduct in a short span of time. Barring those experiments, the authors might want to confine their conclusions to refer to "SMAD4-dependent TGFβ and BMP signaling" as opposed to "TGFβ and BMP signaling". For instance, would it be more accurate to specifically say that SMAD4-dependent TGFβ signaling (as opposed to TGFβ signaling altogether) is transient and adaptive?

We agree with the reviewer that it would be of interest to perform live cell imaging with Smad1 and Smad2, but that these experiments are beyond the scope of the present paper. We note that immunofluorescence data for the activity of these signal transducers are presented in the supplementary data. As suggested, we also now explicitly clarify at the beginning of the results that we are referring to SMAD4-dependent TGFb and BMP signaling, but as this language is cumbersome, we use the shorthand “TGFb and BMP signaling” throughout. Finally, we believe that the discussion we incorporated on the dynamics of Smad2 vs Smad4 (see previous answer) also highlights this important point.

2) SMAD4 signal-to-noise ratio. Their live imaging data suggests that upon ligand stimulation, the increase in nuclear SMAD4 is modest (2-fold or less; Figure 1D). This is rather surprising, because one might anticipate a massive change in SMAD4 levels is needed to drive large transcriptional changes. The modest SMAD4 signal-to-noise ratio also raises questions about the authors' subsequent and quite detailed analysis of temporal dynamics, adaptation and so on. Why is the signal-to-noise ratio so modest? One possibility is that there are high nuclear SMAD4 levels in the steady-state (i.e., there is high SMAD4 background), perhaps because the authors' hESCs are grown in mTeSR1 (Materials and methods), which contains 2 ng/mL of TGFβ1 (Ludwig et al., 2006; Nature Cell Biology), and would be expected to lead to some level of basal nuclear SMAD4. What if the authors transiently pulsed their hESCs with a TGFβ inhibitor before performing the various ACTIVIN and BMP stimulations, would this reduce their background and thus enhance the signal-to-noise ratio?

We agree with the reviewer that the measured changes in Smad4 are likely smaller than the increase in the active Smad complexes which mediate transcription. Background stimulation from the TGFb1 in mTeSR1 does not explain the small change as adding the Activin inhibitor reduced the baseline Smad2 levels but not those of Smad4. This can be seen in Figure 1—figure supplement 1G, where the data are scaled so that the levels with SB are 1. The higher levels of SMAD2 compared to SMAD4 prior to Activin treatment show that SB reduces these baseline levels. We believe that the relatively small increase in nuclear Smad4 results from the fact that Smads are present in the nucleus even in the absence of stimulation, as under all conditions Smads shuttle between the nucleus and cytoplasm. So that even when the pathway is inactive, there is some Smad in the nucleus (see Schmierer, 2005). Also, background fluorescence which is present in both the nucleus and the cytoplasm will tend to bring the ratios closer to 1 both before and after stimulation, resulting in reduced fold changes. In any event, the nuclear to cytoplasmic ratio will be a monotonically increasing function of the total Smad4 in the nucleus, and so is a suitable reporter for Smad levels. We chose to use this ratio instead of absolute levels as it tends to remove noise associated with variations in Smad expression, cell size, and imaging artifacts, and we have found empirically that it is a more reliable and reproducible metric for signaling. We point out that the peak to baseline signaling ratio is not the same thing as signal to noise, which would be either of those levels relative to the standard deviation around that level. The significance of the difference between the two signaling states would be best captured by a statistical test such as Kolmogorov-Smirnov and the differences before and after stimulation are highly significant.

Minor critiques:

1) Figure 1H-J should be more clearly labeled and explained.

We have substantially expanded the discussion of the FRAP experiments and split these off into a separate figure in order to clarify.

2) "Error bars" erroneously merged into 1 word (Figure Legends).

This has been corrected.

3) Data in the Nemashkalo et al., 2017, paper (by the same authors) seems to already suggest that BMP4 leads to sustained SMAD4 signaling, therefore it seems appropriate for Heemskerk et al. to cite this earlier work, as their Figure 1 seems to be confirming-and extending beyond-this earlier paper.

We have added this citation.

4) There seem to be some strange artifacts in certain subpanels. For instance why does the orange line in Figure 3B (reflecting nuclear SMAD4 localization after sudden exposure to high BMP) jitter at regular intervals? It makes sense that the blue lines in Figure 3B,E would jitter (as they might correlate with each step of the ligand concentration increase) but not the orange line in Figure 3B.

We believe that these artifacts result from slight perturbations to the imaging due to the media replacement. Although we could remove these by smoothing, we prefer to present the data as is without manipulation. We have added a note to figure caption explaining these artifacts to avoid confusion.

5) Similarly, in Figure 4B (red line: pulses), there seem to be discontinuous spikes in SMAD4 whenever the media is changed.

See answer to point 4.

6) The authors conclude that BMP4 induces sustained, whereas ACTIVIN induces transient, SMAD4 signaling, yet in Figure 4 they suggest that BMP4 might induce endogenous TGFβ signaling. One alternate interpretation of the authors' data is that maybe BMP4 itself induces a transient SMAD4 response, but that by secondarily inducing TGFβ signaling, this latter wave of TGFβ might lead to sustained SMAD4.

To exclude this possibility we have examined treatment with BMP together with a TGFb inhibitor and Activin together with a BMP inhibitor. These data support our conclusion that BMP signaling is sustained while Activin signaling is transient (Figure 1G,H).

7) Another oddity is in Figure 3C. If hESCs sense the absolute concentration of BMP4 ligand, one would expect that at low levels of BMP, they would respond weakly and at higher levels of BMP, they would respond strongly. Yet in Figure 3C, at the 2-3 hr timepoint, the "ramp" cells (blue line, exposed to low BMP4 at that point in time) express HAND1 mRNA almost as highly as the "step" cells (orange line, exposed to high BMP4 at that point in time). The only other interpretation I can think of is that maybe the low BMP4 levels in the blue line (at 2-3 hrs) are super-threshold already.

We thank the reviewer for pointing this out, and we discovered an error in our previous data – the earlier ramp was performed in steps of 1 ng/ml to 10 ng/ml, whereas we intended to use data from a ramp to 3 ng/ml in steps of 0.3 ng/ml. This slower ramp does show a slightly slower increase, particularly in the SMAD4 signaling. Nonetheless, the reviewer is correct that the ramp reaches the levels of the step quickly. The explanation for this is that the BMP dose response (now in Figure 1E) shows that there is a relatively sharp threshold in BMP signaling over which full activation is achieved. This agrees with our previous data (Nemashkalo et al., 2017), as well as from other labs who have suggested that BMP responses are close to binary (Etoc et al. Dev Cell 2016). Thus, it only takes a few steps in the ramp to reach full activation of the pathway. We note that the dose of BMP will have a graded effect on integrated signaling due to the decay in time of signaling at low doses (Figure 1E). We have revised our discussion of this experiment to reflect this, and of the BMP dose-response to better highlight the sharp transition between off and on.

Reviewer #2:

[…] Overall, this study brings a quantitative perspective to a central developmental pathway, and should help to further strengthen a growing sense that dynamics are critically important in the behavior of many developmental signaling pathways. It suggests that controlling ligand dynamics could facilitate directed differentiation protocols. It should therefore be of broad interest, and we support publication.

Below, we describe several issues that could be addressed to strengthen the work in a revised manuscript:

Major comments:

Figures 1 and Figure 1—figure supplement 1F together show that BMP4 triggers a more stable nuclear localization of BMP-responsive Smad1, while activin results in adaptive nuclear localization of Nodal-responsive Smad2/3. However, the cross-pathway dynamics are unclear. Does BMP4 affect Smad2/3 and does Activin/Nodal affect Smad1, and if so, with what dynamics? This control is important to support the authors claim that the difference in dynamics results from distinct pathway-specific transducers.

We thank the reviewer for this suggestion. To address this issue, we have performed additional experiments in which each pathway is stimulated while the other is inhibited. The results show that BMP has no effect on the dynamics downstream of Activin while Activin/Nodal signaling has a small negative effect on the response to BMP so that inhibiting Activin during BMP signaling gives a stronger response. These data are now included as Figure 1G,H.

Figure 2, B and D, present gene sets that show qualitatively different dynamic response to Activin. The authors suggest the interesting possibility that these sets have different sensitivities to Smad4. However, an alternative explanation could be a difference in RNA half-lives, with the monotonic responses in Figure 2D resulting from stable mRNAs that accumulate over time, effectively integrating the response. Is it possible to discriminate between these explanations (e.g. by examining nascent or intronic RNA), or at least discuss these two possibilities?

To examine this possibility, we have performed gene expression experiments in which we inhibited Activin/Nodal signaling and examined target gene expression in time. If the differences in gene expression dynamics were due to differences in mRNA half-life, we would expect that transcripts that are stably induced in response to Activin stimulation would persist following pathway inhibition. However, we find that targets including Lefty, Nodal, and Nanog which are activated in a sustained manner rapidly decline upon inhibition of the pathway with the small molecule inhibitor SB431542. Thus, different mRNA stabilities cannot explain our gene expression results. This data is now included in Figure 3—figure supplement 1G and discussed in the text. In addition, while not impossible, it would be surprising if the different transcriptional dynamics of HAND1 in response to BMP4 vs Activin in Figure 3C and of Noggin in the supplementary data for Figure 3 were due to differential modulation of RNA life time as it is the same transcript being induced. The more natural interpretation is that this directly reflects SMAD4 dependent transcriptional activity.

Figure 3 presents an elegant comparison between ramps (staircases) and steps of BMP and Activin. The results are striking. However, two issues complicate interpretation. First, the ramp and step differ in total integrated signaling over time, making it difficult, from this single experiment, to disentangle the roles of ligand concentration and dynamics. The simplest way to address this might be to also include a step of half the amplitude, or to do a ramp at double the slope.

We now include data on a faster ramp, as suggested, and find that the response is still strongly reduced from step stimulation, however it is stronger than the slower ramp. These data support our conclusion that the time-derivative of concentration is the key determinant of the response rather than the integrated signal and are now presented in Figure 4G,H. We also point out that signaling responses are coincident by the end of the experiment, so that integrated signaling response to the ramp will never catch up with integrated signaling response to step, even if it were allowed to run longer to get to the same integrated ligand exposure.

Second, based on Figure 1 and Figure 1—figure supplement 1B-C, which shows dose-dependent responses to BMP concentration, we would expect the ramp response in Figure 3B to depend on ligand dose and therefore to grow slowly over time. Why does it step up so rapidly?

The BMP dose response (now in Figure 1E) shows that there is a relatively sharp threshold in BMP signaling over which full activation is achieved. This agrees with our previous data (Nemashkalo et al., 2017), as well as from other labs who have suggested that BMP responses are close to binary (Etoc et al., Dev Cell 2016). Thus, it only takes a few steps in the ramp to reach full activation of the pathway. We have revised our discussion of this experiment to reflect this, and of the BMP dose-response to better highlight the sharp transition between off and on.

In Figure 3, rate responsive behavior is demonstrated by comparing a step to a ramp. Does the Nodal pathway act as a quantitative sensor of the rate of increase of Activin, or does it exhibit some kind of rate threshold? This could be answered by showing how response depends quantitatively on the ramp slope. A mathematical model of the adaptation process could help to provide insight into the expected behavior here.

We now include data showing intermediate rates of increase which yield intermediate responses in the Smad4 signaling (Figure 4G,H). We also note that the supplement includes a model of the adaptation process that generates behavior consistent with these experiments (Supplemental text).

Rate-responsiveness in other systems. Work in bacteria (Young et al., 2013) showed similar rate-responsive behavior in the sigB stress response system, where two kinds of stresses activate sigB in different ways, with one input showing rate-responsive behavior, assayed by temporal "ramps" of different slopes. In that system the mechanism of rate responsiveness has been determined. It would be interesting to discuss the similarity (and differences) in signal processing behaviors between these very different biological systems.

We now cite this paper in the discussion, however, as our focus is on adaptation in developmental signaling, we feel a detailed discussion of this paper is beyond the scope of our manuscript.

In Figure 4, the authors suggest that during development the moving Smad4 signaling front is responsible for the activation of Brachyury expression. Although strong evidence for adaptative and rate sensing responses to Nodal/ActivinA ligands is provided for sparse hESCs cultures, the micropattern experiment cannot distinguish between a rate sensing mechanism and a pure concentration sensing mechanism. Even though Brachyury responds with higher expression to Activin pulses in sparse hESC cultures, these two experimental conditions differ in both the time scales and the number of pulses cells could perceive (i.e. in sparse culture cells are exposed to multiple Activin pulses separated by ~4.5 hours, whereas in micropattern the data is consistent with only one wave of pSmad2/3, one wave of pSmad1 and one front of Smad4 moving at 1 cell diameter (10µm) per hour, and traveling at least 150 µm to fully cover the spatial domain of brachyury expression. Since no data on dynamics of endogenous Nodal expression is provided, it remains unknown if, over time, ligand concentration is locally increasing (or oscillating) on a timescale relevant to Brachyury induction. As result, in the micropattern, in which responses to BMP and endogenous Nodal co-occur, it is unclear if cells are sensing absolute ligand concentration and/or a rapid increase in ligand concentration. Experimentally determining Nodal dynamics in micropatterns remains difficult, and perturbing its dynamics would require significant cell line engineering. Thus, an alternative avenue to support the claim that in the gastrulating human embryo expression of fate determinants results from rate sensing of ligand increments could be to computationally compare whether patterning triggered by adaptation to Nodal achieved by rate sensing rate qualitatively or quantitatively differs from patterning the mesoderm through sensing absolute ligand concentrations. What new or different emergent properties result from sensing rates as opposed to just concentrations?

While we agree in principle that these data in isolation are also consistent with sensing ligand concentration, we note that this experiment does provide an important test, and could have resulted in data disproving our model. The most straightforward hypothesis would be that a static gradient of Nodal signaling is formed, peaks in the region of future mesendoderm, and cells read this gradient. Had this been true, it would have shown that the rate of increase is not relevant to the patterning that occurs. Our data disproved this hypothesis. The fact that there is a rapidly moving wave of signal activation in micropatterned colonies together with our experiments in Figures 1 – 3 (now 1 – 5) that demonstrate that cells respond to the rate of Activin/Nodal increase strongly suggest that cells in micropatterning colonies are also responding to the rate of increase. We agree that direct demonstration would require measurement of Nodal protein dynamics but believe this to be beyond the scope of this manuscript. We also note that the consequences of reading rate of increase during patterning were modeled theoretically in Sorre et al., 2014.

Relation to other dynamic encoding systems: This study provides an important advance in understanding the role of dynamics in core mammalian signaling systems. But the relationship to other work on dynamic signal encoding in signaling pathways, as reviewed for example in Purvis and Lahav (2013), Behar and Hoffmann (2010), etc. is minimal. Specifically, it would be useful to know which aspects of the Nodal system in hES resemble or differ from phenomena observed in other pathways or other cell types. For example, different target programs are activated by sustained vs pulsatile activation of p53 and Notch, among other pathways. Also, a recent report (Frick et al., 2017 PNAS) studied the response of C2C12 cells to TGFb ligands, which belong to the Activin/Nodal family and are also transduced by Smads2/3, and showed that the pathway senses fold changes. How does this fold change detection relate to the adaptive rate-responsive behavior observed here? Discussion of these issues would help to connect this study to the broader literature and unify results from multiple studies.

We now cite and briefly discuss the work on p53 and NFkB mentioned above, however, as our focus is on adaptation in developmental signaling, we feel a detailed discussion is beyond the scope of our manuscript. Our results on ramps and pulses are not consistent with fold change detection and so we do not discuss this in the paper.

Minor Comments:

In Figure 2, please include the base of the log on the y axis: is it 2, E, or something else?

It is base 2. We have now clarified this in the figure legend.

In Figure 1E there is an orange line near the y-axis. It remains unclear if this line represents the fact that after photobleaching Smad4 fluorescence does not recover to the original amplitude.

The orange line represents data from the photobleaching experiment. The drop near the y-axis is the actual bleaching event. We have now clarified this point in the figure legend.

Reviewer #3:

[…] Perhaps the biggest suggestion I can make with regards to the work is that it is a bit too dense and hard to grasp. Overall the work follows three parts: 1) An analysis of mechanism for adaptation using innovative combination of FRAP and parameter identification of ODE model. 2) analysis of gene expression dynamics response to dynamic SMAD4, and 3) implication of the differences in response dynamics on overall self-organization patterns in the in-vitro gastrulation model. I will refer to these as parts 1-3 below for brevity.

I am not sure why the authors restricted themselves to only four figures, although I can make an educated guess;). To make the paper easier to follow I strongly recommend to take advantage of eLife lack of restriction on length and number of figures and use the space to elaborate much more about the work.

In part 1 the author only use two paragraphs on main text to explain a lot of stuff. To understand the adaptation and the FRAP experiments one has to think about two different math models that use a parameter fitting procedure to identify underlying mechanisms. The material is mostly in the supp text, but supp text should not replace main text rather provide more details for completeness. I recommend to move some of the material from supp mat to the main text, explain the FRAP procedure, show the alternative fits for different parameters etc. This part can easily take 2-3 figures.

We thank the reviewer for these suggestions. We have now expanded the paper to 6 figures, including expanding the discussion in the FRAP section with material from the supplement.

Similarly, part II seemed rushed. The authors show how some genes follow SMAD4 dynamics and some do not. But the sparse text doesn't allow them to go into details, what is different about these genes? Is it RNA halflife? Even without additional experiments, few simple models would be very helpful to understand Figure 2A-D. I am not sure I got the point of panels 2EF, so either explain that in more details or remove those. The different gene expression response to different stimuli regimes (Figure 3) is very nice with its demonstration of the role of SMAD4 dynamics.

We have now included additional data on gene expression that rule out differences in RNA half-life (see response to reviewer 2 above). We have also expanded the discussion to better explain the gene expression data in Figures 2EF (now Figure 3EF).

Finally, the last part is very short and descriptive, The author measure different SMADs over space and time, show that there are different responses and that these change gene expression. I was left a bit confused as to why these differences arise? What aspect of self-organization is in play to cause these patterns. The fact that all these results are described in a single paragraphs of 322 words might be part of the problem.

We have now split this into two figures and expanded the explanation. We hope it is clearer now.

Minor Comments:

1) The use of SMAD4 signaling in figure titles should be replaced with y-axis label saying SMAD4_{nuc}

We now use the label Smad4 (N:C) throughout to indicate that what we are measuring the ratio of nuclear to cytoplasmic Smad4.

2) The cartoon in Figure 2 is confusing, why don't all the no sequestered dots recover?

Fluorescent proteins which have been bleached never recover regardless of sequestration. In the cartoon, the non-sequestered fluorescent dots simply re-equilibrate between the nucleus and cytoplasm while the sequestered ones remain in place. We have expanded the discussion of this figure and hope it is now clearer.

https://doi.org/10.7554/eLife.40526.032

Article and author information

Author details

  1. Idse Heemskerk

    Department of Biosciences, Rice University, Houston, United States
    Present address
    Department of Cell and Developmental Biology, University of Michigan MedicalSchool, Ann Arbor, United States
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Supervision, Funding acquisition, Investigation, Visualization, Writing—original draft, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-8861-7712
  2. Kari Burt

    Department of Biosciences, Rice University, Houston, United States
    Contribution
    Investigation, Writing—review and editing
    Competing interests
    No competing interests declared
  3. Matthew Miller

    Department of Biosciences, Rice University, Houston, United States
    Contribution
    Investigation, Visualization
    Competing interests
    No competing interests declared
  4. Sapna Chhabra

    Systems, Synthetic and Physical Biology Program, Rice University, Houston, United States
    Contribution
    Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-2479-0911
  5. M Cecilia Guerra

    Department of Biosciences, Rice University, Houston, United States
    Contribution
    Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  6. Lizhong Liu

    Department of Biosciences, Rice University, Houston, United States
    Contribution
    Investigation, Writing—review and editing
    Competing interests
    No competing interests declared
  7. Aryeh Warmflash

    1. Department of Biosciences, Rice University, Houston, United States
    2. Department of Bioengineering, Rice University, Houston, United States
    Contribution
    Conceptualization, Supervision, Funding acquisition, Investigation, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    Aryeh.Warmflash@rice.edu
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-5679-2268

Funding

National Science Foundation (MCB-1553228)

  • Aryeh Warmflash

Cancer Prevention and Research Institute of Texas (RR140073)

  • Aryeh Warmflash

Gillson Longenbaugh Foundation (15-0217)

  • Aryeh Warmflash

The Branco Weiss Fellowship – Society in Science

  • Idse Heemskerk

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

We thank Elena Camacho Aguilar for careful reading of the manuscript and Anastasiia Nemashkalo for providing data for Figure 1h. This work was supported by the National Science Foundation (grant MCB-1553228), the Cancer Prevention Research Institute of Texas (grant RR140073), the Longenbaugh Foundation, Rice University, and the Society in Science - Branco Weiss fellowship, administered by ETH Zurich.

Senior Editor

  1. Naama Barkai, Weizmann Institute of Science, Israel

Reviewing Editor

  1. Roel Nusse, Stanford University, United States

Reviewers

  1. Kyle M Loh, Stanford University, United States
  2. Roy Wollman, University of California, San Diego, United States

Publication history

  1. Received: October 1, 2018
  2. Accepted: February 14, 2019
  3. Version of Record published: March 4, 2019 (version 1)

Copyright

© 2019, Heemskerk et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 870
    Page views
  • 151
    Downloads
  • 0
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Cell Biology
    2. Computational and Systems Biology
    Nona Farbehi et al.
    Research Article
    1. Computational and Systems Biology
    2. Developmental Biology
    Erika Tsingos et al.
    Research Article