GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL WITH IMMIGRATION OF DIFFERENT COMPARTMENTS*

GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL WITH IMMIGRATION OF DIFFERENT COMPARTMENTS*

Available online at www.sciencedirect.com sCIENcE@DiREcT* Acta Mathematica Scientia 2006,26B(3):551-567 www.wipm.ac.cn/publish/ GLOBAL DYNAMICS OF A...

955KB Sizes 3 Downloads 75 Views

Available online at www.sciencedirect.com sCIENcE@DiREcT*

Acta Mathematica Scientia 2006,26B(3):551-567 www.wipm.ac.cn/publish/

GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL WITH IMMIGRATION OF DIFFERENT COMPARTMENTS* Li Jzanquan ( +&& ) Ma Zhien ( $b& ) Department of Mathematics, Xi'an Jiaotong University, Xi'an 710049, China E-mail: lwdzjOZ63.net; lianq- [email protected]; zhmaOmail.xjtu. edu. c n

Zhang Juan ( %4?, )

Abstract The SEIR epidemic model studied here includes constant inflows of new susceptibles, exposeds, infectives, and recovereds. This model also incorporates a population size dependent contact rate and a disease-related death. As the infected fraction cannot be eliminated from the population, this kind of model has only the unique endemic equilibrium that is globally asymptotically stable. Under the special case where the new members of immigration are all susceptible, the model considered here shows a threshold phenomenon and a sharp threshold has been obtained. In order to prove the global asymptotical stability of the endemic equilibrium, the authors introduce the change of variable, which can reduce our four-dimensional system t o a three-dimensional asymptotical autonomous system with limit equation. Key words SEIR model, population size dependent contact rate, compartment, infected individual, compound matrix 2000 MR Subject Classification

1

92D30, 34D05

Introduction

The incidence of a disease is the number of new cases per unit time, and it plays an important role in the study of mathematical epidemiology. Thieme and Castillo-Chavez [l] argued that the general form of a population-size-dependent incidence should be written as P X ( N ) $ I , where S is the number of susceptibles at time t , 1 the number of infectives at time t , and N the total population size at time t , then S I 5 N . ,l? is the probability per unit time of transmitting a disease between two individuals in contact. And X ( N ) is the probability for an individual to take part in a contact. In many articles, X(N) is also called the contact rate, and

+

A

PX(N) = C ( N ) ,which is the average number of adequate contacts of an individual per unit time, is said to be an adequate contact rate. An adequate contact is a contact which is sufficient for transmission of the infection from an infective to a susceptible. Many contact rate forms are used in the incidence term in deterministic epidemic models described by differential equations. 'Received October 31, 2004. This research is supported by the NNSF of China (19971066)

552

Vo1.26 Ser.B

ACTA MATHEMATICA SCIENTIA

For example, the standard incidence starts with the assumption that the adequate contact rate is a constant A. The bilinear incidence PSI = P N S I implies that the adequate contact rate is PN that is linearly proportional to the total population size N . This would be realistic when the total population size N is not too large because the number of adequate contacts made by an individual per unit time should increase as the total population size N increases. On the contrary, if the population size N is quite large, because the number of adequate contacts made by an infective per unit time should be limited, or grows less rapidly as N increases, the linear contact rate PN is not available and the constant adequate contact rate A may be more realistic. Hence the two above-mentioned adequate contact rates are actually two extreme cases for the total population size N being very small and very large, respectively. Other expressions have also been used for the incidence term. For human diseases, Anderson [2] argued that if C ( N ) = x N 6 is an appropriate contact rate, then the incidence is X N 6 S I / N . The data in Anderson’s article [2] for five human diseases in communities with population sizes from 1 000 to 400 000 imply that b is between 0.03 and 0.07. This also shows that the incidence given by C ( N ) S I / N with a contact rate such as C ( N ) = XN0.05 would be even better. Using a handling time argument in modeling the contact rate for venereal diseases (in analogy to Holling’s [3] derivation of a predator’s functional response to the amount of prey) suggests that the contact rate function C ( N ) should be Modeling the formation of short-time social complexes, Heesterbeek and Metz [4] derived an effective contact rate of the form C ( N ) = All above-mentioned adequate contact rate forms satisfy the following two assumptions

s.

I+bNyTm.

(HI) C ( N )is a nonnegative nondecreasing continuous function as N 2 0 and is continuously differentiable as N > 0.

D(N)= is a nonincreasing continuously differentiable function as N > 0. D ( 0 ) # 0. And C ’ ( N ) ID’(N)I # 0. Here we consider the SEIR model with the incidence term of the form C ( N ) S I / N ,where C ( N )is a general population-size-dependent function satisfying the assumptions (HI) and (Hz). Then the five above-mentioned adequate contact rates are all included in our general form. Note that the incidence C ( N ) S I / N becomes the standard incidence when C ( N )is a constant X and becomes the bilinear incidence PSI when C ( N ) = P N . And the incidence C ( N ) S I / N is a saturating form when C ( N ) = I + b N +BJbmN or C ( N )= as these two adequate coqtact rates grow linearly for small N and approach a constant for large N . This general adequate contact form C ( N ) has been considered by Brauer \5-81, Pugliese [9, 101, Thieme [ll,121, Zhou ahd Hethcote [13], and Zhou [14]. Population-size-dependent contact rates have been also used in AIDS models by Castillo-Chavez et al. [15] and Thieme and Castillo-Chavez [I]. But it has not been used in the SEIR model. Here we consider the SEIR model with the general population-size-dependent adequate contact rate. Let E represent the number of the exposed compartment (or class) of individuals in the latent period, the exposed is infected but not yet infectious and enters the class (or compartment) 1 of infectives after the latent period ends, infectives are infectious in the sense that they are capable of transmitting the infection. An infected individual is an individual who is in the exposed class E or in the infective class I . Communicable diseases may be introduced into a (H2)

+

9

g,

N0.3

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

553

population by the arrival of infected individuals from outside of the population. For example, travelers may return home from a foreign trip with an infection acquired abroad. So in this manuscript, we construct the model for disease transmission that includes the immigration of distinct compartments, that is, there are susceptibles, exposeds, infectives, and recovereds in the new number of immigration. We assume a constant flow A of new members into the population per unit time, of which a fraction q is exposed, a fraction p is infective, and a fraction b is recovered, so the fraction 1 - p - q - b is susceptible. Here q , p , and b are nonnegative constants and 0 5 q p b 5 1. Assuming that natural deaths occur at a rate proportional to the population size N , then the natural death rate term is p N , here p is the natural death rate constant. Under our assumptions, the infection cannot be eliminated from the population because there is a constant flow of new infected individuals when q f p > 0. In order to eradicate the disease, it would be necessary to isolate the fraction of arriving infected individuals. Our model has some features in common with models that include vertical transmission, but vertical transmission models normally include a flow of new infectives proportional to the number of infectives who are already in the population and thus may have a disease-free equilibrium [16, 22, 241. Greenhalgh [ 171 considered SEIR models that incorporate density dependence in the death rate. Cooke and van den Driessche [18] introduced and studied SEIRS models with two delays. Recently, Greenhalgh [19]studied Hopf bifurcations in models of the SEIRS type, with densitydependent contact rate and death rate. Li and Muldowney [20] and Li et al. [21], studied the global dynamics of the SEIR models with a nonlinear incidence, rate and with a standard incidence, respectively. Li et al. [22], analyzed the global dynamics of an SEIR model with vertical transmission and a bilinear incidence. Research on epidemic models of SEIR or SEIRS types with the general population-size-dependent adequate contact rate C(N ) satisfying the assumptions (HI) and (Hz) are scarce in the published reports. In the present article, we construct and study the SEIR epidemic model with a general contact rate and immigration of distinct compartments, give the rigorous treatment of the global stability of the endemic equilibrium, and establish the complete global dynamics. Especially, to establish the global stability of the endemic equilibrium, we introduce the change of variables by which our four-dimensional SEIR model can be reduced to a three-dimensional asymptotical autonomous differential system with a limit system.

+ +

2

Model Formulation

The total population size N ( t ) is divided into four distinct epidemic subclasses (compartments) of individuals which are susceptible, exposed, infectious, and recovered (with permanent immunity), with the denoted by S ( t ) ,E ( t ) , I @ ) ,and R(t),respectively. Our assumptions on the dynamical transfer of the population are demonstrated in the following diagram: PA

IJs

1

I

554

ACTA MATHEMATICA SCIENTIA

Vo1.26 Ser.B

Here q, p,and b are all nonnegative constants and 0 5 q + p + b 5 1. The parameters y are all positive constants. a is a nonnegative constant and represents the death rate because of disease (the disease-related death rate). p is the rate constant for natural death, y is the rate constant for recovery, and E is the rate constant at which the exposed individuals become infective, so that $ is the mean latent period. The recovered individuals are assumed to acquire permanent immunity, so there is no transfer from class R back to class S. The positive constant A is the constant recruitment rate into the population, so that A / p represents a carrying capacity, or maximum possible population size, rather than the population size N . C ( N )satisfying assumptions (HI) and (Hz) is the adequate contact rate. p,

E,

Then the SEIR model with the adequate contact rate C ( N ) and immigration of different compartment individuals is derived on the basis of the basic assumptions, by using the transfer diagram, and is described by the following system of differential equations

S'

S N

= (1 - p - q - b)A - PX(N)--I

-

pS,

S E' = qA + PX(N)--I - p E - E E , N I' = p A + t E - p I - CYI- 71, R' = bA

+y l

-

+

pR,

+

N(t) = S(t) E(t) I(t)

+ R(t).

Rom biological considerations, we study (1) in the closed set

l-

where

=

{ ( S , E , I , R )E R$

:OsS +

E + I + R = N 5 Alp},

R$ denotes the nonnegative cone of R4 including its lower dimensional faces. It can be

verified that r is positively invariant with respect to (1). d r and the interior of J? respectively.

0

r denote the boundary and

The continuity of the right side of (1) on N > 0 implies that solutions exist in and because solutions remain in r, they are continuous for all t 2 0 (see [28]). If D ( N ) E C'[O, $1

w),

(where D ( N ) = = then the right side of (1) is globally Lipschitzian in I? so that the initial value problem has a unique solution for all t 2 0 that depends continuously on the data. If D'(0) is unbounded, then D ' ( N ) is still bounded on [b,$1. Then in the interior of the region r, the right side of (1) is locally Lipschitzian, as a Lipschitz constant exists for every closed bounded subset of the interior of r. Thus for any point in the interior of I?, the initial value problem has a unique solution for all t 2 0. If 0 < q p b < 1, the solutions starting at points on the boundary of r enter the interior of r. If q p = 0, solutions starting on the S axis approach PO= ($, O,O, 0), but all the other solutions starting at points on the boundary of I? enter the interior of I?. If q p = 1 - 6, for all solutions starting at points on the boundary of I?, they enter the interior of r or stay on d r . So the initial value problem is well posed in the closed set r.

+ + +

+

No.3

3

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

555

Equilibria

Equilibria for system (1) can be found by setting the right sides of the four differential equations of (1) equal to zero, giving the algebraic system, (1 - q - p -

b ) A = S I D ( N )+pS,

+ D ( N ) S I= (p + E)E, P A + E E= ( p + a + y)I, qA

b A +?I = p R . Adding all equations in (2), we obtain 1

I = -(A - p N ) . a Let 6 = p

+ a + y and w = p +

E

, then R and E

(3)

can be expressed in terms of N as

Thus, S can also be denoted in terms of N as

Substituting (4) and (5) into the second equation in

+

( F ( N ) @ D ( N ) ) -1( A t a where

6W F ( N ) = --C(N) a€

(-

(a),we have

-

6W 0.5

-

1

A 6W + b) --D(N) - -. P E

(7)

The existence and number of equilibria correspond to those of the roots of the equation (6) in the interval [0,A / p ] . Let

Then, Ro is the number of secondary infections made by a topical newly infective individual entering the population a t the disease-free equilibrium during his or her average infective period. So ROis the basic reproduction number when new members of immigration are all susceptible, or p + q b = 0. Theorem 3.1 Suppose p q b = 0 or new members of immigration are all susceptible. Then, (i) The point PO= ( A / p , 0 , 0,O) E I?, the disease-free equilibrium of the system (1) and it exists for all nonnegative values of its parameters. PO is globally asymptotically stable when Ro 5 1 and it is unstable when Ro > 1. (ii) When Ro > 1, the solutions of the system (1) starting sufficiently close to POin r move away from POexcept those starting on the invariant S-axis which approach POalong this axis.

+

+ +

556

ACTA MATHEMATICA SCIENTIA

Vo1.26 Ser.B

(iii) When Ro > 1, the system (1) has a unique interior (endemic) equilibrium PI = ( 4 E, l , I I , R1) and its coordinates satisfy (3) - (5) when p = 0, q = 0, b = 0, where N = N1 is the unique root of the equation F ( N ) = 0 in the interval (0, A / p ) . Proof When p q b = 0, that is, p = q = b = 0, the equation (6) becomes F ( N ) ( Ap N ) = 0. From this, we see that either A - p N = 0 (then I = 0, which corresponds to the disease-free equilibrium) or F ( N ) = 0 (the endemic equilibrium). If I = 0, it is easy to deduce that S = E = R = 0; thus, the point PO= ( A / @O, , O , 0) E l? is the disease-free equilibrium of the system (1) and it exists for all nonnegative values of its parameters. When Ro 5 1, consider the function L = d3 + w I . Its derivative along the solutions of (1)

+ +

2,

is

The inequality holds because S 5 N 5 $, C ( N ) is nondecreasing, and Ro 5 1. Furthermore, the Liapunov-Lasalle theorem (cf. [28], pp.296-297) implies that all paths in r, approach the largest positively invariant subset of the set A4 where L' = 0. Here, M is the set where I = 0. On the boundary where I = 0 of ,'I we have E = 0, R' = -pR so R = R(0)e-Pt t 0 as t + +00, S' = A - ,US so S = f [S(O)-+ A as t t o o . Hence, all solution paths P in r approach the disease-free equilibrium Po. When & > 1, we can rewrite L' as

+

---f

If > 1, then L' > 0 for S is sufficiently close to $ except the invariant S-axis. Solutions starting sufficiently close to PO move away from the neighborhood of PO,except those on the invariant S-axis, on which (1) reduces to S' = A - pS and S ( t ) -+ $ as t -+ f00. Thus, all of these establish the first two conclusions of the theorem. From the assumptions ( H I ) and (H2), the function F ( N ) is increasing, so that the equation F ( N ) = 0 has at most a unique root. And because

6w 6w F ( A / p )= C ( A / p )- - = -(Ro r

bW lim F ( N ) = -C(0)

N-0

(Yt

E

6~

- QE - -(YE

A

- l),

6W lirn D ( N ) - -.

p N-0



Thus, when & > 1, F ( A / p ) > 0. From the assumptions (HI) and (Hz), by noting that b w - a ~= ( a + ~ + p ) ( ~ + p ) - a t> 0 and D ( N ) = if C(0)= 0, then either lim D ( N ) >_ 0 or

9,

N-+O

lirn D ( N ) = +oo so that either lim F ( N ) < 0 or lim F ( N ) = -00, respectively. If C(0)> 0,

N -0

N-0

then lim D ( N ) N -0

=

+00

N-0

so that lim F ( N ) = -m. Anyway when Ro > 1, in either case of N-+O

C(0) = 0 or C(0) > 0, the equation F ( N ) = 0 has only a positive root N1 in the interval (0, A / p ) , that is, the system (1) has a unique endemic equilibrium PI = (S1,El, I I ,R1) except disease-free equilibrium PO,where SI = a P + D ( NCX1A) ( A - P N 1 ) , El = =6 ( A - p N l ) , I I = ~ ( A - ~ N I ) , R1 = & ( A - p N 1 ) , that is, coordinates of PI satisfy (3) - ( 5 ) when p = 0,q = 0, b = 0, where

N

= N1.

+ +

Theorem 3.2 When p q b = 1, that is, all of new members of immigration are not susceptible, system (1) has a unique boundary equilibrium P 2 = (0, 6A 7A(wp+tq) P P6W and it is globally asymptotically stable.

e,v ,+

N0.3

557

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

+ +

Proof Let p q b = 1 in system (l),it is easy to see that P 2 is unique equilibrium of the system (1). The global asymptotical stability of point Pz can also be proved by using the function L = S. Theorem 3.3 When 0 < p+q < 1- b, system (1) only has equilibrium P* = (S*, E * ,I * , R*) and P* is in the interior of the feasible region I?. Proof When 0 < p q < 1 - b, system (1) has no disease-free equilibrium, the equation (6) becomes

+

Substituting ( 7 ) into (9), we obtain A 6W G ( N )= -C(N) Q€

so

(-

6W

-1

QE

6w

G ' ( N ) = -C'(N) a€

+ b PP -)E -

(-6w - 1 + b

=o. + QEA( (Aq-+PwNp)) - -6w E

A --D(N) P

(10)

QPA(EQ+ U P ) € ( A- P N ) ~.

-

QE

Note that

+ + P)(E + P )

(a y

- QE

- QPP

Q€

+b

From the assumptions (HI) and (Hz), the function G ( N ) is increasing, so that the equation (10) has at most a unique root. And we have

6W lim G ( N ) = -C(O) -

a€

N+O

(-

6W

6W

1+ b -

-

QE

6w + PP - ) C ( A / p ) - - + - A ( q + W P ) lim N-$ + p + y and w = + 1-1, Q

lim G ( N )= (1 - b N+AlP

E

By noting that S = Q

E

E

1

-= +m. A

-

pN

E

6w

1

-(w +u p ) - E = E (Q€(P+ q ) + a p p - 6,) E CY

1 =E

(

-

a41 - P - 4 ) - Q P U

-

P) - ( P + y ) ( P + 4 )

< 0. By the assumptions of (HI) and (H2), if C(0)= 0, then either lim C(N) 2 0 or lim N-0

+m, so that either lim G ( N ) < 0 or lirn G ( N ) = -m. If C ( 0 ) > 0, then lirn

N-+O

N -0

N-0

C")

=

N-0

= +m

so that lim G ( N ) = -m. Thus, in either case C(0)= 0 or C(0) > 0, the equation G ( N ) = 0 N-0

exists only one positive root N* in the interval (0, A / p ) . by notingthat 0 < p + q

2,

( -6W- l + b - - ) - DPP( NAo ) - T f - ( q6w + w p 6)

G(No)= -(l--)-D(No)6w ~p A QE

6

P

= (1 - p - b)-D(N0) A P

< 1 - b , then

QE

+ P-64

> 0,

f

P

EP

558

ACTA MATHEMATICA SCIENTIA

and so the unique root N* of the equation G ( N ) = 0 satisfies N * < NO5

E* = -6( A - p N * ) at

Vo1.26 Ser.B

$, thus

6P - PA - = -(No - N * ) > 0, E

ac

and

S=

a(1-p-q-b)A ~ r p D ( N * ) ( A- p N " )

+

> 0.

Thus, system (1) has the unique equilibrium P* = ( S * ,E*,I * , R*). The following results, which may be found in [23], are used in the next section for the study of the global stability of the endemic equilibrium P' and PI. Consider the following systems

where f and g are continuous and locally Lipschitz in x E R" and solutions exist for any positive time. System (12) is called asymptotically autonomous with limit system (13) if f ( t , x ) + g(z) as t + 00 uniformly for 2 E R".

Lemma 3.1 Let e be a locally asymptotically stable equilibrium of (13) and w be the w-limit set of a forward bounded solution ~ ( tof) (12). If w contains a point yo such that the solution of (13), with y(0) = yo converges to e as t 00, then w = { e } , that is, z ( t ) + e as ---f

t

+ 00.

Corollary 3.1 If solutions of system (12) are bounded and the equilibrium e of the limit system (13) is globally asymptotically stable, then any solution x ( t ) of system (12) satisfies x(t) --t e as t + m.

4

Global Stability of the Equilibria

In this section, we establish that all solutions of (1) in the interior of r converge to P* if 0 < p + q < 1 - b and t o PI if p + q + b = 0 and Ro > 1, respectively. This implies that PI and 0 P*, when they exist, are globally asymptotically stable in r, respectively. For this purpose, we first introduce the change of variable that can reduce the four-dimensional autonomous system (1) into a three-dimensional asymptotical autonomous system with a limit system, then prove that the equilibria of this limit system corresponding t o PI and P* are globally asymptotically stable, by using the geometric approach to global stability problems in Li and Muldowney [27]. This implies the main conclusion in Corollary 3.1 of the third section. For system ( l ) ,make the change of variable

M=-Y (S+E+Z)+- a+TR, a a

No.3

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

559

then (1) is equivalent to the following system

S’ = (1 - 9 - p - b)A - D ( N ) S I - pS,

+ D ( N ) S I- W E , I‘ = PA+ EE- 61, E‘ = qA

MI= (:+b)A-pM, N = &(S+ E + I + M ) , where w = p

+ E and 6 = p + + y. The fourth equation of the system (15) implies that Q

Substituting (16) into the first three equations of the system (15) gives the following threedimensional nonautonomous differential system

As N

+

& S + E + I + (2+ b) $

(

1” =

3 as t

+

+m, the system (17) is a three-dimensional

asymptotically autonomous differential system with limit system

s’ = (1 - p - 4 - b ) ~ D- ( ~ ) S-I ps, E’ = qA + D ( 3 ) S I - W E , I t = P A+ E E- 6 I , S+E+I+(2+6)$ It can be verified that the region

{

T = ( s , E , I ) E R ~ :+O ~ S + E + I I A / ~ } is positively invariant with respect to (18), where R$ denotes the nonnegative cone of R3 including its lower dimensional faces. Thus, system (18) is bounded. Lemma 4.1 Suppose that p q b = 0. Then, (i) If RO 5 1, then QO = (A/p,O,O)is the only equilibrium in T and it is globally asymptotically stable. If Ro > 1, then QObecomes unstable whereas Q1 = ( & , E I , I l ) , which is corresponding to the equilibrium PI of system (l),and emerges as a unique equilibrium in the

+ +

interior of T , where ROis defined by (8), S1 = p + D (AN 1 ) I 1, El = ~ I II I, = & ( A- p N l ) , and NI is the unique root of the equation F ( N ) = 0, where F ( N ) is defined by (7). (ii) When Ro > 1, the solutions of system (18) starting sufficiently close to QOin T move away from QOexcept those starting on the invariant S-axis approach QOalong this axis.

560

Vo1.26 Ser.B

ACTA MATHEMATICA SCIENTIA

Lemma 4.2

When 0 < p

+q < 1

-

b, the system (18) has the only equilibrium Q* =

(S*, E*,I*),which is corresponding t o the equilibrium P* of system (l),and Q* is in the interior of the feasible region T , where S* = u l r a+(D1 -( pN-.q)-(bA > pAN . ) , E' = &(No 6 - N * ) ,I* = ; ( A - / I N * ) , and N * is the unique root of the equation G ( N ) = 0 where G ( N ) is defined by (10). Lemma 4.3 The system (18) is uniformly persistent if either p q b = 0 and Ro > 1 or 0 < p i q < 1 - b holds.

+ +

+

Proof When 0 < p q < 1 - b , it is easy t o see that the vector field of system (18) is transversal to the boundary of T on all its faces. When p q + b = 0, the vector field of the system (18) is transversal to the boundary of T on all its faces except the S-axis, which is invariant with respect to (18) and on the S-axis the equation for S is = A - pS, which implies that S ( t ) + as t -++m. Therefore, when p q b = 0, QOis the only w-limit point on the boundary of T . The system (18) is said to be uniformly persistent (see [23, 251) if there exists a constant

+

4

+ +

c,O

< c < 1, such that any solution ( S ( t ) E, ( t ) ,I ( t ) ) with (S(O),E(O),I(O))& satisfies liminf i ( S ( t ) ,E ( t ) ,I ( t ) ) l 2 c t-+x

The conclusion of the lemma follows from Theorem 4.3 in [26]. As the maximal invariant set on the boundary d T of T is an empty set when 0 < p q < 1 - h and the singleton { Q o } and QO is isolated when p q b = 0 and Ro > 1, thus, the hypothesis ( H ) of [26] holds for (18). The condition for uniform persistence in Theorem 4.3 of [26] is equivalent to QObeing unstable when p q b = 0 and & > 1. Thus, when the condition p q b = 0 or Ro > 1 and the 0 condition 0 < p + q < 1 - b holds, system (18) is uniformly persistent in T .

+ +

+ +

+

+ +

The main aim of this section is to prove the following Theorem. Theorem 4.1 (i) When p + q + b = 0 and Ro > 1, the unique endemic equilibrium PI of the system (1) is globally asymptotically stable in the interior of r. Moreover, PI attracts all trajectories in r except those on the invariant S-axis that converge to PO along this axis. (ii) When 0 < p + q < 1 - b, the unique endemic equilibrium P" of system (1) is globally asymptotically stable in the interior of r. By inspecting the vector field given by ( l ) ,we see that, when 0 < p

+q

0

< 1 - b, all

trajectories starting from the boundary dT of I? enter the interior r of r, and when p + q + b = 0, they do except those solutions on the S-axis that converge to PO,along this invariant axis. Thus 0 we need only to prove that r is the attractive region of P' when 0 < p q < 1- b and PI when p q b = 0, respectively. In accordance with Corollary 3.1, it is sufficient to prove that the equilibrium Q* of (18) when 0 < p + q < 1 - b and Q1 of (18) when p + q b = 0, which are respectively corresponding to the equilibrium P* of (1) when 0 < p q < 1 - b and PI of (1) when p q b = 0, are globally asymptotically stable in the interior of the region T whenever they exist. For this purpose, we utilize a geometric approach t o the global stability problem developed in Li and Muldowney [27] and Smith [29]. We briefly outline a general mathematical framework of this geometric approach for proving global stability in the appendix. Here we apply the theory, in particular Theorem A.2 in the appendix, to prove the following theorem. In addition, it is remarked in [27] that, under the assumptions of Theorem A.2 in the appendix, the condition 42 < 0 also implies the local stability of the equilibrium 3, as, assuming the

+

+ +

+ +

+

+

N0.3

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

56 1

contrary, 5 is both the alpha and omega limit point of a homoclinic orbit that is ruled out by the condition < 0. From Theorem 4.1, together with Theorem 3.1, Theorem 3.2 and Theorem 3.3, we can make the following three conclusions: (1) when new members of immigration are all susceptible (that is, p q b = 0),there exists threshold behavior and the basic reproduction number Ro is a sharp threshold parameter, and the global dynamical behaviors of system (1) and the outcome of the disease are completely determined by Ro. In other words, under the condition p q b = 0, if Ro 5 1, the disease-free equilibrium POis globally asymptotically stable so that the disease dies out, while if Ro > 1, the endemic equilibrium P* is globally stable so that the disease remains endemic. (2) When all the new members of immigration are not susceptible ( p q = 1 - b), the disease remains endemic and all individuals will be infected. (3) When new members of immigration include all different compartment (susceptible, exposed, infective, recovered compartment), that is, 0 < p q < 1 - b, the infection cannot be eliminated from the population and the disease always remains endemic. Lemma 4.4 (i) When 0 < p q < 1 - b, Q’ is globally asymptotical stable. (ii) When p q b = 0 and Ro > 1, Q1 is globally asymptotically stable. Proof By Lemma 4.3, either p q b = 0 and Ro > 1 or 0 < p q < 1- b, there exists a compact set K in the interior of T that is absorbing (18). Thus, in the closed set T the system (18) satisfies the assumptions ( H 3 ) , (H4), and (Hs) in the appendix. Let z = (S,E , I ) and f(z) denote the vector field of (18). The Jacobian matrix J = associated with a general solution z ( t ) of (18) is

+ +

+ +

+

+

+ +

+

+ +

all

+

a12

a13

-6 where

The second compound matrix J[’1 of J can be calculated as follows:

where

562

ACTA MATHEMATICA SCIENTIA

vo1.26 Ser.B

The proof of the theorem consists of choosing a suitable vector norm I . I in R3 and a 3 x 3 matrix-valued function A ( s ) , such that the quantity 42 defined by (A.2) in the appendix is negative. We set A as the following diagonal matrix

(

A ( S , E , I )= diag 1, 7

.

-

7

7

)

Then A is C’ and nonsingular in the interior of T . Then

A f A-’

= diag

(0,

();

Therefore, the matrix B = A j A - l

,

();

)

= diag ( 0 , 2 - ,:

-

):

+ AJ[21A-1can be written in the block form

with

where

El I‘ a c11 = - - - - j l - 6 - I D ( f i ) - D / ( I q S ) , E I Q+Y E’ I’ a c22 = - - - - 6 - w + -D / ( N ) S I . E I a+r

+

Choosing the vector norm I . I in R (23 )- R 3 ~

The Lozinskii measure p ( B ) with respect t o 1 . I can be estimated as follows (see [32]):

where 91 = Pl(B11)

+ IB121,

Q2

=

IB211+

Pl(B2Z),

and lB2ll are matrix norms with respect to 11 vector norm, and p l denotes the Lozinskii measure with respect t o 11 vector norm. More specifically, 1B2ll = and noting the assumptions (HI), (Hz), and S, 1 5 Tfi,then lBl2l

D($)

+ - - aD ’ ( f i ) I a+-r

1 D ( k ) + D’(fi)f= i (rnD(I?,)’

= C’(fi) 2 0,

No.3

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

563

Thus

Noting that B 1 1 is a scalar, its Lozinskii measure with respect to any vector norm in p 1 is equal to B11. To calculate p1(B22), we add the absolute value of the off-diagonal elements to the diagonal one in each column of B22, and then take the maximum of two sums. Using D ' ( f i ) 5 0, (21), and w = p E , we have

+

E' I' Pl(B22) = - - - - 6 - p, E I Therefore,

+

D ( 8 ) + --D'(fi)I) ff+T EE E' I' g2 = - + - - - - p - 6. I E I Rewriting the last two equations of (18), we have g 1 = -p - w - I D ( f i )

ff

(

E' 9A .. S I -E+ w - - = DE( N ) - , E I' pA I +6--=---. 1

SI E ,

(22)

(24)

EE 1

Substituting (24) into (22) and (25) into (23), together with Ill(&) I 0 and D ( f i ) 2 0, we obtain g l I - E' - p - - - 9 A g 2 = E' - - p - - PA

E

E'

E

I '

Thus,

E'

P(B) I - P .

Along each solution ( S ( t )E, ( t ) ,I ( t ) )to (18) with (S(O),E(O),I ( 0 ) ) in the compact absorbing set K , which is in the interior of T , we have

which implies @ = limsup sup t-m

5

$ s," p ( B ( s ( s ,so)))

d s 5 -p.

Z ~ E K

Conclusion

This article has considered an SEIR model that incorporates the constant inflow of new individuals and the exponential natural death, as well as disease-related death, so that the total population size may vary in time. One of features of the SEIR model considered here is that the adequate contact rate has been extended to a very general form C ( N ) . As described in the first section, possible choices of the adequate contact function C ( N ) are a constant A, an increasing function such as a linear function P N , the Holling type I1 functional response and the Heesterbeek and Metz function l + b N +a JNm .

%,

564

ACTA MATHEMATICA SCIENTIA

Vo1.26 Ser.B

The other feature of the SEIR model considered here is that the new members into population include different compartment individuals which are susceptible, exposed, infective, and recovered. Under our assumptions, the infection cannot be eliminated from the population because of the constant inflow of new infected individuals; In order to eradicate the disease it would be necessary to isolate the fraction of arriving infected individuals. Our model has some features in common with models that include vertical transmission, but vertical transmission models normally include a flow of new infectives proportional to the number of infectives already in the population and thus may have a disease-free equilibrium [16, 22, 241. Especially, in the case where p + Q + b = 0, there exists a threshold phenomenon and the threshold is the basic reproduction number Ro of the system (1). The global dynamical behaviors of system (1) and the outcome of the disease are completely determined by &. And RO is the product of the constant & and the value C ( A / p ) of the adequate contact rate function C ( N ) a t N = ,; which is the population size a t the disease-free equilibrium, where the constant $ = is the average infectious period and 2 = E is the fraction of exposed fife individuals surviving the latent class E. That is, from expression (8) of Ro, it is easy to see the contribution of the contact rate to the basic reproduction number. If let C ( N )be the constant A€ of the SEIR A, then & defined by (8) becomes the basic reproduction number (fi+c)(/J+a+r) model with the standard incidence in [21]. If let C ( N )be the linear function P N , then Ro defined by (8) becomes the basic reproduction number I A of the SEIR model with the bilinear incidence PSI i n [20].

&

(fi+c)E+a+7)

References

8

9

10 11

Thieme H R, Castillo-Chavez C. On the role of variable infectivity in the dynamics of the human immunodeficiency virus epidemic. In: Castillo-Chavez C, Ed. Mathematical and Statistical Approaches to AIDS Epidemiology (Lect Notes Biomath vo1.83). Berlin, Heidelberg, New York: Springer, 1989. 157-177 Anderson R M . Transmission dynamics and control of infectious diseases. In: Anderson R M, May M R, Ed. Population Biology of Infectious Diseases, Life Sciences Research Report 25, Dahlem Conference. Berlin: Springer, 1982. 149-176 Dietz K , Overall population patterns in the transmission cycle of infectious disease agents. In: Anderson R M, May R M, Ed. Population Biology of Infectious Diseases. Belin, Heidelberg, New York: Springer, 1982 Heesterbeek J A P, Metz J A J. The saturating contact rate in marriage and epidemic models. J Math Biol, 1993, 31: 529-539 Brauer F, van den Driessche P. Models for transmission of disease with immigration of infectives. Math Biosci, 2001, 171: 143-154 Brauer F. Models for universally fatal diseases 11. In: Busenberg S, Martelli M, Ed. Differential Equations Models in Biology, Epidemiology and Ecology. Lecture Notes in Biomathematics 92. Belin, Heidelberg, New York: Springer, 1991. 5 7 4 9 Brauer F. Models for the spread of universally fatal biseases. J Math Biol, 1990, 28: 451-462 Brauer F. Epidemic models in populations of varying size. In: CastillO-Chavez C, Levin S A, Shoemaker C, Ed. Mathematical Approaches to Problems in Resource Management and Epidemiology. Lecture Notes in Biomathematics 81. Belin, Heidelberg, New York: Springer, 1989. 109-123 Pugliese A. An SEI epidemic model with varying population Size. In: Busenberg S, Martelli M, Ed. Differential Equation Models in Biology, Epidemiology and Ecology. Lecture Notes in Biomathematics 92. Belin, Heidelberg, New York: Springer, 1991. 121-138 Pugliese A. Population models for diseases with no recovery. J Math biol, 1990, 28: 65-82 Thieme H R. Persistence under relaxed point-dissipativity (with Application to an Endemic Model). SIAM J Math Anal, 1993, 24: 407-435

No.3

565

Zhang et al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

12 Thieme H R. Epidemic and demographic interaction in the spread.of potentially fatal diseases in growing populations. Math Biosci, 1992, 111: 99-130 13 Zhou J, Hethcote H W. Population size dependent incidence in models for diseases without immunity. J Math Biol, 1994, 32: 809-834 14 Zhou J. An epidemiological model with population size dependent incidence. Rocky Mt J Math, 1993, 24: 85-101 15 Castillo-Chavez C, Cooke K , Huang W, Levin S A. On the role of long incubation periods in the dynamics of AIDS I: Single Population Models. J Math Biol, 1989, 27: 373-398 16 Busenberg S, Cooke K L. Verticl Transmitted Diseases: Models and Dynamics. Biomathematics, Vol 53. Berlin: Springer, 1993 17 Greenhalgh D. Some results for a SEIR epidemic model with density dependence in the death rate. IMA J Math Appl Medicine & Biology, 1992, 9: 67-106 18 Cooke K, van den Driessche P. Analysis of an SEIRS epidemic model with two delays. J Math Biol, 1996, 35: 240-260 19 Greenhalgh D, Hopf bifurcation in epidemic models with a latent period and nonpermanent immunity. Math Comput Modelling, 1997, 25(2): 85-107 20 Li M Y, Muldowney J S. Global stability for the SEIR model in epidemiology. Math Biosci, 1995, 125: 155-164 21 Li M Y, Graef John R, Wang L. Global dynamics of a SEIR model with varying total population size. Math Biosci, 1999, 160: 191-213 22 Li M Y, Smith H L, Wang L. Global dynamics of an SEIR epidemic model with vertical transmission. SIAM J Appl Math, 2001, 62(1): 58-69 23 Thieme H R. Convergence results and a Poicare-Bendixson trichotomy for asymptotically automous differential equations. J Math Biol, 1992, 30: 755-763 24 Brauer F. Models for Diseases with Vertical Transmission and Nonlinear Population Dynamics. Math Biosci, 1995, 128: 13-24 25 Butler G, Waltman P. Persistence in dynamical systems. Proc Amer Math SOC,1986, 96: 425-430 26 Freedman H I, Tang M X, Ruan S G. Uniform persitence and flows near a closed positively invariant set. J Dynam Diff Equat 1994, 6: 583-600 27 Li M Y, Muldowney J S. A geometric approach to global stability problems. SIAM J Math Anal, 1996, 27(4): 1070-1083 28 Hale J K. Ordinary Differential Equations. New York: Wiley-Interscience, 1969 29 Smith R A. Some applications of Hausdorff dimension inequalities for ordinary differential equations. Proc Roy SOCEdinburgh Sect A, 1986, 104: 235-259 30 Pugh C C, Robinson C. The C’ closing lemma, including hamiltonians. Ergod Theor Dynam Sys, 1983, 3: 261-313 31 Pugh C C. An improved closing lemma and the general density Theorem. Amer J Math, 1976, 89: 1010-1021 32 Martin R H, JR. Logarithmic norms and projections applied to linear differential systems. J Math Anal Appl, 1974, 45: 432-454 33 Li Y, Muldowney J S. On Bendixson’s criterion. J Differential Equations, 1994, 106: 27-39

Appendix Geometric Approach to Global-Stability Problem Let the mapping x H f (x)from an open subset T c Rn to Rn be such that each solution z ( t ) to the differential equation 5’

= f(x)

(-4.1)

is uniquely determined by its initial value x ( 0 ) = 20,and denote this solution by z(t,zo). Assume that (H3) T is simply connected; (H4) there is a compact absorbing set K c T ;

566

ACTA MATHEMATICA SCIENTIA

Vo1.26 Ser.B

(H5) 5 is the only equilibrium of (A.l) in T. The Global Stability Problem arising in [27] is to find conditions under which the global stability of IC with respect to T is implied by its local stability. In [27], the main idea of this geometric approach to the global stability problem is as follows: Assume that (A.l) satisfies a condition in T , which precludes the existence of periodic solutions and suppose that this condition is robust, in the sense that it is also satisfied by ordinary differential equations that are C'-close to (A.1); Then every nonwandering point of (A.l) is an equilibrium, as otherwise, by the C1 closing lemma of Pugh [30, 311, we can perturb (A.l) near such a nonequilibrium nonwandering point to get a periodic solution. As a special case, every omega limit point of (A.l) is an equilibrium. This implies that 5 attracts points in T . As a consequence, its global stability is implied by the local stability. For this purpose, we introduce the notion of robustness of a Bendixson criterion under local C1perturbations of f . Definition 1 A point zo E T is wandering for (A.l) if there exists a neighborhood U of 50 and t* > 0 such that U n x ( t , 17) is empty for all t > t' Thus, for example, any equilibrium, alpha limit point, or omega limit point is nonwandering. Definition 2 A function g E C'(T -+ R") is called a C1 local E-perturbation of f at xo E T if there exists an open neighborhood U of xo in T such that the support supp(f -9) c U and If - g(cl < E , where

and 1 . 1 denotes a vector norm on R" and the operator norm that it induces for linear mappings from R" to R". Definition 3 A Bendixson Criterion for (A.1) is a condition satisfied by f , which precludes the existence of nonconstant periodic solutions to (A.1). Definition 4 A Bendixson criterion is said to be robust under C' local perturbations of f at zo if, for each sufficiently small E > 0 and neighborhood U of 2 0 , it is also satisfied by each C' local €-perturbations g such that supp(f - g) c U . The following is the local version of the global-stability principle proved by [27]. Theorem A . l Suppose that assumptions (H4) and (H5) hold and that (A.l) satisfies a Bendixson criterion that is robust under C1 local perturbations of f at all nonequilibriurn nonwandering points for (A.l), then 3 is globally asymptotically stable with respect to T , provided it is stable. The following is to introduce the quantity q 2 given in [27]. Assume that (A.l) has a compact absorbing set Ii' c T.Then every solution z ( ~ , z oof ) (A.l) exists for any t > 0. The quantity 9 2 is well defined:

where

af [21 A-'

B = AfA-' + A -

BX

and x I+ A ( x ) is a ( g ) x ( g ) matrix-valued function, which is C' in T and A f = ( D A ) ( f ) or, equivalently, A f is the matrix obtained by replacing each entry aij in A by its directional

N0.3

Zhang e t al: GLOBAL DYNAMICS OF AN SEIR EPIDEMIC MODEL

567

derivative in the direction of f, !?.ff2] az is the second additive compound matrix of the Jacobian of f. If 1 . 1 is a vector norm on R(y),then p ( B ) is the Lozinskii measure with respect matrix to I .I, which is defined by ( I hBI - 1 p ( B ) = lim h ' h+O+

+

Under assumptions (H3) and (H4), [27] proved that q 2 < 0 is a Bendixson criterion for (A.l) and it is robust under C1 local perturbations of f at all nonequilibrium nonwandering points for (A.l) by means of the local version of C1 closing lemma of Pugh [30, 311. Then we have the following theorem Theorem A.2 Under assumptions (H3), (H4), and (Hs), the unique equilibrium a: is globally asymptotically stable in T if ij2 < 0. The criterion & < 0 provides the flexibility of a choice of (y) x (y) arbitrary function in addition to the choice of vector norms 1. 1 in deriving suitable conditions. It has been remarked in [27] that, under the assumptions of Theorem A.2, the classical result of Lyapunov, the classical Bendixson-Dulac condition, and the criterion in [33] are obtained as special cases (cf.[27]). In addition, it has also been stated in [27] that, under the assumptions of Theorem A.2 in the appendix, the condition 172 < 0 also implies the local stability of the equilibrium 5, because, assuming the contrary, 5 is both the alpha and the omega limit point of a homoclinic orbit that is ruled out by the condition 4 2 < 0.