Available online at www.sciencedirect.com
Applied Mathematics and Computation 199 (2008) 285–291 www.elsevier.com/locate/amc
Global stability of an SIR epidemic model with constant infectious period Fengpan Zhang a, Zi-zhen Li
a,b,*
, Feng Zhang
c
a
b
Institute of Bioinformatics, School of Mathematics and Statistics, Lanzhou University, Lanzhou 730000, China Key Laboratory of Arid Agroecology under the Ministry of Education, Lanzhou University, Lanzhou 730000, China c Department of Mathematics, College of Science, Gansu Agricultural University, Lanzhou 730070, China
Abstract In this paper, we derive and study an SIR epidemic model with constant infectious period which is incorporated as a time delay. Both trivial and endemic equilibrium are found, and their stability is investigated. Using Lyapunov functional approach, sufficient conditions for global stability of endemic equilibrium is obtained. Ó 2007 Elsevier Inc. All rights reserved. Keywords: SIR epidemic model; Time delay; Constant infectious period; Global asymptotic stability; Lyapunov functional
1. Introduction One of the main issues in the study of behavior of epidemics is the analysis of steady states of the model and their stability. Generally, a model contains a disease-free equilibrium and one or more endemic equilibria. The stability of a disease-free steady state as well as the existence of other non-trivial equilibria can be determined using the so-called basic reproduction ratio, which quantifies how many secondary infections appear from a single infected put in a population of susceptible [1]. When the basic reproduction number is less than unity, the disease-free equilibrium is locally asymptotically stable, and therefore, the disease dies out after some period of time. Similarly, when the endemic equilibrium is a global attractor, epidemiologically this means that the disease will prevail and persist in a population. The stability of epidemic models has been studied in many papers. But most of them are concerned with local stability of equilibria. The fraction of papers that obtain global stability of these models is relatively low, especially the models with time delays. Beretta and Tackeuchi [2,3] studied the global stability of an SIR epidemic model with time delays. Recently, Takeuchi, Ma and Beretta [4] considered a delayed SIR
* Corresponding author. Address: Key Laboratory of Arid Agroecology under the Ministry of Education, Lanzhou University, Lanzhou 730000, China. E-mail addresses:
[email protected] (F. Zhang),
[email protected] (Z.-z. Li).
0096-3003/$ - see front matter Ó 2007 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2007.09.053
286
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
epidemic model with finite incubation time. They proved that the endemic equilibrium is globally stable if the length of incubation time is small. In this paper, we discuss the equilibrium and their stability of SIR epidemic model with constant infectious period which is incorporated as a time delay. Most of the models used in studies of the dynamics of epidemics have employed a simple description of the disease process. One particular assumption made is that the time for which individuals remain infectious can be described by an exponential distribution. This distribution is biologically unrealistic, however, because it corresponds to the assumption that the chance of recovery in a given time interval is independent of the time since infection. This leads to the distribution of infectious periods being too dispersed. In reality, infectious periods are fairly closely centered about the mean duration of infection. A constant infectious period and a gamma distributed infectious period are more realistic assumption to real world. The general infectious period distribution is described by its probability function. Let P ðtÞ be the probability of remaining infectious t units after becoming infectious. Assume that P ðtÞ is a non-increasing function with P ð0Þ ¼ 1 and P ð1Þ ¼ 0. These conditions allow many different P ðtÞ such as those corresponding to a constant infectious period, an exponentially distributed infectious period and a gamma distributed infectious period. Let QðtÞ be the probability of being alive at time t0 þ t given that an individual is alive at time t0 . We use QðtÞ ¼ elt since QðtÞ is the probability which is independent of the age of the individual. The function elt is the only QðtÞ for which the model is translation invariant, i.e., a semiflow [5]. Assume that the population size is constant and that the population is uniform and homogeneously mixing. In the SIR model, SðtÞ, IðtÞ, RðtÞ denote the faction of the population that are susceptible, infectious and removed, respectively. The constant contact rate b is the average number of contacts (sufficient for transmission) of an infective per unit time. Let the initial susceptible and removed fractions be S 0 > 0 and R0 > 0 and let I 0 ðtÞelt be the fraction of the population that was initially infectious and is still alive and infectious at time t. The function I 0 ðtÞ is a non-increasing function, I 0 ðtÞ > 0 and I 0 ðtÞ 6 I 0 ð0ÞP R 1 ðtÞ soltthat I 0 ðtÞ approaches 0 as t approaches 1. Because of deaths the average effective infectious period P ðtÞe dt is slightly less than the 0 R1 average infectious period s ¼ 0 P ðtÞdt. The integral equation for IðtÞ is Z t lt bSðxÞIðxÞP ðt xÞelðtxÞ dt; IðtÞ ¼ I 0 ðtÞe þ 0
where the second term is the sum of those who become infectious in the time interval ½0; t and are still alive and infectious at time t. The removed fraction satisfies Z t RðtÞ ¼ R0 ðtÞelt þ ½I 0 ð0Þ I 0 ðtÞelt þ bSðxÞIðxÞ½1 P ðt xÞelðtxÞ dt: 0
The second term represents those initial infectives who have recovered and are still alive at time t, and the third term represents those who became infectious in the interval ½0; t and are still alive at time t, but are no longer infectious. When the waiting time in the infective class is exponentially distributed, then P ðtÞ ¼ et=s and I 0 ðtÞ ¼ I 0 ð0Þet=s . The system reduces to the standard SIR model [6]. Their global stability results have been obtained [6]. When the infectious period is constant, then P ðtÞ is 1 on ½0; s and is 0 otherwise. Then the system reduces to a system of ordinary differential equations on ½0; s, and on ½s; 1Þ, it reduces to the system of delay-differential equations dSðtÞ ¼ l lSðtÞ bSðtÞIðtÞ dt dIðtÞ ¼ bSðtÞIðtÞ bels Sðt sÞIðt sÞ lIðtÞ dt dRðtÞ ¼ bels Sðt sÞIðt sÞ lRðtÞ dt
ð1Þ
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
287
where it is assumed that all newborns are susceptible and infection confers permanent immunity. Parameters in the system are as follows: l is a natural death and birth rate; b is the average number of adequate contacts of an infectious individual per unit time; s is the length of the infectious period. The term bels Sðt sÞIðt sÞ reflects the fact that an individual has recovered from infection and still are alive after infectious period s. SðtÞ þ IðtÞ þ RðtÞ ¼ N ðtÞ denotes the number of a population at time t. Under the assumption that birth and death rates are the same, the total population N ðtÞ evolves according to dNdtðtÞ ¼ lð1 N ðtÞÞ, and N ðtÞ ! 1 as t ! 1. The first two equations in system (1) do not depend on the third equation, and therefore this equation can be omitted without loss of generality. Hence, system (1) can be rewritten as dSðtÞ ¼ l lSðtÞ bSðtÞIðtÞ dt dIðtÞ ¼ bSðtÞIðtÞ bels Sðt sÞIðt sÞ lIðtÞ dt
ð2Þ
2. The stability analysis for equilibrium 2.1. Disease-free equilibrium and its stability Now we analysis system (2) by finding its equilibria and studying their stability. Steady states of system satisfy the following system of equations l lS bSI ¼ 0 bSI bels SI lIðtÞ ¼ 0
ð3Þ
It is easy to check that the trivial equilibrium E0 ¼ ð1; 0Þ and one non-trivial equilibrium for certain parameter values. We start with analyzing the behavior of the original system (2) near E0 . The eigenvalues of the linearization of system (2) near the steady state E0 are k1 ¼ l and k2 ¼ bð1 els Þ l. All parameters of the model are assumed to be positive. Therefore, for k1 ; k2 to be negative, i.e. for a disease-free equilibrium to be locally asymptotically stable, the following condition has to be required bð1 els Þ l < 0: Let us define the basic reproduction number of the infection as R0 ¼
bð1 els Þ : l
Using R0 we can get the following theorem indicating the stability of E0 ¼ ð1; 0Þ. Theorem 1. The disease-free equilibrium E0 is locally asymptotically stable if R0 < 1 and unstable if R0 > 1 2.2. Endemic equilibrium and its stability From the Theorem 1, when the basic reproduction number R0 > 1, the trivial equilibrium E0 of system (2) becomes unstable. Under this condition, there is a non-trivial equilibrium of system (2), i.e. endemic equilibrium, denoted by E ¼ ðS s ; I s Þ. It is easy to express this endemic equilibrium using R0 1 lðR0 1Þ E ¼ ; : R0 b The linearization matrix of system (2) near the steady state E ¼ ðS s ; I s Þ has the following characteristic equation ðk þ lÞ2 ðk þ lÞbðS s I s Þ þ ðk þ lÞbS s eðkþlÞs ¼ 0: We obtain the two eigenvalues of the linearization matrix of system (2) k1 ¼ l;
288
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
and ðk þ lÞ bðS s I s Þ þ bS s eðkþlÞs ¼ 0; let x ¼ ks: Rewrite it as x þ p qex ¼ 0; where p ¼ s½bðS s I s Þ l and q ¼ sbS s els . Hays (1950) showed that a necessary and sufficient condition for every root of the above characteristic equation having negative real part is 1) p < 1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2) p < q < a21 þ p2 where a1 is a root of a ¼ pðtanðaÞÞ; 0 < a < p Since p q ¼ s½bðS s I s Þ l ðsbS s els Þ ¼ s½bS s ð1 els Þ bI s l ¼ sbI s < 0: Thus we have Theorem 2. As long as R0 > 1 holds, the endemic equilibrium of system (2) is feasible and stays locally asymptotically stable. This local stability results can also be obtained by Grossman [7]. 3. The global asymptotic stability of the endemic equilibrium The objective of this section is to establish sufficient conditions under which the endemic equilibrium is globally stable. Let us begin by proving the following lemma which will be used in our further calculations. Lemma 1. Let the initial data be SðsÞ ¼ Sð0Þ > 0; IðsÞ ¼ Ið0Þ P 0 for all s 2 ½s; 0Þ and Rð0Þ > 0. Then solutions SðtÞ; IðtÞ and RðtÞ of system (1) are positive for all t > 0 Proof. We assume for contradiction that there exists the first time t1 such that Iðt1 ÞSðt1 Þ ¼ 0. Assume that Sðt1 Þ ¼ 0, then IðtÞ P 0 for all t 2 ½0; t1 . We get from equation for SðtÞ of system (1) dSðtÞ ¼ l > 0: dt t¼t1 j < 0, which is a contradiction. Since Sð0Þ > 0, for Sðt1 Þ ¼ 0 we must have dSðtÞ dt t¼t1 Next, we show that IðtÞ is positive for all t > 0. Similarly let t1 be the first time when Iðt1 ÞSðt1 Þ ¼ 0. Assume that Iðt1 Þ ¼ 0, then SðtÞ P 0 for all t 2 ½0; t1 . Then we have Z t1 elðt1 sÞ IðsÞSðsÞds: Iðt1 Þ ¼ b t1 s
Since SðsÞ > 0; IðsÞ > 0 for all s 2 ½0; t1 Þ, the right hand side of above equation is greater than zero; while the left hand side, Iðt1 Þ, is zero, which is a contradiction. Finally, we show the positivity of RðtÞ. Assume that t1 is the first time such that Rðt1 Þ ¼ 0, then dRðtÞ ¼ bels Iðt sÞSðt sÞ: dt t¼t1 Because it is established that SðtÞ and IðtÞ are positive for all t > 0, therefore dRðtÞ j is also positive. Since dt t¼t1 Rð0Þ > 0, for Rðt1 Þ ¼ 0 it is required to have dRðtÞ j < 0, which is a contradiction. This completes the dt t¼t1 proof. h
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
289
Lemma 2. Let the initial data for system (1) be SðsÞ ¼ Sð0Þ ¼ S 0 > 0; IðsÞ ¼ Ið0Þ ¼ I 0 > 0 for all s 2 ½s; 0Þ and Rð0Þ ¼ R0 > 0. Then SðtÞ 6 max f1; S 0 þ I 0 þ R0 g ¼ M Proof. Let N ðtÞ ¼ SðtÞ þ IðtÞ þ RðtÞ. From system (1) we know that dNdtðtÞ ¼ lð1 N ðtÞÞ. So N ðtÞ is a monotone function and N ðtÞ ! 1 as t ! 1. Suppose that N ð0Þ 6 1, then N ðtÞ 6 1. From Lemma 1 the solutions of system (1) are positive. It follows that SðtÞ 6 1 for all t 6 0. On the contrary, if N ð0Þ > 1 then N ðtÞ < N ð0Þ and hence, SðtÞ < N ð0Þ for all t > 0. The proof is complete. h Introducing new variables as u1 ¼ S S s ; u2 ¼ I I s ; u3 ¼ R Rs : After substituting these variables, the original system can be rewritten the following form du1 ¼ lu1 bSu2 bI s u1 dt du2 ð4Þ ¼ bSu2 þ bI s u1 bels I s u1 ðt sÞ bels Su2 ðt sÞ lu2 dt du3 ¼ bels I s u1 ðt sÞ þ bels Su2 ðt sÞ lu3 dt Now, proving that a trivial solution of system(4) is globally asymptotically stable, is equivalent to the fact that the endemic equilibrium Es of system (1) is globally asymptotically stable. We shall employ Lyapunov functional technique to prove it. Now let us introduce the following functional 1 1 2 V ðuÞ ¼ xðu1 þ u2 Þ þ ðu22 þ u23 Þ; 2 2 where x > 0 is an arbitrary real constant. The derivative of V is V 0 ðuÞ ¼ xðu1 þ u2 Þðu01 þ u02 Þ þ u2 u02 þ u3 u03 ¼ xðu1 þ u2 Þ½lu1 bSu2 bI s u1 þ bSu2 þ bI s u1 bels I s u1 ðt sÞ bels Su2 ðt sÞ lu2 þ u2 ½bSu2 þ bI s u1 bels I s u1 ðt sÞ bels Su2 ðt sÞ lu2 þ u3 ½bels I s u1 ðt sÞ þ bels Su2 ðt sÞ lu3 ¼ xlu21 ðxl þ l bSÞu22 lu23 þ ð2xl þ bI s Þu1 u2 xbels I s u1 u1 ðt sÞ xbels I s u2 u1 ðt sÞ xbels Su1 u2 ðt sÞ xbels Su2 u2 ðt sÞ bels I s u2 u1 ðt sÞ bels Su2 u2 ðt sÞ þ bels I s u3 u1 ðt sÞ þ bels Su3 u2 ðt sÞ: Choosing x as follows: x¼
bI s ; 2l
and applying Cauchy–Chwartz inequality to all ui uj -type terms, we obtain the following expression: 1 1 V 0 ðuÞ 6 xlu21 ðxl þ l bMÞu22 lu23 þ xbels I s ðu21 þ u21 ðt sÞÞ þ ðx þ 1Þbels I s ðu22 þ u21 ðt sÞÞ 2 2 1 1 1 ls 2 2 ls 2 2 þ xbe Mðu1 þ u2 ðt sÞÞ þ ðx þ 1Þbe Mðu2 þ u2 ðt sÞÞ þ bels I s ðu23 þ u21 ðt sÞÞ 2 2 2 1 ls þ be Mðu23 þ u22 ðt sÞÞ 2 1 1 ¼ xl xbels ðM þ I s Þ u21 xl þ l bM ðx þ 1Þbels ðM þ I s Þ u22 2 2 1 l bels ðM þ I s Þ u23 þ ðx þ 1Þbels I s u21 ðt sÞ þ ðx þ 1Þbels Mu22 ðt sÞ: 2
290
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
We choose Lyapunov function to be the form Z t Z u21 ðhÞdh þ ðx þ 1Þbels M U ðut Þ ¼ V ðuÞ þ ðx þ 1Þbels I s ts
t
u22 ðhÞdh;
ts
and hence, U 0 ðut Þ ¼ V 0 ðuÞ þ ðx þ 1Þbels I s u21 ðtÞ ðx þ 1Þbels I s u21 ðt sÞ þ ðx þ 1Þbels Mu22 ðtÞ ðx þ 1Þbels Mu22 ðt sÞ: Substituting inequality for V 0 ðuÞ, we get 1 1 U 0 ðut Þ 6 xl xbels ðM þ I s Þ u21 xl þ l bM ðx þ 1Þbels ðM þ I s Þ u22 2 2 1 ls l be ðM þ I s Þ u23 þ ðx þ 1Þbels I s u21 ðt sÞ þ ðx þ 1Þbels Mu22 ðt sÞ 2 þ ðx þ 1Þbels I s u21 ðtÞ ðx þ 1Þbels I s u21 ðt sÞ þ ðx þ 1Þbels Mu22 ðtÞ ðx þ 1Þbels Mu22 ðt sÞ 1 1 ls 2 ls ¼ xl xbe ðM þ I s Þ u1 xl þ l bM ðx þ 1Þbe ðM þ I s Þ u22 2 2 1 ls l be ðM þ I s Þ u23 þ ðx þ 1Þbels I s u21 ðtÞ þ ðx þ 1Þbels Mu22 ðtÞ 2 1 ls ls ¼ xl xbe ðM þ I s Þ ðx þ 1Þbe I s u21 2 1 1 ls ls ls 2 ½xl þ l bM ðx þ 1Þbe ðM þ I s Þ ðx þ 1Þbe M u2 l be ðM þ I s Þ u23 : 2 2 Therefore, 1 1 U 0 ðut Þ 6 ½xl bels ðxM þ 3xI s þ 2I s Þu21 ½xl þ l bM bels ð3ðx þ 1ÞM þ ðx þ 1ÞI s Þu22 2 2 1 ls ½l be ðM þ I s Þu23 : 2 The right-hand expression of the above inequality is always negative provided that 1 xM þ 3xI s þ 2I s 1 3ðx þ 1ÞM þ ðx þ 1ÞI s 1 M þ I s ln ln ln ; ; s > max : b b b 2xl 2xl þ l bM 2l A direct application of the Lyapunov–LaSalle type theorem (Theorem 2.5.3 of Kuang [8, p. 30] shows that limt!1 ui ðtÞ ¼ 0; i ¼ 1; 2; 3. We have proved the following theorem. Theorem 3. Let the initial conditions for system (1) be SðsÞ ¼ S 0 ðsÞ P 0; IðsÞ ¼ I 0 ðsÞ P 0 for s 2 ½s; 0Þ with ls Þ > 1. Then for any infectious period s S 0 ð0Þ > 0; I 0 ð0Þ > 0 and Rð0Þ > 0. Assume further that R0 ¼ bð1e l satisfying 1 xM þ 3xI s þ 2I s 1 3ðx þ 1ÞM þ ðx þ 1ÞI s 1 M þ I s ln ln ln ; ; s > max ; b b b 2xl 2xl þ l bM 2l where M ¼ maxf1; Sð0Þ þ Ið0Þ þ Rð0Þg, the endemic equilibrium e ¼ fR ; I g is globally asymptotically stable. 4. Conclusions This paper has been concerned with a disease model with constant infectious period, which is regarded as a type of time delay. Some previous efforts on incorporating delays in epidemic models have been mainly concerned on inclusion of latency periods which assumes that the force of infection at a present time is determined
F. Zhang et al. / Applied Mathematics and Computation 199 (2008) 285–291
291
by the number of infectives in the past. Others studied temporal immunity that allow individuals to have immunity for some time after they recover from infection [9]. However, it is a new idea to model constant infectious period which is incorporated to simple epidemic model SIR. It is epidemiologically realistic model that epidemic model SIR incorporates constant infectious period in this study. More realistic assumption is that both the incubation and infectious period are constant rather than the standard exponential distribution that is commonly implemented. Keeling and Grenfell [10] explored simple epidemic model SIR with constant incubation and infectious period by numerical simulation. Their results agree well with observations of disease transmission within households [11]. However, it is difficult to study analytically this more realistic model. We have analytically studied modified SIR model incorporating constant infectious period. Our results show that the introduction of distributed delays for non-cyclic infectious disease models do not change local asymptotic behavior of the models; that is, distributed delays can not lead to periodic solutions. Sufficient conditions have been given ensuring the existence of the endemic equilibrium for system and stability of the endemic equilibrium is investigated. By using Lyapunov functional technique, we have been able to show that under certain restrictions on the parameter values and infectious period, the endemic equilibrium is globally asymptotically stable, epidemiologically this means that the disease will prevail and persist in a population. Acknowledgements This work was supported by the National Social Science Foundation of China (No. 04AJL007) and National Natural Science Foundation of China (No. 03470298), and Specialized Research Fund for the Doctoral program of Higher Education of China (No. 200220730017). We are grateful to Prof. Shigui Ruan working at Miami University in USA for his constructive comments. References [1] Y. Takeuchi, W. Ma, E. Beretta, Global asymptotic properties of a delay SIR epidemic model with finite incubation times, Nonlinear Anal. 42 (2000) 931–947. [2] E. Beretta, Y. Takeuchi, Global stability of an SIR epidemic mode1 with time delays, J. Math. Biol. 83 (1995) 250–260. [3] E. Beretta, Y. Takeuchi, Convergence results in SIR epidemic mode1 with varying population sizes, Nonlinear Anal. 28 (1997) 1909– 1921. [4] O. Diekmann, J.A.P. Heesterbeek, Mathematical Epidemiology of Infectious Diseases: Model Building, Analysis and Interpretation, Wiley, NewYork, 2000. [5] Herbert W. Hethcote, David W. Tudor, Integral equation models for endemic infectious diseases, J. Math. Biol. 9 (1980) 37–47. [6] R.M. Anderson, R.M. May, Infectious Diseases of Humans: Dynamics and Control, Oxford University Press, 1991. [7] H.W. Grossman, Osillatory phenomena in a model of infectious diseases, Theor. Popul. Biol. 18 (1980) 204–243. [8] Y. Kuang, Delay-Differential Equations with Application in Population Biology, Academic Press, New York, 1993. [9] Yuliya N. Kyrychko, Konstantin B. Blyuss, Global properties of a delayed SIR model with temporary immunity and nonlinear incidence rate, Nonlinear Anal. (2005) 495–507. [10] M.J. Keeling, B.T. Grenfell, Disease extinction and community size: modelling the persistence of measles, Science 275 (1997) 65–67. [11] R.E. Hope-Simpson, Infectiousness of communicable diseases in the household, Lancet 2 (1952) 549–554.