# Stability and Hopf Bifurcation Analysis of an Epidemic Model by Using the Method of Multiple Scales.

1. IntroductionMathematical models describing the transmission of infectious diseases have played an important role in understanding the mechanism of disease transmission and controlling the spread of infectious diseases. In the literatures, many classical epidemic models have been proposed and studied, and many authors also attempt to develop more realistic mathematical models. In 1927, by using "compartment modeling," Kermack and Mckendrick [1] described an epidemic model and computed the theoretical number of infectious individuals. From then on, "compartment modeling" is used until now. Compartmental models are also called SIR models or SIRS models, where S, I, and R denote the number of susceptible, infectious, and recovered individuals. Some SIR and SIRS models have been analyzed in detail [2-7]. However, many diseases incubate inside the hosts for a period of time before the hosts become infectious. Therefore, the SEIR epidemic models, where E denotes the number of individuals who are infected but not yet infectious, are developed to investigate the role of the incubation period in disease transmission. Some authors have studied the SEIR epidemic models [8-12]. For some diseases, it is necessary that some infectious individuals are quarantined. Then, the SIQR and SEIQR models, where Q denotes the number of infectious individuals that have been quarantined and controlled, are studied [13-16].

Incidence rate plays an important role in the transmission of disease. In earlier models [17], based on the law of mass action, the bilinear incidence rate [beta]SI, where [beta] is the average number of contacts per infectious individual per day is used. Capasso and Serio [18] introduced a saturated incidence rate g(I)S into epidemic models, where g(I) = [beta]l/(1 + al),where [beta]l measures the infection force of the disease and 1/(1 + al) describes the psychological effect or inhibition effect from the behavioral change of the susceptible individuals with the increase of the infectious individuals.

Based on the work of Capasso and Serio [18], in this paper, we investigate an SIQR epidemic model with a saturated incidence rate g(x)S, where g(x) depends on the ratio of I and S. Therefore, we have

g = g (I/S) = [beta](I/S)1 + I/S = [beta]I/S + I. (1)

Then, the saturated incidence rate is [beta]SI/(S + I). Considering the effects of time delay, we obtain the following SIQR epidemic model:

[mathematical expression not reproducible], (2)

where S(t), I(t), Q(t), and R(t) denote the number of susceptible, infectious, quarantined, and recovered individuals at time t. A is the recruitment rate of the population. d is the natural death rate of the population. [[mu].sub.1] and [[mu].sub.2] are the extra disease-related death rate in classes I and Q, respectively. k is the rate constant for individuals leaving infectious compartment I for quarantine compartment Q x [c.sub.1] and [c.sub.2] are the removal rate constants from classes I and Q, respectively. In this model, we assume that the susceptible individuals were infected before time delay t which is the latent period. Then, the infected individuals become infectious individuals and some of infectious individuals are quarantined. In classes I and Q, some individuals are cured and removed.

The outline of this paper is as follows. In Section 2, the equilibria and their stability are analyzed. The critical value of time delay at which Hopf bifurcation occurs is obtained. In Section 3, by using an extended method of multiple scales, we obtained the normal form of the Hopf bifurcation. In Section 4, The above theoretical results are validated by numerical simulations with the help of dynamical software WinPP. In Section 5, conclusions are given.

2. Stability of Equilibria and Existence of Hopf Bifurcation

In this paper, we will study system (2) which is modeled by delayed differential equations. Noticing that ([2.sub.1]) and ([2.sub.2]) are uncoupled with ([2.sub.3]) and ([2.sub.4]), we will study the system consisting of ([2.sub.1]) and ([2.sub.2]):

[mathematical expression not reproducible]. (3)

System (3) always has a disease-free equilibrium [E.sub.1] = ([LAMBDA]/d,0). If [beta] > [c.sub.1] + k + d + [[mu].sub.1], system (3) has an unique endemic equilibrium [E.sub.2]([S.sup.*], [I.sup.*]), where

[mathematical expression not reproducible]. (4)

In order to investigate the stability of the two equilibria [E.sub.1] and [E.sub.2] simultaneously, we assume that E is one of the two equilibria. Then, the Jacobian matrix of system (3) at point E is given by

[mathematical expression not reproducible]; (5)

that is,

[mathematical expression not reproducible], (6)

where [mathematical expression not reproducible]. Equation (6) can be rewritten as

[mathematical expression not reproducible], (7)

where A = f + d, B = p - q, C = fp - qd, D = fd.

For [tau] = 0, the characteristic equation becomes

[[lambda].sup.2] + (A + B) [lambda] + C +D = 0. (8)

Remark 1. For r = 0, if C + D > 0 and A + B > 0, then the roots of (8) are real and negative, and then the equilibrium E is asymptotically stable.

For [tau] [not equal to] 0, if [lambda] = i[omega] ([omega] > 0) is a root of (7), then we have

[mathematical expression not reproducible]. (9)

Separating the real and imaginary parts, we have

B[omega] sin [omega][tau] + C cos [omega][tau] = [[omega].sup.2] - D, B[omega] cos [omega][tau] - C sin [omega][tau] = -[omega]A, (10)

which leads to the following fourth-degree polynomial equation:

[[omega].sup.4] + ([A.sup.2] - [B.sup.2] - 2D) [[omega].sup.2] + [D.sup.2] - [C.sup.2] = 0. (11)

Then, by discussing the distribution of the solutions of (11), we analyze the stability of the equilibrium and the existence of Hopf bifurcation. The distribution of the solutions of (11) can be divided into the following three cases.

Case 1. If

[A.sup.2] -[B.sup.2] - 2D > 0, [D.sup.2] -[C.sup.2] > 0, (12)

or

[DELTA] = [([A.sup.2] -[B.sup.2] - 2D).sup.2] - 4 ([D.sup.2] - [C.sup.2]) < 0, (13)

then (7) has no positive root.

Case 2. If

[D.sup.2] -[C.sup.2] < 0, (14)

then (7) has only one positive root [[omega].sub.+].

Case 3. If

[A.sup.2] - [B.sup.2] - 2D < 0, [D.sup.2] - [C.sup.2] > 0, [DELTA] > 0, (15)

then (7) has two positive roots:

[[omega].sub.[+ or -]] = [square root of 2]/2 [square root of [B.sup.2] -2D) [+ or -] [square root of [DELTA]]]], (16)

and no such solutions otherwise.

Noticing that [[omega].sub.+] in Case 2 is the same as [[omega].sub.+] in Case 3, we suppose that (7) has two positive roots w[+ or -]. Then, from (10), we can determine that

[mathematical expression not reproducible], (17)

at which (7) has a pair of purely imaginary roots [+ or -] i[[omega].sub.[+ or - ]]. Supposing that

[lambda]([tau]) = d([tau]) + i[omega]([tau]), (18)

then

d([[tau].sup.[+ or -].sub.j]) = 0, [omega]([[tau].sup.[+ or -].sub.j]) = [[omega].sub.[+ or -]]. (19)

Substituting [lambda]([tau]) into (7) and taking the derivative with respect to t, we have

[mathematical expression not reproducible], (20)

which, together with (10) and (16), leads to

[mathematical expression not reproducible] (21)

Thus, if [DELTA] [not equal to] 0, we have

[mathematical expression not reproducible]. (22)

Therefore, we can obtain the following theorem.

Theorem 2. Let [[omega].sub.[+ or -]], [[tau].sup.[+ or -].sub.j] (j = 0,1,...) be defined by (16) and (17), respectively. For system (3),

(i) if [A.sup.2] - [B.sup.2] - 2D > 0 and [D.sup.2] - [C.sup.2] > 0 or [DELTA] < 0, then equilibrium E of system (3) is asymptotically stable for all [tau] > 0;

(ii) if [D.sup.2] - [C.sup.2] < 0, then equilibrium E of system (3) is asymptotically stable when [tau] [member of] [0, [[tau].sup.+.sub.0]) and unstable when [tau] > [[tau].sup.+.sub.0]; system (3) undergoes a Hopf bifurcation at E when [tau] = [[tau].sup.+.sub.0];

(iii) if [A.sup.2] - [B.sup.2] - 2D < 0, [D.sup.2] - [C.sup.2] > 0, and [DELTA] > 0, then there is positive integer k, such that equilibrium E switches k times from stability to instability to stability; that is, E is asymptotically stable when

[mathematical expression not reproducible], (23)

and unstable when

[mathematical expression not reproducible]. (24)

System (3) undergoes a Hopf bifurcation at E when [tau] = [[tau].sup.[+ or - ].sub.j].

3. The Normal Form of Hopf Bifurcation

In this section, we suppose that system (3) undergoes a Hopf bifurcation from equilibrium [E.sup.*]([S.sup.*],[I.sup.*]). Time delay t is chosen as the bifurcation parameter and its critical value is [[tau].sub.c]. By using the method of multiple scales [19], we will derive the normal form of the Hopf bifurcation. By translating the equilibrium to the origin, system (3) is rewritten as the following equation:

[??] = F (x, [x.sub.[tau]]), (25)

where x = [(S(t),I(t)).sup.T] [member of] [R.sup.2] is the state variable vector and [x.sub.[tau]] = [(S(t-T[tau], I(t-T).sup.T].

According to the MMS, a monoparametric family of solution of the type is as follows:

[mathematical expression not reproducible], (26)

where [T.sub.k] = [[epsilon].sup.k]t (k = 0,2) and [epsilon] is a nondimensional parameter. The solution does not depend on slow scale [T.sub.1] = et because secular terms first appear at 0([[epsilon].sup.3]). Therefore, we assume a two-scale expansion of the solution of (25).

Bifurcation parameter t is ordered as

[tau] = [[tau].sub.c] + [[epsilon].sup.2][[tau].sub.[epsilon]]. (27)

The delay term in (25) can be further expanded as

[mathematical expression not reproducible], (28)

where [D.sub.k] = [partial derivative]/[partial derivative][T.sub.k].

By substituting (26), (27), and (28) into (25), expanding F, and balancing the coefficients of [[epsilon].sup.j] (j = 1,2,3.), a set of perturbative equations are obtained.

First, for [epsilon]-order terms, we have

[mathematical expression not reproducible], (29)

where

[mathematical expression not reproducible]. (30)

Equation (29) has the following general solution:

[mathematical expression not reproducible], (31)

where A is complex constant and u is the right eigenvector given by

[mathematical expression not reproducible], (32)

then

[mathematical expression not reproducible]. (33)

Next, for [[epsilon].sup.2]-order terms, we have

[mathematical expression not reproducible], (34)

where [mathematical expression not reproducible] are the second-order partial derivatives of F(x, [x.sub.[tau]]) with respect to x and [x.sub.[tau]].

Substituting (31) into (34), we obtain

[mathematical expression not reproducible]. (35)

Solving (35), it yields

[mathematical expression not reproducible], (36)

where vectors [z.sub.11] and [z.sub.11] [member of] [C.sup.2] are obtained by solving the following equations:

[mathematical expression not reproducible]; (37)

that is,

[mathematical expression not reproducible], (38)

where [mathematical expression not reproducible]. Then, we can obtain

[mathematical expression not reproducible], (39)

Finally, for [[epsilon].sup.3]-order terms, we get

[mathematical expression not reproducible], (40)

where [mathematical expression not reproducible] are the third-order partial derivatives of F(x, [x.sub.[tau]]) with respect to x and [x.sub.[tau]].

Substituting (31) and (36) into (40) and eliminating the secular terms, we can obtain the equations including [D.sup.2]A. Eliminating the coefficients of [D.sup.2]A by using the left eigenvectors and reabsorbing parameter e [20], the normal form is determined by

[mathematical expression not reproducible], (41)

where the expressions of coefficients [C.sub.[tau]][[tau].sub.[epsilon]] and [C.sub.111] are given by

[mathematical expression not reproducible], (42)

where [mathematical expression not reproducible] is the left eigenvector obtained by solving the following equations:

[mathematical expression not reproducible], (43)

where H denotes the transpose conjugate and T denotes the transpose.

To express the normal form in real form, a polar form representation for the complex amplitude is introduced:

A = 1/2 a[e.sup.i[theta]]. (44)

Substituting (44) into (41) and separating the real and imaginary parts in (41), the generalized amplitude and phase modulation equations are drawn:

[mathematical expression not reproducible], (45)

where [mathematical expression not reproducible].

From (45), we can get that the amplitude of the periodic solution is

a = 2 [square root of [R.sub.1]/[R.sub.111]. (46)

The stability of the periodic solution is determined by the sign of the eigenvalue of the Jacobian matrix of (451). The eigenvalue is

[lambda] = -2[R.sub.1]. (47)

Then, we can obtain the following theorem.

Theorem 3. The amplitude of the bifurcating periodic solution of system (3) is 2 [square root of -[R.sub.1]/[R.sub.111]; the stability of the bifurcating periodic solution is determined by [lambda]: the bifurcating periodic solution is stable (unstable), if [lambda] < 0 ([lambda] > 0).

4. An Example

In this section, we perform some numerical simulations to verify our analysis. For the following parameter values [LAMBDA] = 10, [beta] = 0.6, d = 0.02, [c.sub.1] = 0.08, k = 0.01, [[mu].sub.1] = 0.03, we can see that endemic equilibria [E.sup.*] (20.8333, 68.4524) of system (3) exist. Then, a periodic solution with frequency [omega] = 0.32395 bifurcates from the endemic equilibria at critical value [tau] = [[tau].sub.c] = 4.94223. By Theorem 2, endemic equilibria [E.sup.*] are asymptotically stable when [tau] [member of] [0, [[tau].sub.c]) and unstable when [tau] > [[tau].sub.c].

By using the dynamical software WinPP, some numerical simulations are given. First, by Theorem 2, endemic equilibria [E.sup.*] are asymptotically stable when [tau] [member of] [0, [[tau].sub.c]) as shown in Figure 1 and the time histories of the state variables Q and R are also given in Figure 1.

Endemic equilibria [E.sup.*] are unstable when [tau] > [[tau].sub.c]. When [tau] crosses critical value [[tau].sub.c], a stable periodic solution yields as shown in Figure 2.

By using the method developed in this paper, the amplitude of the periodic solution can be obtained. The relationship between the amplitude and bifurcation parameter [tau] is given in (45). Then, the analysis results are compared with that of numerical simulations as shown in Figure 3.

5. Conclusions

In this paper, the stability and Hopf bifurcation of an epidemic model with time delay are studied. The stability of disease-free equilibrium and endemic equilibrium is studied and the conditions of the stability and unstability of the two equilibria are obtained. Then, the critical value of [tau] at which a Hopf bifurcation occurs is obtained. By using an extended method of multiple scales, the normal form of Hopf bifurcation is obtained. In order to verify the theoretical results, the Runge-Kutta scheme is adopted to produce the numerical results. The comparison between analytical predictions and numerical results of the amplitude of Hopf bifurcation reveals a qualitatively and quantificationally excellent agreement.

In order to obtain the normal form of Hopf bifurcation occurring in delayed differential equations, we develop an extended method of multiple scales. Choosing time delay as a bifurcation parameter and using the method, one can get the normal form of Hopf bifurcation without the tedious computation encountered in the center manifold reduction. Moreover, when we develop the method, we choose n-dimensional delayed differential equation. Therefore, the method can be used for studying the Hopf bifurcation occurring in a general n-dimensional delayed differential equation.

Through observing the results, there are two contributions in this paper. First, we improve an epidemic model and introduce a nonlinear incidence rate which depends on the ratio of the numbers of susceptible and infectious individuals. Second, we extend the method of multiple scales. By using this method, the normal form of Hopf bifurcation occurring in n-dimensional delayed differential equations can be obtained.

http://dx.doi.org/10.1155/2016/2034136

Competing Interests

The authors declare that there are no competing interests regarding the publication of this article.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (11302072) and the Doctoral Research Fund of Henan Institute of Engineering ([D.sup.2]012021).

References

[1] W. O. Kermack and A. G. Mckendrick, "A contribution to the mathematical theory of epidemics," Proceedings of the Royal Society of London A, vol. 115, no. 772, pp. 700-721, 1927.

[2] Y. Takeuchi, W. Ma, and E. Beretta, "Global asymptotic properties of a delay SIR epidemic model with finite incubation times," Nonlinear Analysis: Theory, Methods & Applications, vol. 42, no. 6, pp. 931-947, 2000.

[3] T. Kuniya, "Global stability analysis with a discretization approach for an age-structured multigroup SIR epidemic model," Nonlinear Analysis: Real World Applications, vol. 12, no. 5, pp. 2640-2655, 2011.

[4] R. Xu and Z. E. Ma, "Global stability of a SIR epidemic model with nonlinear incidence rate and time delay," Nonlinear Analysis: Real World Applications, vol. 10, no. 5, pp. 3175-3189, 2009.

[5] R. Xu, Z. E. Ma, and Z. P. Wang, "Global stability of a delayed SIRS epidemic model with saturation incidence and temporary immunity," Computers & Mathematics with Applications, vol. 59, no. 9, pp. 3211-3221, 2010.

[6] C. C. McCluskey, "Global stability for an SIR epidemic model with delay and nonlinear incidence," Nonlinear Analysis: Real World Applications, vol. 11, no. 4, pp. 3106-3109, 2010.

[7] Y. Muroya, Y. Enatsu, and T. Kuniya, "Global stability for a multi-group SIRS epidemic model with varying population sizes," Nonlinear Analysis: Real World Applications, vol. 14, no. 3, pp. 1693-1704, 2013.

[8] H. Y. Shu, D. J. Fan, and J. 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.

[9] C. J. Sun and Y.-H. Hsieh, "Global analysis of an SEIR model with varying population size and vaccination," Applied Mathematical Modelling, vol. 34, no. 10, pp. 2685-2697, 2010.

[10] N. Yi, Q. L. Zhang, K. Mao, D. M. Yang, and Q. Li, "Analysis and control of an SEIR epidemic system with nonlinear transmission rate," Mathematical and Computer Modelling, vol. 50, no. 9-10, pp. 1498-1513, 2009.

[11] J. Zhang, J. Q. Li, and Z. E. Ma, "Global dynamics of an SEIR epidemic model with immigration of different compartments," Acta Mathematica Scientia, vol. 26, no. 3, pp. 551-567, 2006.

[12] X. Y. Zhou and J. A. Cui, "Analysis of stability and bifurcation for an SEIR epidemic model with saturated recovery rate," Communications in Nonlinear Science and Numerical Simulation, vol. 16, no. 11, pp. 4438-4450, 2011.

[13] L.-I. Wu and Z. L. Feng, "Homoclinic bifurcation in an SIQR model for childhood diseases," Journal of Differential Equations, vol. 168, no. 1, pp. 150-167, 2000.

[14] H. Hethcote, M. Zhien, and L. Shengbing, "Effects of quarantine in six endemic models for infectious diseases," Mathematical Biosciences, vol. 180, pp. 141-160, 2002.

[15] S. A. A. El-Marouf and S. M. Alihaby, "Global analysis of an epidemic model with non-linear incidence rate," Journal of Mathematics and Statistics, vol. 7, no. 4, pp. 319-325, 2011.

[16] Y. Z. Pei, S. Y. Liu, S. J. Gao, S. P. Li, and C. G. Li, "A delayed SEIQR epidemic model with pulse vaccination and the quarantine measure," Computers & Mathematics with Applications, vol. 58, no. 1, pp. 135-145, 2009.

[17] K. L. Cooke, "Stability analysis for a vector disease model," The Rocky Mountain Journal of Mathematics, vol. 9, no. 1, pp. 31-42, 1979.

[18] V. Capasso and G. Serio, "A generalization of the Kermack-McKendrick deterministic epidemic model," Mathematical Biosciences, vol. 42, no. 1-2, pp. 43-61, 1978.

[19] A. H. Nayfeh and D. T. Mook, Nonlinear Oscillations, John Wiley & Sons, New York, NY, USA, 1979.

[20] A. Luongo, A. Paolone, and A. Di Egidio, "Multiple timescales analysis for 1:2 and 1:3 resonant Hopf bifurcations," Nonlinear Dynamics, vol. 34, no. 3-4, pp. 269-291, 2003.

Wanyong Wang and Lijuan Chen

College of Science, Henan University of Engineering, Zhengzhou 451191, China

Correspondence should be addressed to Wanyong Wang; wangwanyong630@163.com

Received 7 May 2016; Accepted 2 August 2016

Academic Editor: Oleg V. Gendelman

Caption: Figure 1: The response of state variables of system (2) with [LAMBDA] = 10, [beta] = 0.6, d = 0.02, [c.sub.1] = 0.08, k = 0.01, = 0.03, [c.sub.2] = 0.02, [[mu].sub.2] = 0.04, and [tau] = 4.8 < [[tau].sub.c] = 4.94223.

Caption: Figure 2: The response of state variables of system (2) with [LAMBDA] = 10, [beta] = 0.6, d = 0.02, [c.sub.1] = 0.08, k = 0.01, = 0.03, [c.sub.2] = 0.02, = 0.04, and [tau] = 5.1 > [[tau].sub.c] = 4.94223.

Caption: Figure 3: Comparison between the theoretical solution (solid) and numerical simulation (dot) for the amplitude of periodic solution.

Printer friendly Cite/link Email Feedback | |

Title Annotation: | Research Article |
---|---|

Author: | Wang, Wanyong; Chen, Lijuan |

Publication: | Mathematical Problems in Engineering |

Article Type: | Report |

Date: | Jan 1, 2016 |

Words: | 3429 |

Previous Article: | Hybrid Models Based on Singular Values and Autoregressive Methods for Multistep Ahead Forecasting of Traffic Accidents. |

Next Article: | Disturbance Compensation Based Finite-Time Tracking Control of Rigid Manipulator. |

Topics: |