Global stability analysis for a generalized delayed SIR model with vaccination and treatment

In this work, we investigate the stability of an SIR epidemic model with a generalized nonlinear incidence rate and distributed delay. The model also includes vaccination term and general treatment function, which are the two principal control measurements to reduce the disease burden. Using the Lyapunov functions, we show that the disease-free equilibrium state is globally asymptotically stable if \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathcal{R}_{0} \leq 1 $\end{document}R0≤1, where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathcal{R}_{0} $\end{document}R0 is the basic reproduction number. On the other hand, the disease-endemic equilibrium is globally asymptotically stable when \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathcal{R}_{0} > 1 $\end{document}R0>1. For a specific type of treatment and incidence functions, our analysis shows the success of the vaccination strategy, as well as the treatment depends on the initial size of the susceptible population. Moreover, we discuss, numerically, the behavior of the basic reproduction number with respect to vaccination and treatment parameters.


Introduction
Mathematical modeling has become a powerful and important tool to understand infectious disease dynamic behavior and to improve control of the disease in a population. These models are often described by many forms such as: SI, SIS, SIR, or SIRS models, where S stands for susceptible subpopulation, I is infected subpopulation, and R is recovered subpopulation. The progress of a disease in a population is dictated by the nature and the mode of transmission between infected and susceptible individuals. The mode of transmission is the method of transfer by which the infection moves or carries from one place to another to reach the new host (for example airborne, saliva, vector-borne, and bodily fluids). Hence, it is natural to adapt these models to the concerned disease by choosing the right incidence function. It is known that the function forms of the incidence rate of the infection have a crucial role in the modeling of the infection dynamics, many forms of incidence function have been considered by the researchers in mathematical epidemiology, for example, the bilinear incidence rate βSI, where β is the transmission rate of infection, the saturated incidence rate βSI 1+αI , with α defined as the inhibitory coefficient, and many other forms (see [1][2][3][4][5][6][7]). To make a model more realistic, the introduction of the time delay is more interesting, and considerable attention has been paid by several authors to studying the dynamics of epidemic models with discrete or distributed time delay (see [3,4,[8][9][10][11]). Vaccination and treatment are the two main public health control strategies that help to minimize the burden of an infectious disease spread and to delay a possible outbreak. Vaccination has the role of preventing healthy people from getting infected by a disease, while treatment cures a disease and can also be used as a prophylactic. These control strategies are usually used together to contain the disease spread (see [12] in the context of influenza). Tulu et al., in [13], developed a mathematical model to study the effect of both vaccination and quarantine on the spread of Ebola virus, they applied the vaccination strategy to the susceptible individuals. However, in [14], the authors studied the global dynamics of an SEIRS epidemic model with preventive vaccination applied to the newborns. Various vaccination policies were studied in different mathematical models (see [8,[15][16][17][18]). It is well known in classical epidemic models that the recovery rate due to treatment is proportional to the number of the infected individuals. However, this proportionality is not satisfied in the reality because of limited medical facilities (see [19]). In order to include the limited capacity of medical resources, Chauhan et al.,in [20], considered the piecewise linear treatment function of the form where I 0 is the capacity of treatment. Recently, Li introduced the following saturated treatment function [21]: where a represents the maximal medical resources supplied per united time and is half-saturation constant, which measures the effect of being delayed for treatment. Other works have investigated the effects of the treatment on an epidemic (see [19][20][21][22][23][24][25]) and also its optimal control (see [26]). The motivation of this work comes from [10,11], where the authors studied an SIR epidemic model with nonlinear incidence function, and from [19][20][21], where the authors considered a special type of treatment function. The present work would be a continuation and generalization of the above cited works. It is concerned with a generalized SIR epidemic model with distributed delay, vaccination, and treatment. This model incorporates distributed delay, general incidence function, vaccination, and general function treatment. In fact, we apply the vaccination to both susceptible and newborn individuals. On the newborn individuals, we apply the mechanism of "all-or-nothing" vaccine. Recall that an "all-or-nothing" vaccine offers complete protection to a subset of the vaccinated individuals, but the remainder of them stays susceptible to catching the disease. Second, we consider a class of treatment functions satisfying suitable conditions, and it is more general than the one given by (1) or (2). Moreover, it is necessary to point out that the delay in this model represents the incubation time taken to become infectious. This model can be applied to investigate the impact of the vaccination and the treatment in containing the spread of infections which have an incubation time to become infectious, for example, SARS-CoV(see [27,28]). Our purpose in this work is to investigate the impact of the combined vaccination and treatment strategies on the dynamic behavior of the considered model. We prove that the basic reproduction number R 0 depends explicitly on the vaccination parameters and the general treatment function T(I). Moreover, we discuss the global stability of the model near equilibria (the disease-free equilibrium E 0 and the disease-endemic equilibrium E * ) by means of R 0 and Lyapunov's method. Furthermore, to verify the theoretical results, numerical simulations are performed for special treatment and incidence functions. For illustration, we give some numerical results on the behavior of the basic reproduction number R 0 with respect to vaccination and treatment parameters.
The paper is organized as follows. We give a mathematical model formulation in Sect. 2. In Sect. 3, we propose a mathematical analysis of the considered model. More precisely, we calculate the basic reproduction number R 0 , and we determine the disease-free equilibrium E 0 and the endemic equilibrium E * . Moreover, we prove the local stability of the disease-free equilibrium and the global stability of E 0 and E * . In Sect. 4, we give some numerical examples with an incidence and treatment functions satisfying assumptions presented in the previous sections. We finish the paper, in the last section, by providing some concluding remarks.

Mathematical model and preliminaries
In this work, we are interested in a general SIR epidemic model with distributed delay, vaccination, and treatment. The dynamics are governed by the diagram in Fig. 1.
From Fig. 1, we have the following SIR model: where S(t), I(t), and R(t) denote the numbers of susceptible, infective, and recovered individuals at time t respectively. The susceptibles are augmented by the birth of newborns.
Here, we assume that the birth rate b and death rate μ are not the same. The parameter p Figure 1 Flow diagram of the disease transmission is the fraction of the vaccinated newborns. A fraction ∈ [0.1) (the all-or-nothing parameter) of the vaccinated newborns exhibits an unsuccessful vaccination and passes directly to the susceptible class. Our vaccine has an efficacy of 1 -(see [27][28][29][30][31][32][33]). For simplicity, we assume that the recovered class stands also for the vaccinated state. Hence, susceptible individuals get vaccinated with rate d.
The nonlinear incidence rate and distributed delay are considered to represent wide class epidemic model similarly as in [10,11]. More precisely, by taking β the disease transmission coefficient, individuals leave the susceptible class at a rate where h represents the maximum time taken to become infectious. The function g that satisfies h 0 g(τ ) dτ = 1 is assumed to be nonnegative. The function f : R 2 + → R 2 + is assumed to be continuously differentiable in the interior of We also assume that the disease causes death with rate c and γ is the natural recovery rate of the infected individuals.
The function T : R + → R + represents the treatment function which we assume to be continuously differentiable and concave down satisfying the following hypotheses: The treatment rate T(I) I is monotone increasing. The assumption of the concavity of the treatment function refers to the fact that the supply of the treatment drugs increases as the disease kicks off in the population until it reaches a maximum level, then the treatment drug stocks start going down due to the exhaustive consumption.
Hypothesis (T 1 ) means that there is no treatment if there is no infection, while hypothesis (T 2 ) reflects the increasing effort needed from the public health authorities to provide treatment during the time of the infections.

Proposition 2.1 The solution of (3), with initial condition
Proof We prove, by contradiction, that the solution (S(t), I(t), R(t)) is positive. Let t 1 = min{t ≥ 0 : S(t)I(t) = 0}, and we assume that S(t 1 ) = 0, which implies that, for all 0 ≤ t ≤ t 1 , It follows that This contradicts S(t 1 ) = 0. Using a similar argument, we can prove that S(t) > 0 and I(t) > 0 for all t ≥ 0. The positivity of R follows from the inequality which implies that For the boundedness, we note that It follows that lim t→+∞ n(t) = b μ , which completes the proof.
The local existence and boundedness of the solution of (3) imply the global existence of the solution.
As the variable R does not appear in the first two equations for system (3), we focus our analysis on the reduced system 3 Analysis of the model

Existence of equilibria points
System (5) has a disease-free equilibrium On the other hand, using the next generation method [36], the basic reproduction number should be as follows.

Lemma 3.1 The basic reproduction number is
Note that S 0 depends on the vaccination of susceptible population and the treatment terms.
Proof Let X = (I, S) T , then it follows from system (5) that The Jacobian of matrices F and ν at the disease-free equilibrium E 0 is given by where f 2 (E 0 ) is the derivative of f with respect to I at E 0 . The inverse of V is given by Thus, the next generation matrix for system (5) is Since R 0 is the spectral radius of the matrix FV -1 , it follows that the basic reproduction number is To prove the existence of an endemic equilibrium, we need the following lemma. Proof Let K be the function defined on R + by We have Since K is continuous, the equation K(u) = 0 has a unique positive solution in the interval (0, b a ).
Next result shows the existence of the endemic equilibrium.
Let K be the function defined for R + \ {0} to R by By hypotheses (H 2 ) and (T 2 ), K is strictly monotone decreasing on R + \ {0} satisfying Moreover, by Lemma 3.2, there exists a unique solution I 0 of the following equation: and then Hence, there exists a unique positive real I * such that 0 < I * < I 0 and K I * = 0, which allows us to conclude that E * = (S * , I * ) is the unique endemic equilibrium of system (5).

Local stability analysis
In this section, we discuss the local stability of the disease-free equilibrium of system (5).
We have the following result.
Then the characteristic equation of system (8) at E 0 is of the form It is clear that λ = -(μ + d) is a root of (9). All other roots λ of (9) are determined by the following equation: Then by separating real ( ) and imaginary ( ) parts, we derive Using the first equation of the above system, we obtain We suppose, by contradiction, that there exists λ ∈ C such that (λ) ≥ 0, and it satisfies equality (10). Then Since the function T is concave down, it follows that T (0) ≥ 0. Moreover, we know that f 2 (E 0 ) > 0, which implies If R 0 < 1, then βf 2 (E 0 ) < μ + c + γ + T (0) and which gives a contradiction with inequality (11). Then the real parts of all the eigenvalues of (9) are negative. Therefore, if R 0 < 1, the disease-free equilibrium E 0 of system (5) is locally asymptotically stable. Now, let From the fact that P(0) = (μ+c+γ +T (0))(1-R 0 ) < 0 if R 0 > 1 and lim λ− →+∞ P(λ) = +∞, we conclude that there is at least one positive root of (9). Hence, if R 0 > 1, E 0 is unstable.

Global stability of the disease-free equilibrium
The next result gives the condition of the global asymptotic stability of the disease-free equilibrium E 0 of system (5).
Proof To prove this result, we consider the following Lyapunov function: where σ = μ + c + γ , and From hypothesis (T 2 ), it follows that Hypothesis (H 1 ) implies that and hypothesis (H 2 ) gives that Hence, Then the condition R 0 ≤ 1 implies that d dt V (t) = 0}. By the Lyapunov-LaSalle principle, we conclude that the disease-free equilibrium E 0 of system (5) is globally asymptotically stable.

Global stability of the endemic equilibrium
In this section, we aim to show the global asymptotic stability of the endemic equilibrium E * of system (5) via a Lyapunov stability approach.
Proof Let G be the function defined from R + to R by It is clear that G(x) ≥ 0 if x > 0 and G(x) = 0 if x = 1. Let us consider the following Lyapunov function: Moreover, we have Moreover, hypothesis (H 1 ) implies that and hypothesis (T 2 ) gives Hence, d dt U(t) ≤ 0. We conclude that the endemic equilibrium of system (5) is globally asymptotically stable.

Numerical results
In this section, we present the numerical simulation of the model by considering the following delayed SIR epidemic model with vaccination, treatment, and distributed time delay: The function g is chosen, as in [37], in the following form: On the other hand, the treatment function T, similarly to [23], is defined by The reproduction number R 0 is given by For our system (12) without vaccination and treatment, the reproduction number is given by Hence R 0 can be rewritten as If R 0 ≤ 1, then the disease will die out (the disease-free equilibrium E 0 is globally asymptotically stable) without any control measures.
However, if R 0 > 1, then where S 0 is given in (6). Similarly, This shows that, during the epidemic R 0 > 1, if the number of susceptible population is below the thresholdS, then the disease can be controlled by vaccination and treatment. However, if the number of susceptible population is above the thresholdS, then the disease will persist in the population.
To make sense of our simulation, we will focus on the case of R 0 > 1, and we choose the parameters p and a to guarantee the clearance of the disease from the population by these two public health control measures.
All the numerical simulations are performed using the explicit Runge-Kutta-like method (dde45) [38].
First, we start with the case of no vaccination and no treatment (p = 0 and a = 0). In this situation our model is similar to that of Enatsu et al. [10], in which the authors claim that when the basic reproduction number, denoted by R 0 , is greater than one (R 0 > 1), the disease persists. However, our numerical simulation ( Fig. 3 and Fig. 4) shows that the disease will die out even if R 0 > 1.
Next, we consider the case with vaccination and no treatment, with R 0 > 1 and R 0 < 1. As shown in Fig. 3, the disease dies out, which corresponds to our theoretical result.
Finally, we give the simulation of the case with vaccination and treatment, with R 0 > 1 and R 0 < 1. As shown in Fig. 4, the treatment with vaccination helps the eradication of the infection from the population.
For more illustration, it is very interesting to discuss the behavior of the basic reproduction number R 0 with respect to vaccination and treatment parameters. Namely, the  In Fig. 5 we show the effect of vaccination and treatment parameters on the dynamic of R 0 . We notice that the critical valuesp = 0.875,d = 1.02, andā = 1.157 represent separated values between the endemic state and the disease-free state for (j), (k), and (l) respectively (it means the cases R 0 < 1 and R 0 > 1).

Conclusion
In this work, we analyzed a delayed SIR model with generalized incidence function and distributed delay as the contact between infected individuals and healthy ones does not result in an immediate infection. The delay, presented in this work, reflects the time that it takes to have an infection after the contact. The model also included the two main types of disease control measures: vaccine and treatment. The question that arises in using these two measures is how each vaccination should depend on the treatment. In fact, as the treatment is the first control measure to be taken either as a prophylactic or antiviral, the vaccination implementation should take into consideration the effect of the treatment on the disease infectiousness. Moreover, the function treatment was chosen to reflect the reality of drug stock supply during the time of the infections. Our analysis showed that when R 0 ≤ 1, the disease-free equilibrium is globally asymptotically stable, and when R 0 > 1, then there is a unique disease-endemic equilibrium, which is globally asymptotically stable. To put this result in context, we chose the treatment function T(I) = aI(t) 1+ξ I(t) (see [23]).
In our analysis we showed that when the disease is endemic, in the absence of the vaccination and treatment, then there are two possible scenarios: (a) if the number of susceptible population is below the thresholdS, then the disease can be controlled by vaccination and treatment; (b) if the susceptible population is above the thresholdS, then the disease will persist in the population. This finding reflects the limited capability of the control measure to eradicate the disease if the population is too large.