Asymptotic stability of a two-group stochastic SEIR model with infinite delays

Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

Contents lists available at ScienceDirect

Commun Nonlinear Sci Numer Simulat journal homepage: www.elsevier.com/locate/cnsns

Short communication

Asymptotic stability of a two-group stochastic SEIR model with infinite delays Meng Liu a, Chuanzhi Bai a,⇑, Ke Wang b a b

School of Mathematical Science, Huaiyin Normal University, Huaian 223300, PR China Department of Mathematics, Harbin Institute of Technology, Harbin 15000, PR China

a r t i c l e

i n f o

Article history: Received 18 June 2013 Received in revised form 16 February 2014 Accepted 17 February 2014 Available online 28 February 2014

a b s t r a c t A two-group stochastic SEIR epidemic model with infinite delays is proposed and investigated. Sufficient conditions for asymptotic stability are established. Some simulation figures are introduced to support the results. Ó 2014 Elsevier B.V. All rights reserved.

Keywords: Multi-group epidemic system Random perturbations Itô formula

1. Introduction In order to describe the spread of many infectious diseases in heterogeneous populations, such as HIV/AIDS, measles, gonorrhea and mumps, multi-group epidemic systems have received great attention recently. According to modes of contact patterns, transmission, or different geographic distributions, a host population can be divided into several groups. Therefore inter-group and within-group interactions can be modeled separately. Lajmanovich and Yorke did pioneering work in [1]. They proposed a class of SIS multi-group models for the transmission dynamics of gonorrhea, and showed the global stability of a unique endemic equilibrium by constructing quadratic Lyapunov functions. From then on, various forms of multi-group epidemic models have been developed and investigated, see e.g. [2–8]. In the analysis of multi-group epidemic models, one of main mathematical challenges is the global stability of the endemic equilibrium [2–5]. A complete resolution of this problem has been given by Li and his co-workers until recently. In [6–9], for a class of multi-group SEIR systems, the authors established the global stability of a unique endemic equilibrium by using the graph-theoretic approach to the construction of Lyapunov functions. Especially, Li et al. [9] proposed the following general multi-group epidemic model with infinite delays 8 n X R þ1 > S > > dSdtk ¼ Kk  bkj Sk 0 hj ðsÞEj ðt  sÞds  dk Sk ; > > > j¼1 > > > n < X R þ1 E dEk ¼ bkj Sk 0 hj ðsÞEj ðt  sÞds  ðdk þ ek ÞEk ; ð1Þ dt > > j¼1 > > > I dIk > > ¼ ek Ek  ðdk þ ck ÞIk ; > dt > : dR R k ¼ ck Ik  dk Rk dt ⇑ Corresponding author. Tel.: +86 051784183732. E-mail address: [email protected] (C. Bai). http://dx.doi.org/10.1016/j.cnsns.2014.02.025 1007-5704/Ó 2014 Elsevier B.V. All rights reserved.

3445

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

where k ¼ 1; . . . ; n. Sk ; Ek ; Ik and Rk stand for the susceptible, infected but non-infectious, infectious, and recovered populations in the kth group, respectively; Kk denotes the influx of individuals into the kth group; bkj represents the transmission R þ1 S E I R coefficient between compartments Sk and Ej ; hj ðsÞ is a nonnegative function satisfying 0 hj ðsÞds ¼ 1; dk ; dk ; dk and dk represent the death rates of S, E, I and R populations in the kth group, respectively; ek is the rate of becoming infectious in the kth group; ck denotes the recovery rate of infectious individuals in the kth group. All parameter values are assumed S

E

I

R

to be nonnegative and Kk ; dk ; dk ; dk ; dk > 0 for all k ¼ 1; . . . ; n. Since the variables Ik and Rk do not appear in the first two equations of (1), Li, Shuai and Wang [9] considered the following reduced model for Sk and Ek

8 n X R þ1 S > dSk > ¼ Kk  bkj Sk 0 hj ðsÞEj ðt  sÞds  dk Sk ; > dt > < j¼1

ð2Þ

n X > R þ1 > E dEk > > bkj Sk 0 hj ðsÞEj ðt  sÞds  ðdk þ ek ÞEk ; : dt ¼ j¼1

It has been shown [9] that system (2) has two possible equilibria: the infection-free equilibrium P0 ¼ ðS01 ; 0; . . . ; S0n ; 0Þ, where 

S

E

S0k ¼ KdSk , and the endemic equilibrium P  ¼ ðS1 ; E1 ; . . . ; Sn ; En Þ. Denote dk ¼ minfdk ; dk þ ek g, and define k

( 0 6 Sk 6

H¼ (

U¼

) Kk Kk ; 0 6 S þ E ð0Þ 6 ; E ðsÞ P 0; s 2 ð1; 0; k ¼ 1; . . . ; n : k k k  S dk dk

) Kk 0 < Sk < S ; 0 < Sk þ Ek ð0Þ <  ; Ek ðsÞ > 0; s 2 ð1; 0; k ¼ 1; . . . ; n : dk dk

Kk

Li et al. [9] have proved that Lemma 1. Assume that the transmission matrix ðbij Þnn is irreducible. Define M 0 ¼ the spectral radius.



bkj S0k

dEk þ k

e



and R0 ¼ qðM 0 Þ, where q denotes

nn

(i) If R0  1, then P 0 is the only equilibrium for system (2) in H. Moreover, P0 is globally asymptotically stable in H. (ii) If R0 > 1, then system (2) has a unique endemic equilibrium P  in U. Moreover,P is globally asymptotically stable in U. On the other hand, epidemic systems are often subject to environmental noise. Therefore lots of stochastic epidemic models have been developed and studied, see [10–20]. Beretta et al. [10] investigated the stability of endemic equilibrium of a stochastic delay SIR epidemic model. Carletti [11] considered the stability of a stochastic model for phage–bacteria interaction. Tornatore et al. [12] studied the stability of disease free equilibrium of a stochastic SIR model. Dalal et al. [13,14] proved that their stochastic epidemic systems have nonnegative solutions and studied the asymptotic stability of the solutions. Yu et al. [15] studied the stability of endemic equilibrium of a stochastic two-group SIR model. Ji et al. [20] considered a multigroup stochastic SIR epidemic model. The authors investigated the stability of disease free equilibrium and the persistence of endemic equilibrium. However, to the best of our knowledge, no results related to multi-group stochastic epidemic model with time delay have been reported. Motivated by these, in this paper, based on system (1), we propose a two-group stochastic SEIR model with infinite delay in Section 2. In Section 3, we establish the sufficient conditions for asymptotic stability of endemic equilibrium. Some simulation figures are introduced to illustrate our main result in Section 4. Section 5 gives the conclusions and future directions. 2. Stochastic model derivation In this paper, we only consider the case of n ¼ 2 in system (1):

8 2 X > R þ1 S > dSk > ¼ Kk  bkj Sk 0 hj ðsÞEj ðt  sÞds  dk Sk ; > dt > > > j¼1 > > > > 2 < X R þ1 E dEk ¼ bkj Sk 0 hj ðsÞEj ðt  sÞds  ðdk þ ek ÞEk ; dt > > j¼1 > > > > I > > dIdtk ¼ ek Ek  ðdk þ ck ÞIk ; > > > : dRk R ¼ ck Ik  dk Rk dt

ð3Þ

3446

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

where k ¼ 1; 2. Here we suppose that ðbij Þ22 is irreducible, i.e. b12 > 0 and b21 > 0, where i; j ¼ 1; 2. It follows from Lemma 1 0 1 0 0 b11 S1

B dE þ e that if R0 ¼ q@ b1 S01 21 2 dE2 þe2

P ¼

b12 S1

dE1 þe1 b22 S02

C A > 1, then system (3) has a unique endemic equilibrium

dE2 þe2

ðS1 ; E1 ; I1 ; R1 ; S2 ; E2 ; I2 ; R2 Þ:

We assume stochastic perturbations are of white noise type, which are directly proportional to distances SðtÞ; EðtÞ; IðtÞ and RðtÞ from values of S ; E ; I and R , influence on the becomes

dSk dt

;

dEk dt

;

dIk dt

k and dR , respectively (see e.g. [10,15,18]). So system (3) dt

" # 8 2 X R þ1 > S  > > dS ¼ K  b S h ðsÞE ðt  sÞds  d S j j k > k kj k 0 k k dt þ r1k ðSk  Sk ÞdB1k ; > > > j¼1 > > " # > > 2 X > R þ1 < E dEk ¼ bkj Sk 0 hj ðsÞEj ðt  sÞds  ðdk þ ek ÞEk dt þ r2k ðEk  Ek ÞdB2k ; j¼1 > > > h i > > I > > dIk ¼ ek Ek  ðdk þ ck ÞIk dt þ r3k ðIk  Ik ÞdB3k ; > > > h i > > R : dRk ¼ ck Ik  dk Rk dt þ r4k ðRk  Rk ÞdB4k :

ð4Þ

with initial conditions

Sk ð0Þ > 0;

Ek ðhÞ ¼ uk 2 BCðð1; 0; ð0; þ1ÞÞ;

Ik ð0Þ P 0; Rk ð0Þ P 0;

ð5Þ

where Bik ; i ¼ 1; 2; 3; 4; k ¼ 1; 2, are independent standard Brownian motions defined on a complete probability space ðX; F ; PÞ and r2ik represent the intensities of Bik ; BCðð1; 0; ð0; þ1ÞÞ is the space of bounded and continuous functions from ð1; 0 to ð0; þ1Þ with the norm jjujj ¼ suph0 juðhÞj. For more biological motivation on this type of modelling in epidemic model we refer the reader to [10,15,18]. By the work of Wei and Wang [21], for each initial condition (5), system (4) has a global solution. It is easy to see that stochastic system (4) has the same equilibrium points as system (3). In the next section, we shall study the stability of the equilibrium P of (4) by constructing some Lyapunov functions. 3. Stability of P  Consider the n-dimensional stochastic functional differential equation

dXðtÞ ¼ FðX t ; tÞdt þ GðX t ; tÞdBðtÞ;

X 0 ¼ u 2 BCðð1; 0; Rn Þ

n

ð6Þ n

where X t ¼ Xðt þ sÞ; s  0; R Þ is the space of bounded and continuous functions from ð1; 0 to R with the norm jjujj ¼ suph0 juðhÞj. Suppose that Eq. (6) has a zero solution. Let C 2;1 ðRn  ½0; þ1Þ; ð0; þ1ÞÞ be the family of all nonnegative functions Vðx; tÞ defined on Rn  ½0; þ1Þ such that they are continuously twice differentiable in x and once in t. For Vðx; tÞ 2 C 2;1 ðRn  ½0; þ1Þ; ð0; þ1ÞÞ function, define the following operator LV of Eq. (6)

h i LVðx; tÞ ¼ V t ðx; tÞ þ V x ðx; tÞFðxt ; tÞ þ 0:5 trace GT ðxt ÞV xx ðx; tÞGðxt Þ :

Definition 1 (See e.g. Mao [23]). (i) The zero solution of Eq. (6) is said to be stochastically stable or stable in probability if for every pair of r > 0, there exists a d > 0 such that if jjujj < d then

e 2 ð0; 1Þ and

P fjXðt; uÞj < r for all t P 0g P 1  e: (ii) The zero solution of Eq. (6) is said to be stochastically asymptotically stable if it is stochastically stable and, moreover, for every e 2 ð0; 1Þ, there exists a d > 0 such that if jjujj < d then

 P

 lim Xðt; uÞ ¼ 0 P 1  e:

t!þ1

By the change of variables

xk ¼ Sk  Sk ;

yk ¼ Ek  Ek ;

zk ¼ Ik  Ik ;

wk ¼ Rk  Rk ;

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

3447

we obtain the following system h i 8 R þ1 R þ1 R þ1 R þ1 S > dxk ¼  ðdk þ bk1 E1 þ bk2 E2 Þxk þ bk1 xk 0 h1 ðsÞy1 ðt  sÞds þ bk2 xk 0 h2 ðsÞy2 ðt  sÞds þ bk1 Sk 0 h1 ðsÞy1 ðt  sÞds þ bk2 Sk 0 h2 ðsÞy2 ðt  sÞds dt þ r1k xk dB1k ; > > > > h i > R þ1 R þ1 R þ1 R þ1 > > < dyk ¼ ðbk1 E1 þ bk2 E2 Þxk þ bk1 xk 0 h1 ðsÞy1 ðt  sÞds þ bk2 xk 0 h2 ðsÞy2 ðt  sÞds þ bk1 Sk 0 h1 ðsÞy1 ðt  sÞds þ bk2 Sk 0 h2 ðsÞy2 ðt  sÞds  ðdEk þ ek Þyk dt þ r2k yk dB2k ; h i I > > > > dzk ¼ ek yk  ðdk þ ck Þzk dt þ r3k zk dB3k ; > > h i > > : dwk ¼ c zk  dR wk dt þ r4k wk dB4k : k k

ð7Þ 

It is easy to see that the stability of P of (4) is equivalent to the stability of zero solution of system (7). Now we are in the position to give our main result. Theorem 2. Assume that B ¼ ðbij Þ22 is irreducible and R0 > 1. If E

r21k < 2dSk ; r22k <

2ðbk1 E1 þ bk2 E2 Þðdk þ ek Þ S dk

þ

E dk

þ ek þ

bk1 E1

þ bk2 E2

r23k < 2dIk þ 2ck ; r24k < 2dRk ;

;

ð8Þ

then the endemic equilibrium P  is stochastically asymptotically stable. Proof. From the stability theory of stochastic functional differential equations, we only need to find a Lyapunov function VðnÞ such that LVðnÞ 6 0 for sufficiently small d > 0 and the identity holds if and only if n ¼ 0 (see e.g. [22,23]). P Define V 1 ¼ 0:5 2k¼1 ak y2k , where ak ; k ¼ 1; 2, are positive constants to be chosen later. Making use of Itô’s formula leads to  Z þ1 Z þ1 2 X LV 1 ¼ ak yk ðbk1 E1 þ bk2 E2 Þxk þ bk1 xk h1 ðsÞy1 ðt  sÞds þ bk2 xk h2 ðsÞy2 ðt  sÞds k¼1

þ bk1 Sk

0

Z

þ1

0

h1 ðsÞy1 ðt  sÞds þ bk2 Sk

0

Z

þ1

0

 Z 2 2 h i X X E ¼ ak ðdk þ ek Þy2k þ ak bk1 Sk k¼1

þ1

0

k¼1

 2 X E h2 ðsÞy2 ðt  sÞds  ðdk þ ek Þyk þ 0:5 ak r22k y2k h1 ðsÞy1 ðt  sÞds þ bk2 Sk

 Z 2 2 X  X ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 þ k¼1

k¼1

þ1

  h2 ðsÞy2 ðt  sÞds yk

0

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

Z

k¼1

0

0

k¼1

0

0

þ1

 h2 ðsÞy2 ðt  sÞds xk yk

 2 # X    Z þ1 Z þ1 2 y y ðt  sÞ y ðt  sÞ y E þ ¼ ak Ek ðdk þ ek ÞEk k ak Ek bk1 Sk E1 h1 ðsÞ 1  ds þ bk2 Sk E2 h2 ðsÞ 2  ds k Ek Ek E1 E2 0 0 k¼1 k¼1  Z þ1  Z þ1 2 2 X  X ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 h1 ðsÞy1 ðt  sÞds þ bk2 h2 ðsÞy2 ðt  sÞds xk yk þ 2 X

"

Z

k¼1

ð9Þ

Set

k1 ¼ b S E ; b k1 k 1

k2 ¼ b S E ; b k2 k 2

k ¼ 1; 2:

Thus

k1 þ b k2 ¼ ðdE þ ek ÞE ; b k k

k ¼ 1; 2:

Then it follows from (9) that LV 1 ¼

"  2 # X  Z þ1   Z þ1 2 y ðt  sÞ y ðt  sÞ y k1 þ b k2 Þ yk k1 k2 ak Ek ðb ak Ek b h1 ðsÞ 1  ds þ b h2 ðsÞ 2  ds k þ E Ek E E 0 0 k 1 2 k¼1 k¼1  Z þ1  Z þ1 2 2 X  X ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 h1 ðsÞy1 ðt  sÞds þ bk2 h2 ðsÞy2 ðt  sÞds xk yk þ

2 X

k¼1

k¼1

0

0

k¼1

0

0

k¼1

0

0

" "Z "Z  2 # 2  2 # 2  2 # 2 2 2 þ1 þ1 X X X y ðt  sÞ y y ðt  sÞ y k1 þ b k2 Þ yk k1 k2 þ 0:5 ak Ek b þ 0:5 ak Ek b 6 ak Ek ðb h1 ðsÞ 1  ds þ k h2 ðsÞ 2  ds þ k E E E E E 0 0 k 1 k 2 k k¼1 k¼1 k¼1  Z þ1  Z þ1 2 2 X X  þ ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 h1 ðsÞy1 ðt  sÞds þ bk2 h2 ðsÞy2 ðt  sÞds xk yk k¼1

" Z  2 #  2  2 # Z þ1 2 þ1 X y ðt  sÞ y ðt  sÞ k1 þ b k2 Þ yk k1 k2 þ 0:5 ak Ek b 6 0:5 ak Ek ðb h1 ðsÞ 1  ds þ b h2 ðsÞ 2  ds Ek E1 E2 0 0 k¼1 k¼1  Z þ1  Z þ1 2 2 X  X ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 h1 ðsÞy1 ðt  sÞds þ bk2 h2 ðsÞy2 ðt  sÞds xk yk þ 2 X

k¼1

"

3448

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

where the last inequality follows from the famous Hölder’s inequality. Let

21 ; a1 E1 ¼ b

12 ; a2 E2 ¼ b

That is to say, choose

a1 ¼ b21 S2 ;

a2 ¼ b12 S1

Then we have

" "  2 #  2 # 21 ðb 11 þ b 12 Þ y1 12 ðb 21 þ b 22 Þ y2 LV 1  0:5b þ 0:5b E1 E2 "  2  2 # Z þ1 Z þ1 y1 ðt  sÞ y2 ðt  sÞ    þ 0:5b21 b11 h1 ðsÞ ds þ b12 h2 ðsÞ ds E1 E2 0 0 "  2  2 # Z þ1 Z þ1 y1 ðt  sÞ y2 ðt  sÞ    þ 0:5b12 b21 h1 ðsÞ ds þ b22 h2 ðsÞ ds E1 E2 0 0 þ

2 X  ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk k¼1

 Z 2 X ak bk1 þ

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

0

k¼1

Z

þ1

 h2 ðsÞy2 ðt  sÞds xk yk

ð10Þ

0

Define

" Z 21 b 11 V 2 ¼ 0:5b "

þ1

0

12 b 21 þ 0:5b

#  2 2 Z þ1 Z t  y 1 ð sÞ y2 ðsÞ  b d s ds þ h ðsÞ d s ds 12 2 E1 E2 ts 0 ts #  2 2 Z t Z þ1 Z t  y 1 ð sÞ y2 ðsÞ  h1 ðsÞ dsds þ b22 h2 ðsÞ dsds : E1 E2 ts 0 ts

h1 ðsÞ Z

þ1

0

Z

t

Then

"   "    2 #  2 # 2 2 y1 y2 y1 y       LV 2 ¼ 0:5b21 b11  þ b12  þ 0:5b12 b21  þ b22 2 E1 E2 E1 E2 "  2  2 # Z þ1 Z þ1 y1 ðt  sÞ y2 ðt  sÞ     0:5b21 b11 h1 ðsÞ ds þ b12 h2 ðsÞ ds E1 E2 0 0 "  2  2 # Z þ1 Z þ1 y ðt  sÞ y ðt  sÞ 12 b 21 22  0:5b h1 ðsÞ 1  ds þ b h2 ðsÞ 2  ds E1 E2 0 0 " "  2 #  2 # y1 y       þ 0:5b12 ðb21 þ b22 Þ 2 ¼ 0:5b21 ðb11 þ b12 Þ  E1 E1 "  2  2 # Z þ1 Z þ1 y1 ðt  sÞ y2 ðt  sÞ     0:5b21 b11 h1 ðsÞ ds þ b12 h2 ðsÞ ds E1 E2 0 0 "  2  2 # Z þ1 Z þ1 y ðt  sÞ y ðt  sÞ 12 b 21 22  0:5b h1 ðsÞ 1  ds þ b h2 ðsÞ 2  ds E1 E2 0 0

ð11Þ

Define V 3 ¼ V 1 þ V 2 , then it follows from (10) and (11) that

LV 3 ¼ LV 1 þ LV 2 6

 Z 2 2 X  X ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk þ ak bk1 k¼1

0

k¼1

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

Define 2 X V 4 ¼ 0:5 bk ðxk þ yk Þ2 ; k¼1

V 5 ¼ 0:5

2 X ck z2k ; k¼1

2 X V 6 ¼ 0:5 dk w2k k¼1

where bk ; ck and dk are positive constants to be chosen later. Using Itô’s formula gives

Z 0

þ1

 h2 ðsÞy2 ðt  sÞds xk yk :

ð12Þ

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

LV 4 ¼

3449

2 2 h i X X S E bk ðxk þ yk Þ dk xk  ðdk þ ek Þyk þ 0:5 bk ðr21k x2k þ r22k y2k Þ k¼1

k¼1

2 h i X S E S E ¼ bk ðdk  0:5r21k Þx2k  ðdk þ ek  0:5r22k Þy2k  ðdk þ dk þ ek Þxk yk

ð13Þ

k¼1

and

LV 5 ¼

2 2 2 h i h i X X X I I ck zk ek yk  ðdk þ ck Þzk þ 0:5 ck r23k z2k ¼ ck ðdk þ ck  0:5r23k Þz2k þ ek yk zk ; k¼1

k¼1

ð14Þ

k¼1

as well as

LV 6 ¼

2 2 2 h i h i X X X R R dk wk ck zk  dk wk þ 0:5 dk r24k w2k ¼ dk ðdk  0:5r24k Þw2k þ ck wk zk : k¼1

k¼1

ð15Þ

k¼1

Define V ¼ V 3 þ V 4 þ V 5 þ V 6 , then it follows from (12)–(15) that

LV ¼ LV 3 þ LV 4 þ LV 5 þ LV 6 6

2 X  ak 0:5r22k y2k þ ðbk1 E1 þ bk2 E2 Þxk yk k¼1

 Z 2 X ak bk1 þ

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

Z

0

k¼1

þ1

 h2 ðsÞy2 ðt  sÞds xk yk

0

2 h i X S E S E þ bk ðdk  0:5r21k Þx2k  ðdk þ ek  0:5r22k Þy2k  ðdk þ dk þ ek Þxk yk k¼1 2 2 h i X h i X I R ck ðdk þ ck  0:5r23k Þz2k þ ek yk zk þ dk ðdk  0:5r24k Þw2k þ ck wk zk þ k¼1

k¼1

2 n h i X S E ¼ bk ðdk  0:5r21k Þx2k  bk ðdk þ ek  0:5r22k Þ  0:5ak r22k y2k k¼1

h i o S E I R þ ak ðbk1 E1 þ bk2 E2 Þ  bk ðdk þ dk þ ek Þ xk yk þ ck ek yk zk  ck ðdk þ ck  0:5r23k Þz2k  dk ðdk  0:5r24k Þw2k þ dk ck wk zk  Z þ1  Z þ1 2 X ak bk1 h1 ðsÞy1 ðt  sÞds þ bk2 h2 ðsÞy2 ðt  sÞds xk yk þ 0

k¼1

0

 Z 2 X ¼: L0 V þ ak bk1 k¼1

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

Z

0

þ1

 h2 ðsÞy2 ðt  sÞds xk yk ;

0

where

L0 V ¼

2 n X

h i h i S E S E bk ðdk  0:5r21k Þx2k  bk ðdk þ ek  0:5r22k Þ  0:5ak r22k y2k þ ak ðbk1 E1 þ bk2 E2 Þ  bk ðdk þ dk þ ek Þ xk yk

k¼1

o I R þ ck ek yk zk  ck ðdk þ ck  0:5r23k Þz2k  dk ðdk  0:5r24k Þw2k þ dk ck wk zk : Choose

b1 ¼

b2 ¼

a1 ðb11 E1 þ b12 E2 Þ S d1

þ

E d1

þ e1

a2 ðb21 E1 þ b22 E2 Þ S d2

þ

E d2

þ e2

¼

¼

b21 S2 ðb11 E1 þ b12 E2 Þ S

b12 S1 ðb21 E1 þ b22 E2 Þ S

E

d2 þ d2 þ e2

Moreover, using Cauchy inequality to

1 4

ek yk zk 6 ðdIk þ ck  0:5r23k Þz2k þ 1 4

E

d1 þ d1 þ e1

ck wk zk 6 ðdRk  0:5r24k Þw2k þ

:

ek yk zk and ck wk zk , we can see that I dk

e2k y2k ; þ ck  0:5r23k

c2k z2k R dk

;

 0:5r24k

;

Substituting the above inequalities into (16) results in

ð16Þ

3450

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

L0 V 6

( " 2 X S E bk ðdk  0:5r21k Þx2k  bk ðdk þ ek  0:5r22k Þ  0:5ak r22k  k¼1

"

#

dk c2k

I

 0:75ck ðdk þ ck  0:5r23k Þ 

#

ck e2k I

dk þ ck  0:5r23k )

y2k

R

R

dk  0:5r24k

z2k  0:75dk ðdk  0:5r24k Þw2k :

2 X

2  Ak xk þ Bk y2k þ C k z2k þ Dk w2k ; ¼:  k¼1

where S

ck e2k

E

Ak ¼ bk ðdk  0:5r21k Þ;

Bk ¼ bk ðdk þ ek  0:5r22k Þ  0:5ak r22k  dk c2k

I

C k ¼ 0:75ck ðdk þ ck  0:5r23k Þ 

R dk

 0:5r24k

;

I dk

þ ck  0:5r23k

;

R

Dk ¼ 0:75dk ðdk  0:5r24k Þ:

If (8) holds, then Ak > 0 and E

E

bk ðdk þ ek  0:5r22k Þ  0:5ak r22k ¼ bk ðdk þ ek Þ  0:5ðak þ bk Þr22k > 0; I

R

dk þ ck > 0:5r23k ;

dk > 0:5r24k :

Now choose ck such that

dk þ ck  0:5r23k h I

0 < ck <

e

2 k

i E bk ðdk þ ek Þ  0:5ðak þ bk Þr22k :

Then Bk > 0. Choose dk such that R

0 < dk <

dk  0:5r24k

c2k

I

0:75ck ðdk þ ck  0:5r23k Þ:

Thus C k > 0 and Dk > 0. Consequently

 Z 2 2 X

2  X LV 6  Ak xk þ Bk y2k þ C k z2k þ Dk w2k þ ak bk1 k¼1

k¼1

þ1

h1 ðsÞy1 ðt  sÞds þ bk2

0

Z

þ1

 h2 ðsÞy2 ðt  sÞds xk yk ;

0

Let us suppose that Pfjy1 ðsÞj 6 d; jy2 ðsÞj 6 dg ¼ 1, then

 Z 2 X ak bk1

þ1 0

k¼1

h1 ðsÞy1 ðt  sÞds þ bk2

Z

þ1 0

 2 X h2 ðsÞy2 ðt  sÞds xk yk 6 d ak ðbk1 þ bk2 Þjxk yk j k¼1

2 X  0:5d ak ðbk1 þ bk2 Þðx2k þ y2k Þ: k¼1

Therefore 2 X

 LV 6  ½Ak  0:5dak ðbk1 þ bk2 Þx2k þ ½Bk  0:5dak ðbk1 þ bk2 Þy2k þ C k z2k þ Dk w2k : k¼1

Consequently, for sufficiently small d > 0 we have LV < 0 except the zero point. Then the required assertion follows immediately. Remark 1. By comparing Lemma 1 with Theorem 2 we can see that if the positive equilibrium of the deterministic model is stable, then the stochastic system will keep this nice property provided the noise is sufficiently small.

4. Numerical simulations In this section let use the Milstein method mentioned in Higham [24] to substantiate this result. For convenience, set hj ðsÞ ¼ es ; s P 0, and the initial condition is Ej ðhÞ ¼ tj eh ; h  0, where tj > 0. Then Eq. (4) becomes

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

" # 8 2 2 X X Rt s > S  et > t > dS ¼ K  t b S  e b S e E ðsÞds  d S j kj k j k k > kj k 0 k k dt þ r1k ðSk  Sk ÞdB1k ; 2 > > > j¼1 j¼1 > > " # > > 2 2 > X X Rt s < E et t dEk ¼ 2 tj bkj Sk þ e bkj Sk 0 e Ej ðsÞds  ðdk þ ek ÞEk dt þ r2k ðEk  Ek ÞdB2k ; j¼1 j¼1 > > > h i > > I > > dIk ¼ ek Ek  ðdk þ ck ÞIk dt þ r3k ðIk  Ik ÞdB3k ; > > > h i > > : dR ¼ c I  dR R dt þ r ðR  R ÞdB : k 4k k 4k k k k k k

3451

ð17Þ

Here, Sk ; Ek ; Ik ; Rk ; Sk ; Ek ; Ik ; Rk and v j have units of ‘‘population’’; Kk has units of population/time; bkj has units of 1/(popS E I R ulation  time); dk ; dk ; dk ; dk ; ek ; ck and rik has units of 1/time, k; j ¼ 1; 2; i ¼ 1; 2; 3; 4. Consider the discretization equation 8 " # 2 2 i X X X pffiffiffiffiffiffi ðiÞ r2 ðiÞ ðiÞ > > S ðiÞ ðiþ1Þ ðiÞ ðiÞ ðiÞ ðiÞ ðiÞ 2 eiDt iDt s Dt ðsÞ > bkj tj Sk  e bkj Sk e Ej ðsÞDt  dk Sk Dt þ r1k ðSk  Sk Þ Dt n1k þ 21k Sk ðSk  Sk Þ½ðn1k Þ  1Dt; > Sk  Sk ¼ Kk  2 > > > s¼1 j¼1 j¼1 > > > " # > 2 2 i > X X X pffiffiffiffiffiffi ðiÞ r2 ðiÞ ðiÞ > ðiþ1Þ i D t < E ðiÞ ðiÞ ðiÞ ðiÞ ðiÞ ðiÞ 2 e iDt s Dt ðsÞ bkj tj Sk þ e bkj Sk e Ej ðsÞDt  ðdk þ ek ÞEk Dt þ r2k ðEk  Ek Þ Dt n2k þ 22k Ek ðEk  Ek Þ½ðn2k Þ  1Dt; Ek  Ek ¼ 2 s¼1 j¼1 j¼1 > > > h i > pffiffiffiffiffiffi 2 > > Iðiþ1Þ  IðiÞ ¼ ek EðiÞ  ðdI þ c ÞIðiÞ Dt þ r3k ðIðiÞ  I Þ Dt nðiÞ þ r3k IðiÞ ðIðiÞ  I Þ½ðnðiÞ Þ2  1Dt; > k k > k k k k k k k 3k k 3k 2 k > > > h i pffiffiffiffiffiffi 2 > 2 > : Rðiþ1Þ  RðiÞ ¼ c IðiÞ  dR RðiÞ Dt þ r4k ðRðiÞ  R Þ Dt nðiÞ þ r4k RðiÞ ðRðiÞ  R Þ½ðnðiÞ Þ  1Dt: k k k k k k k k k 4k k k 4k 2

ðiÞ

where njk ; i ¼ 1; 2; . . . ; n; j ¼ 1; 2; 3; 4; k ¼ 1; 2, are the Gaussian random variables which follow Nð0; 1Þ. In this section, S E I R  for simplicity, we always assume that dk ¼ dk ¼ dk ¼ dk ¼ d and t1 ¼ 2:5; t2 ¼ 3:5; S1 ð0Þ ¼ 4; S2 ð0Þ ¼ 2:5, k I1 ð0Þ ¼ 2:2; I2 ð0Þ ¼ 2:2; R1 ð0Þ ¼ 2:5; R2 ð0Þ ¼ 3. 1 ¼ 0:4; e1 ¼ 0:9; c ¼ 0:8; K2 ¼ 5:75; b In the above equation, we choose K1 ¼ 4:6; b11 ¼ 0:135; b12 ¼ 0:1; d 21 1 2 ¼ 0:5; e2 ¼ 0:9; c ¼ 1:75. Then S ¼ 5; E ¼ 2; I ¼ 1:5; R ¼ 3 and S ¼ 4:5; E ¼ 2:5; I ¼ 1; ¼ 0:0889; b22 ¼ 0:24; d 2 1 1 1 2 2 1 2 R2 ¼ 3:5. In Fig. 1 and Fig. 2, the unit of ‘‘population’’ is thousands of people and the unit of ‘‘time’’ is days. For example, the unit of b11 is 1/(thousands of people  days). The only difference between conditions of Fig. 1 and Fig. 2 is that the values of rjk ; j ¼ 1; 2; 3; 4; k ¼ 1; 2, are different. In Fig. 1, we choose rjk ¼ 0; j ¼ 1; 2; 3; 4; k ¼ 1; 2. In view of Lemma 1, the equilibrium position P is asymptotically stable. Fig. 1 confirms this. In Fig. 2, we choose r11 ¼ 0:8; r21 ¼ 0:6; r31 ¼ 1:1; r41 ¼ 0:75; r12 ¼ 0:9; r22 ¼ 0:65; r32 ¼ 2; r42 ¼ 0:95. Thus (8) holds. Then it follows from Theorem 2 that the equilibrium position P is asymptotically stable. See Fig. 2 (in order to avoid confusion, we plot S1 ; E1 ; I1 ; R1 and S2 ; E2 ; I2 ; R2 in Fig. 2(a) and (b), respectively). By comparing Fig. 1 with Fig. 2, we can see that if the positive equilibrium of the deterministic model is asymptotically stable, then the stochastic system will keep this nice property provided the noise is sufficiently small.

5.5

S1

5

S2

4.5

E1 E2

4

I1 I2

3.5

R1

3

R2 2.5 2 1.5 1 0.5

0

5

10

15

20

25

Time

1 ¼ 0:4; e1 ¼ 0:9; c ¼ 0:8; K2 ¼ 5:75; b ¼ 0:0889; Fig. 1. Solutions of system (17) for t1 ¼ 2:5; t2 ¼ 3:5; K1 ¼ 4:6; b11 ¼ 0:135; b12 ¼ 0:1; d 21 1 2 ¼ 0:5; e2 ¼ 0:9; c ¼ 1:75; r ¼ 0; j ¼ 1; 2; 3; 4; k ¼ 1; 2; Dt ¼ 0:01; S1 ð0Þ ¼ 4; S2 ð0Þ ¼ 2:5; I1 ð0Þ ¼ 2:2; I2 ð0Þ ¼ 2:2; R1 ð0Þ ¼ 2:5; R2 ð0Þ ¼ 3. b22 ¼ 0:24; d jk 2

3452

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453 6

S1(t) E1(t)

5

I1(t) R1(t)

4

3

2

1

(a) 0

0

5

10

15

20

25

Time 7 S2(t) E2(t)

6

I2(t) R2(t)

5

4

3

2

1

0

(b) 0

5

10

15

20

25

Time

1 ¼ 0:4; e1 ¼ 0:9; c ¼ 0:8; K2 ¼ 5:75; b ¼ 0:0889; b ¼ Fig. 2. Solutions of system (17) for t1 ¼ 2:5; t2 ¼ 3:5; K1 ¼ 4:6; b11 ¼ 0:135; b12 ¼ 0:1; d 21 22 1  ¼ 0:5; e ¼ 0:9; c ¼ 1:75; r ¼ 0:8; r ¼ 0:6; r ¼ 1:1; r ¼ 0:75; r ¼ 0:9; r ¼ 0:65; r ¼ 2; r ¼ 0:95; Dt ¼ 0:01; S ð0Þ ¼ 4; S ð0Þ ¼ 2:5; 0:24; d 2 2 11 21 31 41 12 22 32 42 1 2 2 I1 ð0Þ ¼ 2:2; I2 ð0Þ ¼ 2:2; R1 ð0Þ ¼ 2:5; R2 ð0Þ ¼ 3. (a) is the figure of S1 ; E1 ; I1 ; R1 ; (b) is the figure of S2 ; E2 ; I2 ; R2 .

5. Concluding remarks A two-group stochastic SEIR epidemic model with infinite delays was proposed and investigated. Sufficient conditions for asymptotic stability were established. The results and simulation figures show that if the positive equilibrium of the deterministic system is stable, then the stochastic model will preserve this nice property provided the noise is sufficiently small. Owing to its theoretical and practical significance, stability of multi-group system with time delays has deserved a lot of attention, but mainly in deterministic case (see e.g. [5,8,9]). The present paper is the first attempt, up to our knowledge, to study the stochastic case. Some interesting questions deserve further investigation. In this paper, we consider the two group case. It is interesting to study the general n group case. In fact, we also attempted to investigate this problem. Unfortunately, there are some technical obstacles that can not be overcome at present stage. Moreover, in this paper, we suppose that the stochastic perturbations are proportional to SðtÞ  S ; EðtÞ  E ; IðtÞ  I and RðtÞ  R , it is interesting to study other types of stochastic perturbations (see e.g. [16–20,25–27]).

Acknowledgments The authors thank the editor and reviewers for their important and valuable comments. This research is supported by the Natural Science Foundation of PR China (Nos. 11301207, 11171081, and 11171056), Natural Science Foundation of Jiangsu Province (Nos. BK2011407 and BK20130411), Natural Science Research Project of Ordinary Universities in Jiangsu Province (No. 13KJB110002).

M. Liu et al. / Commun Nonlinear Sci Numer Simulat 19 (2014) 3444–3453

3453

References [1] Lajmanovich A, Yorke JA. A deterministic model for gonorrhea in a nonhomogeneous population. Math Biosci 1976;28:221–36. [2] Hethcote HW, Thieme HR. Stability of the endemic equilibrium in epidemic models with subpopulations. Math Biosci 1985;75:205–27. [3] Beretta E, Capasso V. Global stability results for a multigroup SIR epidemic model. In: Hallam TG, Gross LJ, Levin SA, editors. Mathematical ecology. Teaneck, NJ: World Scientific; 1988. p. 317–42. [4] Huang W, Cooke KL, Castillo-Chavez C. Stability and bifurcation for a multiple-group model for the dynamics of HIV/AIDS transmission. SIAM J Appl Math 1992;52:835–54. [5] Thieme HR. Mathematics in population biology. Princeton: Princeton University Press; 2003. [6] Guo H, Li MY, Shuai Z. Global stability of the endemic equilibrium of multigroup SIR epidemic models. Can Appl Math Q 2006;14:259–84. [7] Guo H, Li MY, Shuai Z. A graph-theoretic approach to the method of global Lyapunov functions. Proc Amer Math Soc 2008;136:2793–802. [8] Li MY, Shuai Z. Global-stability problem for coupled systems of differential equations on networks. J Differ Equ 2010;248:1–20. [9] Li MY, Shuai Z, Wang C. Global stability of multi-group epidemic models with distributed delays. J Math Anal Appl 2010;361:38–47. [10] Beretta E, Kolmanovskii V, Shaikhet L. Stability of epidemic model with time delays influenced by stochastic perturbations. Math Comput Simul 1998;45:269–77. [11] Carletti M. On the stability properties of a stochastic model for phage–bacteria interaction in open marine environment. Math Biosci 2002;175:117–31. [12] Tornatore E, Buccellato SM, Vetro P. Stability of a stochastic SIR system. Phys A 2005;354:111–26. [13] Dalal N, Greenhalgh D, Mao XR. A stochastic model of AIDS and condom use. J Math Anal Appl 2007;325:36–53. [14] Dalal N, Greenhalgh D, Mao XR. A stochastic model for internal HIV dynamics. J Math Anal Appl 2008;341:1084–101. [15] Yu J, Jiang D, Shi N. Global stability of two-group SIR model with random perturbation. J Math Anal Appl 2009;360:235–44. [16] Lu Q. Stability of SIRS system with random perturbations. Phys A 2009;388:3677–86. [17] Chen G, Li TC. Stability of stochastic delayed SIR model. Stochastics Dyn 2009;9:231–52. [18] Jiang D, Ji C, Shi N, Yu J. The long time behavior of DI SIR epidemic model with stochastic perturbation. J Math Anal Appl 2010;372:162–80. [19] Meng X. Stability of a novel stochastic epidemic model with double epidemic hypothesis. Appl Math Comput 2010;217:506–15. [20] Ji C, Jiang D, Shi N. Multigroup SIR epidemic model with stochastic perturbation. Phys A 2011;390:1747–62. [21] Wei F, Wang K. The existence and uniqueness of the solution for stochastic functional differential equations with infinite delay. J Math Anal Appl 2007;331:516–31. [22] Hasminskij RZ. Stochastic stability of differential equations. The Netherlands: Sijthoof & Noordhoof; 1980. [23] Mao XR. Stochastic differential equations and their applications. Chichester: Horwood publishing; 1997. [24] Higham D. An algorithmic introduction to numerical simulation of stochastic differential equations. SIAM Rev 2001;43:525–46. [25] Liu M, Wang K. Analysis of a stochastic autonomous mutualism model. J Math Anal Appl 2013;402:392–403. [26] Liu M, Wang K. Dynamics of a Leslie–Gower Holling-type II predator–prey system with Lévy jumps. Nonlinear Anal 2013;85:204–13. [27] Liu M, Wang K. Stochastic Lotka–Volterra systems with Lévy noise. J Math Anal Appl 2014;410:750–63.