Stationary distribution and extinction of a stochastic tuberculosis model

Journal Pre-proof Stationary distribution and extinction of a stochastic tuberculosis model Ming-Zhen Xin, Bin-Guo Wang

PII: DOI: Reference:

S0378-4371(19)32085-0 https://doi.org/10.1016/j.physa.2019.123741 PHYSA 123741

To appear in:

Physica A

Received date : 4 August 2019 Revised date : 11 November 2019 Please cite this article as: M.-Z. Xin and B.-G. Wang, Stationary distribution and extinction of a stochastic tuberculosis model, Physica A (2019), doi: https://doi.org/10.1016/j.physa.2019.123741. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier B.V.

Journal Pre-proof

of

Stationary distribution and extinction of a stochastic tuberculosis model

p ro

Ming-Zhen Xin, Bin-Guo Wang∗†

School of Mathematics and Statistics, Lanzhou University, Lanzhou, Gansu 730000, People’s Republic of China.

Pr e-

Abstract

In this paper, we focus on a stochastic tuberculosis model. We first obtain the existence of a stationary distribution to the positive solutions by stochastic Lyapunov function method. Then we establish sufficient conditions for extinction of the disease. Finally, some examples and numerical simulations are provided to illustrate our theoretical results. Keywords: Stochastic tuberculosis model; extinction; stationary distribution; Lyapunov function

1

Introduction

Jo

urn

al

Tuberculosis (TB) is a widespread infectious disease usually caused by Mycobacterium tuberculosis. Tuberculosis frequently attacks the lungs, but can also affect other parts of the body. According to the report of WHO[1], the best estimate is that 10.0 million people (range, 9.0 to 11.1 million) developed TB disease in 2017: 5.8 million men, 3.2 million women and 1.0 million children. There were cases in all countries and age groups, but overall 90% were adults (aged 15 years), 9% were people living with HIV (72% in Africa) and two thirds were in eight countries: India (27%), China (9%), Indonesia (8%), the Philippines (6%), Pakistan (5%), Nigeria (4%), Bangladesh (4%) and South Africa (3%). These and 22 other countries in WHOs list of 30 high TB burden countries accounted for 87% of the worlds cases. Mathematical models which describe the dynamics of infectious diseases have become useful tools to analyze the spread and control of Tuberculosis and to evaluate the effectiveness of various control strategies in the past few decades[2, 3, 4, 5]. The fast and slow progression was considered earlier by some authors to study the transmission of TB [19, 20, 21]. In 2010, Liu et al. [6] introduced the fast and slow progression based on the real situation of tuberculosis disease. The total population is divided into four distinct epidemiological subclasses of individuals which are the susceptible(S), the exposed(E), the ∗

The research is supported in part by NSF of China (11501269, 11731005) and Natural Science Foundation of Gansu Province 1506RJYA207. † Corresponding author ([email protected])

1

Journal Pre-proof

2

Xin and Wang

p ro

of

infected(I) and the recovered(R). The model is given by a system of ordinary differential equations β(t)SI dS − µS, dt = Λ − N β(t)SI dE = (1 − q) − (k(t) + µ)E, dt N β(t)SI dI (1.1) + k(t)E − (r + d + µ)I, dt = q N dR dt = rI − µR, N = S + E + I + R,

urn

al

Pr e-

where Λ is the recruitment rate, µ is the natural death rate, q is the fraction of fast developing infectious cases, r is the treatment/recovery rate, and d represents the diseaseinduced death rate. All parameters are assumed to be positive. These parameters are positive constants and independent of time t, and q < 1. β(t) and k(t) denote the periodic infection rate and reactivation rate. In the special case of β(t) ≡ β, k(t) ≡ k, R0 = β k (d+r+µ) (q + µ+k (1 − q)) denotes the basic reproduction number. The disease persists in the population and tends to a positive steady state when R0 > 1, the disease dies out eventually when R0 < 1. In the real world, population system is inevitably affected by the environmental white noise, which affects the population dynamics significantly [22, 23]. As a matter of fact, due to the continuous fluctuations in the environment, the parameters involved in the system do not persist at such a steady state and always fluctuate around some average values. Deterministic models ignore the stochastic nature of the rare events that lead to initial resistance generation and spread while the stochastic approach can capture the intrinsic fluctuations when the size of population is small. Stochastic differential equation models play an important role in population dynamics as they can provide some additional degree of realism compared to their corresponding deterministic counterparts. Many authors have considered stochastic biological models and stochastic epidemic models [7, 8, 9, 10, 11]. Motivated by the referred works, in this paper, we adopt the approach used by Imhof and Walcher, and assume that the environmental white noise is proportional to variables S(t), E(t), I(t) and R(t). Now we consider the following stochastic model

Jo

dS = [Λ − βSI N − µS]dt + σ1 SdB1 (t), dE = [(1 − q) βSI N − (k + µ)E]dt + σ2 EdB2 (t), βSI dI = [q N + kE − (r + d + µ)I]dt + σ3 IdB3 (t), dR = [rI − µR]dt + σ4 RdB4 (t), N = S + E + I + R,

(1.2)

where Bi (t)(i = 1, 2, 3, 4) are mutually independent standard Brownian motions defined on a complete probability space (Ω, F, {Ft }t≥0 , P) with a filtration {Ft }t≥0 satisfying the usual conditions (i.e. it is increasing and right continuous while F0 contains all P-null sets), σi2 > 0(i = 1, 2, 3, 4) denote the intensities of white noise, the other parameters in system (1.2) have the same meaning as in system (1.1). In the special case of q = 0, Liu et al. [13] have established conditions for the existence of ergodic stationary distribution and extinction to the model. Throughout this paper, we let R4+ = {(x1 , x2 , x3 , x4 ) ∈ R4 ; xi > 0, i = 1, 2, 3, 4}.

Journal Pre-proof

A stochastic tuberculosis model

3

2

p ro

of

This paper is organized as follows. In Section 2, we show that there exists a unique global positive solution of system (1.2) with initial value (S0 , E0 , I0 , R0 ) ∈ R4+ . In Section 3, by constructing a suitable stochastic Lyapunov function and a rectangular set, we verify the existence of an ergodic stationary distribution of the system (1.2). In Section 4, we establish sufficient conditions for the extinction of the disease. In Section 5, some examples and numerical simulations are provided to illustrate our theoretical results. Finally, some conclusions are presented to close the paper.

Existence and uniqueness of the positive solution

Pr e-

To investigate the dynamical behavior of an epidemic model, the first concerning thing is whether the solution is global and positive. In this section, motivated by the methods mentioned in [15] and we show that there is a unique global positive solution of system (1.2). Consider an l-dimensional stochastic differential equation dz(t) = f (z(t), t)dt + g(z(t), t)dB(t) on t ≥ t0

(2.1)

with initial value z(t0 ) = z0 ∈ Rl . Define the differential operator L associated with Eq.(2.1) by L=

l l X ∂ ∂ 1 X T ∂2 + fk (z, t) + [g (z, t)g(z, t)]kj . ∂t ∂zk 2 ∂zk ∂zj k,j=1

al

k=1

For simplicity, we define a ∨ b = max{a, b} and a ∧ b = min{a, b}.

urn

Theorem 2.1 For any initial value (S0 , E0 , I0 , R0 ) ∈ R4+ , there is a unique solution (S(t), E(t), I(t), R(t)) of system (1.2) on R4 and the solution will remain in R4+ with probability one, that is to say, (S(t), E(t), I(t), R(t)) ∈ R4+ for all t ∈ R almost surely (a.s.).

Jo

Proof. Since the coefficients of system (1.2) satisfy the local Lipschitz condition, then for any initial value (S0 , E0 , I0 , R0 ) ∈ R4+ , there exists a unique local solution (S(t), E(t), I(t), R(t)) on t ∈ [0, τe ), where τe is the explosion time. To show this solution is global, we only need to prove τe = ∞ a.s. Firstly, we show that S(t), E(t), I(t) and R(t) do not explode to infinity in a finite time. To this end, let k0 > 0 be sufficiently large such that S(0), E(0), I(0) and R(0) lie within the interval [ k10 , k0 ]. For each integer k ≥ k0 , define the following stopping time τk = inf{t ∈ [0, τe ) : min{(S(t), E(t), I(t), R(t)} ≤

1 or max{S(t), E(t), I(t), R(t)} ≥ k}. k

Throughout this paper we set inf ∅ = ∞ (as usual ∅ denotes the empty set). Clearly, τk is increasing as k → ∞. Set τ∞ = limk→∞ τk a.s. If this assertion is false, then there is

Journal Pre-proof

4

Xin and Wang

of

a pair of constants T > 0 and  ∈ (0, 1) such that P{τ∞ ≥ T } > . Hence there exists an integer k1 ≥ k0 such that P{τk ≤ T } ≥ , ∀ k ≥ k1 . (2.2) Define a C 2 −function V : R4+ → R+ by

V (S, E, I, R) = (S − 1 − ln S) + (E − 1 − ln E) + (I − 1 − ln I) + (R − 1 − ln R),

p ro

where a is a positive constant to be determined later. The nonnegativity of this function can be obtained from u − 1 − ln u ≥ 0, for any u > 0. Let k ≥ k0 and T > 0 be arbitrary. Applying Itb o, s formula to V (S, E, I, R), we have dV (S, E, I, R) = LV (S, E, I, R)dt + σ1 (S − 1)dB1 (t) + σ2 (E − 1)dB2 (t) + σ3 (I − 1)dB3 (t) + σ4 (R − 1)dB4 (t),

Pr e-

where LV : R4+ → R is defined by

1 βSI 1 1 βSI )[Λ − − µS] + σ12 + (1 − )[(1 − q) − (k + µ)E] S N 2 E N 1 1 βSI 1 + σ22 + (1 − )[ + kE − (r + d + µ)I] + σ32 2 I N 2 1 1 2 + (1 − )[rI − µR] + σ4 R 2 βSI βSI βSI =Λ− − µS + (1 − q) − (µ + k)E + q N N N 1 2 1 2 1 2 1 2 + kE − (µ + d + r)I + rI − µR + aσ1 + σ2 + σ3 + σ4 2 2 2 2 Λ I (1 − q)βSI − −β +µ− + (µ + k) S N EN qβS E I − − k + (µ + d + r) − r + µ N I R (1 − q)βSI ≤ Λ − µN − dI − + 3µ + d + r + k + β EN 1 1 1 1 + σ12 + σ22 + σ32 + σ42 . 2 2 2 2

urn

al

LV (S, E, I, R) = (1 −

Let K := Λ + 3µ + d + r + k + β 12 σ12 + 21 σ22 + 12 σ32 + 21 σ42 , then, we can obtain that

Jo

dV (S, E, I, R) ≤ Kdt + σ1 (S − 1)dB1 (t) + σ2 (E − 1)dB2 (t) + σ3 (I − 1)dB3 (t) + σ4 (R − 1)dB4 (t).

(2.3)

Integrating both side (2.3) from 0 to τk ∧ T and then taking the expectations yield EV (S(τk ∧ T ), E(τk ∧ T ), I(τk ∧ T ), R(τk ∧ T )) ≤ V (S0 , E0 , I0 , R0 ) + KE(τk ∧ T ).

Therefore, EV (S(τk ∧ T ), E(τk ∧ T ), I(τk ∧ T ), R(τk ∧ T )) ≤ V (S0 , E0 , I0 , R0 ) + KT.

(2.4)

Journal Pre-proof

A stochastic tuberculosis model

5

of

Set Ωk = {τk ≤ T } for k ≥ k1 and by virtue of (2.2), we obtain P(Ωk ) ≥ . Note that for every ω ∈ Ωk , there is S(τk , ω),I(τk , ω),E(τk , ω) or R(τk , ω) equals either k1 or k. Consequently, V (S(τk , ω), E(τk , ω), I(τk , ω), R(τk , ω)) is no less than either k − 1 − ln k or 1 1 1 k − 1 − ln k = k − 1 + ln k. Hence, we can get

p ro

1 (S(τk , ω), E(τk , ω), I(τk , ω), R(τk , ω)) ≤ [k − 1 − ln k] ∧ [ − 1 + ln k]. k It follows from (2.4) that

V (S0 , E0 , I0 , R0 ) + KT ≥ E[IΩk (ω)V (S(τk , ω), E(τk , ω), I(τk , ω), R(τk , ω))] 1 ≥ [k − 1 − ln k] ∧ [ − 1 + ln k], k where IΩk is the indicator function of Ωk . Letting k → ∞, then we have

(2.5)

Pr e-

∞ > V (S0 , E0 , I0 , R0 )) + KT = ∞,

which yields the contradiction. Thus, we get τ∞ = ∞. This means that S(t), E(t), I(t) and R(t) will not explode in a finite time almost surely. 

3

Stationary distribution

urn

al

When considering epidemic dynamical systems, we are also interested to know when the disease will persist and prevail in a population. In the deterministic models, it can be solved by proving that the endemic equilibrium of the corresponding model is a global attractor or globally asymptotically stable. But for system (1.2), there exists no endemic equilibrium. In this section, based on the theory of Hasminskii [18], we prove that there is a stationary distribution which reveals that the disease will be persistent in the mean. Here we present some theory about the stationary distribution which is introduced in [18]. Let X(t) be a homogeneous Markov process in Ed (Ed denotes d-dimensional Euclidean space), and be described by the following stochastic differential equation dX(t) = b(X)dt +

The diffusion matrix is defined as follows

k X

gr (X)dBr (t).

r=1

Jo

A(x) = (ai,j (x)), ai,j (x) =

k X

gri (x)grj (x).

r=1

Lemma 3.1 The Markov process X(t) has a unique ergodic stationary distribution µ(·) if there exists a bound D ⊂ E d with regular P boundary Γ and the following conditions: (1) there is a positive number M such that di,j=1 ai,j (x)ξi ξj ≥ M kξk2 , x ∈ D, ξ ∈ Rd . (2) there exists a nonnegative C 2 -function V such that LV is negative for any Ed \D. Then Z Z 1 T Px { lim f (X(t))dt = f (x)µ(dx)} = 1 T →+∞ T 0 Ed for all x ∈ Ed , where f is a function integrable with respect to the measure µ.

Journal Pre-proof

6

Xin and Wang Define a parameter

(µ +

σ12 2 )(µ

βk(1 − q)µ +k+

σ22 2 )(µ

+d+r+

σ32 2 )

qβµ

+ (µ +

σ12 2 )(µ

+d+r+

of

R0s =

σ32 2 )

.

p ro

4 , system (1.2) admits a unique Theorem 3.2 Assume that R0s > 1, (S0 , E0 , I0 , R0 ) ∈ R+ stationary distribution µ(·) and it has the ergodic property.

Choose M = 4 X

i,j=1

min

(S,E,I,R)∈ Dδ ⊂R4+

Pr e-

Proof. Theorem 2.1 tells us that for any initial value (S0 , E0 , I0 , R0 ) ∈ R4+ , there is a unique global solution (S(t), E(t), I(t), R(t)) ∈ R4+ . In order to prove this Theorem, we only need to validate condition (1) and (2) in Lemma 3.1 . Now we verify the condition (1). The diffusion matrix of system (1.2) is given by   2 2 σ1 S 0 0 0  0 σ22 E 2 0 0  . (3.1) Hx =  2 2  0 0 σ3 I 0  0 0 0 σ42 R2 {σ12 S 2 , σ22 E 2 , σ32 I 2 , σ42 R2 }, we can get that

aij (S, E, I, R)ξi ξj = σ12 S 2 ξ12 + σ22 E 2 ξ22 + σ32 I 2 ξ32 + σ42 R2 ξ42 ≥ M kξk2 , f or any (S, E, I, R) ∈ Dδ , ξ = (ξ1 , ξ2 , ξ3 , ξ4 ) ∈ R4 .

al

Then the condition (1) in Lemma 3.1 is satisfied. Construct a C 2 -function Q : R4+ → R in the following form Q(S, E, I, R) = M [(−c1 ln S − c2 ln E − c3 ln I + c4 (S + E + I + R) − c5 ln S

where

βk(1 − q)µΛ

(µ +

σ12 2 2 ) (µ

+k+

Jo

c1 =

urn

+ c6 (S + E + I + R)] − ln S − ln E − ln R + (S + E + I + R) 1 + (S + E + I + R)(θ+1) , θ+1

σ22 2 )

, c2 =

c5 =

βk(1 − q)µΛ (µ + k + qβµΛ

(µ +

σ12 2 2 )

σ22 2 2 ) (µ

, c6 =

+

σ12 2 )

qβµ µ+

σ12 2

, c3 = Λ, c4 =

βk(1 − q)µ (µ + k +

σ22 2 )

.

M > 0 satisfies the following condition −M Λ(µ + d + r + where

σ32 )(R0s − 1) + G < −2. 2

(3.2)

Journal Pre-proof

A stochastic tuberculosis model

sup 4 (S,E,I,R)∈R+

{3µ +

σ12 σ2 σ2 + 3µ + k + Λ + 2 + 4 + B 2 2 2

of

G=

7

It is easy to check that lim inf

k→∞,(S,E,I,R)∈R4+ \Uk

p ro

1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 )} < ∞. 2 2

Q(S, E, I, R) = ∞.

Pr e-

where Uk = ( k1 , k) × ( k1 , k) × ( k1 , k) × ( k1 , k). Furthermore, Q(S, E, I, R) is a continuous function. Hence, Q(S, E, I, R) must have a minimum point (S0 , E0 , I0 , R0 ) in the interior of R4+ . Then we define a nonnegative C 2 -function V : R4+ → R4+ as follows V (S, E, I, R) = Q(S, E, I, R) − Q(S0 , E0 , I0 , R0 ). Let V1 (S, E, I, R) = −c1 ln S−c2 ln E−c3 ln I +c4 (S+E+I +R)−c5 ln S+c6 (S+E+I +R), then c1 Λ c1 I σ2 βc2 (1 − q)SI σ2 +β + c1 (µ + 1 ) − + c2 (µ + k + 2 ) S N 2 EN 2 qβS E σ32 Λ − c3 − c3 + c3 (µ + d + r + ) + c4 Λ − c4 µN − c4 dI − c5 N I 2 S I σ12 + βc5 + c5 (µ + ) + c6 Λ − c6 µN − c6 dI N 2 Λ β(1 − q)SI E Λ qβS = −c1 − c2 − c3 − c4 µN − c5 − c6 µN − c3 S EN I S N σ12 σ22 σ32 + c1 (µ + ) + c2 (µ + k + ) + c3 (µ + d + r + ) + c4 Λ 2 2 2 I σ12 + c5 (µ + ) + c6 Λ + β(c1 + c5 ) 2 N 1 1 σ2 ≤ −4(c1 βc2 (1 − q)c3 c4 Λµ) 4 − 3(c6 µc5 Λc3 qβ) 3 + c1 (µ + 1 ) 2 2 2 2 σ σ σ + c2 (µ + k + 2 ) + c3 (µ + d + r + 3 ) + c4 Λ + c5 (µ + 1 ) + c6 Λ 2 2 2 βk(1 − q)µΛ qβΛ σ32 I =− − + (µ + d + r + ) + β(c + c ) 1 5 2 2 2 σ σ σ 2 N (µ + 1 )(µ + k + 2 ) µ + 1

Jo

urn

al

LV1 (S, E, I, R) = −

2

2

2

σ2 I = −Λ(µ + d + r + 3 )(R0s − 1) + β(c1 + c5 ) . 2 N

Let V2 (S, E, I, R) = − ln S − ln E − ln R, then LV2 (S, E, I, R) = −

I σ 2 (1 − q)βSI σ2 I σ2 Λ +β +µ+ 1 − +µ+k + 2 −r +µ+ 4. S N 2 EN 2 R 2

Journal Pre-proof

8

Xin and Wang

Let V3 (S, E, I, R) = S + E + I + R, then

Let V4 (S, E, I, R) =

1 θ+1 (S

+ E + I + R)θ+1 , then

of

LV3 (S, E, I, R) = Λ − µ(S + E + I + R) − dI.

where B=

sup 4 (S,E,I,R)∈R+

Pr e-

p ro

LV4 (S, E, I, R) = (S + E + I + R)θ [Λ − µ(S + E + I + R) − dI] 1 + θ(S + E + I + R)(θ−1) × (σ12 S 2 + σ22 E 2 + σ32 I 2 + σ42 R2 ) 2 ≤ (S + E + I + R)θ [Λ − µ(S + E + I + R)] 1 + θ(S + E + I + R)θ+1 (σ12 ∨ σ22 ∨ σ32 ∨ σ42 ) 2 1 = Λ(S + E + I + R)θ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S + E + I + R)θ+1 2 1 1 2 2 2 ≤ B − [µ − θ(σ1 ∨ σ2 ∨ σ3 ∨ σ42 )](S + E + I + R)θ+1 2 2 1 1 ≤ B − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 ) 2 2 {Λ(S + E + I + R)θ

1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )(S + E + I + R)θ+1 ]} < ∞. 2 2

al

Thus,

I Λ I σ2 σ32 )(R0s − 1) + β(c1 + c5 ) ] − + β + µ + 1 2 N S N 2 (1 − q)βSI σ42 I σ32 − +µ+k+ −r +µ+ + Λ − µ(S + E + I + R) − dI EN 2 R 2 1 1 + B − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 ) 2 2 σ2 I ≤ −M Λ(µ + d + r + 3 )(R0s − 1) + β[M (c1 + c5 ) + 1] 2 N µ(1 − q)βSI 1 − 2( )2 E 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 ) 2 2 I Λ σ2 σ2 σ2 − r − − dI + µ + 1 + 3µ + k + Λ + 2 + 4 + B R S 2 2 2

Jo

urn

LV (S, E, I, R) ≤ M [−Λ(µ + d + r +

Now we are in the position to construct a compact subset D such that the condition (2) holds. Define the bounded closed set D = {1 ≤ S ≤

1 1 1 1 , 2 ≤ I ≤ , 3 ≤ E ≤ , 4 ≤ R ≤ }, 1 2 3 4

Journal Pre-proof

A stochastic tuberculosis model

9

For the convenience, we divide R4+ \D into eight domains, D3 = {(S, E, I, R) ∈ R4+ : S ≥ 1 , I ≥ 2 , 0 < E < 3 },

of

D1 = {(S, E, I, R) ∈ R4+ : 0 < S < 1 }, D2 = {(S, E, I, R) ∈ R4+ : 0 < I < 2 , S ≥ 1 },

p ro

D4 = {(S, E, I, R) ∈ R4+ : 0 < R < 4 , I ≥ 2 }, 1 1 D5 = {(S, E, I, R) ∈ R4+ : S > }, D6 = {(S, E, I, R) ∈ R4+ : I > }, 1 2 1 1 D7 = {(S, E, I, R) ∈ R4+ : E > }, D8 = {(S, E, I, R) ∈ R4+ : R > }. 3 4

where i > 0(i = 1, 2, 3, 4) are sufficiently small constants satisfying the following conditions

2 = 21 , −M Λ(µ + d + r +

Λ + F ≤ −1. 1

σ32 )(R0s − 1) + β[M (c1 + c3 ) + 1]1 + G < −1. 2

3 = 41 , −2(

(3.5)

1 + F < −1. 1

(3.6)

(3.7)

1 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 2(θ+1) + L < −1. 4 2 

(3.8)

1 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 3(θ+1) + L < −1. 4 2 

(3.9)

1 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 4(θ+1) + L < −1. 4 2 

(3.10)

urn

al

1 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] θ+1 + L < −1. 4 2 1

Jo F =

(3.4)

µ(1 − q)β ) + F < −1. 1

4 = 31 , −γ

where

(3.3)

Pr e-

−

sup

4 (S,E,I,R)∈R+

1

1

1

{β[M (c1 + c5 ) + 1]

I σ2 σ2 σ2 + µ + 1 + 3µ + k + Λ + 2 + 4 + B N 2 2 2

1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 )} < ∞. 2 2

where L=

sup 4 (S,E,I,R)∈R+

{β[M (c1 + c5 ) + 1]

I σ2 σ2 σ2 + µ + 1 + 3µ + k + Λ + 2 + 4 + B.} N 2 2 2

Journal Pre-proof

10

Xin and Wang

of

Next, we will show that LV (S, E, I, R) ≤ −1 on R4+ \D, which is equivalent to proving it on the above eight domains. Case 1. If (S, E, I, R) ∈ D1 , then we can get that Λ I σ2 σ2 σ2 + β[M (c1 + c5 ) + 1] + µ + 1 + 3µ + k + Λ + 2 + 4 + B S N 2 2 2 1 1 2 2 2 2 θ+1 θ+1 θ+1 θ+1 − [µ − θ(σ1 ∨ σ2 ∨ σ3 ∨ σ4 )](S +E +I +R ) 2 2 Λ Λ ≤ − + F ≤ − + F. S 1

p ro

LV (S, E, I, R) ≤ −

According to (3.3)

LV (S, E, I, R) ≤ −1 f or any (S, E, I, R) ∈ D1 .

Pr e-

Case 2. If (S, E, I, R) ∈ D2 , then we can get that LV (S, E, I, R) ≤ −M Λ(µ + d + r +

σ32 I σ2 )(R0s − 1) + β[M (c1 + c3 ) + 1] + 3µ + 1 + 3µ 2 N 2

σ22 σ42 + +B 2 2 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 )} 2 2 σ2 I ≤ −M Λ(µ + d + r + 3 )(R0s − 1) + β[M (c1 + c3 ) + 1] + G 2 S σ32 2 ≤ −M Λ(µ + d + r + )(R0s − 1) + β[M (c1 + c3 ) + 1] + G. 2 1

According to (3.4)

al

+Λ+

urn

LV (S, E, I, R) ≤ −M (R0s − 1) + β[M (c1 + c3 ) + 1]2 + G ≤ −1 f or any (S, E, I, R) ∈ D2 . Case 3. If (S, E, I, R) ∈ D3 , then we can get that µ(1 − q)βSI 1 I ) 2 + β[M (c1 + c5 ) + 1] E N 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 ) 2 2 σ12 σ2 σ2 +µ+ + 3µ + k + Λ + 2 + 4 + B 2 2 2 µ(1 − q)βSI 1 µ(1 − q)β1 2 1 ≤ −2( ) 2 + F ≤ −2( ) 2 + F. E 3

Jo

LV (S, E, I, R) ≤ −2(

According to (3.5) LV (S, E, I, R) ≤ −2(

µ(1 − q)β 1 ) 2 + F < −1 f or any (S, E, I, R) ∈ D3 . 1

Journal Pre-proof

A stochastic tuberculosis model Case 4. If (S, E, I, R) ∈ D4 , then we can get that I I µ(1 − q)βSI 1 β[M (c1 + c5 ) + 1] − 2( )2 R N E 1 1 − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )](S θ+1 + E θ+1 + I θ+1 + Rθ+1 ) 2 2 σ12 σ2 σ2 +µ+ + 3µ + k + Λ + 2 + 4 + B 2 2 2 2 I ≤ −r + F ≤ −γ + F. R 4

According to (3.6)

γ + F < −1 f or any (S, E, I, R) ∈ D4 . ε1

Pr e-

LV (S, E, I, R) ≤ −

p ro

of

LV (S, E, I, R) ≤ −r

Case 5. If (S, E, I, R) ∈ D5 , then we can get that

According to (3.7)

al

1 1 LV( S, E, I, R) ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )]S θ+1 2 2 I σ2 σ2 σ2 + β[M (c1 + c5 ) + 1] + µ + 1 + 3µ + k + Λ + 2 + 4 + B N 2 2 2 1 1 2 2 2 2 θ+1 +L ≤ − [µ − θ(σ1 ∨ σ2 ∨ σ3 ∨ σ4 )]S 2 2 1 1 1 ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] θ+1 + L. 2 2 1

urn

LV (S, E, I, R) < −1 f or any (S, E, I, R) ∈ D5 .

Case 6. If (S, E, I, R) ∈ D6 , then we can get that

Jo

1 1 LV( S, E, I, R) ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )]E θ+1 2 2 σ2 σ2 σ2 I + β[M (c1 + c5 ) + 1] + µ + 1 + 3µ + k + Λ + 2 + 4 + B N 2 2 2 1 1 ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )]E θ+1 + L 2 2 1 1 1 ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 2(θ+1) + L. 2 2  1

According to (3.8) LV (S, E, I, R) < −1 f or any (S, E, I, R) ∈ D6 .

11

Journal Pre-proof

12

Xin and Wang

Case 7. If (S, E, I, R) ∈ D7 , then we can get that

p ro

of

1 1 LV (S, E, I, R) ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )]I θ+1 2 2 I σ2 σ2 σ2 + β[M (c1 + c5 ) + 1] + µ + 1 + 3µ + k + Λ + 2 + 4 + B N 2 2 2 1 1 2 2 2 2 θ+1 +L ≤ − [µ − θ(σ1 ∨ σ2 ∨ σ3 ∨ σ4 )]I 2 2 1 1 1 ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 3(θ+1) + L. 2 2  1

According to (3.9)

Pr e-

LV (S, E, I, R) < −1 f or any (S, E, I, R) ∈ D7 . Case 8. If (S, E, I, R) ∈ D8 , then we can get that

1

urn

According to (3.10)

al

1 1 LV (S, E, I, R) ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )]Rθ+1 2 2 I σ2 σ2 σ2 + β[M (c1 + c5 ) + 1] + µ + 1 + 3µ + k + Λ + 2 + 4 + B N 2 2 2 1 1 2 2 2 2 θ+1 ≤ − [µ − θ(σ1 ∨ σ2 ∨ σ3 ∨ σ4 )]R +L 2 2 1 1 1 ≤ − [µ − θ(σ12 ∨ σ22 ∨ σ32 ∨ σ42 )] 4(θ+1) + L. 2 2 

LV (S, E, I, R) < −1 f or any (S, E, I, R) ∈ D8 .

Thus, we can see that for a sufficiently small i , LV (S, E, I, R) ≤ −1, ∀(S, E, I, R) ∈ R4+ \D.

Jo

Consequently, condition (2) in Lemma 3.1 is satisfied. According to Lemma 3.1, we can obtain that system (1.2) is ergodic and has a unique stationary distribution.  Remark 3.3 Theorem 3.2 reveals that system (1.2) has a unique ergodic stationary disβk(1−q)µ qβµ tribution µ(·), if R0s = + > 1. The expression σ2 σ2 σ2 σ2 σ2 R0s

(µ+

1 2

)(µ+k+

2 2

)(µ+d+r+

3 2

)

(µ+

1 2

)(µ+d+r+

3 2

)

of is coincide with the number the threshold R0 of the deterministic autonomous system (1.1) if we don’t consider the white nosie. This shows that we generalize the results of the deterministic autonomous system.

Journal Pre-proof

A stochastic tuberculosis model

4

13

Extinction of the disease

β(qµ + k) (µ +

k)(( β(qµ+k) (µ+k)

+

σ22 2 )

∧

( µ+d+r 2

β(qµ+k) 2(µ+k) )

+

∧ (µ + d + r +

σ32 2 ))

.

p ro

R0E =

of

In this section, we shall consider the extinction of the disease. For convenience and simplicity in the following analysis, define the following notations:

Lemma 4.1 Let (S(t), E(t), I(t), R(t)) be the solution of system (1.2) with any initial value (S0 , E0 , I0 , R0 ) ∈ R4+ . Then

Moreover, if µ >

S(t) E(t) I(t) R(t) = 0, lim = 0, lim = 0, lim = 0 a.s. t→∞ t→∞ t→∞ t t t t

σ12 ∨σ22 ∨σ32 ∨σ42 , 2

1 lim t→∞ t 1 t→∞ t lim

Z

t

0

Z

then

Pr e-

lim

t→∞

1 S(u)dB1 (u) = 0, lim t→∞ t

t

1 t→∞ t

I(u)dB3 (u) = 0, lim

0

Z

t

E(u)dB2 (u) = 0 a.s.

0

Z

t

R(u)dB4 (u) = 0 a.s.

0

The proof is similar to [14, Lemma 1 and Lemma 2] and hence is omitted.

al

Lemma 4.2 [15, Theorem 1.3.4] Let M = {Mt }t≤0 be a real-valued continuous local martingale vanishing at t = 0.Then Mt = 0 a.s. t→∞ hM, M it

lim hM, M it = ∞ a.s. ⇒ lim

and also

urn

t→∞

lim sup

Mt hM, M it < ∞ a.s. ⇒ lim = 0 a.s. t→∞ t t

Jo

Theorem 4.3 Let (S(t), E(t), I(t), R(t)) be the solution of system (1.2) with any initial σ 2 ∨σ 2 ∨σ 2 ∨σ 2 value (S0 , E0 , I0 , R0 ) ∈ R4+ . If R0E < 1 and µ > 1 2 2 3 4 , then the solution (S(t), E(t), I(t), R(t)) of the system (1.2) satisfies. lim sup t→∞

ln(kE(t) + (µ + k)I(t)) ≤ K(R0E − 1) < 0 a.s. t lim R(t) = 0 a.s.

t→∞

and the distribution of S(t)converges weakly to the measure which has the density f ∗ (x) = Cx−2 x

−2µ 2 σ1

−

e

2Λ 2x σ1

,

Journal Pre-proof

14

Xin and Wang

where Γ(

2µ + σ12 −1 )] , σ12

β(qµ + k) σ22 µ + d + r β(qµ + k) σ2 + )∧( + ) ∧ (µ + d + r + 3 )). (µ + k) 2 2 2(µ + k) 2

p ro

K = ((

2 2µ+σ1 2 σ1

of

2Λ − C = [( 2 ) σ1

Proof. Consider the following stochastic differential equation du = [Λ − µu]dt + σ1 udB1 (t).

(4.1)

It is easy to check that Eq. (4.1) has a stationary solution u(t) (see[16] )which has the Cx−2 x

density = theorem it follows that

−2µ 2 σ1

e

−

2Λ 2x σ1

1 lim t→∞ t

2 2µ+σ1 2 σ1

Pr e-

f ∗ (x)

,where C =

Z

t

S(s)ds =

0

[( 2Λ ) σ12

Z

∞

−

Γ(

2µ+σ12 −1 )] . σ12

xf ∗ (x)dx a.s.

From the ergodic

(4.2)

0

From SDE(4.1), by direct calculation, we obtain ∞

0

xf ∗ (x)dx = Cσ12

urn

∞

x

− 2µ2 −1 − σ1

e

2Λ 2x σ1

dx

0

Z

∞

2Λ − σ2µ2 −1 σ2µ2 −1 t 2Λ ) 1 t 1 e ( 2 )dt σ12 σ1 0 2µ 2 2µ σ 2 = Cσ12 ( 1 ) σ1 Γ( 2 ) 2Λ σ1 =

Jo

Z

al

Z

Cσ12

(

2µ 2Λ Γ( σ12 ) = 2 2µ σ1 Γ( 2 ) + 1 σ1

σ12

2Λ σ12 2µ Λ = . µ =

Then, let X(t) be the solution of SDE Eq. (4.1) with the initial value u(0) = S(0) > 0, then using the comparison theorem of 1−dimensional stochastic differential equation (see[17]), we have S(t) ≤ u(t) a.s.

(4.3)

Journal Pre-proof

A stochastic tuberculosis model

15

On the other hand, let C(t) = kE(t) + (µ + k)I(t), By Itˆ o’s formula, we get SI kβ SI σ 2 k 2 E 2 + σ32 (µ + k)2 I 2 N + µqβ N − (µ + k)(µ + d + r) − 2 kE + (µ + k)I 2(kE + (µ + k)I)2 σ2 kE σ3 (µ + k)I + dB2 (t) + dB3 (t) kE + (µ + k)I kE + (µ + k)I

−

kβ(qµ+k) µ+k E

kβI + µqβI + σ22 k 2 E 2

σ32 (µ

k)2 I 2

+ + 2(kE + (µ + k)I)2

−

kβ(qµ+k) µ+k E

− (µ + k)(µ + d + r)I

p ro

≤

of

d ln C(t) =

kE + (µ + k)I σ2 kE σ3 (µ + k)I + dB2 (t) + dB3 (t) kE + (µ + k)I kE + (µ + k)I

kβ(qu+k) β(qµ + k) ( µ+k E + (µ + k)(µ + d + r)I) k 2 σ22 E 2 + (k + µ)2 σ32 I 2 − − µ+k kE + (µ + K)I 2(kE + (µ + k)I)2 σ2 kE σ3 (µ + k)I + dB2 (t) + dB3 (t) kE + (µ + k)I kE + (µ + k)I

=−

k 2 ( β(qµ+k) µ+k +

Pr e-

=

σ22 2 2 )E

+ 2k(µ + k)( µ+d+r + 2

β(qµ+k) 2(µ+k) )EI

(kE + (µ + k)I)2 σ2

al

(µ + k)2 (µ + d + r + 23 )I 2 β(qµ + k) σ2 kE − + + dB2 (t) 2 (kE + (µ + k)I) µ+k kE + (µ + k)I σ3 (µ + k)I + dB3 (t) kE + (µ + k)I β(qµ + k) β(qµ + k) σ22 µ + d + r β(qµ + k) ≤ −( + )∧( + ) µ+k (µ + k) 2 2 2(µ + k) σ2 σ2 kE σ3 (µ + k)I ∧ (µ + d + r + 3 ) + dB2 (t) + dB3 (t). 2 kE + (µ + k)I kE + (µ + k)I

urn

(4.4)

Integrating (4.4) from 0 to t and then dividing by t on both sides, we can get

Jo

ln C(t) β(qµ + k) σ22 µ + d + r β(qµ + k) σ2 ≤ −( + )∧( + ) ∧ (µ + d + r + 3 ) t (µ + k) 2 2 2(µ + k) 2 Z t Z t β(qµ + k) σ2 kE σ3 (µ + k)I + + dB2 (s) + dB3 (s). µ+k 0 kE + (µ + k)I 0 kE + (µ + k)I Rt σ2 kE(s) Let M (t) = 0 kE(s)+(µ+k)I(s) dB2 (s). Obviously, hM, M it lim sup = lim sup t t→∞ t→∞

By using the Lemma 4.2, we have Rt lim

t→∞

Rt

σ22 k2 E(s)2 0 (kE(s)+(µ+k)I(s))2 ds

t

σ2 kE(s) 0 kE(s)+(µ+k)I(s) dB2 (s)

t

= 0.

(4.5)

≤ σ22 a.s.

(4.6)

Journal Pre-proof

16

Xin and Wang

Similarly, we can get σ3 (µ+k)I(s) 0 kE(s)+(µ+k)I(s) dB3 (s)

t

= 0.

of

lim

t→∞

Rt

lim sup t→∞

p ro

From (4.5), (4.6) and (4.7), we can get

ln C(t) β(qµ + k) β(qµ + k) σ22 µ + d + r β(qµ + k) ≤ −( + )∧( + ) t µ+k (µ + k) 2 2 2(µ + k) σ2 ∧ (µ + d + r + 3 ) 2 = K(R0E − 1).

Pr e-

Since R0E < 1, from (4.8), we can get

(4.7)

(4.8)

lim E(t) = 0 a.s. , lim I(t) = 0 a.s.

t→∞

t→∞

From system (1.2) we can see that limt→∞ R(t) = 0 a.s. when limt→∞ I(t) = 0 a.s. Thus for any small  > 0 there exists t0 and a set Ω ⊂ Ω such that P(Ω ) > 1 −  and βSI ≤ S for t ≥ t0 and ω ∈ Ω . S − µS)dt + σ1 SdB1 (t) ≤ dS ≤ (Λ − µS)dt + σ1 SdB1 (t). S+E+I +R

(4.9)

al

(Λ −

urn

According to (4.2), (4.3) and (4.9), the distribution of the process S(t) converges to the measure with the density f ∗ .  Remark 4.4 If we take q = 0, then R0E = ˆs = et al. [13] point out that if R 0

βk σ2 σ2 βk βk (µ+k)[( µ+k + 22 )∧( µ+d+r + 2(µ+k) )∧(µ+d+r+ 23 )] 2

βk(µ+k) σ2 k2 σ 2 1 [(µ+k)2 (d+µ+r+ 23 )∧ 2 2 ] 2

. Liu

< 1, the disease will die out. It

ˆ s ≤ RE , we improve the result in the this case. is easy to check that R 0 0

5

Jo

β(qµ+k) In the determine model, if R0 = (µ+k)(µ+d+r) > 1(that is β(qµ+k) > (mu + d + r)), µ+k E this disease will persist. But from this expression of R0 , we can find that if the density of noise is enough large that the disease maybe die out.

Numerical simulations

In this section,we give some examples to illustrate the obtained theoretical results. We illustrate our findings by the method developed in [12]. Consider the corresponding

Journal Pre-proof

A stochastic tuberculosis model

Values 1531845 0.001167 0.05 0.005 0.425 0.00017088 728430000 52336000 231976 60458024 1312480000

Reference [6] [6] [6] [6] [6] [6] [6] [6] [6] [6] [6]

p ro

of

Table 1: List of parameters Interpretation Recruitment Natural death rate the fraction of fast developing infectious cases disease-induced death cases Removal rate of infectious cases progression rate Initial value of susceptible individuals Initial value of exposed individuals Initial value of infectious individuals Initial value of removed individuals Initial value of the total population

discretization equations:

Pr e-

Symbol Λ µ q d r k S(0) E(0) I(0) R(0) N (0)

17

√ Sj+1 = Sj + [Λ − βSj Ij /(Sj + Ej + Ij + Rj ) − µSj ]∆t + σ1 Sj ∆tc1,j

σ12 Sj 2 (c1,j − 1)∆t, 2 √ = Ej + [(1 − q)βSj Ij /(Sj + Ej + Ij + Rj ) − (k + µ)Ej ]∆t + σ2 Ej ∆tc2,j

+ Ej+1

σ22 Ej 2 (c2,j − 1)∆t, 2 √ = Ij + [qβSj Ij /(Sj + Ej + Ij + Rj ) + kEj − (d + r + µ)Ij ]∆t + σ3 Ij ∆tc3,j

+

+

σ32 Ij 2 (c3,j − 1)∆t, 2

al

Ij+1

urn

√ σ 2 Rj Rj+1 = Rj + [rIj − µRj ]∆t + σ4 Rj ∆tc4,j + 4 (c24,j − 1)∆t, 2

Jo

where the time increment ∆t > 0, σi2 > 0, (i = 1, 2, 3, 4) denote the intensities of white noise, i,j (i = 1, 2, 3, 4) are the Gaussian random variables which follows the distribution N (0, 1). We choose β = 4.455(see[6]), σ12 = 0.0008, σ22 = 0.0008, σ32 = 0.0008, σ42 = 0.0008, and other parameter values see Table 1. In this case, by numerical simulations , we can see the exposed component E(t) and the infectious component I(t) of system (1.2) have positive values (see Fig.1). This means that the disease of system (1.2) remains in the population. In this situation, R0S = 1.1024 > 1, µ > 12 (σ12 ∨ σ22 ∨ σ32 ∨ σ42 ). According to Theorem 3.2, the system (1.2) has a unique stationary distribution. The numerical simulations are consistent with the result. If we increase the recruitment of susceptible σ22 to 0.0048, and other parameters consistent with Fig.1. In this case, by simulations , we can see that the exposed component E(t) and the infectious component I(t) of system (1.1) have positive values. The exposed component E(t) and the infectious component I(t) of system (1.2) tend to 0 as t approaches to +∞ . This implies that white noise makes the disease die out (see Fig.2).

Journal Pre-proof

Xin and Wang

Pr e-

p ro

of

18

Figure 1: k = 4.455, σ12 = 0.0008, σ22 = 0.0008, σ32 = 0.0008, σ42 = 0.0008, and other parameter values see Table 1. The blue lines represent the solution of system (1.2),and the yellow lines represent the solution of the corresponding undisturbed system (1.1). The right column displays the histogram of the probability density functions of S(t), E(t), I(t) and R(t).

Jo

urn

al

Now, if we decrease infection rate β to 3.7125, and other parameters consistent with Fig.1. In this case, by simulations , we can get that exposed component E(t) and the infectious component I(t) of system (1.2) tend to 0 as t approaches to +∞. This implies that the disease of system (1.2) dies out (see Fig 3). Compare Fig.3 with Fig.1, we can concluded that the infectious component of system (1.2) will die out by decreasing the number of infection rate. From Fig.3, we can find that the infectious component of system (1.1) have positive value , but the infectious component of systems (1.2) tend to 0 as t approaches to +∞. We conclude that the white noise is beneficial to the extinction of the disease. In this paper, we deal with a stochastic tuberculosis model . In fact, conditions for the disease die out are rough to some extent due to mathematical methods. Besides, there is a region for parameters between extinct conditions and a stationary distribution to the positive solutions, such as Fig.2 and Fig.3. Although, by numerical simulations, we can see the value of the infectious component as t approaches to +∞ , we want to explain it in mathematical theory.

Journal Pre-proof

19

Pr e-

p ro

of

A stochastic tuberculosis model

Jo

urn

al

Figure 2: k = 4.455, σ12 = 0.0008, σ22 = 0.0048, σ32 = 0.0008, σ42 = 0.0008, and other parameter values see Table 1.The blue lines represent the solution of system (1.2),and the yellow lines represent the solution of the corresponding undisturbed system (1.1). The right column displays the histogram of the probability density functions of S(t), E(t), I(t) and R(t).

Figure 3: k = 3.725, σ12 = 0.0008, σ22 = 0.0008, σ32 = 0.0008, σ42 = 0.0008, The blue lines represent the solution of system(1.2), and the yellow lines represent the solution of the corresponding undisturbed system (1.1). The right column displays the histogram of the probability density functions of S(t), E(t), I(t) and R(t).

Journal Pre-proof

20

Xin and Wang

References

of

[1] World Health Organization. Global Tuberculosis Report 2018. World Health Organization, Geneva (2018).

p ro

[2] Z. Feng, M. Iannelli, F. A. Milner, A two-strain tuberculosis model with age of infection, SIAM J. Appl. Math., 62, 1634-1635 (2002). [3] Y. Zhou, K. Khan, Z. Feng, J. Wu, Projection of tuberculosis incidence with increasing immigration trends, J. Theoret. Biol., 254, 215-228 (2008). [4] L. Liu, J. Wu, X.-Q. Zhao, The impact of migrant workers on the tuberculosis transmission: General models and a case study for China, Math. Biosci. Eng., 9, 785-807 (2012).

Pr e-

[5] D. P. Moualeu, M. Weiser, R. Ehrig, P. Deuflhard, Optimal control for a tuberculosis model with undetected cases in Cameroon. Commun. Nonlinear Sci. Numer. Simul., 20, 986-1003 (2015). [6] L. Liu, X.-Q. Zhao, Y. Zhou, A tuberculosis model with seasonality, Bull. Math. Biol., 72, 931-952 (2010). [7] X. Fan, Z. Wang, Stability analysis of an SEIR epidemic model with stochastic perturbation and numerical simulation, Int. J. Nonlinear Sci. Numer. Simul., 14, 113–121 (2013).

al

[8] Q. Liu, D. Jiang, N. Shi, T. Hayat, B. Ahmad, Stationary distribution and extinction of stochastic SEIR epidemic model with standard incidence, Physica A, 476, 58–69 (2017).

urn

[9] Q. Liu, D. Jiang, N. Shi, T. Hayat, A. Alsaedi, Asymptotic behavior of a stochastic delayed SEIR epidemic model with nonlinear incidence, Physica A, 462, 870–882 (2016) . [10] A. Julia, L.-H. Mariajesus, Cumulative and maximum epidemic sizes for a nonlinear SEIR stochastic model with limited resources, Discrete Contin. Dyn. Syst. Ser. B, 23, 3137–3151 (2018).

Jo

[11] Q. Yang, X. Mao, Extinction and recurrence of multi-group SEIR epidemic models with stochastic perturbations, Nonlinear Anal. Real World Appl., 14, 1434–1456 (2013). [12] D.J. Higham, An algorithmic introduction to numerical simulation of stochastic differential equation, SIAM Rev., 43, 525–546 (2001). [13] Q. Liu, D. Jiang, N. Shi, T. Hayat, B. Ahmad, Stationary distribution and extinction of astochastic SEIR epidemic model with standard incidence, Physica A, 476, 58-69 (2017).

Journal Pre-proof

A stochastic tuberculosis model

21

[14] Q. Liu, Q. Chen, D. Jiang, The threshold of a stochastic delayed SIR epidemic model with temporary immunity, Physica A, 450, 115–125 (2016).

of

[15] X. Mao, Stochastic Differential Equations and Their Application, Horwood, Chichester, (1997).

p ro

[16] Y.A. Kutoyants, Statistical Inference for Ergodic Diffusion Processes, Springer, London, (2003). [17] I. Nobuyuki, S. Watanabe, A comparison theorem for solutions of stochastic differential equations and its applications, Osaka J. Math. 14, 619–633 (1977). [18] R.Z. Khasminskii, Stochastic Stability of Differential Equations, Sijthoff & Noordhoff. (1980).

Pr e-

[19] S.M. Blower, The intrinsic transmission dynamics of tuberculosis epidemics, Nat. Med., 1, 815–821 (1995). [20] S.M. Blower, T. Chou, Modeling the emergence of the hot zones: tuberculosis and the amplification dynamics of drug resistance, Nat. Med., 10, 1111-1116 (1996). [21] E. Ziv, C.L. Daley, S.M. Blower, Eearly therapy for latent tuberculosis infection, Am. J. Epidemiol., 153, 381-385 (2001). [22] T.C. Gard, Persistence in stochastic food web models,Bull. Math. Biol., 46, 357-370 (1984).

Jo

urn

al

[23] T.C. Gard, Stability for multispecies population models in random environments, Nonlinear Anal., 10, 1411-1419 (1986).

Journal Pre-proof Dear Editors: On behalf of my co-authors, we declare that we have no conflict of interest.

of

Sincerely yours,

Jo

urn

al

Pr e-

p ro

Bin-Guo Wang.

Journal Pre-proof

Dear Reviewers: The following is the highlights of this paper: 1. Since Tuberculosis (TB) is a widespread infectious disease, we consider a stochastic tuberculosis model. 2.A stationary distribution and the extinction of the disease are showed.

of

3. We consider different number of the intensities of the white noise. The result of numerical simulations shows that white noise is beneficial to the

p ro

extinction of the disease.

Thank you very much for your attention and time! Sincerely yours,

Pr e-

Bin-Guo Wang School of Mathematics and Statistics Lanzhou University

Lanzhou, 730000, People’s Republic of China TEL:9318912481 FAX:9318912481

Jo

urn

al

E-mail: [email protected]

1