Stability Analysis of a Multigroup Epidemic Model with General Exposed Distribution and Nonlinear Incidence Rates Academic research paper on "Mathematics"

Hindawi Publishing Corporation Abstract and Applied Analysis Volume 2013, Article ID 354287,11 pages http://dx.doi.org/10.1155/2013/354287

Research Article

Stability Analysis of a Multigroup Epidemic Model with General Exposed Distribution and Nonlinear Incidence Rates

Ling Zhang,1 Jingmei Pang,2 and Jinliang Wang2

1 School of Science, Department of Fundamental Mathematics, Jiamusi University, Jiamusi 154007, China

2 School of Mathematical Science, Heilongjiang University, Harbin 150080, China

Correspondence should be addressed to Jinliang Wang; jinliangwang@hit.edu.cn Received 26 January 2013; Revised 16 June 2013; Accepted 22 July 2013 Academic Editor: Pagavathi Balasubramaniam

Copyright © 2013 Ling Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

We investigate a class of multigroup epidemic models with general exposed distribution and nonlinear incidence rates. For a simpler case that assumes an identical natural death rate for all groups, and with a gamma distribution for exposed distribution is considered. Some sufficient conditions are obtained to ensure that the global dynamics are completely determined by the basic production number R0. The proofs of the main results exploit the method of constructing Lyapunov functionals and a graph-theoretical technique in estimating the derivatives of Lyapunov functionals.

1. Introduction

Multigroup epidemic models have been used in the literature to describe the transmission dynamics of many different infectious diseases such as mumps, measles, gonorrhea, HIV/AIDS and vector borne diseases such as Malaria [1]. In the models, heterogeneous host population can be divided into several homogeneous groups according to modes of transmission, contact patterns, or geographic distributions, so that within-group and intergroup interactions can be modeled separately. It is well known that global dynamics of multigroup models with higher dimensions, especially the global stability of the endemic equilibrium, are a very challenging problem. Guo et al. [2] proposed a graph-theoretic approach to the method of global Lyapunov functions and used it to resolve the open problem on the uniqueness and global stability of the endemic equilibrium of a multigroup SIR model with varying subpopulation sizes. Subsequently, a series of studies on the global stability of multigroup epidemic models were produced in the literature (see e.g., [2-5]).

In the present paper, a more general multigroup epidemic model is proposed and studied to describe the disease spread in a heterogeneous host population with general exposed distribution and nonlinear incidence rate. The host population is divided into m distinct groups (m > 1). For 1 < i < m,

the (th group is further partitioned into four disjoint classes: the susceptible individuals, exposed individuals, infectious individuals, and recovered individuals, whose numbers of individuals at time t are denotedby St(t), Et(t), It(t), and Rt(t), respectively. Susceptible individuals infected with the disease but not yet infective are in the exposed (latent) class.

It is pointed in [6] that a fixed latent period can be considered as an approximation of the mean latent period, and this would be appropriate for those diseases whose latent periods vary only relatively slightly. For example, poliomyelitis has a latent period of 1-3 days (comparing to its much longer infectious period of 14-20 days). However disease such as tuberculosis, including bovine tuberculosis (a disease spread from animal to animal mainly by direct contact), may take months to develop to the infectious stage and also can relapse. Since the time it takes from the moment of new infection to the moment of becoming infectious may differ from disease to disease, even for the same disease, it differs from individual to individual, and it is indeed a random variable. It is thus of interest from both mathematical and biological viewpoints to investigate whether sustained oscillations are the result of general exposed distribution.

Following the method of [6], we also assume that the disease does not cause deaths during the latent period, taking the natural death rate into consideration. Let P(t) denote

the probability that an exposed individual remains in the time t after entering the exposed class. For 1 < i, j < m, Pij > 0 denotes the coefficient of transmission between compartments Si and Ij. It is assumed that m-square matrix (Pij)\<i,j<m is irreducible [7]. So the proportion of exposed individuals can be expressed by the integral

where the sum takes into account cross-infections from all groups. Integrals in (1) are in the Riemann-Stieltjes sense. Pj(t) satisfies the following reasonable properties:

(A)Pi : [0,œ)

[0,1] is nonincreasing, piecewise

continuous with possibly finitely many jumps and satisfies Pj(0+) = 1, and limt^mPj(t) = 0 with

J™ Pj(t)dt is positive and finite. Differentiating (1) gives

E', (t) = Zfiij ftj (St (t),Ij (t))

the model to be studied takes the following differential and integral equations form:

S't (t) = cp, (S, (t))-fß,j f,j (S, (t),Ij (t)), j=i

<f) = pu Î fa (,s' M ' h M) e-S,(t-u)p> « - u) du> (i) E'(t) = I^f'j(S'(t) ' h(t))

Zß,j I f,j (S, (u),Ij (u))

j=l J0

X Pj (t -u)du- SiEi (t).

e-S1(t-u) (2)

e-S1(t-u)

Xßj I fj (S, (u),Ij M)

j=l J0

x hj (t -u)du- SiEi (t),

' (t) = Tßj J0 fj (s> M'lj wy8'™

x hj (t -u)du- (St + £, + y,) Ij (t), R't (t) = y,I, (t)-S,R, (t).

Since the variables Ei and Ri do not appear in the first and third equations of model (5), Ei(t) and Ri(t), i = 1,...,m, can be decoupled from the Si(t) and Ii(t) equations; we only need to consider the subsystem of (5) consisting of only the Si and Ii equations:

The first term on the right hand side in (2) is the rate at which new infected individuals come into the exposed class, and the last term explains the natural deaths. The second term accounts for the rate at which the individuals move to the infectious class (noting that P'j(t - u) < 0 due to the aformentioned property) from the exposed class; hence

t «)= - Iß,J I fj (Si M.Ij (u))e~S>™

j=i J0

x p'. (t - u) du - (S, + e, + y,) I, (t).

Let hj(t) = -Pj(t) be the probability density function for the time (a random variable) it takes for an infected individual in the (th group to become infectious. Then (4) becomes

Ï (t) = Ißij I f,j (S, (u),Ij (u))

i=i J0

e-S1(t-u)

x h: (t -u)du- (St + £, + yt ) Ij (t).

Within the ith group, fi(Si) denotes the growth rate of Si, which includes both the production and the natural death of susceptible individuals. Therefore, under the assumptions,

s' (t) = <Pi (S, (t))-Iß,j f,j (s, (t),Ij (t)),

t (t) = Ißij I fj (Si (u).Ij (»))

e-Sj(t-u) (6)

x hj (t -u)du- (Si + £i + yi) Ij (t),

where Si denotes the natural death rates of Ii compartments in the (th group, £i is the death rate caused by disease in the (th group, and yi is the rate of recovery of infectious individuals in the (th group. In what follows we investigate the global stability of system (5).

When m = 1, P(t) = eet, and with bilinear incidence rate, system (5) will reduce to the standard SEIR ordinary differential equation (ODE) model studied in [8, 9], and if we further assume that P(t) is a step function, system (5) becomes the SEIR model with a discrete delay studied in [10]. Recently, a model of this type, but including the possibility of disease relapse, has been proposed in [11,12] to investigate the transmission of herpes, and its global dynamics have been completely investigated in [5,13].

To express the main idea and the approaches more clearly, we consider a simpler case in which all groups share the same natural death rate: Sj = 8 for j = 1,2, ...,m. Further, we assume that the functions hj(u) are disease specific only,

implying that hj(u) = h(u) for j = 1,2,... ,m. We choose the gamma distribution:

k (u) = Kb (u) =

(n- \)!bn

where b > 0 is a real number and n > 1 is an integer, which is widely used and can approximate several frequently used distributions. For example, when b ^ 0+, hnb(s) will approach the Dirac delta function, and when n = 1, hnb(s) is an exponentially decaying function.

The main object of this paper is to carry out the well-known "linear chain trick" to system (6), transfer system into an equivalent ordinary differential equations system, and establish its global dynamics. We derive the basic reproductive number R0 and show that R0 completely determines the global dynamics of system (6). More specifically, if R0 < 1,the disease-free equilibrium is globally asymptotically stable and the disease dies out; if R0 > 1, a unique endemic equilibrium exists and is globally asymptotically stable, and the disease persists at the endemic equilibrium. The global stability of P* rules out any possibility for Hopf bifurcations and the existence of sustained oscillations. We should point out here that this work is motivated by Yuan and Zou [11,12,14]. In the proof we demonstrate that the graph-theoretic approach developed in [2, 3] can be successfully applied to construct suitable Lyapunov functionals and thus prove the global stability of the endemic equilibrium for model (6) with general exposed distribution and nonlinear incidence rate. Our work is also based on a recent work by Sun and Shi [15], which resolved the dynamics of multigroup SEIR epidemic models with nonlinear incidence of infection and nonlinear removal functions between compartments.

In Section 2, we first give the model, preliminaries and the basic reproduction number R0. The global stability of the corresponding equilibria for R0 <1 and R0 > 1 is shown, respectively, in Section 3—the key results of this paper. And in Section 4, some numerical simulations are shown to illustrate the effectiveness of the proposed result.

2. Preliminaries

We make the following basic assumptions for the intrinsic growth rate of susceptible individuals in the ith group pi(Si) and the transmission functions f(St, Ij).

(A^ p{ are C1 non-increasing functions on [0, rn) with pt(0) > 0, and there is a unique positive solution % = S° for the equation p{(%) = 0. Pj(S) > 0 for 0<S<S0, and (pt (S) <0 for S > S°; that is

h (s,)-<p, (s?)](s, -s0)<o

for S, = S0, i = 1,2,..., m.

(A 2)ftj (S,,Ij)<clj (S,)Ij for all I j >0.

(A 3) cj (S,) < Cj (S0), 0<S, < S0, i,j=1,..., m.

Following the technique and method in [14], define

which can absorb the exponential term e Su into the delay kernel. The second equation in (6) can be rewritten as

i; «)=!

=i(1 + 8b) Jo - iß + h + Yt) 1,.

ffij (S,,ij)K:b (t-u)du

For I = 1,..., n,let

m ß- b c*

^ (f) = güW I fj (SIj) kl'b (f - u) du

i = 1,2, ...,m.

Thus, for I e {2,..., n}, we obtain

y,; = ki (0)I

-f (s„ Ij)

=1(1+8b) ßijb (* (l-1)(t-u)1-2

j=l(1+8b)n J-ra (l-1)®

e-(t-u)/b fij (s,, Ij) du

( - u)

j=1(1 + 8b)n J-m (l-1)\bl+1

= [yg-i - y,,i ] b .

For I = 1,we have

e-^fj (S,, Ij) du

= y ßj b t e-{t-u)ih yi'1 = y (1 + 8b)n \

\ (Si, Ij) du,

=1(1 +8b) b (13)

i = 1,...,m.

It follows that

*= Z (1 + lb)nf (s',Ij")

ßj t e-(t-u)/h J=1 (1 + 8b)n l-m $

fj (S,, Ij) du (14)

m ß.. 1 = y777ß^ fj (si,Ij)-Tyi,i, i=1,...,'

1 (1 + 8b)n

Thus the integro-differential system (6) is equivalent to the ordinary differential equations

s', (t) = 9i (Si (t))-lPtj fij (St (t),Ij (t)), j=i

1 m 1 & (f) = Ifrj ftj (Si (t), Ij (tj) - ^ (t),

(1+8b)npi y',,2 (t) = l(yt,i (t)-yt,2 (t)), i = 1,2,... ,m,

y'hn (V = 1 (yi,n-i (t)-yt,n (t)), b

I', (t) = 1y,,n (t)-(S + ei + Yl)h (t). b

For initial condition

(Si (0),yhi (0),...,yhn (0),li (0),

S2 (0),y2ii (0),...,yhn (0),I2 (0),..., (16)

Sm (0),ym,i (0),...,ym,n (0)Jm (0))e Rm(n+2\

the existence, uniqueness, and continuity of the solution (S,, yi,i ,yi,2,..., yi,n, h) of system (15) follow from the standard theory of Volterra integro-differential equation [16]. It can also be verified that every solution of (15) with nonnegative initial condition remains nonnegative.

It follows from (A i) and the first equation of (15) that limsupt^TO Sj(t) < S° for all i = 1,2,..., m. Let be the maximum of the function fi on R+ and let q be a positive real number such that q > bNf . Denote by Yi the (th tube for system (15); that is,

Yi = (Si,yii,yia,...,yi,n,Ii).

It follows from a similar argument to that in [14] that we can show that the set De defined by

-*m(n+2) +

r = {{Si,yi,i,yi,2,...,yi,n,Ii) 6 Rm

St <S° + e, Si + (1 + Sb)nyiA <q + S°, q + S® + le

y,,i <

(1+Sb)n ' q + S° + (n + 1) e

b(1 + Sb)n (S + £i + Vi) i=1,2,...,m, I = 2,3,... ,n}

is a forward invariant compact absorbing set for system (15) for e > 0 and that the set ro (i.e., when e = 0) is a forward invariant compact set.

Under the assumption (Ai), we know that system (15) always has the disease-free equilibrium

Po = (s0,0,...,0,lO,sO ,0,...,0,I°o,...,S°m,0,...,0,I0m )

An equilibrium P* of (6) has the form P* = (S*, I** ,S*2,r2,..., S*m, I*) e R2m with S** > 0, I* > 0, i = 1,...,m. Translating to the equivalent system (15), P* is corresponding to

p* = (s*i,y*^,i,..., y*,n, Ii*,s*, yh,..., y*,n, I**,...,

„* « « J* ) e Rm(n+2)

m' sm,V ' ' ' ' sm,n' m)

P in the interior of ro is called an endemic equilibrium, and it satisfies the following equilibrium equations:

0 = 9, (SD-Zfrjf,j (s* J]),

1 m 1 0=77Tx^f'j ft,1!)-^,

(1 + Sb)nj=l

0 = 1 (yl - y*,i), b

0 = 1 (fi,n-i - fi,n) , b

0=\fi,n -(s + £i +Yi )l*i. b

The basic reproduction number R0 is defined as the expected number of secondary cases produced by single infectious individual during its entire period of infectiousness in a completely susceptible population. For system (15), we can calculate it as the spectral radius of a matrix called the next generation matrix. Let

(1 + Sb)n

Cml (S°m)ß ml

(1 + Sb)n

lm l m

(1 + Sb)n

mm\ m)rmm

(1 + Sb)n

v = diag (S + ei + yt)

/ S + el +Yl 0

0 5 + £2 + y2

8+£m + Ym

Then the next generation matrix is

/ <11 ($)Pii

,(s°1)ßlr

(1 + Sb)n (S + ei +Yl ) (1 + Sb)n (S + em + ym)

fv-1 =

. Cm1 ^mt Cmm

( (1+Sb)n (S + et + ft) ''' (1+Sb)" (S + em + Ym) )

and hence, the basic reproduction number R0 is

R0 = p(fv-1) = max {|A|;A e o(fv-1)\, (24)

where p(-) and a(-) denote the spectral radius and the set of eigenvalues of a matrix, respectively. Since it can be verified that system (15) satisfies conditions (A^-(A5) of Theorem 2 of [17], we have the following proposition.

Lemma 1. For system (15), the disease-free equilibrium P0 is locally asymptotically stable if R0 < 1, while it is unstable if

Ro >1.

Following the method of [2], one defines a matrix

M = V-iF

Ii (Si)ßn

(1 + Sb)n (S + ei +Yl) (1 + Sb)n (S + ei + Yl)

I (Sl)ß m

(1+Sb)n (S + em + Ym)

^mm (S1 ) ßmm

(1+SbT (S + em + Ym) )

whose spectral radius has a similar threshold property to that of R0, since both of the nonnegative matrices fv-1 and M0 are irreducible, and hence from the Perron-Frobenius theorem [7] that their spectral radii are given by each of their simple eigenvalues. Thus, we

p(v-1f) = p(M0). Then the following lemma immediately follows.

Lemma 2. p(M0) < 1 ifandonlyif R0 < 1. 3. Main Results

The following main theorems are summarized in terms of system (15).

Theorem 3. Assume that the functions pt and fj satisfy assumptions (A 1)-(A3), and the matrix B = (pj)mxm is irreducible and R0 is defined by (24).

(i) If R0 < 1, then P0 is the unique equilibrium of system (15), and P0 is globally asymptotically stable in r0.

(ii) If R0 > 1, then P0 is unstable, and system (15) is uniformly persistent in r0.

Proof. Let us define M(S) = (frfr(Si)/(1 + Sb)n

(S + £i + Y,))mxm, where S = (Sv S2,..., Sm)T. Note that M(S0) = M0. Since B = (pij)mxm is irreducible, the matrix M0 is also irreducible.

First we claim that there does not exist any endemic equilibrium P in Q. Suppose that S = S0. Then we have

0 < M(S) < M0. Since nonnegative matrix M(S) + M0 is irreducible, it follows from the Perron-Frobenius theorem (see Corollary 2.1.5 of [7]) that p(M(S)) < p(M0) < 1. This implies that M(S)I = I has only the trivial solution

1 = 0, where I = (I1,..., Im)T. Hence the claim is true. Next we claim that the disease-free equilibrium P0 is globally asymptotically stable in r0. From the Perron-Frobenius (see Theorem 2.1.4 of [7]), the nonnegative irreducible matrix M0 has a strictly positive left eigenvector (w1,w2,... ,wm) associated with the eigenvalue p(M0) such that

(Wi,W2, ...,(*m)p (M0) = (wi,U>2,..., Wm) M°. (26)

Using this positive eigenvector, we construct the following Lyapunov function:

Computing the derivative of VDFE along the solutions of (15) in r0, we get

V fe = I

=1{\ + Sb)n (S + t, + y,i) ^ißiicij (sd

h (Si' Ij) -

1 (1 + Sb)n(S + £i + y)Ij -WiIi

I UißijCjj(S0) . . hd + Sbr (S+e, + yi) j Wii

= (w1,...,wm) [M0I-I] = [p(M°)-l](w1,...,wm)L

Thus, under the assumption R0 = p(M0) < 1, VDfe < 0, and Vd fe = 0 if and only if I = 0 and S = S0 = (S0V S02,...,S0m). Suppose that p(M0) = 1. Then it follows from the previous that Vpfe = 0 implies

((V1,...,(vm)M°I = (œ1,...,œm)I.

Hence, if S = S0, then (wj.....wm)M(S) < (wj,..., wm)M0 =

p(M0)(wJ,..., wm) = (wj,..., wm) and thus I = 0 is the only solution of (29). Summarizing the statements, we see that VFE = 0 if and only if I = 0 or S = S0, which implies that the compact invariant subset of the set where VFE = 0 is only the singleton P0. Thus, by LaSalle's invariance principle

[18], it follows that the disease-free equilibrium E0 is globally asymptotically stable in r0. If R0 = p(M0) > 1, then

= [p(M0)-1](cVi,CV2,...,Um)>0,

and then, by continuity, we can obtain

Vd fe = («i ,...,Wm)[M0I-l]>0, (31)

in a neighborhood of P0 in r0; then P0 is unstable.

Assume R0 = p(M0) > 1. By the uniform persistence result from [19] and a similar argument as in the proof of [2], the instability of P0 implies the uniform persistence of (15). This together with the dissipativity of (15) resulted from the forward invariant and compact property of r0 stated previously, implies which that (15) has an equilibrium in r0, denoted by P (see, e.g., Theorem D.3 in [20]). □

In what follows we prove that the endemic equilibrium P of system (15) is globally asymptotically stable when R0 >1. Throughout the paper, we denote

H(z) = z - 1- lnz.

Then H(z) > 0 for z > 0 and has global minimum at z = 1.

For convenience of notations, set fij = P*/*(S* ,1**), 1 < i, ] < m, and

Ißll - ß 21 1=1_

- ßl2 Iß2

-ßlm -ß 2.

ßml -ßm2

Iß ml

Then, B is also irreducible. One knows that the solution space of the linear system

has dimension 1 and

(v_, V2,... Vm) = (Cn,...,Cmm)

gives a base of this space, where Ckk > 0, k = 1,2, ...,m, is the cofactor of the kth diagonal entry of B. To get the global stability of P , the following assumptions in [15] are proposed:

(A4) : [^i(Si) - fi(S*)](Si -S*)<0 for S, = 81,8, e [0,S?j,

(A5) : ForS, = S*, [%(Si)-fi(S*)]-[fii(St,I^)-ftt(S*,1*)] < 0.

(A6) : For S{,Ij >0,

(fa (s;,infij (Si ,ij) - f,, (s.iDfj (s;j;)) (fa(s;,infij(s,,ij) fi,(s^)fj(s;,i;) ( b - b .

< 0. (36)

A difficult mathematical question for system (15) is that of

whether the endemic equilibrium P is unique when R0 >1 — *

and whether P is globally asymptotically stable when it is unique. Our main global stability result is given.

Theorem 4. Consider system (15). Assume that (Aa)-(A6) hold and the matrix B = (fij)mxm is irreducible. If R0 > 1,

then there is a unique endemic equilibrium P for system (15), and P is globally asymptotically stable in r0.

Proof. We show that P is globally asymptotically stable in r0, which implies that there exists a unique endemic equilibrium. Consider a Lyapunov function as

Vw = S - f

ee , '',

■i-fu (s;,inl

+ (1 + 8b)n

$ fi, (W )

1 (y,,j - y*i,j -yijln

••J

+I - j* - j* in —

• • , T*

This function has a linear part Vee expressed by

I(y,,j -yO) + i, j=_

L ee = Si + (1 + 8b)' First, calculating the derivatives of LEE, we obtain

Lee = Vi (Si)-(l + Sb)n (S + e, + y1)l1. (39)

Calculating the time derivative of VEE along the solutions of system (15) and using equilibrium equation (21), we have

vEe = 4e -ffjffy S, + (1+Sb)n = % (S,)-(1 + Sb)n (S + e, + Yl)l,

n V*- T*

1 fi» + T1'

Figure 1: Trajectories of S1 (t),I1(t),S2(t),and I2(t) for R0 = 0.051 < 1,and P0 = (3,0,0,0,3,0,0,0) is globally stable. S1(t), S2(t), I1(t),and I2(t) versus t are illustrated by (a), (b), (c), and (d). Initial values are S1(0) = 9,S2(0) = 1,yu(0) = 2, yia(0) = 2,y21(0) = 0, y2a(0) = 0,I1(0) = 6, and I2(0) = 2.

+ (1 + 8b)n

ßafn (s„Ij )y*i

=1 (1+8b)"y,,i

y,,i 1 ■

T + fyy*\ y

b bk=2 \ yi,k

yi,k-i

yini { ^ \ T*

- (8 + e, + Yi) Ii

= <Pi (Si )[1-

fii (S*,I*) fii (Si, I*)

f (S* r) y*ifij (Si,Ij) - (1 + 8b)n y ?ikyt,k-i

(1 + 8b)n * (1 + 8b)n * yhnI*

+ —?—ny*--?—y*

i'n y* Ii

+ (1 + 8b)n (8 + ei +yi) I*

fj (Si,Ij)fu (s*,In

-ijjij [s* ,I*)'

j=i fj (s* ,I* )fii (Si, I*)

-(1 + 8b)n (8 + e, +y,)I,

fii (s^I*

= (<Pi (St )-9i (S*)) (1 -

fii (Si, I*

n + 2 -

fn (s;,i;) y y**,ky,k-i

fn M) hy^tk-i

ytJ* I, fj (si,Ij)yh

y*J, I* fj (s*,I;)yhi

f,j (s, ,Ij)fi, (s;jn

ftj (s*,I;)fii (Si, I*)

0 h_ 0

Zu (t) (a)

0 h_ 0

yia (t) (b)

/2.1 (t)

0 hz 0

/2.2(f) (d)

Figure 2: Numerical simulation of (45) with R0 = 0.051 < 1; hence P0 = (3,0,0,0,3,0,0,0) is globally stable. Graphs (a) and (b) illustrate that Sj (t), yi i (t), yi 2(t) and Ii (t) will eventually towards to steady state. Graphs (c) and (d) illustrate that S2(t), y2 i (t), y2 2(t), and I2(t) will eventually towards to steady state. Initial values are Si (0) = 9,S2(0) = 1,^(0) = 2,yi2(0) = 2,y2i(0) = 0,y22 (0) = 0,Ii(0) = 6, and I2 (0) = 2.

It follows from the assumptions (A4)-(A5) that V^E can be estimated by

^EE < Iß.j i.j=l

G, (I,)-Gj (lj) + H

; j ; )

fii (s-Jj fii (S,,I*

'fa (^yji )

/ * / Wi

k=2 V yi,kyi,k-l v- r

S i,n i

ijfi, (s„i:)ftj (s*,i*) j*; fa (shWfij fab).

2.5 2 1.5 1

Figure 3: Trajectories of S1(t),I1(t),S2(t), and I2(t) for R0 = 1.67355 > 1, and P = (0.347644, 0.0760948, 0.0760948, 4.51674, 0.330353, 0.0765909, 0.0765909, 4.4678) is globally stable. S1(t), S2(t), I1(t), and I2(t) versus t are illustrated by (a), (b), (c), and (d). Initial values are S1(0) = 6,S2(0) = 2,yu1(0) = 3,y12(0) = 3,y21(0) = 0.1,y22(0) = 0.1,^(0) = 1.5,and/2(0) = 0.5.

fn (siInfn (s„Ij) 1 fii (s,,I* )f,j (s;,I*)

Ijfii (St,I*)ftj (S*,I**) I*fii (s;,I:)f,j (s,,Ij)

From the assumption (A6) and (32), we know that

vEe < m {Gt (I,)-Gj (Ij)\,

where G,(I,) = -I,/1* + ln(I,/I*).

Obviously, the equalities in (41) and (42) hold if and only

f„ (s*,Ii

fii (s,, I*

fn (s;,i,

fii (Si, I*

Wi (S,)-9i (S*)]=0,

fj (Si,Ij)fu (S*,I*) fa (s*,I;)fii (Si, I*)

fn (s,,I^)f,j (s;j;)Ij

fü (st,I*)ftj (S,,Ij)I*

That is, Sj = S*, It = I*, i = 1,2, ...,m. We can show that VEE and pj satisfy the assumptions of Theorem 3.1 and Corollary 3.3 in [21]. Therefore, the function

l=xv,ve]

is a Lyapunov function for system (15); namely, L'|(15) < 0 for P e r0. One can only show that the largest invariant subset, where L'|(15) = 0, is the singleton {P } by the same

/2.1 (t) (c)

/2.2(f) (d)

Figure 4: Numerical simulation of (45) with R^ = 1.67355 > 1; hence P' = (0.347644, 0.0760948, 0.0760948, 4.51674, 0.330353, 0.0765909, 0.0765909, 4.4678) is globally stable. Graphs (a) and (b) illustrate that Si (t), yi i (t), yi 2(t), and Ii (t) will eventually towards to steady state. Graphs (c) and (d) illustrate that S2(t), y2i{t), y22(t), and I2(t) will eventually towards to steady state. Initial values are Si(0) = 6, S2(0) = 2,yu(0) = 3,yu(0) = 3,y2i(0) = 0.1,y22'(0) = 0.1,Ij(0) = 1.5,andl2(0) = 0.5.

argument as in [2-5,13, 21]. By LaSalle's invariance principle, P is globally asymptotically stable in r0. This completes the proof of Theorem 4. □

Remark 5. We show a complete proof for global asymptotic stability of unique endemic equilibrium of system (15). In the case of fj(S,, Ij) = S,Ij, system (15) will reduce to the system studied in [14,22]. Here Theorem 4 extends related results in [14,22] to a result to a more general case allowing a nonlinear incidence rate. Our result also cover the related results of single group model in [13] for the case of f(S, I) = f(S)I.

4. Numerical Example

Consider the system (15) when m = 2, n = 2, fi(Si(t)) = 3 - S,, and fj(S,, Ij) = S,Ij, i, j = 1,2. One then has a two-group model as follows:

S[ (t) = 3-Si - [fiiiSi (t) Ii (t) + pnSi (t) I2 (i)], y'u (t) = -^gn [PiiSi (t) ii (t) + pnSi (t) 12 (i)] 1

-~yi,i (t), b

y[a (t) = \(yi,i (t)-yi,2 (t)), b

I'l (t) = iyU2 (t)-(8 + ei +Yi)Ii (t), b

S2 (t) = 3-S2 - [ß2iS2 (t) Ii (t) + ß22S2 (t) I2 (t)] 1

yi (t) =

(1+8b) 1

n [ß2iS2 (t)Ii (t)+ß22S2 (t)I2 (t)]

- $y2,1 (t) , b

y'2,2 (t) = 1 (y2,1 (t) - y2,2 (t)) , b

i2 (t) = 1y2,2 (t) -(8+ £2 + Y2) I2 (t) ■ b

If we choose parameters as ß11 = 5/24, ß12 = 1, ß21 = 1/36, ß22 = 1/2, 8 = 0.8, e1 = 2, e2 = 2, y1 = 1/4, and Y2 = 1/4, we can compute R0 = 0.051 < 1, and hence P0 = (3, 0, 0,0, 3, 0,0,0) is the unique equilibrium of system (45) and it is globally stable from Theorem 4 (see Figures 1 and 2).

On the other hand, if ßj are chosen as ß11 = 0.7, ß12 = 1, ß21 = 0.8, ß22 = 1,8 = 0.5, e1 = 0.02, e2 = 0.03, y1 = 0.05,

and y2 = 0,05, we can compute R0 = 1.67355 > 1, and hence

—* /

P = (0.347644, 0.0760948, 0.0760948, 4.51674, 0.330353, 0.0765909, 0.0765909, 4.4678) is the unique equilibrium of system (45) and it is globally stable from Theorem 4 (see Figures 3 and 4).

Acknowledgments

The authors wish to thank the reviewers for their valuable comments and suggestions that led to truly significant improvement of the paper. J. Wang is supported by the Science and Technology Research Project of the Department ofEducation ofHeilongjiang Province (no. 12531495).

References

[1] H. R. Thieme, Mathematics in Population Biology, Princeton University Press, Princeton, NJ, USA, 2003.

[2] H. Guo, M. Y. Li, and Z. Shuai, "Global stability of the endemic equilibrium of multigroup SIR epidemic models," Canadian Applied Mathematics Quarterly,vol. 14, no. 3,pp. 259-284,2006.

[3] H. Guo, M. Y. Li, and Z. Shuai, "A graph-theoretic approach to the method of global Lyapunov functions," Proceedings of the American Mathematical Society, vol. 136, no. 8, pp. 2793-2802, 2008.

[4] H. Shu, D. Fan, and J. Wei, "Global stability of multi-group SEIR epidemic models with distributed delays and nonlinear transmission," Nonlinear Analysis: Real World Applications, vol. 13, no. 4, pp. 1581-1592, 2012.

[5] J. Wang, Y. Takeuchi, and S. Liu, "A multi-group SVEIR epidemic model with distributed delay and vaccination," International Journal of Biomathematics, vol. 5, no. 3, Article ID 1260001, 18 pages, 2012.

[6] P. van den Driessche, L. Wang, and X. Zou, "Modeling diseases with latency and relapse," Mathematical Biosciences and Engineering, vol. 4, no. 2, pp. 205-219, 2007.

[7] A. Berman and R. J. Plemmons, Nonnegative Matrices in the Mathematical Sciences, Academic Press, New York, NY, USA, 1979.

[8] A. Korobeinikov and P. K. Maini, "A Lyapunov function and global properties for SIR and SEIR epidemiological models with nonlinear incidence," Mathematical Biosciences and Engineering, vol. 1, no. 1, pp. 57-60, 2004.

[9] M. Y. Li and J. S. Muldowney, "Global stability for the SEIR model in epidemiology," Mathematical Biosciences, vol. 125, no. 2, pp. 155-164, 1995.

[10] G. Huang and Y. Takeuchi, "Global analysis on delay epidemi-ological dynamic models with nonlinear incidence," Journal of Mathematical Biology, vol. 63, no. 1, pp. 125-139, 2011.

[11] P. van den Driessche and X. Zou, "Modeling relapse in infectious diseases," Mathematical Biosciences, vol. 207, no. 1, pp. 89103, 2007.

[12] P. van den Driessche, L. Wang, and X. Zou, "Impact of group mixing on disease dynamics," Mathematical Biosciences, vol. 228, no. 1, pp. 71-77, 2010.

[13] S. Liu, S. Wang, and L. Wang, "Global dynamics of delay epidemic models with nonlinear incidence rate and relapse," Nonlinear Analysis: Real World Applications, vol. 12, no. 1, pp. 119-127, 2011.

[14] Z. Yuan and X. Zou, "Global threshold property in an epidemic model for disease with latency spreading in a heterogeneous host population," Nonlinear Analysis: Real World Applications, vol. 11, no. 5, pp. 3479-3490, 2010.

[15] R. Sun and J. Shi, "Global stability of multigroup epidemic model with group mixing and nonlinear incidence rates," Applied Mathematics and Computation, vol. 218, no. 2, pp. 280286, 2011.

[16] R. K. Miller, Nonlinear Volterra Integral Equations, W. A. Benjamin, New York, NY, USA, 1971.

[17] P. van den Driessche and J. Watmough, "Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission," Mathematical Biosciences, vol. 180, no. 1-2, pp. 29-48, 2002.

[18] J. P. Lasalle, The Stability of Dynamical Systems, Regional Conference Series in Applied Mathematics, SIAM, Philadelphia, Pa, USA, 1976.

[19] H. I. Freedman, M. X. Tang, and S. G. Ruan, "Uniform persistence and flows near a closed positively invariant set," Journal of Dynamics and Differential Equations, vol. 6, no. 4, pp. 583-600, 1994.

[20] H. L. Smith and P. Waltman, The Theory of the Chemostat: Dynamics of Microbial Competition, Cambridge University Press, Cambridge, UK, 1995.

[21] M. Y. Li and Z. S. Shuai, "Global-stability problem for coupled systems of differential equations on networks," Journal of Differential Equations, vol. 248, no. 1, pp. 1-20, 2010.

[22] Z. Shuai and P. van den Driessche, "Impact of heterogeneity on the dynamics of an SEIR epidemic model," Mathematical Biosciences and Engineering, vol. 9, no. 2, pp. 393-411, 2012.

Copyright of Abstract & Applied Analysis is the property of Hindawi Publishing Corporation and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.