A stochastic mathematical model of two different infectious epidemic under vertical transmission

1.   Introduction

In all periods of the development of human society, there are arduous struggles against various infectious diseases ^[1]. What makes matters worse is multiple infectious diseases often exist on human individuals at the same time, and the coordinated and cross-infection between infectious diseases makes the course of the disease more complicated and difficult to deal with. We often refer to this situation as the parallel development of multiple infectious diseases ^[2]. The probability of several infectious diseases in a patient at the same time is related to the environment, susceptible population and human immunity. If the sanitary environment is poor, and there are sewage, excreta, animal and plant residues, mosquitoes and mice everywhere, the probability of $N$ kinds of infectious diseases taking place in parallel will be very high. If people with bad living habits, sanitation workers, medical staff in the infection department and other personnel are frequently exposed to a large number of pathogenic microorganisms, the probability of suffering from $N$ kinds of infection rate is also very high. In addition, people with immune system defects or injuries have a high rate of $N$ infectious diseases ^[3]. Therefore, it is very common for multiple infectious diseases to occur in the same individual. For example, patients with advanced AIDS have a $100\%$ chance of contracting $N$ diversified infectious diseases when exposed to other infection environment besides AIDS ^[4,5].

Recently, the researches on modelling two types of infectious diseases have drawn many applied mathematicians' attention ^{[6,7,8,9,10,11]}. In ^[6], the authors discussed an epidemic model with double hypothesis combined two different transmission mechanisms. In ^[7], a stochastic SIRS epidemic model with two different viruses is built to study the effect of intensities from white noise on each individual, and further discussed the dynamics of threshold around disease-free equilibrium as well as endemic equilibrium. In ^[8], they proposed a new mathematical model with nonlinear incidence rate and two fundamental epidemic hypothesis.

Besides, some researchers have also noticed the phenomena of infectious diseases under vertical transmission and they feel interested in modelling the epidemic with this phenomena ^{[11,12,13,14,15]}. For instance, Ackleh A.S. et al. studied the epidemic model with multi-disease and vertical transmission ^[10,11]. What is vertical transmission infection? It refers to the phenomenon that some infectious diseases pass from the infected mother to the offspring through the blood system and reproductive system, so that the child becomes an infectious disease patient at birth ^[15].

Although the literatures above-mentioned studied the very important phenomenon of multiple infectious diseases taking place in parallel with vertical infection, we note that most of them established deterministic models and did not take into account the impact of various stochastic factors on the dynamic behavior of the system, which is not objective enough. In view of this, stochastic models and stochastic analysis methods are necessary to describe a large number of biological and social phenomena ^{[16,17,18,19]}. In the process of infectious disease modeling, the contact infection incidence of infectious disease is the most noteworthy for two reasons, one is that it can describe the rate at which susceptible becomes infectious, i.e., the incidence can outline the transmission intensity of infectious diseases to some extent; the other is that it can be easily affected by many subjective and objective factors, such as weather, mood, habits and so on. Therefore, it is more objective to establish a stochastic model to describe the infection process of some infectious diseases. In recent decades, many applied mathematicians and epidemiologists have established a series of stochastic epidemic models and achieved great results by introducing random perturbation into deterministic model ^{[16,17,18,19,20,21,22,23,24,25]}. Among these stochastic modeling and analysis literatures, Hattaf K. and coauthors investigated stochastic model with specific functional response and temporary immunity in ^[20], and proposed a delayed stochastic SIR epidemic model with a general incidence rate in ^[21]. Din A. and his coauthors recently also carried out a lot of research on dengue ^[24] and hepatitis B epidemic ^[22,25] as well as some social behaviours including drinking alcohol ^[23] by constructing corresponding stochastic models, they studied the stochastic dynamics such as stationary distribution of extinction and persistence ^{[22,23,24,25]} and also considered the optimal control on hepatitis B epidemic stochastic model ^[22].

The reference ^[11] have discussed the dynamic properties of the system under some different parameters and initial values, on this basis as well as the foregoing mentioned reasons and also motivated by related references mentioned before, in this paper, we will study a stochastic epidemic model with two types of infectious diseases and vertical transmission. By constructing suitable $Lyapunov$ function and using $It\hat{o}$'s formula and $Chebyshev$'s inequality, the dynamic behaviors of our model is analyzed. Next, numerical simulations are carried out to testify our obtained theoretical conclusions. Finally, we propose two interesting related problems to be discussed and further studied in the last section.

2.   Model formulation

In reference ^[11], the authors studied the following differential equations:

where $S(t)$ represents the number of susceptible units at time $t$; $I_{1}(t)$ and $I_{2}(t)$ that represents the number of infected units at time $t$ and $R(t)$ that represents the number of recovered units at time $t$, the initial condition is noted as $S(0) > 0, I_{1}(0) > 0, I_{2}(0) > 0, R(0) > 0$, the total population is thus $N(t) = S(t)+I_{1}(t)+I_{2}(t)+R(t)$. It is assumed that there is a constant quantity of populations entering to the deterministic model with recruitment rate $\Lambda > 0$. There is a vertical transmission of both the two diseases; specifically, the infectious individual gives new birth to a infected individual at the rates $0\leq p_{1}\leq 1$ and $0\leq p_{2}\leq 1$ for the disease $I_{1}$ and $I_{2}$, respectively. Consequently, $p_{1}I_{1}$ and $p_{2}I_{2}$ units will enter into infected compartments $I_{1}$ and $I_{2}$, respectively, and the same quantities are removing from recruitment in the susceptible compartment. The diseases are both transmitted by contact between the individuals in the $S$ compartment and those in $I_{i}(i = 1, 2)$ compartments with nonlinear incidence rate for $I_{1}$ that is described by $\frac{\beta_{1}SI_{1}}{1+I_{1}}$, where $\beta_{1} > 0$ represents the infection power rate, and linear incidence rate for $I_{2}$ that is similarly described by $\beta_{2}SI_{2}$, in which $\beta_{2} > 0$ represents the infection rate. The units in the $I_{1}$ compartment are experiencing death because of the disease at infection death rate $\alpha_{1}\geq0$. They will also recover from the disease and get immunity at a recovered rate $\delta > 0$. The units in the $I_{2}$ compartment are also confronting death for the disease at infection death rate $\alpha_{2}\geq0$. They will also recover from the disease and return back to be susceptible at recovery rate $\gamma > 0$. Furthermore, the individuals in the $R$ compartment will lose immunity rate $0\leq\eta < 1$. There is a natural death rate $\mu > 0$ for every member in the total population. A particularly noteworthy assumption is that both the two diseases cannot be transmitted to the same individual at the same time. Moreover, to insure that the recruitment $\Lambda$ in the susceptible compartment is always positive, the following hypotheses are assumed to be hold always:

The model (2.1) can usually be applied in the environment where the people will probably suffer from two different infectious epidemic, but the same individual will not suffer from these two diseases at the same time, because in one way people around are highly alert to the patients with these two diseases, in another way patients themselves will consciously minimize contact with the outside world. The environment where the two infectious diseases coexist is quite complex. There are some uncertain factors in the process of disease infection and transmission. These factors often play a key role in the dynamics of infectious disease infection and transmission, which can not be ignored. Therefore, in this paper, taking the effect of stochastically fluctuating environment into account, we can embody the fluctuations from the environment in the parameter $\beta_{1}, \beta_{2}$ as follows

where $B_{i}(t)(i = 1, 2)$ is standard Brownian motions with $B(0) = 0$, and with intensity of white noise $\sigma_{i}^{2} > 0(i = 1, 2)$. Therefore, based on deterministic model (2.1), we can derive the corresponding stochastic model which takes the following form:

Next, we introduce some fundamental theory in stochastic differential equation (please see ^[16]).

Let $(\Omega, \mathcal{F}, \mathbf{P})$ be a complete probability space with a filtration $\{\mathcal{F}_{t}\}_{t\geq0}$ which satisfies the usual conditions(i.e., it is right continuous and $\mathcal{F}_{0}$ contains all $\mathbf{P}-null$ sets). Let $B(t)$ be an n-dimensional standard Brownian motion defined on this space.

We consider the n-dimensional stochastic differential equation of $It\hat{o}$ type

with initial value $x(t_{0}) = x_{0}\in\mathbb{R}_{+}^{n} = \{x\in\mathbb{R}_{+}^{n}: x_{i} > 0, \; 1\leq i\leq n\}$. And then we define the differential operator $\mathbf{L}$ along with system (2.4) as

We act the operator $\mathbf{L}$ on a function $V\in C^{2, 1}(\mathbb{R}^{n}\times \mathbb{R}_{+}; \mathbb{R}_{+})$, thus

where

Lemma 1. ($It\hat{0}'s\, \, formula$) (^[20]) Assume $x(t)$ to be an $It\hat{o}$ process on $t\geq 0$ of system (2.4), $V\in C^{2, 1}(\mathbb{R}^{n}\times \mathbb{R}_{+}; \mathbb{R}_{+})$, then the function $V(x(t), t)$ is still an $It\hat{o}$ process with the stochastic differential given by

Lemma 2. ($Chebyshev's \, \, inequality$) (^[16])

3.   Existence as well as uniqueness of positive solution

In this part, in order to show our model is well-posed, by using the Lyapunov analysis method and $It\hat{o}'s$ formula, we will prove that there is a unique global positive solution of system (2.3).

Theorem 3.1. For any initial value $(S(0), I_{1}(0), I_{2}(0), R(0))\in\mathbb{R}_{+}^{4}$, there exists a unique solution $(S(t), I_{1}(t), I_{2}(t), R(t))$ of system (2.3) for $t\geq0$, and the solution will remain in $\mathbb{R}_{+}^{4}$ with probability 1, i.e., $(S(t), I_{1}(t), I_{2}(t), R(t))\in\mathbb{R}_{+}^{4}$ for all $t\geq0$ almost surely.

proof. It can be easily found that the coefficients of the system (2.3) are locally Lipschitz continuous (which means for SDE $dX_{t} = b(t, X_{t})dt+\sigma(t, X_{t})dB_{t}$, $b(t, X_{t})$ and $\sigma(t, X_{t})$ satisfy $|b(t, x)-b(t, y)|+|\sigma(t, x)-\sigma(t, y)|\leq D |x-y|$ for $\forall x, y$ and $t > 0$, here $D > 0$ is a constant ^[20]) for all given initial values $(S(0), I_{1}(0), I_{2}(0), R(0))\in\mathbb{R}_{+}^{4}$, therefore, there must exist a unique local solution $(S(t), I_{1}(t), I_{2}(t), R(t))$ on $t\in[0, \tau_{e})$, where $\tau_{e}$ is called explosion time ^[19]. Next, in order to show this solution is global, we need to prove that $\tau_{e} = \infty$ a.s.. Let $k_{0}\geq0$ be sufficiently large so that all the initial values $S(0), I_{1}(0), I_{2}(0)$ and $R(0)$ can lie within the interval $[\frac{1}{k_{0}}, k_{0}]$. For each integer $k\geq k_{0}$, we define the stopping time

and then we set $\inf\emptyset = \infty$(as usual $\emptyset$ denotes the empty set). Obviously, $\tau_{k}$ is increasing when $k\rightarrow \infty$. We set $\tau_{\infty} = \lim\limits_{k\rightarrow \infty}{\tau_{k}}$, hence $\tau_{\infty}\leq\tau_{e}$ a.s.. To finish the proof, we need to show that $\tau_{\infty} = \infty$ a.s. in the end. We select the proof by contradiction. If this statement is not ture, then there exist two constants $T > 0$ and $\varepsilon\in(0, 1)$ such that

Hence, there is an integer $k_{1}\geq k_{0}$ such that for any $k\geq k_{1}$

We define a $C^{2}$-function V: $\mathbb{R}_{+}^{4}\rightarrow \mathbb{R}_{+}$ by

The nonnegativity of this function can be easily derived from $u-1-\ln u\geq0, \; u > 0$.

Let $k\geq k_{0}$ and $T > 0$ be arbitrary. By virtue of the $It\hat{o}'s$ formula, we obtain

where

where $H$ is obviously a positive constant. Hence,

We integrate both sides of the above inequality from 0 to $\tau_{k}\wedge T$, then we can derive

where $T\wedge \tau_{k} = \min\{\tau_{k}, T\}.$ And then we take the expectations along the two sides of the above equality, we derive the further relation equality as

Let $\Omega_{k} = \{\tau_{k}\leq T\}$ for $k\geq k_{1}$ and from Eq (3.1), we have $\mathbf{P}(\Omega_{k})\geq\varepsilon$. For every $\nu\in\Omega_{k}$, $S(\tau_{k}, \nu)$, $I_{1}(\tau_{k}, \nu)$, $I_{2}(\tau_{k}, \nu)$, $R(\tau_{k}, \nu)$ equals either $k$ or $\frac{1}{k}$; hence $V(S(\tau_{k}, \nu)$, $I_{1}(\tau_{k}, \nu)$, $I_{2}(\tau_{k}, \nu)$, $R(\tau_{k}, \nu))$ is no less than $\min\{k-1-\ln k, \frac{1}{k}-1-\ln \frac{1}{k}\}$.

Then we obtain

where $1_{\Omega_{k}(\nu)}$ is the indicator function of $\Omega_{k}$.

Taking limit process $n\rightarrow \infty$, the contradiction $\infty = V(S(0), I_{1}(0), I_{2}(0), R(0))+HT < \infty$ appears.

Thus, we complete the proof.

4.   Stochastic ultimate boundedness

In this section, we will study the stochastic permanence of system (2.3).

Definition 4.1. (^[20]) The solution $X(t) = (S(t), I_{1}(t), I_{2}(t), R(t))$ of system (2.3) is called stochastically permanent, if for all $\varepsilon\in(0, 1)$ and all initial value $(S(0), I_{1}(0), I_{2}(0), R(0))\in\mathbb{R}_{+}^{4}$, two positive constants $\lambda = \lambda(\varepsilon)$ and $\chi = \chi(\varepsilon)$ can be found to make the solution of system (2.3) $X(t)$ satisfies the properties

Theorem 4.1. For any initial value $(S(0), I_{1}(0), I_{2}(0), R(0))\in\mathbb{R}_{+}^{4}$, the solution $X(t) = (S(t), I_{1}(t), I_{2}(t), R(t))$ satisfies

where $\vartheta$ is an arbitrary positive constant satisfying

in which

proof. Firstly, we define a function $V(S, I_{1}, I_{2}, R) = \frac{1}{S+I_{1}+I_{2}+R}$, $(S(t), I_{1}(t), I_{2}(t), R(t))\in\mathbb{R}_{+}^{4}$; applying $It\hat{o}'s$ formula, we get

Selecting a positive constant $\vartheta$ that satisfies Eq (4.2) and using $It\hat{o}'s$ formula, we can obtain

where $G = \mu V+V^{2}(\alpha_{1}I_{1}+\alpha_{2}I_{2}-\Lambda)+2V^{3}(\frac{\sigma_{1}^{2}S^{2}I_{1}^{2}}{(1+I_{1})^{2}}+\sigma_{2}^{2}S^{2}I_{2}^{2})$.

Since

Therefore,

Let $k > 0$ be small enough and satisfy Eq (4.4), by virtue of $It\hat{o}$'s formula

where

and the definition of $C$ can be found in the statement of the theorem.

Consequently,

We again integrate both sides of the above inequality from $0$ to $t$, and then we can get

Hence, in the end, we obtain

For $(S, I_{1}, I_{2}, R)\in\mathbb{R}_{+}^{4}$, we know that

Thus, we come to the final conclusion

which completes the proof.

Theorem 4.2. The solutions of system (2.3) are stochastically permanent on condition that

proof. From Theorem 4.1, we have $\mathbf{P}\{|X(t)| > \chi\} < \varepsilon$ which implies $\mathbf{P}\{|X(t)|\leq\chi\}\leq1-\varepsilon$. This follows that $\liminf\limits_{t\rightarrow \infty}\mathbf{P}\{|X(t)|\leq\chi\}\geq1-\varepsilon$.

By Theorem 4.1, we get $\limsup\limits_{t\rightarrow \infty}\mathbf{E}(\frac{1}{|X(t)|^{\vartheta}})\leq Q$.

For any $\varepsilon > 0$, let $\lambda = \frac{\varepsilon^{\frac{1}{\vartheta}}}{Q^{\frac{1}{\vartheta}}}$; and then by virtue of $Chebyshev$'s inequality,

Hence,

which means the following assertion is correct.

The proof is thus completed.

Compared with the existence, uniqueness and boundedness of positive solutions, the stochastic permanence condition of the model solution is much more stringent. Specifically, some main parameters, i.e., the recrument rate $\Lambda$, disease-related mortality $\alpha_{1}, \alpha_{2}$, natural mortality $\mu$ and all the perturbed intensities $\sigma_{1}, \sigma_{2}$ should satisfy $\max\{\sigma_{1}^{2}, \sigma_{2}^{2}\} < \Lambda-\frac{\mu+\max\{1, \alpha_{1}, \alpha_{2}\}}{2}$. From a biological point of view, stochastic permanence means that after a long enough time, the probability of population size of all compartment in model (2.3) being too large or too small is very low. Stochastic permanence is equivalent to both random upper bounded and random lower bounded. It implies that the population size in all compartments is relatively stable. It should be further pointed out that when the conditions of the theorem cannot be met, the model will die out, that is, part of the compartments population will eventually die out.

5.   Numerical simulations

In this part, we will use the Milstein method proposed in Higham ^[21] to illustrate our main assertions for the stochastic system (2.3). The Milstein method has been developed from the deterministic Euler method with strong order of convergence $\gamma = \frac{1}{2}$. By adding a correction to the stochastic increment to raise the strong order of Euler method to 1, making Milstein¡¯s method. Therefore, the order of convergence of Milstein method in this paper is 1.

To begin with, let we consider the following discretization equations:

where $\xi_{k}(k = 1, \cdots, n)$ is called the Guassian random variables which Obey Gaussian distribution $N(0, 1)$.

To testify our main results in this paper, according to the biological meanings of relative parameters and initial values, we select the rational parameters values as $\Lambda = 2, \beta_{1} = 0.75, \beta_{2} = 0.3, \gamma = 0.75, p_{1} = 0.01, p_{2} = 0.05, \mu = 0.3, \eta = 0.5, \delta = 0.7, \alpha_{1} = 0.1, \alpha_{2} = 0.6$, and also proper initial value $(S(0), I_{1}(0)$, $I_{2}(0), R(0))$ = (5, 4, 2, 3). Then the related pathwise estimation of the solutions of system (2.3) are presented in Figure 1. Since the intensity of white noise $\sigma_{1}$ and $\sigma_{2}$ are both difficult to measure in actual situations, therefore, we estimate the values of them according to the conditions presented in our theorems with the method in ^[22]. Let $\sigma_{1} = 0.04, \sigma_{2} = 0.06$, the solutions of system (2.3) are stochastically permanent(Figure 1(a)). Let $\sigma_{1} = 0.2, \sigma_{2} = 0.1$ and the condition of Theorem 4.2 is satisfied, we can see that the larger intensity of the white noise will weaken the stability of the system (Figure 1(b)).

When $\sigma_{1} = \sigma_{2} = 0.0$, system (2.3) will be deterministic and its time-series plots shown by Figure 2(a). We select $\sigma_{1} = 1.0, \sigma_{2} = 1.3$ not to satisfy the condition of Theorem 4.2, then the noise can make the population become largely fluctuating. In this condition, the solutions of system (2.3) are not stochastically permanent(please see Figure 2(b)).

6.   Conclusions

In this paper, in order to investigate the effect of stochastic perturbation from the environment on dynamical behaviours of two different infectious epidemic with vertical transmission, we derived a stochastic model based on the deterministic model from literature ^[11]. For the established stochastic model, using the methods and techniques of stochastic analysis, we systematically study the existence and uniqueness of positive solutions, stochastic ultimate boundedness, stochastic persistence and other related issues of the system. We draw the following conclusions:

1) No matter which positive number the initial value takes and independent on the selection of parameters, stochastic system (2.3) always exists a unique positive solution and never explode(i.e., the positive solution is ultimate bounded). It shows our stochastic model is well-posed and biologically significant.

2) Further, when the intensities of stochastic perturbation from environment $\sigma_{1}^{2}, \sigma_{2}^{2}$ are both not rather large with a upper bound described by $\max\{\sigma_{1}^{2}, \sigma_{2}^{2}\} < \Lambda-\frac{\mu+\max\{1, \alpha_{1}, \alpha_{2}\}}{2}$, in any limited time, the population of all compartments in the model will not be very large or very small, that is, it will maintain the so-called stochastic permanence. This means that when the stochastic disturbance intensity is large enough, for example, they are both greater than the population recruitment rate $\Lambda$, the population of at least one compartment will become extinct.

3) Compared to the corresponding deterministic model in ^[11], which showed that both the two disease-free equilibrium and endemic equilibrium point are not only locally asymptotically stable but also global asymptotically stable provided that the parameters in the model satisfied some quantitative relationships, the stochastic model no longer has stability, we can only discuss rather weaker properties, such as permanence and boundedness. This also shows the destructive impact from stochastic disturbance on the properties of the system.

7.   Discussion

In this paper, although we prove the existence and uniqueness, ultimate boundedness and stochastic persistence of the positive solution to the established stochastic model, there are still some more detailed and interesting problems that have not been solved. We will discuss them here in order to stimulate readers' interest, or we will continue to study them in the near future.

In the process of proof, we find that the system satisfies stochastic permanence only when the random disturbance intensity of the system is small enough, in other words, when the random disturbance intensity is large enough, the system will go to extinct. A very important and interesting question is how to control the random disturbance intensity or clarify the relationship between disturbance intensity and system parameters to make the populations in the two disease compartment be extinct while other compartments be stochastic permanent? If this problem can be solved, it will provide a valuable control strategy for effectively controlling infectious diseases in the model.

On the other hand, we also noticed a very interesting phenomenon that individuals will acquire knowledge about this disease when a disease spreads within a community, these knowledge or cognition will form memory and have a subtle impact on people's behavior. Therefore, it will be interesting to study the memory effect on the dynamics of our model by using the new generalized fractional derivative method that we can refer to the document ^[29], where the author put forward a novel definition of fractional derivative taking non-singular kernel in the sense of Caputo. The new generalized fractional derivative method is recently rather popular, we can also refer to the publications ^[30,31,32], using the methods and skills in them, we are sure to depict the memory effect in our model.

Acknowledgments

We would like to thank the funding support from the Natural Science Foundation of Gansu Province, China (No.21JR7RA217, No.20JR5RA460), and National Natural Science Foundation of China (No. 1216053). Besides, we are also very grateful to the two anonymous reviewers since their professional opinions help us improve the quality of the manuscript.

Conflict of interest

We declare that there is no conflict of interest.