Correlation between external regulators governs the mean-noise relationship in stochastic gene expression

1.   Introduction

Gene expression is a complex stochastic process where transcription as a key step occurs often in a bursty manner, creating cell-to-cell variation in mRNA and further protein abundance ^[1,2,3], In unicellular organisms, this variation improves fitness by generating phenotypic differences within a population of genetically identical cells, thus enabling a rapid response to fluctuating environments ^[4,5]. In multi-cellular organisms, the variation plays an important role in development since it allows identical progenitor cells to acquire very different fates ^[6,7]. Owing to this functional importance of variation, an important task in the post-genome era is to identify and dissect the molecular mechanisms that generate and control variation.

Recent studies ^{[8,9,10,11,12,13,14,15,16,17,18]} have focused on the relationship between the cell-to-cell variance ($\sigma^{2}$) of mRNAs or proteins and its mean expression level ($\mu$). This relationship is of particular importance for probing the evolutionary origins and consequences of noise in gene expression since it allows biologists to correctly normalize mRNA or protein variations to identify which mRNAs or proteins display unexpectedly high single-cell variances around their mean levels. Experimental data have indicated that the variance is linearly related to the mean when plotted on a log-log scale, suggesting a power-law like relationship (${{\sigma}^{2}}\propto{{\mu}^{k}}$).The dissection and interpretation of this relationship have been aided by analyzing theoretical or experimental models of gene expression ^{[19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36]}, but conflicting views or evidence have also appeared: (1) The existence of a mean-noise trend line, which is followed by a large number of genes, seems to be in conflict with the idea that gene-specific mechanisms of transcriptional regulation (for which there is an abundance of in vitro biochemical evidence ^[24]) can lead to a characteristic noise signature, quantified by the dependence of variance on mean; (2) Several experimental studies ^[8,9] implies that the slope of the log-log line for the relationship between noise (${{{\sigma }^{2}}}/{{{\mu }^{2}}}$) and mean ($\mu$) is always negative whereas another experimental study ^[10] indicates that this slope may be positive or negative (referring to Figure 1A).

Given these conflicting views or evidences, two important questions then arise: what molecular mechanisms underlying expression noise govern the mean-noise relationship? To what degree is this relationship universal? In this article, we address these issues and reveal that the sign of the slope is controlled completely by the correlation coefficient between external regulators, but independent of promoter structure.

2.   Materials and methods

Since the size of relative fluctuations (i.e., the noise intensity) defined as the ratio of variance to mean can quantify the noise level more reasonably than variance ^[25], we mathematically formulate and dissect the mean-noise relationship as

where $\alpha$ and $\beta$ are two parameters that represent the level of internal noise, constant $C$ represents the level of extrinsic noise such as cell cycle and DNA replication. It was noted that when the average expression level of molecular was low, the uncertainty of molecular fluctuation caused the internal noise of cell expression to be much higher than the external noise. Therefore, we ignore the extrinsic noise $C$ and rewrite (2.1) as follow

In principle, these parameters can be easily estimated from experimental data (Figure 1A), as done in the literature ^{[8,9,10,11,12,13,14,15,16]}, but molecular mechanisms influencing them are unclear. Note that in contrast to $\alpha$, $\beta$ is a more important parameter since its sign determines the qualitative dependence of the noise level on the mean whereas its size determines the degree of the dependence of the former on the latter. We use gene expression model to study mean-noise relationship. In fact, the degradation of mRNA is much faster than that of protein ^[37]. The two stage process (transcription and translation) can be integrated into one single step process with translational bursting ^[38]. Meanwhile, all the models used in this paper assume proteins are produced one molecule at a time (i.e., without translational bursting), although it is a common mode in eukaryotes and prokaryotes ^[21,22]. These conditions are sufficient to ensure that protein production can be done in one step ^[39], which has been applied in some theoretical models ^[40]. Furthermore, stationary mRNA levels also can well explain protein levels in many conditions ^[39]. Therefore, it is reasonable to directly regard the expression of mRNA as protein abundance in the gene expression model.

2.1.   Exploratory analysis of the mean-noise relationship with arbitrary promoter structure

First, we want to explore whether the architecture of promoters affects the mean-noise relationship. We introduce the multi-state promoter gene expression model without regulation ^[24], which is an extension of the common ON-OFF model ^{[26,27,28,29]}. The advantage of this model is that it can have any number of ON and OFF states instead of the original two limited states, and the states can be switched arbitrarily between each other, which is more in line with the actual complex promoter in biology.

We assume that the promoter of a gene contains a total of $N$ states, including $L$ active states and $K = N-L$ inactive states. Figure 1B shows a five-state gene expression model containing 3 inactive states and 2 active states. In order to simplify the model, we further assume that the promoter will translate only in the active state, regardless of the leakage of mRNA in the inactive state and their transition rate among states accords with the Markov process, which is a constant. Then for the degradation rate of protein, we assume it obeys the linear hypothesis. To sum up, the biochemical reaction of the promoter multi-state model can be written as

where ${{G}_{i}}$ and ${{G}_{j}}$ represent two different states of the promoter, respectively. ${{X}}$ represents the protein. Let ${{P}_{k}}\left(m; t \right)$ denote the probability distribution of the number of protein $m$ molecules in the ${{k}^{th}}$ promoter state at $t$ time, then $P = {{\left({{P}_{1}}, \cdots, {{P}_{N}} \right)}^{T}}$ is a column vector composed of $N$ probability distributions. The transition rate between promoters of each state is described by a $N\times N$ transition matrix $A = \left({{k}_{ij}} \right)$, where ${{k}_{ij}}$ represents the transition rate from state $i$ to state $j$. We use two diagonal matrices $\Lambda = diag\left({{\lambda }_{1}}, \cdots, {{\lambda }_{N}} \right)$ and $\delta = diag\left({{\delta }_{1}}, \cdots, {{\delta }_{N}} \right)$ to describe the rate of translation and degradation of protein, where ${{\lambda }_{i}}$ is the translation rate of promoter state $i$ (when it is in inactive state, ${{\lambda }_{i}} = 0$), ${{\delta }_{i}}$ is the degradation rate of protein produced by promoter state $i$. Therefore, under this multi-state promoter model, the chemical master equation can be described in the following form:

where $E$ is the step operators and ${{E}^{-1}}$ is its inverse, and $I$ is the identity operator. For this chemical master equation, we can use the convergent moment approach to reconstruct the stationary distribution ^[42]. Denoted by ${{b}_{1}}$ and ${{b}_{2}}$ the first- and second-order binomial moments, the mean expression level ${{\mu }_{p}}$ and the intensity of the noise $\eta _{p}^{2}$ in protein can be computed as follow:

Next, we present numerical results. We computationally find that promoter sequence influences $\alpha$ more remarkably than $\beta$. In particular, the sign of $\beta$ is always negative and independent of promoter structure (Figure 1C). This indicates that promoter structure primarily affects the intercept rather than the slope in the log-log line described by (2.2).

2.2.   External regulators control mean-noise relationship

The expression of a gene is inevitably regulated by intracellular or/and extracellular factors. These factors would be so complex that any existing theoretical models of gene expression are flawed. Here we particularly assume that for gene regulation, except for their own product proteins to regulate themselves, all other regulations belong to external regulation. So the transcription factor is also one of the external regulators (ERs). In order to reveal the essential molecular mechanism of how ERs regulate the relationship between mean and noise, here we introduce a simple yet representative gene model (Figure 2A), which captures main molecular events occurring in gene expression. This model assumes that the gene has one active (ON) and one inactive (OFF) states with bidirectional transitions between them. This model is a special case of (2.3) mentioned above. Thus, the gene model used is described by

where ${{G}_{ON}}$ and ${{G}_{OFF}}$ represent the active and inactive states of the promoter, respectively. $X$ represents the protein. ${{k}_{OFF}}$ and ${{k}_{ON}}$ are transition rates from inactive and active states and vice versa. $\lambda$ is the rate of translation when the gene is in the active state and $\delta$ is the rate of protein decay. Note that if ${{k}_{OFF}}$ and ${{k}_{ON}}$ are constants, the corresponding gene expression model is just the common ON-OFF model. If any transition rate between ON and OFF is not a constant, this actually corresponds to the translational regulation in biology. Therefore, we consider two ERs regulate the gene expression by adjusting two transition rates between ON and OFF states. This setting is used to model either the effect of regulation from the other genes inside the cell or uncertainty of extracellular environments or both. In this situation, ${{k}_{OFF}}$ and ${{k}_{ON}}$ are functions of concentration of ERs.

For clear demonstration, let us first consider two extreme cases: (1) Two distinct ERs independently regulate two transition rates between promoter activity states (Figure 2A) with the settings ${{k}_{ON}} = {{f}_{1}}\left(E{{R}_{1}} \right) = {1}/{\left(1+ER_{1}^{{{n}_{1}}} \right)}$ and ${{k}_{OFF}} = {{f}_{2}}\left(E{{R}_{2}} \right) = {1}/{\left(1+ER_{2}^{{{n}_{2}}} \right)}$; (2) A common ER simultaneously regulates the transition rates (Figure 2B) with the settings ${{k}_{ON}} = {{f}_{1}}\left(ER \right) = {1}/{\left(1+E{{R}^{{{n}_{1}}}} \right)}$ and ${{k}_{OFF}} = {{f}_{2}}\left(ER \right) = {1}/{\left(1+E{{R}^{{{n}_{2}}}} \right)}$. In each case, ERs take values randomly in a finite range. In Figure 2, $E{{R}_{1}}$ and $E{{R}_{2}}$ are generated in the interval $\left(0.001, 1000 \right)$. We find that the mean-noise relationship is characterized by a negative slope in case (1) (Figure 2E), and by a negative or positive slope in case (2) (Figure 2G). These numerical results are in accord with experimental observations ^[8,10]. In both cases, the dependence of the proteins mean levels on the inducers concentration exhibits different characteristics (comparing Figure 2D with Figure 2F).

Then, consider a more general case: two ERs regulate two transition rates between the promoter activity states in a correlated manner with a correlation coefficient $r$ (Figure 2H). In this case, two regulated parameters are set as ${{k}_{ON}} = 4+2ER_{1}^{{{n}_{1}}}$ and ${{k}_{OFF}} = 1+ER_{2}^{{{n}_{2}}}$, where $E{{R}_{1}}$ and $E{{R}_{2}}$ are generated randomly and with a correlation coefficient $r$ in the interval $\left(0.001, 1000 \right)$. We can see that small $r$'s can correspond to the model in Figure 2A whereas large $r$'s to the model in Figure 2B. To reveal how regulation affects the mean-noise relationships, we study the dependence of the slope on the correlation coefficient. Figure 2I shows that the sign of $\beta$ is controlled completely by correlation: small or moderate $r$'s correspond to negative slopes whereas $r$'s close to 1 correspond to positive slopes. We choose four representative points on the $\beta$ vs $r$ curve, three of which are located at the shadowed part with negative or zero slopes and the one is located at the outside of the shadowed part with a positive slope, to demonstrate different characteristics of the ER-regulated mean-noise relationship (Figure 2J). Interestingly, this relationship is randomly distributed if the slope is zero (Figure 2J).

2.3.   Theoretical analysis of the mean-noise relationship

Next, we perform theoretical analysis, focusing on the essential mechanism governing the mean-noise relationship. In general, the protein noise can be divided into two parts

where $\eta _{int}^{2}$ represents noise generated in the process of synthesis and degradation of protein, whereas $\eta _{prom}^{2}$ represents the noise comes from the promoter architecture.

To reveal a general mechanism underlying the mean-noise relationship, we consider a gene model as which is different from the model shown in Figure 1B:

where OFF-state and ON-state lifetimes follow general distributions with probability density functions denoted by ${{f}_{off}}\left(t \right)$ and ${{f}_{on}}\left(t \right)$ ^[30,35]. Assume that translation and protein degradation are described as first-order processes with rate constants $\mu$ and $\delta$ respectively ^{[29,33,34,35,36]}. We point out that the corresponding gene model can be taken as a $G{{I}^{X}}/M/\infty$ system known in the queueing literature (see ^[43] for details). To determine the protein mean and noise and further the mean-noise relationship in this model, the key is to derive analytical expressions for the first two moments of the protein number distribution.

Suppose that there are $N$ proteins at initial time $t = 0$. At time $t$, every protein with the molecule number denoted by ${{X}_{i}}$ ($1\le i\le N$) has a survival probability $P = {{e}^{-\delta t}}$, and the total number of protein molecules is given by $S = {{X}_{1}}+\ldots +{{X}_{N}}$. Denote by $W\left(z; t \right)$ the moment-generating function (MGF) of the protein distribution at time $t$. In particular, denote ${{W}_{init}}\left(z \right)\text{ = }W\left(z; 0 \right)$. Assume that these random variables are independent of one another, with each following a Bernoulli distribution with MGF given by $1+p\left({{e}^{z}}-1 \right)$. Applying the theorem of total expectation in probability theory, we find that $W\left(z; t \right)$ can be expressed by $W\left(z; 0 \right)$ ^[29,30]. That is,

With (2.9) combined with a boundary condition for one cycle of OFF and ON states, we can further arrive at the following integral equation of ${{W}_{init}}\left(z \right)$

where $\rho \left(u \right) = {\lambda \left(1-{{e}^{-\delta u}} \right)}/{\delta }$ represents the transient average for a birth-death process associated with creation and degradation of protein. This equation is a pivot for calculating the first two moments of the steady-state protein distribution.

Let ${{\langle {{m}^{k}} \rangle }_{s}}$ and ${{\langle {{m}^{k}} \rangle }_{u}}$ be the ${{k}^{th}}$ moment of the protein copy number distribution at time $t = s$ during the OFF-state and at time $t = {{\tau }_{off}}+u$ during the ON-state, respectively. Then, moments of the steady-state distribution of the total protein number are calculated according to

where $\langle {{\tau }_{off}} \rangle = \int_{0}^{\infty }{t{{f}_{off}}\left(t \right)dt}$ and $\langle {{\tau }_{on}} \rangle = \int_{0}^{\infty }{t{{f}_{on}}\left(t \right)dt}$ are the mean times that the gene dwells at OFF and ON states respectively, whereas $\int_{s}^{\infty }{{{f}_{off}}\left(t \right)dt}$ and $\int_{u}^{\infty }{{{f}_{on}}\left(t \right)dt}$ represent survival probabilities (which are actually the cumulative distribution functions of the distributions of off-state and on-state lifetimes, respectively). By (2.11), we find that the exact mean and noise intensity for protein can be expressed separately as (see Appendix)

and

where ${{L}_{off}}\left(s \right)$ and ${{L}_{on}}\left(s \right)$ are Laplace transforms of ${{f}_{off}}\left(t \right)$ and ${{f}_{on}}\left(t \right)$, respectively. The (2.13) highlights how ON and OFF lifetime distributions contribute to the promoter noise and further the protein noise.

While (2.13) is valid for general waiting-time distributions, it is of interest to consider specific examples. Consider the common ON-OFF model of gene expression at the translation level, where transition rates between promoter activity states are assumed as two constants. In this case, ${{f}_{off}}\left(t \right) = \left({1}/{\langle {{\tau }_{off}} \rangle }\; \right){{e}^{-\left({1}/{\langle {{\tau }_{off}} \rangle }\; \right)t}}$ and ${{f}_{on}}\left(t \right) = \left({1}/{\langle {{\tau }_{on}} \rangle }\; \right){{e}^{-\left({1}/{\langle {{\tau }_{on}} \rangle }\; \right)t}}$. Therefore, the protein mean is given by $\langle m \rangle = {\left(\lambda \langle {{\tau }_{on}} \rangle \right)}/{\left(\langle {{\tau }_{off}} \rangle +\langle {{\tau }_{on}} \rangle \right)}$, which is a known result ^[24], whereas the promoter noise is given by $\eta _{\text{promoter}}^{2} = {{{\langle {{\tau }_{off}} \rangle }^{2}}}/{\left(\langle {{\tau }_{on}} \rangle +\langle {{\tau }_{off}} \rangle +\langle {{\tau }_{on}} \rangle \langle {{\tau }_{off}} \rangle \right)}\;$, which is also a known result ^[24].

After having derived the analytical expressions for protein mean and noise, we next determine two parameters $\alpha$ and $\beta$ in the mean-noise relationship described by (2.2). For this, if denote $x = \log \left(\langle m \rangle \right)$ and $y = \log \left(\eta _{\operatorname{int}}^{2}+\eta _{\text{prom}}^{2} \right) = \log \left(\eta _{m}^{2} \right)$, then we set

where $\operatorname{cov}\left(x, y \right)$ represents the covariance between $x$ and $y$, and $\sigma _{x}^{2}$ represents the variance of $x$ (similarly, we let $\sigma _{y}^{2}$ represent the variance of $y$). For a gene model with general waiting-time distributions, it is difficult to determine the sign of $\beta$. For the preceding ON-OFF model, however, it is possible. In fact, we first observe that the protein noise is necessarily correlated to the protein expression level through two random variables ${{\tau }_{on}}$ and ${{\tau }_{off}}$. Then, we note that stochastic variable $y$ is inversely correlated to stochastic variable $x$ because of $y = -x+\log \left(1+{\eta _{ext}^{2}}/{\eta _{\operatorname{int}}^{2}}\; \right)$. Thus, from the expression $\beta = \frac{{{\sigma }_{y}}}{{{\sigma }_{x}}}\frac{\operatorname{cov}\left(x, y \right)}{{{\sigma }_{x}}{{\sigma }_{y}}} = r\frac{{{\sigma }_{y}}}{{{\sigma }_{x}}}$, we know $\beta < 0$, where $r$ represents the correlation coefficient between $x$ and $y$. Such analysis indicates that in the unregulated case, $\beta$ is always negative, in accordance with the results obtained by analysis of experimental data as well as with the above numerical results. Moreover, promoter structure has little effect on the size of $\beta$, as shown in the above numerical simulation. Note that we demonstrate theoretically that $\beta$ is always negative in the unregulated case, but cannot give a theoretical result in correlation regulated case because of mathematical difficulty. Nevertheless, we have numerically show that mean-noise relationship is controlled by correlation between transcription factors.

3.   Conclusions and discussion

In summary, we have revealed the essential molecular mechanism that controls two parameters $\alpha$ and $\beta$ including the sign of $\beta$ in the empirical mean-noise relationship formulated by $\log \left({{\eta }^{2}} \right) = \beta \log \mu +\log \alpha$. This mechanism revealed not only can well interpret diverse experimental phenomena in the existing literature but also can be used in analysis of sources of phenotypic heterogeneity in isogenic populations. The settings of $\alpha$ and $\beta$ can provide insight into how detailed biochemical processes of gene expression impact the dependence of expression noise on mean. In addition, these settings can be used in the analysis of inverse problems, e.g., using experimental measurements of expression noise to determine parameters of the underlying kinetic models. Such efforts, in turn, can lead to further insights into cellular factors that impact gene regulation, based on experimental observations of expression noise ^{[44,45,46,47]}.

However, there are still some limitations in our model. First, we use a model without translation bursts, though they are common in prokaryotic and eukaryotic cells ^[21,22]. Hence we ignore the translational bursting noise. Second, in our article, we do not consider the external noise. External noise plays an important role in biological functions, such as cell cycle and DNA replication ^[31,48,49]. The defects of these models can be explored as a clue in future work.

This study would have a number of potential applications. For example, from a pharmaceutical science and drug screening perspective, simply checking the mean-noise dependence line in the log space (i.e., so-called noise screening) presents an orthogonal axis to detect synergistic drug combinations ^[14,50]. Compared with random synergy screening, noise screening requires substantially fewer tests. In fact, blind synergy searches for pairwise combinations of $M$ compounds require $\sim {{M}^{2}}$ tests whereas noise screening needs only $\sim M$ tests. In addition, noise screening might help identify compounds for manipulating other fate-switching phenotypes such as cellular reprogramming, metastasis, and bacterial persistence.

Acknowledgments

This work was supported by grants 91530320, 11775314, 11475273, and 11631005 from Natural Science Foundation of P. R. China; 2014CB964703 from Science and Technology Department, P. R. China; 2020B1212060032 from Guangdong Province Key Laboratory of Computational Science at the Sun Yat-sen University.

Conflict of interests

The authors have declared that no competing interest exists.

Appendix

As is well known, transcription and translation occur often in a burst fashion, leading to variation in protein. Translational bursting kinetics can be characterized by burst frequency and burst size. Suppose the time intervals between bursts are independent of the burst size. These intervals are described by the distribution of the off-state life times. This description corresponds to a $G/M/\infty$ queue with batch arrival in queuing theory [1].

Consider a gene model [2, 3], where the promoter is assumed to have one on-state and one off-state (i.e., the so-called on-off model) but off and on state life lengths are assumed to follow general distributions, denoted respectively by $f\left({{\tau }_{off}} \right)$ and $g\left({{\tau }_{on}} \right)$. For simplicity, we consider only zero-order translation and first-order protein degradation kinetics, with constant translation and degradation rates denoted respectively by $\lambda$ and $\delta$.

Denote the moment-generating function of the protein distribution at time $t$ as $W\left(z; t \right)$. In addition, denote ${{W}_{init}}\left(z \right) = W\left(z; 0 \right)$. Since each protein molecule has a survival probability of ${{e}^{-\delta t}}$, the moment-generating function at time $t$ can be expressed as a function of ${{W}_{init}}\left(z \right)$ [3–5], i.e.,

In particular, at the end of an off-state, $W\left(z; t \right)$ can be expressed as

Except that an on-period degradation of the protein molecules that have been present at the beginning of the burst continues as described through the function $1+{{e}^{-\delta t}}\left({{e}^{z}}-1 \right)$ with ${{W}_{\deg }}\left(z; {{\tau }_{off}}+u \right) = {{W}_{init}}\left(\log \left(1+{{e}^{-\delta \left(u+{{\tau }_{off}} \right)}}\left({{e}^{z}}-1 \right) \right) \right)$, proteins are also created and degraded according to a birth-death process with exponential waiting times. Moreover, the transient distribution for this birth-death process is a Poisson distribution with average ${{\lambda }_{u}} = {\lambda \left(1-{{e}^{-\delta u}} \right)}/{\delta }\;$ and moment-generating function ${{W}_{effect}}\left(z; {{\tau }_{off}}+u \right) = {{e}^{{{\lambda }_{u}}\left({{e}^{z}}-1 \right)}}$ [6]. The combination of them produces an effective burst size. During an on-state, the probability distribution of the protein number is given by the convolution of the distribution of the number of those molecules that are still present from previous bursts and the effective burst-size distribution. Thus, we have

Note that complete of one cycle (off- and on-state) defines a boundary condition:

With (5.4), we obtain the following integral equation with respect to ${{W}_{init}}\left(z \right)$

In general, solving (5.5) is very difficult for general waiting time distributions $f\left(t \right)$ and $g\left(t \right)$. However, we can derive the formal expressions for first two moments of the protein distribution at time $t = 0$. In fact, differentiating both sides of (5.5) at $z = 0$ yields

where ${{\rho }_{u}} = {\lambda \left(1-{{e}^{-\delta u}} \right)}/{\delta}$. According to the definition of moment-generating function, we know that ${{W}_{init}}\left(0 \right) = 1$ but ${{{W}'}_{init}}\left(0 \right)$ is equal to the protein mean at time $t = 0$, denoted by ${{\langle m \rangle }_{init}}$. Thus, we obtain

Similarly, we have the expression for the second-order moment at time $t = 0$

Finally, using these expressions for moments at the initial moment, we can calculate first two steady-state moments of protein. First, note that the moments of the full protein distribution can be obtained by averaging over the moments at times weighted according to the probability that the next burst has not yet occurred at time $t$, i.e., the survival probability ${{P}_{s}}\left(t \right) = 1-F\left(t \right)$ with $F\left(t \right)$ being the cumulative distribution function of the distribution of off-state life times defined as $F\left(t \right) = \int_{0}^{t}{f\left(s \right)ds}$. Similarly, we define $G\left(t \right) = \int_{0}^{t}{g\left(s \right)ds}$. Then, to obtain the moments of the total protein number distribution, we need to average over all times according to their probabilities. This can be expressed as

with $k = 1, 2$, where $F\left(s \right)$ and $G\left(s \right)$ are the cumulative distribution functions of the duration distributions of OFF- and ON-states, respectively. In (5.9), ${{\langle {{m}^{k}} \rangle }_{s}}$ is the ${{k}^{th}}$ moment of the protein copy number distribution at time $S$ during the off-state, which is given by the ${k}$ order derivative of (5.1) at $z = 0$, ${{\langle {{m}^{k}} \rangle }_{u}}$ is the ${{k}^{th}}$ moment of the protein copy number distribution at time $t = \langle {{\tau }_{off}} \rangle +u$ during the on-state, which is given by the ${{k}}$ order derivative of (5.3) at $z = 0$. By calculation, we find these quantities are

First, we derive the analytical expression for the protein mean. Note that

Because of

we have

Therefore,

which is the same as in the case of no regulation.

Then, we derive the analytical expression for the mRNA noise intensity. The key is to give the second-order moment $\langle {{m}^{2}} \rangle$. Note that

Therefore

Using (5.8), we obtain

By calculation, we know

and

Substituting them into (5.13) yields

Note that

and

Thus, we have

Using (5.7), we thus obtain the analytical expression for the protein noise intensity

If $f\left(s \right)$ and $g\left(u \right)$ follow exponential distribution, i.e., $f\left(s \right) = \left({1}/{\langle {{\tau }_{off}} \rangle }\; \right){{e}^{-{s}/{\langle {{\tau }_{off}} \rangle }\; }}$ and $g\left(u \right) = \left({1}/{\langle {{\tau }_{on}} \rangle }\; \right){{e}^{-{s}/{\langle {{\tau }_{on}} \rangle }\; }}$, then we have

and

Therefore, we finally obtain

which is a know result [7].

References

1. L. Liu, B. R. K. Kashyap, J. G. C. Templeton, On the $G{{I}^{X}}/G/\infty$ system, J. Appl. Probab., 27 (1990), 671–683.

2. T. Jia, R. V. Kulkarni, Intrinsic noise in stochastic models of gene expression with molecular memory and bursting, Phys. Rev. Lett., 106 (2011), 058102.

3. J. J. Zhang, T. S. Zhou, Stationary moments, distribution conjugation and phenotypic regions in stochastic gene transcription, Math. Biosci. Eng., 16 (2019), 6134–6166.

4. E. Cinlar, Introduction to stochastic processes, in Prentice Hall: Englewood Cliffs, (1975).

5. A. Schwabe, K. N. Rybakova, F. J. Bruggeman, Transcription stochasticity of complex gene regulation models, Biophys. J., 103 (2012), 1152–1161.

6. M. Hemberg, H. Barahona, Perfect sampling of the master equation for gene regulatory networks, Biophys. J., 93 (2007), 401–410.

7. J. J. Zhang, L. N. Chen, T. S. Zhou, Analytical distribution and tunability of noise in a model of promoter progress, Biophys. J., 102 (2012), 1247–1257.