 Research article
 Open Access
 Published:
Stability and bifurcations in a discretetime epidemic model with vaccination and vital dynamics
BMC Bioinformatics volume 21, Article number: 525 (2020)
Abstract
Background
The spread of infectious diseases is so important that changes the demography of the population. Therefore, prevention and intervention measures are essential to control and eliminate the disease. Among the drug and nondrug interventions, vaccination is a powerful strategy to preserve the population from infection. Mathematical models are useful to study the behavior of an infection when it enters a population and to investigate under which conditions it will be wiped out or continued.
Results
A discretetime SIS epidemic model is introduced that includes a vaccination program. Some basic properties of this model are obtained; such as the equilibria and the basic reproduction number \(\mathcal {R}_0\). Then the stability of the equilibria is given in terms of \(\mathcal {R}_0\), and the bifurcations of the model are studied. By applying the forward Euler method on the continuous version of the model, a discretized model is obtained and analyzed.
Conclusion
It is proven that the diseasefree equilibrium and endemic equilibrium are stable if \(\mathcal {R}_0<1\) and \(\mathcal {R}_0>1\), respectively. Also, the diseasefree equilibrium is globally stable when \(\mathcal {R}_0\le 1\). The system has a transcritical bifurcation when \(\mathcal {R}_0=1\) and it might also have perioddoubling bifurcation. The sufficient conditions for the stability of equilibria in the discretized model are established. The numerical discussions verify the theoretical results.
Background
The spread of infectious diseases in populations and how to control and eliminate them from the population are important and necessary subjects. Mathematical models are introduced to study what happens when an infection enters in a population, and under which conditions the disease will be wiped out from population or persists in population. The literature about mathematical epidemic models that have been constructed and analyzed for various types of diseases is very rich; see, for example, [1, 2]. Such models can be formulated either as continuoustime models by differential equations or as discretetime ones by difference equations. Recently, discrete models have gotten more interest because epidemic data are collected in discrete time intervals and numerical schemes also use discretization for solving differential equations. Moreover, the discretetime epidemic model exhibits more complex dynamics. Allen [3] studied some discretetime SI, SIR, and SIS epidemic models and showed that the simple discretetime SI and SIR epidemic models without births or deaths mimic the behavior of the continuoustime models, while the behavior in the discretetime SI, SIR, and SIS models with recovery or births differ from their continuous analogs. CastilloChavez and Yakubu [4] have investigated a discrete SIS model with complex dynamics. Brauer et al. [5] introduced a discreteepidemic framework and highlight, emphasizing the final size of the epidemic, the similarities between singleoutbreak comparable classical continuoustime epidemic models, and the discretetime models. Farnoosh and Parsamanesh [6] investigated the stability and bifurcation in a discrete SIS model with bilinear incidence. The vaccination program implements not only on susceptible individuals but also on newcomers to the population. Parsamanesh and Mehrshad [7] performed a similar investigation on an SIS model with a temporary vaccination program with standard incidence. A discrete SIRS epidemic model with vaccination and a general infection probability function is investigated by Xiang et al. [8]. The vaccination performs only on susceptible individuals but not on the recruitment of the population. The local and global stabilities of diseasefree equilibrium were derived as well as the local stability of the endemic equilibrium.
The discrete models have usually constructed either directly from properties of the disease and the population in which the propagation occurs, or by discretizing continuous model by a method based on numerical analysis such as the forward/backward Euler method, Mickens’ nonstandard method, or a mixed type formula which uses implicit and explicit Euler methods. For example, Roeger and Barnard [9] applied the central difference method to a continuous SIR epidemic model and showed the local stability of the equilibria. Liu et al. [10] presented four discrete epidemic models with the nonlinear incidence rate, using the forward Euler and backward Euler methods. They discussed the effect of two discretizations on the stability of the endemic equilibrium for these models. They concluded the forward Euler method makes the models have much richer dynamical behavior than the continuous models while the backward Euler method preserves the global asymptotic stability of endemic equilibria. Aranda et al. [11] provided a discrete model by discretization of the continuoustime model for transmission of Babesiosis disease in bovine and tick populations. They proved local stability of the equilibria and global stability of the diseasefree equilibrium and obtained similar conclusions as they got in the continuous case. In particular, Mickens [12] consider the nonstandard discretization for numerical methods which ensures positivity of the solutions of the difference equations. Recently, Izzo and Vecchio [13] obtained a discrete epidemic model by using a mixed type formula, which uses implicit and explicit Euler methods for discretization. Their discretization is similar to Mickens’ nonstandard one. But their discretization showed some good properties for positivity, boundedness, and global behavior of the solution. There are some works devoted to the local stability of a discrete model and few others investigated the global stability of the model. Hu et al. [14] established the criteria on the local stability of the diseasefree equilibrium and endemic equilibrium for a class of SIRS epidemic model with a nonlinear incidence, by using the linearization method and the expression of roots of the cubic polynomial equation. Ma et al. [15] studied the global stability of the endemic equilibrium of a discrete SIR model and get sufficient conditions by using the comparison principle. Cui and Zhang [16] gave sufficient conditions for the global stability of a discrete SIR model derived by employing a nonstandard finite difference scheme. Van den Driessche and Yakubu [17] proved the global stability and persistence for some discrete epidemic models that are formulated for SEIR infections, cholera in humans, and anthrax in animals, under some demographic assumptions. Among these models, the susceptible infected susceptible (SIS) epidemic models are one of the wellknown types of epidemic models. To consider the effect of vaccination as an efficient strategy to control and eliminate infections, it is possible to add a compartment for the vaccinated individuals to the SIS model and obtain the SIS epidemic model with vaccination, namely SIVS epidemic model [18, 19]. These models may be deterministic [6] or stochastic [20], with constant [21] or variable [19] population size, and with standard [22] or bilinear [6, 20] incidence rate. Motivated by the above studies, in this paper, we consider the discretetime SIS epidemic model presented in [7], but with a different recruitment rate and without diseasecaused death. Indeed, a model with the standard incidence in which a (perfect but temporary) vaccination program has been included. For this model using a similar approach used in [7], we investigate the local stability and bifurcations theoretically and numerically. Moreover, the global stability of the diseasefree equilibrium of the model is proven. Furthermore, a generalization of the model is given by forwarding Euler discretization, and then its stability is studied. The effect of the vaccination in the model for controlling and eliminating the disease will also be shown.
The organization of the paper reads as follows: In the next section, the model is introduced, and equilibria of the model and its basic reproduction number are obtained. Two next subsections are devoted to studying the stability of the equilibria and bifurcations of the model, respectively. Then by using the forward Euler method, a discretetime model is obtained from a continuous version of the model, and the stability of its equilibria is analyzed. After a numerical discussion, finally we summarize the results.
Results
The model
Suppose that the individuals in a population are partitioned into susceptible individuals, infected individuals, and vaccinated individuals. Also, consider \(\Delta t\) as the appropriate time increment such that the changes in the model may take place at times \(0, \Delta t, 2\Delta t, 3\Delta t, \dots\). The number of total individuals at time \(t=n\Delta t\), for some n, is denoted by \(N_t\) and numbers of individuals in other compartments in the same time are as \(S_t\), \(I_t\), and \(V_t\).
All possible changes in the model and transmissions between its subpopulations together with their transmission rates have been shown in Fig. 1. Here, all parameters are assumed to be nonnegative, and N and \(\mu\) are positive. Also, \(\mu\) is the natural death rate, \(\beta\) is the contact rate, \(\gamma\) is the cure rate, \(\epsilon\) is the rate of losing immunity, while q and p are the vaccination rate in newcomers and susceptible individuals, respectively. The model can be illustrated by the following system of difference equations:
The susceptible individuals become infected at standard incidence rate \(\beta S_t I_t/N_t\). Moreover, summing equations in system (1), we see that \(N_{t+1}=N_t\), and then the population size will remain a constant value. Thus, by letting \(V_t=NS_tI_t\), the corresponding difference equation is deleted and the following system of two difference equations is obtained:
System (2) is considered under the following conditions, which are sufficient but not necessary for the nonnegativity of solutions.
These conditions are extracted from system (2) and they are natural requirements for the model. The first, state the rate at which the susceptible individuals who die or get infected or become vaccinated is less than one within a unit time. The second, says the rate at which the infected people who die or get recovered is less than one within a unit time.
The equilibria of the model are solutions of the following system:
From the first equation, we must have either \(\bar{I}=0\) or \(\beta \bar{S} /N  (\mu + \gamma )=0\). When \(\bar{I}=0\), the equilibrium is named the diseasefree equilibrium and is written as
while if \(\beta \bar{S} /N  (\mu + \gamma )=0\), we obtain
This equilibrium in which \(\bar{I}\ne 0\), is called the endemic equilibrium and is written as
Notice that \(I^*>0\) if and only if \([(1q)\mu +\epsilon ]\beta (\mu +p+\epsilon )(\mu +\gamma )>0\) if and only if
The quantity \(\mathcal {R}_0\) is referred to as the basic reproduction number of model (2) and is interpreted as the number of individuals who become infected by entering one infected individual into a fully susceptible population; see [23]. Thus it is reasonable intuitively, as we showed mathematically, that when \(\mathcal {R}_0<1\), each infected individual transmits the infection to less than one other individual and it is expected that the infection vanishes. While in the case \(\mathcal {R}_0>1\), each infected individual transmits the infection to more than one other individual and so the infection does not wipe out. Here, the assumption of \(\mathcal {R}_0\) was made directly from the expression of \(I^*\) and condition of its positivity in endemic equilibrium. Also, there are some methods to calculate \(\mathcal {R}_0\) for discrete models; for example, see [24]. We see that \(\mathcal {R}_0\) is independent of the total population size N. Also
and
Therefore, we can state the following lemma about the existence of equilibria of the model.
Lemma 1
For SIVS epidemic model (2), the diseasefree equilibrium \(Q^0\) always exists and the endemic equilibrium \(Q^*\) also exists if \(\mathcal {R}_0>1\).
Stability of the equilibria
We study stability of the system at an equilibrium by considering eigenvalues of the corresponding Jacobian matrix at that equilibrium. When eigenvalues are less than one, the system is stable.
Theorem 2
The diseasefree equilibrium is stable if and only if \(\mathcal {R}_0<1\).
Proof
The Jacobian matrix of model (2) at (I, S) is
Therefore, the Jacobian matrix at \(Q^0\) is given by
The eigenvalues of \(J(Q^0)\) are \(\lambda _1=1(\mu +\gamma )+(\mu +\gamma )\mathcal {R}_0\) and \(\lambda _2=1(\mu +p+\epsilon )\). obviously, \(\vert \lambda _2\vert <1\) by assumptions (3) and \(\vert \lambda _1\vert <1\) if and only if \(\mathcal {R}_0<1\). \(\square\)
Theorem 3
When \(\mathcal {R}_0>1\), the endemic equilibrium \(Q^*\) is stable and otherwise is unstable.
Proof
At \(Q^*\), we have \(\beta S^*/N=(\mu +\gamma )\) and so
Thus we get
and by assuming
we can rewrite them as
The characteristic equation of \(J^*\) is of the form \(P(\lambda )=\lambda ^2tr(J^*)\lambda +det(J^*)\) and according to the Jury conditions, all eigenvalues of \(J^*\) are from module less than one if and only if (see [25])
First, \(1+det(J^*)<2\) holds if and only if \(b_1+b_2<0\). Besides, \(\beta I^*/N>(\mu +\epsilon )\beta I^*/N\) and so \((\mu +p+\epsilon )+\beta I^*/N>(\mu +\epsilon )\beta I^*/N\), that is, \(b_1>b_2\) and thus the condition \(1+det(J^*)<2\) holds.
Second, for \(tr(J^*)>0\) we must show that \(tr(J^*) <1+det(J^*)\), which holds since it is equivalent to \(b_2>0\). If \(tr(J^*)<0\), then we have to prove \(tr(J^*) <1+det(J^*)\), which holds if and only if \(42b_1+b_2>0\). Indeed
since \(\mu +p+\epsilon <1\) and \(\beta I^*/N<1\). Therefore, when \(\mathcal {R}_0>1\), the Jury conditions are satisfied and the proof is completed. Now, in the next theorem we prove the global stability of the diseasefree equilibrium. \(\square\)
Theorem 4
The diseasefree equilibrium of model (2) is globally asymptotically stable if \(\mathcal {R}_0\le 1\).
Proof
First we divide the equations in system (2) by total population size N and letting \(s_t=S_t/N\) and \(i_t=I_t/N\), we obtain the following equivalent system
The basic reproduction number of this system is \(\mathcal {R}_0^{eqv}=\frac{\beta }{\mu +\gamma }\). Now, considering the function \(W(i_t)=i_t\), we see that \(W>0\), and \(W=0\) if and only if \(i=0\), that is at diseasefree equilibrium of (6). On the other hand, we have
and thus when \(\mathcal {R}_0^{eqv}\le 1\), we have \(W(i_{t+1})W(i_t)\le 0\). Therefore W is a Lyapunov function [26] and the diseasefree equilibrium of system (6) and as a result \(Q^0\) is globally asymptotically stable. \(\square\)
Bifurcations of the model
In a discretetime system, bifurcations occur at the equilibria of the given system when there exist some eigenvalues of the Jacobian matrix with module one. Indeed, for an eigenvalue \(\lambda\), if \(\lambda =1\), then a transcritical bifurcation (or flip bifurcation) occurs and when \(\lambda =1\), a perioddoubling bifurcation occurs; see [2, 25]. While a Neimark–Sacker bifurcation, which is the same as the Hopf bifurcation in continuous systems [27], occurs if there is a pair of conjugate complex eigenvalues with module one, \(\vert \lambda \vert =1\).
As we saw, the eigenvalues of \(J(Q^0)\) are \(\lambda _1=1(\mu +\gamma )+(\mu +\gamma )\mathcal {R}_0\) and \(\lambda _2=1(\mu +p+\epsilon )\). Also, \(\lambda _1=1\) if and only if \(\mathcal {R}_0=1\) and thus a transcritical bifurcation occurs at \(Q^0\) when \(\mathcal {R}_0=1\). On the other hand, \(\lambda _1=1\) if and only if \(\mathcal {R}_0=1\frac{2}{\mu +\gamma }\), but this is impossible because \(\mu +\gamma <1\) and \(\mathcal {R}_0\) becomes a negative value. This shows that a perioddoubling bifurcation does not occur at \(Q^0\). Also, the eigenvalues of \(J(Q^0)\) are both real, and therefore a NeimarkSacker bifurcation does not take place, too. Thus we have the following theorem.
Theorem 5
At diseasefree equilibrium, \(Q^0\) of SIVS epidemic model (2), transcritical bifurcation happens if \(\mathcal {R}_0=1\) while a perioddoubling bifurcation and a NeimarkSacker bifurcation do not take place.
Now, we consider bifurcations at the endemic state. The following theorem is devoted to this purpose.
Theorem 6
None of transcritical, perioddoubling, and NeimarkSacker bifurcations occur at endemic equilibrium \(Q^*\) for SIVS model (2).
Proof
We have \(\lambda =1\) is an eigenvalue of Jacobian matrix \(J(Q^*)\) if it is a root of the corresponding characteristic equation, \(1tr(J^*)+det(J^*)=0\). This holds if and only if \(b_2=0\), if and only if \(\beta I^*/N=0\), if and only if \(\mathcal {R}_0=1\), since
But by Lemma 1, \(Q^*\) exists when \(\mathcal {R}_0>1\). However, \(\lambda =1\) is an eigenvalue of \(J(Q^*)\) if \(P(1)=0\). This is satisfied if and only if \(42b_1+b_2=0\) that can be written as
or equivalently
Now, notice that as we concluded previously, \(P(1)>0\) when \(\mathcal {R}_0>1\). Also, \(\mathcal {R}_0=1\) implies \(\beta I^*/N=0\) and this results in \(2(\mu +p+\epsilon )=0\) which is impossible. These discussions state that a perioddoubling bifurcation does not happen at \(Q^*\).
If we write the characteristic equation of \(J^*\) as \(P(\lambda )=\lambda ^2+a_1\lambda +a_2\), we see that
Hence, the roots of \(P(\lambda )\) are both real and thus a NeimarkSacker bifurcation cannot appear at \(Q^*\). \(\square\)
Remark 1
If we omit the restriction \(\beta <1\) from the system and allow \(\beta\) to take values greater than or equal to one, then from (7), we get
Thus a perioddoubling bifurcation occurs at \(Q^*\) for \(\beta \ge 1\) if
The model obtained by the forward Euler discretization
The model described in Fig. 1 can be stated as a continuoustime model by the following system of ordinary differential equations (see [28, 29]):
We see that \(\dot{N}=dN/dt=0\) and therefore the population size is constant. Similar to the discretetime model, we get the following twodimensional system by substituting \(V=NSI\) and omitting variable V from the system:
Now in this section, we discretize and analyze model (9) by using the forward Euler method. Substituting \(\dot{S}=(S_{t+1}S_t)/\Delta\) and \(\dot{I}=(I_{t+1}I_t)/\Delta\), where \(\Delta\) is the fixed step size of the discretization, we obtain the discrete version of model as follows:
It can be seen that the equilibria of this model and the corresponding basic reproduction number are similar to model (2). The diseasefree equilibrium of discretized model, \(Q_d^0\), always exists while, its endemic equilibrium, \(Q_d^*\), exists only when \(\mathcal {R}_0>1\).
Theorem 7
When \(\mathcal {R}_0<1\), the diseasefree equilibrium of model (10) is stable if
Proof
The Jacobian matrix of the model at (I, S) is given by
and at the diseasefree equilibrium, it is
Thus the eigenvalues of \(J(Q_d^0)\) are \(\lambda _1=1+\Delta (\mu +\gamma )(\mathcal {R}_01)\) and \(\lambda _2=1\Delta (\mu +p+\epsilon )\). Therefore, \(\vert \lambda _1\vert <1\) if and only if \(\Delta <\frac{2}{(\mu +\gamma )(1\mathcal {R}_0)}\), and \(\vert \lambda _2\vert <1\) if and only if \(\Delta <\frac{2}{(\mu +p+\epsilon )}\). \(\square\)
Theorem 8
When \(\mathcal {R}_0>1\), the endemic equilibrium of model (10) is stable if \(\Delta <\Delta ^*\), where \(\Delta ^*\) is the smallest root of \(b_2x^22b_1x+4\), in which \(b_1=(\mu +p+\epsilon )+\beta I^*/N\) and \(b_2=(\mu +\epsilon )\beta I^*/N\).
Proof
The Jacobian matrix at endemic equilibrium is
According to Jury conditions (Schur–Cohn criterion), the matrix \(J_d^*\) is stable (i.e., the roots of its characteristic equation \(P_d(\lambda )=\lambda ^2+a_1\lambda +a_2\) lie inside the unit disk) if and only if the following conditions hold (see [30]):

(i)
\(1  a_2 > 0\),

(ii)
\(P(1) = 1 + a_1 + a_2 > 0\),

(iii)
\(P(1) = 1  a_1 + a_2 > 0\).
Here, \(a_1=tr(J_d^*)\) and \(a_2=det(J_d^*)\). We see, \(tr(J_d^*)=2\Delta b_1\) and \(det(J_d^*)=1\Delta b_1+\Delta ^2 b_2\). Thus condition (i) holds if and only if \(\Delta b_1\Delta ^2 b_2>0\), or equivalently \(\Delta <\frac{b_1}{b_2}\). Condition (ii), \(P(1)=1tr(J_d^*)+det(J_d^*)>0\), holds if and only if \(b_2 \Delta ^2>0\), which holds because \(b_2=(\mu +\epsilon )\beta I^*/N>0\). Condition (iii), \(P(1)=1+tr(J_d^*)+det(J_d^*)>0\), holds if and only if \(b_2 \Delta ^22b_1\Delta +4>0\). Since
thus \(b_2 \Delta ^22b_1\Delta +4>0\) has two roots of the form \(\frac{b_1\pm \sqrt{b_1^24b_2}}{b_2}\). Now, if we denote two roots as \(r_1\) and \(r_2\) (suppose \(r_1<r_2\)), then \(P(1)\) is positive when \(\Delta <r_1\) or \(\Delta >r_2\), since \(b_2>0\). Moreover, we can easily see that \(b_1\sqrt{b_1^24b_2}>0\), and thus \(r_1>0\). Therefore we can state, conditions (i)–(iii) hold if \(\Delta <r_1\), because we also have \(r_1<\frac{b_1}{b_2}\). \(\square\)
Remark 2
Model (2) was formulated straightly by considering a population and its transmissions. Also, the model can be concluded from discretized model (10) for \(\Delta =1\) and assumptions (3).
Discussion
In this section, we consider numerically theoretical results obtained in the paper. For this purpose assume that the parameters in the model are as \(q=0.4\), \(p=0.2\), \(\gamma =0.15\), \(\mu =0.1\), and \(\epsilon =0.25\). Moreover, consider units of time and population as one day and one million individuals, respectively. Let the number of initial individuals in each subpopulations be as \(I_0=0.4\), \(S_0=0.8\), and \(V_0=0.5\). We take the contact rate \(\beta\) as the bifurcation parameter and get the bifurcation diagram as it is shown in Fig. 2. We see that at \(\beta =0.443\), the dynamic of the system changes: The diseasefree equilibrium that was stable for values \(\beta <0.443\) becomes unstable and instead the endemic equilibrium becomes stable. Indeed, at \(\beta =0.4435\), we have \(\mathcal {R}_0=1\) and a transcritical bifurcation occurs. Moreover, it is seen that at \(\beta =2.428\), the endemic equilibrium becomes unstable and a perioddoubling bifurcation happens and after that, the system remains unstable. This value for \(\beta\) also is obtained according to Remark 1 as \(\beta =2.4279\). Figure 3 shows the Lyapunov exponents of the Jacobian matrix for the same values of \(\beta\). Here also it is observable that for values \(\beta = 0.443, 2.428\), and 3.235, the Lyapunov exponent is positive as seen in the bifurcation diagram. Figure 4 presents solutions of the system for various values of \(\beta\) and the behavior of the solutions is the same as we expect from the bifurcation diagram and the Lyapunov exponents. For \(\beta =0.4\), we have \(\mathcal {R}_0=0.9018<1\) and as we expect from Theorem 2, the disease will vanish. While, for values \(\beta =0.55\), \(\beta =2.45\), and \(\beta =2.5\), we have \(\mathcal {R}_0=1.2400\), \(\mathcal {R}_0=5.5236\), and \(\mathcal {R}_0=5.6364\), respectively, that all are greater than one and according to Theorem 3, the infection remains at a positive level. In addition, Fig. 5 displays some parts of solutions of infected population \(I_t\) for different values of \(\beta\). It is observable that the behavior of solutions corresponds to those are in the bifurcation diagram. The effect of the discretization of the continuoustime model (9) by applying the forward Euler method has been considered in Fig. 6. The bifurcation diagram shows the dynamics of the infected population in the discretized model (10) when the step size \(\Delta\) varies. The contact rate has been supposed as \(\beta =2.7\) and other parameters are the same as preceding simulations. As it was established in Theorem 8, the endemic equilibrium \(Q_d^*\) is stable for \(\Delta <\Delta ^*=0.8946\) while for greater values of \(\Delta\), the endemic equilibrium becomes unstable.
Conclusions
In this paper, we introduced and studied an SIS epidemic model that includes a vaccination program. The equilibria of the model were detected: The diseasefree equilibrium \(Q^0\) in which the infection will be extinct, and the endemic equilibrium \(Q^*\) in which the disease will persist in the population. It was proved that under some assumptions on parameters for positivity of solutions, \(Q^0\) and \(Q^*\) are stable if \(\mathcal {R}_0<1\) and \(\mathcal {R}_0>1\), respectively. It was proven \(Q^0\) is also globally asymptotically stable when \(\mathcal {R}_0\le 1\). Thus the basic reproduction number \(\mathcal {R}_0\) plays important role in determining the dynamics of the model. To clear the effect of vaccination, a sensitivity analysis is performed on \(\mathcal {R}_0\) to determine how sensitive the model is to changes in the value of parameters that are related to vaccination. The normalized forward sensitivity index of \(\mathcal {R}_0\) concerning parameter x has been defined in [31] as \(\Psi _x^{R_0} = \frac{x}{R_0} \times \frac{\partial {R_0}}{\partial x}\). The sensitivity indices for vaccination proportions p and q are calculated as \(\Psi _q^{R_0} = \frac{  q\mu }{(1  q)\mu + \varepsilon } < 0\) and \(\Psi _p^{R_0} = \frac{  p}{\mu + p + \varepsilon } < 0\). These indices are both negative meanings \(\mathcal {R}_0\) is inversely related with p and q; an increase in parameters will cause a decrease in value of \(\mathcal {R}_0\), that is more vaccinated individuals, the infection becomes more controllable. On the other hands, The basic reproduction number for the corresponding model without vaccination (i.e. \(p=q=0\)) is \(\mathcal {R}_0=\frac{\beta }{\mu +\gamma }\). Thus \(\mathcal {R}_0=\Big (\frac{(1q)\mu +\varepsilon }{\mu +p+\varepsilon }\Big ){\tilde{\mathcal {R}_0}}\) and \(\mathcal {R}_0<\tilde{\mathcal {R}_0}\). When \(\tilde{\mathcal {R}_0}>1\), then the vaccination must be used such that \(\mathcal {R}_0<1\) to eliminate an infection.
Furthermore, the bifurcations of the model were investigated and it was proved that when \(\mathcal {R}_0=1\) system has a transcritical bifurcation and although the Neimark–Sacker bifurcation does not appear, it may have a perioddoubling bifurcation if we ignore the restriction \(\beta <1\). To study the discretization of the continuous version of the model, we applied the forward Euler method and analyzed the effect of step size of the discretization on the dynamics of the model. We established the sufficient condition for stability of diseasefree equilibrium \(Q_d^0\) and endemic equilibrium \(Q_d^*\) in the discretized model. Finally, we examine the results obtained in the paper in numerical example by considering the bifurcation diagram, the Lyapunov exponents of the Jacobian matrix, and graphs of solutions for values of \(\beta\) and \(\Delta\). It was observed that the numerical discussions verify the theoretical results.
Availability of data and materials
Not applicable.
Abbreviations
 SI(S):

Susceptibleinfected(susceptible)
 SIR(S):

Susceptibleinfectedremoved(susceptible)
 SIVS:

Susceptibleinfectedsusceptible epidemic model with vaccination
References
 1.
Brauer F, CastilloChavez C. Mathematical models in population biology and epidemiology, vol. 1. Berlin: Springer; 2001.
 2.
Allen LJ. Introduction to mathematical biology. Upper Saddle River: Pearson/Prentice Hall; 2007.
 3.
Allen LJ. Some discretetime SI, SIR, and SIS epidemic models. Math Biosci. 1994;124(1):83–105.
 4.
CastilloChavez C, Yakubu A. Discretetime SIS models with complex dynamics. Nonlinear Anal Theory Methods Appl. 2001;47(7):4753–62.
 5.
Brauer F, Feng Z, CastilloChavez C. Discrete epidemic models. Math Biosci Eng. 2010;7(1):1–15.
 6.
Farnoosh R, Parsamanesh M. Disease extinction and persistence in a discretetime SIS epidemic model with vaccination and varying population size. Filomat. 2017;31(15):4735–47.
 7.
Parsamanesh M, Mehrshad S. Stability of the equilibria in a discretetime SIVS epidemic model with standard incidence. Filomat. 2019;33(8):2393–408.
 8.
Xiang L, Zhang Y, Huang J. Stability analysis of a discrete SIRS epidemic model with vaccination. J Differ Equ Appl. 2020;26(3):309–27.
 9.
Roeger LW, Barnard RW. Preservation of local dynamics when applying central difference methods: application to SIR model. J Differ Equ Appl. 2007;13(4):333–40.
 10.
Liu J, Peng B, Zhang T. Effect of discretization on dynamical behavior of SEIR and SIR models with nonlinear incidence. Appl Math Lett. 2015;39:60–6.
 11.
Aranda DF, Trejos DY, Valverde JC. A discrete epidemic model for bovine Babesiosis disease and tick populations. Open Phys. 2017;15(1):360–9.
 12.
Mickens RE. Discretizations of nonlinear differential equations using explicit nonstandard methods. J Comput Appl Math. 1999;110(1):181–5.
 13.
Izzo G, Vecchio A. A discrete time version for models of population dynamics in the presence of an infection. J Comput Appl Math. 2007;210(1–2):210–21.
 14.
Hu Z, Teng Z, Jiang H. Stability analysis in a class of discrete SIRS epidemic models. Nonlinear Anal Real World Appl. 2012;13(5):2017–33.
 15.
Ma X, Zhou Y, Cao H. Global stability of the endemic equilibrium of a discrete SIR epidemic model. Adv Differ Equ. 2013;1:42.
 16.
Cui Q, Zhang Q. Global stability of a discrete SIR epidemic model with vaccination and treatment. J Differ Equ Appl. 2015;21(2):111–7.
 17.
Van den Driessche P, Yakubu A. Disease extinction versus persistence in discretetime epidemic models. Bull Math Biol. 2019;81(11):4412–46.
 18.
KribsZaleta CM, VelascoHernández JX. A simple vaccination model with multiple endemic states. Math Biosci. 2000;164(2):183–201.
 19.
Yang W, Sun C, Arino J. Global analysis for a general epidemiological model with vaccination and varying population. J Math Anal Appl. 2010;372(1):208–23.
 20.
Farnoosh R, Parsamanesh M. Stochastic differential equation systems for an SIS epidemic model with vaccination and immigration. Commun Stat Theory Methods. 2017;46(17):8723–36.
 21.
Parsamanesh M. Global stability analysis of a VEISV model for network worm attack. Univ Politeh Buchar Sci Bull Ser A Appl Math Phys. 2017;79(4):179–88.
 22.
Li J, Ma Z. Qualitative analyses of SIS epidemic model with vaccination and varying total population size. Math Comput Modell. 2002;35(11):1235–43.
 23.
Hethcote HW. The mathematics of infectious diseases. SIAM Rev. 2000;42(4):599–653.
 24.
Allen LJ, Van den Driessche P. The basic reproduction number in some discretetime epidemic models. J Differ Equ Appl. 2008;14(10–11):1127–47.
 25.
Elaydi S. An introduction to difference equations. Berlin: Springer; 2005.
 26.
Stuart A, Humphries AR. Dynamical systems and numerical analysis, vol. 2. Cambridge: Cambridge University Press; 1998.
 27.
Kuznetsov YA. Elements of applied bifurcation theory, vol. 112. Cambridge: Springer; 2013.
 28.
Parsamanesh M, Erfanian M. Global dynamics of an epidemic model with standard incidence rate and vaccination strategy. Chaos Solitons Fractals. 2018;117(1):192–9.
 29.
Parsamanesh M, Farnoosh R. On the global stability of the endemic state in an epidemic model with vaccination. Math Sci. 2018;12(14):313–20.
 30.
Jury EI. Theory and application of the ztransform method. Hoboken: Wiley; 1964.
 31.
Chitnis N, Hyman JM, Cushing JM. Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model. Bull Math Biol. 2008;70(5):1272–96.
Acknowledgements
We thank the reviewers’ valuable comments for improving the quality of this work.
Funding
No funding was obtained for this study.
Author information
Affiliations
Contributions
MP and ME introduced the model. MP, ME, and SM obtained and proved the mathematical results. MP and ME performed the numerical analysis. MP, ME, and SM drafted the manuscript. All the authors have read and approved the final manuscript.
Corresponding author
Ethics declarations
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare that they have no competing interests.
Additional information
Publisher's Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.
About this article
Cite this article
Parsamanesh, M., Erfanian, M. & Mehrshad, S. Stability and bifurcations in a discretetime epidemic model with vaccination and vital dynamics. BMC Bioinformatics 21, 525 (2020). https://doi.org/10.1186/s12859020038391
Received:
Accepted:
Published:
Keywords
 SIS epidemic model
 Discretetime system
 Stability
 Lyapunov exponent
 Bifurcation