Dynamics of the stochastic chemostat with Monod-Haldane response function

The stochastic chemostat model with Monod-Haldane response function is perturbed by environmental white noise. This model has a global positive solution. We demonstrate that there is a stationary distribution of the stochastic model and the system is ergodic under appropriate conditions, on the basis of Khasminskii’s theory on ergodicity. Sufficient criteria for extinction of the microbial population in the stochastic system are established. These conditions depend strongly on the Brownian motion. We find that even small scale white noise can promote the survival of microorganism populations, while large scale noise can lead to extinction. Numerical simulations are carried out to illustrate our theoretical results.

The chemostat, known as a continuous stir tank reactor (CSTR) in the engineering literature, is a basic piece of laboratory apparatus used for the continuous culture of microorganisms. It occupies a central place in mathematical ecology and has played an important role in many fields. In ecology it is often viewed as a model of simple lake or an ocean system. It can also model the wastewater treatment process 1 or biological waste decomposition 2 . In some fermentation processes, the chemostat plays a central role in the commercial production of genetically altered organisms. The theoretical investigation of the dynamics of microbial interaction in the chemostat was initiated by Monod 3 and Novick et al. 4 in 1950. There is also an extensive literature on the chemostat (for example, see refs [5][6][7][8][9][10][11]  where S(t), x(t) stand for the concentrations of nutrient and microbial population at time t respectively; S 0 denotes the imput concentration of nutrient and Q represents the volumetric flow rate of the mixture of nutrient and microorganism; the coefficient δ is the ratio of the biomass of the microbial population produced by the nutrient consumed. The growth rate of the microbial population is represented by the function μ(S) (μ ≤ < ∞ S m ( ) ), which is generally assumed to be non-negative. That is μ = (0) 0, μ(S) > 0 for S > 0. Some experiments and observations indicate that not only insufficient nutrient but also excessive nutrient may inhibit the growth of a microbial population in the chemostat. To model such growth, Andrews 12 suggested a non-monotonic response function: , (1 2) which is called the Monod-Haldane growth rate (inhibition rate). Parameter m is the maximum growth rate of the microorganism, a is the Michaelis-Menten (or half-saturation) constant. The term KS 2 models the inhibitory effect of the substrate at high concentrations. The chemostat model with Monod-Haldane response function has been studied by many researchers (see e.g. Wang et al. 13 , Pang et al. 14 , Baek et al. 15 2 It is not difficult to verify that one of the two interior equilibria is a stable nodal point and the other is a saddle point. In mathematical biology, when the orbits tend to the stable nodal point (S * , x * ), the continuous culture of microorganism is considered to be successful. We refer the reader to Chen 16 for more details.
System (1.3) is deterministic with parameters assumed to be constant. Environmental fluctuations are ignored, which, from the biological point of view, is unrealistic. Ecosystem dynamics are inevitably affected by environmental noise and it is more realistic to include the effect of stochasticity rather than to study models that are entirely deterministic. With respect to chemostat models, even though the experimental results that are observed in well-controlled laboratory conditions have been shown to match closely with the prediction of deterministic models involving ordinary differential equations, we cannot ignore the difference that may occur in operational conditions. Imhof et al. 17 derived a stochastic chemostat model by considering a discrete-time Markov process with jumps corresponding to the addition of a centered Gaussian term to the deterministic model. With the time step converging to zero, they reported that the stochastic model may lead to extinction in some cases in which the deterministic model predicts persistence. In 18 , Campillo et al. considered a set of stochastic chemostat models that are valid on different scales. Xu and Yuan 19 dealt with the stochastic chemostat in which the maximal growth rate m is perturbed by white noise and obtained a new break-even concentration λ ∼ which completely determines the persistence or extinction of the microbial population. Zhao and Yuan 20 further computed the probability for extinction and persistence in mean of the microbial population in 19 using stochastic calculus. Wang et al. 21 investigated the periodic solutions for the stochastic chemostat model with periodic washout rate, on the basis of Khasminskii's theory (see, ref. 22 Chapter 3) for periodic Markov processes.
In this paper, the dynamical behaviour of a two-dimensional chemostat model with Monod-Haldane response function under stochastic perturbation is investigated. The white noise is incorporated in stochastic system (1.3) to model the effect of a randomly fluctuating environment. The substrate and microorganism population are usually estimated by an average value plus errors (i.e. noise intensities) which are usually normally distributed. The standard deviations of the errors may depend on the population sizes. Utilizing the approach used in 17 to include stochastic effects (readers can also refer to 23 , Appendix A to see the construction of this kind of stochastic model), the stochastic chemostat model with Monod-Haldane response function takes the form: where B 1 (t), B 2 (t) are independent standard (i.e. with variance t) Brownian motion defined in a complete probability space are their intensities. Since an ODE model is never exact, but only an approximation of reality; adding a linear perturbation using Brownian motion (terms σ  SB 1 1 , σ  xB 2 2 ) as in model (1.5) helps to determine how the dynamics can possibly change by approximating the error in the ODE model, and the error is normally distributed, σ σ > ( 0) i i 2 , i = 1, 2 are constants that reflect the sizes of error (stochastic effects).
Very few other studies have appeared on the stochastic chemostat model with the Monod-Haldane functional response. Campillo et al. 24 considered a stochastic model of the chemostat, with both non-inhibitory (Monod) and inhibitory (Haldane) growth functions, as a diffusion process and used a finite difference scheme to approximate the solutions of the associated Fokker-Planck equation. They considered the stochastic model from a 'demographic noise' point of view. Instead of using S and x in the stochastic terms in model (1.5), they used S and x , respectively.
In this paper, our aim is to reveal how the environmental noise affects a microbial population in the chemostat with Monod-Haldane response function. First, in section 2 we prove that there is a unique positive solution for model (1.5). Then in section 3 we show that for any initial value  ∈ + S x ( (0), (0)) 2 , there is a stationary distribution for system (1.5) and it is ergodic under appropriate conditions. Sufficient conditions for extinction of the microorganism are presented in section 4. Finally, to illustrate our main conclusions, examples and numerical simulations are given in section 5.

Existence and uniqueness of the positive solution
Before investigating the dynamical behavior of system (1.5), the existence of a global positive solution is proved. In order to ensure the stochastic chemostat make sense, we need to show at least not only that this SDE model has a unique global solution but also that the solution will remain in  + 2 whenever it starts from there. First we give the definition of a solution of stochastic differential equation (see 25 , Chapter 2 for more details). Throughout this paper, unless otherwise specified, let Ω ≥ P ( , , { } , ) satisfying the usual conditions (i.e. it is right continuous and  0 contains all P-null sets). Let l l m 0 be both Borel measurable. Consider the l-dimensional stochastic differential equation of It ô type with initial value X(t 0 ) = X 0 . SDE (2.1) is equivalent to the following stochastic integral equation: is called a solution of equation (2.1) if it has the following properties: where P denotes the probability of an event. By utilizing the methods described in 25 , the coefficients of the equations would be required to satisfy a local Lipschitz condition and a linear growth condition. However, the Haldane function μ → = + + S S ( ) mS a S KS 2 is nonlinear, coefficients of system (1.5) do not satisfy the linear growth condition. Thus in this section, using the Lyapunov analysis method 26 , we prove that the solution of the system (1.5) is positive and global.
Proof: Consider the diffusion process as follows Since the coefficients of system (2.3) are locally Lipschitz continuous, there is a unique local solution for system (2.3). Let S = e μ , x = e μ , Itô's formula (given in Section 3) implies that system (1.5) has a unique local positive solution. Hence it suffices to prove that this unique local positive solution of system (1.5) is global.
On the basis of the discussion above, we know that there is a unique local solution (S(t), where τ e is the blow up time. If we can prove τ e = ∞ a.s., then the solution will be global. . For each integer ≥ n n 0 , define the stopping time: Throughout this paper, we set φ = ∞ in f (as usual, φ denotes the empty set). Clearly, τ n is increasing as for all ≥ t 0. In other words, to complete the proof all we need to show is that τ = ∞ . . ∞ a s If not, there exists a pair of constants T > 0 and ε ∈ (0, 1) such that Hence there is an integer n 1 ≥ n 0 such that Obviously, this function is non-negative for all S, x ≥ 0. Using Itô's formula (see Theorem 3.1), we get where  is the generating operator of system (1.5) and Taking expectation of both sides implies that Note that for every ω ∈ Ω n , there exists at least one τ ω S( , ) n or τ ω x( , ) n that equals either n or n 1 , and then It then follows from (2.4) and (2.5) that

Stationary distribution and ergodicity
In the study of the dynamics for stochastic systems, ergodicity is one of the most important and significant characteristics. The ergodic property for chemostat implies that the stochastic model has a unique stationary distribution which predicts the survival of the microbial population in the future. In addition, the ergodic property gives a reason why the integral average of a population system converges to a fixed point whilst the population system may fluctuate around as time goes by. Thus research on the ergodicity of population system is essential from a biological perspective (see 27 ). In this section, before the proof of the ergodicity, several auxiliary results are derived for the stationary distribution. For more details, we refer the reader to 28 . First we introduce the multi-dimensional Itô formula.
Theorem 3.1 (The multi-dimensional Itô formula) (Chapter 1) 25 Let X(t) be a l-dimensional It ô process on t ≥ 0 with the stochastic differential . Then V(X(t), t) is an Itô process with the stochastic differential given by Process in  l ( l denotes the l-dimensional Euclidean space) described by the stochastic differential equation

Remark 3.1
The proof is given in 22 . To show the existence of the stationary distribution μ ⋅ ( ) of system (1.5), it is enough for us to take  + 2 as the whole space. To validate (B.1), we need to prove that for any bounded domain A regular process X(t) described by (3.2) with nonsingular diffusion matrix (i.e., the smallest eigenvalue of Λ(x) is bounded away from zero in every bounded domain in  l .) is said to be recurrent if there exists a bounded domain U such that for all

Remark 3.2
The definition of regularity and recurrence come from 32 . Let X(t) be a regular temporally homogeneous Markov process in  l . If X(t) is recurrent with respect to some bounded domain U, then it is recurrent with respect to any nonempty domain in  l . Since the existence of the positive solution for model (1.5) has been obtained by Theorem 2.1, it is enough for us to take  + 2 as the whole space. Before proving the ergodicity, we define the notation = . Define a continuous Theorem 3.2 If σ 1 , σ 2 satisfy: x ( (0), (0)) 2 , there is a stationary distribution μ ⋅ ( ) for system (1.5) and the system is ergodic.
Proof: For simplicity let S, x stand for S(t), x(t), respectively. We will verify that (B.1) and (B.2) hold under condition (3.3-3.4), according to Lemma 3.1 on the existence of the stationary distribution. System (1.5) can be written as the following system: hence the diffusion matrix is Let D to be any bounded domain in  + 2 , then there exists a positive constant This implies that condition (B.1) is satisfied. Next we will construct a nonnegative C 2 -function V(S, x) and a closed set where R is a positive constant. This assures that condition (B.2) is satisfied. Choose a positive constant M big enough such that and the function Φ S ( ) will be given later in (B). Define a C 2 -function: We derive the unique solution S 0 of (3.5) from the monotonicity property of the left part. Thus the following equation  , then we get  Substituting (3.7-3.9) into (3.6), The above cases lead to < − V 1  , respectively. Moreover, by reviewing condition (A) we obtain Take ε small enough, and let

Extinction
In this section, we discuss conditions that predict the failure of the continuous culture of microorganisms both in the case of environmental noise is ignored and in the case of big intensities of white noises, i.e. noise may lead to extinction of the microorganism in the reactor.
. Thus we conclude that that is, x(t) tends to zero exponentially with probability one, .
The proof is complete. ☐

Remark 4.1
We refer to the condition (ii) in Theorem 4.1, which tells us that the microorganism species may die out when dilution rate Q and white noise are not large. While if the strength of white noise is large enough such that (i) holds, then the microorganism population will also become extinct, which never happens in the deterministic system (1.3) without environmental perturbations. Moreover, if the dilution rate Q is big enough which leads to the extinction of the species in the deterministic chemostat while the noise is not big (condition (i)), the microorganism species in stochastic chemostat (1.5) will also die out.  .

→∞
x t a s lim ( ) 0, t Therefore the microorganism species x will go to extinction exponentially almost surely. In other words, the microorganism population will die out at an exponential rate with probability one.

Examples and numerical simulations
Generally, nonlinear stochastic differential equations, such as system (1.5), are usually too complex to be solved exactly, we can only theoretically prove the existence and uniqueness of the positive solution (see Theorem 2.1). In order to illustrate the analytical results in Section 3 and 4, we need to obtain approximate solutions of stochastic system (1.5) with given initial values and parameter values, via numerical methods and algorithms that can be implemented by Matlab. Due to Brownian motion and noise terms σ SdB t ( ) 1 1 , σ xdB t ( ) 2 2 , we use Milstein's higher order method 33 to obtain the approximate solutions of stochastic chemostat model (1.5). J. Higham has showed the convergence of Milstein's numerical method for stochastic differential equations (see Section 6, ref. 33 ). We consider the following discretization equation in Milstein's type to iteratively calculate the approximate solutions of stochastic system (1.5) in Matlab programs: 2 where ξ i 1, and ξ i 2, are N(0, 1)-distributed independent Gaussian random variables, σ σ , 1 2 are intensities of white noise and time increment ∆ > t 0. In this way, the approximate solutions (performance results of (5.1) from Matlab) will converge to the explicit solutions of stochastic system (1.5) very fast.
Example 5.1 From Theorem 3.2, we expect that under some appropriate conditions, the stochastic system (1.5) has a stationary distribution. Choose the initial value = .
. S x ( (0), (0)) (1 2, 0 8) and assume that the parameters in (1. in Theorem 3.2 are all satisfied. Then we conclude that there is a stationary distribution (see the right two subgraphs in Fig. 1) of the stochastic system (1.5) and it is ergodic. Numerical simulations in Fig. 1 support this conclusion clearly, illustrating that the standard deviation σ keeps processes S(t), x(t) moving around the solution of the corresponding deterministic chemostat model.  Figure 2 illustrates that the white noise still keeps processes S t ( ), x t ( ) moving around the solution of the deterministic chemostat model, while the random oscillations become stronger.

Example 5.3
In order to show how the dilution rate and the random perturbation influence the extinction of the microbial population in system (1.5), we choose the initial value = . S x ( (0), (0)) (3, 0 5) and assume that = .
Q 2 2, σ = . We therefore conclude that by Theorem 4.1 (i), the solution of (1.5) obeys  = . 0 1, 02 1 2 . The stochastic chemostat model still has a stationary distribution, even though the amplitude oscillations are stronger. and so x t ( ) tends to zero exponentially with probability one. These results are illustrated in the numerical simulations in Fig. 4. We can observe from Figs 3-4 that under the same random noises, the extinction time of the microbial population occurs later when = .
Q 1 6 than for = . Q 2 2 for both the deterministic system (1.3) and the stochastic system (1.5).  . In the second subgraph, the large washout rate leads to the extinction of x(t) both in the deterministic and stochastic chemostat models. . ⁎ E (0 2715, 2 7657) and E* is a globally asymptotically stable node. The numercial simulations in Fig. 5 confirm the extinction of x t ( ) in the stochastic system (1.5) and the stability of E* in the deterministic system (1.3), and it is clear that with the increasing of the intensity of white noise, the tendency for exponential extinction increases.

Discussion and Conclusion
Understanding the effects of environmental stochasticity on the survival of microbial populations is of theoretical and practical importance in population biology. In this paper we consider a stochastic chemostat model with the Monod-Haldane response function. We first determine when this system has a unique global positive solution. Then we show that the stochastic chemostat model admits a unique stationary distribution which is ergodic if the scale of the random perturbations is relatively small. Sufficient criteria are provided for the extinction of the population of microorganisms. We also find that both strong random noise and a large dilution rate can lead to extinction. Our conclusions are all expressed in terms of the system parameters and the intensity of the Brownian motion. This means that white noise can have a major impact on the survival or extinction of microorganism populations.
How does environmental stochasticity change the predicted outcome for the microbial population in (1.3)? Actually, the stochastic chemostat (1.5) has no equilibrium. The stationary distribution shows that the solutions of the stochastic system fluctuate in a neighborhood of the positive equilibrium of the corresponding deterministic system, which can be regarded as weak stability. In Theorem 3.2, we define a new dilution rate  Q which is in terms of the original dilution rate Q and the random noise on the microbial population σ 2 , i.e. σ = +  Q Q   Table 1. Effects of the dilution rate and stochasticity on the survival-extinction of the microorganism.
has a unique stationary distribution and it is ergodic, which means the microbial population survives for all time regardless of the initial conditions. This result illustrates that random noise in low levels is advantageous for species survival. In Figs 1 and 2, numerical simulations show the dynamics of model (1.5) with different σ i 2 , = i 1,2. Comparing Fig. 1 with Fig. 2 we see that with increasing σ i 2 , the oscillations become stronger, but both figures reveal that the microbial population persists and demonstrate the existence of the unique stationary distribution.
With increasing Q, Theorem 4.1 shows that x t ( ) tends to extinction exponentially. Furthermore, it can be observed from Figs 3-5 that the extinction of the microbial population occurs more quickly for system (1.5) with σ = .
Q 2 2 tends to extinction faster than that with = . Q 1 6. That is, if the volumetric flow rate Q is increased, the microbial population tends to extinction faster. Thus from the biological point of view, a well-controlled volumetric flow rate in the chemostat is also the key to getting a successful continuous culture of microorganisms.
We have already shown that if the densities of the environmental noise are small, the deterministic and stochastic systems have similar dynamical behavior, both in the case of persistence and extinction. Taking the dilution rate = .
Q 1 1 in system (1.3) with σ = 0 i , = i 1, 2, there exists an asymptotically stable interior equilibrium. Increasing σ 2 , the stochastic trajectory x t ( ) will first fluctuate around this equilibrium, then the phenomenon of noise-induced extinction occurs when σ 2 is big enough. In this stochastic case, the microorganism tends to extinction exponentially even faster than in the case of large dilution rate = .
Q 2 2. These results are summarized in the following Table 1.
To complement our theoretical result we study the global dynamics of the deterministic chemostat system (1.3) and the stochastic chemostat system (1.5) in cases that are not counter by our theorems. We set parametric values δ = .
= . . Then according to the mathematical results in 16 , the deterministic chemostat model (3) has two interior equilibrium points, one an asymptotically stable node = .
E (2 4, 0) 0 that is also an asymptotically stable node. The topological structure of the deterministic chemostat model (1.3) in Fig. 6(a) shows that E 0 and E 1 are both stable. The microbial population will survive if and only if the orbits tend to the interior stable node E 1 . This means in order for the population to survive, we must choose appropriate initial concentrations of nutrient and microbial population.
0 05 1 , σ = . 0 05 2 . With the same initial conditions and parameter values, Fig. 6(b) suggests that processes S t ( ), x t ( ) move around the trajectories of the deterministic chemostat model (1.3) with small fluctuations. They have similar properties to that of the deterministic model. Note that the existence conditions for a stationary distribution are only sufficient, but not necessary. In case of = .
Q 0 7, these parameters do not satisfy the hypotheses of Theorem 3.2, and so we cannot estimate whether the stochastic chemostat is ergodic or not. This is an open problem.