Deviant effects in molecular reaction pathways

Extinction scenarios in evolutionary processes: a multinomial Wright-Fisher approach

We study a discrete-time multi-type Wright-Fisher population process. The mean-field dynamics of the stochastic process is induced by a general replicator difference equation. We prove several results regarding the asymptotic behavior of the model, focusing on the impact of the mean-field dynamics on it. One of the results is a limit theorem that describes sufficient conditions for an almost certain path to extinction, first eliminating the type which is the least fit at the mean-field equilibrium. The effect is explained by the metastability of the stochastic system, which under the conditions of the theorem spends almost all time before the extinction event in a neighborhood of the equilibrium. In addition to the limit theorems, we propose a maximization principle for a general deterministic replicator dynamics and study its implications for the stochastic model.

Out-of-equilibrium gene expression fluctuations in the presence of extrinsic noise

Cell-to-cell variability in protein concentrations is strongly affected by extrinsic noise, especially for highly expressed genes. Extrinsic noise can be due to fluctuations of several possible cellular factors connected to cell physiology and to the level of key enzymes in the expression process. However, how to identify the predominant sources of extrinsic noise in a biological system is still an open question. This work considers a general stochastic model of gene expression with extrinsic noise represented as fluctuations of the different model rates, and focuses on the out-of-equilibrium expression dynamics. Combining analytical calculations with stochastic simulations, we characterize how extrinsic noise shapes the protein variability during gene activation or inactivation, depending on the prevailing source of extrinsic variability, on its intensity and timescale. In particular, we show that qualitatively different noise profiles can be identified depending on which are the fluctuating parameters. This indicates an experimentally accessible way to pinpoint the dominant sources of extrinsic noise using time-coarse experiments.

A Comparison of Weighted Stochastic Simulation Methods for the Analysis of Genetic Circuits

Rare events are of particular interest in synthetic biology because rare biochemical events may be catastrophic to a biological system by, for example, triggering irreversible events such as off-target drug delivery. To estimate the probability of rare events efficiently, several weighted stochastic simulation methods have been developed. Under optimal parameters and model conditions, these methods can greatly improve simulation efficiency in comparison to traditional stochastic simulation. Unfortunately, the optimal parameters and conditions cannot be deduced a priori. This paper presents a critical survey of weighted stochastic simulation methods. It shows that the methods considered here cannot consistently, efficiently, and exactly accomplish the task of rare event simulation without resorting to a computationally expensive calibration procedure, which undermines their overall efficiency. The results suggest that further development is needed before these methods can be deployed for general use in biological simulations.

Operant conditioning of stochastic chemical reaction networks

Adapting one's behavior to environmental conditions and past experience is a key trait of living systems. In the biological world, there is evidence for adaptive behaviors such as learning even in naturally occurring, non-neural, single-celled organisms. In the bioengineered world, advances in synthetic cell engineering and biorobotics have created the possibility of implementing lifelike systems engineered from the bottom up. This will require the development of programmable control circuitry for such biomimetic systems that is capable of realizing such non-trivial and adaptive behavior, including modification of subsequent behavior in response to environmental feedback. To this end, we report the design of novel stochastic chemical reaction networks capable of probabilistic decision-making in response to stimuli. We show that a simple chemical reaction network motif can be tuned to produce arbitrary decision probabilities when choosing between two or more responses to a stimulus signal. We further show that simple feedback mechanisms from the environment can modify these probabilities over time, enabling the system to adapt its behavior dynamically in response to positive or negative reinforcement based on its decisions. This system thus acts as a form of operant conditioning of the chemical circuit, in the sense that feedback provided based on decisions taken by the circuit form the basis of the learning process. Our work thus demonstrates that simple chemical systems can be used to implement lifelike behavior in engineered biomimetic systems.

Stochastic self-assembly of reaction-diffusion patterns in synaptic membranes

Synaptic receptor and scaffold molecules self-assemble into membrane protein domains, which play an important role in signal transmission across chemical synapses. Experiment and theory have shown that the formation of receptor-scaffold domains of the characteristic size observed in nerve cells can be understood from the receptor and scaffold reaction and diffusion processes suggested by experiments. We employ here kinetic Monte Carlo (KMC) simulations to explore the self-assembly of synaptic receptor-scaffold domains in a stochastic lattice model of receptor and scaffold reaction-diffusion dynamics. For reaction and diffusion rates within the ranges of values suggested by experiments we find, in agreement with previous mean-field calculations, self-assembly of receptor-scaffold domains of a size similar to that observed in experiments. Comparisons between the results of our KMC simulations and mean-field solutions suggest that the intrinsic noise associated with receptor and scaffold reaction and diffusion processes accelerates the self-assembly of receptor-scaffold domains, and confers increased robustness to domain formation. In agreement with experimental observations, our KMC simulations yield a prevalence of scaffolds over receptors in receptor-scaffold domains. Our KMC simulations show that receptor and scaffold reaction-diffusion dynamics can inherently give rise to plasticity in the overall properties of receptor-scaffold domains, which may be utilized by nerve cells to regulate the receptor number at chemical synapses.

Quantifying the roles of space and stochasticity in computer simulations for cell biology and cellular biochemistry

Most of the fascinating phenomena studied in cell biology emerge from interactions among highly organized multimolecular structures embedded into complex and frequently dynamic cellular morphologies. For the exploration of such systems, computer simulation has proved to be an invaluable tool, and many researchers in this field have developed sophisticated computational models for application to specific cell biological questions. However, it is often difficult to reconcile conflicting computational results that use different approaches to describe the same phenomenon. To address this issue systematically, we have defined a series of computational test cases ranging from very simple to moderately complex, varying key features of dimensionality, reaction type, reaction speed, crowding, and cell size. We then quantified how explicit spatial and/or stochastic implementations alter outcomes, even when all methods use the same reaction network, rates, and concentrations. For simple cases, we generally find minor differences in solutions of the same problem. However, we observe increasing discordance as the effects of localization, dimensionality reduction, and irreversible enzymatic reactions are combined. We discuss the strengths and limitations of commonly used computational approaches for exploring cell biological questions and provide a framework for decision making by researchers developing new models. As computational power and speed continue to increase at a remarkable rate, the dream of a fully comprehensive computational model of a living cell may be drawing closer to reality, but our analysis demonstrates that it will be crucial to evaluate the accuracy of such models critically and systematically.

Stochastic fluctuations in apoptotic threshold of tumour cells can enhance apoptosis and combat fractional killing

Fractional killing, which is a significant impediment to successful chemotherapy, is observed even in a population of genetically identical cancer cells exposed to apoptosis-inducing agents. This phenomenon arises not from genetic mutation but from cell-to-cell variation in the activation timing and level of the proteins that regulates apoptosis. To understand the mechanism behind the phenomenon, we formulate complex fractional killing processes as a first-passage time (FPT) problem with a stochastically fluctuating boundary. Analytical calculations are performed for the FPT distribution in a toy model of stochastic p53 gene expression, where the cancer cell is killed only when the p53 expression level crosses an active apoptotic threshold. Counterintuitively, we find that threshold fluctuations can effectively enhance cellular killing by significantly decreasing the mean time that the p53 protein reaches the threshold level for the first time. Moreover, faster fluctuations lead to the killing of more cells. These qualitative results imply that fluctuations in threshold are a non-negligible stochastic source, and can be taken as a strategy for combating fractional killing of cancer cells.

Chemical Systems Involving Two Competitive Self-Catalytic Reactions

Self-catalytic reactions are chemical phenomena, in which a product catalyzes the reactions of substrates further to yield products. A significant amplification of product concentration occurs during the reactions in a dilute solution, which exhibit notable properties such as sigmoidal kinetics, seeding effects, and thermal hysteresis. Chemical systems involving two competitive self-catalytic reactions can be considered, in which the competitive formation of two products occurs, which is affected by environmental changes, subtle perturbations, and fluctuations, and notable chemical phenomena appear such as formation of different structures in response to slow/fast temperature changes, chiral symmetry breaking, shortcut in reaction time, homogeneous-heterogeneous transitions, and mechanical responses. Studies on such chemical systems provide understanding on biological systems and can also be extended to the development of novel functional materials.

Chemical Kinetics Roots and Methods to Obtain the Probability Distribution Function Evolution of Reactants and Products in Chemical Networks Governed by a Master Equation

In this paper first, we review the physical root bases of chemical reaction networks as a Markov process in multidimensional vector space. Then we study the chemical reactions from a microscopic point of view, to obtain the expression for the propensities for the different reactions that can happen in the network. These chemical propensities, at a given time, depend on the system state at that time, and do not depend on the state at an earlier time indicating that we are dealing with Markov processes. Then the Chemical Master Equation (CME) is deduced for an arbitrary chemical network from a probability balance and it is expressed in terms of the reaction propensities. This CME governs the dynamics of the chemical system. Due to the difficulty to solve this equation two methods are studied, the first one is the probability generating function method or z-transform, which permits to obtain the evolution of the factorial moment of the system with time in an easiest way or after some manipulation the evolution of the polynomial moments. The second method studied is the expansion of the CME in terms of an order parameter (system volume). In this case we study first the expansion of the CME using the propensities obtained previously and splitting the molecular concentration into a deterministic part and a random part. An expression in terms of multinomial coefficients is obtained for the evolution of the probability of the random part. Then we study how to reconstruct the probability distribution from the moments using the maximum entropy principle. Finally, the previous methods are applied to simple chemical networks and the consistency of these methods is studied.

Reducing the extinction risk of stochastic populations via nondemographic noise

We consider nondemographic noise in the form of uncertainty in the reaction step size and reveal a dramatic effect this noise may have on the stability of self-regulating populations. Employing the reaction scheme mA→kA but allowing, e.g., the product number k to be a priori unknown and sampled from a given distribution, we show that such nondemographic noise can greatly reduce the population's extinction risk compared to the fixed k case. Our analysis is tested against numerical simulations, and by using empirical data of different species, we argue that certain distributions may be more evolutionary beneficial than others.

Stochastic lattice model of synaptic membrane protein domains

Neurotransmitter receptor molecules, concentrated in synaptic membrane domains along with scaffolds and other kinds of proteins, are crucial for signal transmission across chemical synapses. In common with other membrane protein domains, synaptic domains are characterized by low protein copy numbers and protein crowding, with rapid stochastic turnover of individual molecules. We study here in detail a stochastic lattice model of the receptor-scaffold reaction-diffusion dynamics at synaptic domains that was found previously to capture, at the mean-field level, the self-assembly, stability, and characteristic size of synaptic domains observed in experiments. We show that our stochastic lattice model yields quantitative agreement with mean-field models of nonlinear diffusion in crowded membranes. Through a combination of analytic and numerical solutions of the master equation governing the reaction dynamics at synaptic domains, together with kinetic Monte Carlo simulations, we find substantial discrepancies between mean-field and stochastic models for the reaction dynamics at synaptic domains. Based on the reaction and diffusion properties of synaptic receptors and scaffolds suggested by previous experiments and mean-field calculations, we show that the stochastic reaction-diffusion dynamics of synaptic receptors and scaffolds provide a simple physical mechanism for collective fluctuations in synaptic domains, the molecular turnover observed at synaptic domains, key features of the observed single-molecule trajectories, and spatial heterogeneity in the effective rates at which receptors and scaffolds are recycled at the cell membrane. Our work sheds light on the physical mechanisms and principles linking the collective properties of membrane protein domains to the stochastic dynamics that rule their molecular components.

An open-loop approach to calculate noise-induced transitions

Bistability permits the co-existence of two distinct cell fates in a population of genetically identical cells. Noise induced transitions between two fates of a bistable system are difficult to calculate due to the intricate interplay between nonlinear dynamics and noise in bistable positive feedback loops. Here we opened multivariable feedback loops at the slowest variable to obtain the open-loop function and the fluctuations in the open-loop output. By the subsequent reclosing of the loop, we calculated the mean first passage time (MFPT) using the Fokker-Planck equation in good agreement with the exact stochastic simulation. When an external component interacts with a feedback component, it amplifies the extrinsic noise in the loop. Consequently, the open-loop function is shifted and the transition rates between the two states in the closed loop are increased. Despite this shift, the open-loop output reflects the system faithfully to predict the MFPT in the feedback loop. Therefore, the open-loop approach can help theoretical analysis. Furthermore, the measurement of the mean value, variance, and the reaction time-scale of the open-loop output permits the prediction of MFPT simply from experimental data, which underscores the practical value of the stochastic open-loop approach.

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.

Roles of factorial noise in inducing bimodal gene expression

Some gene regulatory systems can exhibit bimodal distributions of mRNA or protein although the deterministic counterparts are monostable. This noise-induced bimodality is an interesting phenomenon and has important biological implications, but it is unclear how different sources of expression noise (each source creates so-called factorial noise that is defined as a component of the total noise) contribute separately to this stochastic bimodality. Here we consider a minimal model of gene regulation, which is monostable in the deterministic case. Although simple, this system contains factorial noise of two main kinds: promoter noise due to switching between gene states and transcriptional (or translational) noise due to synthesis and degradation of mRNA (or protein). To better trace the roles of factorial noise in inducing bimodality, we also analyze two limit models, continuous and adiabatic approximations, apart from the exact model. We show that in the case of slow gene switching, the continuous model where only promoter noise is considered can exhibit bimodality; in the case of fast switching, the adiabatic model where only transcriptional or translational noise is considered can also exhibit bimodality but the exact model cannot; and in other cases, both promoter noise and transcriptional or translational noise can cooperatively induce bimodality. Since slow gene switching and large protein copy numbers are characteristics of eukaryotic cells, whereas fast gene switching and small protein copy numbers are characteristics of prokaryotic cells, we infer that eukaryotic stochastic bimodality is induced mainly by promoter noise, whereas prokaryotic stochastic bimodality is induced primarily by transcriptional or translational noise.

Analytic approaches to stochastic gene expression in multicellular systems

Deterministic thermodynamic models of the complex systems, which control gene expression in metazoa, are helping researchers identify fundamental themes in the regulation of transcription. However, quantitative single cell studies are increasingly identifying regulatory mechanisms that control variability in expression. Such behaviors cannot be captured by deterministic models and are poorly suited to contemporary stochastic approaches that rely on continuum approximations, such as Langevin methods. Fortunately, theoretical advances in the modeling of transcription have assembled some general results that can be readily applied to systems being explored only through a deterministic approach. Here, I review some of the recent experimental evidence for the importance of genetically regulating stochastic effects during embryonic development and discuss key results from Markov theory that can be used to model this regulation. I then discuss several pairs of regulatory mechanisms recently investigated through a Markov approach. In each case, a deterministic treatment predicts no difference between the mechanisms, but the statistical treatment reveals the potential for substantially different distributions of transcriptional activity. In this light, features of gene regulation that seemed needlessly complex evolutionary baggage may be appreciated for their key contributions to reliability and precision of gene expression.

Role of intracellular stochasticity in biofilm growth. Insights from population balance modeling

There is increasing recognition that stochasticity involved in gene regulatory processes may help cells enhance the signal or synchronize expression for a group of genes. Thus the validity of the traditional deterministic approach to modeling the foregoing processes cannot be without exception. In this study, we identify a frequently encountered situation, i.e., the biofilm, which has in the past been persistently investigated with intracellular deterministic models in the literature. We show in this paper circumstances in which use of the intracellular deterministic model appears distinctly inappropriate. In Enterococcus faecalis, the horizontal gene transfer of plasmid spreads drug resistance. The induction of conjugation in planktonic and biofilm circumstances is examined here with stochastic as well as deterministic models. The stochastic model is formulated with the Chemical Master Equation (CME) for planktonic cells and Reaction-Diffusion Master Equation (RDME) for biofilm. The results show that although the deterministic model works well for the perfectly-mixed planktonic circumstance, it fails to predict the averaged behavior in the biofilm, a behavior that has come to be known as stochastic focusing. A notable finding from this work is that the interception of antagonistic feedback loops to signaling, accentuates stochastic focusing. Moreover, interestingly, increasing particle number of a control variable could lead to an even larger deviation. Intracellular stochasticity plays an important role in biofilm and we surmise by implications from the model, that cell populations may use it to minimize the influence from environmental fluctuation.

Stochastic effects as a force to increase the complexity of signaling networks

Cellular signaling networks are complex and appear to include many nonfunctional elements. Recently, it was suggested that nonfunctional interactions of proteins cause signaling noise, which, perhaps, shapes the signal transduction mechanism. However, the conditions under which molecular noise influences cellular information processing remain unclear. Here, we explore a large number of simple biological models of varying network sizes to understand the architectural conditions under which the interactions of signaling proteins can exhibit specific stochastic effects—called deviant effects—in which the average behavior of a biological system is substantially altered in the presence of molecular noise. We find that a small fraction of these networks does exhibit deviant effects and shares a common architectural feature whereas most of the networks show only insignificant levels of deviations. Interestingly, addition of seemingly unimportant interactions into protein networks gives rise to deviant effects.

Modeling heterogeneous responsiveness of intrinsic apoptosis pathway

BACKGROUND

Apoptosis is a cell suicide mechanism that enables multicellular organisms to maintain homeostasis and to eliminate individual cells that threaten the organism's survival. Dependent on the type of stimulus, apoptosis can be propagated by extrinsic pathway or intrinsic pathway. The comprehensive understanding of the molecular mechanism of apoptotic signaling allows for development of mathematical models, aiming to elucidate dynamical and systems properties of apoptotic signaling networks. There have been extensive efforts in modeling deterministic apoptosis network accounting for average behavior of a population of cells. Cellular networks, however, are inherently stochastic and significant cell-to-cell variability in apoptosis response has been observed at single cell level.

RESULTS

To address the inevitable randomness in the intrinsic apoptosis mechanism, we develop a theoretical and computational modeling framework of intrinsic apoptosis pathway at single-cell level, accounting for both deterministic and stochastic behavior. Our deterministic model, adapted from the well-accepted Fussenegger model, shows that an additional positive feedback between the executioner caspase and the initiator caspase plays a fundamental role in yielding the desired property of bistability. We then examine the impact of intrinsic fluctuations of biochemical reactions, viewed as intrinsic noise, and natural variation of protein concentrations, viewed as extrinsic noise, on behavior of the intrinsic apoptosis network. Histograms of the steady-state output at varying input levels show that the intrinsic noise could elicit a wider region of bistability over that of the deterministic model. However, the system stochasticity due to intrinsic fluctuations, such as the noise of steady-state response and the randomness of response delay, shows that the intrinsic noise in general is insufficient to produce significant cell-to-cell variations at physiologically relevant level of molecular numbers. Furthermore, the extrinsic noise represented by random variations of two key apoptotic proteins, namely Cytochrome C and inhibitor of apoptosis proteins (IAP), is modeled separately or in combination with intrinsic noise. The resultant stochasticity in the timing of intrinsic apoptosis response shows that the fluctuating protein variations can induce cell-to-cell stochastic variability at a quantitative level agreeing with experiments. Finally, simulations illustrate that the mean abundance of fluctuating IAP protein is positively correlated with the degree of cellular stochasticity of the intrinsic apoptosis pathway.

CONCLUSIONS

Our theoretical and computational study shows that the pronounced non-genetic heterogeneity in intrinsic apoptosis responses among individual cells plausibly arises from extrinsic rather than intrinsic origin of fluctuations. In addition, it predicts that the IAP protein could serve as a potential therapeutic target for suppression of the cell-to-cell variation in the intrinsic apoptosis responsiveness.

Looking for chemical reaction networks exhibiting a drift along a manifold of marginally stable states

I recently reported some examples of mass-action equations that have a continuous manifold of marginally stable stationary states [Brogioli, D., 2010. Marginally stable chemical systems as precursors of life. Phys. Rev. Lett. 105, 058102; Brogioli, D., 2011. Marginal stability in chemical systems and its relevance in the origin of life. Phys. Rev. E 84, 031931]. The corresponding chemical reaction networks show nonclassical effects, i.e. a violation of the mass-action equations, under the effect of the concentration fluctuations: the chemical system drifts along the marginally stable states. I proposed that this effect is potentially involved in abiogenesis. In the present paper, I analyze the mathematical properties of mass-action equations of marginally stable chemical reaction networks. The marginal stability implies that the mass-action equations obey some conservation law; I show that the mathematical properties of the conserved quantity characterize the motion along the marginally stable stationary state manifold, i.e. they allow to predict if the fluctuations give rise to a random walk or a drift under the effect of concentration fluctuations. Moreover, I show that the presence of the drift along the manifold of marginally stable stationary-states is a critical property, i.e. at least one of the reaction constants must be fine tuned in order to obtain the drift.

The Relevance of Topology in Parallel Simulation of Biological Networks

Important achievements in traditional biology has deepened the knowledge about living systems leading to an extensive identification of parts-list of the cell as well as of the interactions among biochemical species responsible for cell's regulation. Such an expanding knowledge also introduces new issues. For example the increasing comprehension of the inter- dependencies between pathways (pathways cross-talk) has resulted, on one hand, in the growth of informational complexity, on the other, in a strong lack of information coherence. The overall grand challenge remains unchanged: to be able to assemble the knowledge of every 'piece' of a system in order to figure out the behavior of the whole (integrative approach). In light of these considerations high performance computing plays a fundamental role in the context of in-silico biology. Stochastic simulation is a renowned analysis tool, which, although widely used, is subject to stringent computational requirements, in particular when dealing with heterogeneous and high dimensional systems. Here we introduce and discuss a methodology aimed at alleviating the burden of simulating complex biological networks. Such a method, which springs from graph theory, is based on the principle of fragmenting the computational space of a simulation trace and delegating the computation of fragments to a number of parallel processes.

Discreteness-induced concentration inversion in mesoscopic chemical systems

Molecular discreteness is apparent in small-volume chemical systems, such as biological cells, leading to stochastic kinetics. Here we present a theoretical framework to understand the effects of discreteness on the steady state of a monostable chemical reaction network. We consider independent realizations of the same chemical system in compartments of different volumes. Rate equations ignore molecular discreteness and predict the same average steady-state concentrations in all compartments. However, our theory predicts that the average steady state of the system varies with volume: if a species is more abundant than another for large volumes, then the reverse occurs for volumes below a critical value, leading to a concentration inversion effect. The addition of extrinsic noise increases the size of the critical volume. We theoretically predict the critical volumes and verify, by exact stochastic simulations, that rate equations are qualitatively incorrect in sub-critical volumes.

Stochastic steady state gain in a gene expression process with mRNA degradation control

Recent analyses with high-resolution single-molecule experimental methods have shown highly irregular and variable bursting of mRNA in a wide range of organisms. Noise in gene expression is thought to be beneficial in cell fate specifications, as it can lay a foundation for phenotypic diversification of isogenetic cells in the homogeneous environment. However, because the stability of proteins is, in many cases, higher than that of mRNAs, noise from transcriptional bursting can be considerably buffered at the protein level, limiting the effect of noisy mRNAs at a more global regulation level. This raises a question as to what constructive role noisy mRNAs can play in the system-level dynamics. In this study, we have addressed this question using the computational models that extend the conventional transcriptional bursting model with a post-transcriptional regulation step. Surprisingly, by comparing this stochastic model with the corresponding deterministic model, we find that intrinsic fluctuations can substantially increase the expected mRNA level. Because effects of a higher mRNA level can be transmitted to the protein level even with slow protein degradation rates, this finding suggests that an increase in the protein level is another potential effect of transcriptional bursting. Here, we show that this striking steady state increase is caused by the asynchronous nature of molecular reactions, which allows the transcriptional regulation model to create additional modes of qualitatively distinct dynamics. Our results illustrate non-intuitive effects of reaction asynchronicity on system dynamics that cannot be captured by the traditional deterministic framework. Because molecular reactions are intrinsically stochastic and asynchronous, these findings may have broad implications in modelling and understanding complex biological systems.

Deterministic Function Computation with Chemical Reaction Networks

Chemical reaction networks (CRNs) formally model chemistry in a well-mixed solution. CRNs are widely used to describe information processing occurring in natural cellular regulatory networks, and with upcoming advances in synthetic biology, CRNs are a promising language for the design of artificial molecular control circuitry. Nonetheless, despite the widespread use of CRNs in the natural sciences, the range of computational behaviors exhibited by CRNs is not well understood. CRNs have been shown to be efficiently Turing-universal (i.e., able to simulate arbitrary algorithms) when allowing for a small probability of error. CRNs that are guaranteed to converge on a correct answer, on the other hand, have been shown to decide only the semilinear predicates (a multi-dimensional generalization of "eventually periodic" sets). We introduce the notion of function, rather than predicate, computation by representing the output of a function f : ℕ k → ℕ l by a count of some molecular species, i.e., if the CRN starts with x1, …, xk molecules of some "input" species X1, …, Xk , the CRN is guaranteed to converge to having f(x1, …, xk ) molecules of the "output" species Y1, …, Yl . We show that a function f : ℕ k → ℕ l is deterministically computed by a CRN if and only if its graph {(x, y) ∈ ℕ k × ℕ l ∣ f(x) = y} is a semilinear set. Finally, we show that each semilinear function f (a function whose graph is a semilinear set) can be computed by a CRN on input x in expected time O(polylog ∥x∥1).

Intrinsic noise alters the frequency spectrum of mesoscopic oscillatory chemical reaction systems

Mesoscopic oscillatory reaction systems, for example in cell biology, can exhibit stochastic oscillations in the form of cyclic random walks even if the corresponding macroscopic system does not oscillate. We study how the intrinsic noise from molecular discreteness influences the frequency spectrum of mesoscopic oscillators using as a model system a cascade of coupled Brusselators away from the Hopf bifurcation. The results show that the spectrum of an oscillator depends on the level of noise. In particular, the peak frequency of the oscillator is reduced by increasing noise, and the bandwidth increased. Along a cascade of coupled oscillators, the peak frequency is further reduced with every stage and also the bandwidth is reduced. These effects can help understand the role of noise in chemical oscillators and provide fingerprints for more reliable parameter identification and volume measurement from experimental spectra.

Construction and accuracy of partial differential equation approximations to the chemical master equation

The mesoscopic description of chemical kinetics, the chemical master equation, can be exactly solved in only a few simple cases. The analytical intractability stems from the discrete character of the equation, and hence considerable effort has been invested in the development of Fokker-Planck equations, second-order partial differential equation approximations to the master equation. We here consider two different types of higher-order partial differential approximations, one derived from the system-size expansion and the other from the Kramers-Moyal expansion, and derive the accuracy of their predictions for chemical reactive networks composed of arbitrary numbers of unimolecular and bimolecular reactions. In particular, we show that the partial differential equation approximation of order Q from the Kramers-Moyal expansion leads to estimates of the mean number of molecules accurate to order Ω(-(2Q-3)/2), of the variance of the fluctuations in the number of molecules accurate to order Ω(-(2Q-5)/2), and of skewness accurate to order Ω(-(Q-2)). We also show that for large Q, the accuracy in the estimates can be matched only by a partial differential equation approximation from the system-size expansion of approximate order 2Q. Hence, we conclude that partial differential approximations based on the Kramers-Moyal expansion generally lead to considerably more accurate estimates in the mean, variance, and skewness than approximations of the same order derived from the system-size expansion.

Computation of steady-state probability distributions in stochastic models of cellular networks

Cellular processes are "noisy". In each cell, concentrations of molecules are subject to random fluctuations due to the small numbers of these molecules and to environmental perturbations. While noise varies with time, it is often measured at steady state, for example by flow cytometry. When interrogating aspects of a cellular network by such steady-state measurements of network components, a key need is to develop efficient methods to simulate and compute these distributions. We describe innovations in stochastic modeling coupled with approaches to this computational challenge: first, an approach to modeling intrinsic noise via solution of the chemical master equation, and second, a convolution technique to account for contributions of extrinsic noise. We show how these techniques can be combined in a streamlined procedure for evaluation of different sources of variability in a biochemical network. Evaluation and illustrations are given in analysis of two well-characterized synthetic gene circuits, as well as a signaling network underlying the mammalian cell cycle entry.

Marginal stability in chemical systems and its relevance in the origin of life

Concentration fluctuations are always present in solutions; it has been noticed that, in chemical systems, they can lead to deviations from what is expected from mass-action equations. I recently described the class of the "marginally stable" chemical systems; namely, a system that have an infinity of stationary states forming a continuous curve, and I showed that they present such deviations, which appear as a drift along the stationary-state curve [Phys. Rev. Lett. 105, 058102 (2010)]. Here I describe various marginally stable chemical reaction networks, including replicating molecules, and I present numerical calculations based on reaction-diffusion master equations, showing that the thermodynamic fluctuations induce a drift. This drift can be interpreted in terms of evolution toward a more efficiently replicating system and is analogous to a Darwinian evolution. The concentration fluctuations observed during the drift are scale invariant. Relevance of this phenomenon to the origin of life is discussed. I propose that marginal stability is the mathematical property defining chemical reaction networks potentially involved in the origin of life.

Connection between noise-induced symmetry breaking and an information-decoding function for intracellular networks

The biological function of noise-induced symmetry breaking (NISB) is still unclear even though it may potentially occur in noisy intracellular systems. In this work, I demonstrate that information decoding from a noisy signal is a potential biological function of NISB by revealing that NISB naturally emerges from an optimal information-decoding dynamics and that several intracellular networks can be identified with the information-decoding dynamics. I also propose a mean first passage time profile as a way to experimentally identify NISB.

Transcriptional regulation: effects of promoter proximal pausing on speed, synchrony and reliability

Recent whole genome polymerase binding assays in the Drosophila embryo have shown that a substantial proportion of uninduced genes have pre-assembled RNA polymerase-II transcription initiation complex (PIC) bound to their promoters. These constitute a subset of promoter proximally paused genes for which mRNA elongation instead of promoter access is regulated. This difference can be described as a rearrangement of the regulatory topology to control the downstream transcriptional process of elongation rather than the upstream transcriptional initiation event. It has been shown experimentally that genes with the former mode of regulation tend to induce faster and more synchronously, and that promoter-proximal pausing is observed mainly in metazoans, in accord with a posited impact on synchrony. However, it has not been shown whether or not it is the change in the regulated step per se that is causal. We investigate this question by proposing and analyzing a continuous-time Markov chain model of PIC assembly regulated at one of two steps: initial polymerase association with DNA, or release from a paused, transcribing state. Our analysis demonstrates that, over a wide range of physical parameters, increased speed and synchrony are functional consequences of elongation control. Further, we make new predictions about the effect of elongation regulation on the consistent control of total transcript number between cells. We also identify which elements in the transcription induction pathway are most sensitive to molecular noise and thus possibly the most evolutionarily constrained. Our methods produce symbolic expressions for quantities of interest with reasonable computational effort and they can be used to explore the interplay between interaction topology and molecular noise in a broader class of biochemical networks. We provide general-purpose code implementing these methods.

Gene activation is an inherently random process because numerous diffusing proteins and DNA must first interact by random association before transcription can begin. For many genes the necessary protein–DNA associations only begin after activation, but it has recently been noted that a large class of genes in multicellular organisms can assemble the initiation complex of proteins on the core promoter prior to activation. For these genes, activation merely releases polymerase from the preassembled complex to transcribe the gene. It has been proposed on the basis of experiments that such a mechanism, while possibly costly, increases both the speed and the synchrony of the process of gene transcription. We study a realistic model of gene transcription, and show that this conclusion holds for all but a tiny fraction of the space of physical rate parameters that govern the process. The improved control of cell-to-cell variations afforded by regulation through a paused polymerase may help multicellular organisms achieve the high degree of coordination required for development. Our approach has also generated tools with which one can study the effects of analogous changes in other molecular networks and determine the relative importance of various molecular binding rates to particular system properties.

The role of incoherent microRNA-mediated feedforward loops in noise buffering

MicroRNAs are endogenous non-coding RNAs which negatively regulate the expression of protein-coding genes in plants and animals. They are known to play an important role in several biological processes and, together with transcription factors, form a complex and highly interconnected regulatory network. Looking at the structure of this network, it is possible to recognize a few overrepresented motifs which are expected to perform important elementary regulatory functions. Among them, a special role is played by the microRNA-mediated feedforward loop in which a master transcription factor regulates a microRNA and, together with it, a set of target genes. In this paper we show analytically and through simulations that the incoherent version of this motif can couple the fine-tuning of a target protein level with an efficient noise control, thus conferring precision and stability to the overall gene expression program, especially in the presence of fluctuations in upstream regulators. Among the other results, a nontrivial prediction of our model is that the optimal attenuation of fluctuations coincides with a modest repression of the target expression. This feature is coherent with the expected fine-tuning function and in agreement with experimental observations of the actual impact of a wide class of microRNAs on the protein output of their targets. Finally, we describe the impact on noise-buffering efficiency of the cross-talk between microRNA targets that can naturally arise if the microRNA-mediated circuit is not considered as isolated, but embedded in a larger network of regulations.

Automatic design of digital synthetic gene circuits

De novo computational design of synthetic gene circuits that achieve well-defined target functions is a hard task. Existing, brute-force approaches run optimization algorithms on the structure and on the kinetic parameter values of the network. However, more direct rational methods for automatic circuit design are lacking. Focusing on digital synthetic gene circuits, we developed a methodology and a corresponding tool for in silico automatic design. For a given truth table that specifies a circuit's input–output relations, our algorithm generates and ranks several possible circuit schemes without the need for any optimization. Logic behavior is reproduced by the action of regulatory factors and chemicals on the promoters and on the ribosome binding sites of biological Boolean gates. Simulations of circuits with up to four inputs show a faithful and unequivocal truth table representation, even under parametric perturbations and stochastic noise. A comparison with already implemented circuits, in addition, reveals the potential for simpler designs with the same function. Therefore, we expect the method to help both in devising new circuits and in simplifying existing solutions.

Synthetic Biology is a novel discipline that aims at the construction of new biological systems able to perform specific tasks. Following the example of electrical engineering, most of the synthetic systems so far realized look like circuits where smaller DNA-encoded components are interconnected by the exchange of different kinds of molecules. According to this modular approach, we developed, in a previous work, a tool for the visual design of new genetic circuits whose components are displayed on the computer screen and connected through hypothetical wires where molecules flow. Here, we present an extension of this tool that automatically computes the structure of a digital gene circuit–where the inputs and the output take only 0/1 values–by applying procedures commonly used in electrical engineering to biology. In this way, our method generalizes and simplifies the design of genetic circuits far more complex than the ones so far realized. Moreover, different from other currently used methods, our approach limits the use of optimization procedures and drastically reduces the computational time necessary to derive the circuit structure. Future improvements can be achieved by exploiting some more biological mechanisms able to mimic Boolean behavior, without a substantial growth of the algorithmic complexity.

Non-equilibrium phase transition in mesoscopic biochemical systems: from stochastic to nonlinear dynamics and beyond

A theory for an non-equilibrium phase transition in a driven biochemical network is presented. The theory is based on the chemical master equation (CME) formulation of mesoscopic biochemical reactions and the mathematical method of large deviations. The large deviations theory provides an analytical tool connecting the macroscopic multi-stability of an open chemical system with the multi-scale dynamics of its mesoscopic counterpart. It shows a corresponding non-equilibrium phase transition among multiple stochastic attractors. As an example, in the canonical phosphorylation-dephosphorylation system with feedback that exhibits bistability, we show that the non-equilibrium steady-state (NESS) phase transition has all the characteristics of classic equilibrium phase transition: Maxwell construction, a discontinuous first-derivative of the 'free energy function', Lee-Yang's zero for a generating function and a critical point that matches the cusp in nonlinear bifurcation theory. To the biochemical system, the mathematical analysis suggests three distinct timescales and needed levels of description. They are (i) molecular signalling, (ii) biochemical network nonlinear dynamics, and (iii) cellular evolution. For finite mesoscopic systems such as a cell, motions associated with (i) and (iii) are stochastic while that with (ii) is deterministic. Both (ii) and (iii) are emergent properties of a dynamic biochemical network.

Discrete diffusion models to study the effects of Mg2+ concentration on the PhoPQ signal transduction system

Background

The challenge today is to develop a modeling and simulation paradigm that integrates structural, molecular and genetic data for a quantitative understanding of physiology and behavior of biological processes at multiple scales. This modeling method requires techniques that maintain a reasonable accuracy of the biological process and also reduces the computational overhead. This objective motivates the use of new methods that can transform the problem from energy and affinity based modeling to information theory based modeling. To achieve this, we transform all dynamics within the cell into a random event time, which is specified through an information domain measure like probability distribution. This allows us to use the “in silico” stochastic event based modeling approach to find the molecular dynamics of the system.

Results

In this paper, we present the discrete event simulation concept using the example of the signal transduction cascade triggered by extra-cellular Mg²⁺ concentration in the two component PhoPQ regulatory system of Salmonella Typhimurium. We also present a model to compute the information domain measure of the molecular transport process by estimating the statistical parameters of inter-arrival time between molecules/ions coming to a cell receptor as external signal. This model transforms the diffusion process into the information theory measure of stochastic event completion time to get the distribution of the Mg²⁺ departure events. Using these molecular transport models, we next study the in-silico effects of this external trigger on the PhoPQ system.

Conclusions

Our results illustrate the accuracy of the proposed diffusion models in explaining the molecular/ionic transport processes inside the cell. Also, the proposed simulation framework can incorporate the stochasticity in cellular environments to a certain degree of accuracy. We expect that this scalable simulation platform will be able to model more complex biological systems with reasonable accuracy to understand their temporal dynamics.

HIV promoter integration site primarily modulates transcriptional burst size rather than frequency

Mammalian gene expression patterns, and their variability across populations of cells, are regulated by factors specific to each gene in concert with its surrounding cellular and genomic environment. Lentiviruses such as HIV integrate their genomes into semi-random genomic locations in the cells they infect, and the resulting viral gene expression provides a natural system to dissect the contributions of genomic environment to transcriptional regulation. Previously, we showed that expression heterogeneity and its modulation by specific host factors at HIV integration sites are key determinants of infected-cell fate and a possible source of latent infections. Here, we assess the integration context dependence of expression heterogeneity from diverse single integrations of a HIV-promoter/GFP-reporter cassette in Jurkat T-cells. Systematically fitting a stochastic model of gene expression to our data reveals an underlying transcriptional dynamic, by which multiple transcripts are produced during short, infrequent bursts, that quantitatively accounts for the wide, highly skewed protein expression distributions observed in each of our clonal cell populations. Interestingly, we find that the size of transcriptional bursts is the primary systematic covariate over integration sites, varying from a few to tens of transcripts across integration sites, and correlating well with mean expression. In contrast, burst frequencies are scattered about a typical value of several per cell-division time and demonstrate little correlation with the clonal means. This pattern of modulation generates consistently noisy distributions over the sampled integration positions, with large expression variability relative to the mean maintained even for the most productive integrations, and could contribute to specifying heterogeneous, integration-site-dependent viral production patterns in HIV-infected cells. Genomic environment thus emerges as a significant control parameter for gene expression variation that may contribute to structuring mammalian genomes, as well as be exploited for survival by integrating viruses.

Marginally stable chemical systems as precursors of life

Current research on the origin of life aims at finding the simplest entity that can undergo spontaneous Darwinian evolution toward increasing replication efficiency. Here I consider some of the models of self-replicating molecular systems, and I show that they exhibit a distinct feature, namely, an infinity of stationary states forming a continuous curve; i.e., they are only marginally stable. I show that, in marginally stable chemical systems, thermodynamic fluctuations induce a drift directed toward increasing replication efficiency. This drift represents a form of evolution, taking place slowly, cooperatively, in macroscopic volumes of water.

Counter-intuitive stochastic behavior of simple gene circuits with negative feedback

It has often been taken for granted that negative feedback loops in gene regulation work as homeostatic control mechanisms. If one increases the regulation strength a less noisy signal is to be expected. However, recent theoretical studies have reported the exact contrary, counter-intuitive observation, which has left a question mark over the relationship between negative feedback loops and noise. We explore and systematically analyze several minimal models of gene regulation, where a transcriptional repressor negatively regulates its own expression. For models including a quasi-steady-state assumption, we identify processes that buffer noise change (RNA polymerase binding) or accentuate it (repressor dimerization) alongside increasing feedback strength. Moreover, we show that lumping together transcription and translation in simplified models clearly underestimates the impact of negative feedback strength on the system's noise. In contrast, in systems without a quasi-steady-state assumption, noise always increases with negative feedback strength. Hence, subtle mathematical properties and model assumptions yield different types of noise profiles and, by consequence, previous studies have simultaneously reported decrease, increase or persistence of noise levels with increasing feedback. We discuss our findings in terms of separation of timescales and time correlations between molecular species distributions, extending current theoretical findings on the topic and allowing us to propose what we believe new ways to better characterize noise.

Equilibriumlike behavior in chemical reaction networks far from equilibrium

In an equilibrium chemical reaction mixture, the number of molecules present obeys a Poisson distribution. We report that, surprisingly, the same is true of a large class of nonequilibrium reaction networks. In particular, we show that certain topological features imply a Poisson distribution, whatever the reaction rates. Such driven systems also obey an analog of the fluctuation-dissipation theorem. Our results shed light on the fundamental question of when equilibrium concepts might apply to nonequilibrium systems and may have applications to models of noise in biochemical networks.

Temperature control of fimbriation circuit switch in uropathogenic Escherichia coli: quantitative analysis via automated model abstraction

Uropathogenic Escherichia coli (UPEC) represent the predominant cause of urinary tract infections (UTIs). A key UPEC molecular virulence mechanism is type 1 fimbriae, whose expression is controlled by the orientation of an invertible chromosomal DNA element-the fim switch. Temperature has been shown to act as a major regulator of fim switching behavior and is overall an important indicator as well as functional feature of many urologic diseases, including UPEC host-pathogen interaction dynamics. Given this panoptic physiological role of temperature during UTI progression and notable empirical challenges to its direct in vivo studies, in silico modeling of corresponding biochemical and biophysical mechanisms essential to UPEC pathogenicity may significantly aid our understanding of the underlying disease processes. However, rigorous computational analysis of biological systems, such as fim switch temperature control circuit, has hereto presented a notoriously demanding problem due to both the substantial complexity of the gene regulatory networks involved as well as their often characteristically discrete and stochastic dynamics. To address these issues, we have developed an approach that enables automated multiscale abstraction of biological system descriptions based on reaction kinetics. Implemented as a computational tool, this method has allowed us to efficiently analyze the modular organization and behavior of the E. coli fimbriation switch circuit at different temperature settings, thus facilitating new insights into this mode of UPEC molecular virulence regulation. In particular, our results suggest that, with respect to its role in shutting down fimbriae expression, the primary function of FimB recombinase may be to effect a controlled down-regulation (rather than increase) of the ON-to-OFF fim switching rate via temperature-dependent suppression of competing dynamics mediated by recombinase FimE. Our computational analysis further implies that this down-regulation mechanism could be particularly significant inside the host environment, thus potentially contributing further understanding toward the development of novel therapeutic approaches to UPEC-caused UTIs.

Detailed simulations of cell biology with Smoldyn 2.1

Most cellular processes depend on intracellular locations and random collisions of individual protein molecules. To model these processes, we developed algorithms to simulate the diffusion, membrane interactions, and reactions of individual molecules, and implemented these in the Smoldyn program. Compared to the popular MCell and ChemCell simulators, we found that Smoldyn was in many cases more accurate, more computationally efficient, and easier to use. Using Smoldyn, we modeled pheromone response system signaling among yeast cells of opposite mating type. This model showed that secreted Bar1 protease might help a cell identify the fittest mating partner by sharpening the pheromone concentration gradient. This model involved about 200,000 protein molecules, about 7000 cubic microns of volume, and about 75 minutes of simulated time; it took about 10 hours to run. Over the next several years, as faster computers become available, Smoldyn will allow researchers to model and explore systems the size of entire bacterial and smaller eukaryotic cells.