Efficient analysis of stochastic gene dynamics in the non-adiabatic regime using piecewise deterministic Markov processes

Abstract

Single-cell experiments show that gene expression is stochastic and bursty, a feature that can emerge from slow switching between promoter states with different activities. In addition to slow chromatin and/or DNA looping dynamics, one source of long-lived promoter states is the slow binding and unbinding kinetics of transcription factors to promoters, i.e. the non-adiabatic binding regime. Here, we introduce a simple analytical framework, known as a piecewise deterministic Markov process (PDMP), that accurately describes the stochastic dynamics of gene expression in the non-adiabatic regime. We illustrate the utility of the PDMP on a non-trivial dynamical system by analysing the properties of a titration-based oscillator in the non-adiabatic limit. We first show how to transform the underlying chemical master equation into a PDMP where the slow transitions between promoter states are stochastic, but whose rates depend upon the faster deterministic dynamics of the transcription factors regulated by these promoters. We show that the PDMP accurately describes the observed periods of stochastic cycles in activator and repressor-based titration oscillators. We then generalize our PDMP analysis to more complicated versions of titration-based oscillators to explain how multiple binding sites lengthen the period and improve coherence. Last, we show how noise-induced oscillation previously observed in a titration-based oscillator arises from non-adiabatic and discrete binding events at the promoter site.

Keywords: circadian rhythm; gene expression; intrinsic noise; mathematical model; stochastic cycles.

Figure 1.

Schematic diagrams of the idealized (a) activator-titration circuit (ATC) and (b) repressor-titration circuit (RTC). Protein X is a transcription factor and Y is an inhibitor that can irreversibly associate with X to form an inactive complex. In the ATC, X is an activator that can sequentially bind multiple DNA sites in the promoter of gene Y and increase transcription of the inhibitor. In the RTC, X is a repressor that can bind its own promoter and repress transcription. (Online version in colour.)

Figure 2.

Sample paths of the full CME of the ATC and RTC in the (a) adiabatic regime (λ = 1000) and (b) non-adiabatic regime (λ = 1) for a single binding site (). (c) Sample paths of the constructed piecewise deterministic Markov process when λ = 1 (§2.4). (d) The alternative deterministic limit of the processes (§2.8). (e) Quantification of the periods of the stochastic cycles. (Online version in colour.)

Figure 3.

Schematic diagrams of the derived piecewise deterministic Markov process (PDMP) for (a) idealized activator-titration circuit (ATC) and (b) idealized repressor-titration circuit (RTC). Both models have a single promoter site . The linearized PDMP for ATC and RTC are shown in (c) and (d), respectively, where the green circular arrows indicates the direction of the emergent stochastic cycles. Dark blue and light red boxes denote promoter states with different production rates where and , respectively. The transitions between the two regimes and () are due to titration. (Online version in colour.)

Figure 4.

The waiting time distributions of the next dissociation and binding event in the linear PDMP. The waiting time distribution of the dissociation event is exponential. The waiting time distribution of the binding event is derived from the survival function (2.11), and depends on the initial concentration of x. (Online version in colour.)

Figure 5.

Schematic diagrams of the KB model [15] and the VKBL model [28] of the activator-titration circuit (ATC). The parametrizations were adopted from the original papers. (Online version in colour.)

Figure 6.

Schematic diagram of the linearized PDMP describing the idealized ATC with multiple binding sites () in the non-adiabatic regime. The green (circular) arrows indicates the emergent cycles which are predominantly observed in the simulations. The path of the dotted arrow is also observed, but less frequently (figure 11b). (Online version in colour.)

Figure 7.

The linearized PDMP for the KB model [15]. To make the expressions compact, we define and . Green (circular) arrows indicate the direction of the emergent cycle which are frequently observed. (Online version in colour.)

Figure 8.

The reduced PDMP approximating the KB model [15]. The waiting time of the deterministic activation is computed as the first moment of the cumulative distribution equation (E 2). (Online version in colour.)

Figure 9.

Numerically measured data of the KB model [15]. (a) A sample path of the full CME, (b) a sample path of the linear PDMP (figure 7). The stochastic cycles of the model were measured using the protocol provided in appendix B. The full CME exhibited a minor fraction of short period cycles in the genetic states. To quantify the predominant longer periods (100 ≲ period), we discarded any periods less than 50 to generate panels (c,d). Panel (c) presents the coefficient of variation (CV) of the stochastic periods measured in the full CME as a function of the scaling factor ϖ. Each point was computed from 10⁴ stochastic cycles. The larger the ϖ, the more short-lived is the mRNA. The long-lived mRNA introduces a delay in protein production and improves the CV. Panel (d) presents the probability distribution of the stochastic periods when the scaling factor ϖ = 1, as measured from 10⁵ stochastic cycles of the full CME, PDMP and reduced PDMP (figure 8). (Online version in colour.)

Figure 10.

Slow fluctuations in the random binding and unbinding events induce the excitable mode of the VKBL model [28]. Column (a) is the full CME, (b) is the PDMP and (c) is the ADL of the VKBL model. The first row are the dynamics of the mRNAs and the second are the dynamics of the activator (A), inhibitor (R) and heterodimer complex (C). Insets are longer time series. (Online version in colour.)

Figure 11.

Sample paths of the full CME model of the ATC and RTC in the (a) adiabatic regime (λ = 1000) and (b) non-adiabatic regime (λ = 1) for multiple binding sites (. (c) Sample paths of the constructed piecewise deterministic Markov process when λ = 1 (§2.4). (d) The alternative deterministic limit of the processes (§2.8). (e) Quantification of the periods of the stochastic cycles. (Online version in colour.)