Absolute rate theories of epigenetic stability

Irreversibility in bacterial regulatory networks

Irreversibility, in which a transient perturbation leaves a system in a new state, is an emergent property in systems of interacting entities. This property has well-established implications in statistical physics but remains underexplored in biological networks, especially for bacteria and other prokaryotes whose regulation of gene expression occurs predominantly at the transcriptional level. Focusing on the reconstructed regulatory network of Escherichia coli, we examine network responses to transient single-gene perturbations. We predict irreversibility in numerous cases and find that the incidence of irreversibility increases with the proximity of the perturbed gene to positive circuits in the network. Comparison with experimental data suggests a connection between the predicted irreversibility to transient perturbations and the evolutionary response to permanent perturbations.

Switching off: The phenotypic transition to the uninduced state of the lactose uptake pathway

The lactose uptake pathway of E. coli is a paradigmatic example of multistability in gene regulatory circuits. In the induced state of the lac pathway, the genes comprising the lac operon are transcribed, leading to the production of proteins that import and metabolize lactose. In the uninduced state, a stable repressor-DNA loop frequently blocks the transcription of the lac genes. Transitions from one phenotypic state to the other are driven by fluctuations, which arise from the random timing of the binding of ligands and proteins. This stochasticity affects transcription and translation, and ultimately molecular copy numbers. Our aim is to understand the transition from the induced to the uninduced state of the lac operon. We use a detailed computational model to show that repressor-operator binding and unbinding, fluctuations in the total number of repressors, and inducer-repressor binding and unbinding all play a role in this transition. Based on the timescales on which these processes operate, we construct a minimal model of the transition to the uninduced state and compare the results with simulations and experimental observations. The induced state turns out to be very stable, with a transition rate to the uninduced state lower than 2×10^-9 per minute. In contrast to the transition to the induced state, the transition to the uninduced state is well described in terms of a 2D diffusive system crossing a barrier, with the diffusion rates emerging from a model of repressor unbinding.

Mechanism of life-long maintenance of neuron identity despite molecular fluctuations

Cell fate is maintained over long timescales, yet molecular fluctuations can lead to spontaneous loss of this differentiated state. Our simulations identified a possible mechanism that explains life-long maintenance of ASE neuron fate in Caenorhabditis elegans by the terminal selector transcription factor CHE-1. Here, fluctuations in CHE-1 level are buffered by the reservoir of CHE-1 bound at its target promoters, which ensures continued che-1 expression by preferentially binding the che-1 promoter. We provide experimental evidence for this mechanism by showing that che-1 expression was resilient to induced transient CHE-1 depletion, while both expression of CHE-1 targets and ASE function were lost. We identified a 130 bp che-1 promoter fragment responsible for this resilience, with deletion of a homeodomain binding site in this fragment causing stochastic loss of ASE identity long after its determination. Because network architectures that support this mechanism are highly conserved in cell differentiation, it may explain stable cell fate maintenance in many systems.

Stochastic epigenetic dynamics of gene switching

Epigenetic modifications of histones crucially affect eukaryotic gene activity, while the epigenetic histone state is largely determined by the binding of specific factors such as the transcription factors (TFs) to DNA. Here, the way in which the TFs and the histone state are dynamically correlated is not obvious when the TF synthesis is regulated by the histone state. This type of feedback regulatory relation is ubiquitous in gene networks to determine cell fate in differentiation and other cell transformations. To gain insights into such dynamical feedback regulations, we theoretically analyze a model of epigenetic gene switching by extending the Doi-Peliti operator formalism of reaction kinetics to the problem of coupled molecular processes. Spin-1 and spin-1/2 coherent-state representations are introduced to describe stochastic reactions of histones and binding or unbinding of TFs in a unified way, which provides a concise view of the effects of timescale difference among these molecular processes; even in the case that binding or unbinding of TFs to or from DNA is adiabatically fast, the slow nonadiabatic histone dynamics gives rise to a distinct circular flow of the probability flux around basins in the landscape of the gene state distribution, which leads to hysteresis in gene switching. In contrast to the general belief that the change in the amount of TF precedes the histone state change, flux drives histones to be modified prior to the change in the amount of TF in self-regulating circuits. Flux-landscape analyses shed light on the nonlinear nonadiabatic mechanism of epigenetic cell fate decision making.

Reconstructing an epigenetic landscape using a genetic pulling approach

Cells use genetic switches to shift between alternate stable gene expression states, e.g., to adapt to new environments or to follow a developmental pathway. Conceptually, these stable phenotypes can be considered as attractive states on an epigenetic landscape with phenotypic changes being transitions between states. Measuring these transitions is challenging because they are both very rare in the absence of appropriate signals and very fast. As such, it has proved difficult to experimentally map the epigenetic landscapes that are widely believed to underly developmental networks. Here, we introduce a nonequilibrium perturbation method to help reconstruct a regulatory network's epigenetic landscape. We derive the mathematical theory needed and then use the method on simulated data to reconstruct the landscapes. Our results show that with a relatively small number of perturbation experiments it is possible to recover an accurate representation of the true epigenetic landscape. We propose that our theory provides a general method by which epigenetic landscapes can be studied. Finally, our theory suggests that the total perturbation impulse required to induce a switch between metastable states is a fundamental quantity in developmental dynamics.

Quantifying the Stability of Coupled Genetic and Epigenetic Switches With Variational Methods

The Waddington landscape provides an intuitive metaphor to view development as a ball rolling down the hill, with distinct phenotypes as basins and differentiation pathways as valleys. Since, at a molecular level, cell differentiation arises from interactions among the genes, a mathematical definition for the Waddington landscape can, in principle, be obtained by studying the gene regulatory networks. For eukaryotes, gene regulation is inextricably and intimately linked to histone modifications. However, the impact of such modifications on both landscape topography and stability of attractor states is not fully understood. In this work, we introduced a minimal kinetic model for gene regulation that combines the impact of both histone modifications and transcription factors. We further developed an approximation scheme based on variational principles to solve the corresponding master equation in a second quantized framework. By analyzing the steady-state solutions at various parameter regimes, we found that histone modification kinetics can significantly alter the behavior of a genetic network, resulting in qualitative changes in gene expression profiles. The emerging epigenetic landscape captures the delicate interplay between transcription factors and histone modifications in driving cell-fate decisions.

Dissipation in Non-Steady State Regulatory Circuits

In order to respond to environmental signals, cells often use small molecular circuits to transmit information about their surroundings. Recently, motivated by specific examples in signaling and gene regulation, a body of work has focused on the properties of circuits that function out of equilibrium and dissipate energy. We briefly review the probabilistic measures of information and dissipation and use simple models to discuss and illustrate trade-offs between information and dissipation in biological circuits. We find that circuits with non-steady state initial conditions can transmit more information at small readout delays than steady state circuits. The dissipative cost of this additional information proves marginal compared to the steady state dissipation. Feedback does not significantly increase the transmitted information for out of steady state circuits but does decrease dissipative costs. Lastly, we discuss the case of bursty gene regulatory circuits that, even in the fast switching limit, function out of equilibrium.

Stochastic analysis of genetic feedback controllers to reprogram a pluripotency gene regulatory network

Cellular reprogramming is traditionally accomplished through an open loop (OL) control approach, wherein key transcription factors (TFs) are injected in cells to steer the state of the pluripotency (PL) gene regulatory network (GRN), as encoded by TFs concentrations, to the pluripotent state. Due to the OL nature of this approach, the concentration of TFs cannot be accurately controlled. Recently, a closed loop (CL) feedback control strategy was proposed to overcome this problem with promising theoretical results. However, previous analyses of the controller were based on deterministic models. It is well known that cellular systems are characterized by substantial stochasticity, especially when molecules are in low copy number as it is the case in reprogramming problems wherein the gene copy number is usually one or two. Hence, in this paper, we analyze the Chemical Master Equation (CME) for the reaction model of the PL GRN with and without the feedback controller. We computationally and analytically investigate the performance of the controller in biologically relevant parameter regimes where stochastic effects dictate system dynamics. Our results indicate that the feedback control approach still ensures reprogramming even when both the PL GRN and the controller are stochastic.

Relatively slow stochastic gene-state switching in the presence of positive feedback significantly broadens the region of bimodality through stabilizing the uninduced phenotypic state

Within an isogenic population, even in the same extracellular environment, individual cells can exhibit various phenotypic states. The exact role of stochastic gene-state switching regulating the transition among these phenotypic states in a single cell is not fully understood, especially in the presence of positive feedback. Recent high-precision single-cell measurements showed that, at least in bacteria, switching in gene states is slow relative to the typical rates of active transcription and translation. Hence using the lac operon as an archetype, in such a region of operon-state switching, we present a fluctuating-rate model for this classical gene regulation module, incorporating the more realistic operon-state switching mechanism that was recently elucidated. We found that the positive feedback mechanism induces bistability (referred to as deterministic bistability), and that the parameter range for its occurrence is significantly broadened by stochastic operon-state switching. We further show that in the absence of positive feedback, operon-state switching must be extremely slow to trigger bistability by itself. However, in the presence of positive feedback, which stabilizes the induced state, the relatively slow operon-state switching kinetics within the physiological region are sufficient to stabilize the uninduced state, together generating a broadened parameter region of bistability (referred to as stochastic bistability). We illustrate the opposite phenotype-transition rate dependence upon the operon-state switching rates in the two types of bistability, with the aid of a recently proposed rate formula for fluctuating-rate models. The rate formula also predicts a maximal transition rate in the intermediate region of operon-state switching, which is validated by numerical simulations in our model. Overall, our findings suggest a biological function of transcriptional “variations” among genetically identical cells, for the emergence of bistability and transition between phenotypic states.

Identifying the mechanism underlying the coexistence of multiple stable phenotypic states has been a challenging scientific problem for more than half a century, and an appropriate mathematical model at the single-cell level is also in high demand. Single-cell measurements conducted in the past ten years have shown that gene-state switching is slow relative to the typical rates of active transcription and translation; hence the recently proposed fluctuating-rate model is a good candidate for describing the single-cell dynamics. We use the classic gene regulation module of the lac operon as an archetype and build a specific fluctuating-rate model based on the recently identified operon-state switching mechanism. This model is analyzed to dissect the interplay between positive feedback and the stochastic switching of gene states in the emergence of bistability/multistablity and the transition between phenotypic states. We show that relatively slow operon-state switching stabilizes the uninduced state and that the positive feedback stabilizes the induced state. Thus, the parameter range for bistability is significantly broadened. In addition, recently proposed landscape theory and rate formula predict opposite phenotype-transition rate dependence on operon-state switching rates for the two types of bistability.

A stochastic and dynamical view of pluripotency in mouse embryonic stem cells

Pluripotent embryonic stem cells are of paramount importance for biomedical sciences because of their innate ability for self-renewal and differentiation into all major cell lines. The fateful decision to exit or remain in the pluripotent state is regulated by complex genetic regulatory networks. The rapid growth of single-cell sequencing data has greatly stimulated applications of statistical and machine learning methods for inferring topologies of pluripotency regulating genetic networks. The inferred network topologies, however, often only encode Boolean information while remaining silent about the roles of dynamics and molecular stochasticity inherent in gene expression. Herein we develop a framework for systematically extending Boolean-level network topologies into higher resolution models of networks which explicitly account for the promoter architectures and gene state switching dynamics. We show the framework to be useful for disentangling the various contributions that gene switching, external signaling, and network topology make to the global heterogeneity and dynamics of transcription factor populations. We find the pluripotent state of the network to be a steady state which is robust to global variations of gene switching rates which we argue are a good proxy for epigenetic states of individual promoters. The temporal dynamics of exiting the pluripotent state, on the other hand, is significantly influenced by the rates of genetic switching which makes cells more responsive to changes in extracellular signals.

Stochastic timing in gene expression for simple regulatory strategies

Timing is essential for many cellular processes, from cellular responses to external stimuli to the cell cycle and circadian clocks. Many of these processes are based on gene expression. For example, an activated gene may be required to reach in a precise time a threshold level of expression that triggers a specific downstream process. However, gene expression is subject to stochastic fluctuations, naturally inducing an uncertainty in this threshold-crossing time with potential consequences on biological functions and phenotypes. Here, we consider such ‘timing fluctuations’ and we ask how they can be controlled. Our analytical estimates and simulations show that, for an induced gene, timing variability is minimal if the threshold level of expression is approximately half of the steady-state level. Timing fluctuations can be reduced by increasing the transcription rate, while they are insensitive to the translation rate. In presence of self-regulatory strategies, we show that self-repression reduces timing noise for threshold levels that have to be reached quickly, while self-activation is optimal at long times. These results lay a framework for understanding stochasticity of endogenous systems such as the cell cycle, as well as for the design of synthetic trigger circuits.

Markov State Models of gene regulatory networks

Background

Gene regulatory networks with dynamics characterized by multiple stable states underlie cell fate-decisions. Quantitative models that can link molecular-level knowledge of gene regulation to a global understanding of network dynamics have the potential to guide cell-reprogramming strategies. Networks are often modeled by the stochastic Chemical Master Equation, but methods for systematic identification of key properties of the global dynamics are currently lacking.

Results

The method identifies the number, phenotypes, and lifetimes of long-lived states for a set of common gene regulatory network models. Application of transition path theory to the constructed Markov State Model decomposes global dynamics into a set of dominant transition paths and associated relative probabilities for stochastic state-switching.

Conclusions

In this proof-of-concept study, we found that the Markov State Model provides a general framework for analyzing and visualizing stochastic multistability and state-transitions in gene networks. Our results suggest that this framework—adopted from the field of atomistic Molecular Dynamics—can be a useful tool for quantitative Systems Biology at the network scale.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-017-0394-4) contains supplementary material, which is available to authorized users.

A computational method for the investigation of multistable systems and its application to genetic switches

Background

Genetic switches exhibit multistability, form the basis of epigenetic memory, and are found in natural decision making systems, such as cell fate determination in developmental pathways. Synthetic genetic switches can be used for recording the presence of different environmental signals, for changing phenotype using synthetic inputs and as building blocks for higher-level sequential logic circuits. Understanding how multistable switches can be constructed and how they function within larger biological systems is therefore key to synthetic biology.

Results

Here we present a new computational tool, called StabilityFinder, that takes advantage of sequential Monte Carlo methods to identify regions of parameter space capable of producing multistable behaviour, while handling uncertainty in biochemical rate constants and initial conditions. The algorithm works by clustering trajectories in phase space, and iteratively minimizing a distance metric. Here we examine a collection of models of genetic switches, ranging from the deterministic Gardner toggle switch to stochastic models containing different positive feedback connections. We uncover the design principles behind making bistable, tristable and quadristable switches, and find that rate of gene expression is a key parameter. We demonstrate the ability of the framework to examine more complex systems and examine the design principles of a three gene switch. Our framework allows us to relax the assumptions that are often used in genetic switch models and we show that more complex abstractions are still capable of multistable behaviour.

Conclusions

Our results suggest many ways in which genetic switches can be enhanced and offer designs for the construction of novel switches. Our analysis also highlights subtle changes in correlation of experimentally tunable parameters that can lead to bifurcations in deterministic and stochastic systems. Overall we demonstrate that StabilityFinder will be a valuable tool in the future design and construction of novel gene networks.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-016-0375-z) contains supplementary material, which is available to authorized users.

Multiple coupled landscapes and non-adiabatic dynamics with applications to self-activating genes

Many physical, chemical and biochemical systems (e.g. electronic dynamics and gene regulatory networks) are governed by continuous stochastic processes (e.g. electron dynamics on a particular electronic energy surface and protein (gene product) synthesis) coupled with discrete processes (e.g. hopping among different electronic energy surfaces and on and off switching of genes). One can also think of the underlying dynamics as the continuous motion on a particular landscape and discrete hoppings among different landscapes. The main difference of such systems from the intra-landscape dynamics alone is the emergence of the timescale involved in transitions among different landscapes in addition to the timescale involved in a particular landscape. The adiabatic limit when inter-landscape hoppings are fast compared to continuous intra-landscape dynamics has been studied both analytically and numerically, but the analytical treatment of the non-adiabatic regime where the inter-landscape hoppings are slow or comparable to continuous intra-landscape dynamics remains challenging. In this study, we show that there exists mathematical mapping of the dynamics on 2(N) discretely coupled N continuous dimensional landscapes onto one single landscape in 2N dimensional extended continuous space. On this 2N dimensional landscape, eddy current emerges as a sign of non-equilibrium non-adiabatic dynamics and plays an important role in system evolution. Many interesting physical effects such as the enhancement of fluctuations, irreversibility, dissipation and optimal kinetics emerge due to non-adiabaticity manifested by the eddy current illustrated for an N = 1 self-activator. We further generalize our theory to the N-gene network with multiple binding sites and multiple synthesis rates for discretely coupled non-equilibrium stochastic physical and biological systems.

Contribution of Bistability and Noise to Cell Fate Transitions Determined by Feedback Opening

Alternative cell fates represent a form of non-genetic diversity, which can promote adaptation and functional specialization. It is difficult to predict the rate of the transition between two cell fates due to the strong effect of noise on feedback loops and missing parameters. We opened synthetic positive feedback loops experimentally to obtain open-loop functions. These functions allowed us to identify a deterministic model of bistability by bypassing noise and the requirement to resolve individual processes in the loop. Combining the open-loop function with kinetic measurements and reintroducing the measured noise, we were able to predict the transition rates for the feedback systems without parameter tuning. Noise in gene expression was the key determinant of the transition rates inside the bistable range. Transitions between two cell fates were also observed outside of the bistable range, evidenced by bimodality and hysteresis. In this case, a slow transient process was the rate-limiting step in the transitions. Thus, feedback opening is an effective approach to identify the determinants of cell fate transitions and to predict their rates.

Dichotomous noise models of gene switches

Molecular noise in gene regulatory networks has two intrinsic components, one part being due to fluctuations caused by the birth and death of protein or mRNA molecules which are often present in small numbers and the other part arising from gene state switching, a single molecule event. Stochastic dynamics of gene regulatory circuits appears to be largely responsible for bifurcations into a set of multi-attractor states that encode different cell phenotypes. The interplay of dichotomous single molecule gene noise with the nonlinear architecture of genetic networks generates rich and complex phenomena. In this paper, we elaborate on an approximate framework that leads to simple hybrid multi-scale schemes well suited for the quantitative exploration of the steady state properties of large-scale cellular genetic circuits. Through a path sum based analysis of trajectory statistics, we elucidate the connection of these hybrid schemes to the underlying master equation and provide a rigorous justification for using dichotomous noise based models to study genetic networks. Numerical simulations of circuit models reveal that the contribution of the genetic noise of single molecule origin to the total noise is significant for a wide range of kinetic regimes.

Extending the dynamic range of transcription factor action by translational regulation

A crucial step in the regulation of gene expression is binding of transcription factor (TF) proteins to regulatory sites along the DNA. But transcription factors act at nanomolar concentrations, and noise due to random arrival of these molecules at their binding sites can severely limit the precision of regulation. Recent work on the optimization of information flow through regulatory networks indicates that the lower end of the dynamic range of concentrations is simply inaccessible, overwhelmed by the impact of this noise. Motivated by the behavior of homeodomain proteins, such as the maternal morphogen Bicoid in the fruit fly embryo, we suggest a scheme in which transcription factors also act as indirect translational regulators, binding to the mRNA of other regulatory proteins. Intuitively, each mRNA molecule acts as an independent sensor of the input concentration, and averaging over these multiple sensors reduces the noise. We analyze information flow through this scheme and identify conditions under which it outperforms direct transcriptional regulation. Our results suggest that the dual role of homeodomain proteins is not just a historical accident, but a solution to a crucial physics problem in the regulation of gene expression.

Role of DNA binding sites and slow unbinding kinetics in titration-based oscillators

Genetic oscillators, such as circadian clocks, are constantly perturbed by molecular noise arising from the small number of molecules involved in gene regulation. One of the strongest sources of stochasticity is the binary noise that arises from the binding of a regulatory protein to a promoter in the chromosomal DNA. In this study, we focus on two minimal oscillators based on activator titration and repressor titration to understand the key parameters that are important for oscillations and for overcoming binary noise. We show that the rate of unbinding from the DNA, despite traditionally being considered a fast parameter, needs to be slow to broaden the space of oscillatory solutions. The addition of multiple, independent DNA binding sites further expands the oscillatory parameter space for the repressor-titration oscillator and lengthens the period of both oscillators. This effect is a combination of increased effective delay of the unbinding kinetics due to multiple binding sites and increased promoter ultrasensitivity that is specific for repression. We then use stochastic simulation to show that multiple binding sites increase the coherence of oscillations by mitigating the binary noise. Slow values of DNA unbinding rate are also effective in alleviating molecular noise due to the increased distance from the bifurcation point. Our work demonstrates how the number of DNA binding sites and slow unbinding kinetics, which are often omitted in biophysical models of gene circuits, can have a significant impact on the temporal and stochastic dynamics of genetic oscillators.

Effects of Collective Histone State Dynamics on Epigenetic Landscape and Kinetics of Cell Reprogramming

Cell reprogramming is a process of transitions from differentiated to pluripotent cell states via transient intermediate states. Within the epigenetic landscape framework, such a process is regarded as a sequence of transitions among basins on the landscape; therefore, theoretical construction of a model landscape which exhibits experimentally consistent dynamics can provide clues to understanding epigenetic mechanism of reprogramming. We propose a minimal gene-network model of the landscape, in which each gene is regulated by an integrated mechanism of transcription-factor binding/unbinding and the collective chemical modification of histones. We show that the slow collective variation of many histones around each gene locus alters topology of the landscape and significantly affects transition dynamics between basins. Differentiation and reprogramming follow different transition pathways on the calculated landscape, which should be verified experimentally via single-cell pursuit of the reprogramming process. Effects of modulation in collective histone state kinetics on transition dynamics and pathway are examined in search for an efficient protocol of reprogramming.

DNA-Binding Kinetics Determines the Mechanism of Noise-Induced Switching in Gene Networks

Gene regulatory networks are multistable dynamical systems in which attractor states represent cell phenotypes. Spontaneous, noise-induced transitions between these states are thought to underlie critical cellular processes, including cell developmental fate decisions, phenotypic plasticity in fluctuating environments, and carcinogenesis. As such, there is increasing interest in the development of theoretical and computational approaches that can shed light on the dynamics of these stochastic state transitions in multistable gene networks. We applied a numerical rare-event sampling algorithm to study transition paths of spontaneous noise-induced switching for a ubiquitous gene regulatory network motif, the bistable toggle switch, in which two mutually repressive genes compete for dominant expression. We find that the method can efficiently uncover detailed switching mechanisms that involve fluctuations both in occupancies of DNA regulatory sites and copy numbers of protein products. In addition, we show that the rate parameters governing binding and unbinding of regulatory proteins to DNA strongly influence the switching mechanism. In a regime of slow DNA-binding/unbinding kinetics, spontaneous switching occurs relatively frequently and is driven primarily by fluctuations in DNA-site occupancies. In contrast, in a regime of fast DNA-binding/unbinding kinetics, switching occurs rarely and is driven by fluctuations in levels of expressed protein. Our results demonstrate how spontaneous cell phenotype transitions involve collective behavior of both regulatory proteins and DNA. Computational approaches capable of simulating dynamics over many system variables are thus well suited to exploring dynamic mechanisms in gene networks.

Distinguishing the rates of gene activation from phenotypic variations

Background

Stochastic genetic switching driven by intrinsic noise is an important process in gene expression. When the rates of gene activation/inactivation are relatively slow, fast, or medium compared with the synthesis/degradation rates of mRNAs and proteins, the variability of protein and mRNA levels may exhibit very different dynamical patterns. It is desirable to provide a systematic approach to identify their key dynamical features in different regimes, aiming at distinguishing which regime a considered gene regulatory network is in from their phenotypic variations.

Results

We studied a gene expression model with positive feedbacks when genetic switching rates vary over a wide range. With the goal of providing a method to distinguish the regime of the switching rates, we first focus on understanding the essential dynamics of gene expression system in different cases. In the regime of slow switching rates, we found that the effective dynamics can be reduced to independent evolutions on two separate layers corresponding to gene activation and inactivation states, and the transitions between two layers are rare events, after which the system goes mainly along deterministic ODE trajectories on a particular layer to reach new steady states. The energy landscape in this regime can be well approximated by using Gaussian mixture model. In the regime of intermediate switching rates, we analyzed the mean switching time to investigate the stability of the system in different parameter ranges. We also discussed the case of fast switching rates from the viewpoint of transition state theory. Based on the obtained results, we made a proposal to distinguish these three regimes in a simulation experiment. We identified the intermediate regime from the fact that the strength of cellular memory is lower than the other two cases, and the fast and slow regimes can be distinguished by their different perturbation-response behavior with respect to the switching rates perturbations.

Conclusions

We proposed a simulation experiment to distinguish the slow, intermediate and fast regimes, which is the main point of our paper. In order to achieve this goal, we systematically studied the essential dynamics of gene expression system when the switching rates are in different regimes. Our theoretical understanding provides new insights on the gene expression experiments.

Electronic supplementary material

The online version of this article (doi:10.1186/s12918-015-0172-0) contains supplementary material, which is available to authorized users.

Stochastic phenotype transition of a single cell in an intermediate region of gene state switching

Multiple phenotypic states often arise in a single cell with different gene-expression states that undergo transcription regulation with positive feedback. Recent experiments show that, at least in E.coli, the gene state switching can be neither extremely slow nor exceedingly rapid as many previous theoretical treatments assumed. Rather, it is in the intermediate region which is difficult to handle mathematically. Under this condition, from a full chemical-master-equation description we derive a model in which the protein copy number, for a given gene state, follows a deterministic mean-field description while the protein-synthesis rates fluctuate due to stochastic gene state switching. The simplified kinetics yields a nonequilibrium landscape function, which, similar to the energy function for equilibrium fluctuation, provides the leading orders of fluctuations around each phenotypic state, as well as the transition rates between the two phenotypic states. This rate formula is analogous to Kramers' theory for chemical reactions. The resulting behaviors are significantly different from the two limiting cases studied previously.

What makes the lac-pathway switch: identifying the fluctuations that trigger phenotype switching in gene regulatory systems

Multistable gene regulatory systems sustain different levels of gene expression under identical external conditions. Such multistability is used to encode phenotypic states in processes including nutrient uptake and persistence in bacteria, fate selection in viral infection, cell-cycle control and development. Stochastic switching between different phenotypes can occur as the result of random fluctuations in molecular copy numbers of mRNA and proteins arising in transcription, translation, transport and binding. However, which component of a pathway triggers such a transition is generally not known. By linking single-cell experiments on the lactose-uptake pathway in E. coli to molecular simulations, we devise a general method to pinpoint the particular fluctuation driving phenotype switching and apply this method to the transition between the uninduced and induced states of the lac-genes. We find that the transition to the induced state is not caused only by the single event of lac-repressor unbinding, but depends crucially on the time period over which the repressor remains unbound from the lac-operon. We confirm this notion in strains with a high expression level of the lac-repressor (leading to shorter periods over which the lac-operon remains unbound), which show a reduced switching rate. Our techniques apply to multistable gene regulatory systems in general and allow to identify the molecular mechanisms behind stochastic transitions in gene regulatory circuits.

Construction of an effective landscape for multistate genetic switches

Multistate genetic switches play a crucial role during embryonic development and tumorigenesis. An archetypical example is the three-way switch regulating epithelial-hybrid-mesenchymal transitions. We devise a special WKB-based approach to investigate white Gaussian and shot noise effects on three-way switches, and construct an effective landscape in good quantitative agreement with stochastic simulations. This approach allows efficient analytical or numerical calculation of the landscape contours, the optimal path, and the state relative stability for general multicomponent multistate switches.

Stem cell differentiation as a many-body problem

Stem cell differentiation has been viewed as coming from transitions between attractors on an epigenetic landscape that governs the dynamics of a regulatory network involving many genes. Rigorous definition of such a landscape is made possible by the realization that gene regulation is stochastic, owing to the small copy number of the transcription factors that regulate gene expression and because of the single-molecule nature of the gene itself. We develop an approximation that allows the quantitative construction of the epigenetic landscape for large realistic model networks. Applying this approach to the network for embryonic stem cell development explains many experimental observations, including the heterogeneous distribution of the transcription factor Nanog and its role in safeguarding the stem cell pluripotency, which can be understood by finding stable steady-state attractors and the most probable transition paths between those attractors. We also demonstrate that the switching rate between attractors can be significantly influenced by the gene expression noise arising from the fluctuations of DNA occupancy when binding to a specific DNA site is slow.

Epigenetic state network approach for describing cell phenotypic transitions

Recent breakthroughs of cell phenotype reprogramming impose theoretical challenges on unravelling the complexity of large circuits maintaining cell phenotypes coupled at many different epigenetic and gene regulation levels, and quantitatively describing the phenotypic transition dynamics. A popular picture proposed by Waddington views cell differentiation as a ball sliding down a landscape with valleys corresponding to different cell types separated by ridges. Based on theories of dynamical systems, we establish a novel 'epigenetic state network' framework that captures the global architecture of cell phenotypes, which allows us to translate the metaphorical low-dimensional Waddington epigenetic landscape concept into a simple-yet-predictive rigorous mathematical framework of cell phenotypic transitions. Specifically, we simplify a high-dimensional epigenetic landscape into a collection of discrete states corresponding to stable cell phenotypes connected by optimal transition pathways among them. We then apply the approach to the phenotypic transition processes among fibroblasts (FBs), pluripotent stem cells (PSCs) and cardiomyocytes (CMs). The epigenetic state network for this case predicts three major transition pathways connecting FBs and CMs. One goes by way of PSCs. The other two pathways involve transdifferentiation either indirectly through cardiac progenitor cells or directly from FB to CM. The predicted pathways and multiple intermediate states are supported by existing microarray data and other experiments. Our approach provides a theoretical framework for studying cell phenotypic transitions. Future studies at single-cell levels can directly test the model predictions.

Non-equilibrium transition state rate theory

Transition state or Kramers' rate theory has been used to quantify the kinetic speed of many chemical, physical and biological equilibrium processes successfully.

Transcription factor binding kinetics constrain noise suppression via negative feedback

Negative autoregulation, where a transcription factor regulates its own expression by preventing transcription, is commonly used to suppress fluctuations in gene expression. Recent single molecule in vivo imaging has shown that it takes significant time for a transcription factor molecule to bind its chromosomal binding site. Given the slow association kinetics, transcription factor mediated feedback cannot at the same time be fast and strong. Here we show that with a limited association rate follows an optimal transcription factor binding strength where noise is maximally suppressed. At the optimal binding strength the binding site is free a fixed fraction of the time independent of the transcription factor concentration. One consequence is that high-copy number transcription factors should bind weakly to their operators, which is observed for transcription factors in Escherichia coli. The results demonstrate that a binding site's strength may be uncorrelated to its functional importance.

Time scales in epigenetic dynamics and phenotypic heterogeneity of embryonic stem cells

A remarkable feature of the self-renewing population of embryonic stem cells (ESCs) is their phenotypic heterogeneity: Nanog and other marker proteins of ESCs show large cell-to-cell variation in their expression level, which should significantly influence the differentiation process of individual cells. The molecular mechanism and biological implication of this heterogeneity, however, still remain elusive. We address this problem by constructing a model of the core gene-network of mouse ESCs. The model takes account of processes of binding/unbinding of transcription factors, formation/dissolution of transcription apparatus, and modification of histone code at each locus of genes in the network. These processes are hierarchically interrelated to each other forming the dynamical feedback loops. By simulating stochastic dynamics of this model, we show that the phenotypic heterogeneity of ESCs can be explained when the chromatin at the Nanog locus undergoes the large scale reorganization in formation/dissolution of transcription apparatus, which should have the timescale similar to the cell cycle period. With this slow transcriptional switching of Nanog, the simulated ESCs fluctuate among multiple transient states, which can trigger the differentiation into the lineage-specific cell states. From the simulated transitions among cell states, the epigenetic landscape underlying transitions is calculated. The slow Nanog switching gives rise to the wide basin of ESC states in the landscape. The bimodal Nanog distribution arising from the kinetic flow running through this ESC basin prevents transdifferentiation and promotes the definite decision of the cell fate. These results show that the distribution of timescales of the regulatory processes is decisively important to characterize the fluctuation of cells and their differentiation process. The analyses through the epigenetic landscape and the kinetic flow on the landscape should provide a guideline to engineer cell differentiation.

Embryonic stem cells (ESCs) can proliferate indefinitely by keeping pluripotency, i.e., the ability to differentiate into any cell-lineage. ESCs, therefore, have been the focus of intense biological and medical interests. A remarkable feature of ESCs is their phenotypic heterogeneity: ESCs show large cell-to-cell fluctuation in the expression level of Nanog, which is a key factor to maintain pluripotency. Since Nanog regulates many genes in ESCs, this fluctuation should seriously affect individual cells when they start differentiation. In this paper we analyze this phenotypic fluctuation by simulating the stochastic dynamics of gene network in ESCs. The model takes account of the mutually interrelated processes of gene regulation such as binding/unbinding of transcription factors, formation/dissolution of transcription apparatus, and histone-code modification. We show the distribution of timescales of these processes is decisively important to characterize the dynamical behavior of the gene network, and that the slow formation/dissolution of transcription apparatus at the Nanog locus explains the observed large fluctuation of ESCs. The epigenetic landscapes are calculated based on the stochastic simulation, and the role of the phenotypic fluctuation in the differentiation process is analyzed through the landscape picture.

Quantifying Waddington landscapes and paths of non-adiabatic cell fate decisions for differentiation, reprogramming and transdifferentiation

Cellular differentiation, reprogramming and transdifferentiation are determined by underlying gene regulatory networks. Non-adiabatic regulation via slow binding/unbinding to the gene can be important in these cell fate decision-making processes. Based on a stem cell core gene network, we uncovered the stem cell developmental landscape. As the binding/unbinding speed decreases, the landscape topography changes from bistable attractors of stem and differentiated states to more attractors of stem and other different cell states as well as substates. Non-adiabaticity leads to more differentiated cell types and provides a natural explanation for the heterogeneity observed in the experiments. We quantified Waddington landscapes with two possible cell fate decision mechanisms by changing the regulation strength or regulation timescale (non-adiabaticity). Transition rates correlate with landscape topography through barrier heights between different states and quantitatively determine global stability. We found the optimal speeds of these cell fate decision-making processes. We quantified biological paths and predict that differentiation and reprogramming go through an intermediate state (IM1), whereas transdifferentiation goes through another intermediate state (IM2). Some predictions are confirmed by recent experimental studies.

Metastable behavior in Markov processes with internal states

A perturbation framework is developed to analyze metastable behavior in stochastic processes with random internal and external states. The process is assumed to be under weak noise conditions, and the case where the deterministic limit is bistable is considered. A general analytical approximation is derived for the stationary probability density and the mean switching time between metastable states, which includes the pre exponential factor. The results are illustrated with a model of gene expression that displays bistable switching. In this model, the external state represents the number of protein molecules produced by a hypothetical gene. Once produced, a protein is eventually degraded. The internal state represents the activated or unactivated state of the gene; in the activated state the gene produces protein more rapidly than the unactivated state. The gene is activated by a dimer of the protein it produces so that the activation rate depends on the current protein level. This is a well studied model, and several model reductions and diffusion approximation methods are available to analyze its behavior. However, it is unclear if these methods accurately approximate long-time metastable behavior (i.e., mean switching time between metastable states of the bistable system). Diffusion approximations are generally known to fail in this regard.

Eddy current and coupled landscapes for nonadiabatic and nonequilibrium complex system dynamics

Physical and biological systems are often involved with coupled processes of different time scales. In the system with electronic and atomic motions, for example, the interplay between the atomic motion along the same energy landscape and the electronic hopping between different landscapes is critical: the system behavior largely depends on whether the intralandscape motion is slower (adiabatic) or faster (nonadiabatic) than the interlandscape hopping. For general nonequilibrium dynamics where Hamiltonian or energy function is unknown a priori, the challenge is how to extend the concepts of the intra- and interlandscape dynamics. In this paper we establish a theoretical framework for describing global nonequilibrium and nonadiabatic complex system dynamics by transforming the coupled landscapes into a single landscape but with additional dimensions. On this single landscape, dynamics is driven by gradient of the potential landscape, which is closely related to the steady-state probability distribution of the enlarged dimensions, and the probability flux, which has a curl nature. Through an example of a self-regulating gene circuit, we show that the curl flux has dramatic effects on gene regulatory dynamics. The curl flux and landscape framework developed here are easy to visualize and can be used to guide further investigation of physical and biological nonequilibrium systems.

Extrinsic noise driven phenotype switching in a self-regulating gene

Analysis of complex gene regulation networks gives rise to a landscape of metastable phenotypic states for cells. Heterogeneity within a population arises due to infrequent noise-driven transitions of individual cells between nearby metastable states. While most previous works have focused on the role of intrinsic fluctuations in driving such transitions, in this Letter we investigate the role of extrinsic fluctuations. First, we develop an analytical framework to study the combined effect of intrinsic and extrinsic noise on a toy model of a Boolean regulated genetic switch. We then extend these ideas to a more biologically relevant model with a Hill-like regulatory function. Employing our theory and Monte Carlo simulations, we show that extrinsic noise can significantly alter the lifetimes of the phenotypic states and may fundamentally change the escape mechanism. Finally, our theory can be readily generalized to more complex decision making networks in biology.

Positive feedback produces broad distributions in maximum activation attained within a narrow time window in stochastic biochemical reactions

How do single cell fate decisions induced by activation of key signaling proteins above threshold concentrations within a time interval are affected by stochastic fluctuations in biochemical reactions? We address this question using minimal models of stochastic chemical reactions commonly found in cell signaling and gene regulatory systems. Employing exact solutions and semi-analytical methods we calculate distributions of the maximum value (N) of activated species concentrations (P(max)(N)) and the time (t) taken to reach the maximum value (P(max)(t)) within a time interval in the minimal models. We find, the presence of positive feedback interactions make P(max)(N) more spread out with a higher "peakedness" in P(max)(t). Thus positive feedback interactions may help single cells to respond sensitively to a stimulus when cell decision processes require upregulation of activated forms of key proteins to a threshold number within a time window.

MiRNAs confer phenotypic robustness to gene networks by suppressing biological noise

miRNAs are small non-coding RNAs able to modulate target gene expression. It has been postulated that miRNAs confer robustness to biological processes, but clear experimental evidence is still missing. Here, using a synthetic biological approach, we demonstrate that microRNAs provide phenotypic robustness to transcriptional regulatory networks by buffering fluctuations in protein levels. We construct a network motif in mammalian cells exhibiting a 'toggle-switch' phenotype in which two alternative protein expression levels define its ON and OFF states. The motif consists of an inducible transcription factor that self-regulates its own transcription and that of a miRNA against the transcription factor itself. We confirm, using mathematical modelling and experimental approaches, that the microRNA confers robustness to the toggle-switch by enabling the cell to maintain and transmit its state. When absent, a dramatic increase in protein noise level occurs, causing the cell to randomly switch between the two states.

Predicting rates of cell state change caused by stochastic fluctuations using a data-driven landscape model

We develop a potential landscape approach to quantitatively describe experimental data from a fibroblast cell line that exhibits a wide range of GFP expression levels under the control of the promoter for tenascin-C. Time-lapse live-cell microscopy provides data about short-term fluctuations in promoter activity, and flow cytometry measurements provide data about the long-term kinetics, because isolated subpopulations of cells relax from a relatively narrow distribution of GFP expression back to the original broad distribution of responses. The landscape is obtained from the steady state distribution of GFP expression and connected to a potential-like function using a stochastic differential equation description (Langevin/Fokker-Planck). The range of cell states is constrained by a force that is proportional to the gradient of the potential, and biochemical noise causes movement of cells within the landscape. Analyzing the mean square displacement of GFP intensity changes in live cells indicates that these fluctuations are described by a single diffusion constant in log GFP space. This finding allows application of the Kramers' model to calculate rates of switching between two attractor states and enables an accurate simulation of the dynamics of relaxation back to the steady state with no adjustable parameters. With this approach, it is possible to use the steady state distribution of phenotypes and a quantitative description of the short-term fluctuations in individual cells to accurately predict the rates at which different phenotypes will arise from an isolated subpopulation of cells.

Type of noise defines global attractors in bistable molecular regulatory systems

The aim of this study is to demonstrate that in molecular dynamical systems with the underlying bi- or multistability, the type of noise determines the most strongly attracting steady state or stochastic attractor. As an example we consider a simple stochastic model of autoregulatory gene with a nonlinear positive feedback, which in the deterministic approximation has two stable steady state solutions. Three types of noise are considered: transcriptional and translational - due to the small number of gene product molecules and the gene switching noise - due to gene activation and inactivation transitions. We demonstrate that the type of noise in addition to the noise magnitude dictates the allocation of probability mass between the two stable steady states. In particular, we found that when the gene switching noise dominates over the transcriptional and translational noise (which is characteristic of eukaryotes), the gene preferentially activates, while in the opposite case, when the transcriptional noise dominates (which is characteristic of prokaryotes) the gene preferentially remains inactive. Moreover, even in the zero-noise limit, when the probability mass generically concentrates in the vicinity of one of two steady states, the choice of the most strongly attracting steady state is noise type-dependent. Although the epigenetic attractors are defined with the aid of the deterministic approximation of the stochastic regulatory process, their relative attractivity is controlled by the type of noise, in addition to noise magnitude. Since noise characteristics vary during the cell cycle and development, such mode of regulation can be potentially employed by cells to switch between alternative epigenetic attractors.

Influence of decoys on the noise and dynamics of gene expression

Many transcription factors bind to DNA with a remarkable lack of specificity, so that regulatory binding sites compete with an enormous number of nonregulatory "decoy" sites. For an autoregulated gene, we show decoy sites decrease noise in the number of unbound proteins to a Poisson limit that results from binding and unbinding. This noise buffering is optimized for a given protein concentration when decoys have a 1/2 probability of being occupied. Decoys linearly increase the time to approach steady state and exponentially increase the time to switch epigenetically between bistable states.

A new mechanism of stem cell differentiation through slow binding/unbinding of regulators to genes

Understanding differentiation, a biological process from a multipotent stem or progenitor state to a mature cell is critically important. We developed a theoretical framework to quantify the underlying potential landscape and pathways for cell development and differentiation. We proposed a new mechanism of differentiation and found the differentiated states can emerge from the slow binding/unbinding of regulatory proteins to gene promoters. With slow promoter binding/unbinding, we found multiple meta-stable differentiated states, which can explain the origin of multiple states observed in recent experiments. The kinetic time for the differentiation and reprogramming strongly depends on the time scale of the promoter binding/unbinding processes. We discovered an optimal speed for differentiation for certain promoter binding/unbinding rates. Future experiments might be able to tell if cells differentiate at that optimal speed. We also quantified irreversible kinetic pathways for the differentiation and reprogramming, which captures the non-equilibrium dynamics in multipotent stem or progenitor cells.

Biomolecular dynamics: order-disorder transitions and energy landscapes

While the energy landscape theory of protein folding is now a widely accepted view for understanding how relatively weak molecular interactions lead to rapid and cooperative protein folding, such a framework must be extended to describe the large-scale functional motions observed in molecular machines. In this review, we discuss (1) the development of the energy landscape theory of biomolecular folding, (2) recent advances toward establishing a consistent understanding of folding and function and (3) emerging themes in the functional motions of enzymes, biomolecular motors and other biomolecular machines. Recent theoretical, computational and experimental lines of investigation have provided a very dynamic picture of biomolecular motion. In contrast to earlier ideas, where molecular machines were thought to function similarly to macroscopic machines, with rigid components that move along a few degrees of freedom in a deterministic fashion, biomolecular complexes are only marginally stable. Since the stabilizing contribution of each atomic interaction is on the order of the thermal fluctuations in solution, the rigid body description of molecular function must be revisited. An emerging theme is that functional motions encompass order-disorder transitions and structural flexibility provides significant contributions to the free energy. In this review, we describe the biological importance of order-disorder transitions and discuss the statistical-mechanical foundation of theoretical approaches that can characterize such transitions.

Landscape and global stability of nonadiabatic and adiabatic oscillations in a gene network

We quantify the potential landscape to determine the global stability and coherence of biological oscillations. We explore a gene network motif in our experimental synthetic biology studies of two genes that mutually repress and activate each other with self-activation and self-repression. We find that in addition to intrinsic molecular number fluctuations, there is another type of fluctuation crucial for biological function: the fluctuation due to the slow binding/unbinding of protein regulators to gene promoters. We find that coherent limit cycle oscillations emerge in two regimes: an adiabatic regime with fast binding/unbinding and a nonadiabatic regime with slow binding/unbinding relative to protein synthesis/degradation. This leads to two mechanisms of producing the stable oscillations: the effective interactions from averaging the gene states in the adiabatic regime; and the time delays due to slow binding/unbinding to promoters in the nonadiabatic regime, which can be tested by forthcoming experiments. In both regimes, the landscape has a topological shape of the Mexican hat in protein concentrations that quantitatively determines the global stability of limit cycle dynamics. The oscillation coherence is shown to be correlated with the shape of the Mexican hat characterized by the height from the oscillation ring to the central top. The oscillation period can be tuned in a wide range by changing the binding/unbinding rate without changing the amplitude much, which is important for the functionality of a biological clock. A negative feedback loop with time delays due to slow binding/unbinding can also generate oscillations. Although positive feedback is not necessary for generating oscillations, it can make the oscillations more robust.

Protein fluxes along the filopodium as a framework for understanding the growth-retraction dynamics: the interplay between diffusion and active transport

We present a picture of filopodial growth and retraction from physics perspective, where we emphasize the significance of the role played by protein fluxes due to spatially extended nature of the filopodium. We review a series of works, which used stochastic simulations and mean field analytical modeling to find the concentration profile of G-actin inside a filopodium, which, in turn, determines the stationary filopodial length. In addition to extensively reviewing the prior works, we also report some new results on the role of active transport in regulating the length of filopodia. We model a filopodium where delivery of actin monomers towards the tip can occur both through passive diffusion and active transport by myosin motors. We found that the concentration profile of G-actin along the filopodium is rather non-trivial, containing a narrow minimum near the base followed by a broad maximum. For efficient enough actin transport, this non-monotonous shape is expected to occur under a broad set of conditions. We also raise the issue of slow approach to the stationary length and the possibility of multiple steady state solutions.

Stability and multiattractor dynamics of a toggle switch based on a two-stage model of stochastic gene expression

A toggle switch consists of two genes that mutually repress each other. This regulatory motif is active during cell differentiation and is thought to act as a memory device, being able to choose and maintain cell fate decisions. Commonly, this switch has been modeled in a deterministic framework where transcription and translation are lumped together. In this description, bistability occurs for transcription factor cooperativity, whereas autoactivation leads to a tristable system with an additional undecided state. In this contribution, we study the stability and dynamics of a two-stage gene expression switch within a probabilistic framework inspired by the properties of the Pu/Gata toggle switch in myeloid progenitor cells. We focus on low mRNA numbers, high protein abundance, and monomeric transcription-factor binding. Contrary to the expectation from a deterministic description, this switch shows complex multiattractor dynamics without autoactivation and cooperativity. Most importantly, the four attractors of the system, which only emerge in a probabilistic two-stage description, can be identified with committed and primed states in cell differentiation. To begin, we study the dynamics of the system and infer the mechanisms that move the system between attractors using both the quasipotential and the probability flux of the system. Next, we show that the residence times of the system in one of the committed attractors are geometrically distributed. We derive an analytical expression for the parameter of the geometric distribution, therefore completely describing the statistics of the switching process and elucidate the influence of the system parameters on the residence time. Moreover, we find that the mean residence time increases linearly with the mean protein level. This scaling also holds for a one-stage scenario and for autoactivation. Finally, we study the implications of this distribution for the stability of a switch and discuss the influence of the stability on a specific cell differentiation mechanism. Our model explains lineage priming and proposes the need of either high protein numbers or long-term modifications such as chromatin remodeling to achieve stable cell fate decisions. Notably, we present a system with high protein abundance that nevertheless requires a probabilistic description to exhibit multistability, complex switching dynamics, and lineage priming.