Modeling stochasticity in gene regulation: characterization in the terms of the underlying distribution function

Statistics or biology: the zero-inflation controversy about scRNA-seq data

Researchers view vast zeros in single-cell RNA-seq data differently: some regard zeros as biological signals representing no or low gene expression, while others regard zeros as missing data to be corrected. To help address the controversy, here we discuss the sources of biological and non-biological zeros; introduce five mechanisms of adding non-biological zeros in computational benchmarking; evaluate the impacts of non-biological zeros on data analysis; benchmark three input data types: observed counts, imputed counts, and binarized counts; discuss the open questions regarding non-biological zeros; and advocate the importance of transparent analysis.

Queuing Models of Gene Expression: Analytical Distributions and Beyond

Activation of a gene is a multistep biochemical process, involving recruitments of transcription factors and histone kinases as well as modification of histones. Many of these intermediate reaction steps would have been unspecified by experiments. Therefore, classical two-state models of gene expression established based on the memoryless (or Markovian) assumption would not well describe the reality in gene expression. Recent experimental data have indicated that the inactive phases of gene promoters are differently distributed, showing strong memory. Here, we use a nonexponential waiting-time distribution to model the complex activation process of a gene, and then analyze a queuing model of stochastic transcription. We successfully derive the analytical expression of the stationary mRNA distribution, which provides insight into the effect of molecular memory created by complex activating events on the mRNA expression. We find that the reduction in the waiting-time noise may result in the increase in the mRNA noise, contrary to the previous conclusion. Based on the derived distribution, we also develop a method to infer the waiting-time distribution from a known mRNA distribution. Data analysis on a realistic example verifies the validity of this method.

Gene-Specific Linear Trends Constrain Transcriptional Variability of the Toll-like Receptor Signaling

Single-cell gene expression is inherently variable, but how this variability is controlled in response to stimulation remains unclear. Here, we use single-cell RNA-seq and single-molecule mRNA counting (smFISH) to study inducible gene expression in the immune toll-like receptor system. We show that mRNA counts of tumor necrosis factor α conform to a standard stochastic switch model, while transcription of interleukin-1β involves an additional regulatory step resulting in increased heterogeneity. Despite different modes of regulation, systematic analysis of single-cell data for a range of genes demonstrates that the variability in transcript count is linearly constrained by the mean response over a range of conditions. Mathematical modeling of smFISH counts and experimental perturbation of chromatin state demonstrates that linear constraints emerge through modulation of transcriptional bursting along with gene-specific relationships. Overall, our analyses demonstrate that the variability of the inducible single-cell mRNA response is constrained by transcriptional bursting.

•

Single-cell TNF-α and IL-1β mRNA responses are differentially controlled

•

Variability of TLR-induced responses scale linearly with mean mRNA counts

•

Gene-specific constraints emerge via modulation of transcriptional bursting

•

Chromatin state regulates transcriptional bursting of IL-1β

Hybrid deterministic/stochastic simulation of complex biochemical systems

In a biological cell, cellular functions and the genetic regulatory apparatus are implemented and controlled by complex networks of chemical reactions involving genes, proteins, and enzymes. Accurate computational models are indispensable means for understanding the mechanisms behind the evolution of a complex system, not always explored with wet lab experiments. To serve their purpose, computational models, however, should be able to describe and simulate the complexity of a biological system in many of its aspects. Moreover, it should be implemented by efficient algorithms requiring the shortest possible execution time, to avoid enlarging excessively the time elapsing between data analysis and any subsequent experiment. Besides the features of their topological structure, the complexity of biological networks also refers to their dynamics, that is often non-linear and stiff. The stiffness is due to the presence of molecular species whose abundance fluctuates by many orders of magnitude. A fully stochastic simulation of a stiff system is computationally time-expensive. On the other hand, continuous models are less costly, but they fail to capture the stochastic behaviour of small populations of molecular species. We introduce a new efficient hybrid stochastic-deterministic computational model and the software tool MoBioS (MOlecular Biology Simulator) implementing it. The mathematical model of MoBioS uses continuous differential equations to describe the deterministic reactions and a Gillespie-like algorithm to describe the stochastic ones. Unlike the majority of current hybrid methods, the MoBioS algorithm divides the reactions' set into fast reactions, moderate reactions, and slow reactions and implements a hysteresis switching between the stochastic model and the deterministic model. Fast reactions are approximated as continuous-deterministic processes and modelled by deterministic rate equations. Moderate reactions are those whose reaction waiting time is greater than the fast reaction waiting time but smaller than the slow reaction waiting time. A moderate reaction is approximated as a stochastic (deterministic) process if it was classified as a stochastic (deterministic) process at the time at which it crosses the threshold of low (high) waiting time. A Gillespie First Reaction Method is implemented to select and execute the slow reactions. The performances of MoBios were tested on a typical example of hybrid dynamics: that is the DNA transcription regulation. The simulated dynamic profile of the reagents' abundance and the estimate of the error introduced by the fully deterministic approach were used to evaluate the consistency of the computational model and that of the software tool.

Stochastic models of gene transcription with upstream drives: exact solution and sample path characterization

Gene transcription is a highly stochastic and dynamic process. As a result, the mRNA copy number of a given gene is heterogeneous both between cells and across time. We present a framework to model gene transcription in populations of cells with time-varying (stochastic or deterministic) transcription and degradation rates. Such rates can be understood as upstream cellular drives representing the effect of different aspects of the cellular environment. We show that the full solution of the master equation contains two components: a model-specific, upstream effective drive, which encapsulates the effect of cellular drives (e.g. entrainment, periodicity or promoter randomness) and a downstream transcriptional Poissonian part, which is common to all models. Our analytical framework treats cell-to-cell and dynamic variability consistently, unifying several approaches in the literature. We apply the obtained solution to characterize different models of experimental relevance, and to explain the influence on gene transcription of synchrony, stationarity, ergodicity, as well as the effect of time scales and other dynamic characteristics of drives. We also show how the solution can be applied to the analysis of noise sources in single-cell data, and to reduce the computational cost of stochastic simulations.

Modeling dynamics of cell-to-cell variability in TRAIL-induced apoptosis explains fractional killing and predicts reversible resistance

Isogenic cells sensing identical external signals can take markedly different decisions. Such decisions often correlate with pre-existing cell-to-cell differences in protein levels. When not neglected in signal transduction models, these differences are accounted for in a static manner, by assuming randomly distributed initial protein levels. However, this approach ignores the a priori non-trivial interplay between signal transduction and the source of this cell-to-cell variability: temporal fluctuations of protein levels in individual cells, driven by noisy synthesis and degradation. Thus, modeling protein fluctuations, rather than their consequences on the initial population heterogeneity, would set the quantitative analysis of signal transduction on firmer grounds. Adopting this dynamical view on cell-to-cell differences amounts to recast extrinsic variability into intrinsic noise. Here, we propose a generic approach to merge, in a systematic and principled manner, signal transduction models with stochastic protein turnover models. When applied to an established kinetic model of TRAIL-induced apoptosis, our approach markedly increased model prediction capabilities. One obtains a mechanistic explanation of yet-unexplained observations on fractional killing and non-trivial robust predictions of the temporal evolution of cell resistance to TRAIL in HeLa cells. Our results provide an alternative explanation to survival via induction of survival pathways since no TRAIL-induced regulations are needed and suggest that short-lived anti-apoptotic protein Mcl1 exhibit large and rare fluctuations. More generally, our results highlight the importance of accounting for stochastic protein turnover to quantitatively understand signal transduction over extended durations, and imply that fluctuations of short-lived proteins deserve particular attention.

Decomposing noise in biochemical signaling systems highlights the role of protein degradation

Stochasticity is an essential aspect of biochemical processes at the cellular level. We now know that living cells take advantage of stochasticity in some cases and counteract stochastic effects in others. Here we propose a method that allows us to calculate contributions of individual reactions to the total variability of a system's output. We demonstrate that reactions differ significantly in their relative impact on the total noise and we illustrate the importance of protein degradation on the overall variability for a range of molecular processes and signaling systems. With our flexible and generally applicable noise decomposition method, we are able to shed new, to our knowledge, light on the sources and propagation of noise in biochemical reaction networks; in particular, we are able to show how regulated protein degradation can be employed to reduce the noise in biochemical systems.

Self-regulatory gene: an exact solution for the gene gate model

The stochastic dynamics of gene expression is often described by highly abstract models involving only the key molecular actors DNA, RNA, and protein, neglecting all further details of the transcription and translation processes. One example of such models is the "gene gate model," which contains a minimal set of actors and kinetic parameters, which allows us to describe the regulation of a gene by both repression and activation. Based on this approach, we formulate a master equation for the case of a single gene regulated by its own product-a transcription factor-and solve it exactly. The obtained gene product distributions display features of mono- and bimodality, depending on the choice of parameters. We discuss our model in the perspective of other models in the literature.

Gene expression in self-repressing system with multiple gene copies

We analyze a simple model of a self-repressing system with multiple gene copies. Protein molecules may bound to DNA promoters and block their own transcription. We derive analytical expressions for the variance of the number of protein molecules in the stationary state in the self-consistent mean-field approximation. We show that the Fano factor (the variance divided by the mean value) is bigger for the one-gene case than for two gene copies and the difference decreases to zero as frequencies of binding and unbinding increase to infinity.

Multiscale stochastic modelling of gene expression

Stochastic phenomena in gene regulatory networks can be modelled by the chemical master equation for gene products such as mRNA and proteins. If some of these elements are present in significantly higher amounts than the rest, or if some of the reactions between these elements are substantially faster than others, it is often possible to reduce the master equation to a simpler problem using asymptotic methods. We present examples of such a procedure and analyse the relationship between the reduced models and the original.

Exact and approximate distributions of protein and mRNA levels in the low-copy regime of gene expression

Gene expression at the single-cell level incorporates reaction mechanisms which are intrinsically stochastic as they involve molecular species present at low copy numbers. The dynamics of these mechanisms can be described quantitatively using stochastic master-equation modelling; in this paper we study a generic gene-expression model of this kind which explicitly includes the representations of the processes of transcription and translation. For this model we determine the generating function of the steady-state distribution of mRNA and protein counts and characterise the underlying probability law using a combination of analytic, asymptotic and numerical approaches, finding that the distribution may assume a number of qualitatively distinct forms. The results of the analysis are suitable for comparison with single-molecule resolution gene-expression data emerging from recent experimental studies.

On the spontaneous stochastic dynamics of a single gene: complexity of the molecular interplay at the promoter

BACKGROUND

Gene promoters can be in various epigenetic states and undergo interactions with many molecules in a highly transient, probabilistic and combinatorial way, resulting in a complex global dynamics as observed experimentally. However, models of stochastic gene expression commonly consider promoter activity as a two-state on/off system. We consider here a model of single-gene stochastic expression that can represent arbitrary prokaryotic or eukaryotic promoters, based on the combinatorial interplay between molecules and epigenetic factors, including energy-dependent remodeling and enzymatic activities.

RESULTS

We show that, considering the mere molecular interplay at the promoter, a single-gene can demonstrate an elaborate spontaneous stochastic activity (eg. multi-periodic multi-relaxation dynamics), similar to what is known to occur at the gene-network level. Characterizing this generic model with indicators of dynamic and steady-state properties (including power spectra and distributions), we reveal the potential activity of any promoter and its influence on gene expression. In particular, we can reproduce, based on biologically relevant mechanisms, the strongly periodic patterns of promoter occupancy by transcription factors (TF) and chromatin remodeling as observed experimentally on eukaryotic promoters. Moreover, we link several of its characteristics to properties of the underlying biochemical system. The model can also be used to identify behaviors of interest (eg. stochasticity induced by high TF concentration) on minimal systems and to test their relevance in larger and more realistic systems. We finally show that TF concentrations can regulate many aspects of the stochastic activity with a considerable flexibility and complexity.

CONCLUSIONS

This tight promoter-mediated control of stochasticity may constitute a powerful asset for the cell. Remarkably, a strongly periodic activity that demonstrates a complex TF concentration-dependent control is obtained when molecular interactions have typical characteristics observed on eukaryotic promoters (high mobility, functional redundancy, many alternate states/pathways). We also show that this regime results in a direct and indirect energetic cost. Finally, this model can constitute a framework for unifying various experimental approaches. Collectively, our results show that a gene - the basic building block of complex regulatory networks - can itself demonstrate a significantly complex behavior.

Mathematical modeling as a tool for planning anticancer therapy

We review a large volume of literature concerning mathematical models of cancer therapy, oriented towards optimization of treatment protocols. The review, although partly idiosyncratic, covers such major areas of therapy optimization as phase-specific chemotherapy, antiangiogenic therapy and therapy under drug resistance. We start from early cell cycle progression models, very simple but admitting explicit mathematical solutions, based on methods of control theory. We continue with more complex models involving evolution of drug resistance and pharmacokinetic and pharmacodynamic effects. Then, we consider two more recent areas: angiogenesis of tumors and molecular signaling within and among cells. We discuss biological background and mathematical techniques of this field, which has a large although only partly realized potential for contributing to cancer treatment.

Investigating the two-moment characterisation of subcellular biochemical networks

While ordinary differential equations (ODEs) form the conceptual framework for modelling many cellular processes, specific situations demand stochastic models to capture the influence of noise. The most common formulation of stochastic models for biochemical networks is the chemical master equation (CME). While stochastic simulations are a practical way to realise the CME, analytical approximations offer more insight into the influence of noise. Towards that end, the two-moment approximation (2MA) is a promising addition to the established analytical approaches including the chemical Langevin equation (CLE) and the related linear noise approximation (LNA). The 2MA approach directly tracks the mean and (co)variance which are coupled in general. This coupling is not obvious in CME and CLE and ignored by LNA and conventional ODE models. We extend previous derivations of 2MA by allowing (a) non-elementary reactions and (b) relative concentrations. Often, several elementary reactions are approximated by a single step. Furthermore, practical situations often require the use of relative concentrations. We investigate the applicability of the 2MA approach to the well-established fission yeast cell cycle model. Our analytical model reproduces the clustering of cycle times observed in experiments. This is explained through multiple resettings of M-phase promoting factor (MPF), caused by the coupling between mean and (co)variance, near the G2/M transition.

Translational repression contributes greater noise to gene expression than transcriptional repression

Stochastic effects in gene expression may result in different physiological states of individual cells, with consequences for pathogen survival and artificial gene network design. We studied the contributions of a regulatory factor to gene expression noise in four basic mechanisms of negative gene expression control: 1), transcriptional regulation by a protein repressor, 2), translational repression by a protein; 3), transcriptional repression by RNA; and 4), RNA interference with the translation. We investigated a general model of a two-gene network, using the chemical master equation and a moment generating function approach. We compared the expression noise of genes with the same effective transcription and translation initiation rates resulting from the action of different repressors, whereas previous studies compared the noise of genes with the same mean expression level but different initiation rates. Our results show that translational repression results in a higher noise than repression on the promoter level, and that this relationship does not depend on quantitative parameter values. We also show that regulation of protein degradation contributes more noise than regulated degradation of mRNA. These are unexpected results, because previous investigations suggested that translational regulation is more accurate. The relative magnitude of the noise introduced by protein and RNA repressors depends on the protein and mRNA degradation rates, and we derived expressions for the threshold below which the noise introduced by a protein repressor is higher than the noise introduced by an RNA repressor.

Oscillations and bistability in the stochastic model of p53 regulation

The p53 regulatory pathway controls cell responses, which include cell cycle arrest, DNA repair, apoptosis and cellular senescence. We propose a stochastic model of p53 regulation, which is based on two feedback loops: the negative, coupling p53 with its immediate downregulator Mdm2, and the positive, which involves PTEN, PIP3 and Akt. Existence of the negative feedback assures homeostasis of healthy cells and oscillatory responses of DNA-damaged cells, which are persistent when DNA repair is inefficient and the positive feedback loop is broken. The positive feedback destroys the negative coupling between Mdm2 and p53 by sequestering most of Mdm2 in cytoplasm, so it may no longer prime the nuclear p53 for degradation. It works as a clock, giving the cell some time for DNA repair. However, when DNA repair is inefficient, the active p53 rises to a high level and triggers transcription of proapoptotic genes. As a result, small DNA damage may be repaired and the cell may return to its initial "healthy" state, while the extended damage results in apoptosis. The stochasticity of p53 regulation, introduced at the levels of gene expression, DNA damage and repair, leads to high heterogeneity of cell responses and causes cell population split after irradiation into subpopulations of apoptotic and surviving cells, with fraction of apoptotic cells growing with the irradiation dose.

Family tree of Markov models in systems biology

Motivated by applications in systems biology, a probabilistic framework based on Markov processes is proposed to represent intracellular processes. The formal relationships between different stochastic models referred to in the systems biology literature are reviewed. As part of this review, a novel derivation of the differential Chapman-Kolmogorov equation for a general multidimensional Markov process made up of both continuous and jump processes, is presented. First, the definition of a time-derivative for a probability density is focused, but placing no restrictions on the probability distribution, in particular, it is not assumed to be to be confined to a region that has a surface (on which the probability is zero). In this derivation, the master equation gives the jump part of the Markov process and the Fokker-Planck equation gives the continuous part. As a result, a 'family tree' for stochastic models in systems biology is sketched, providing explicit derivations of their formal relationship and clarifying assumptions involved.