Stochastic dynamic analysis of a chemostat model of intestinal microbes with migratory effect

1.   Introduction

We are not individuals but a collection of living things ^[1]. Up to 10 million microbes inhabit every square centimeter of skin ^[2]. Surface-dwelling microbes even play a crucial role in maintaining the homeostasis of the host microbial community ^[3]. At the same time, our teeth, throat, and esophagus are even more infested with microbes, which accumulate thousands of times more than the surface of the skin. The human body has a large number of microorganisms living in the gut, and the vast majority of them are beneficial to us. A correct view of their value and continuous research and discovery is the proper attitude of human beings towards the microbes in the body. Researchers have long applied various methods to study the role of the microbiome in living organisms and achieved significant results. In addition to the application of biological methods ^[4], the application of mathematical methods is undoubtedly a great help for this research direction ^[5,6].

The flow characteristics of chemostat model are very suitable for simulating the flow ecological environment. And the culture of microorganisms in the chemostat can be regulated by controlling the concentration of nutrients. Combined with the features of intestinal circulation, the chemostat model can be used to further investigate the microbial system in the intestinal tract. Presently, research on chemostat models has been further improved ^[9,10,11], and researchers have studied ecosystems with different types of functional responses through chemostat models, such as Monod-Haldane ^[12,13,23], Lotka-Volterra ^[14], Holling II ^[16], etc. In spite of this, few papers have applied the chemostat model to specific ecological environment. This paper proposes that, it is an innovation to apply the chemostat model to the characterization of intestinal microorganisms. In addition, based on the study of reaction function in the above literature and the characteristics of intestinal microorganisms, we adopted Holling II functional response function. In the gut, microbes circulate with material and may attach to the folds in the intestinal wall as they flow, or they may fall off the wall. It is an innovative point of this paper to consider and study the impact of the migration between the flowing microorganisms and the microorganisms on the intestinal wall as the microorganisms on the two plates. This phenomenon can be described by adding a migration term to the model. In ^[7], the authors considered the variable flow rate, microbial decomposition, and other factors to construct a detailed intestinal chemostat model.

Inevitably, microbial survival in the gut is affected by uncertainties such as body temperature, water supply and pH-value.

In stochastic biological models, scholars use Brownian motion to describe the perturbation of the system caused by uncertain factors and have made important achievements ^{[8,17,18,21,22]}. On the basis of deterministic models, many authors consider the linear disturbance of mortality, which further enrich the content of the models ^{[15,25,27,29]}. In addition, bilinear perturbation and nonlinear perturbation have been considered in ^[24,28], respectively. In ^[26], stochastic models perturbed by Lévy noise are concerned. From what has been discussed above, we assume that uncertain factors can cause an impact on microbial growth, that is, environmental disturbance, mainly interferes with the consumption function. Then the following chemostat model with nonlinear perturbation, intraspecific competition and migration is established.

Let $(\Omega, \mathcal{F}, \mathbb{P})$ be a complete probability space with a filtration $\{\mathcal{F}_t\}_{t\geq0}$ satisfying the usual conditions (i.e., it is increasing and right continuous while $\mathcal{F}_0$ contains all $\mathbb{P}$-null sets). At the same time, there are the corrective coefficients of Wong-Zakai in the model, and the same calculation results can be obtained by using Itô integral and Stratonovich integral ^[30]. The function $B(t)$ is independent standard Brownian motion defined on this complete probability space with $B(0) = 0$ and $\sigma > 0$ represents for the intensity of the white noise. And other parameters are defined in the Table 1.

The characteristics of the proposed model are the application of the chemostat model to simulate the intestinal environment, and the introduction of microbial migration, interspecific competition and nonlinear disturbance on the basis of the classical chemostat model. Compared with ^[7], we consider the effect of environmental disturbance on the microbial response function based on the corresponding autonomous model.

The organizational structure of the paper is as follows. In Section 2, we provide the lemmas and partial proofs used in the argument of this paper. In Section 3, we verify the existence and uniqueness of global positive solution and the boundedness of solution of the system ($1.1$). Section 4 presents all the research results for the proposed system. In the first subsection, sufficient conditions for microbial extinction are proved. Then, the condition of persistence in mean of the system is given. Next, the existence of ergodic stationary distribution is established. Finally, the behavior of system ($1.1$) is simulated numerically in Section 5, and the conclusions are given.

2.   Preliminaries

We let $\mathbb{R}_+ = \{x\in \mathbb{R}: x > 0\}$. If $f(t)$ is an intergrable function on $[0, +\infty)$, let $\left \langle f(t)\right \rangle = \frac{1}{t}\int_0^tf(\theta) \text{d}\theta$. Here we give the common lemmas of stochastic differential equations which will be used in the proof below.

Lemma 2.1 (Stationary distribution ^[20]). Let $X(t)$ be a homogenous Markov process in $E_l$ ($E_l$ denotes $l$-dimensional Euclidean space), and it is described by the following stochastic differential equation

The diffusion matrix is defined as follows

If there is a bounded open set $U\subset R^{n}$ with regular that satisfies the following conditions.

(i) For any $x\in U, \varepsilon \in R^{n}$, there is a constant $\varsigma > 0$, satisfying that:

(ii) For any $x\in R^{n}\setminus U$, there is a $C^{1, 2}$ function $V$, such that $LV < 0$.

Then, the Markov process $X(t)$ exists a stationary distribution $\psi (\cdot)$.

Lemma 2.2 (Existence and uniqueness of global positive solutions). For any initial value $(S(0), x_1(0), x_2(0))\in \mathbb{R}^3_{+}$, there is a unique positive solution $(S(t), x_1(t), x_2(t))$ of system $(1.1)$ on $t\geq 0$ and the solution will remain in $\mathbb{R}^3_{+}$ with probability one.

Proof. Since the coefficients of the system (1.1) satisfy the local Lipschitz conditions, then for any given initial value $(S(0), x_1(0), x_2(0))\in \mathbb{R}^3_{+}$, there exists a unique local solution $(S(t), x_1(t), x_2(t))$ on $t\in[0, \tau_e)$, where $\tau_e$ is the explosion time. Now we just need to prove $\tau_e = \infty$ a.s..

Let $k_0 > 0$ be so sufficiently large that any initial value $S(0)$, $x_1(0)$ and $x_2(0)$ lying within the interval $\left[\frac{1}{k_0}, k_0\right]$. For any integer $k$ satisfying $k\geq k_0$, define the stopping time as follows:

with $\inf\emptyset = \infty$ (where $\emptyset$ denotes the empty set). It is easy to get that $\tau_k$ is increasing as $k \rightarrow \infty$. Let $\tau_\infty = \lim_{k\rightarrow \infty}\tau_k$, hence $\tau_\infty\leq\tau_k$ a.s.. Next, we only need to verify $\tau_\infty = \infty$ a.s.. If this statement is false, then there exist two constants $T > 0$ and $\epsilon\in(0, 1)$ such that

Thus there is an integer $k_1\geq k_0$ such that

Define a $C^2$-function $I$: $\mathbb{R}_{+}^3\rightarrow \mathbb{R_+}$ as follows,

Clearly,

Applying Itô's formula yields

where

then

where $\kappa$ is a positive constant.

The proof is completed.

Lemma 2.3 (Boundedness). For any initial value $(S(0), x_1(0), x_2(0))\in\mathbb{R}^3_{+}$, the solution $(S(t), x_1(t), x_2(t))$ of the stochastic system $(1.1)$ is bounded. And it satisfies

where $\omega = \max\left\{\frac{DS_0}{\mu}, S(0)+\frac{a}{b}x_1(0)+\frac{a}{b}x_2(0)\right\}$.

Proof. Denote

then from the system $(1.1)$, we obtain

where $\mu = \min \left\{D, v\right\}$, $Y(t)$ is the solution of $(2.4)$. And it has the following form

Furthermore,

The proof is completed.

3.   Main results

3.1.   Extinction

In this subsection, we establish threshold condition for intestinal microbial extinction.

Define a parameter

Theorem 3.1. (1) If $\frac{bS_0}{m+S_0} < \frac{2}{\sigma^2}$ and $R_{1}^{*} < 0$, then microbials tend to extinction almost surely.

(2) If $\frac{bS_0}{m+S_0}\geq\frac{2}{\sigma^2}$ and $R_{2}^{*} < 0$, then microbials tend to extinction almost surely.

Proof. Through the system, we can get,

Applying Itô's formula gives

Consider the function

Obviously, when $\frac{bS}{m+S} = \frac{2}{\sigma^2}$, $g(S)$ gets the maximum. Otherwise, $S\leq S_0$. Therefore, the maximum value of $g(S)$ is discussed in the following two cases.

Case 1. If $\frac{bS_0}{m+S_0} < \frac{2}{\sigma^2}$,

Combining ($3.1$), we get

Integrating from $0$ to $t$ and dividing by $t$ on both sides of ($3.2$) give

According to the strong law of large numbers ^[19], we obtain

So we can get

which means

Case 2. If $\frac{bS_0}{m+S_0}\geq \frac{2}{\sigma^2}$,

Combining ($3.1$), we get

Similar to Case 1, we can obtain

which means

3.2.   Persistence in mean

In this subsection, we examine the conditions sufficient for microbes to persist over time.

Define

Theorem 3.2. If $R^{*} > 0$, then for any initial value $(S(t), x_1(t), x_2(t))\in R_{+}^{3}$, the system $(1.1)$ is persistent in mean.

Proof. Integrating the equations of the system ($1.1$) from $0$ to $t$ and dividing $t$ on both sides yield

Then we can get

Applying Itô's formula gives

from which we can obtain

From the strong law of large number of martingales, we have

Since

When $R^{*} > 0$ we can get

where

This completes the proof.

3.3.   Existence of ergodic stationary distribution of system (1.1)

In this subsection, we explore sufficient conditions for the existence of a unique stationary ergodic distribution.

Define

Theorem 3.3. Assume $R > 1$, then system $(1.1)$ exists a unique ergodic stationary distribution.

Proof. In order to prove this theorem, it is necessary to prove that system ($1.1$) satisfies the conditions in Lemma 2.1. Define

$V_1(S, x_1, x_2)$ is a continuous function, so it exists a minimum $V_{1 min}$. We can get a $C^2$-function $V:\mathbb{R}_+^3\rightarrow \mathbb{R}_+:$

Applying Itô's formula gives

where $c_1 = S_0$, $c_2 = \frac{DS_0}{D+v+r_1+\frac{\sigma^2b^2}{2}}$ and $c_3 = \frac{DS_0}{D+v+r_2+\frac{\sigma^2b^2}{2}}$, so that

In the same way, we can get

Then,

where $M$ is a sufficiently large constant such that

where $A, B, C, E, F$ are defined in the following.

Define a compact bounded subset $U$,

and in the set $\mathbb{R}_+^3\setminus U$, choosing $\epsilon$ small enough such that

Next, six domains are given in the following,

Case 1. If $(S, x_{1}, x_{2})\in U_{1}$,

where

Together with $(3.5)$ and $(3.6)$, we get $L V\leq-1$.

Case 2. If $(S, x_{1}, x_{2})\in U_{2}$,

where

On account of $(3.5)$ and $(3.7)$, we get $LV\leq-1$.

Case 3. If $(S, x_{1}, x_{2})\in U_{3}$,

where

In view of $(3.5)$ and $(3.8)$, we get $LV\leq-1$.

Case 4. If $(S, x_{1}, x_{2})\in U_{4}$,

Considering $(3.5)$ and $(3.9)$, we get $L V\leq-1$.

Case 5. If $(S, x_{1}, x_{2})\in U_{5}$,

where

Combining with $(3.5)$ and $(3.10)$, we get $L V\leq-1$.

Case 6. If $(S, x_{1}, x_{2})\in U_{6}$,

where

Together with $(3.5)$ and $(3.11)$, we get $L V\leq-1$.

Evidently, there exists a sufficiently small $\varepsilon$ such that $LV\leq -1,$ for any $(S, x_1, x_2)\in R_{+}^{3}\setminus U$.

In addition, the diffusion matrix of the system $(1.1)$ is

Then there exists a positive $N = \underset{(S, x_1, x_2)\in U_{\varepsilon }}{\min}\left\{ \frac{\sigma^2a^2S^2}{(m+S)^2}(x_1+x_2)^2, \frac{\sigma^2b^2S^2x_1^2}{(m+S)^2}, \frac{\sigma^2b^2S^2x_2^2}{(m+S)^2}\right\}$ such that

where $(S, x_1, x_2)\in U$, $\xi = (\xi_1, \xi_2, \xi_3) \in R_{+}^{3}$. Obviously this satisfies the conditions in the Lemma 2.1.

Hence, the system ($1.1$) exists a unique stationary distribution.

4.   Numerical simulations

For the theoretical results obtained above, we will conduct further verification through numerical simulation.

Example 1. We use numerical simulations to verify the conditions under which microbes would go extinct. For details about parameter values, see Tables 2 and 3.

And the results are as follows.

The changes in the density of microorganisms and nutrient are shown in Figure 1, the red and green tracks gradually approach zero. Under those parameter settings, $R_1^* = -0.0511 < 0$ and $R_2^* = -0.25 < 0$, the threshold conditions of microbial extinction are met. Therefore, the theoretical results are consistent with the numerical simulation results.

From the point of view of extinction threshold, it is easy to analyze that large noise leads to microbial extinction. Next, we further analyze the effect of noise on the survival of gut microbes. In Example 2, we demonstrate that the enteric microorganism could still survive for a long time under the influence of certain noise.

Example 2. The parameters of composite number are selected to meet the conditions of microbial persistence, that is, $R^* > 0$. The correctness of the conclusion is verified by the following numerical simulation. Parameter settings are shown in Table 4, at this point, $R^* = 0.1096 > 0$, and simulation results are shown in Figure 2.

In addition to the trajectories indicated by the legend in Figure 2, the black dashed line, the black solid line and the purple solid line represent the trajectories of $S$, $x_1$ and $x_2$ of the deterministic model, respectively.

As we can see from the figure, when it is not affected by noise, that is, $\sigma = 0$, the microorganism will survive for a long time. When affected by a certain noise, the trajectories fluctuate near the trajectories of the corresponding deterministic model, and all of them are above the X-axis. In other words, microbes are still persistent.

Example 3. On the basis of Example 2, with other parameters unchanged and $\sigma = 0.7$, the numerical simulation results are shown in Figure 3.

Under the influence of high noise, the originally persistent microorganisms tend to become extinct. Therefore, it can be concluded that the strong noise has adverse effects on the survival of microorganisms. This conclusion is reflected in the threshold values of microbial persistence and extinction.

Example 4. Finally, we verify the existence of stationary ergodic distribution through numerical simulation. The density distributions of microorganisms and nutrient are shown in Figure 4. Set the parameters as Table 5.

In this case, $R = 2.0535 > 1$, and the theoretical results are consistent with the numerical simulation.

5.   Conclusions

This paper proposes a stochastic chemostat model for intestinal microorganisms that considers migration and intraspecific competition. The dynamic behavior of the model is analyzed by ordinary differential theory. By constructing the Lyapunov functions, the threshold conditions of intestinal microbial extinction and persistence in mean are obtained, and the condition for the existence of stationary ergodic distribution is analyzed. Analysis shows that high noise increases the risk of extinction of gut microbes. At the same time, too fast intestinal velocity will also affect the survival of microorganisms. Now, we sum up the main results as follows.

$(1)$ Microorganisms tend to become extinct when the following conditions are met.

$(i)$ If $\frac{bS_0}{m+S_0} < \frac{2}{\sigma^2}$ and $R_{1}^{*} < 0$, then microbials tend to extinction almost surely.

$(ii)$ If $\frac{bS_0}{m+S_0}\geq\frac{2}{\sigma^2}$ and $R_{2}^{*} < 0$, then microbials tend to extinction almost surely.

$(2)$ If $R^* > 0$, then the microorganism is persistent in mean. In addtion, we get

$(3)$ Assume $R > 1$, then system $(1.1)$ exists a unique ergodic stationary distribution.

By overcoming the difficulties in Lyapunov function construction and inequality reduction, the sufficient conditions for microbial extinction, persistence and stationary ergodic distribution were obtained. The threshold of microbial persistence or extinction is influenced by environmental disturbance (noise intensity). The results of numerical simulation also show that a large degree of noise will interfere with the survival of microorganisms, so that the microorganisms that can live for a long time tend to become extinct.

Acknowledgments

This work is supported by the Shandong Provincial Natural Science Foundation of China (No.ZR2019MA003), the Research Fund for the Taishan Scholar Project of Shandong Province of China and the SDUST Innovation Fund for Graduate Students (No.YC20210291).

Conflict of interest

The authors declare no conflicts of interest.