Time-dependent propagators for stochastic models of gene expression: an analytical method (Q1659718)

The authors of this very interesting paper develop an analytical method for the efficient evaluation of time-dependent propagators in stochastic gene expression models, for arbitrary values of the model parameters. Two different stochastic models for gene expression are given. The first one, referred to as ``model \textbf{A}'' studied by \textit{S. Iyer-Biswas} and \textit{C. Jayaprakash} [``Mixed Poisson distributions in exact solutions of stochastic autoregulation models'', Phys. Rev. E (3) 90, No. 5, Article ID 052712, 5 p. (2014; \url{doi:10.1103/PhysRevE.90.052712})], incorporates autoregulation. The second model known as ``model \textbf{B}'', studied by \textit{V. Shahrezaei} and \textit{P. S. Swain} [``The stochastic nature of biochemical networks'', Curr. Opinion Biotechnol. 19, No. 4, 369--374 (2008; \url{doi:10.1016/j.copbio.2008.06.011})], describes both mRNA and protein explicitly. Both models are formulated in terms of the so called chemical master equation (CME for short). The basic stochastic model for gene expression is represented by the known reaction scheme containing three stages: a) \(\{\)DNA switching\(\}\), \(\{\)Production of protein\(\}\) and \(\{\)Decay of protein\(\}\). The gene can hence switch between the inactive state \(D\) and the active state \(D^{\ast }\), with switching rates \(c_f\) and \(c_b\), respectively. The active gene produces protein \(P\) with rate \(p_b\), while protein decays with rate \(p_d\). The autoregulatory part of this model is implemented as either positive or negative feedback: autoactivation and autorepression. When there exists autoactivation the protein induces activation of the gene with some activation rate \(a\), thereby it accelerates its own production. If there exists autorepression, then the protein deactivates the active gene with repression rate \(r\), impeding its own production. The authors investigate the CME system that is associated to the above stated reaction, with autoactivation that is modeled by two differential equations having the form: \[ dP_n^{(0)}/dt=-\tilde{k}P_n^{(0)}+ k_bP_n^{(1)}+\tilde{P}^{(0)} \] and \[ dP_n^{(1)}/dt=\tilde{k}P_n^{(0)}- k_bP_n^{(1)}+\tilde{P}^{(1)}+\lambda [P_{n-1}^{(1)}-P_n^{(1)}], \] where \(\tilde{k} = k_f+an/p_d\), \(\tilde{P}^{(j)}=(n+1)P_{n+1}^{j}-nP_{n}^{j}\) (\(j=0,1\)), \(P_{n}^{j}=P_{n}^{j}(t)\) (\(j=0,1\)) is the probability of \(n\) protein at time \(t\) while the gene is either inactive (\(0\)) or active (\(1\)). The time variable is nondimensionalised by the protein decay rate \(p_d\); \(k_f=c_f/p_d\), \(k_b=c_b/p_d\), \(\lambda =p_b/p_d\) are some parameters. Analogously, the CME system for the case of autorepression contains also two similar differential equations. To analyze above stated systems of differential equations the authors introduce the probability-generating function used for the first time by \textit{C. Gardiner} [Stochastic methods. A handbook for the natural and social sciences. 4th revised and augmented ed. Berlin: Springer (2009; Zbl 1181.60001)]. The second model \textbf{B} is a stochastic gene expression model presented by Shahrezaei and Swain [loc. cit.] which explicitly incorporates the transcription stage in the expression process, as well as DNA switching. The following reaction scheme shows the consecutive operations: \( \{\text{DNA switching}\}\), \(\{\text{transcription of DNA to mRNA}\}\), \(\{\text{translation of mRNA to protein}\}\), \(\{\text{decay of mRNA}\}\) and \(\{\text{decay of protein}\}\). As in the model \textbf{A}, autoregulatory terms can be added to the above stated core reaction scheme: \(\{\text{autoactivation through mRNA}\}\), \(\{\text{autorepression through mRNA}\}\), \(\{\text{autoactivation through protein}\}\) and \(\{\text{autorepression through protein}\}\). The autoactivation can be achieved either by mRNA or by protein, with certain rates. Similarly, autorepression can occur either through mRNA or through protein, with respective rates. Here, it is important to note that it is not possible to obtain exact solutions to the CME in explicit form. Therefore, the authors propose an original analytical method for the efficient approximation of the propagator which describe the evolution of such network in time. It is defined as a solution of the corresponding CME. The following procedure is applied to the models under consideration: (i) a probability-generating function is defined; it transforms the CME into a system of partial differential equations; actually they are coupled, linear, first-order, hyperbolic partial differential equations; (ii) application of the method of characteristics then yields a dynamical system containing ordinary differential equations by which obtain the generating function; (iii) the propagator probabilities is reconstructed numerically via the Cauchy integral formula.