# A Stochastic Analysis of Autoregulation of Gene Expression

Renaud Dessalles, Vincent Fromion, Philippe Robert

► To cite this version:

Renaud Dessalles, Vincent Fromion, Philippe Robert. A Stochastic Analysis of Autoregulation of Gene Expression. Journal of Mathematical Biology, 2017, 75 (5), pp.1253-1283. 10.1007/s00285-017-1116-7. hal-01203076

HAL Id: hal-01203076

https://inria.hal.science/hal-01203076v1

Submitted on 6 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution - ShareAlike 4.0 International License
# A STOCHASTIC ANALYSIS OF AUTOREGULATION OF GENE EXPRESSION

RENAUD DESSALLES, VINCENT FROMION, AND PHILIPPE ROBERT

**ABSTRACT.** This paper analyzes, in the context of a prokaryotic cell, the stochastic variability of the number of proteins when there is a control of gene expression by an autoregulation scheme. The goal of this work is to estimate the efficiency of the regulation to limit the fluctuations of the number of copies of a given protein. The autoregulation considered in this paper relies mainly on a negative feedback: the proteins are repressors of their own gene expression. The efficiency of a production process without feedback control is compared to a production process with an autoregulation of the gene expression assuming that both of them produce the same average number of proteins. The main characteristic used for the comparison is the standard deviation of the number of proteins at equilibrium. With a Markovian representation and a simple model of repression, we prove that, under a scaling regime, the repression mechanism follows a Hill repression scheme with an hyperbolic control. An explicit asymptotic expression of the variance of the number of proteins under this regulation mechanism is obtained. Simulations are used to study other aspects of autoregulation such as the rate of convergence to equilibrium of the production process and the case where the control of the production process of proteins is achieved via the inhibition of mRNAs.

## CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Stochastic Models of Protein Production</td><td>5</td></tr><tr><td>3. A Scaling Analysis</td><td>7</td></tr><tr><td>4. Fluctuations of the Number of Proteins</td><td>11</td></tr><tr><td>5. Discussion</td><td>13</td></tr><tr><td>References</td><td>19</td></tr><tr><td>Appendix A. Convergence Results</td><td>22</td></tr></table>

## 1. INTRODUCTION

1.1. **Biological Context.** The gene expression is the process by which genetic information is used to produce functional products of gene expression: proteins and non-coding RNAs. This paper concerns itself with the production of proteins. The information flow from DNA genes to proteins is a fundamental process. It is composed of three main steps: *Gene Activation*, *transcription* and *translation*.

(1) The initiation of transcription is strongly regulated. Schematically the gene is said to be in “inactive state” if a repressor is bound on the gene’s promoter preventing the RNA polymerase from binding and is in “active state” otherwise.

Date: October 8, 2018.
(2) When the gene is in active state, the RNA polymerase binds and initiates transcription that leads to the creation of a mRNA, a copy of a specific DNA sequence.

(3) The translation of the messenger into a protein is achieved by a large com-
plex molecule: the *ribosome*. A ribosome binds to an active mRNA, initi-
ates the translation and proceeds to protein elongation. Once the elongation
terminates, the protein is released in the medium and the ribosome is anew
available for any another translation.

The production of proteins is the most important cellular activity, both for the
functional role and the high associated cost in terms of resources. In a *E. Coli*
bacterium for example there are about 3.6 × 10⁶ proteins of approximately 2000
different types with a large variability in concentration, depending on their types:
from a few dozen up to 10⁵. The gene expression is additionally a highly stochas-
tic process and results from the realization of a very large number of elementary
stochastic processes of different nature. The three main steps are the results of a
large number of encounters of macromolecules following random motions, due in
particular to thermal excitation, in the viscous fluid of the cytoplasm. One of the
key problems is to understand the basic mechanisms which allow a cell to produce
a large number of proteins with very different concentrations and in a random con-
text. This can be seen as a problem of minimization of the variance of the number
of proteins of each type.

To study this problem, one can take a simple stochastic model, with a limited
set *S* of parameters preferably, describing the three steps of the production of a
given type of protein. Once a closed form expression of the variance of the number
of proteins is obtained, it is natural to find the parameters of the set *S* which
minimizes the variance with the constraint that the mean number of proteins is
fixed. See the survey Paulsson [25].

A more effective way to regulate the number of proteins can be of using a direct
feedback control, an *autoregulation* mechanism, so that the production of proteins
is either sped up or slowed down depending on the current number of proteins.
It should be noted that the feedback control loop can involve other intermediate
proteins to achieve this goal, like the classical lac operon, but it is not considered
here. See Yildirim and MacKey [39] for example.

The protein can regulate the gene activation simply, for example by being a
repressor and tend to bind on his own gene’s promoter. This is the *autogenous regulation* scheme. See Goldberger [11] and Maloy and Stewart [21]. See also Thattai
and van Oudenaarden [37]. Other autoregulation mechanisms are possible in cells,
such as an autoregulation on the mRNAs where a protein inhibits its own trans-
lation initiation by binding to the translation initiation region of its own mRNAs.
It occurs for example in the production of ribosomal proteins, see Kaczanowska
and Rydén-Aulin [16]. The idea being that a feedback mechanism may reduce
significantly the number of large excursions from the mean. In this paper, the
mathematical analysis will mainly focus on a negative autogenous feedback, when
the rate of inactivation of the gene expression grows with the number of proteins.

1.2. **Literature.** The classical results concerning the mathematical analysis of the variance of the number of proteins has been investigated in Berg [3] and Rigney [27, 28] and reviewed more recently by Paulsson [25], see also Raj and van Oudenaarden [26] for the biological aspects. These references use the three
stage model, the state of the system is given by three variables: the state of the promoter, the number of mRNAs and the number of proteins. Mathematically, the techniques used rely on the Fokker-Planck equations of the associated three dimensional Markov process and the observation that at equilibrium, a recurrence on the moments of the number of proteins holds. Fromion et al. [9] investigates a more general model (elongation times are not necessarily exponentially distributed in particular) and an alternative technique to a Markovian approach is introduced.

Concerning the evaluation of autoregulation, most of mathematical models use a continuous state, the rate of production of proteins depends linearly on the number of mRNAs and the rate of production of mRNAs is a function $k(p)$ exhibiting a non-linear dependence on the current number $p$ of proteins. In Rosenfeld et al. [31] and Becskei and Serrano [2], based on experiments the constant $k(p)$ is taken a *Hill repression function*, i.e. $k(p) = a/(b + p^n)$ for some constants $a$ and $b$ and $n \ge 1$ is the Hill coefficient. See also Thattai and van Oudenaarden [37]. Related models in a similar framework with further results are presented in Bokes et al. [5] and Yvinec et al. [40]. For most of these models the state of the promoter, active or inactive, which is a source of variability is not taken into account, it is in some way encapsulated in the constant $k(p)$ whose representation is rarely discussed. In Hornos et al. [14] the state of the gene expression, on or off, is taken into account but not the number of mRNAs and therefore the fluctuations generated by transcription. The parameter of activation $k(p)$ is of course crucial in our case since autogenous regulation rely on the state of the promoter which can be inactivated by proteins. Our model includes it. See also Fournier et al. [8] for some simulations of these stochastic models of autoregulation as well as some experiments.

**1.3. Results of the Paper.** The main goal of this paper is to estimate the possible benefit of the autogenous regulation to control the fluctuations of the number of copies of a given protein. The efficiency of a production process without feedback control is compared to a production process with an autoregulation of the gene expression, assuming that both of them produce the same average number proteins. The main characteristic used for the comparison is the standard deviation of the number of proteins at equilibrium. For this purpose, two approaches are used.

**MATHEMATICAL ANALYSIS.** One first studies the distribution of the number of proteins via a stochastic model. When there is no regulation, the corresponding classical mathematical model has been investigated in detail for some time now. In particular, the standard deviation of the number of proteins at equilibrium has a closed form expression in terms of the basic parameters of the production process. See for example the survey Paulsson [25], and also Fromion et al. [9].

To represent the negative feedback of the autogenous regulation, a simple model is used: each protein can be bound, at some rate and for some random duration of time, on its own gene expression. In this situation the gene expression is inactive and the transcription is not possible during that time. This amounts to say that the gene expression is deactivated at a rate proportional to the number of proteins. The activation rate is constant.

As will be seen, the mathematical model of the autogenous regulation is more complicated, in particular there is no recurrence relationship between the moments of the number of proteins at equilibrium as in the classical model of protein production process. For this reason, a limiting procedure is used, it amounts to assume
that the dynamics of the activation of the gene expression and of the evolution of mRNAs occur on a much faster time scale than the dynamics of the proteins. The values of the key parameters are presented in Section 5.1. The scaling parameter is the multiplicative factor describing the difference of speed of these two time scales. The main convergence result is Theorem 2. The assumption of a fast time scale for gene expression activation and mRNAs is quite common in the literature, see Bokes et al. [5] and Yvinec et al. [40]. The techniques used in these references rely on singular perturbation methods to deal with the two time scales. In our setting, a probabilistic approach is used, as will be seen, it gives precise results on the asymptotic stochastic evolution of the number of proteins.

Under this limiting regime it is shown that, asymptotically, the protein produc-
tion process can be described as a birth and death process. See Keilson [17] for
example. In state $x \in \mathbb{N}$, the birth rate is given $a/(b+x)$ for some constants $a$ and $b$. This is a contribution of the paper that, with a simple model of the autoregula-
tion, one can show that the repression mechanism follows indeed a Hill repression
scheme with an hyperbolic control, i.e. with Hill coefficient 1. The death rate is not
changed by the limiting procedure, it is proportional to $x$. Consequently, one can
get an asymptotic closed form expression of the standard deviation of the number
of proteins by using the explicit representation of the equilibrium of this birth and
death process. See Corollary 4.1. It is shown that, in this limiting regime, the
standard deviation is reduced by 30%. The corresponding results are presented in
Section 3 and Section 4 and in Appendix A. The mathematical results are obtained
via convergence theorems for sequence of Markov process, the proof of a stochastic
averaging principle and a saddle point approximation result.

SIMULATIONS. We also analyze, via simulations, autogenous regulation but also other aspects related to the regulation of protein production. This is presented in Section 5. Simulations are used mainly because of the complexity of the mathematical models of some aspects of the autogenous regulation. By using plausible biological parameters, one gets an improvement of 15% for the standard deviation of the number of proteins can be expected. This is significantly less than the performances of the limiting mathematical model studied in Section 3. The main reason seems to be that the scaling parameter is not, in some cases, sufficiently large to have a reasonable accuracy with the limit given by the convergence result of Theorem 2.

Via simulations, one also investigates the case when the regulation is not on the gene expression but on the corresponding mRNAs: a protein can block an mRNA for some time. In this situation, it could be expected that the production process is modulated more smoothly by playing on the inactivation of a fraction of the mRNAs and not on the rough on-off control of the gene expression. It is shown that the improvement is real but not that big (less than 10%). It is nevertheless remarkable that if the average life time of mRNAs is significantly increased, our experiments show that the benefit of such regulation can be of the order of more than 30% on the standard deviation of the number of proteins.

Coming back to regulation on the gene expression. Our experiments show that, despite the impact of autogenous regulation on fluctuations of the number of proteins can be limited, it has nevertheless a very interesting property. Starting with a number of proteins significantly less (or greater) than the average number of proteins at equilibrium, the autogenous regulation returns to the “correct” number of
<figure><img src="image_2.png" /><figcaption>FIGURE 1. Classical Three Stage Model for Protein Production.</figcaption></figure>

proteins much faster than the classical production process without regulation. This is a clear advantage of this mechanisms to adapt quickly when biological conditions change due to an external stress for example. See Section 5.6. This phenomenon has been observed, via experiments, in Rosenfeld et al. [31]. See also Camas et al. [6]. Finally Section 5.5 investigates the comparison of production processes with and without a feedback on the gene expression through the estimation of their respective power spectral density.

## Acknowledgments

The authors are very grateful to an anonymous referee for the important work done on an earlier version of the paper.

## 2. STOCHASTIC MODELS OF PROTEIN PRODUCTION

We present the stochastic models used to investigate the protein production process. We will use the three step model describing the activation-deactivation of the gene, the transcription phase and the translation phase. Like in most of the literature, it is assumed that the various events, like the encounter of two macromolecules, occurring within the cell have a duration with an exponential distribution. We start with the classical model used in this domain since the late 70's by Berg [3] and Rigney [27, 28]. See also Thattai and van Oudenaarden [37] and Paulsson [25].

### 2.1. The Classical Model of Protein Production.

(1) The gene is activated at rate $λ_1^+$ and inactivated at rate $λ_1^-$.

(2) If the gene is active, an mRNA is produced at rate $λ_2$. An mRNA is degraded at rate $μ_2$.

(3) Given M mRNAs at some moment, a protein is produced at rate $λ_3 M$. Each protein is degraded at rate $μ_3$.

The stochastic processes describing the protein production process are: $I(t)$ the state of the gene at time $t$ which is 0 if it is inactive and 1 otherwise. The number of mRNA at time $t$ is $M(t)$ and $P(t)$ denotes the number of proteins at that moment. The process $(I(t), M(t), P(t))$ is Markovian with state space

$$S \stackrel{\text{def.}}{=} \{0, 1\} \times \mathbb{N}^2,$$
its transition rates are given by, if $(I(t), M(t), P(t))=(i, m, p) \in S$,

$$ \left\{ \begin{array}{ll} (0, m, p) \rightarrow (1, m, p) & \text{at rate } \lambda_1^+, \qquad (1, m, p) \rightarrow (0, m, p) \text{ at rate } \lambda_1^-, \\ (i, m, p) \rightarrow (i, m+1, p) & \lambda_2 i, \qquad (i, m, p) \rightarrow (i, m-1, p) \quad \mu_2 m, \\ (i, m, p) \rightarrow (i, m, p+1) & \lambda_3 m, \qquad (i, m, p) \rightarrow (i, m, p-1) \quad \mu_3 p. \end{array} \right. $$

See Figure 1. This Markov process has a unique invariant distribution. An explicit expression of the distribution of $P$ at equilibrium is not known but, due to the linear transition rates, the moments of $P$ can be calculated recursively. In the following $(I, M, P)$ will denote random variables whose law is invariant for $(I(t), M(t), P(t))$.

**Proposition 1.** *At equilibrium, the two first moments of $P$ can be expressed by*

$$ (1) \qquad \mathbb{E}(P) = \frac{\lambda_1^+ \lambda_2 \lambda_3}{\lambda_1^+ + \lambda_1^- \mu_2 \mu_3} $$

$$ (2) \quad \mathrm{var}(P) = \mathbb{E}(P) \left( 1 + \frac{\lambda_3}{\mu_2 + \mu_3} + \frac{\lambda_1^- \lambda_2 \lambda_3 (\lambda_1^+ + \lambda_1^- + \mu_2 + \mu_3)}{(\lambda_1^+ + \lambda_1^-) (\mu_2 + \mu_3) (\lambda_1^+ + \lambda_1^- + \mu_2) (\lambda_1^+ + \lambda_1^- + \mu_3)} \right). $$

See Paulsson [25], Shahrezaei and Swain [34], Swain et al. [35] and Fromion et al. [9] for example.

**2.2. A Stochastic Model of Protein Production with Autogenous Regulation.** The regulation is done via proteins which can inactivate the gene corresponding to the protein. If there are $P$ proteins at some moment then the gene is activated at a rate proportional to $P$. Compared to the above model, only the first step changes.

(1) The inactive gene is activated at rate $\lambda_1^+$ and inactivated at rate $\lambda_1^- P$ otherwise.

See Figure 2. For the sake of simplicity, we use the same notations $\lambda_1^+$ and $\lambda_1^-$ as for the classical model of protein production instead of $\lambda_{F,1}^+$ and $\lambda_{F,1}^-$ for example. It should be noted that in our comparisons in Section 5, these quantities are not necessarily the same for these two models.

The corresponding Markov process is denoted as $(I_F(t), M_F(t), P_F(t))$, its transitions have the same rate as $(I(t), M(t), P(t))$ except for those concerning the first coordinate.

$$ \left\{ (0, m, p) \rightarrow (1, m, p) \quad \text{at rate } \lambda_1^+, \qquad (1, m, p) \rightarrow (0, m, p) \quad \text{at rate } \lambda_1^- p. \right. $$

As before, $(I_F, M_F, P_F)$ will denote random variables whose law is the invariant distribution of the Markov process $(I_F(t), M_F(t), P_F(t))$. The following proposition is the analogue of Proposition 1 for the feedback model but with unknown quantities related to the activity of the gene, $\mathbb{E}(I_F)$, and the correlation of the activity of the gene and the number of mRNAs, $\mathbb{E}(I_F M_F)$.

**Proposition 2.** *At equilibrium, the first moment of $P_F$ can be expressed by*

$$ (3) \qquad \mathbb{E}(P_F) = \mathbb{E}(I_F) \frac{\lambda_2 \lambda_3}{\mu_2 \mu_3}. $$
<figure><img src="image_3.png" /><figcaption>FIGURE 2. Three Stage Model for Protein Production with Auto-genous Regulation</figcaption></figure>

Proof. By equality of input and output for $(M(t))$ and $(P(t))$ at equilibrium, one gets the relations

$$ \lambda_2 \mathbb{E}(I_F) = \mu_2 \mathbb{E}(M_F), \quad \lambda_3 \mathbb{E}(M_F) = \mu_3 \mathbb{E}(P_F), $$

and therefore Relation (3). □

<figure><img src="image_4.png" /></figure>

It does not seem that an expression for $\mathbb{E}(I_F)$ can be obtained, the relation $\lambda_1^- \mathbb{E}(I_F P_F) = \lambda_1^+(1-\mathbb{E}(I_F))$ of equality of flows for activation/deactivation process introduces the correlation between $I_F$ and $P_F$. This is in fact the main obstacle to get more insight on the fluctuations of the number of proteins. The next section investigates a scaling where the activation/deactivation phase is much more rapid than the production process of proteins.

### 3. A SCALING ANALYSIS

It has been seen in the previous section that, for the feedback mechanism, an explicit representation of the variance of the number of proteins at equilibrium seems to be difficult to derive. In this section we use the fact that the time scale of the first two steps, activation/deactivation of the gene and production of mRNAs is more rapid than the time scale of protein production. This is illustrated by the fact that the lifetime of an mRNA is of the order of 2mn. whereas the doubling time of a bacteria is around 40mn giving a lifetime of a protein of the order of one hour. See Taniguchi et al. [36], Li and Elf [19] and Hammar et al. [13]. As will be seen, this assumption simplifies the analysis of the feedback mechanism. We will be able to get an asymptotic explicit expression for the distribution of the number of proteins at equilibrium.

A (large) scaling parameter $N$ is used to stress the difference of time scale. When there is a feedback control, an upper index $N$ is added to the variables so that the corresponding Markov process is denoted as $(X_F^N(t)) = (I_F^N(t), M_F^N(t), P_F^N(t))$ on the state space $S = \{0,1\} \times \mathbb{N}^2$. The transition rates of the Markov process are given by

$$ (4) \quad \begin{cases} (0, m, p) \to (1, m, p) \text{ at rate } \lambda_1^+ N, & (1, m, p) \to (0, m, p) \text{ at rate } \lambda_1^- Np, \\ (i, m, p) \to (i, m+1, p) & i\lambda_2 N, & (i, m, p) \to (i, m-1, p) & \mu_2 mN, \\ (i, m, p) \to (i, m, p+1) & \lambda_3 m, & (i, m, p) \to (i, m, p-1) & \mu_3 p. \end{cases} $$

The initial state is constant with $N$ given by $X_F^N(0) = (i_0, m_0, p_0) \in S$.
The aim of this section is of proving that the non-Markovian process ($P_F^N(t)$) converges in distribution to a limiting Markov process ($\overline{P}_F(t)$). As will be seen, an averaging principle, proved in the appendix, holds: locally the “fast” process ($I_F^N(t)$, $M_F^N(t)$) reaches very quickly some equilibrium depending on the current value of the “slow” variable $P_F^N(t)$. It turns out that the equilibrium of this limiting process ($\overline{P}_F(t)$) can be analyzed in detail. The proof of the averaging principle relies on stochastic calculus applied to Markov processes in the same spirit as in Papanicolau et al. [24] in a Brownian setting, see also Kurtz [18].

**Notations.** Throughout the rest of this paper, we will use the following notations $\rho_1 = \lambda_1^+/\lambda_1^-$ and, for $i = 1, 2$, $\rho_i = \lambda_i/\mu_i$.

**3.1. Scaling of the Classical Model of Protein Production.** One first states a scaling result for the classical model of protein production. The result being much simpler to prove than the corresponding result, Theorem 2, for the feedback process, its proof is skipped. One denotes by ($X^N(t)$) = ($I^N(t)$, $M^N(t)$, $P^N(t)$) the corresponding Markov process, its transition rates are the same as for feedback in Relation (4) except for deactivation:

$$ (1, m, p) \rightarrow (0, m, p) \text{ at rate } \lambda_1^- N. $$

The following result shows that, in the limit, the evolution of the number of proteins converges to the time evolution of an $M/M/\infty$ queue. See Chapter 6 of Robert [29] for example.

**Theorem 1.** If $X^N(0) = (i_0, m_0, p_0) \in S$, the sequence of processes ($P^N(t)$) converges in distribution to a birth and death process ($\overline{P}(t)$) on $\mathbb{N}$ whose respective birth and death rates ($\beta_x$) and ($\delta_x$) are given by

$$ \beta_x = \frac{\lambda_3 \rho_2 \rho_1}{\rho_1 + 1} \quad \text{and} \quad \delta_x = \mu_3 x. $$

The equilibrium distribution of ($\overline{P}(t)$) is a Poisson distribution with parameter $\rho_1\rho_2\rho_3/(1 + \rho_1)$.

*Proof.* The intuition of this result can be described quickly as follows. The processes ($I^N(t)$, $M^N(t)$) live on a much faster time scale than ($P^N(t)$) and therefore reach quickly the equilibrium. When $N$ gets large, the process ($M^N(t)$) is an $M/M/\infty$ queue with arrival rate $\lambda_2\lambda_1^+/(\lambda_1^+ + \lambda_1^-)$ and service rate $\mu_2$. See Chapter 6 of Robert [29] for example. Its equilibrium distribution is therefore Poisson with parameter $\rho_2\rho_1/(1 + \rho_1)$. The process ($P^N(t)$) can then be seen as an $M/M/\infty$ queue with arrival rate $\lambda_3\rho_2\rho_1/(1 + \rho_1)$ and service rate $\mu_3$, i.e. a birth and death process with the transition rates of the theorem. Its equilibrium is Poisson with parameter $\rho_1\rho_2\rho_3/(1 + \rho_1)$.

The proof of a corresponding result in a more complicated setting, for the production process with feedback, is done below. For this reason the proof of this result is skipped. □

**3.2. Scaling of the Production Process with Feedback.** The following theorem is the main result of this section. As in the case of the classical model of protein production, it relies on the fact that, due to the scaling, the activation/deactivation of the gene and the production of mRNAs occur on a fast time scale so that an
<figure><img src="image_5.png" /><figcaption>FIGURE 3. Feedback Model with Scaling Parameter N</figcaption></figure>

averaging principle holds. See below. Some of the technical results used to establish the following theorem are presented in the Appendix.

**Theorem 2 (Hill Repression Scheme).** If $X_F^N(0) = (i_0, m_0, p_0) \in S$, the sequence of processes $(P_F^N(t))$ converges in distribution to a birth and death process $(\bar{P}_F(t))$ on $\mathbb{N}$ whose respective birth and death rates $(\beta_x)$ and $(\delta_x)$ are given by

$$ \beta_x = \frac{\lambda_3 \rho_2 \rho_1}{\rho_1 + x} \quad \text{and} \quad \delta_x = \mu_3 x, $$

with $\rho_1 = \lambda_1^+/\lambda_1^-$ and $\rho_2 = \lambda_2/\mu_2$.

*Proof.* If $f$ is a function on $\mathbb{N}$ with finite support then

$$ V_f^N(t) \stackrel{\text{def.}}{=} f(P_F^N(t)) - f(p_0) - \int_0^t \lambda_3 M_F^N(u) \Delta^+(f)(P_F^N(u)) \, du - \int_0^t \mu_3 P_F^N(u) \Delta^-(f)(P_F^N(u)) \, du, $$

is a local martingale. See Rogers and Williams [30] for example. The operators $\Delta^+$ and $\Delta^-$ are defined as follows, for a real-valued function $f$ on $\mathbb{N}$,

$$ \Delta^+(f)(x) = f(x+1) - f(x) \quad \text{and} \quad \Delta^-(f)(x) = f(x-1) - f(x), \quad x \in \mathbb{N}. $$

With a similar method as in the proof of Assertion 1) of Lemma 1 in the appendix and by using the criterion of the modulus of continuity, see Theorem 7.2 page 81 of Billingsley [4], it is easy to show that the two processes

$$ \left( \int_0^t \lambda_3 M_F^N(u) \Delta^+(f)(P_F^N(u)) \, du \right) \quad \text{and} \quad \left( \int_0^t \mu_3 P_F^N(u) \Delta^-(f)(P_F^N(u)) \, du \right) $$

are tight. Because of the tightness of $(P_F^N(t))$ of Proposition 5 of the appendix, one can take $(N_k)$ a subsequence such that the process

$$ \left( P_F^{N_k}(t), \int_0^t \lambda_3 M_F^{N_k}(u) \Delta^+(f)(P_F^{N_k}(u)) \, du, \int_0^t \mu_3 P_F^{N_k}(u) \Delta^-(f)(P_F^{N_k}(u)) \, du \right) $$

converges in distribution.

Let $(\bar{P}_F(t))$ be a possible limit of $(P_F^{N_k}(t))$, then by continuity of the mapping

$$ (z(t)) \mapsto \left( \int_0^t z(u) \Delta^-(f)(z(u)) \, du \right) $$
on $\mathcal{D}([0, T])$ endowed with the Skorohod topology then, for the convergence in distribution

$$ \lim_{k \to +\infty} \left( P_F^{N_k}(t), \int_0^t P_F^{N_k}(u) \Delta^-(f)(P_F^{N_k}(u)) \, du \right) = \left( \bar{P}_F(t), \int_0^t \bar{P}_F(u) \Delta^-(f)(\bar{P}_F(u)) \, du \right). $$

For $t \le T$, by using the definition of $\Lambda^N$ and of $\mathcal{E}_T$ in Section A.2 of the Appendix, one has the relation

$$ \int_0^t M_F^{N_k}(u) \Delta^+(f)(P_F^{N_k}(u)) \, du = \int_{\mathcal{E}_T} m \Delta^+(f)(p) \mathbb{1}_{[0,t]}(u) \Lambda^{N_k}(dz), $$

hence, by Proposition 6 of Appendix, for the convergence in distribution

$$ \begin{aligned} & \lim_{k \to +\infty} \int_{\mathcal{E}_T} m \Delta^+(f)(p) \mathbb{1}_{[0,t]}(u) \Lambda^{N_k}(dz), \\ &= \int_0^t \sum_{p \in \mathbb{N}} \Delta^+(f)(p) \sum_{(i,m) \in \{0,1\} \times \mathbb{N}} m l_u(i, m, p) \, du, \\ &= \int_0^t \sum_{p \in \mathbb{N}} \Delta^+(f)(p) \frac{\lambda_1^+}{\lambda_1^+ + \lambda_1^- p} \frac{\lambda_2}{\mu_2} \nu_u(p) \, du \end{aligned} $$

by Relation (17) of Proposition 6 of the Appendix. By convergence of the sequence $(\Lambda^{N_k})$ this last expression can be expressed as

$$ \begin{aligned} & \left( \int_0^t \sum_{p \in \mathbb{N}} \Delta^+(f)(p) \frac{\lambda_1^+}{\lambda_1^+ + \lambda_1^- p} \nu_u(p) \, du \right) \\ &= \lim_{k \to +\infty} \left( \int_0^t \Delta^+(f)(P_F^{N_k}(u)) \frac{\lambda_1^+}{\lambda_1^+ + \lambda_1^- P_F^{N_k}(u)} \, du \right) \\ &\stackrel{\text{dist.}}{=} \left( \int_0^t \Delta^+(f)(\bar{P}_F(u)) \frac{\lambda_1^+}{\lambda_1^+ + \lambda_1^- \bar{P}_F(u)} \, du \right) \end{aligned} $$

for the convergence in distribution.

For $0 \le s \le t$, the characterization of a Markov process as the solution of a martingale problem gives the relation

$$ \mathbb{E} \left( f(P_F^N(t)) - f(P_F^N(s)) - \int_s^t \lambda_3 M_F^N(u) \Delta^+(f)(P_F^N(u)) \, du - \int_s^t \mu_3 P_F^N(u) \Delta^-(f)(P_F^N(u)) \, du \mid \mathcal{F}_s \right) = 0, $$

from which we deduce the identity

$$ \mathbb{E} \left( f(\bar{P}_F(t)) - f(\bar{P}_F(s)) - \int_s^t \lambda_3 \frac{\rho_1 \rho_2}{\rho_1 + \bar{P}_F(u)} \Delta^+(f)(\bar{P}_F(u)) - \int_s^t \mu_3 \bar{P}_F(u) \Delta^-(f)(\bar{P}_F(u)) \, du \mid \mathcal{F}_s \right) = 0. $$
See Theorem II.2.42 of Jacod and Shiryaev [15]. Consequently, a possible limit is the solution of the martingale problem associated to the birth and death process with birth rate ($\beta_x$) and death rate ($\delta_x$) and with initial state in $p_0$. One gets therefore the desired convergence in distribution of ($P_F^N(t)$). The theorem is proved. □

There exist cases where the autoregulation is not achieved by the regulated protein but by a complex of this protein, e.g by a dimer (2 copies of the protein) or a tetramer (4 copies) to cite few examples. In order to handle such cases, it is necessary to add to the gene expression model, a preliminary step describing the reaction scheme of the complex formation based on the law of mass action. In general, the dynamics involved in the reaction scheme are (very) rapid compared to the other processes of the gene expression and leads, by a singular perturbation like argument, to represent in case of deterministic model the rate of production of mRNAs as a non-linear function of protein concentration. Furthermore, when the reaction scheme possesses suitable properties, a Hill like repression function could also be obtained. See Weiss [38] for details. In the stochastic context, that leads to introduce a suitable scaling factor in the dynamics of the complex formation and to extend the previous derivation in the previous theorem to Hill functions, $x \mapsto a/(b+x^n)$, with order n greater than 1.

The next section analyzes, in this limiting regime, the fluctuations of the number of proteins at equilibrium.

## 4. FLUCTUATIONS OF THE NUMBER OF PROTEINS

This section is devoted to the analysis of the equilibrium of the asymptotic process ($\bar{P}_F(t)$) of Theorem 2 describing the evolution of the number of proteins with feedback. We start with a classical result for birth and death processes.

**Proposition 3.** The invariant distribution $\pi_F$ of the birth and death process ($\bar{P}_F(t)$) of Theorem 2 is given by

$$ \pi_F(x) = \frac{1}{Z} \frac{(\rho_2 \rho_3)^x}{x!} \prod_{i=0}^{x-1} \frac{\rho_1}{\rho_1 + i}, \quad x \in \mathbb{N}, $$

where $\rho_1 = \lambda_1^+/\lambda_1^-$, $\rho_i = \lambda_i/\mu_i$ for $i = 1, 2$ and $Z$ is the normalization constant.

The expression of $\pi_f$ is explicit but with a normalization constant which is not simple. The constant $Z$ can be expressed in terms of hypergeometric functions. See Abramowitz and Stegun [1] for example. Even if we can get a numerical evaluation of the average and of the variance of $\pi_F$, it is much more awkward to get some insight on the dependence of these quantities with respect to some of the parameters like $\rho_2$ or $\rho_3$ for example. In the following we give an asymptotic description of the ratio of the variance and the mean of the number of proteins at equilibrium when the value of the quantity $\rho_1\rho_2\rho_3$ is large. In a biological context the numerical value of this parameter is not always large but this limit results sheds some light on the qualitative behavior of the auto-regulation mechanism. See Corollary 4.1 for example. A Laplace method is in particular used to investigate the asymptotic behavior of the first two moments of $\pi_F$.

Theorem 1 shows that the distribution of the process ($P(t)$) at equilibrium is Poisson with parameter $E(P(t)) = x_\rho = \rho_1\rho_2\rho_3/(1 + \rho_1)$. In particular, one has the relation $\text{var}(P(t)) = E(P(t))$. In the rest of this section, we will be interested in the corresponding quantity for the feedback process.
For $\eta > 0$ and $\rho > 0$, denote by $\nu_\rho$ the probability distribution on $\mathbb{N}$ defined by

$$ (5) \qquad \nu_{\rho,\eta}(k) = \frac{1}{Z_{\rho}} \frac{\rho^k}{k!} \prod_{i=1}^{k} \frac{1}{\eta + i} = \frac{1}{Z_{\rho}} \exp \left( \sum_{i=1}^{k} \log \left( \frac{x}{i(\eta + i)} \right) \right), $$

where $Z_\rho$ is the normalization constant. It is easily seen that $\pi_F$ is $\nu_{\rho,\eta}$ with $\rho = \rho_1\rho_2\rho_3$ and $\eta = \rho_1 - 1$.

**Proposition 4.** If, for $\rho > 0$ and $\eta > -1$, $A_\rho$ is a random variable with distribution $\nu_{\rho,\eta}$ defined by Equation (5), then for the convergence in distribution

$$ \lim_{\rho \to +\infty} \frac{A_\rho - a_\rho}{\sqrt{a_\rho}} = \mathcal{N}(0, 1/\sqrt{2}), $$

where $a_\rho = (\sqrt{\eta^2 + 4\rho} - \eta)/2$ and $\mathcal{N}(0, 1/\sqrt{2})$, is a centered Gaussian random variable with variance $1/2$. In particular, for the convergence in distribution,

$$ \lim_{\rho \to +\infty} \frac{A_\rho}{\sqrt{\rho}} = 1. $$

*Proof.* If $\phi$ is a bounded function on $\mathbb{R}$, denote

$$ \Delta_{\rho}(\phi) \stackrel{\text{def.}}{=} \frac{1}{\sqrt{a_{\rho}}} \sum_{k=0}^{+\infty} \phi\left(\frac{k - \lceil a_{\rho} \rceil}{\sqrt{a_{\rho}}}\right) \exp\left(\sum_{i=\lceil a_{\rho} \rceil}^{k} \log\left(\frac{\rho}{i(\eta + i)}\right)\right), $$

the definition of $\nu_{\rho,\eta}$ gives that

$$ (6) \qquad \mathbb{E} \left( \phi \left( \frac{A_{\rho} - \lceil a_{\rho} \rceil}{\sqrt{a_{\rho}}} \right) \right) = \frac{\Delta_{\rho}(\phi)}{\Delta_{\rho}(1)} $$

Fix $\phi$ some continuous function with compact support on $[-K_0, K_0]$ for some $K_0>0$. Since $a_\rho$ is the solution of the equation $a_\rho(\eta+a_\rho)=\rho$, a change of variable gives the relation

$$ \Delta_{\rho}(\phi) = \frac{1}{\sqrt{a_{\rho}}} \sum_{k=-\lceil K_0 \sqrt{a_{\rho}} \rceil}^{\lceil K_0 \sqrt{a_{\rho}} \rceil} \phi\left(\frac{k}{\sqrt{a_{\rho}}}\right) \exp\left(\sum_{i=0}^{k} \log\left(\frac{a_{\rho}(\eta + a_{\rho})}{(i + \lceil a_{\rho} \rceil)(\eta + \lceil a_{\rho} \rceil + i)}\right)\right). $$

The uniform estimation

$$ \sum_{i=0}^{k} \log \left( \frac{a_{\rho}(\eta + a_{\rho})}{(i + \lceil a_{\rho} \rceil)(\eta + \lceil a_{\rho} \rceil + i)} \right) = \int_{0}^{k} \log \left( \frac{a_{\rho}(\eta + a_{\rho})}{(u + \lceil a_{\rho} \rceil)(\eta + \lceil a_{\rho} \rceil + u)} \right) du + O\left(\frac{1}{\sqrt{a_{\rho}}}\right) $$
for all $k \in \mathbb{Z}$ with $|k| \le K_0 \sqrt{a_\rho}$ and the fact that $\phi$ has a compact support give that the quantity $\Delta_\rho(\phi)$ is equivalent to

$$ \frac{1}{\sqrt{a_\rho}} \sum_{k=-\lfloor K_0 \sqrt{a_\rho} \rfloor}^{\lfloor K_0 \sqrt{a_\rho} \rfloor} \phi\left(\frac{k}{\sqrt{a_\rho}}\right) \exp\left(\int_0^k \log\left(\frac{a_\rho(\eta+a_\rho)}{(u+\lfloor a_\rho \rfloor)(\eta+\lfloor a_\rho \rfloor+u)}\right) du\right) \\ = \frac{1}{\sqrt{a_\rho}} \sum_{k=-\lfloor K_0 \sqrt{a_\rho} \rfloor}^{\lfloor K_0 \sqrt{a_\rho} \rfloor} \phi\left(\frac{k}{\sqrt{a_\rho}}\right) \\ \times \exp\left(\int_0^{k/\sqrt{a_\rho}} \sqrt{a_\rho} \log\left(\frac{a_\rho(\eta+a_\rho)}{(u\sqrt{a_\rho}+\lfloor a_\rho \rfloor)(\eta+\lfloor a_\rho \rfloor+u\sqrt{a_\rho})}\right) du\right). $$

Again, with the uniform estimation

$$ \sqrt{a_\rho} \log \left( \frac{a_\rho(\eta + a_\rho)}{(u\sqrt{a_\rho} + \lfloor a_\rho \rfloor)(\eta + \lfloor a_\rho \rfloor + u\sqrt{a_\rho})} \right) = -2u + O\left(\frac{1}{\sqrt{a_\rho}}\right), $$

for $u$ in some fixed finite interval, one gets that

$$ \Delta_\rho(\phi) \sim \frac{1}{\sqrt{a_\rho}} \sum_{k=-\lfloor K_0 \sqrt{a_\rho} \rfloor}^{\lfloor K_0 \sqrt{a_\rho} \rfloor} \phi\left(\frac{k}{\sqrt{a_\rho}}\right) \exp\left(-2 \int_0^{k/\sqrt{a_\rho}} u \, du\right) \sim \int_{-\infty}^{+\infty} \phi(v) e^{-v^2} \, dv. $$

With similar estimations for $\Delta_\rho(1)$ (which imply in fact the tightness of the random variables $(A_\rho - \lfloor a_\rho \rfloor)/\sqrt{a_\rho}$) and Relation (6), the proposition is proved. $\square$

**Corollary 4.1** (Asymptotic Number of Proteins with Regulation). If $\bar{P}_F$ is a random variable with distribution $\pi_F$ then, for the convergence in distribution

$$ (7) \qquad \lim_{\rho_2 \rho_3 \to +\infty} \frac{\mathbb{E}(\bar{P}_F)}{\sqrt{\rho_1 \rho_2 \rho_3}} = 1 \quad \text{and} \quad \lim_{\rho_2 \rho_3 \to +\infty} \frac{\text{var}(\bar{P}_F)}{\mathbb{E}(\bar{P}_F)} = \frac{1}{2}. $$

Furthermore, for the convergence in distribution,

$$ \lim_{\rho_2 \rho_3 \to +\infty} \frac{\bar{P}_F - a_\rho}{\sqrt{a_\rho}} = \mathcal{N}\left(0, 1/\sqrt{2}\right), $$

$$ \text{where } a_\rho = \left( \sqrt{(\rho_1 - 1)^2 + 4\rho_1\rho_2\rho_3} - \rho_1 + 1 \right) / 2. $$

The equivalent of Relation (7) for the scaling of the classical model of protein production is

$$ \mathbb{E}(P) = \frac{\rho_1}{1 + \rho_1} \rho_2 \rho_3 \text{ and } \frac{\text{var}(P)}{\mathbb{E}(P)} = 1. $$

by Theorem 1. it shows that a feedback mechanism reduces the variance of the number of proteins in this limiting regime by a factor 2 for the ratio of the second moment and the first moment.

## 5. DISCUSSION

In this section, other aspects of regulation of protein production are discussed via simulations in a plausible biological context whose parameters are going to be defined. Simulation follows the models in Section 2.2 and simulates the variables $I_F$, $M_F$, and $P_F$, not their scaling limits.
**5.1. Numerical Values of Biological Parameters.** For the model with feedback, there are six parameters to determine. By using the literature one can estimate the common orders of magnitude of these parameters in a biological context. We therefore propose a set of parameters corresponding to an “ordinary” gene.

(1) Gene regulation. The parameter $\lambda_1^-$ gives the rate at which a given protein reaches its own promoter. It has been shown that this motion combines a three-dimensional diffusion in the cytoplasm and one-dimensional sliding along the DNA, see Halford [12].

Experiments on the lac repressor, using live-cell single-molecule imaging techniques, show that this time is of the order of 5 min, see Li and Elf [20] and Hammar et al. [13]. For this reason we will take $\lambda_1^- = 3.3 \times 10^{-3} s^{-1}$.

The parameter $\lambda_1^+$ can be quite variable, depending on the affinity of the protein to the DNA sequence, we set $\lambda_1^+ = 1 s^{-1}$.

(2) mRNAs. The lifetime of an mRNA is $\mu_2^{-1} \approx 4$ min, see Taniguchi et al. [36]. When the gene expression is always active (corresponding to the case where our variable $I$ remains equals to 1), there is an average of 2 messengers, that is to say $\lambda_2^{-1} = \mu_2^{-1}/m = 120 s$ which gives $\lambda_2 = 8.3 \times 10^{-3} s^{-1}$.

(3) Proteins. A doubling time for the cell of $t_{1/2} \approx 40$ min gives a protein decay of around one hour. For this reason one takes $\mu_3 = \log 2/t_{1/2} = 2.8 \times 10^{-4} s^{-1}$ for the rate of protein decay. It is assumed that a give type of protein that is produced in $p = 300$ copies when the gene expression is always active. From one messenger, a protein should be produced in a duration of time of the order of $\lambda_3^{-1} = m \times \mu_3^{-1}/p$ which gives $\lambda_3 = 4 \times 10^{-2} s^{-1}$.

These parameters may correspond to an “ordinary bacterial” gene: in a E. Coli genome of 4300 genes, there are around $3.6 \times 10^6$ proteins and $1.4 \times 10^3$ mRNAs per gene, see Table 1 of Chapter 3 of Neidhardt [23], the number of messengers and proteins is of the order of magnitude of our numerical estimation of the parameters.

**5.2. Impact of Autogenous Regulation on Gene Expression.** We have compared two mechanisms: the classical model without regulation and the autogenous regulation process. The mean number of proteins is the same as well as the mean number of mRNAs produced $\mathbb{E}(M) = \mathbb{E}(M_F)$. Parameters $\lambda_1^+$ and $\lambda_1^-$ are adapted in the classical model to fulfill these conditions. The other parameters are as defined in the previous section.

The comparison is shown in Figure 4. The mean number of proteins is 178, as can be seen that the curve for the autogenous regulation is slightly more concentrated around the mean but not that much. The values of the corresponding standard deviations are not really different $\sqrt{\text{var}(P)}=42.2$ and $\sqrt{\text{var}(P_F)}=35.8$. The impact of the autogenous regulation on the variability of the number of proteins is non-trivial but not really spectacular for the set of parameters associated to a “typical” gene. This is significantly less than the performances of the limiting mathematical model studied in Section 3. The main reason seems to be that the scaling parameter is not, in some cases, sufficiently large to have a reasonable accuracy with the limit given by the convergence result of Theorem 2.

**5.3. The Limiting Scaling Regime as a Lower Bound.** Roughly speaking, Theorem 1 and Corollary 4.1 give that for $N$ and $\rho_2\rho_3$ large, then the ratio $\text{var}(P_F^N)/\mathbb{E}(P_F^N)$ converges to $1/2$. In Figure 5, one considered a simulation with
<figure><img src="image_6.png" /><figcaption>FIGURE 4. Simulations: Protein distribution with and without autogenous regulation with a fixed mean number of proteins of 178.</figcaption></figure>

fixed product $\rho_2\rho_3 = 71.43$ with $N$ varying. The interesting feature is that the ratio is decreasing with $N$, this suggests that the variance of the limit of the scaling procedure should provide a lower bound for the variance of the real model. We have not been able to show rigorously this phenomenon. For $N = 250$, the value of the ratio $\text{var}(P_F^N)/\text{E}(P_F^N) = .7964$ which is quite far from its limiting value 1/2 given by Corollary 4.1. This can be explained by the fact that the quantities $N$ and $\rho_2\rho_3$ are not very large.

5.4. **Regulation of the Production Process on mRNAs.** The regulation on the gene has the effect of an ON/OFF mechanism. When the gene is active, it is producing mRNAs at full speed and no mRNA is produced when it is inactive. This suggests that the production of proteins follows roughly the same pattern: steady production rate at some instants and little is produced otherwise. This scheme can consequently increase the variability of the production process of proteins. A possible idea to reduce the variance due to the activation/inactivation of the gene is to transfer the activation/inactivation process at the level the mRNAs. This possibility is investigated in this section. Each mRNA can be inactivated by a protein at rate $\lambda_2^-$, in this state it cannot produce proteins. An inactivated mRNAs becomes active at rate $\lambda_2^+$. In this way the production process can, hopefully, be modulated more smoothly by playing on the inactivation of a fraction of the mRNAs. In this way at time $t$, if the number of active [resp. inactive] mRNAs is $M(t)$ [resp. $M^*(t)$], the process $(M(t), M^*(t), P(t))$ is Markov with transition
<figure><img src="image_7.png" /><figcaption>FIGURE 5. Simulations: Evolution of the ratio var($P_F^N$)/E($P_F^N$) as a function of N.</figcaption></figure>

rates, for $(m, m^*, p) \in \mathbb{N}^3$,

$$ \begin{cases} (m, m^*, p) \to (m + 1, m^*, p) \text{ at rate } \lambda_2, \\ (m, m^*, p) \to (m - 1, m^* + 1, p) \text{ at rate } \lambda_2^+ mp, \\ (m, m^*, p) \to (m + 1, m^* - 1, p) \text{ at rate } \lambda_2^+ m^*, \end{cases} $$

the other transitions are as before, active of inactive mRNAs die at rate $\mu_2$ and proteins are produced at rate $\lambda_3 m$ and die at rate $\mu_3$.

To compare the two regulation processes, either on the gene or on mRNAs, simulations have been done with the following constraints: the average number of proteins is fixed around 1400. To have a fair comparison, we add the constraint that the number of mRNAs produced should be the same in all simulations. The numerical values have been estimated by using similar methods as in section 3 but for this setting. Experiment (3) considers the case of an average lifetime of an mRNA of 40mn, if this is far from a “normal” biological setting, as it will be seen, this scenario has the advantage of stressing the importance of this parameter in this configuration.

## Numerical Values of Parameters.

(1) Regulation on the gene.

<table><thead><tr><th>λ<sub>1</sub><sup>+</sup></th><th>λ<sub>1</sub><sup>-</sup></th><th>λ<sub>2</sub></th><th>μ<sub>2</sub></th><th>λ<sub>3</sub></th><th>μ<sub>3</sub></th></tr></thead><tbody><tr><td>0.21"</td><td>5'</td><td>12"</td><td>4'</td><td>25"</td><td>1h.</td></tr></tbody></table>

(2) Regulation on mRNAs (I).

For this experiment, the expected lifetime of an mRNA is twice the corresponding value of case (1).

<table><thead><tr><th>λ<sub>2</sub></th><th>λ<sub>2</sub><sup>+</sup></th><th>λ<sub>2</sub><sup>-</sup></th><th>μ<sub>2</sub></th><th>λ<sub>3</sub></th><th>μ<sub>3</sub></th></tr></thead><tbody><tr><td>23"</td><td>2"</td><td>45'</td><td>8'</td><td>25"</td><td>1h.</td></tr></tbody></table>
(3) Regulation on mRNAs (II).

For this second experiment on the regulation of mRNAs, the expected lifetime of an mRNA is 10 times than in case (1).

<table><thead><tr><th>λ<sub>2</sub></th><th>λ<sub>2</sub><sup>+</sup></th><th>λ<sub>2</sub><sup>-</sup></th><th>μ<sub>2</sub></th><th>λ<sub>3</sub></th><th>μ<sub>3</sub></th></tr></thead><tbody><tr><td>23.8"</td><td>2"</td><td>45'</td><td>40'</td><td>25"</td><td>1h.</td></tr></tbody></table>

**Results of the Experiments.** Table 1 shows that the mean number of mRNAs produced per unit of time is essentially the same in all experiments as well as the mean number of active mRNAs. It should be noted the impact of regulation on mRNAs for the standard deviation of the number of proteins when the mean life time is 8mn is not really significant (10% gain) than the regulation on the gene. When the mean lifetime is 40mn the improvement, 36%, of the standard deviation becomes significant, showing that in this case the production process is “smoothed” by this mechanism. The three distributions of the number of proteins of these experiments are presented in Figure 6.

<table><caption>TABLE 1. Comparison of Regulation Processes on Gene or on mRNAs with Different Lifetimes</caption><thead><tr><th>Regulation on</th><th>Gene</th><th>mRNAs/8mn</th><th>mRNAs/40mn</th></tr></thead><tbody><tr><td>Mean Nb of mRNAs</td><td>10.33</td><td>19.74</td><td>99.04</td></tr><tr><td>Mean Nb of Active mRNAs</td><td>10.33</td><td>9.77</td><td>9.81</td></tr><tr><td>Mean Nb of Proteins</td><td>1403.63</td><td>1400.29</td><td>1403.36</td></tr><tr><td>Standard Deviation of Nb of Proteins</td><td>92.66</td><td>84.22</td><td>59.04</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>FIGURE 6. Simulations: Probability Distribution of the Number of Proteins with Regulation on Gene or on mRNAs, $\mu_2^{-1}$ is the average lifetime of an mRNA. The average number of proteins is 1400.</figcaption></figure>
5.5. **Impact of Feedback on Frequency.** In this section, we study the nature of the fluctuations of the number of proteins at equilibrium from the point of view of signal processing or automatic control. The aim of a feedback is often of changing the nature of the signal, attenuating disturbances by reducing, for instance, high frequencies. In these cases, spectral analysis gives a characterization of the nature of changes.

By analogy, we consider our model as a system that has to achieve a command (the production of a given mean number of proteins) and where the resulting signal $P(t)$ or ($P_F(t)$ is altered by some noise. In this framework, one can study if the effect of the feedback has an impact on the signal, by rejection of some frequency ranges.

To do so, consider the signals ($P(t)$) and ($P_F(t)$) of two simulations with or without autogenous regulation. The analysis of these signals is done by estimating the power spectral density, that describes the spectral characteristics of stochastic process. We estimate the power spectral density for each signal, using classical estimator of smoothed periodogram. See George et al. [10] and Chapter 10 of Miller et al. [22] for example.

The result is shown in Figure 7. Both spectra seem to represent a low-pass filter with a cut off frequency in the order of magnitude of the dilution factor $\mu_3 = 2.8 \times 10^{-4} s^{-1}$. The two power spectral densities do not seem to exhibit significant differences. The feedback has therefore no noticeable effect in terms of reduction of frequency disturbances.

<figure><img src="image_9.png" /><figcaption>FIGURE 7. Power spectral density estimation of signals with and without regulation</figcaption></figure>

5.6. **Versatility of the Protein Production Process.** This section is devoted to the impact of autogenous regulation on another aspect of protein production. Up to now, we have considered the production process of proteins at equilibrium, by assuming that the production rate of a given protein has to be fixed. It may happen nevertheless that, due to an external stress, such as antibiotics, DNA damage by UV, see Camas et al. [6], or nutriment absorption, see Schleif [33], the cell has to
change rapidly its production rate to quickly produce a large amount of proteins for example. The affinity of the transcription factor for the promoter of the gene can be adapted for that purpose. Conversely, when the external stress disappears, the production of the protein has to be quickly reduced to minimize the consumption of resources.

We consider the situation when the two production processes, with and without autogenous regulation, give the same average output of proteins at equilibrium. Two cases are investigated: when the initial number of proteins is below the value equilibrium, see Figure 8, or above this value, see Figure 9. As it can be seen, the autogenous production process converges more rapidly to equilibrium in both cases. Our simulations show that when the initial value is 290, the autogenous production process is 40% faster than the process without feedback to reach the level 1300 (the equilibrium is at 1400 in this case). A similar result holds in the other case.

These interesting properties are related to the modulation of the gene activity. In the experiment of Figure 8, for the autogenous process the rate of activity of the gene is of the order of 50% at the beginning and it is only of the order of 0.1 later at equilibrium. Without regulation this rate is constant throughout the simulation. This explains the “fast start” of the autogenous process. An analogous explanation holds for the experiment of Figure 8, in the autogenous process. The gene is rapidly switched off due to the large number of proteins, thereby decreasing rapidly the number of proteins. This is consistent with experiments described in Camas et al. [6] and especially Rosenfeld et al. [31] where the improvement has been estimated at 80% in some cases.

<figure><img src="image_10.png" /><figcaption>FIGURE 8. Simulations: Evolution of the Mean Number of Proteins: Initial Point at 290, equilibrium at 1400. Time scale: ×10⁴ sec.</figcaption></figure>

REFERENCES

[1] Milton Abramowitz and Irene A. Stegun, *Handbook of mathematical functions with formulas, graphs, and mathematical tables*, National Bureau of Standards Applied Mathematics Series,
<figure><img src="image_11.png" /><figcaption>FIGURE 9. Simulations: Evolution of the Mean Number of Proteins: Initial Point at 1400, equilibrium at 290. Time scale: ×10⁴ sec.</figcaption></figure>

vol. 55, For sale by the Superintendent of Documents, U.S. Government Printing Office, Washington, D.C., 1964. MR 0167642 (29 #4914)

[2] Attila Becskei and Luis Serrano, *Engineering stability in gene networks by autoregulation*, Nature **405** (2000), 590–593.

[3] O. G. Berg, *A model for the statistical fluctuations of protein numbers in a microbial population*, Journal of theoretical biology **71** (1978), no. 4, 587–603.

[4] P. Billingsley, *Convergence of probability measures*, second ed., Wiley Series in Probability and Statistics: Probability and Statistics, John Wiley & Sons Inc., New York, 1999, A Wiley-Interscience Publication.

[5] Pavol Bokes, John R. King, A.T. Wood, and Matthew Loose, *Multiscale stochastic modelling of gene expression*, Journal of Mathematical Biology **65** (2012), no. 3, 493–520 (English).

[6] Francisco M. Camas, Jesús Blázquez, and Juan F. Poyatos, *Autogenous and nonautogenous control of response in a genetic network*, Proceedings of the National Academy of Sciences **103** (2006), no. 34, 12718–12723.

[7] Donald A. Dawson, *Measure-valued Markov processes*, École d'Été de Probabilités de Saint-Flour XXI—1991, Lecture Notes in Math., vol. 1541, Springer, Berlin, 1993, pp. 1–260.

[8] T. Fournier, J.P. Gabriel, C. Mazza, J. Pasquier, J.L. Galbete, and N. Mermod, *Steady-state expression of self-regulated genes*, Bioinformatics **23** (2007), no. 23, 3185–3192.

[9] V. Fromion, E. Leoncini, and P. Robert, *Stochastic gene expression in cells: A point process approach*, SIAM Journal on Applied Mathematics **73** (2013), no. 1, 195–211.

[10] William K. Jr. George, Paul D. Beuther, and John L. Lumley, *Processing of Random Signals*, Proceedings of the Dynamic Flow Conference 1978 on Dynamic Measurements in Unsteady Flows (Bengt Wulff Hansen M.Sc, ed.), Springer Netherlands, 1978, pp. 757–800 (en).

[11] Robert F. Goldberger, *Strategies of genetic regulation in prokaryotes*, Biological Regulation and Development (Robert F. Goldberger, ed.), Biological Regulation and Development, vol. 1, Springer US, 1979, pp. 1–18 (English).

[12] Stephen E. Halford, *An end to 40 years of mistakes in DNA-protein association kinetics?*, Biochemical Society Transactions **37** (2009), no. 2, 343 (en).

[13] Petter Hammar, Prune Leroy, Anel Mahmutovic, Erik G. Marklund, Otto G. Berg, and Johan Elf, *The lac repressor displays facilitated diffusion in living cells*, Science **336** (2012), no. 6088, 1595–1598.

[14] J. E. M. Hornos, D. Schultz, G. C. P. Innocentini, J. Wang, A. M. Walczak, J. N. Onuchic, and P. G. Wolynes, *Self-regulating gene: An exact solution*, Phys. Rev. E **72** (2005), 051907.
[15] Jean Jacod and Albert N. Shiryaev, *Limit theorems for stochastic processes*, Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of Mathematical Sciences], vol. 288, Springer-Verlag, Berlin, 1987.

[16] Magdalena Kaczanowska and Monica Rydén-Aulin, *Ribosome Biogenesis and the Translation Process in Escherichia coli*, Microbiology and Molecular Biology Reviews **71** (2007), no. 3, 477–494 (en).

[17] J. Keilson, *Markov chains models-rarity and exponentiality*, Applied Mathematical sciences, vol. 28, Springer Verlag, New York, 1979.

[18] T.G. Kurtz, *Averaging for martingale problems and stochastic approximation*, Applied Stochastic Analysis, US-French Workshop, Lecture notes in Control and Information sciences, vol. 177, Springer Verlag, 1992, pp. 186–209.

[19] Gene-Wei Li and Johan Elf, *Single molecule approaches to transcription factor kinetics in living cells*, FEBS Letters **583** (2009), no. 24, 3979–3983.

[20] __________, *Single molecule approaches to transcription factor kinetics in living cells*, FEBS Letters **583** (2009), no. 24, 3979–3983.

[21] S. Maloy and V. Stewart, *Autogenous regulation of gene expression*, Journal of Bacteriology **175** (1993), no. 2, 307–316.

[22] Scott Miller and Donald Childers, *Probability and Random Processes: With Applications to Signal Processing and Communications*, Academic Press, January 2012 (en).

[23] Frederick C. Neidhardt (ed.), *Escherichia coli and salmonella: Cellular and molecular biology*, 2 edition ed., ASM Press, Washington, D.C, May 1996 (English).

[24] G. C. Papanicolaou, D. Stroock, and S. R. S. Varadhan, *Martingale approach to some limit theorems*, Papers from the Duke Turbulence Conference (Duke Univ., Durham, N.C., 1976), Paper No. 6, Duke Univ., Durham, N.C., 1977, pp. ii+120 pp. Duke Univ. Math. Ser., Vol. III.

[25] J. Paulsson, *Models of stochastic gene expression*, Physics of Life Reviews **2** (2005), no. 2, 157–175.

[26] A. Raj and A. van Oudenaarden, *Nature, nurture or chance: Stochastic gene expression and its consequences*, Cell **135** (2008), no. 2, 216–226.

[27] D. Rigney, *Stochastic model of constitutive protein levels in growing and dividing bacterial cells*, Journal of Theoretical Biology **76** (1979), no. 4, 453–480.

[28] D. Rigney and W. Schieve, *Stochastic model of linear, continuous protein synthesis in bacterial populations*, Journal of Theoretical Biology **69** (1977), no. 4, 761–766.

[29] Philippe Robert, *Stochastic networks and queues*, Stochastic Modelling and Applied Probability Series, vol. 52, Springer, New-York, June 2003.

[30] L. C. G. Rogers and David Williams, *Diffusions, Markov processes, and martingales. Vol. 2: Itô calculus*, John Wiley & Sons Inc., New York, 1987.

[31] Nitzan Rosenfeld, Michael B Elowitz, and Uri Alon, *Negative autoregulation speeds the response times of transcription networks*, Journal of Molecular Biology **323** (2002), no. 5, 785–793.

[32] Walter Rudin, *Real and complex analysis*, third ed., McGraw-Hill Book Co., New York, 1987.

[33] R Schleif, *Regulation of the L-arabinose operon of Escherichia coli*, Trends in genetics: TIG **16** (2000), no. 12, 559–565.

[34] Vahid Shahrezaei and Peter S. Swain, *Analytical distributions for stochastic gene expression*, Proceedings of the National Academy of Sciences **105** (2008), no. 45, 17256–17261.

[35] Peter S. Swain, Michael B. Elowitz, and Eric D. Siggia, *Intrinsic and extrinsic contributions to stochasticity in gene expression*, Proceedings of the National Academy of Sciences **99** (2002), no. 20, 12795–12800.

[36] Yuichi Taniguchi, Paul J. Choi, Gene-Wei Li, Huiyi Chen, Mohan Babu, Jeremy Hearn, Andrew Emili, and X. Sunney Xie, *Quantifying E. coli proteome and transcriptome with single-molecule sensitivity in single cells*, Science **329** (2010), no. 5991, 533–538.

[37] M. Thattai and A. van Oudenaarden, *Intrinsic noise in gene regulatory networks*, Proceedings of the National Academy of Sciences **98** (2001), no. 15, 8614–8619.

[38] J. N. Weiss, *The Hill equation revisited: uses and misuses*, The FASEB Journal **11** (1997), no. 11, 835–841.

[39] N. Yildirim and M. C. Mackey, *Feedback regulation in the lactose operon: A mathematical modeling study and comparison with experimental data*, Biophysical Journal **84** (2003), no. 5, 2841–2851.
[40] Romain Yvinec, Changjing Zhuge, Jinzhi Lei, and Michael C. Mackey, *Adiabatic reduction of a model of stochastic gene expression with jump Markov process*, Journal of Mathematical Biology **68** (2014), no. 5, 1051–1070 (English).

## APPENDIX A. CONVERGENCE RESULTS

We first introduce some notations that will be used throughout this section.

A.1. **Evolution Equations.** We will use the Skorohod's topology for convergence in distribution in the space $D([0, T], \mathbb{R}_+)$ of càdlàg processes. See Chapter 3 of Billingsley [4] for example. To simplify the presentation, all our processes will be defined on the same probability space in the following way.

Let $N_i^+$, $N_i^-$, $i=1, 2, 3$ be independent Poisson processes on $\mathbb{R}_+^2$ with rate 1 defined on a probability space $(\Omega, \mathcal{F}, \mathbb{P})$. If $A \in \mathcal{B}(\mathbb{R}_+^2)$ is a Borelian subset of $\mathbb{R}_+^2$ and $(i, c) \in \{1, 2, 3\} \times \{+, -\}$, $N_i^c(A)$ denotes the number of points of the process $N_i^c$ in the subset $A$. For $t \ge 0$, one denotes by $\mathcal{F}_t$ the $\sigma$-field generated by the random variables

$$ N_i^c(B \times [0, t]) \text{ for } B \in \mathcal{B}(\mathbb{R}_+) \text{ and } (i, c) \in \{1, 2, 3\} \times \{+, -\}. $$

It is easily seen that the process $(X_F^N(t))$ has the same distribution as the solution of the following stochastic differential equations (SDE)

$$ (8) \quad dI_F^N(t) = 1_{\{I_F^N(t-)=0\}} N_1^+([0, \lambda_1^+ N] \times [dt]) \\ \qquad - 1_{\{I_F^N(t-)=1\}} N_1^-([0, \lambda_1^- N P_F^N(t-) \times [dt]) $$

$$ (9) \quad dM_F^N(t) = 1_{\{I_F^N(t-)=1\}} N_2^+([0, \lambda_2 N] \times [dt]) - N_2^-([0, \mu_2 N M_F^N(t-) \times [dt])) $$

$$ (10) \quad dP_F^N(t) = N_3^+([0, \lambda_3 M_F^N(t-) \times [dt])) - N_3^-([0, \mu_3 P_F^N(t-) \times [dt])) $$

with the same initial condition. For any $N \ge 1$, $(X_F^N(t))$ is a Markov process adapted to the filtration $(\mathcal{F}_t)$. These SDE can be rewritten as, for some function $f$ with finite support on $S$,

$$ \begin{aligned} (11) \quad f(X_F^N(t)) = f(X_F^N(0)) &+ \int_0^t \lambda_1^+ N(1 - I_F^N(u)) \Delta_1(f)(X_F^N(u)) \, du \\ &+ \int_0^t \lambda_1^- N P_F^N(u) I_F^N(u) \Delta_1(f)(X_F^N(u)) \, du \\ &+ \int_0^t \lambda_2 N I_F^N(u) \Delta_2^+(f)(X_F^N(u)) \, du + \int_0^t \mu_2 N M_F^N(u) \Delta_2^-(f)(X_F^N(u)) \, du \\ &+ \int_0^t \lambda_3 N M_F^N(u) \Delta_3^+(f)(X_F^N(u)) \, du + \int_0^t \mu_3 P_F^N(u) \Delta_3^-(f)(X_F^N(u)) \, du + W_f^N(t) \end{aligned} $$

where, for $x = (i, m, p) \in S$, the operators $\Delta^{+/-}$ are defined by

$$ \left\{ \begin{aligned} \Delta_1(f)(x) &= f(1-i, m, p) - f(x) \\ \Delta_2^+(f)(x) &= f(i, m+1, p) - f(x), \quad \Delta_2^-(f)(x) = f(i, m-1, p) - f(x) \\ \Delta_3^+(f)(x) &= f(i, m, p+1) - f(x), \quad \Delta_3^-(f)(x) = f(i, m, p-1) - f(x), \end{aligned} \right. $$
and ($W_f^N(t)$) is a local martingale whose previsible increasing process is given by

$$ (12) \quad \begin{aligned} \langle W_f^N \rangle (t) = & \int_0^t [\lambda_1^+ N(1-I_F^N(u)) + \lambda_1^- N P_F^N(u) I_F^N(u)] [\Delta_1(f)(X_F^N(u))]^2 \, du \\ & + \int_0^t \lambda_2 N I_F^N(u) [\Delta_2^+(f)(X_F^N(u))]^2 \, du + \int_0^t \mu_2 N M_F^N(u) [\Delta_2^-(f)(X_F^N(u))]^2 \, du \\ & + \int_0^t \lambda_3 M_F^N(u) [\Delta_3^+(f)(X_F^N(u))]^2 \, du + \int_0^t \mu_3 P_F^N(u) [\Delta_3^-(f)(X_F^N(u))]^2 \, du. \end{aligned} $$

See Rogers and Williams [30] for example.

**Definition 1.** Let $(\overline{M}^N(t), \overline{P}^N(t))$ be the Markov process with transition rates given by

$$ (13) \quad \begin{cases} (m, p) \to (m+1, p) \text{ at rate } \lambda_2 N, & (m, p) \to (m-1, p) \quad \mu_2 m N, \\ (m, p) \to (m, p+1) \quad \lambda_3 m, & (m, p) \to (m, p-1) \quad \mu_3 p \end{cases} $$

and initial state $(\overline{M}^N(0), \overline{P}^N(0)) = (m_0, p_0)$.

The process $(\overline{M}^N(t), \overline{P}^N(t))$ is simply the analogue of our process $(M_F^N(t), P_F^N(t))$ when the gene is always active.

**Lemma 1.** (1) For the convergence in distribution for the uniform norm on compact sets

$$ \lim_{N \to +\infty} \left( \int_0^t \overline{M}^N(u) \, du \right) = (\rho_2 t). $$

(2) For $T > 0$,

$$ \sup_{N \ge 1} \mathbb{E} \left( \sup_{0 \le t \le T} \overline{P}^N(t) \right) < +\infty. $$

*Proof.* From Relations (13), it is easily seen that the process $(\overline{M}^N(t))$ can be expressed ($L_1(Nt)$) where ($L_1(t)$) is an $M/M/\infty$ queue with arrival rate $\lambda_2$ and service rate $\mu_2$ with $L_1(0) = m_0$. See Chapter 6 of Robert [29] for example. Elementary stochastic calculus gives, for $t > 0$,

$$ (14) \quad L_1(Nt) = m_0 + \lambda_2 Nt - \mu_2 \int_0^{Nt} L_1(u) \, du + \mathcal{M}_1^N(t), $$

where $(\mathcal{M}_1^N(t))$ is a local martingale whose previsible increasing process is given by

$$ \langle \mathcal{M}_1^N \rangle (t) = \lambda_2 Nt + \mu_2 \int_0^{Nt} L_1(u) \, du. $$

Doobs'Inequality shows that the process $(\mathcal{M}_1^N(t)/N)$ vanishes for the convergence in distribution as $N$ gets large.

For $\varepsilon > 0$ and $x \in \mathbb{N}$, if

$$ T_x = \inf\{t \ge 0 : L_1(u) \ge x\}, $$

Proposition 6.10 of Robert [29] shows the convergence in distribution

$$ \lim_{x \to +\infty} \frac{\rho_2^x}{(x-1)!} T_x = E_0 $$
where $E_0$ is an exponential random variable with parameter $\mu_2 \exp(-\rho_2)$. This shows in particular the process $(L_1(Nt)/N)$ converges in distribution to 0 for the uniform convergence on compact intervals since

$$ \mathbb{P}\left(\sum_{0 \le t \le T} \frac{L_1(Nt)}{N} \ge \varepsilon\right) \le \mathbb{P}(T_{[\varepsilon N]}) \le NT. $$

From Equation (14), one gets

$$ \int_{0}^{t} \bar{M}^{N}(u) \, du = \frac{1}{N} \int_{0}^{Nt} L_{1}(u) \, du = \rho_{2}t + \frac{1}{\mu_{2}} \left( \frac{m_{0}}{N} - \frac{L_{1}(Nt)}{N} + \frac{M_{1}^{N}(t)}{N} \right) $$

and therefore assertion 1) of the lemma.

For the last assertion, the method is similar: one first write the evolution equation

$$ \bar{P}^N(t) = p_0 + \lambda_3 \int_0^t \bar{M}^N(u) \, du - \mu_3 \int_0^t \bar{P}^N(u) \, du + M_2^N(t), $$

where $(M_2^N(t))$ is a local martingale whose previsible increasing process is given by

$$ \langle \mathcal{M}_2^N \rangle (t) = \lambda_3 \int_0^t \bar{M}^N(u) \, du + \mu_3 \int_0^t \bar{P}^N(u) \, du. $$

Define $\bar{P}_*^N(t) = \sup\{\bar{P}^N(u) : 0 \le u \le t\}$, then for $0 \le t \le T$

$$ (15) \quad \mathbb{E} \left( \bar{P}_*^N(t) \right) \le p_0 + \lambda_3 \mathbb{E} \left( \int_0^T \bar{M}^N(u) \, du \right) \\ \qquad + \mathbb{E} \left( \sup_{0 \le u \le t} |\mathcal{M}_2^N(u)| \right) + \mu_3 \int_0^t \mathbb{E} \left( \bar{P}_*^N(u) \right) \, du. $$

Doob's Inequality gives, for $t \le T$,

$$ \mathbb{E} \left( \sup_{0 \le u \le t} |\mathcal{M}_2^N(u)| \right) \le 2\lambda_3 \int_0^T \mathbb{E} \left( \bar{M}^N(u) \right) du + 2\mu_3 \int_0^T \mathbb{E} \left( \bar{P}_*^N(u) \right) du, $$

and from the ergodic theorem for $(L_1(t))$ (recall that $\bar{M}^N(t) = L_1(Nt)$) one gets

$$ \lim_{N \to +\infty} \mathbb{E} \left( \int_{0}^{T} \bar{M}^{N}(u) \, du \right) = \rho_{3}T. $$

One concludes by using Equation (15) and Gronwall's Lemma. $\square$

<figure><img src="image_12.png" /></figure>

**Proposition 5.** The sequence $(P_F^N(t))$ is tight for the convergence in distribution of càdlàg processes.

*Proof.* Aldous' criterion for tightness is used. See Theorem 4.5 page 320 of Jacod and Shiryaev [15] for example. For $T > 0$, one denotes by $\mathcal{T}_T$ the set of stopping times associated to the filtration $(\mathcal{F}_t)$ which are bounded by $T$. For $\eta > 0$, let $\tau_1$, $\tau_2 \in \mathcal{T}_T$ be such that $\tau_1 \le \tau_2 \le \tau_1 + \eta$. The respective probabilities that, on the time interval $[\tau_1, \tau_2]$, no protein is made or that no protein is degraded are respectively given by

$$ \mathbb{E} \left( \exp \left( -\lambda_3 \int_{\tau_1}^{\tau_2} M_F^N(u) \, du \right) \right) \text{ and } \mathbb{E} \left( \exp \left( -\mu_3 \int_{\tau_1}^{\tau_2} P_F^N(u) \, du \right) \right) $$
By using the strong Markov property, one gets the relation

$$ \mathbb{P}(|P_F^N(\tau_1) - P_F^N(\tau_2)| \ge 1) \le 1 - \mathbb{E} \left( \exp \left( -\lambda_3 \int_{\tau_1}^{\tau_2} M_F^N(u) \, du \right) \right) \\ + 1 - \mathbb{E} \left( \exp \left( -\mu_3 \int_{\tau_1}^{\tau_2} P_F^N(u) \, du \right) \right). $$

With a simple coupling using the same Poisson processes $N_{2/3}^{+/-}$ of Equations (9) and (10) gives a process as in Definition 1 on the same probability space such that the relations $M_F^N(t) \le \bar{M}^N(t)$ and $P_F^N(t) \le \bar{P}^N(t)$ hold almost surely for all $t \ge 0$. From the last relation, one gets the inequality

$$ \mathbb{P}(|P_F^N(\tau_1) - P_F^N(\tau_2)| \ge 1) \le 1 - \mathbb{E} \left( \exp \left( -\lambda_3 \int_{\tau_1}^{\tau_1+\eta} \bar{M}^N(u) \, du \right) \right) \\ + 1 - \mathbb{E} \left( \exp \left( -\mu_3 \eta \sup_{0 \le t \le T} \bar{P}^N(t) \right) \right) \\ \le 1 - \mathbb{E} \left( \exp \left( -\lambda_3 \sup_{0 \le t \le T} \int_t^{t+\eta} \bar{M}^N(u) \, du \right) \right) \\ + 1 - \mathbb{E} \left( \exp \left( -\mu_3 \eta \sup_{0 \le t \le T} \bar{P}^N(t) \right) \right). $$

Lemma 1 gives the relation

$$ \lim_{N \to +\infty} \sup_{\tau_1 \in \mathcal{T}_T} \mathbb{E} \left( \exp \left( -\lambda_3 \int_{\tau_1}^{\tau_1+\eta} \bar{M}^N(u) \, du \right) \right) = e^{-\lambda_3 \rho_2 \eta} $$

and, for $\varepsilon > 0$, the existence of $K > 0$ such that

$$ \sup_{N \ge 1} \mathbb{P} \left( \sup_{0 \le t \le T} \bar{P}^N(t) \ge K \right) \le \varepsilon. $$

Consequently

$$ \lim_{\eta \to 0} \lim_{N \to +\infty} \sup_{\substack{\tau_1, \tau_2 \in \mathcal{T}_T \\ \tau_1 \le \tau_2 \le \tau_1 + \eta}} \mathbb{P} \left( |P_F^N(\tau_1) - P_F^N(\tau_2)| \ge 1 \right) = 0, $$

hence, by Aldous' criterion, the tightness of the sequence $(P_F^N(t))$ is established. The proposition is proved.

□

**A.2. Convergence of Occupation Measures.** For $N \ge 1$ and $T > 0$, one defines the random measure $\Lambda^N$ on $\mathcal{E}_T \stackrel{\text{def.}}{=} \{0, 1\} \times \mathbb{N}^2 \times [0, T]$ as follows, for a non-negative Borelian function $G$ on $\mathcal{E}_T$,

$$ \Lambda^N(G) = \int_0^T G(X_F^N(u), u) \, du. $$

If $A$ is a Borelian subset of $\mathcal{E}_T$, $\Lambda^N(A)$ denotes $\Lambda^N(1_A)$.

**Proposition 6.** The sequence $\Lambda^N$ of random measures is tight and any of its limiting points $\Lambda$ can be written as

$$ \Lambda(F) = \sum_{(i,m,p) \in \mathcal{S}} \int_0^T G(i, m, p, u) \pi_p(i, m) \nu_u(p) \, du. $$
where, for any $u \le T$, $\nu_u$ is a positive measure on $\mathbb{N}$ such that, almost surely,

$$ \int_0^t \nu_u(\mathbb{N}) \, du = t, \quad \forall t \le T, $$

and, for $p \in \mathbb{N}$, $\pi_p$ is the invariant distribution of the Markov process on $\{0,1\} \times \mathbb{N}$ whose transition rates are given by, for $(i,m) \in \{0,1\} \times \mathbb{N}$,

$$ (16) \qquad \begin{cases} (i,m) \to (1-i,m) & \text{at rate } \lambda_1^+ i + \lambda_1^- p(1-i), \\ (i,m) \to (i,m+1) & \lambda_2 i, \\ (i,m) \to (i,m-1) & \mu_2 m. \end{cases} $$

Additionally, one has

$$ (17) \qquad \sum_{(i,m) \in \{0,1\} \times \mathbb{N}} m \pi_p(i,m) = \frac{\lambda_1^+}{\lambda_1^+ + \lambda_1^- p} \frac{\lambda_2}{\mu_2}. $$

*Proof.* For $K>0$, if $K_K$ is the compact subset $\{0,1\} \times [0,K]^2 \times [0,T]$ of $\mathcal{E}_T$, then

$$ \mathbb{E} (\Lambda^N (\mathcal{E}_T \setminus K_K)) \le \int_0^T \mathbb{P} (M_F^N(u) \ge K) \, du + T \mathbb{P} \left( \sup_{0 \le u \le T} P_F^N(u) \le K \right). $$

By using the same coupling as in the proof of Proposition 5, one gets that

$$ \mathbb{E} (\Lambda^N (\mathcal{E}_T \setminus K_K)) \le \int_0^T \mathbb{P} (\bar{M}^N(u) \ge K) \, du + T \mathbb{P} \left( \sup_{0 \le u \le T} \bar{P}^N(u) \le K \right). $$

By Lemma 1, for $\varepsilon > 0$, there exists some $K$ such that

$$ \sup_{N \ge 1} \mathbb{E} (\Lambda^N (\mathcal{E}_T \setminus K_K)) \le \varepsilon. $$

Consequently, the sequence $(\Lambda^N)$ of random Radon measures on $\mathcal{E}_T$ is tight. See Dawson [7, Lemma 3.28, page 44] for example.

Let $\Lambda$ be a limiting point of some subsequence $(\Lambda^{N_k}(\cdot))$. By using Radon-Nikodym's Theorem, see Chapter 8 of Rudin [32] for example, it is not difficult to see that there exists some non-negative random variables $(\ell_u(x)(\omega), (\omega, x, u) \in \Omega \times S \times [0, T])$ such that $(\omega, x, u) \mapsto \ell_u(x)(\omega)$ is measurable and $\Lambda$ can be expressed as

$$ \Lambda(G) = \sum_{x \in S} \int_0^T G(x, u) \ell_u(x) \, du. $$

From the domination relation of Lemma 1, one gets that, almost surely, there is no loss of mass, i.e.

$$ (18) \qquad \int_0^t \ell_u(S) \, du = t, \quad \forall t \le T, $$

holds almost surely. Now take a function $f$ with bounded support on $S$, by using Equation (12), it is not difficult to show that the process $(\langle W_f^N \rangle(t))$ satisfies the relation

$$ \lim_{N \to +\infty} \frac{1}{N^2} \mathbb{E} (\langle W_f^N \rangle (T)) = 0, $$

by Doob's Inequality this implies that the martingale $(W_f^N(t)/N)$ converges in distribution to 0 for the uniform norm on $[0, T]$.
By dividing Relation (11) by $N$, one gets that, for the convergence in distribution, the relation

$$ \lim_{N \to +\infty} \left( \int_0^t \lambda_1^+ (1-I_F^N(u)) \Delta_1(f)(X_F^N(u)) \, du + \int_0^t \lambda_1 P_F^N(u) I_F^N(u) \Delta_1(f)(X_F^N(u)) \, du + \int_0^t \lambda_2 I_F^N(u) \Delta_2^+(f)(X_F^N(u)) \, du + \int_0^t \mu_2 M_F^N(u) \Delta_2^-(f)(X_F^N(u)) \, du \right) = 0. $$

holds. The convergence of the sequence $(\Lambda^{N_k})$ gives that the relation

$$ \sum_{x=(i,m,p) \in S} \int_0^t \lambda_1^+(1-i) \Delta_1(f)(x) \ell_u(x) \, du + \int_0^t \lambda_1^- ip \Delta_1(f)(x) \ell_u(x) \, du + \int_0^t \lambda_2 i \Delta_2^+(f)(x) \ell_u(x) \, du + \int_0^t \mu_2 m \Delta_2^-(f)(x) \ell_u(x) \, du = 0 $$

holds almost surely for all $0 \le t \le T$ and for all indicator functions of elements $S$. Now, for $p \in \mathbb{N}$ and $g$ a function with finite support on $\{0, 1\} \times \mathbb{N}$, define $f(i, m, p) = g(i, m)$, the above relation gives

$$ (19) \quad \sum_{x=(i,m,p) \in S} \left( \begin{aligned} & \ell_u(i,m,p) \lambda_1^+(1-i) \Delta_1(g)(i,m) + \ell_u(i,m,p) \lambda_1^- ip \Delta_1(g)(i,m) \\ & + \ell_u(i,m,p) \lambda_2 i \Delta_2^+(g)(i,m) + \ell_u(i,m,p) \mu_2 m \Delta_2^-(g)(i,m) \end{aligned} \right) = 0 $$

holds almost surely for all $u \in A \subset [0, T]$ and $[0, T] - A$ is negligible for Lebesgue measure. Relation (19) shows that for $u \in A$, the vector $(\ell_u(i, m, p))$ is proportional to the invariant distribution $\pi_p$ of the Markov process on $\{0, 1\} \times \mathbb{N}$ whose transition rates are given by Relations (16).

One gets therefore the existence of a constant $\nu_u(p)$ such that $\ell_u(i, m, p) = \nu_u(p)\pi_p(i, m)$ for all $(i, m, p) \in S$. Equation (18) gives the relation

$$ \int_0^t \nu_u(\mathbb{N}) \, du = t, \quad \forall t \le T. $$

Hence one has $\nu_u(\mathbb{N}) = 1$ almost surely for all $u \in A_1 \subset [0, T]$ and $[0, T] - A_1$ is negligible for Lebesgue measure.

Straightforward calculations as in the proof of Proposition 2 complete the proof of the proposition to give Relation (17). □

E-mail address: Renaud.Dessalles@inria.fr

E-mail address: Vincent.Fromion@jouy.inra.fr

(R. Dessalles, V. Fromion) MAIAGE, INRA, DOMAINE DE VILVERT, 78350 JOUY-EN-JOSAS, FRANCE

E-mail address: Philippe.Robert@inria.fr

URL: http://www-rocq.inria.fr/~robert

(R. Dessalles, Ph. Robert) INRIA PARIS—ROCQUENCOURT, DOMAINE DE VOLUMEAU, 78153 LE CHESNAY, FRANCE