Threshold behaviour of a stochastic SIRS $\mathrm {L\acute{e}vy}$ jump model with saturated incidence and vaccination

1.   Introduction

The method of mathematical modelling has been an important way to study the spreading and controlling of infectious disease. Since the pioneer work of Kermack and McKendrick ^[1], many deterministic mathematical models have been proposed to help us understand the spreading and controlling of infectious disease ^[2,3,4,5]. In ^[6], Lahrouz et al. put forward and analyzed a deterministic SIRS model with general incidence and vaccination. The authors derived the basic reproduction number and provided sufficient conditions for the persistence and extinction of disease. However, the model does not consider the environmental fluctuations. Many researchers have demonstrated that environmental fluctuation has largely affected the spreading of infectious disease ^[7,8]. For human disease, due to the nondeterminacy of person-to-person contacts, the growth and spreading of infectious disease are inherently random.

There are several ways to model the influence of environmental fluctuations by stochastic differential equations ^[9,10]. One of the most important way is to consider system perturbation which arises from the approach in ^[11]. On the basis of the work of Lahrouz et al. ^[6], the SIRS model incorporated by system perturbation can be described by

where the numbers of the susceptible, infectious and recovered individuals at time $t$ are represented by $S(t)$, $I(t)$ and $R(t)$, respectively. $S(t^-)$, $I(t^-)$ and $R(t^-)$ are the left limits of $S(t)$, $I(t)$ and $R(t)$, respectively. $Z(t) = (Z_{1}(t)$, $Z_{2}(t)$, $Z_{2}(t))$ is a 3-dimensional stochastic process which models the random perturbation of the system. The parameters are all assumed to be positive. The biological meanings of the parameters are summarised in the following Table 1.

When $Z_{i}(t) = 0$, system (1.1) is reduced into a deterministic model irrespective of the effects of environmental fluctuations. In fact, the spreading of disease is inevitably disturbed by environmental fluctuations. Usually, white noise ^[12,13] is used to model the impact of environmental noise. However, when there exist large occasionally environmental shocks, such as earthquakes and floods, a stochastic model only with white noise cannot explain these discontinuous perturbations. Since Bao et al. ^[14] firstly proposed that $\mathrm {L\acute{e}vy}$ process should be suitable to describe large occasionally fluctuations, some researchers use non-Gaussian $\mathrm {L\acute{e}vy}$ jump noise to model these discontinuous phenomena ^[15,16]. For the sake of modelling the influence of environmental fluctuations, the form of $Z_{i}(t)$ is represented by

and $Z_{i}(t)$ is a $\mathrm {L\acute{e}vy}$ process and has the following $\mathrm {L\acute{e}vy}$-Khintchine representation

where $u = (u_{1}, u_{2}, u_{3})\in\mathbb{R}^3$ and A is a $3\times3$ diagonal matrix and its diagonal elements are $\sqrt{\sigma_i}$.

In (1.2), defined on a complete probability space $(\Omega, \mathcal{F}, P)$, $B_i(t)$ $(i = 1, 2, 3)$ are standard Brownian motions are independent with $B_i(0) = 0$. The intensities of the white noise is denoted by $\sigma_i$ $(i = 1, 2, 3)$. $\widetilde{N}(\mathrm{d}t, \mathrm{d}z): = N(\mathrm{d}t, \mathrm{d}z)-v(\mathrm{d}z)\mathrm{d}t$ represents the compensated Poisson random measure, where $N(\mathrm{d}t, \mathrm{d}z)$ is the Poisson random measure and its characteristic measure is $v(\mathrm{d}z)$ which is finite $\mathrm {L\acute{e}vy}$ measure. The $\mathrm {L\acute{e}vy}$ measure satisfies $\int_{\mathbb{Y}}\mathrm{min}(|\eta_i(z)|^2, 1)v(\mathrm{d}z) < \infty \quad (i = 1, 2, 3)$, see Theorem 1.2.14 in ^[16]. The effects of jumps is represented by $\eta_i:\mathbb{Y}\to\mathbb{R}$ $(i = 1, 2, 3)$ which is supposed to be continuously differentiable and bounded.

On the basis of the above discussion, the model (1.1) incorporated by linear system perturbation can be described by

The spreading and extinction of infectious disease are two important and interesting topics in the control of infectious disease owing to their theoretical and practical meanings. The aim of the paper is to study the extinction and persistence of stochastic model (1.3). In order to investigate the extinction and persistence of stochastic model (1.3), it is required to estimate the solution of stochastic system (1.3). Usually, for continuous martingales, Zhou and Zhang ^[17] have used Burkholder-Davis-Gundy (BDG) inequality to prove the asymptotic pathwise estimation (see Lemmas 2.1–2.2 in ^[17]) of the solution. However, in this paper, due to the introduction of $\mathrm {L\acute{e}vy}$ jump, the jump stochastic integral process is a discontinuous martingale. Therefore, the Kunita's inequality which is suitable for discontinuous martingales is used to estimate the asymptotic pathwise of the solution for the proposed model, instead of BDG inequality.

This paper is organized as follow: In Section 2, the asymptotic pathwise of the solution is estimated for stochastic model (1.3). In Section 3, sufficient conditions are provided for persistence and extinction of stochastic model (1.3). In Section 4, the paper ends with some discussions and numerical simulations.

2.   Existence and uniqueness of the positive solution

For the jump diffusion coefficients of stochastic model (1.3), we suppose that

$\bf{(H_1)}$ For each $n > 0$, there exists a constant $L_n > 0$ such that $\int_{\mathbb{Y}}|H_i(x, z)-H_i(y, z)|^2v(\mathrm{d}z)\leq L_n|x-y|^2 \ (i = 1, 2, 3)$, where $H_1(x, z) = \eta_1(z)S(t^-), H_2(x, z) = \eta_2(z)I(t^-), H_3(x, z) = \eta_3(z)R(t^-)$ with $|x| \vee |y|\leq n$;

$\bf{(H_2)}$ $\int_{\mathbb{Y}}|\eta_i(z)-\mathrm{ln}(1+\eta_i(z))|v(\mathrm{d}z) < \infty$ for $\eta_i(z) > -1\ (i = 1, 2, 3)$.

Theorem 2.1. Assume $\bf{(H_1)}$ and $\bf{(H_2)}$ hold. Then the stochastic system (1.3) has a unique global solution $(S(t), I(t), R(t))\in{\mathbb{R}_{+}^3}$ for all $t\geq0$ and any given initial value $(S(0), I(0), R(0))\in{\mathbb{R}_{+}^3}$ a.s..

Proof. According to the assumption $\bf{(H_1)}$, the coefficients of the model (1.3) satisfy local Lipschitz conditions. By Theorem 2.1 of ^[14], for any $t\in[0, \tau_{e})$, there is a unique local solution $(S(t), I(t), R(t))\in{\mathbb{R}_{+}^3}$ and $\tau_{e}$ is the explosion time. To illustrate the solution is global, it is necessary to show that $\tau_{e} = +\infty$ a.s.. Define the stopping time as

Set $\tau_{\infty} = \lim\limits_{m\to\infty}\tau_{m}\leq\tau_{e}$ a.s.. It is sufficient to prove $\tau_{\infty} = +\infty$ a.s.. If it does not hold, then there exist constants $T > 0$ and $\epsilon\in(0, 1)$ such that

Consequently, for all $m\ge m_1$, there is an integer $m_1$ such that $P(\tau_{m}\le T)\ge\epsilon$. Define

where $0 < b < \frac{\mu+c}{\beta}$. By Itô's formula with $\mathrm {L\acute{e}vy}$ jump process, we have

where $LV$ is the generating operator of stochastic system (1.3):

Applying the inequality $\eta_i(z)-\mathrm{ln}(1+\eta_i(z))\geq 0$ for $\eta_i(z) > -1$ and assumption $\bf{(H_2)}$, we obtain

where $K$ is a positive constant.

The remaining steps are similar to ^[8], so they are omitted.

3.   Asymptotic pathwise estimation of the stochastic solution

Firstly, we provide some preliminary results for the asymptotic pathwise estimation of the solution of stochastic system (1.3).

For continuous martingales, Zhou and Zhang ^[17] used BDG inequality to prove asymptotic pathwise estimation of the solution. However, in this paper, for discontinuous martingales, Kunita's inequality is needed to prove asymptotic pathwise estimation of the solution and it is different from the BDG inequality for continuous martingales used by Zhou and Zhang ^[17].

Let us first state the Kunita inequality (Theorem 4.4.23 of ^[16]). Consider the stochastic integral for jump process

of the predictable integrand $(H(s, z))_{(s, z)\in\mathbb{R}_{+}\times\mathbb{Y}}$.

Lemma 3.1. Assume that $Y(t)$ is a semi-martingale represented by (3.1) ^[16]. Then there exists a constant $C_p > 0$ such that for any $p\ge2$,

For convenience, denote

Now, let us state the results on the asymptotic pathwise estimation.

Theorem 3.1. Suppose $\bf{(H_1)}$–$\bf{(H_2)}$ hold. If $\mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$ for some $p > 1$, then

Moreover,

Proof. The proof is enlightened by ^[18]. Let $X(t) = S(t)+I(t)+R(t)$. Define

where $p > 0$ is a constant to be determined later. It follows that

where

Set $\rho = \mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda$. It follows from $\rho = \mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$ that

For any $\kappa\in\mathbb{R}$, direct calculation yields that

Integrating from $0$ to $t$, we have

Taking expectations of (3.5) leads to

According to the inequality (3.3), for any $\kappa < \rho p$, we have

where

Therefore, for any $\kappa\in(0, \rho p)$, according to (3.6) and (3.7), we obtain

Consequently,

which indicates that we have $M > 0$ such that

For $k = 1, 2, \ldots$ and sufficiently small $\theta > 0$, by (3.2) and (3.3) that

Direct calculation yields that

Taking expectations of the above inequality, it follows that

where

Now let us estimate $I_1$, $I_2$ and $I_3$. It is easy to see that there exists $c_{11} > 0$ and $c_{12} > 0$ such that

then we have

where $c_1 = p(c_{11}+c_{12})$ is a positive constant. By using BDG inequality, it follows that

where $c_2$ is positive constant. Since $\widetilde{N}(\mathrm{d}t, \mathrm{d}z): = N(\mathrm{d}t, \mathrm{d}z)-v(\mathrm{d}z)\mathrm{d}t$, then,

Therefore, it follows that

Moreover, it is easy to see that we can choose $\theta > 0$ small enough such that

It follows from (3.8) and (3.9) that

Set arbitrary $\epsilon > 0$. For all $k\geq1$, applying Chebyshev's inequality yields

By Borel-Cantelli Lemma (^[19]),

holds for all but finitely many $k$. Then, there exists $k_0(\omega)$ such that whenever $k\geq k_0$, then

Therefore, we have

Letting $\epsilon\to0$ yields

For $p > 1$,

i.e., for any $b\in(0, 1-\frac{1}{p})$, we can find a finite random time $T(\omega)$ such that

It implies

Thus, we get

The positivity of the solution and the equality (3.12) imply

In the following, Kunita's inequality is used to estimate the asymptotic pathwise of the jump stochastic integral process instead of BDG inequality.

Theorem 3.2. Assume $\bf{(H_1)}$–$\bf{(H_2)}$ hold. If there exists some $p > 2$ such that $\mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$, then

Proof. Denote

By Kunita's inequality, there exists $C_{p} > 0$ for any $p\geq2$ such that

where the inequality (3.8) has been used in the last inequality. By the above inequality and (3.10), then

Let $\epsilon > 0$ be arbitrary. Then, by Doob's martingale inequality,

By Borel-Cantelli Lemma (^[19]),

holds for all but finitely many $k$. Therefore, there exists a positive $k_0(\omega)$ such that for all $k > k_0$, we have

Therefore, we obtain

Letting $\epsilon\to0$, we have

Similar to the proof of Theorem 3.1, it follows that

Together with $\liminf\limits_{t\to\infty}\frac{|X_{1}(t)|}{t}\geq0$, this yields

We also obtain

Theorem 3.3. Assume $\bf{(H_1)}$–$\bf{(H_2)}$ hold. If there exists some $p > 1$ such that $\mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$, then

The proof is omitted here because it is similar to that of Lemma 2.2 in ^[18].

4.   Threshold behaviour

When $Z_i(t) = 0$ $(i = 1, 2, 3)$, the stochastic system (1.3) is reduced into a deterministic system investigated by Lahrouz et al. ^[6]. In ^[6], the basic reproduction number $R_0$ is derived represented by

The purpose of this paper is to provide sufficient conditions for the persistence and extinction of disease for stochastic model (1.3). Define

In the following, we will show that $R_{0}^{s}$ completely determines the persistence and extinction of the disease. Throughout this section, always assume

For convenience, denote $\langle x(t)\rangle = \frac{1}{t}\int_{0}^{t}x(s)\mathrm{d}s$.

Theorem 4.1. Suppose $\bf{(H_1)}$–$\bf{(H_3)}$ hold. If $R_{0}^{s} < 1$ and there exists some $p > 2$ such that $\mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$, then

i.e., the disease will die out a.s..

Proof. From stochastic model (1.3),

Integrating the above equality from $0$ to $t$ and dividing both sides by $t$,

where

Similarly, for the first equation of system (1.3),

where

For the third equation in system (1.3), we obtain

where

According to Theorem 3.1–3.3, then

By Itô's formula with jump process,

It follows that

Substituting (4.1)–(4.3) into (4.5) leads to

Besides, define

By assumption $\bf{(H_3)}$, we obtain

Law of large numbers (Theorem 1 of ^[20]) yields that

Employing the law of large number yields $\lim\limits_{t\to\infty}\frac{B_{2}(t)}{t} = 0$, a.s.. This with (4.4) and (4.7) yields

Then

Moreover, from (4.1), (4.3) and (4.4) we obtain

Now, sufficient condition for the persistence of the disease is provided.

Theorem 4.2. Suppose $\bf{(H_1)}$–$\bf{(H_3)}$ hold. If $R_{0}^{s} > 1$ and there exists some $p > 2$ such that $\mu-\frac{p-1}{2}\sigma^2-\frac{1}{p}\lambda > 0$, then

i.e., the disease will persist a.s..

Proof. According to (4.6), we have

where

Theorem 3.1–3.3, the law of large number, (4.4) and (4.7) yield that $\lim\limits_{t\to\infty}\Psi(t) = 0$. It follows from Lemma 2 of Liu and Wang ^[21] that

Meanwhile, from (4.1) and (4.4) we obtain

Consequently, it follows that

By virtue of (4.3), we have

5.   Conclusions

In the paper, we analyze the asymptotic behaviour of a stochastic SIRS system with $\mathrm {L\acute{e}vy}$ process. Due to the introduction of $\mathrm {L\acute{e}vy}$ jump, the jump stochastic integral process is a discontinuous martingale. Thus the Kunita's inequality is used to estimate the asymptotic pathwise of solution for stochastic system (1.3). On this basis, the basic reproduction number for stochastic system (1.3) is defined:

If $R_{0}^{s} < 1$ and some other conditions hold, the disease will die out. If $R_{0}^{s} > 1$ and some other conditions hold, the disease will persist. By comparing $R_{0} = \frac{\beta \Lambda((1-\nu)\mu+\delta)}{\mu(\mu+\delta)(\mu+c+\xi)}$, it is easy to see that $R_{0}^{s} < R_{0}$. This indicates that $\mathrm {L\acute{e}vy}$ noise can suppress the outbreak of infectious disease.

Next, let us make numerical simulations to verify the theoretical results. The initial value is given by $(S(0), I(0), R(0)) = (80, 8, 5)$. The parameters $\nu$, $\Lambda$, $\mu$, $\beta$, $a$, $\delta$, $c$ and $\xi$ are chosen as $0.7346$, $6$, $0.04$, $0.02$, $1$, $0.005$, $0.01$ and $0.8$, respectively. By calculation, we can obtain $R_0 = 1.1023$. This implies that the disease will persist when the system does not have environmental perturbation.

To consider the effect of $\mathrm {L\acute{e}vy}$ process in the spreading of disease, the values of $\sigma_{2}$ are both taken as $0.4$. We first choose $\eta_{2}$ as $0.08$. By calculating the value of $R_{0}^{s}$, we get the value is $1.0046$. It follows from Theorem 4.2 that the disease will persist. The simulated result is shown in Figure 1(a). Then we choose $\eta_{2}$ as $0.2$. The value of $R_{0}^{s}$ is $0.9874$ in this case. By Theorem 4.1, the disease will die out. The simulation result is presented in Figure 1(b). From Figure 1, it can be seen that the larger jump noise induces the extinction of disease, i.e., the $\mathrm {L\acute{e}vy}$ noise suppress the outbreak of infectious disease while the corresponding deterministic model is persist.

Acknowledgments

Y. Zhu is supported by the Postgraduate Research & Practice Innovation Program of Jiangsu Province (No. KYCX22_0391) and the Fundamental Research Funds for the Central Universities (No. 30922010813). L. Wang is supported by the National Natural Science Foundation of China (No. 12001271), the Natural Science Foundation of Jiangsu Province (No. BK20200484) and the Fundamental Research Funds for the Central Universities (No. 30922010813). Z. Qiu is supported by the National Natural Science Foundation of China (No. 12071217).

Conflict of interest

The authors declare that there is no conflict of interest.