Stochasticity in the miR-9/Hes1 oscillatory network can account for clonal heterogeneity in the timing of differentiation

Abstract
Introduction
Results
Discussion
Materials and methods
Appendix
References
Article and author information
Metrics

Abstract

Recent studies suggest that cells make stochastic choices with respect to differentiation or division. However, the molecular mechanism underlying such stochasticity is unknown. We previously proposed that the timing of vertebrate neuronal differentiation is regulated by molecular oscillations of a transcriptional repressor, HES1, tuned by a post-transcriptional repressor, miR-9. Here, we computationally model the effects of intrinsic noise on the Hes1/miR-9 oscillator as a consequence of low molecular numbers of interacting species, determined experimentally. We report that increased stochasticity spreads the timing of differentiation in a population, such that initially equivalent cells differentiate over a period of time. Surprisingly, inherent stochasticity also increases the robustness of the progenitor state and lessens the impact of unequal, random distribution of molecules at cell division on the temporal spread of differentiation at the population level. This advantageous use of biological noise contrasts with the view that noise needs to be counteracted.

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

Introduction

The nervous system is a complex organ system with many different cell types. These are generated by actively dividing progenitor cells in a sequential order during embryogenesis, while progenitor cells themselves change their character along the way (Götz and Huttner, 2005). In cases where histogenesis of neural tissues has been examined, such as the mammalian cortex and the vertebrate retina, a temporal order has been found such that certain cell types are generated early while others are generated late (Franco et al., 2012; Bertrand et al., 2002; Boije et al., 2014; Greig et al., 2013; Franco and Müller, 2013; Goetz et al., 2014; Shen et al., 2006; Gaspard et al., 2008). In order to understand neural development it is essential to understand how the timing of differentiation is regulated.

The predominant view is that the timing of differentiation is precise and tightly controlled. Indeed, given that the final developed nervous system is consistently structurally organised, this is an intuitive assumption. Therefore, it was surprising that recent experiments based on single-cell lineage mapping in the zebrafish retina have shown that the timing of differentiation can be a fundamentally stochastic process, although the collective behaviour of many progenitors forms a reproducible retina (He et al., 2012). A similar conclusion was reached in a study where rat retinal progenitor cells were cultured at single-cell density and constant extracellular environment (Gomes et al., 2011). In addition to the retina, lineage tracing has also demonstrated that progenitors produce a temporally ordered but variable number of neurons within the mouse neocortex (Shen et al., 2006; Gaspard et al., 2008) although other studies reported a high degree of determinism in this system (Gao et al., 2014). Nevertheless, inherent stochasticity may be a widespread phenomenon in differentiation.

Previous models aiming to explain this phenomenon, such as the branching models of population dynamics (Rué and Martinez Arias, 2015; Klein and Simons, 2011), take into account 3 types of division: symmetric proliferative, asymmetric and symmetric differentiative. The probabilities of division type are either fixed (e.g. in adult tissue homeostasis; reviewed in Simons and Clevers (Simons and Clevers, 2011)) or change with time (e.g. in developing tissues (Slater et al., 2009)), but are usually fitted to the data (He et al., 2012; Gomes et al., 2011). Such models are useful in understanding population dynamics but lack a molecular mechanistic explanation.

At a molecular level, the balance between differentiation and the maintenance of progenitors in the developing nervous system is regulated by transcriptional repressors such as the Hes gene family, including Hes1 and Hes5 (Kageyama et al., 2007, 2008). When Hes genes are absent, neural progenitors prematurely differentiate and cause a wide range of defects in brain formation (Hatakeyama et al., 2004; Nakamura et al., 2000). Conversely, overexpression of Hes genes leads to inhibition of neurogenesis and over-maintenance of neural progenitors (Ishibashi et al., 1995).

The development of live-cell imaging with unstable Luciferase (LUC) reporters has shown that the dynamics of Hes1 gene expression changes during neural development (Imayoshi et al., 2013). Hes1 is expressed in an oscillatory manner with a period of around 2 hr in neural progenitors (Imayoshi et al., 2013; Shimojo et al., 2008) but is expressed at a low steady state in differentiated neurons (Imayoshi et al., 2013; Sasai et al., 1992). Based on expression dynamics and the functional studies mentioned above, it has been proposed that a Hes1 oscillatory state is necessary for the maintenance of progenitors, while low, non-oscillatory levels are associated with a transition to neuronal differentiation (Kageyama et al., 2008). The most direct evidence for the functional importance of oscillatory dynamics in general, comes from optogenetics studies of the Hes1 targets Ascl1 (Imayoshi et al., 2013) and Dll1 (Shimojo et al., 2016). It was shown that light-induced oscillatory expression of Ascl1 increased the proportion of dividing cells in Ascl1-null neural progenitors, whereas sustained Ascl1 expression increased the efficiency of neuronal differentiation (Imayoshi et al., 2013). It has similarly been shown that light-induced sustained expression of the Delta ligand Dll1 leads to higher levels of the cell cycle inhibitor p21 than oscillatory expression (Shimojo et al., 2016). Together, these suggested that the expression dynamics of Ascl1 and Dll1 encode information for a choice between proliferation and differentiation within neural progenitors. The hypothesis that gene expression dynamics change as cells make cell-state transitions in development is consistent with previous studies in theoretical biology (Furusawa and Kaneko, 2012; Huang, 2011; Garcia-Ojalvo and Martinez Arias, 2012; Rué and Martinez Arias, 2015). In these studies, the interactions of multiple genes in regulatory networks can lead to the emergence of transient stem cell dynamics, which evolve to an attracting stable configuration of gene expression corresponding to distinct cell types.

Experimental and theoretical work has shown that Hes1 oscillates in neural progenitors possibly due to a combination of delayed negative self-repression and relatively fast degradation of Hes1 mRNA and HES1 protein, previously measured in fibroblasts (Jensen et al., 2003; Monk, 2003; Hirata et al., 2002; Momiji and Monk, 2008). Until recently it was not understood how oscillations of Hes1 could be terminated and the timing of differentiation controlled. Recent experimental results have shown that Hes1 is a primary target of the microRNA miR-9, and HES1 also periodically represses the transcription of miR-9, thus forming a double negative feedback loop (Bonev et al., 2012). However, mature miR-9 is very stable and accumulates over time in a gradual manner. It has been proposed that accumulating levels of miR-9 beyond a certain level can cause oscillations of Hes1 to cease, leading to differentiation (Bonev et al., 2012; Tan et al., 2012). Experimentally it has been shown that Hes1 is a target of miR-9, that depleting miR-9 prevents or delays differentiation and that Hes1 changes dynamics of expression as cells differentiate (Bonev et al., 2012; Bonev et al., 2011; Imayoshi et al., 2013; Tan et al., 2012; Coolen et al., 2012). However, a theoretical approach unifying these phenomena was lacking.

The Hes1/miR-9 interaction has been mathematically modelled using deterministic delay differential equations (Goodfellow et al., 2014). This showed that interactions within the Hes1/miR-9 network can lead to the emergence of bistability (alternative cell fates) and Hopf bifurcations (separating stable states from an oscillatory regime). Computationally it was shown that indeed, accumulation of miR-9 can lead to a transition from oscillatory to stable expression of HES1 protein and that the timing of this transition is controlled by the parameters and the initial conditions of the model. Although this mathematical model provided a useful starting point, a deterministic system is not able to explain the stochastic nature of the differentiation process that has been reported in vivo.

Here, we show that high cell-to-cell variability exists in the abundance of Hes1 mRNA and protein in undifferentiated neural cells (progenitor and stem cells) as evidenced by quantitative single-molecule Fluorescent in Situ Hybridisation (smFISH) and Fluorescence Correlation Spectroscopy (FCS), respectively. Furthermore, the copy number of Hes1 mRNA, protein and miR-9 per cell is low, (determined by smFISH, FCS and stem loop qRT-PCR, respectively) supporting the inclusion of 'finite number' intrinsic stochasticity in the computational model of the Hes1/miR-9 interaction. Oscillations produced by a stochastic model match single-cell imaging data well and lead to systematic changes in the dynamics and output of this network. In the deterministic model, equivalent cells (i.e clonal or isogenic) that start from the same initial conditions terminate oscillations i.e. 'differentiate' at exactly the same time. However, in the stochastic framework the time-to-differentiation of such cells becomes distributed, with the spread of the distribution inversely related to the absolute number of interacting molecules. Due to the ability of intrinsic noise to induce oscillatory dynamics, Hes1 oscillations are maintained for a greater range of parameters than in the deterministic system, which can be interpreted as increased robustness of the progenitor state in the stochastic model. Computationally, the average time of differentiation can be shifted by varying the molecular copy number or by varying the initial miR-9 levels. This supports the idea that the number of interacting molecules within a cell is an important parameter that can be tuned to control the onset of differentiation. We validate predictions of our stochastic model by measuring the time distribution of differentiation when miR-9 is added to the system. Finally, we show that the stochastic network is able to better 'absorb' noise attributable to unequal distribution of some molecular components at division, maintaining its output as if the divisions were symmetric and thus giving a robust outcome at a population level.

Together, our work demonstrates that randomness originating in the finite copy number within transcriptional programmes provides a theoretical framework to understand the heterogeneity observed in single-cell time series of neural stem (NS) cells and its consequences for neural development. Developing systems may exploit stochasticity in gene regulatory networks to stagger differentiation time, cope with noise introduced by unequal inheritance of components and safeguard the progenitor state.

Results

Developing a stochastic model for the Hes1/miR-9 oscillator

Previous analysis focussed on a model of the Hes1/miR-9 network using delay differential equations (Goodfellow et al., 2014), which are themselves an extension of earlier models by Jensen et al. (2003) and Monk (2003); Momiji and Monk (2008), as outlined in the Materials and methods.

A deterministic system relies on differential equations, making an assumption that variables such as mRNA and protein can take a continuous range. This is very familiar to the biologist and such equations are commonly used for modelling biological processes (Alon, 2007). However, in reality, biochemical reactions are the result of random encounters between discrete numbers of molecules, and some of these molecules may be present at low numbers. The finite number of molecules interacting within the system leads to inherent randomness known as intrinsic stochasticity (Swain et al., 2002).

Thus, while the units of the deterministic differential equations are concentrations, a stochastic model requires a description in terms of discrete numbers of molecules. To convert from concentration into the number of molecules we use the system size parameter, $Ω$ , which can be seen as a measure of the volume of the biological system (Kaern et al., 2005; Elf and Ehrenberg, 2003). Concentrations and copy numbers (i.e. number of molecules) are then related via

C o p y n u m b e r = Ω \times Concentration

The copy number can increase when the system size (volume) increases and the concentration is constant, or when the system size (volume) is constant and the concentration increases.

In order to investigate the effects of intrinsic noise on the dynamics of differentiation we reformulated the Hes1/miR-9 model in a stochastic framework, which describes the system in terms of its underlying microscopic reactions, namely the production and degradation of Hes1 mRNA, HES1 protein and microRNA miR-9. We use three complementary stochastic modelling approaches: simulation of the Chemical Master Equation (CME) with the Stochastic Simulation Algorithm with delay (dSSA), simulation of the Chemical Langevin Equation (CLE) and the Linear Noise Approximation (LNA) as described fully in the Materials and methods.

Briefly, the CME uses the network of interactions and their associated rates to formulate how each reaction acts to increase or decrease the likelihood of transitioning to a new state. For example, if the number of protein molecules is very high, the next most likely reaction may be a protein degradation event. Other reactions can still be triggered, however, and this leads to the heterogeneity over time.

Single realisations of the CME were simulated computationally using Gillespie’s Stochastic Simulation Algorithm with delay (dSSA) (Gillespie, 1977; Anderson, 2007; Cai, 2007) which is a statistically exact solution of the master equation but can be computationally expensive, particularly at high system size ( $Ω$ ) because every reaction is simulated. Therefore, to investigate different parameters and to understand the behaviour of a population of cells, both of which require a large number of repeats, we used methods that are more computationally efficient such as Chemical Langevin Equations (CLE) and the Linear Noise Approximation (LNA), which add a random fluctuating term by bundling reactions together (Gillespie, 2000; Grima et al., 2011). The Linear Noise Approximation (LNA) is a (semi-)analytical approach that does not require multiple simulations and is hence the fastest of all the methods for scanning parameter space. However, it can be less accurate at low system size (Grima, 2012).

A 'finite number' stochastically simulated network is justified by experimental data

In our conceptual framework, the timing of differentiation corresponds to the time required to reach a low stable HES1 state in the model, brought about by an increase in miR-9 over time (Goodfellow et al., 2014). The increase of miR-9 over time in cultured neural progenitor cells has been shown before (Bonev et al., 2012). Here, using a miR-9 sensor we show that the activity of miR-9 increases specifically in neuronally differentiated cells (Figure 1—figure supplement 1).

In order to justify the finite-number approach in our modelling we used smFISH to determine the absolute copy number of Hes1 mRNA in C17.2 neural progenitor cells and neural stem cells (NS) expressing GFP from the Mapt (tau) locus (Tau-GFP NS cells). Multiple fluorescent probes were designed to target the Hes1 transcript (Figure 1A) and the accuracy of smFISH for absolute quantitation was validated in Figure 1—figure supplement 2. The Hes1 mRNA copy number was highly variable between individual cells grown under identical conditions and also between cell lines (Figure 1B and D), but generally quite low. In C17.2 neural progenitor cells, we detected an average of 72 copies Hes1 mRNA/cell by smFISH (Figure 1C), and an average of 88 copies by qRT-PCR (data not shown), confirming the result. In Tau-GFP NS cells grown under proliferative conditions (GFP negative cells; Figure 1D and E), the average was only 11 copies of Hes1 mRNA/cell by smFISH. The long-tail distribution of mRNA molecules in NS cells may be indicative of slow switching between transcriptionally active and inactive promoter states (Munsky et al., 2012). Other studies using smFISH have reported numbers of mRNA molecules ranging from a few copies to $10^{3}$ (Neuert et al., 2013; Bahar Halpern et al., 2015; Battich et al., 2015). A few spontaneously differentiated cells, identified by high-intensity GFP expression (Figure 1—figure supplement 3), were mostly negative for Hes1 mRNA transcript, indicating that the gene is switched off upon differentiation (average 0.8 copies per cell) (Figure 1D and F). This was confirmed by using an intronic probe for Hes1 which showed transcription in GFP negative cells (Figure 1H, corresponding exonic shown in G, I) but no transcription in GFP positive neuronally differentiated cells (Figure 1J).

Figure 1 with 5 supplements see all

Download asset Open asset

Absolute quantification of mRNA copy number and active transcription sites.

(A) Schematic showing the mouse *Hes1* primary transcript, with 29 smFISH probes targeting exonic sequence (red) and 38 smFISH probes targeting intronic sequence (blue). (B) smFISH for *Hes1* mRNA (white) in C17.2 cells; nuclei are stained with DAPI (blue). (C) Quantification of *Hes1* mRNA number in C17.2 cells (n=103). (D) smFISH for *Hes1* mRNA (white) in Tau-GFP NS cells; nuclei are stained with DAPI (blue), and GFP protein is stained with rabbit anti-GFP > anti-rabbit Alexa Fluor 488 (green). The field contains one GFP positive cell, as determined by GFP immunofluorescence intensity analysis (Figure 1—figure supplement 3). (E) Quantification of *Hes1* mRNA number in GFP negative NS cells (n=253). (F) Quantification of *Hes1* mRNA number in GFP positive neuronal differentiated NS cells (n=49). (G, H) Double smFISH for *Hes1* mRNA (G) and *Hes1* intron (H) in the same Tau-GFP NS cells; arrow indicates an active transcription site, marked by colocalised exonic and intronic *Hes1* smFISH signals in the nucleus. (I, J) Quantification of *Hes1* active transcription sites in GFP negative NS cells (I, n=253) and in GFP positive neuronal differentiated NS cells (J, n=49). Source data contained in Figure 1—source data 1.

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

Figure 1—source data 1 mRNA counts by smFISH in C17.2 and NS cells.: https://doi.org/10.7554/eLife.16118.003
Download elife-16118-fig1-data1-v1.zip
Figure 1—source data 2 Quantification of miR-9 sensor activity.: https://doi.org/10.7554/eLife.16118.004
Download elife-16118-fig1-data2-v1.pzfx
Figure 1—source data 3 Validation of smFISH accuracy.: https://doi.org/10.7554/eLife.16118.005
Download elife-16118-fig1-data3-v1.pzf
Figure 1—source data 4 Quantification of Tau-GFP signal intensity.: https://doi.org/10.7554/eLife.16118.006
Download elife-16118-fig1-data4-v1.pzf
Figure 1—source data 5 Quantification of miR-9 copy number in C17.2 cells by qRT-PCR.: https://doi.org/10.7554/eLife.16118.007
Download elife-16118-fig1-data5-v1.xlsx

For mature miR-9, we used a stem-loop Taqman qRT-PCR and synthetic miR-9 (Figure 1—figure supplement 4) to estimate the absolute copy number of endogenous miR-9, expressed as a population average copy number per cell in C17.2 cells. Mature miR-9 was present in an average copy number of 1100 ± 20 (S.E.M) molecules per cell, significantly lower than other reported microRNAs and surprisingly low considering the hundreds of miR-9 targets (Bonev et al., 2011) (Figure 1—figure supplement 5).

Finally, to quantify the abundance of endogenous HES1 protein per nucleus in NS cells we generated two stable HES1 reporter cell-lines: VENUS:HES1 driven either by the constitutive UbC promoter or the 2.7 kb Hes1 promoter (Figure 2A,B). FCS was used to count the number of VENUS:HES1 fluorophores (see Experimental methods). This gave an estimate of 3660 ± 120 (S.E.M) molecules of VENUS:HES1 per nucleus when driven by the UbC promoter (Figure 2D and Figure 2—figure supplement 1). We next used western blotting to find the ratio of expression between VENUS:HES1 and endogenous HES1 (Figure 2C) and subsequently normalised to the percentage of VENUS:HES1 positive cells in the population. This gave variable ratios between experiments and showed an approximate ratio of 1:0.8 ectopic to endogenous HES1. This suggests an average abundance of endogenous HES1 of 2930 molecules per nucleus which is significantly lower than the median protein abundance of 50,000 copies per cell previously reported from global quantification studies (Schwanhausser et al., 2011). Furthermore, FCS quantification of VENUS:HES1 driven by the Hes1 promoter gave a similar value, that is 2450 ± 190 (S.E.M) molecules per nucleus (Figure 2D). These complementary approaches highlight the low average abundance and high variability of HES1 molecules in the cell.

Figure 2 with 1 supplement see all

Download asset Open asset

Quantifying HES1 protein number with FCS.

Schematic and image of stably infected NS-E cells containing (A) *UbC*-VENUS:HES1 reporter construct and (B) 2.7 kb-*Hes1pr*-VENUS:HES1 reporter construct. VENUS:HES1 signal in yellow. Scale bars show 30 μm. (C) HES1 western blot with LAMIN B1 loading control. Nuclear lysates of NS-E cells stably infected with *UbC*-VENUS:HES1 or control NS-E cells. Both cell-lines treated with DMSO or 10 μM MG132 proteasome inhibitor for 3 hr. MG132 treated cells show increased HES1 and VENUS:HES1. (D) FCS quantification of VENUS:HES1 molecules per nucleus in stable *UbC*-VENUS:HES1 and 2.7 kb-*Hes1pr*-VENUS:HES1 reporter NS-E cell-lines. 236 cells (mean 3660 ± 120 (S.E.M) molecules) and 98 cells (mean 2450 ± 190 (S.E.M) molecules) were analysed in 4 and 2 experiments respectively. Error bars show S.D. FCS data and Western blots at different exposures are contained in Figure 2—source data 1.

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

Figure 2—source data 1 FCS data and Western blots.: https://doi.org/10.7554/eLife.16118.014
Download elife-16118-fig2-data1-v1.xlsx

We then investigated whether a stochastic or deterministic simulation of the Hes1/miR-9 interaction model is a better fit to the experimental data of HES1 protein and mRNA dynamics at the single-cell and population level. We simulated the model of the Hes1/miR-9 network using the default parameters of Goodfellow et al. (Goodfellow et al., 2014) (parameter set 1 – Appendix 1—table 1) in both a deterministic and stochastic framework (using the dSSA algorithm). While oscillations in the deterministic model remain coherent (Figure 3B), the stochastic model shows variability in its dynamics and period (Figure 3C). Single-cell luciferase imaging of primary NS cells containing a Luciferase2-HES1 (LUC2:HES1) fusion BAC reporter (Imayoshi et al., 2013) shows a high degree of noise evident as peak-to-peak variability in amplitude and period. This was also evident in our experiments (Figure 3A and Figure 3—figure supplement 1 and 2 for images and traces of multiple cells). At the population level, analysed by qRT-PCR of serum synchronised C17.2 cells, oscillations in Hes1 mRNA quickly dampen, presumably because cells drift out of synchrony (Figure 3D). The average expression of a population of 2000 stochastic simulations using the dSSA shows faster dampening (Figure 3F) than deterministic simulations (Figure 3E). Thus, the stochastic simulation matches the experimental data better, both at the single-cell and the population level.

Figure 3 with 2 supplements see all

Download asset Open asset

Stochastic simulations of HES1 oscillations match the experimental data better.

(A) Time series of LUC2:HES1 reporter protein expression in a single primary NS cell determined by bioluminescence imaging. (B) Deterministic simulation of *Hes1*/miR-9 network at a single-cell level, plot shows HES1 protein concentration as a function of time. (C) Stochastic simulation of *Hes1*/miR-9 network with the dSSA at a single cell level. (D) q-PCR quantification of endogenous *Hes1* mRNA in C17.2 cells at the population level following serum synchronization. Error bars show S.E.M of 2 technical replicates of 2 biological experiments. (E) Population level deterministic simulation of *Hes1*/miR-9 network through averaging 2000 single-cells, with *Hes1* mRNA shown as a function of time. (F) Population level stochastic simulation of Hes1/miR-9 network with the dSSA and averaging over 2000 cells, with *Hes1* mRNA shown as a function of time. All simulations performed with parameter set 1 (Appendix 1—table 1) with system size $Ω$ = 1 and initial conditions m(0) = p(0) = r(0) = 0. Note that the computational simulations are less damped than the experimental data, presumably because other sources of variability in experimental data lead to greater damping. MATLAB code and experimental data contained in Figure 3—source data1.

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

Figure 3—source data 1 MATLAB code to simulate the Hes1/miR-9 network with a deterministic or stochastic model (with dSSA). Includes time series for experimental data.: https://doi.org/10.7554/eLife.16118.017
Download elife-16118-fig3-data1-v1.zip
Figure 3—source data 2 Multiple time series examples of single-cell experimental data.: https://doi.org/10.7554/eLife.16118.018
Download elife-16118-fig3-data2-v1.pzfx

Stochasticity spreads the timing of differentiation

By measuring the time taken to reach a low HES1 state in multiple simulations we can describe the distribution of timing events expected in a population of cells. The contribution of stochasticity was investigated by varying the system size parameter $Ω$ , which changes the number of molecules in the system. All analysis was performed using parameters for which the amplitude of oscillations is large and the transition to the low HES1 state is very clear, (parameter set 2 – Appendix 1—table 1). Once the HES1 protein expression reaches a low state we observed that it never again exceeds a low threshold of 20x $Ω$ . We therefore defined 'differentiation' as the time taken for HES1 protein to decrease below this threshold.

In deterministic simulations, all cells that start with the same parameters and conditions differentiate at the same time. At a high system size ( $Ω$ =50), corresponding to high molecule number (mRNA∼500; protein and miR-9∼30,000), there is low variability in the stochastic dynamics and oscillations appear regular (Figure 4A). The time taken to reach the low HES1 protein threshold was calculated for 2000 simulations using the dSSA, and the distribution of times to reach the low HES1 state is narrow for $Ω$ =50 (Figure 4B), approaching the deterministic solution.

Figure 4 with 1 supplement see all

Download asset Open asset

The distribution of time-to-differentiation widens at low system sizes.

(A, D, G) Examples of stochastic simulations using the dSSA at $Ω$ =50, 5 and 1, respectively. Blue, mRNA; green, protein; red, miR-9. (B, E, H) The distribution of times to reach the low HES1 protein state of 2000 simulations using the dSSA. (C, F, I) The distribution of times to reach the low HES1 protein state of 2000 simulations using the CLE. Note that even at low system size, the CLE is in excellent agreement with the Gillespie dSSA, which is consistent with previous studies showing that the CLE can be highly accurate when reactions only involve the creation or degradation of one species and act to either increase or decrease the molecule number by one (Grima et al., 2011) and confirms the results. All simulations performed with parameter set 2 (Appendix 1—table 1) and initial conditions m(0) = 20x $Ω$ , p(0) = 400x $Ω$ and r(0) = 0. MATLAB code contained in Figure 4—source data 1.

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

Figure 4—source data 1 MATLAB code for measuring time to differentiation and systematic shift in mean.: https://doi.org/10.7554/eLife.16118.022
Download elife-16118-fig4-data1-v1.zip

As the system size is decreased to $Ω$ =5 (mRNA∼50; protein and miR-9∼3000) and then further to $Ω$ =1 (mRNA∼10; protein and miR-9∼600) molecular fluctuations cause increasingly varying amplitude and period of HES1 oscillations within each simulation (Figures 4D and G, respectively). At the same time, the distribution in the timing of differentiation of 2000 independent simulations becomes increasingly broad as the system size decreases (Figures 4E and H). Similar results were obtained with the CLE (Figures 4C-I). We conclude that the greater random molecular fluctuations caused by decreased molecular copy number acts to increase the uncertainty in the timing to reach a low HES1 protein state in each cell, such that differentiation events become distributed in the population over a longer period of time.

Stochasticity shifts the average onset of differentiation

When we reduced the system size we unexpectedly found that the average of the distribution also became gradually shifted towards earlier times (e.g. compare Figure 4B with Figure 4E and H, where scaling of mean with system size shown in Figure 4—figure supplement 1). This shift indicates that stochasticity also causes a systematic change in timing, whereby the low HES1 state is reached on average earlier at low system sizes, even though the parameters of the model remain exactly the same. The origin of the systematic change in timing can be understood by considering how the production of miR-9 changes in the presence of stochasticity.

Transcription of miR-9 is repressed by HES1 protein, and consequently the production of miR-9 ( $r$ ) increases when HES1 protein levels are low, assuming activation in the absence of repression. The production of miR-9 is therefore intrinsically linked to HES1 protein, and stochastic variability in HES1 protein causes stochasticity in the production of miR-9 (Figure 5A, miR-9 production is blue line). As the repression of miR-9 transcription by HES1 protein is modelled as a nonlinear Hill-type function, the higher incidence of low HES1 troughs at lower system sizes gives rise to large spikes of miR-9 production (compare blue line in Figure 5A with Figure 5B). In turn, this leads to a faster increase of the stable mature miR-9 (red line in Figure 5A and B) and a reduction in the average time to differentiation.

Figure 5

Download asset Open asset

Stochasticity causes systematic changes to the mean time-to-differentiation.

The time course of HES1 protein and miR-9 expression, compared with the production of miR-9 at a system size of (A) $Ω$ =5 or (B) $Ω$ =1 and simulated with the CLE. Green, protein; red, miR-9; blue, production of miR-9. All simulations performed with parameter set 2 (Appendix 1—table 1) and initial conditions m(0) = 20x $Ω$ , p(0) = 400x $Ω$ and r(0) = 0. (C) The distribution of time to reach the low HES1 protein state when $p_{1}$ = 700, $α_{r}$ = 0.25. Black, $Ω$ =50; red, $Ω$ =5, blue, $Ω$ = 1. Note that here differentiation is defined to reach a miR-9 high state (r>170), as the transition to the low protein state is much less well defined. MATLAB code contained in Figure 5—source data 1.

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

Figure 5—source data 1 MATLAB code to simulate Hes1/miR-9 network using the CLE with tracked miR-9 production and noise-induced delayed differentiation.: https://doi.org/10.7554/eLife.16118.025
Download elife-16118-fig5-data1-v1.zip

However, the systematic change in timing caused by stochasticity is not guaranteed to speed up the system. Instead of low troughs in HES1 causing spikes in miR-9 production, in an alternative parameter set the peaks in protein could cause a disproportionate repression in miR-9 production. This occurs, for example, when the transcriptional repression threshold of miR-9 production is roughly twice that of the Hes1 promoter ( $p_{1}$ = 750), as shown in Figure 5C. In this scenario the systematic effect of stochasticity on timing is reversed and the low system size simulations reach the low HES1 state at a later time. Note that the relative effect on timing is smaller, but this may be expected as HES1 is only weakly repressing miR-9 transcription for this parameter setting.

Increased miR-9 concentration shifts timing and reduces the spread of differentiation

The concentration of miR-9 is a key driver of the timing to differentiation of a single cell (single simulation) in the deterministic model (Goodfellow et al., 2014). We therefore tested the effect of altering the concentration of this key component in the stochastic model and in a population-level simulation.

In the stochastic model increasing the amount of miR-9 brought forward and reduced the spread of time to differentiation events in a simulated population (Figure 6A–C). To experimentally validate this theoretical finding we measured the timing and spread of differentiation in a population of C17.2 and murine Tau-GFP NS cells transfected with miR-9 mimic. Increasing miR-9 in C17.2 neural progenitor cells increased the number of Tuj1-positive neurons earlier than in control cultures after several days, without affecting the total cell number (Figure 6—figure supplement 1). The shorter time to differentiate was confirmed by FACS analysis of Tau-GFP NS cells treated again with miR-9 mimic (Figure 6—figure supplement 2). Transfecting the C17.2 cells with a control mimic and measuring the number of Tuj1-positive neuronal cells showed a broad spread in the timing of differentiation (standard deviation = 48 hr, Figure 6D), as expected from a stochastic process. When cells were transfected with 50 nM or 100 nM of miR-9 mimic the average time to differentiation decreased and the spread became narrower in a miR-9 initial concentration dependent way (standard deviation = 24 hr and 17 hr respectively, Figure 6E and F). Thus, both computational prediction and experimental validation suggest that the initial concentration of a single species, miR-9, is sufficient to change the variance in the time to differentiation.

Figure 6 with 2 supplements see all

Download asset Open asset

The spread of time-to-differentiation decreases in both the computational model and experiments with increased initial miR-9.

(A, B, C) The distribution of times to reach the low HES1 state of 3000 simulations using initial conditions of miR-9 r(0) = 0, 40 and 80, respectively. The initial conditions of mRNA and protein were m(0) = 20x $Ω$ , p(0) = 400x $Ω$ . All simulations performed with parameter set 2 (Appendix 1—table 1) with system size $Ω$ = 1 and using the CLE. (D, E, F) Newly differentiated Tuj1 positive C17.2 cells compared to previous time point as a percentage of the population. Cells were transfected with 100 nM control mimic or 50 nM and 100 nM miR-9 mimic, respectively. Cells were cultured in serum-free differentiation conditions for the time shown. 5 images analysed per condition in at least 2 experiments. MATLAB code and experimental data contained in Figure 6—source data 1.

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

Figure 6—source data 1 MATLAB code to simulate timing to differentiation with different initial conditions of miR-9. Includes experimental data for the spread of differentiation of C17.2 cells.: https://doi.org/10.7554/eLife.16118.027
Download elife-16118-fig6-data1-v1.zip
Figure 6—source data 2 Data for the speed up of differentiation in C17.2 cells with added miR-9.: https://doi.org/10.7554/eLife.16118.028
Download elife-16118-fig6-data2-v1.pzfx
Figure 6—source data 3 Data for the speed up of differentiation in NS cells measured with FACS.: https://doi.org/10.7554/eLife.16118.029
Download elife-16118-fig6-data3-v1.pzfx

Stochasticity increases the robustness of the progenitor state

The progenitor state needs to be robust to avoid premature depletion of progenitors during the lengthy course of development. In order to investigate the robustness of the Hes1/miR-9 regulatory network we compared the termination of HES1 oscillations in stochastic and deterministic frameworks. This comparison was carried out systematically in the space of the parameters $r_{1}$ and $m_{1}$ using the linear noise approximation (LNA) (see Appendix for details). These control how much miR-9 is required to repress Hes1 mRNA translation by half, and how nonlinear (sigmoidal) the translational repression is as a function of miR-9. We have previously shown that translational repression of Hes1 mRNA by miR-9 is particularly important in terminating the Hes1 oscillator with a low HES1 protein level (Goodfellow et al., 2014).

Different behaviours between the deterministic and the stochastic model are shown in this parameter space (Figures 7A and B). In the first region, where translational repression is highly nonlinear (high $m_{1}$ ) and translational repression strength is low (high amount of miR-9 required; high $r_{1}$ ), miR-9 never accumulates to the levels required to terminate oscillations (Figure 7A area 1). As a result, deterministic simulations in this parameter region show self-sustained limit cycle oscillations ( $m_{1}$ = 3.5, $r_{1}$ = 900) (Figure 7C). The equivalent stochastic simulations using the dSSA also display sustained periodic behaviour (Figure 7D). This shows that when the translational repressive effect of miR-9 is weak, oscillations are not terminated in either deterministic or stochastic settings.

Figure 7 with 2 supplements see all

Download asset Open asset

Stochasticity increases the parameter space where HES1 oscillates.

(A) The areas of HES1 protein oscillations in the deterministic system, as a function of translation repression Hill co-efficient ( $m_{1}$ ) and translation repression threshold ( $r_{1}$ ). (B) The areas of HES1 protein oscillations in the stochastic system at the steady state as calculated by the LNA. The colour scale represents the coherence of the stochastic oscillations. (C, E, G) Deterministic simulations of the system for point 1 ( $m_{1}$ = 3.5, $r_{1}$ = 900), point 2 ( $m_{1}$ = 3.5, $r_{1}$ = 600) and point 3 ( $m_{1}$ = 1.5, $r_{1}$ = 600), respectively. Green, Protein; red, miR-9. (D, F, H) Stochastic simulations (dSSA) of the system for point 1 ( $m_{1}$ = 3.5, $r_{1}$ = 900), point 2 ( $m_{1}$ = 3.5, $r_{1}$ = 600) and point 3 ( $m_{1}$ = 1.5, $r_{1}$ = 600), respectively. All simulations performed with parameter set 1 (Appendix 1—table 1) and initial conditions m(0) = 20x $Ω$ , p(0) = 400x $Ω$ and r(0) = 0. MATLAB code contained in Figure 4—source data 1.

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

Figure 7—source data 1 MATLAB code to scan parameters for deterministic and stochastic oscillations.: https://doi.org/10.7554/eLife.16118.033
Download elife-16118-fig7-data1-v1.zip
Figure 7—source data 2 MATLAB code to generate power spectra of examples.: https://doi.org/10.7554/eLife.16118.034
Download elife-16118-fig7-data2-v1.zip

In the second region, (Area 2 in the deterministic model; Figure 7A), defined by highly nonlinear translational repression (high $m_{1}$ ) and high translational repression strength (low amount of miR-9 required to repress translation; low $r_{1}$ ), both deterministic and stochastic models produce HES1 protein levels that decrease and reach a steady state (Figure 7E and F, respectively). The power spectrum of the stochastic simulations at $m_{1}$ = 3.5, $r_{1}$ = 600 shows no peak, and hence the coherence (defined in Materials and methods, Figure 7—figure supplement 1) is zero (Figure 7—figure supplement 2A) The resulting behaviour of such a process is aperiodic (Figure 7—figure supplement 2B); therefore even though HES1 protein still shows some fluctuations (Figure 7F), the lack of periodicity allows us to conclude that oscillations are terminated in the stochastic setting. Thus, in this area oscillations are terminated in both models.

Of particular interest is the third region (Area 3 in the stochastic model; Figure 7B), which is characterised by less nonlinear repression (low $m_{1}$ ) and low translational repression (i.e. high amount of miR-9 to repress translation by half, high $r_{1}$ ). In this area of parameter space the periodic behaviour is distinct between the deterministic and stochastic framework. The deterministic solution in this regime exhibits damped protein oscillations that eventually reach a steady state (Figure 7G). In contrast, the time series of a stochastic simulation using the dSSA at this parameter configuration shows that even though miR-9 accumulates over time, it is unable to terminate protein oscillations (Figure 7H). When the system oscillates in the stochastic regime ( $m_{1}$ = 1.5, $r_{1}$ = 600), there is a well-defined peak in the power spectrum (Figure 7—figure supplement 2C). The peak at 0.04 min $^{- 1}$ corresponds to a period of 160 min in the time series (Figure 7—figure supplement 2D). Thus, this region represents an expansion of the range of parameters in which HES1 protein oscillates.

We propose that stochasticity has made the oscillatory state more robust to an increase in miR-9 in the sense that the stochastic HES1 oscillator continues to operate under conditions (i.e. parameter values) where the deterministic equivalent will have stopped. It is still possible for the accumulation of miR-9 to terminate HES1 oscillations in the stochastic model but the parameter space in which this happens is reduced.

Robustness of timing to perturbations from asymmetric inheritance at cell division

Having seen that stochastic oscillations are inherently more robust in that they are observed for an increased parameter space, we also tested whether they would show increased robustness to external perturbations. Such perturbations can arise by the unequal distribution of molecules during cell division and can be modelled by changing the initial conditions of the simulations. We consider two models of cell division: one where molecular constituents are distributed equally between daughter cells (and hence initial conditions are the same at 50% of the parent cell) and another where molecules are randomly inherited.

First, we investigated how Hes1 mRNA is actually partitioned in NS cells that are undergoing fate symmetric divisions (i.e. grown under proliferative divisions). Very few of these divisions partitioned Hes1 mRNA equally and some showed a highly asymmetric distribution (Figure 8A–D). The majority of cell divisions (36/51) were mostly consistent with a binomial mode of division (p-value>0.05, binomial distribution with Bonferroni multiple testing adjusted p-value), with a few outliers. A binomial distribution would be expected if each mRNA molecule is partitioned independently into the two daughter cells and where the probability of inheriting an mRNA is equal to the ratio of the cell areas of the daughter cells. Previous studies in bacteria have also shown that partitioning is approximately binomial (Golding et al., 2005).

Figure 8

Download asset Open asset

*Hes1* mRNA is inherited unequally at cell division.

(A, B, C) smFISH with an exonic *Hes1* probe showing asymmetry in *Hes1* mRNA number (white) between sister NS cells, just prior to their complete separation at cytokinesis. Aurora-B kinase antibody staining (red) marks the midbody at cytokinesis; nuclei are DAPI stained (blue). (D) Quantification of *Hes1* mRNA to cell area ratios between sister cells from 51 NS cell divisions at cytokinesis, from one experiment. Dark blue area represents the 90% confidence bounds that the ratios are drawn from a binomial distribution, where the average amount of mRNA received by each cell is normalised to the relative cell areas. Light blue represents multiple testing Bonferroni corrected 90% confidence bounds that the ratios are drawn from a binomial distribution. mRNA count data contained in Figure 8—source data 1.

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

Figure 8—source data 1 smFISH mRNA counts of cells at cell division.: https://doi.org/10.7554/eLife.16118.038
Download elife-16118-fig8-data1-v1.pzf

Using the binomial distribution we then calculated the impact of partitioning noise into the performance of the HES1 oscillator under deterministic and stochastic conditions. To generate initial conditions we simulated a parent cell deterministically for 12 hr, starting with parameter set 2 and initial conditions mRNA = 20x $Ω$ , protein = 400x $Ω$ , miR-9 = 0. Using these parent cell conditions, we then generated daughter cells assuming either equal 50/50 inheritance (no heterogeneity) or a binomial distribution of inheritance for each species (schematic shown Figure 9A). Figure 9 shows the difference in the time distribution of differentiation assuming equal or binomial partitioning, with deterministic (Figure 9B,C) or stochastic (Figure 9D,E) dynamics.

Figure 9

Download asset Open asset

Stochastic system is robust to the binomial distribution of constituents at cell division.

(A) Schematic to show the difference between symmetric and asymmetric division. (B) The distribution of time to reach the low HES1 state with 50/50 initial conditions using deterministic dynamics as system size of $Ω$ = 50, 1 and 0.2, respectively. (C) The distribution of times to reach the low HES1 state of 3000 simulations with binomial initial conditions using deterministic dynamics at system size of $Ω$ = 50, 1 and 0.2, respectively. (D) The distribution of times to reach the low HES1 state of 3000 simulations with 50/50 initial conditions using the CLE at system size of $Ω$ = 50, 1 and 0.2, respectively. (E) The distribution of times to reach the low HES1 state of 3000 simulations with binomial initial conditions using the CLE at system size of $Ω$ = 50, 1 and 0.2, respectively. All simulations performed with parameter set 2 (Appendix 1—table 1) and initial conditions m(0) = 20x $Ω$ , p(0) = 400x $Ω$ and r(0) = 0. The standard deviation of the distribution is denoted with $σ$ . In the simulations shown, parent cells were simulated deterministically for 12 hr prior to 'division' however, the distributions remained statistically the same (p-value>0.05) when either 11 or 13 hr deterministic simulations prior to division were used (data not shown). MATLAB code contained in Figure 9—source data 1.

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

Figure 9—source data 1 MATLAB code for deterministic and stochastic dynamics with wither 50/50 or binomial division.: https://doi.org/10.7554/eLife.16118.040
Download elife-16118-fig9-data1-v1.zip

With 50/50 inheritance at cell division (i.e. no heterogeneity) and deterministic dynamics there is no randomness in the system and hence all cells 'differentiate' at the same time (Figure 9B). The time is the same for all system sizes, as the number of molecules has no effect on the dynamics. Next, we introduce heterogeneity in the starting conditions of the cell, occurring only through random partitioning and using a binomial distribution. At a high system size, the distribution of time-to-differentiation is very narrow (standard deviation = 8.58) (Figure 9C). This makes sense, because the relative variance of a binomial distribution is lower at high copy number and hence at large molecule number the daughter cells have similar inheritance and therefore similar timing to differentiation. As the system size is lowered ( $Ω$ =1 and 0.2), the variance increases and the variability in the daughter cell initial conditions causes the time-to-differentiation to become more spread out, compared either to the high system size (either mode of partitioning) or the low system size with equal partitioning (standard deviation = 165, p-value<0.05 for Kolmogorov-Smirnov test between $Ω$ = 0.2 and 50) (Figure 9C). This finding denotes a lack of robustness to partitioning errors when one looks at the outcome of the HES1 oscillator at the population level.

We now consider the effect of unequal inheritance in the case of stochastic dynamics. At both high and low system size (low and high intrinsic noise, respectively) there is little difference between the distribution created with 50/50 (constant) or unequal (binomial) initial conditions as indicated by a low value of the Kolmogorov-Smirnov distance between the two distributions and a p-value indicating a non-significant difference (Figure 9D and E).

We conclude that the stochastic system is more robust than the deterministic system to noise introduced at cell division because it produces the same distribution of differentiation at the population level irrespective of whether the cells distribute molecules randomly or equally at cell division, at all system sizes. In other words, the stochastic system shows tolerance to the unequal distribution of species at division and maintains its performance. By contrast, in the deterministic framework the distribution of differentiation events is changed significantly by random partitioning of molecules, and this change increases as the system size is lowered.

Discussion

At the molecular level, gene expression involves discrete reaction steps that result in the birth and death of individual molecules and it is therefore subject to the intrinsic noise generated by the random encounters of a finite number of interacting molecules (Lestas et al., 2010). It is widely assumed that such noise is detrimental to the robust performance of biological systems and indeed, a number of mechanisms to filter noise have been described (Arias and Hayward, 2006; Bahar Halpern et al., 2015). It remains an open question whether biological systems exploit noise, as is the question of whether noise optimization itself is an evolvable trait in development (reviewed in Rao et al., 2002; Kaern et al., 2005; Balázsi et al., 2011). Here, we have considered the influence of intrinsic noise generated by a low molecular number of interacting species, on the output of HES1 protein oscillations that are generated by delayed negative feedback and are tuned by miR-9. Our findings led us to propose that stochasticity generated by intrinsic noise offers a number of advantages to the differentiation process.

Modelling the differentiation process as an exit from the oscillatory state with low HES1 protein (Goodfellow et al., 2014), we found that in the stochastic framework the time taken to exit HES1 oscillations varied from cell to cell, even though they started from identical conditions. Hence, there was a distribution in the timing of differentiation. This is an important feature of the stochastic model, as it predicts that it is possible for an initially equivalent population of cells (e.g. clonally related, with the same transcriptional profile and under identical conditions) to change cellular identity at different times as a result of internal noise alone and not as a consequence of any other external influences; at the population level this would create a distribution in the timing of differentiation (Gomes et al., 2011; Boije et al., 2014). We suggest that it may be advantageous for differentiation events to have a broad distribution, as this could be used to generate cell-type diversity, if coupled with mechanisms of cell-type specification (Qian et al., 2000; Temple, 2001). However, a distribution in timing of differentiation is not only linked to the generation of different cell types, because when one considers the generation of a single particular cell type, there is again a distribution in the timing of birth of this cell type (Gaspard et al., 2008; Gomes et al., 2011). This suggests that a distribution in the timing of differentiation is a widespread property, which is perhaps advantageous in allowing a time window for feedback control such that the right cell numbers are produced.

Our modelling also showed that there is a continuum of behaviours since the spread of the timing distribution is a function of the system size, where high system size gives rise to a narrow distribution and low system size to a broad distribution. We speculate that this may partly explain the existence of differentiation models that are highly deterministic, as for example in Drosophila neuroblasts which generate cells with a rigidly programmed sequence of fates (Pearson and Doe, 2004; Kohwi and Doe, 2013; Rué and Martinez Arias, 2015) and others that are highly stochastic, as the above mentioned example of vertebrate retinal development.

Our findings suggest that the incorporation of stochasticity into the HES1 oscillator model is advantageous because it increases the robustness of the progenitor state to intrinsic and extrinsic noise. The robustness to increased intrinsic noise is demonstrated by the increase in the region of the parameter space within which HES1 oscillates. This happens because in the parameter space near the oscillatory regime (Hopf bifurcation), the deterministic solution would exhibit damped oscillations but noise keeps pushing the system away from the steady state. This is consistent with previous reports that described noise induced oscillations (McKane and Newman, 2005; Barrio et al., 2006; Galla, 2009). However, while previous studies have looked at whether noise can induce oscillations in systems at steady state (statistically stationary), we have used LNA to look at whether initial transient oscillations are eventually terminated. In this way, we showed that miR-9 accumulation can terminate stochastic oscillations if the threshold of miR-9 required for translational repression is sufficiently low, and the Hill coefficient of repression is sufficiently high. This is important, as normal neural development relies on the balance of prolonged progenitor maintenance (i.e. robust HES1 oscillations) and differentiation (i.e. the ability to terminate these).

The outcome of deterministic dynamics is very different in symmetric versus randomly asymmetric divisions, meaning that deterministic dynamics are not robust to heterogeneities introduced at cell division. In particular, the narrow spread of differentiation that is produced by a deterministic molecular network would be impossible to maintain if the distribution of molecules at division is random and the system size low, both of which are experimentally observed. By contrast, the stochastic system has a more robust performance to extrinsic noise introduced by cell-division asymmetries, as the outcome does not change between equal or random division. Interestingly, at low system size the random partitioning of molecules at division produces a distributed time of differentiation, even when oscillations are modelled in a deterministic framework. We interpret this to mean that noise due to cell-division asymmetries has a similar effect to noise due to low system size on the outcome of the HES1 oscillator. This fits well with previous data on stochasticity introduced by partitioning errors (Huh and Paulsson, 2011). Our findings agree that noise introduced in this way is very difficult to distinguish from other sources of noise (Huh and Paulsson, 2011) and both could contribute in the heterogeneity of timing of clonal cells, although the effect does not appear to be additive.

In addition to the spread of differentiation events, the mean time-to-differentiation was also systematically affected in the stochastic framework, which can be explained by the dynamics of HES1 and miR-9 production in conjunction with the parameterisation of the system. This shows that in addition to parameters and initial conditions (Goodfellow et al., 2014), the average time-to-differentiation can also be tuned by the absolute number of molecules reacting in the system through the system size. The mean time-to-differentiation is important as it controls the total cellular number of a developed organ, along with other quantities such as the proliferation rate.

Here, we have assumed that the intrinsic noise is a result of low interacting numbers of molecules in the Hes1/miR-9 network. Indeed, this was experimentally supported by the absolute quantitation of Hes1 mRNA, protein and miR-9. The average but absolute Hes1 mRNA, protein and miR-9 numbers are consistent with a medium to low system size between $Ω$ =5 and $Ω$ =1, with variability observed both between cells and between different cell lines. In particular, while the Hes1 mRNA copy number is close to the reported global median of 17 copies per cell, the HES1 protein abundance is well below the median value of 50,000 copies per cell (Schwanhäusser et al., 2011). The average measured miR-9 copy number (1100 copies) in C17.2 cells is also on the low side of the previously reported range of 10 to more than 30,000 copies per cell (Chen et al., 2005). We consider the abundance of HES1 protein and miR-9 to be low not only compared to average global abundancies, but also considering the number of their predicted molecular targets (3,074 for HES1, (Margolin et al., 2009); and more than 500 for miR-9, (Bonev et al., 2011)). However, the low molecular copy number is not an exclusive source of cell-to-cell heterogeneity. The Hes1 gene is also involved in the Notch signalling pathway of cell-cell communication, and other studies have considered the effect of noise in the Notch pathway (Jenkins et al., 2015). Additionally, spatial models of Notch signalling have considered the effect of stochastic motions of the nucleus within a gradient of Notch signalling molecules (Aggarwal et al., 2016). Such sources of heterogeneity are not mutually exclusive, and future studies may extend the model to include these effects.

In our model there is inherent stochasticity in every reaction controlling the birth and death of each species but the transcription rate can take a continuous range of values. Other studies of the Hes1/Her1 network have considered binary promoter noise, where the transcription rate switches discontinuously between an ON and OFF state (Lewis, 2003; Sturrock et al., 2013; Jenkins et al., 2015). In the promoter ON/OFF model, periods of gene activity are interspersed by transcriptional silence, leading to bursty behaviour (Suter et al., 2011). Bursting models of transcription have been shown to fit the statistics of mRNA distributions measured by smFISH in multiple different systems (Munsky et al., 2012; Neuert et al., 2013; Bahar Halpern et al., 2015), and have also been fitted to live time series of transcription (Harper et al., 2011; Corrigan and Chubb, 2014; Zechner et al., 2014). Other recent studies have explicitly modelled the coupling of the Hes1 system to the cell cycle (Pfeuty, 2015a, 2015b), but carry the disadvantage that they contain many parameters that are not reliably known.

The challenge is to integrate these different phenomena, each potentially contributing a different source of noise, with the ultimate goal of connecting mechanistic explanations of developmental gene expression programmes with models of population dynamics. Our model provided an experimentally grounded and molecularly explicit theoretical framework to understand the underlying heterogeneity in the differentiation of clonal isogenic stem or progenitor cells, which have been uncovered by lineage tracing experiments in vivo, or in vitro cultures under the same conditions (Rapaport et al., 2004). Our work has also highlighted the need for integrated computational modelling and the value of an absolute quantitation in the experimental approach.

Materials and methods

Computational methods

Request a detailed protocol

\begin{aligned} (2) & \frac{d m}{d t} & = G (p (t - τ)) - S (r) m (t), \\ (3) & \frac{d p}{d t} & = F (r) m (t) - \frac{l n (2)}{μ_{p}} p (t), \\ (4) & \frac{d r}{d t} & = G_{r} (p (t)) - \frac{l n (2)}{μ_{r}} r (t), \end{aligned}

where $m$ , $p$ and $r$ are the concentrations of Hes1 mRNA, HES1 protein and miR-9 respectively. The model parameters $μ_{p}$ and $μ_{r}$ are the half-lives of HES1 protein and miR-9 degradation, respectively. The quantities $G (p (t - τ))$ , $S (r (t))$ , $F (r (t))$ and $G_{r} (p (t))$ are Hill-type functions that control the rates of synthesis and degradation. The variable τ indicates the time delay in the repression of Hes1 transcription by HES1 protein, accounting for the total time taken from mRNA transcriptional initiation to HES1 protein returning to the nucleus to repress transcription. For a given set of initial conditions and parameters the solution of these differential equations is always the same, and they are therefore described as deterministic.

HES1 protein represses Hes1 mRNA production through binding to the promoter, and this is modelled with a monotonically decreasing Hill function. In the model it takes the form

G (p (t - τ)) = \frac{α_{m}}{1 + (p (t - τ) / P_{0})^{n_{0}}},

where $α_{m}$ is the basal transcription rate of Hes1 mRNA in the absence of HES1 protein (previously set to 1; [Goodfellow et al., 2014]).

The Hill function is characterized by two parameters: the repression threshold $p_{0}$ , which represents the amount of protein required to repress transcription to half of basal levels, and the Hill coefficient $n_{0}$ , which controls how steep (sigmoidal) the repression is as a function of protein concentration.

The degradation of Hes1 mRNA is dependent on the levels of miR-9 through the Hill function $S (r)$

S (r) = b_{l} + \frac{(b_{u} - b_{l})}{1 + (r / r_{0})^{m_{0}}},

The function $S (r)$ represents an effective degradation rate of Hes1 mRNA, and $b_{l}$ and $b_{u}$ impose lower and upper bounds for Hes1 mRNA half-life, respectively, which are set using half-lives measured in Bonev et al. (Bonev et al., 2012).

The function $F (r)$ represents translation repression of HES1 protein by miR-9, and the functional form is also a Hill-type function, as previously assumed in Osella et al. (Osella et al., 2011).

F (r) = \frac{1}{1 + (r / r_{1})^{m_{1}}},

Finally, the HES1 protein represses miR-9 production through binding to the promoter, and this is again modelled with a Hill function

G_{r} (p) = \frac{1}{1 + (p / p_{1})^{n_{1}}},

All analyses and simulations were performed using MATLAB 7.12 (RRID:SCR_001622)(MathWorks, US), and code is provided as associated figure source files. For deterministic simulations of the model, we used the MATLAB function ‘dde23’ with options ‘RelTol’ = 10 $^{- 5}$ and ‘AbsTol’ = 10 $^{- 8}$ . For history vectors we used a single value repeated over prior time points. Numerical bifurcation analysis was performed using DDE-BIFTOOL (Engelborghs et al., 2002).

Developing a stochastic model for the Hes1/miR-9 oscillator

Request a detailed protocol

In order to investigate the effects of intrinsic noise on the dynamics of differentiation we reformulated the Hes1/miR-9 model in a stochastic framework, which describes the system in terms of its underlying microscopic reactions. The numbers of molecules are increased or decreased through six reactions, namely the production and degradation of mRNA, protein and microRNA as follows:

\emptyset \overset{G}{⟹} M m R N A t r a n s c r i p t i o n

\emptyset \overset{G_{r}}{⟶} R m i R - 9 t r a n s c r i p t i o n

M \overset{F}{⟶} M + P H e s 1 t r a n s l a t i o n

M \overset{S}{⟶} \emptyset m R N A d e g r a d a t i o n

P \overset{\ln (2) / μ_{p}}{⟶} \emptyset p r o t e i n d e g r a d a t i o n

R \overset{\ln (2) / μ_{r}}{⟶} \emptyset m i R - 9 d e g r a d a t i o n

The transcription reaction involves a delay, indicated by a double arrow, while all other reactions are assumed to occur instantaneously with the rates of reaction indicated above the arrows. The stochastic dynamics of the Hes1/miR-9 reaction system were captured by the chemical master equation (CME) for the number $n_{m}$ of mRNA molecules, $n_{p}$ of protein molecules and $n_{r}$ of miR-9 molecules. The CME describes how the probability of the copy number $n = (n_{m}, n_{p}, n_{r})$ of the different molecular species evolves over time.

The six reactions in the system result in a stochastic process described by the CME

\begin{array}{ll} \frac{d}{d t} P (n; t) = & Ω (E_{m} - 1) [S (n_{r} / Ω) (n_{m} / Ω) P (n; t)] + Ω \frac{l n (2)}{μ_{p}} (E_{p} - 1) [(n_{p} / Ω) P (n; t)] \\ + Ω (E_{p}^{- 1} - 1) [F (n_{r} / Ω) (n_{m} / Ω) P (n; t)] + Ω \sum_{n^{'}} G (n_{p}^{'} / Ω) (E_{m}^{- 1} - 1) [P (n; t; n^{'}, t - τ)] \\ + Ω \frac{\ln (2)}{μ_{r}} (E_{r} - 1) [(n_{r} / Ω) P (n; t)] + Ω (E_{r}^{- 1} - 1) [G_{r} (n_{p} / Ω) P (n; t)] \end{array}

where $P (n; t)$ is the probability of finding the system in state $n$ at time $t$ , and $P (n; t; n^{'}, t^{'})$ is the probability for the system to be in state $n$ at $t$ and in state $n^{'}$ at time $t^{'}$ . $E_{m}$ , $E_{p}$ and $E_{r}$ are raising operators acting on functions of $n_{m}$ , $n_{p}$ and $n_{r}$ via $E_{m} g (n_{m}, n_{p}, n_{r}) = g (n_{m} + 1, n_{p}, n_{r})$ , $E_{p} g (n_{m}, n_{p}, n_{r}) = g (n_{m}, n_{p} + 1, n_{r})$ and $E_{r} g (n_{m}, n_{p}, n_{r}) = g (n_{m}, n_{p}, n_{r} + 1)$ . $E^{- 1}$ stands for the inverse operation. In the case of two-time quantities, e.g. $P (n; t; n^{'}, t^{'})$ the raising applies with respect to the second argument $n^{'}$ .

We have also developed a set of Chemical Langevin Equations (CLE) (Brett and Galla, 2014), where $ζ (t)$ is a rapidly fluctuating random term. The CLE is more computationally efficient at large system sizes, because it does not keep track of every reaction, as they are effectively bundled into time steps.

\begin{aligned} (16) & \frac{d m}{d t} & = Ω (G (p (t - τ)) - S (r) m (t)) + Ω^{1 / 2} ζ_{m} (t), \\ (17) & \frac{d p}{d t} & = Ω (F (r) m (t) - \frac{\ln (2)}{μ_{p}} p (t)) + Ω^{1 / 2} ζ_{p} (t), \\ (18) & \frac{d r}{d t} & = Ω (G_{r} (p (t)) - \frac{\ln (2)}{μ_{r}} r (t)) + Ω^{1 / 2} ζ_{r} (t), \end{aligned}

This random term is assumed to have a Gaussian distribution with zero mean, and the fluctuations added at two consecutive time points are uncorrelated (white noise). The variance of the white noise terms are given by

\begin{array}{ll} ⟨ ζ_{m} (t) ζ_{m} (t^{'}) ⟩ & = (G (p (t - τ)) + S (r) m) δ (t - t^{'}), \\ ⟨ ζ_{p} (t) ζ_{p} (t^{'}) ⟩ & = (F (r) m + \frac{l n (2)}{μ_{p}} p) δ (t - t^{'}), \\ ⟨ ζ_{r} (t) ζ_{r} (t^{'}) ⟩ & = (G_{r} (p (t)) + \frac{l n (2)}{μ_{r}} r) δ (t - t^{'}), \\ ⟨ ζ_{m} (t) ζ_{p} (t^{'}) ⟩ & = ⟨ ζ_{p} (t) ζ_{m} (t^{'}) ⟩ = 0, \\ ⟨ ζ_{m} (t) ζ_{r} (t^{'}) ⟩ & = ⟨ ζ_{r} (t) ζ_{m} (t^{'}) ⟩ = 0, \\ ⟨ ζ_{p} (t) ζ_{r} (t^{'}) ⟩ & = ⟨ ζ_{r} (t) ζ_{p} (t^{'}) ⟩ = 0 \end{array}

The CLE was simulated numerically using an Euler-Maruyama scheme with time steps of 0.1 mins. For history vectors we used a single value repeated over prior time points.

To perform parameter scans, which requires a large number of repeats, we used the Linear Noise Approximation (LNA) because it does not require multiple simulations and is therefore much faster. The LNA is obtained through performing a system size expansion on the CME (Galla, 2009). The system size expansion decomposes the time evolution of each species into a deterministic and stochastic contribution, where the relative amplitude of stochastic fluctuations decreases with the inverse of the square root of the system size (i.e. relative fluctuations decrease with higher molecule number)

n_{i} = ϕ_{i} + Ω^{- 1 / 2} ϵ_{i}

Where $n_{i}$ is the number of molecules of species $i$ , $ϕ_{i}$ is the concentration of species $i$ as determined by the deterministic rate equations, and $ϵ_{i}$ is the random fluctuation with continuous degrees of freedom around the deterministic solution.

After performing the system size expansion, the deterministic equations are recovered when the system size is infinite (no stochasticity). The linear noise approximation (LNA) refers to the next level of approximation from the system size expansion, valid in the limit of large (but not infinite) system size. The LNA describes stochastic fluctuations around the deterministic trajectory. However, there is a trade-off with accuracy, and the system size expansion can be less accurate than either the dSSA or CLE at low system sizes (Grima et al., 2011; Thomas et al., 2012). The full detailed description of the system size expansion for the Hes1/miR-9 model is shown in Appendix. Even though the CLE and LNA are less accurate than the CME at low system sizes, good agreement between them was observed.

Finally, in addition to the assumptions of the deterministic model (Goodfellow et al., 2014), which mostly refer to unknown kinetic parameter values, the stochastic model makes several conceptual assumptions. We assume the chemical system is in well-mixed conditions and the system size is a fixed value. We also assume that the cells are of the same size (volume), and hence the system size parameter is the same within a population of cells and after cell division. When cell division is included in the model, the time of cell division is not constrained and can occur at any point of the HES1 oscillation. Transcription of Hes1 and miR-9 is controlled by negative feedback and we assume that it is constitutively active when not inhibited. Note that both Hes1 mRNA and miR-9 are transcribed from DNA, which is assumed to be at constant concentration, and hence both Hes1 mRNA and miR-9 are produced from the void. In order to make direct comparisons with the previous deterministic description, we do not explicitly model the noise associated with bursting transcription (Suter et al., 2011), which is therefore not explicitly distinguished from other sources of inherent noise.

Comparison of population level simulations with different modes of division

Request a detailed protocol

We simulated the time to differentiation in 3000 cells using the CLE, with different system sizes ranging from 0.2 to 50 and looked at the outcome of 50/50 or binomial partitioning at division. The timing distributions using constant or random initial condition models were compared using the Kolmogorov-Smirnov distance, which quantifies the distance between two empirical distributions and can take values between 0 and 1, with lower values indicating higher similarity,

K S = sup_{t} | F_{b} (t) - F_{e} (t) |

where $F_{b} (t)$ is the cumulative distribution function of the timing with binomial initial conditions, and $F_{e} (t)$ is the cumulative distribution function of the timing with equal 50/50 initial conditions.

Evaluation of oscillatory dynamics

Request a detailed protocol

Central to oscillations is the notion of periodicity. In a stochastic setting, it becomes desirable to be able to distinguish oscillations from random fluctuations, which are dynamic in time but aperiodic. A common method for detecting periodicity is the Fourier transform, which decomposes a given signal into simple sinusoidal oscillating functions. The simple sinusoidal oscillating functions are described by their frequency, and so a Fourier transform takes data in the time domain and represents them in the frequency domain. The Fourier transform is then squared to obtain the power spectrum, and an oscillator is now defined as any process with a peak in its power spectrum at a particular (non-zero) frequency. If the peak is high and sharp, the time series will exhibit high quality oscillations at the dominant frequency with little variability in the peak-to-peak in the time series. However, if the peak in the power spectrum is shallow and broad, there can be significant variation in the peak-to-peak times. When there is no peak in the power spectrum at a non-zero frequency, the process will no longer oscillate and will instead consist of aperiodic fluctuations. A continuum of behaviours is possible, and a suitable metric needs to be defined to measure how sharp the peak in the power spectrum is.

One such measure is the coherence, which is defined as the area under the curve within 20% of the peak frequency as a fraction of the total area under the curve (Alonso et al., 2007), as shown in Figure 7—figure supplement 1. The coherence is easily calculated, since the power spectrum of a stochastic process can be approximated efficiently using the LNA. Unlike the CLE, the LNA is a linear stochastic differential equation, and it is therefore possible to directly calculate the power spectrum of the stochastic dynamics, even with delay (as described in Appendix).

Experimental methods

Cell culture

Request a detailed protocol

C17.2 cells (RRID:CVCL_4511)(07062902, Sigma, UK) were grown in DMEM (ThermoFisher Scientific, UK ) with 10% FBS. Primary NS cells were isolated from dissected forebrains of E13.5–15.5 embryos from LUC2:HES1 BAC reporter mice (RRID:IMSR_RBRC06013) and cultured as previously described (Pollard, 2013). LUC2:HES1 BAC reporter mice were provided by the RIKEN BRC through the National Bio-Resource Project of the MEXT, Japan (Imayoshi et al., 2013). Tau-GFP NS cells were generated from ES cells expressing GFP from the Mapt (tau) locus (Bibel et al., 2004) as previously described (Conti et al., 2005). Tau-GFP ES and NS-E (NS derived from E14 ES cells) cells were a gift from Jennifer Nichols (Cambridge Stem Cell Institute, UK). All reported cell lines tested negative for mycoplasma.

Cell transfection, differentiation and generation of lentivirus reporter cell lines

Request a detailed protocol

C17.2 cells were transfected with Lipofectamine 2000 (ThermoFisher Scientific, UK) and Tau-GFP NS cells with Lipofectamine 3000 (ThermoFisher Scientific, UK) as per manufacturers’ instructions. 24 hr following transfection, C17.2 cells were differentiated by culture in DMEM with 0.2% FBS (Lundqvist et al., 2013) whereas Tau-GFP NS were differentiated towards neurons by gradual removal of the EGF and FGF2 growth factors (Pollard, 2013). miRVana miR-9 mimic and control mimic were obtained from ThermoFisher Scientific, UK. pDsRed-miR9 Sensor was a gift from Lynn Hudson (Addgene plasmid # 22742) (Lau et al., 2008) and control sensor was pDsRed (Clontech, France). pEGFP N1 (Clontech, France) was used as a transfection control plasmid. To generate the VENUS:HES1 reporter constructs the mouse Hes1 coding (Accession Number NM008235) and 3’UTR sequence was synthesised with flanking Gateway® att recombination sites and inserted into pUC57 vector, resulting in the pUC57-Hes1 vector. The Hes1 gene was then transferred to the 3rd generation ‘pLNT-Venus:#’ lentiviral transfer by LR clonase recombination reaction (Bagnall et al., 2015). The resultant vector is termed pLNT UbC-VENUS:HES1 and allows for constitutive expression of N-terminally fused HES1 using Ubiquitin-ligase C promoter (UbC). Additionally, the UbC promoter was replaced with a 2.7 kb sequence of the Hes1 promoter (upstream of the start codon) using Pac1 and Nhe1 restriction enzymes; the derived vector is termed pLNT 2.7 kb-Hes1pr-VENUS:HES1. Lentivirus production and transduction protocols were performed as described before (Bagnall et al., 2015). NS-E -UbC-VENUS:HES1 cells were FACS sorted on BD Influx Cell Sorter (BD Biosciences, UK) to maximize the VENUS:HES1 positive population.

qRT-PCR

Request a detailed protocol

For synchronised population qRT-PCR, C17.2 cells were serum-starved for 24 hr. Following activation with 10% FBS and at time-points of 30 min, cells were lysed and RNA extracted with TRizol (ThermoFisher Scientific, UK). cDNA was prepared using Superscript III (Invitrogen, UK) as per manufacturers’ instructions and qRT-PCR for Hes1 and Gapdh was performed with Taqman gene expression assays (Hes1 - Mm01342805 m1, Gapdh-Mm03302249 g1, ThermoFisher, Scientific, UK). In order to obtain an absolute measurement of miR-9 copy numbers we compared miR-9 content of total RNA extractions of C17.2 cells with a standard curve of known miR-9 dilutions. The standard curve was generated using the following synthetic mature miR-9 sequences (Eurogentec); 23nt UCUUUGGUUAUCUAGCUGUAUGA, 22 nt UCUUUGGUUAUCUAGCUGUAUG, 21 nt UCUUUGGUUAUCUAGCUGUAU, 20 nt UCUUUGGUUAUCUAGCUGUA, 19 nt UCUUUGGUUAUCUAGCUGU (in order of abundance according to miRBase). Synthetic microRNA was 10X serially diluted to enable a standard curve of 1x10 $^{5}$ – 10 copies (5 data points). C17.2 cells were grown in four T25 cm $^{2}$ tissue culture flasks (Corning, UK) in DMEM with 10% FBS to approx. 80% confluency ( 6x10 $^{5}$ cells). One flask was trypsinised (Trypsin EDTA Sigma, UK) and total cell number counted. The remaining flasks were lysed and total RNA extracted using Ambion miR-Vana kit (AM1561, ThermoFisher, UK) as per manufacturers’ instructions. The efficiency of RNA extraction was calculated as 67% and this was incorporated into subsequent copy number calculations. The standard curve samples and 50ng of each total RNA sample were reverse transcribed with the Applied Biosystems TaqMan microRNA Reverse Transcription kit (4366596) as per the manufacturers’ protocol. The primers used for reverse transcription and Taqman qRT-PCR were the stem-loop primer set for mouse miR-9-5p (product 4427975, assay ID 001089, ThermoFisher, UK). TaqMan qRT-PCR was performed using TaqMan Fast Advanced master mix reagent (Applied Biosystems 4444557) on Applied Biosystems Step-ONE plus real time PCR system.

Bioluminescence imaging

Request a detailed protocol

200000 primary LUC2:HES1 NS cells up to passage 15 were plated on laminin (Sigma, UK) coated 35 mm glass-bottom dishes (Greiner-Bio One) and imaged in NS proliferation media containing 1 mM D-luciferin (Promega, UK) using a 40x oil objective on an Olympus UK LV200 inverted bioluminescence microscope. 10-minute exposure and 2x2 binning were used. Dishes were maintained at 37°C in 5% CO2. Bioluminescent movies were analysed on Imaris (RRID:SCR_007370)(Version 7.2.2, Bitplane). Images were first subjected to a 3x3 median filter to remove bright spot artefacts from cosmic rays, then individual cells were tracked manually using the 'spots' function.

Single-molecule Fluorescent In Situ Hybridization (smFISH) and immunofluorescence

Request a detailed protocol

A set of 29 Quasar 570 5’ end labeled 20 nt smFISH probes were designed against Hes1 exonic sequence, and 38 Quasar 670 5’ end labeled 18 nt smFISH probes designed against Hes1 intronic sequence (Biosearch Technologies, US). Tau-GFP NS cells were grown on coverslips in NS cell media (Pollard, 2013) for 72 hr at 37°C with 5% CO2. Coverslips were washed with 1XPBS, fixed for 20 min with 4% formaldehyde in 1XPBS, washed with 1XPBS, and permeabilized in 70% EtOH at 4°C for 1 hr. Cells were then washed for 5 min in 2XSSC+10% deionized formamide. smFISH probes were diluted to 100 nM in 200 ul hybridization buffer per coverslip (2XSSC, 10% deionized formamide, 10% dextran sulphate). Probes were hybridized for 5 hr in a humidified chamber at 37°C, then cells washed in 2XSSC+10% deionized formamide at 37°C for 30 min. Immunofluorescence was performed post hybridization, incubating with primary antibodies at 4°C for 15 hr, and with secondary antibody at room temperature for 1 hr. Coverslips were mounted using ProLong diamond antifade mountant with DAPI (P36962, ThermoFisher Scientific, UK). Primary antibodies: anti-Tuj1 (RRID:AB_444319)(ab18207, Abcam, UK), rabbit anti-GFP (RRID:AB_10073917)(A11122, ThermoFisher Scientific, UK) and rabbit anti-Aurora-B kinase (RRID:AB_302923)(ab2254 Abcam, UK). Secondary antibody: goat anti-rabbit Alexa Fluor 488 (RRID:AB_143165)(A11008, ThermoFisher Scientific, UK).

smFISH imaging and analysis

Request a detailed protocol

Z-stack images were acquired at 0.2 μm z-increments by Delta Vision widefield microscopy (Olympus IX71 microscope body, Coolsnap HQ CCD camera (Photometrics, US), API XYZ stage), using a 60X (1.42NA) objective. Image stacks were viewed using Imaris (Version 7.2.2, Bitplane), and cell areas segmented manually using background auto-fluorescence in the Quasar 570 channel on a maximum projection of the image. Hes1 mRNA spots per cell were counted by a semi-automated method, using the inbuilt 'analyze spots' function in Imaris. The following settings were used: estimated XY diameter = 0.3 μM, estimated Z diameter = 0.6 μM, background subtraction = yes, filter type = maximum intensity. The threshold cutoff between signal and background was placed manually at the ‘trough’ in the spots profile for the maximum intensity filter. Active transcription sites were identified manually, defined as Hes1 intronic smFISH spots within the nucleus that colocalised with exonic Hes1 smFISH spots.

Flow cytometry

Request a detailed protocol

Tau-GFP NS cells were harvested at 96 hr, 120 hr, 144 hr and 168 hr of neuronal differentiation following transfection with 40 nM of miR-9 mimic or control mimic. Cells were dissociated using accutase (Sigma, UK) and washed in DMEM-F12 supplemented with 0.2% BSA. Flow cytometry analysis for GFP detection was performed on LSR Fortessa (BD Biosciences, UK) and data were analysed with the FlowJo software (RRID:SCR_008520)(TreeStar, US).

Fluorescence Correlation Spectroscopy (FCS) protein quantitation

Request a detailed protocol

FCS was used to count the number of fluorescent protein molecules expressed in live NS-E cells that were transduced with lentiviral VENUS:HES1 constructs. FCS involves laser excitation and data collection in a tightly focused microscopic region known as the confocal volume. Fluorescent particles diffusing in the confocal volume give rise to fluorescent intensity fluctuations $F (t)$ (example shown in Figure 2—figure supplement 1A) which are analysed over lag time $τ > 0$ to produce the normalised autocorrelation function (Figure 2—figure supplement 1B):

R (τ) = \frac{⟨ F (t + τ) F (t) ⟩}{{⟨ F (t) ⟩}^{2}} - 1

The amplitude of the autocorrelation function is proportional to 1/c where c represents the number of molecules detected in the confocal volume, i.e. molecule concentration. The shape of $R (τ)$ can be explained by models taking into consideration multiple species of fluorescent particles characterised by different diffusion constants. Preliminary model selection indicated that two components were required to explain the variability in our data (results not shown) and a non-fluorescent component was added to account for transitions of fluorescent particles to a triplet state. Normalised autocorrelations were fitted with a two component model with triplet state previously used for live cells (Dross et al., 2009; Bagnall et al., 2015).

FCS experiments were carried out using a Zeiss LSM880 microscope with a Plan-Apochromat 40x, 1.4 NA oil-immersion objective on cells cultured and maintained at 37^oC and $5 %$ CO2. Venus (EYFP) fluorescence was excited with 514 nm laser light and emission collected between 516 and 570 nm. Data from individual cell nuclei was collected using 5 x 5 s runs at 0.2 to 0.3% laser power for each measurement using Zen 2.1 (RRID:SCR_013672)(Zeiss) software. We used custom Matlab (RRID:SCR_001622) routines to fit the autocorrelation data (Figure 2—figure supplement 1B) and obtain protein concentration in the confocal volume. Measurements collected from cells with count per molecule lower than 500 kHz were excluded from the final results. Single-cell data of number of molecules in the cell nucleus was obtained by adjusting concentration to the average volumetric ratio between nuclear volume and confocal volume. NS-E cells have a nuclear volume of mean 721 fL ± 105 fL (S.D.) estimated using nuclear staining and 3D reconstruction from z-stack images in Imaris. The confocal volume had been previously determined with mean 0.57fL ± 11 fL (S.D.) (Bagnall et al., 2015).

Western blotting

Request a detailed protocol

Nuclear extracts of UbC-VENUS:HES1 NS-E cells treated with DMSO or MG132 (MerckMillipore, UK) for 3 hr were generated as described in Schreiber et al. (1989). Western blots were performed using 10% NUPAGEBis-Tris gels (ThermoFisher, UK) and Whatman Protran nitrocellulose membrane (Sigma, UK) according to manufacturer’s instructions. Antibodies used were anti-HES1 (RRID:AB_590682)(D134-3, MBL, Japan), anti-LAMIN B1 (RRID:AB_443298)(ab16048, Abcam, UK), HRP-linked anti-rat IgG (RRID:AB_10694715)(7077, CST, UK) and HRP-linked anti-rabbit IgG (RRID:AB_2099233)(7074, CST, UK).

Statistical analysis

Request a detailed protocol

Data were analysed using pairwise hypothesis testing methods and multiple hypothesis tests as indicated in the figure legends using Prism 7 (RRID:SCR_002798 )(GraphPad, U.S.A). Statistical significance is reported for p-values <0.05 (*), <0.01(**) and <0.001 (***). Errors are reported as standard error of the mean (S.E.M) calculated from pooled data taking into account standard deviation (S.D) and sample size (n) as S.E.M=S.D./sqrt(n).

Appendix

Model parameters

Appendix 1—table 1

Parameters of the model and their values.

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

Parameter	Set 1	Set 2	Interpretation	Reference
$τ$	29 min	29 min	Time delay in HES1 protein production	(Lewis, 2003)
$α_{m}$	1 min⁻¹	1 min⁻¹	Hes1 transcription rate in absence of Hes1 protein
$α_{p}$	1 min⁻¹	2 min⁻¹	Translation rate
$α_{r}$	1 min⁻¹	1 min⁻¹	miR-9 transcription rate in absence of Hes1 protein
$p_{0}$	390	390	Amount of protein required to reduce Hes1 transcription by half	Fitted in (Goodfellow et al., 2014)
$n_{0}$	5	5	Quantifies the step-like nature of $G$	(Monk, 2003)
$r_{0}$	100	80	Amount of miR-9 required to reduce Hes1 degradation by half
$m_{0}$	5	5	Quantifies the step-like nature of $S$
$b_{l}$	$l n (2) / 20$ min⁻¹	$l n (2) / 20$ min⁻¹	Lower bound for Hes1 mRNA half-life	(Bonev et al., 2012)
$b_{u}$	$l n (2) / 35$ min⁻¹	$l n (2) / 35$ min⁻¹	Upper bound for Hes1 mRNA half-life	(Bonev et al., 2012)
$r_{1}$	300	100	Amount of miR-9 required to reduce HES1 protein translation by half
$m_{1}$	5	5	Quantifies the step-like nature of $F$
$n_{1}$	5	5	Quantifies the step-like nature of $G_{r}$
$μ_{p}$	22 min	22 min	Half-life of HES1 protein	(Hirata et al., 2002)
$p_{1}$	260	280	Amount of protein required to reduce miR-9 production rate by half
$μ_{r}$	1000 min	1000 min	Half-life of miR-9

Deterministic dynamics

Starting with the microscopic dynamics of the Hes1/miR-9 described in the Materials and Methods, we apply the system size expansion as described in Brett and Galla (Brett and Galla, 2013) and Galla (Galla, 2009). To lowest order in $Ω^{- 1 / 2}$ , one recovers the deterministic dynamics described in Goodfellow et al. (Goodfellow et al., 2014)

\begin{aligned} (23) & \frac{d m}{d t} & = G (p (t - τ)) - S (r) m, \\ (24) & \frac{d p}{d t} & = F (r) m - \frac{l n (2)}{μ_{p}} p, \\ (25) & \frac{d r}{d t} & = G_{r} (p) - \frac{l n (2)}{μ_{r}} r . \end{aligned}

Linear noise approximation

The next-to-leading order in $Ω^{- 1 / 2}$ represents the LNA and recovers the probability distribution describing random deviations in molecule numbers around the deterministic trajectory. The fluctuations in mRNA, protein and miR-9 molecule number are denoted $ϵ_{m}$ , $ϵ_{p}$ and $ϵ_{r}$ , respectively. In order to account for noise induced dynamics below the Hopf bifurcation, we assume that the deterministic dynamics reach a fixed-point. When the deterministic dynamics have reached a fixed-point (i.e. after large time), the dynamics of the fluctuations can be approximated by the linear Langevin equations

\begin{aligned} (26) & \frac{d ϵ_{m}}{d t} & = G^{'} (p^{*}) ϵ_{p} (t - τ) - S (r^{*}) ϵ_{m} (t) - S^{'} (r^{*}) m^{*} ϵ_{r} (t) + ζ_{m} (t), \\ (27) & \frac{d ϵ_{p}}{d t} & = F (r^{*}) ϵ_{m} (t) + F^{'} (r^{*}) m^{*} ϵ_{r} (t) - \frac{l n (2)}{μ_{p}} ϵ_{p} (t) + ζ_{p} (t), \\ (28) & \frac{d ϵ_{r}}{d t} & = G_{r}^{'} (p^{*}) ϵ_{p} (t) - \frac{\ln (2)}{μ_{r}} ϵ_{r} (t) + ζ_{r} (t), \end{aligned}

where $G^{'} (p^{*})$ , $S^{'} (r^{*})$ , $F^{'} (r^{*})$ and $G_{r}^{'} (p^{*})$ represent are linearisations of the Hill functions $G, S$ and $G_{r}$ around the deterministic fixed point, and where $ζ (t)$ represents white noise. Within the linear-noise approximation this noise is assumed to be additive. At the deterministic fixed point, the noise correlators are given by

\begin{array}{ll} ⟨ ζ_{m} (t) ζ_{m} (t^{'}) ⟩ & = (G (p^{*}) + S (r^{*}) m^{*}) δ (t - t^{'}), \\ ⟨ ζ_{p} (t) ζ_{p} (t^{'}) ⟩ & = (F (r^{*}) m^{*} + \frac{l n (2)}{μ_{p}} p^{*}) δ (t - t^{'}), \\ ⟨ ζ_{r} (t) ζ_{r} (t^{'}) ⟩ & = (G_{r} (p^{*}) + \frac{l n (2)}{μ_{r}} r^{*}) δ (t - t^{'}), \\ ⟨ ζ_{m} (t) ζ_{p} (t^{'}) ⟩ & = ⟨ ζ_{p} (t) (ζ_{m} (t^{'}) ⟩ = 0, \\ ⟨ ζ_{m} (t) ζ_{r} (t^{'}) ⟩ & = ⟨ ζ_{r} (t) ζ_{m} (t^{'}) ⟩ = 0, \\ ⟨ ζ_{p} (t) ζ_{r} (t^{'}) ⟩ & = ⟨ ζ_{r} (t) ζ_{p} (t^{'}) ⟩ = 0 \end{array}

The LNA therefore describes a process with linear delay relaxation dynamics supplemented by additive noise. This constitutes a multivariate Ornstein-Uhlenbeck process (with delay). To analyse this further, we carry out a Fourier transform (with respect to time), and obtain

\begin{aligned} (30) & i ω {\tilde{ϵ}}_{m} (ω) & = G^{'} (p^{*}) {\tilde{ϵ}}_{p} (ω) e^{- i ω t} - S (r^{*}) {\tilde{ϵ}}_{m} (ω) - S^{'} (r^{*}) m^{*} {\tilde{ϵ}}_{r} (ω) + {\tilde{ζ}}_{m} (ω), \\ (31) & i ω {\tilde{ϵ}}_{p} (ω) & = F (r^{*}) {\tilde{ϵ}}_{m} (ω) + F^{'} (r^{*}) m^{*} {\tilde{ϵ}}_{r} (ω) - \frac{l n (2)}{μ_{p}} {\tilde{ϵ}}_{p} (ω) + {\tilde{ζ}}_{p} (ω), \\ (32) & i ω {\tilde{ϵ}}_{r} (ω) & = G_{r}^{'} (p^{*}) {\tilde{ϵ}}_{p} (ω) - \frac{l n (2)}{μ_{r}} {\tilde{ϵ}}_{r} (ω) + {\tilde{ζ}}_{r} (ω) . \end{aligned}

In matrix form this can be represented as

\underline{\underline{M}} (ω) \tilde{ϵ} (ω) = \tilde{ζ} (ω),

where

\underline{\underline{M}} (ω) = (\begin{array}{ccc} i ω + S (r^{*}) & - G^{'} (p^{*}) e^{- i ω t} & S^{'} (r^{*}) m^{*} \\ - F (r^{*}) & i ω + \frac{l n (2)}{μ_{p}} & - F^{'} (r^{*}) m^{*} \\ 0 & - G_{r}^{'} (p^{*}) & i ω + \frac{l n (2)}{μ_{r}} \end{array}) .

It is now useful to introduce the matrix of spectra, $\underline{\underline{S}} (ω)$ , with $S_{i j} (ω) = ⟨ \tilde{ϵ_{i}} (ω) \tilde{ϵ_{j}} (ω)^{†} ⟩$ , where $i, j$ represent the variables $m, p$ and $r$ above. One then has

\underline{\underline{S}} (ω) = \underline{\underline{M}} (ω)^{- 1} ⟨ \tilde{ζ} (ω) \tilde{ζ} (ω)^{†} ⟩ \underline{\underline{M}} (ω)^{† - 1},

where $⟨ \tilde{ζ} (ω) \tilde{ζ} (ω)^{†} ⟩$ is

⟨ \tilde{ζ} (ω) \tilde{ζ} (ω^{'})^{†} ⟩ = (\begin{array}{ccc} (G (p^{*}) + S (r^{*}) m^{*}) & 0 & 0 \\ 0 & (F (r^{*}) m^{*} + \frac{\ln (2)}{μ_{p}} p^{*}) & 0 \\ 0 & 0 & (G_{r} (p^{*}) + \frac{\ln (2)}{μ_{r}} r^{*}) \end{array}) δ (ω + ω^{'}) .

The further analysis is then based on a numerical evaluation of Equation 35.

References

(2016) Concentration sensing by the moving nucleus in cell fate determination: A computational analysis
PLoS One 11:e0149213.

https://doi.org/10.1371/journal.pone.0149213
- Google Scholar
1. Alon U
(2007)
Chapman and Hall/CRC

An introduction to systems biology - design principles of biological circuits, Chapman and Hall/CRC, London, Taylor and Francis group..
- Google Scholar
(2007) Stochastic amplification in epidemics
Journal of the Royal Society Interface 4:575–582.

https://doi.org/10.1098/rsif.2006.0192
- Google Scholar
1. Anderson DF
(2007) A modified next reaction method for simulating chemical systems with time dependent propensities and delays
The Journal of Chemical Physics 127:214107.

https://doi.org/10.1063/1.2799998
- Google Scholar
1. Arias AM
2. Hayward P
(2006) Filtering transcriptional noise during development: concepts and mechanisms
Nature Reviews. Genetics 7:34–44.

https://doi.org/10.1038/nrg1750
- Google Scholar
1. Bagnall J
2. Boddington C
3. Boyd J
4. Brignall R
5. Rowe W
6. Jones NA
7. Schmidt L
8. Spiller DG
9. White MRH
10. Paszek P
(2015) Quantitative dynamic imaging of immune cell signalling using lentiviral gene transfer
Integr. Biol. 7:713–725.

https://doi.org/10.1039/C5IB00067J
- Google Scholar
1. Bahar Halpern K
2. Caspi I
3. Lemze D
4. Levy M
5. Landen S
6. Elinav E
7. Ulitsky I
8. Itzkovitz S
(2015) Nuclear Retention of mRNA in Mammalian Tissues
Cell Reports 13:2653–2662.

https://doi.org/10.1016/j.celrep.2015.11.036
- Google Scholar
1. Bahar Halpern K
2. Tanami S
3. Landen S
4. Chapal M
5. Szlak L
6. Hutzler A
7. Nizhberg A
8. Itzkovitz S
(2015) Bursty gene expression in the intact mammalian liver
Molecular Cell 58:147–156.

https://doi.org/10.1016/j.molcel.2015.01.027
- Google Scholar
(2011) Cellular decision making and biological noise: from microbes to mammals
Cell 144:910–925.

https://doi.org/10.1016/j.cell.2011.01.030
- Google Scholar
1. Barrio M
2. Burrage K
3. Leier A
4. Tian T
(2006) Oscillatory regulation of Hes1: Discrete stochastic delay modelling and simulation
PLoS Computational Biology 2:e117.

https://doi.org/10.1371/journal.pcbi.0020117
- Google Scholar
(2015) Control of Transcript Variability in Single Mammalian Cells
Cell 163:1596–1610.

https://doi.org/10.1016/j.cell.2015.11.018
- Google Scholar
(2002) Proneural genes and the specification of neural cell types
Nature Reviews. Neuroscience 3:517–530.

https://doi.org/10.1038/nrn874
- Google Scholar
1. Bibel M
2. Richter J
3. Schrenk K
4. Tucker KL
5. Staiger V
6. Korte M
7. Goetz M
8. Barde YA
(2004) Differentiation of mouse embryonic stem cells into a defined neuronal lineage
Nature Neuroscience 7:1003–1009.

https://doi.org/10.1038/nn1301
- Google Scholar
(2014) Reconciling competence and transcriptional hierarchies with stochasticity in retinal lineages
Current Opinion in Neurobiology 27:68–74.

https://doi.org/10.1016/j.conb.2014.02.014
- Google Scholar
(2011) MicroRNA-9 reveals regional diversity of neural progenitors along the anterior-posterior axis
Developmental Cell 20:19–32.

https://doi.org/10.1016/j.devcel.2010.11.018
- Google Scholar
(2012) MicroRNA-9 Modulates Hes1 ultradian oscillations by forming a double-negative feedback loop
Cell Reports 2:10–18.

https://doi.org/10.1016/j.celrep.2012.05.017
- Google Scholar
1. Brett T
2. Galla T
(2013) Stochastic processes with distributed delays: Chemical langevin equation and linear-noise approximation
Physical Review Letters 110:250601.

https://doi.org/10.1103/PhysRevLett.110.250601
- Google Scholar
1. Brett T
2. Galla T
(2014) Gaussian approximations for stochastic systems with delay: chemical Langevin equation and application to a Brusselator system
The Journal of Chemical Physics 140:124112.

https://doi.org/10.1063/1.4867786
- Google Scholar
1. Cai X
(2007) Exact stochastic simulation of coupled chemical reactions with delays
The Journal of Chemical Physics 126:124108.

https://doi.org/10.1063/1.2710253
- Google Scholar
1. Chen C
2. Ridzon DA
3. Broomer AJ
4. Zhou Z
5. Lee DH
6. Nguyen JT
7. Barbisin M
8. Xu NL
9. Mahuvakar VR
10. Andersen MR
11. Lao KQ
12. Livak KJ
13. Guegler KJ
(2005) Real-time quantification of microRNAs by stem-loop RT-PCR
Nucleic Acids Research 33:e179.

https://doi.org/10.1093/nar/gni178
- Google Scholar
1. Conti L
2. Pollard SM
3. Gorba T
4. Reitano E
5. Toselli M
6. Biella G
7. Sun Y
8. Sanzone S
9. Ying QL
10. Cattaneo E
11. Smith A
(2005) Niche-independent symmetrical self-renewal of a mammalian tissue stem cell
PLoS Biology 3:e283.

https://doi.org/10.1371/journal.pbio.0030283
- Google Scholar
(2012) miR-9 controls the timing of neurogenesis through the direct inhibition of antagonistic factors
Developmental Cell 22:1052–1064.

https://doi.org/10.1016/j.devcel.2012.03.003
- Google Scholar
1. Corrigan AM
2. Chubb JR
(2014) Regulation of transcriptional bursting by a naturally oscillating signal
Current Biology 24:205–211.

https://doi.org/10.1016/j.cub.2013.12.011
- Google Scholar
1. Dross N
2. Spriet C
3. Zwerger M
4. Müller G
5. Waldeck W
6. Langowski J
(2009) Mapping eGFP oligomer mobility in living cell nuclei
PloS One 4:e5041.

https://doi.org/10.1371/journal.pone.0005041
- Google Scholar
1. Elf J
2. Ehrenberg M
(2003) Fast evaluation of fluctuations in biochemical networks with the linear noise approximation
Genome Research 13:2475–2484.

https://doi.org/10.1101/gr.1196503
- Google Scholar
(2002) Numerical bifurcation analysis of delay differential equations using DDE-BIFTOOL
ACM Transactions on Mathematical Software 28:1–21.

https://doi.org/10.1145/513001.513002
- Google Scholar
(2012) Fate-restricted neural progenitors in the mammalian cerebral cortex
Science 337:746–749.

https://doi.org/10.1126/science.1223616
- Google Scholar
1. Franco SJ
2. Müller U
(2013) Shaping our minds: stem and progenitor cell diversity in the mammalian neocortex
Neuron 77:19–34.

https://doi.org/10.1016/j.neuron.2012.12.022
- Google Scholar
1. Furusawa C
2. Kaneko K
(2012) A dynamical-systems view of stem cell biology
Science 338:215–217.

https://doi.org/10.1126/science.1224311
- Google Scholar
1. Galla T
(2009) Intrinsic fluctuations in stochastic delay systems: Theoretical description and application to a simple model of gene regulation
Physical Review E 80:021909.

https://doi.org/10.1103/PhysRevE.80.021909
- Google Scholar
1. Gao P
2. Postiglione MP
3. Krieger TG
4. Hernandez L
5. Wang C
6. Han Z
7. Streicher C
8. Papusheva E
9. Insolera R
10. Chugh K
11. Kodish O
12. Huang K
13. Simons BD
14. Luo L
15. Hippenmeyer S
16. Shi SH
(2014) Deterministic progenitor behavior and unitary production of neurons in the neocortex
Cell 159:775–788.

https://doi.org/10.1016/j.cell.2014.10.027
- Google Scholar
1. Garcia-Ojalvo J
2. Martinez Arias A
(2012) Towards a statistical mechanics of cell fate decisions
Current Opinion in Genetics & Development 22:619–626.

https://doi.org/10.1016/j.gde.2012.10.004
- Google Scholar
(2008) An intrinsic mechanism of corticogenesis from embryonic stem cells
Nature 455:351–357.

https://doi.org/10.1038/nature07287
- Google Scholar
1. Gillespie DT
(1977) Exact stochastic simulation of coupled chemical reactions
The Journal of Physical Chemistry 81:2340–2361.

https://doi.org/10.1021/j100540a008
- Google Scholar
1. Gillespie DT
(2000) The chemical Langevin equation
The Journal of Chemical Physics 113:297.

https://doi.org/10.1063/1.481811
- Google Scholar
(2014)
International Review of Cell and Molecular Biology

273–321, Making of a Retinal Cell, International Review of Cell and Molecular Biology, Elsevier Inc, 10.1016/b978-0-12-800097-7.00007-5.
- Google Scholar
(2005) Real-time kinetics of gene activity in individual bacteria
Cell 123:1025–1036.

https://doi.org/10.1016/j.cell.2005.09.031
- Google Scholar
1. Gomes FL
2. Zhang G
3. Carbonell F
4. Correa JA
5. Harris WA
6. Simons BD
7. Cayouette M
(2011) Reconstruction of rat retinal progenitor cell lineages in vitro reveals a surprising degree of stochasticity in cell fate decisions
Development 138:227–235.

https://doi.org/10.1242/dev.059683
- Google Scholar
(2014) microRNA input into a neural ultradian oscillator controls emergence and timing of alternative cell states
Nature Communications 5:3399.

https://doi.org/10.1038/ncomms4399
- Google Scholar
(2013) Molecular logic of neocortical projection neuron specification, development and diversity
Nature Reviews. Neuroscience 14:755–769.

https://doi.org/10.1038/nrn3586
- Google Scholar
(2011) How accurate are the nonlinear chemical Fokker-Planck and chemical Langevin equations?
The Journal of Chemical Physics 135:084103.

https://doi.org/10.1063/1.3625958
- Google Scholar
1. Grima R
(2012) A study of the accuracy of moment-closure approximations for stochastic chemical kinetics
The Journal of Chemical Physics 136:154105.

https://doi.org/10.1063/1.3702848
- Google Scholar
1. Götz M
2. Huttner WB
(2005) The cell biology of neurogenesis
Nature Reviews. Molecular Cell Biology 6:777–788.

https://doi.org/10.1038/nrm1739
- Google Scholar
1. Harper CV
2. Finkenstädt B
3. Woodcock DJ
4. Friedrichsen S
5. Semprini S
6. Ashall L
7. Spiller DG
8. Mullins JJ
9. Rand DA
10. Davis JR
11. White MR
(2011) Dynamic analysis of stochastic transcription cycles
PLoS Biology 9:e1000607.

https://doi.org/10.1371/journal.pbio.1000607
- Google Scholar
1. Hatakeyama J
2. Bessho Y
3. Katoh K
4. Ookawara S
5. Fujioka M
6. Guillemot F
7. Kageyama R
(2004) Hes genes regulate size, shape and histogenesis of the nervous system by control of the timing of neural stem cell differentiation
Development 131:5539–5550.

https://doi.org/10.1242/dev.01436
- Google Scholar
1. He J
2. Zhang G
3. Almeida AD
4. Cayouette M
5. Simons BD
6. Harris WA
(2012) How variable clones build an invariant retina
Neuron 75:786–798.

https://doi.org/10.1016/j.neuron.2012.06.033
- Google Scholar
1. Hirata H
2. Yoshiura S
3. Ohtsuka T
4. Bessho Y
5. Harada T
6. Yoshikawa K
7. Kageyama R
(2002) Oscillatory expression of the bHLH factor Hes1 regulated by a negative feedback loop
Science 298:840–843.

https://doi.org/10.1126/science.1074560
- Google Scholar
1. Huang S
(2011) Systems biology of stem cells: three useful perspectives to help overcome the paradigm of linear pathways
Philosophical Transactions of the Royal Society B: Biological Sciences 366:2247–2259.

https://doi.org/10.1098/rstb.2011.0008
- Google Scholar
1. Huh D
2. Paulsson J
(2011) Non-genetic heterogeneity from stochastic partitioning at cell division
Nature Genetics 43:95–100.

https://doi.org/10.1038/ng.729
- Google Scholar
1. Huh D
2. Paulsson J
(2011) Random partitioning of molecules at cell division
PNAS 108:15004–15009.

https://doi.org/10.1073/pnas.1013171108
- Google Scholar
1. Imayoshi I
2. Isomura A
3. Harima Y
4. Kawaguchi K
5. Kori H
6. Miyachi H
7. Fujiwara T
8. Ishidate F
9. Kageyama R
(2013) Oscillatory control of factors determining multipotency and fate in mouse neural progenitors
Science 342:1203–1208.

https://doi.org/10.1126/science.1242366
- Google Scholar
(1995) Targeted disruption of mammalian hairy and Enhancer of split homolog-1 (HES-1) leads to up-regulation of neural helix-loop-helix factors, premature neurogenesis, and severe neural tube defects
Genes & Development 9:3136–3148.

https://doi.org/10.1101/gad.9.24.3136
- Google Scholar
1. Jenkins RP
2. Hanisch A
3. Soza-Ried C
4. Sahai E
5. Lewis J
(2015) Stochastic Regulation of her1/7 Gene Expression Is the Source of Noise in the Zebrafish Somite Clock Counteracted by Notch Signalling
PLoS Computational Biology 11:e1004459.

https://doi.org/10.1371/journal.pcbi.1004459
- Google Scholar
(2003) Sustained oscillations and time delays in gene expression of protein Hes1
FEBS Letters 541:176–177.

https://doi.org/10.1016/S0014-5793(03)00279-5
- Google Scholar
1. Kaern M
2. Elston TC
3. Blake WJ
4. Collins JJ
(2005) Stochasticity in gene expression: from theories to phenotypes
Nature Reviews. Genetics 6:451–464.

https://doi.org/10.1038/nrg1615
- Google Scholar
(2007) The Hes gene family: repressors and oscillators that orchestrate embryogenesis
Development 134:1243–1251.

https://doi.org/10.1242/dev.000786
- Google Scholar
(2008) Dynamic Notch signaling in neural progenitor cells and a revised view of lateral inhibition
Nature Neuroscience 11:1247–1251.

https://doi.org/10.1038/nn.2208
- Google Scholar
1. Klein AM
2. Simons BD
(2011) Universal patterns of stem cell fate in cycling adult tissues
Development 138:3103–3111.

https://doi.org/10.1242/dev.060103
- Google Scholar
1. Kohwi M
2. Doe CQ
(2013) Temporal fate specification and neural progenitor competence during development
Nature Reviews Neuroscience 14:823–838.

https://doi.org/10.1038/nrn3618
- Google Scholar
1. Lau P
2. Verrier JD
3. Nielsen JA
4. Johnson KR
5. Notterpek L
6. Hudson LD
(2008) Identification of dynamically regulated microRNA and mRNA networks in developing oligodendrocytes
The Journal of Neuroscience : The Official Journal of the Society for Neuroscience 28:11720–11730.

https://doi.org/10.1523/JNEUROSCI.1932-08.2008
- Google Scholar
(2010) Fundamental limits on the suppression of molecular fluctuations
Nature 467:174–178.

https://doi.org/10.1038/nature09333
- Google Scholar
1. Lewis J
(2003) Autoinhibition with transcriptional delay: a simple mechanism for the zebrafish somitogenesis oscillator
Current Biology 13:1398–1408.

https://doi.org/10.1016/s0960-9822(03)00534-7
- Google Scholar
(2013) Optimisation of culture conditions for differentiation of C17.2 neural stem cells to be used for in vitro toxicity tests
Toxicology in Vitro 27:1565–1569.

https://doi.org/10.1016/j.tiv.2012.04.020
- Google Scholar
(2009) ChIP-on-chip significance analysis reveals large-scale binding and regulation by human transcription factor oncogenes
PNAS 106:244–249.

https://doi.org/10.1073/pnas.0806445106
- Google Scholar
1. McKane AJ
2. Newman TJ
(2005) Predator-prey cycles from resonant amplification of demographic stochasticity
Physical Review Letters 94:218102.

https://doi.org/10.1103/PhysRevLett.94.218102
- Google Scholar
1. Momiji H
2. Monk NA
(2008) Dissecting the dynamics of the Hes1 genetic oscillator
Journal of Theoretical Biology 254:784–798.

https://doi.org/10.1016/j.jtbi.2008.07.013
- Google Scholar
1. Monk NA
(2003) Oscillatory expression of Hes1, p53, and NF-kappaB driven by transcriptional time delays
Current Biology 13:1409–1413.

https://doi.org/10.1016/s0960-9822(03)00494-9
- Google Scholar
(2012) Using gene expression noise to understand gene regulation
Science 336:183–187.

https://doi.org/10.1126/science.1216379
- Google Scholar
1. Nakamura Y
2. Sakakibara S
3. Miyata T
4. Ogawa M
5. Shimazaki T
6. Weiss S
7. Kageyama R
8. Okano H
(2000)
The bHLH gene hes1 as a repressor of the neuronal commitment of CNS stem cells

Journal of Neuroscience 20:283–293.
- Google Scholar
(2013) Systematic identification of signal-activated stochastic gene regulation
Science 339:584–587.

https://doi.org/10.1126/science.1231456
- Google Scholar
1. Osella M
2. Bosia C
3. Corá D
4. Caselle M
(2011) The role of incoherent microRNA-mediated feedforward loops in noise buffering
PLoS Computational Biology 7:e1001101.

https://doi.org/10.1371/journal.pcbi.1001101
- Google Scholar
1. Pearson BJ
2. Doe CQ
(2004) Specification of temporal identity in the developing nervous system
Annual Review of Cell and Developmental Biology 20:619–647.

https://doi.org/10.1146/annurev.cellbio.19.111301.115142
- Google Scholar
1. Pfeuty B
(2015a) A computational model for the coordination of neural progenitor self-renewal and differentiation through Hes1 dynamics
Development 142:477–485.

https://doi.org/10.1242/dev.112649
- Google Scholar
1. Pfeuty B
(2015b) Neuronal specification exploits the inherent flexibility of cell-cycle gap phases
Neurogenesis 2:e1095694.

https://doi.org/10.1080/23262133.2015.1095694
- Google Scholar
1. Pollard SM
(2013) In vitro expansion of fetal neural progenitors as adherent cell lines
Methods in Molecular Biology 1059:13–24.

https://doi.org/10.1007/978-1-62703-574-3_2
- Google Scholar
1. Qian X
2. Shen Q
3. Goderie SK
4. He W
5. Capela A
6. Davis AA
7. Temple S
(2000) Timing of CNS cell generation: a programmed sequence of neuron and glial cell production from isolated murine cortical stem cells
Neuron 28:69–80.

https://doi.org/10.1016/s0896-6273(00)00086-6
- Google Scholar
1. Rao CV
2. Wolf DM
3. Arkin AP
(2002) Control, exploitation and tolerance of intracellular noise
Nature 420:231–237.

https://doi.org/10.1038/nature01258
- Google Scholar
1. Rapaport DH
2. Wong LL
3. Wood ED
4. Yasumura D
5. LaVail MM
(2004) Timing and topography of cell genesis in the rat retina
The Journal of Comparative Neurology 474:304–324.

https://doi.org/10.1002/cne.20134
- Google Scholar
1. Rué P
2. Martinez Arias A
(2015) Cell dynamics and gene expression control in tissue homeostasis and development
Molecular Systems Biology 11:792.

https://doi.org/10.15252/msb.20145549
- Google Scholar
(1992) Two mammalian helix-loop-helix factors structurally related to Drosophila hairy and Enhancer of split
Genes & Development 6:2620–2634.

https://doi.org/10.1101/gad.6.12b.2620
- Google Scholar
(1989) Rapid detection of octamer binding proteins with 'mini-extracts', prepared from a small number of cells
Nucleic Acids Research 17:6419.

https://doi.org/10.1093/nar/17.15.6419
- Google Scholar
1. Schwanhäusser B
2. Busse D
3. Li N
4. Dittmar G
5. Schuchhardt J
6. Wolf J
7. Chen W
8. Selbach M
(2011) Global quantification of mammalian gene expression control
Nature 473:337–342.

https://doi.org/10.1038/nature10098
- Google Scholar
1. Shen Q
2. Wang Y
3. Dimos JT
4. Fasano CA
5. Phoenix TN
6. Lemischka IR
7. Ivanova NB
8. Stifani S
9. Morrisey EE
10. Temple S
(2006) The timing of cortical neurogenesis is encoded within lineages of individual progenitor cells
Nature Neuroscience 9:743–751.

https://doi.org/10.1038/nn1694
- Google Scholar
1. Shimojo H
2. Isomura A
3. Ohtsuka T
4. Kori H
5. Miyachi H
6. Kageyama R
(2016) Oscillatory control of Delta-like1 in cell interactions regulates dynamic gene expression and tissue morphogenesis
Genes & Development 30:102–116.

https://doi.org/10.1101/gad.270785.115
- Google Scholar
(2008) Oscillations in notch signaling regulate maintenance of neural progenitors
Neuron 58:52–64.

https://doi.org/10.1016/j.neuron.2008.02.014
- Google Scholar
1. Simons BD
2. Clevers H
(2011) Strategies for homeostatic stem cell self-renewal in adult tissues
Cell 145:851–862.

https://doi.org/10.1016/j.cell.2011.05.033
- Google Scholar
1. Slater JL
2. Landman KA
3. Hughes BD
4. Shen Q
5. Temple S
(2009) Cell lineage tree models of neurogenesis
Journal of Theoretical Biology 256:164–179.

https://doi.org/10.1016/j.jtbi.2008.09.034
- Google Scholar
(2013) Spatial stochastic modelling of the Hes1 gene regulatory network: intrinsic noise can explain heterogeneity in embryonic stem cell differentiation
Journal of the Royal Society Interface 10:20120988.

https://doi.org/10.1098/rsif.2012.0988
- Google Scholar
1. Suter DM
2. Molina N
3. Naef F
4. Schibler U
(2011) Origins and consequences of transcriptional discontinuity
Current Opinion in Cell Biology 23:657–662.

https://doi.org/10.1016/j.ceb.2011.09.004
- Google Scholar
(2002) Intrinsic and extrinsic contributions to stochasticity in gene expression
PNAS 99:12795–12800.

https://doi.org/10.1073/pnas.162041399
- Google Scholar
(2012) MicroRNA9 regulates neural stem cell differentiation by controlling Hes1 expression dynamics in the developing brain
Genes to Cells : Devoted to Molecular & Cellular Mechanisms 17:952–961.

https://doi.org/10.1111/gtc.12009
- Google Scholar
1. Temple S
(2001) The development of neural stem cells
Nature 414:112–117.

https://doi.org/10.1038/35102174
- Google Scholar
(2012) Intrinsic noise analyzer: a software package for the exploration of stochastic biochemical kinetics using the system size expansion
PLoS One 7:e38518.

https://doi.org/10.1371/journal.pone.0038518
- Google Scholar
1. Zechner C
2. Unger M
3. Pelet S
4. Peter M
5. Koeppl H
(2014) Scalable inference of heterogeneous reaction kinetics from pooled single-cell recordings
Nature Methods 11:197–202.

https://doi.org/10.1038/nmeth.2794
- Google Scholar

Article and author information

Author details

Nick E Phillips

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
NEP, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Competing interests
The authors declare that no competing interests exist.
Cerys S Manning

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
CSM, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Contributed equally with
Tom Pettini, Veronica Biga, Elli Marinopoulou and Peter Stanley

Competing interests
The authors declare that no competing interests exist.
Tom Pettini

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
TP, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Contributed equally with
Cerys S Manning, Veronica Biga, Elli Marinopoulou and Peter Stanley

Competing interests
The authors declare that no competing interests exist.
Veronica Biga

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
VB, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Contributed equally with
Cerys S Manning, Tom Pettini, Elli Marinopoulou and Peter Stanley

Competing interests
The authors declare that no competing interests exist.
Elli Marinopoulou

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
EM, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Contributed equally with
Cerys S Manning, Tom Pettini, Veronica Biga and Peter Stanley

Competing interests
The authors declare that no competing interests exist.
Peter Stanley

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
PS, Conception and design, Acquisition of data, Analysis and interpretation of data

Contributed equally with
Cerys S Manning, Tom Pettini, Veronica Biga and Elli Marinopoulou

Competing interests
The authors declare that no competing interests exist.
James Boyd

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
JBo, Acquisition of data, Analysis and interpretation of data

Competing interests
The authors declare that no competing interests exist.
James Bagnall

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
JBa, Analysis and interpretation of data, Contributed unpublished essential data or reagents

Competing interests
The authors declare that no competing interests exist.
Pawel Paszek

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
PP, Analysis and interpretation of data, Contributed unpublished essential data or reagents

Competing interests
The authors declare that no competing interests exist.
David G Spiller

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
DGS, provided FCS training, designed experiments and had supervisory input into the FCS work, Analysis and interpretation of data

Competing interests
The authors declare that no competing interests exist.
Michael RH White

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
MRHW, supervisory input into the FCS work, Analysis and interpretation of data

Competing interests
The authors declare that no competing interests exist.

"This ORCID iD identifies the author of this article:" 0000-0002-3617-3232
Marc Goodfellow
1. College of Engineering, Mathematics and Physical Sciences, University of Exeter, Exeter, United Kingdom
2. Centre for Biomedical Modelling and Analysis, University of Exeter, Exeter, United Kingdom
3. EPSRC Centre for Predictive Modelling in Healthcare, University of Exeter, Exeter, United Kingdom
Contribution
MG, Analysis and interpretation of data, Drafting or revising the article

Competing interests
The authors declare that no competing interests exist.

"This ORCID iD identifies the author of this article:" 0000-0002-7282-7280
Tobias Galla

Theoretical Physics, School of Physics and Astronomy, University of Manchester, Manchester, United Kingdom

Contribution
TG, Analysis and interpretation of data, Drafting or revising the article

Competing interests
The authors declare that no competing interests exist.
Magnus Rattray

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
MR, Conception and design, Acquisition of data, Analysis and interpretation of data

Competing interests
The authors declare that no competing interests exist.
Nancy Papalopulu

Faculty of Biology, Medicine and Health, University of Manchester, Manchester, United Kingdom

Contribution
NP, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

For correspondence
nancy.papalopulu@manchester.ac.uk

Competing interests
The authors declare that no competing interests exist.

"This ORCID iD identifies the author of this article:" 0000-0001-6992-6870

Funding

Biotechnology and Biological Sciences Research Council (Doctoral Training Centre in Systems Biology)

Nick E Phillips

Wellcome Trust (Sir Henry Wellcome Postdoctoral Fellowship. 103986/Z/14/Z)

Cerys S Manning

Wellcome Trust (WT ISSF, Biomedical research Consortia grant, 097820/Z/11/B)

Tom Pettini

Medical Research Council (MR/K015885/1)

James Boyd
David G Spiller
Michael RH White

Biotechnology and Biological Sciences Research Council (David Phillips Fellowship, BB/I017976/1)

Pawel Paszek

Biotechnology and Biological Sciences Research Council (BB/K003097/1)

David G Spiller
Michael RH White

Engineering and Physical Sciences Research Council (EP/N014391/1)

Marc Goodfellow

Wellcome Trust (WT105618MA)

Marc Goodfellow

Wellcome Trust (090868/Z/09/Z)

Nancy Papalopulu

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

Acknowledgements

This work was supported by a Wellcome Trust Senior Research Fellowship to NP (090868/Z/09/Z), a Sir Henry Wellcome Fellowship to CM (103986/Z/14/Z), a Wellcome Trust Institutional Strategic Support Award (097820/Z/11/B) and a BBSRC Doctoral Training Centre in Systems Biology studentship to NEP. PP holds a BBSRC David Phillips Research Fellowship (BB/I017976/1). MG gratefully acknowledges the financial support of the EPSRC via grant EP/N014391/1. The contribution of MG was generously supported by a Wellcome Trust Institutional Strategic Support Award (WT105618MA). J Boyd was funded by MRC grant MR/K015885/1. DS and MW’s work is funded by an MRC grant MR/K015885/1 and a BBSRC grant BB/K003097/1. The authors would also like to thank the Biological Services Facility (BSF), the Bioimaging and Flow Cytometry Facilities of the University of Manchester for technical support, in particular to Dr Gareth Howell for expert advice on challenging FACS sorts. We thank Dr. Ximena Soto for advice and discussions and Dr Angelica Santiago-Gomez from Robert B Clarke’s group at the MCRC, Manchester for technical support and advice with western blotting. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Ethics

Animal experimentation: All animal work was performed under regulations set out by the UK Home Office Legislation under the 1986 United Kingdom Animal Scientific Procedures Act.

Copyright

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

3,065

views
582

downloads
44

citations

Views, downloads and citations are aggregated across all versions of this paper published by eLife.

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)

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

Mendeley

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

Nick E Phillips
Cerys S Manning
Tom Pettini
Veronica Biga
Elli Marinopoulou
Peter Stanley
James Boyd
James Bagnall
Pawel Paszek
David G Spiller
Michael RH White
Marc Goodfellow
Tobias Galla
Magnus Rattray
Nancy Papalopulu

(2016)

Stochasticity in the miR-9/Hes1 oscillatory network can account for clonal heterogeneity in the timing of differentiation

eLife 5:e16118.

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

Categories and tags

Research organism

Mouse

Share this article

Cite this article

Absolute quantification of mRNA copy number and active transcription sites.

Figure 1—source data 1

Figure 1—source data 2

Figure 1—source data 3

Figure 1—source data 4

Figure 1—source data 5

Quantifying HES1 protein number with FCS.

Figure 2—source data 1

Stochastic simulations of HES1 oscillations match the experimental data better.

Figure 3—source data 1

Figure 3—source data 2

The distribution of time-to-differentiation widens at low system sizes.

Figure 4—source data 1

Stochasticity causes systematic changes to the mean time-to-differentiation.

Figure 5—source data 1

The spread of time-to-differentiation decreases in both the computational model and experiments with increased initial miR-9.

Figure 6—source data 1

Figure 6—source data 2

Figure 6—source data 3

Stochasticity increases the parameter space where HES1 oscillates.

Figure 7—source data 1

Figure 7—source data 2

Hes1 mRNA is inherited unequally at cell division.

Figure 8—source data 1

Stochastic system is robust to the binomial distribution of constituents at cell division.

Figure 9—source data 1

Author details

Nick E Phillips

Contribution

Competing interests

Cerys S Manning

Contribution

Contributed equally with

Competing interests

Tom Pettini

Contribution

Contributed equally with

Competing interests

Veronica Biga

Contribution

Contributed equally with

Competing interests

Elli Marinopoulou

Contribution

Contributed equally with

Competing interests

Peter Stanley

Contribution

Contributed equally with

Competing interests

James Boyd

Contribution

Competing interests

James Bagnall

Contribution

Competing interests

Pawel Paszek

Contribution

Competing interests

David G Spiller

Contribution

Competing interests

Michael RH White

Contribution

Competing interests

Marc Goodfellow

Contribution

Competing interests

Tobias Galla

Contribution

Competing interests

Magnus Rattray

Contribution

Competing interests

Nancy Papalopulu

Contribution

For correspondence

Competing interests