SIRS Epidemic Models with Delays, Partial and Temporary Immunity and Vaccination

: The basic reproduction, or reproductive number, is a useful index that indicates whether or not there will be an epidemic. However, it is also very important to determine whether an epidemic will eventually decrease and disappear or persist as an endemic. Different infectious diseases have different behaviors and mathematical models used to simulated them should capture the most important processes; however, the models also involve simplifications. Influenza epidemics are usually short-lived and can be modeled with ordinary differential equations without considering demographics. Delays such as the infection time can change the behavior of the solutions. The same is true if there is permanent or temporary immunity, or complete or partial immunity. Vaccination, isolation and the use of antivirals can also change the outcome. In this paper, we introduce several new models and use them to find the effects of all the above factors paying special attention to whether the model can represent an infectious process that eventually disappears. We determine the equilibrium solutions and establish the stability of the disease-free equilibrium using various methods. We also show that many models of influenza or other epidemics with a short duration do not have solutions with a disappearing epidemic. The main objective of the paper is to introduce different ways of modeling immunity in epidemic models. Several scenarios with different immunities are studied since a person may not be re-infected because he/she has total or partial immunity or because there were no close contacts. We show that some relatively small changes, such as in the vaccination rate, can significantly change the dynamics; for example, the existence and number of the disease-free equilibria. We also illustrate that while introducing delays makes the models more realistic, the dynamics have the same qualitative behavior.


Introduction
Influenza, or flu, is a contagious disease spread by contact with an infected individual.Seasonal influenza affects between 3% and 11% of the population in the United States every year [1].The flu peaks from December to February of each year.Although some patients die from seasonal flu, it is not a particularly deadly disease with 0.2 deaths per 100,000 population [2].Other contagious diseases that are short-lived are the common cold, Respiratory Syncytial Virus Disease and Streptococcus throat.There are many papers dealing with mathematical models for epidemics, both via direct contact and vector transmitted.Mathematical modeling of epidemics started in the early 20th century.One of the first efforts is [3] and the seminal paper of Kermack and McKendrick [4] is widely recognized as introducing the first compartmental epidemic model based on ordinary differential equations (ODEs).The use of ODEs is a good modeling assumption if the different populations are well-mixed and there is little or no spatial dependence.Kermack and McKendrick introduced the SIR (susceptible-infected-removed) model, which assumes that the population can be in one of three states: susceptible, infected or infective, and removed (recovered with permanent immunity or deceased).The SIR model and its many variations including the SIS (susceptible-infective-susceptible) and SIRS (susceptible-infective-recovered-susceptible) models have been very influential in understanding epidemics and their control.Among the modifications made to the original SIR model of Kermack and McKendrick, are the introduction of demographics (not important for seasonal flu), vaccinations, different age groups, use of anti-virals and isolation.One important modification is the addition of delays.Delays appear in epidemics due to the time it takes certain processes to occur, like the time it takes a susceptible to become infected after a contact with and infected individual, or the time it takes vaccines to become effective after being administered.The two main ways of introducing delays into a model based on ordinary differential equations are by using delay differential equations (DDEs) or by introducing and additional exposed or latent population.DDEs may have different dynamics than ODEs [5,6].Among the possible differences are the introduction of oscillations, introduction of discontinuities in the time derivatives, non-uniqueness and different stability.While there is a wide range of literature about models of epidemics or disease transmission using delay differential equations, most books on DDEs have sparse coverage of contagious disease spread, and books on epidemics only have little material on DDEs.Two exceptions are the following books [7,8].SIR models consider that after being infected, individuals recover and achieve permanent and total immunity from the disease.For seasonal influenza and other contagious diseases the immunity is temporal which leads to SIRS models in which the recovered individuals lose immunity completely after a certain time period.However, there are other possibilities: the individuals whom recovered lose total immunity but keep partial immunity forever; the individuals whom recovered lose total immunity and achieve partial immunity for a time period; some individuals whom recovered keep total and permanent immunity while others only keep partial permanent immunity, and so on.What happens depends on the particular infection and population.Vaccines may also confer total permanent immunity, or total immunity that wanes to partial immunity, or other scenarios.
One of the most important uses of mathematical models of contagious diseases is to predict whether there will be an epidemic or not.The basic reproduction number R 0 defined as the number of secondary infections produced by a newly infected individual at the start of the infective process is a useful index to determine whether the epidemic will die or not.It is usually calculated using the next generation matrix method [9,10].It also gives the stability of the disease-free equilibrium (DFE) of the model.In case that the number of infectives grows initially, it is also very important to determine if the epidemic will eventually die or not.In this paper, we will develop mathematical models of short-lived contagious diseases liked seasonal influenza that include different assumptions about total, partial, permanent and temporal immunity.We will also add vaccinations with different immunity effects, and also add the use of anti-virals and isolation.The models will include delayed versions, using DDEs and also exposed populations.For these models we will study whether the epidemic persists or not.Ref. [11] compares some basic models for influenza that we will use as a starting point.The rest of the paper is organized as follows.
In Section 2, we present several new models, that are extensions to the SIRS model, with different hypotheses about total and partial immunity and vaccination, determine their equilibria, and analyze the stability of their disease-free equilibria.In Section 3, numerical simulations and results are presented, including the effects of introducing delays.We also show that some relatively small changes in the models can significantly change the dynamics.Finally, in Section 4, the discussions and conclusions about the models and results are presented.

Materials and Methods
In this section, we introduce several new mathematical models of short-term epidemics with different types of immunities.We also add vaccination, the use of anti-virals and isolation or quarantine.To determine if an epidemic develops according to the model used, the equilibrium points need to be found.If there is an equilibrium point with no infectives and susceptibles equal to the total population, and if the number of infectives grows with time, then there will be an epidemic.An equilibrium with zero infectives is known as a disease-free equilibrium (DFE).There are different methods to determine if the number of infectives grows: the differential equation for the infectives has a positive right-hand side; the DFE is locally unstable; using the basic reproduction number R 0 gives R 0 > 1; the infective compartments subsystem of equations is locally unstable; the R 0 index calculated using the next generation matrix method is greater than 1 [10,12].In case there is an epidemic, it is very important to establish if it will eventually die.This can be calculated by determining if there is a DFE for positive time that is locally asymptotically stable.This may happen for a finite time or only as time tends to infinity.This determination will be conducted in the following subsections for different models with no demographics.The first subsection is a review of the susceptible-infective-susceptible model (SIS), the susceptibleinfective-recovered model (SIR) and the susceptible-infective-recovered-susceptible model (SIRS) with different immunities.

SIS, SIR and SIRS Models
A very simple epidemic model, for a total constant population N, consisting of only two population compartments, susceptibles S and infectives I is the SIS model.It assumes that a susceptible gets infected with probability β/N > 0 by a contact with an infective, and then recovers with no immunity and returns to the susceptible compartment.The common cold is an example of a contagious disease that can be modeled by a SIS system of equations.The total population is usually considered constant which assumes that the disease is fast compared to the demographic effects.An almost equivalent assumption, since the death rate is very small compared to the time frame of the disease, is the number of births is equal to the number of deaths with no migration.The SIS model with no demographics is [13,14] where γ > 0 is the recovery rate and N > 0 is the total population.Since the population is constant, S(t) = N − I(t) and we can eliminate one equation The SIS model has two equilibrium states, I = 0 and I = (1 − γ β )N.The last one is physical only for γ β ≤ 1.The basic reproduction number (also called basic reproductive number) R 0 is the number of secondary infectives generated by a single infective at the beginning of the infection process.For the SIS model it can be calculated as the product of the probability of a new infective, β/N, times the number of susceptibles at the beginning of the infective process, N, times the average lifetime of an infection, 1/γ.Therefore R 0 = β/γ.If R 0 > 1 the number of infectives increases away from the DFE, and the other equilibrium makes biological sense since all the states are non-negative.If R 0 < 1 the infection disappears and the DFE equilibrium is the only one that makes biological sense.R 0 > 1 implies that the DFE is locally asymptotically unstable so at least one the eigenvalues of the linearized system about the DFE has positive real part.To determine if the infection is growing it is only necessary to look at the infected compartment, I(t).The next generation matrix method [9,10] is widely used to calculate the basic reproduction number in more complex models in epidemiology.Applying this method to the SIS model there is only one infective compartment, I, so the Jacobian J of this system is a 1 × 1 matrix.Decomposing this matrix as J = T + Σ, with T the transmission and Σ the transition matrices, produces T = βS/N and Σ = −γ.The next generation matrix is K = −TΣ −1 = βS/(γN).At the beginning of the infection S = N − 1, so the spectral radius of K, which is R 0 , can be approximated by R 0 = β/γ.The local stability of the DFE can also be established by finding an index that states when the eigenvalue of the linearized system consisting of the disease compartment evaluated at the DFE has negative real part.A straightforward way of finding whether the number of infectives grows at the start of the epidemic is to the evaluate the differential equation for I at the initial state with S = N: which is positive for β/γ > 1.Note that since at the start of the epidemic I = 1, S = N − 1.However, for large N, S can be taken equal to N. Another way of finding the stability of the DFE at the start of the infection can is by calculating the eigenvalues of the Jacobian matrix corresponding to the full system of equations, but the formulas are usually complicated or cannot be determined analytically.Using the characteristic polynomial and the Routh-Hurwitz criterion is another possibility [13,15].Two delays may be added to the SIS model, τ the time in takes the disease to develop in a susceptible individual after being infected by a contact with an infective individual, and the time it takes for a newly infected individual to recover.This second delay is usually included implicitly since 1/γ is the average recovery time.There are two ways to introduce the first delay.The first way is given by SIS Delay Model I: ( The assumption is that a contact at t − τ produces a new infective at time t.The second way is given by SIS Delay Model II: The assumption now is that an individual susceptible after contact with an infective turns into an infective but takes it a time τ to actually become infected with the disease and be able to infect.Both models make biological sense.Another approach to studying epidemic models is to consider them as chemical reactions, as shown, for example, in [16].As an illustration, the reaction between an infective and a susceptible produces two infectives.
For the reaction to happen both species should be at the same time, which is a caveat for Delay Model II.However, since both ways of introducing the infection delay are used in many papers, we include both to see if they produce the same dynamics qualitatively.The dynamics may be different since their characteristic equations are different.For other epidemic models that introduce the delay in the infection as was performed in Delay Model 1 see, for example, [17][18][19].And for other applications using the delay as in Delay Model 2 see [20,21].
For the SIS model, the equilibrium solutions of both delay models are the same as those of the ODE model.To determine the stability of the DFE for model (3), let I = i and S = N − s.Then the linearized differential equation for i is di dt = βi − γi.By taking i = exp(−rτ) we obtain the characteristic equation This equation is transcendental and there is no exact solution, but r > 0 only for β > γ.So we have the same R 0 as the for ODE model.A similar calculation for model (4) yields the same result.This method has been applied to other models by [20,22].Another way of introducing the time it takes the infection to develop is to add another compartment for the exposed or latent population (E).After a contact with an infective, the susceptible individual goes into the exposed compartment and spends an average time equal to the infection time there.Exposed individuals cannot infect.The SIS model with an exposed population is ( The DFE equilibrium of ( 5) is the same as that of (1).Its stability can be calculated using the next generation matrix method and gives R 0 = β/γ.So the stability is the same as that of (1).Exposed populations have been used in many epidemic models such as those in [23,24].
The SIR model of Kermack and McKendrick [4] considers that and infected recovers with total and permanent immunity.It adds the recovered (R) compartment to the SIS model: Since the total population is constant, one equation can be eliminated using R(t) = N − S(t) − I(t).The new system of equations is The only steady solution of ( 6) is I = 0, S + R = N.At the start of the infection S(0) = N.The second equation of ( 6) states that at the initial time, t = 0, dI dt > 0 for β > γ.So for β/γ > 1 the DFE is unstable, and I(t) has its maximum value when S = Nγ/β.Furthermore, the solution as t tends to infinity is S(0) ) and S(∞) > 0 [13].So if there is an epidemic it does not disappear.
The delay due to the time of infection can be introduced in the two ways as conducted for the SIS model.The equilibrium solution of both delay models is the DFE with the same stability as for the no delay model.The same stability result for the DFE of the no delay ODE model is true for the SIR model with the exposed population.
A model for propagation of many different diseases in a population is the SIRS model with no demographics [13].It consists of three populations: susceptibles (S(t)), infectives or infected (I(t)) and recovered (R(t)).A susceptible turns into an infective with a given probability after a contact with an infective.An infective recovers with a rate γ.A recovered loses immunity at a rate ν and converts into a susceptible.The total population is assumed constant, S(t) + I(t) + R(t) = N, and the disease is assumed to be short enough so births and deaths need not be considered.The infection coefficient will be denoted by β/N.All parameters are positive.The system of ordinary differential equations (ODEs) describing the model is [13,14] Since the total population is constant, one equation can be eliminated using R(t) = N − S(t) − I(t).The new system of equations is The SIRS model has two equilibrium solutions: the DFE, S = N, I = 0, R = 0, and the endemic equilibrium, There are several different methods of introducing delays in the infection into the model: The first one which we will call SIRS Model I [17][18][19] is given by where the hypothesis is that an infective has contact with a susceptible and it takes a time τ for this susceptible to turn infective.It is not counted as a susceptible until that time.A second model is the SIRS Model II [25] described by: where the assumption is that after a contact between a susceptible and an infective, a time τ passes before the susceptible turns into an infective and is then able to infect another susceptible.Instead of considering a recovered population, where individuals spend an average time 1/ν before losing immunity, infective individuals after stopping being infective have to wait a delay time 1/ν before turning infectives.In [26], this is performed in addition to keeping the recovered compartment.
The two models with delays introduced above have the same equilibrium solutions but the actual time dependent solutions will be different.There is also the question of the influence of the initial conditions given in the history interval.For an epidemic they can always be given as the equilibrium solutions with no disease.That is, the natural assumption that the epidemic starts at t = 0 can be used.
A model that takes into account the dependence in the past, caused by the delay in the transmission, by including an exposed population, E, is [27,28] The additional parameter is σ = 1/τ.Models of this form are usually called SEIRS models.
It has the same DFE and R 0 as the non-delay model.

Models with Different Immunities
The three models given in Section 2.1 have different types of immunity.The SIS model has no immunity, the SIR model has total and permanent immunity, and the SIRS model has total but temporal immunity.However, for different infections there may be different types of immunity such as partial and permanent or partial and temporal.Some papers considering immunity for influenza are [29][30][31].Ref. [32] includes a reinfected class with different immunity.Ref. [33] models partial immunity of the recovered population and [34] includes that recovered lose immunity but not all infectives recover.Ref. [35] adds the effect of reinfection of influenza.Ref. [36] also takes into account waning immunity.Ref. [37] models the loss of immunity with time-dependent rates.We will consider several new models with different immunity based on the SIRS model.The first model assumes all infective individuals get partial but permanent immunity.Similar models are given by [34,38].Susceptibles get infected and then move into susceptibles with partial immunity S p who can be reinfected with infection rate β p /N.The equations of the model are Since the total population is constant, one equation can be eliminated using S p (t) = N − S(t) − I(t).The system has two equilibrium solutions: the DFE I = 0, S + S p = N, and S = 0, S p = γN β p , I = N − S p .At the start of the infection, S p = 0 so the equation for I is the same as that in system (6), and R 0 has the same value R 0 = β/γ.The two delay models are constructed as before.The first model is and the second is They have the same equilibrium points as the no-delayed model.The model with the exposed population is constructed similarly and is given by All these three models also have the same R 0 as (12).
The second partial immunity model assumes that all infectives get partial immunity but that it is only temporal.Susceptibles get infected, then go to susceptibles with partial immunity and then back to susceptibles with no immunity.The model is given by the following system of ODEs Again since the total population is constant, one equation can be eliminated.The system has two equilibrium solutions: the DFE I = 0, S = N, S p = 0, and a second equilibrium with S > 0, S p > 0 and I > 0. At the start of the infection, S p = 0 so the equation for I is the same as that in system (6) and R 0 = β/γ.The two delay models are similarly constructed as before.They have the same equilibrium points as the non-delayed model.The model with the exposed population has this population added in the same form.All these three models also have the same R 0 as (12).The third model assumes that all infectives get total immunity but that it is only temporal, and then get permanent partial immunity.Susceptibles get infected, then recover with total temporal immunity, then get partial permanent immunity as susceptibles with partial immunity.The model is given by Again since the total population is constant, one equation can be eliminated.The system has two equilibrium solutions: the DFE I = 0, S + S p = N, and a second equilibrium S = 0, S p = γN β p , I = (β p −γ)Nν β p (γ+ν) and R = N − S p − I.This second equilibrium solution only makes sense for β p > γ.At the start of the infection, S p = 0 so the equation for I is the same as that in system (6) and R 0 = β/γ.However, if β p < γ, the epidemic will eventually disappear.
The fourth model assumes that all infectives get total immunity but that it is only temporal, and then get partial immunity which is also temporal.Susceptibles get infected, then recover with total temporal immunity, then get partial temporal immunity as susceptibles with partial immunity, that then convert to susceptibles.The model is given by Again since the total population is constant, one equation can be eliminated.The system has two equilibrium solutions: the DFE I = 0, S + S p = N and a second endemic equilibrium S > 0, S p > 0, I > 0 and R = N − S p − I.At the start of the infection, S p = 0 so the differential equation for I is the same as that in system ( 6) and R 0 = β/γ.

Models with Loss of Immunity and Two Groups
Different individuals lose immunity at different rates.In this subsection for simplicity we will include two groups of individuals, but the models can be extended to more groups.Ref. [34] has a model with some recovered losing immunity.The first proposed model assumes that all infected get permanent immunity but the two groups lose it at different rates, ν and ν 1 .One group has a stronger immune system due to among other things to previous exposure to a previous variation of influenza or vaccinations on previous seasons [29][30][31].
Here, 0 ≤ f ≤ 1 is the fraction of the population in the first group.The two models with delay in the infection are given by and For both delay models, the differential equation at the beginning of the infection process for the infective compartment is the same as that of the SIS model (1), therefore the same procedure can be used to determine that R 0 = β/γ.Also, for both models the equilibrium points are the same as those of the ODE model ( 19).
The model with an exposed population is It has the same DFE as (19).The stability of the DFE at the beginning of the infection is still given by R 0 = β/γ.In the next model, a fraction f of the infectives gets total and permanent immunity, and the rest gets partial and permanent immunity.The ODEs of the model are The next model assumes that a fraction f gets total and temporal immunity and the rest get partial and temporal immunity.The system of equations is In the next model all infectives get total immunity for a finite time, then lose it different rates ν and ν 1 , and get partial and permanent immunity with different rates of infection β 1 and β 2 .

Models with Vaccination and Immunity Loss
Vaccination is a very effective tool for achieving immunity to seasonal influenza.Usually the effect of the vaccine wanes with time.We consider that the immunity is permanent if it lasts until the next flu season with a new variant of the virus.The vaccinated population will be denoted by V.There are different ways of introducing the vaccination rate.The first way is that it is proportional to the number of susceptibles S [39][40][41].Refs.[42][43][44] also consider the waning of the vaccine.Ref. [45] has a basic vaccination model based on the SIR model.A second way is to assume that only a percentage p r of susceptibles is willing to be vaccinated [46] which is a more realistic assumption [47].This leads to a vaccination rate of the form r(p r S − V).It can also be considered proportional to the total population or just constant [48,49].Ref. [50] considers vaccination, partial immunity and reinfection.Refs.[40,51] also include quarantine.Our first model with vaccination assumes infectives recover and get total and temporal immunity for an average length of time 1/ν.Vaccinated individuals get total immunity but it wanes with time at a rate ω.The model is In this model the rate of vaccination is proportional to the number of susceptibles.Since the equation for the infective population is the same as that for (6), the value of R 0 is still β/γ.If the vaccines does not wane, ω = 0, the system has only one equilibrium point, S = 0, I = 0, R = 0, V = N.If ω ̸ = 0 there are two equilibrium points, I = 0, R = 0, S = (ωN/(ω + r), V = rS/ω and S = γN/β, V = rS/ω, R = γI/ν and I given by the solution to I = N − S − V − R. The first delay model is obtained by replacing S(t) and I(t) in the infection terms by S(t − τ) and I(t − τ).The second one by only replacing the I(t) term.Both have the same equilibrium points and R 0 as (26).The system with an exposed population has the same DFE equilibrium points since I = 0 implies E = 0.The R 0 calculated by the next generation matrix method is still R 0 = β/γ.The second model with vaccination considers that the infective population gets total immunity for a fixed time and then gets partial permanent immunity.The vaccinated population gets total but temporal immunity, and then it wanes to partial permanent immunity.There is a new compartment, S p , representing susceptibles with partial immunity.
This model has a DFE I = 0, S = 0, R = 0, V = 0, S p = N.There is also an endemic equilibrium with S = 0, V = 0. Since at the beginning of the infection, S p = 0, using the differential equation for I, R 0 = β/γ.The two delay models are obtained modifying ( 27) in a similar way as was performed for previous models.The equilibrium states and R 0 are the same as those of ( 27).The corresponding equation with an exposed population also has the same R 0 .
If the vaccination rate is modified to r(p r S − V) instead of rS, to make the model agree with not everybody wanting to be vaccinated, the modified model is This model does not have a steady solution with I = 0 for arbitrary values of the parameters.This an important difference in the dynamics with those of the other models.The next model considers that all vaccinated get total immunity for a time, then a fraction f v keeps the total immunity, while the rest gets partial permanent immunity.All infectives get total temporal immunity after they recover.
If ω = 0 the DFE is I = 0, R = 0, , S p = 0, V = p r S, S = N/(1 + p r ).For ω ̸ = 0 the DFE is I = 0, R = 0, V = 0, S = 0, S p = N.Since S p (0) = 0, we have for both case R 0 = β/γ.The two delay models and the exposed population models are derived from (29) as was performed for the previous models.All three models have the same DFE and R 0 as (29).
As a final and more complete model we will add the use of antivirals to reduce the infection time and isolation to reduce the number of active infectives.There is a new compartment, the quarantined or isolated population Q. Infectives recover with total and temporal immunity, vaccinated get total and temporal immunity and then partial permanent immunity.A fraction f a of the infected gets an antiviral and exit the compartment at a faster rate γ a .A fraction f q of the infected gets quarantined at a rate γ q and leaves the quarantine state at a rate γ q1 .The average time that a quarantined individual spends in the infected and quarantined compartments is the same as the time a non-quarantined individual spends in the infective compartment.Therefore, 1/γ = 1/γ q + 1/γ q1 .The model is The DFE is The delay versions and the model with exposed population have the same DFE and R 0 .

Results
In this section, we present numerical results for the models given in the previous section.The numerical simulations, including the DDE were performed using the julia package DifferentialEqns.jl[52,53].The parameter values used for most simulations are given in Table 1.If a different value is used for a particular simulation, it will be noted in the caption of the corresponding figure.The values of β, γ, ν and τ are based on those given [54] and CDC [1] corresponding to seasonal influenza.All the others are estimated.The value of β has a large variation depending on the conditions: rural, urban, suburban, high density work places, etc.The values of β and γ chosen for most simulations give R 0 = 1.5.Ref. [55] also gives similar parameter values for influenza.The initial conditions for all models are S = 900, I = 100 and all other populations equal to zero.The history for the delay models is S = 1000, and all other delayed variables equal to zero in the interval [−τ, 0].This corresponds to an epidemic starting at t = 0.All simulations, except where specified, have R 0 = 1.5.

SIS, SIR and SIRS Models
The SIS, SIR and SIRS models are widely used for epidemics with by direct contact.In the SIS model, infectives do not get immunity, in the SIR model they get total and permanent immunity and in the SIRS model they get total but temporal immunity.
Figure 1 presents the plots for the SIS model, in the left for the susceptible population and in the right for the infective population.Since R 0 = β/γ > 1 the solution tends to the endemic state.The infective plot for the model with the exposed population has a different limit since there are two infective classes, I and E, and only I is plotted.The solutions have the same form but the one for the ODE model reaches the equilibrium state the fastest.The solutions of the two delay system are all identical since the only difference is when a new infected stops being counted as a susceptible.The model with the exposed population has an average time in the exposed compartment equal to the delay but not the same delay for every individual.Real data is needed to determine which model is the best, since all the models, with no delay, with delay, and with an exposed population, have the same qualitative dynamics.

Models with Different Immunities
In this subsection, we consider models with different immunities: total, partial, permanent and temporal in different combinations.The first model (12) considers that all infectives after recovering get partial and permanent immunity.Figure 4 has the simulated solutions for the infective population for two different values of the infection rate for the population with partial immunity, β p .On the left for β p = 0.25 and the solution tends to I = 0 in spite of having R 0 > 1, and on the right for β p = 0.4 and the solution tends to the endemic state.However, whether the epidemic ends or not depends on the values of S p β p + Sβ − Nγ which appears on the right-hand side of the differential equation for the infective population.The two delay models behave differently at the beginning since in the model with delay I, both S and I depend on the history, while for the model with delay II only I depends on the history.
The next model ( 16) considers that all infectives recover with partial and temporal immunity.Figure 5 has the plots for the infective population, on the left for β p = 0.25, and I tends to a small but non zero value, and on the right for β p = 0.4, and I tends to a much larger value, showing the effect of β.Model (17) is for an infection with total but temporary immunity and then partial and permanent immunity.Figure 6 has the simulated plots for the infective population for two values of β p .The figure on the left is for β p = 0.255 and I tends to the DFE; the figure on the right is for β p = 0.4 and I goes to an endemic equilibrium.For both simulations R 0 > 1 so initially the number of infectives grows.However, whether the epidemic ends or not depends on the values of S p β p + Sβ − Nγ.All models have the same equilibria.The next model (18) supposes that after being infected an individual gets total immunity which is temporal, and then partial immunity that is also temporal, finally losing all immunity.Figure 7 shows the time solution of the infective population.On the left for β p = β/2 and on the right for β p = 0.4.In both cases the epidemic persists.The limit value for β p = 0.4 is slightly larger.Since all individuals eventually lose immunity completely, the epidemic does not end.

Models with Loss of Immunity and Two Groups
Different individuals have different immune system reactions to an infection.For example, all individuals may get total immunity to the infection after recovering but may lose it at different rates.For seasonal influenza some of the reasons are that children have been exposed to fewer variants of the virus, that older individuals have weaker immune systems, or that some individuals may have other risk factors such as smoking.For simplicity we consider just two classes of individuals, but the models can be extended to an arbitrary number.
The first model with two distinct populations (19) considers that the two groups get total and temporal immunity after infection (the usual SIR assumption), but lose the immunity at different rates ν and ν 1 , with fractions f r and (1 − f r ), respectively.The simulations are in Figure 8.Both figures have R 0 > 1, the left one with β = 0.5 and the right one with β = 0.666.The limit value for the right figure is about twice that of the left showing the large influence that the infection rate has.The next model (23) assumes that one class gets total and permanent immunity, while the second gets partial and also permanent immunity.The fraction of individuals in both classes is taken to be f p = 0.5.For the class with partial immunity β p = β/2.Since all individuals get permanent immunity it is expected that the epidemic will not persist.Figure 9 has two simulations for different values of β.The left one for β = 0.5 and the right one for β = 1.5.Even for the fairly large second value of β, that corresponds to R 0 = 4.5, the epidemic dies out.Next we consider that half of the infective populations get total and temporal immunity, and half partial and also temporal immunity (24).Simulations were conducted for two values of β p , the infection rate of the population with partial immunity, β p = β/2 (left plot of Figure 10) and β p = β/8, (right plot).In both cases the epidemic persists and the limit value of I has a relatively small change.The fact that all the individuals get immunity only for a finite time is the reason that the infection persists.Model (25) assumes that all infectives recover with total temporal immunity and then get partial and permanent immunity.A fraction f r = 0.75 has partial immunity with infection rate β 1 = 0.9β and loss of permanent immunity rate ν = 1/30 and the rest has infection rate β 2 = 0.8β and immunity loss rate ν 1 = 1/60 (left plot of Figure 11).The solution tends to the endemic equilibrium.A second simulation was conducted changing the infection rates to β 1 = 0.45β and β 2 = 0.4β (right half of Figure 11).Now the solution tends to the DFE.In the second simulation β 1 /γ < 1 and β 2 /γ < 1 so when most recovered individuals have partial permanent immunity, the differential equation for I shows that I decreases.

Models with Vaccination and Immunity Loss
Vaccines can help control seasonal influenza.Since the influenza virus mutates from season to season, it is hard to predict what is the efficiency of the vaccine.It is also difficult to determine if the vaccine gives total or partial immunity, or if the effect wanes.A vaccinated person may not get influenza due to the vaccine or just because the person was not in close contact with an infective.Similarly, a person becoming infected only once may be due to having acquired total and permanent immunity, or not having close contacts with infectives.In this subsection, we present simulations conducted for different scenarios all with vaccination.
The first scenario represents individuals that get vaccinated at a rate rS, and achieve total immunity.The immunity from becoming infected with the disease is total and disappears at a rate ν, (26).In the left part of Figure 12, the effect of the vaccine wanes at a rate ω = 1/60.In the right part, the effect of the vaccines does not wane, ω = 0.In both cases the solution tends to the DFE.In the first case there is also a second endemic equilibrium.The model in the next scenario is similar to the model given by (26), except that when the vaccine wanes the individual still gets partial permanent immunity with β v = 0.8β. Figure 13 on the left shows the case with ω = 1/60 and the simulated solution almost gets to the DFE but then it rebounds.On the right, the plot corresponds to ω = 0 and the solution trends to the DFE.The figure shows the importance that the waning of the vaccine has.The next simulations were performed with model (28).The vaccination rate is r(p r S − V).For this vaccination rate there is no DFE.The infected get to recovered that then get partial permanent immunity.Vaccinated get total immunity that, in one simulation wanes to partial permanent, and in the second one does not wane.We use p r = 0.4, and the infection rate for the population with partial immunity β v = 0.8β.The left plot of Figure 14 is for the simulation with vaccination waning, and the one of the right for no waning.In both plots the infective solution approaches close to zero and then rebounds.The next simulations are performed with the model given by (29).The vaccination rate is r(p r S − V) with p r = 0.6 and a fraction f v = 0.75 of the vaccinated population gets total and permanent immunity while a fraction 1 − f r gets total immunity that wanes at a rate ω to partial and permanent immunity with infection rate β v = 0.8β.The immunity from an infection is total and temporal.The calculated solution tends to an endemic equilibrium.A second simulation was performed with no waning of the effect of the vaccine ω = 0.The solution tends to the DFE.The plots are in Figure 15.This example shows the importance of the waning or not of the vaccine.The final model considers vaccination, the use of antivirals and isolation.The vaccination rate is r(p r S − V), a fraction f v = 0.75 of the vaccinated keeps the total immunity forever and the rest loses the total immunity and then gets permanent partial immunity with infection rate β v = 0.8β.The infected recover with total and temporal immunity.A fraction f a = 0.2 is given antivirals and is removed from the infective class at a rate γ a = 2γ.A fraction f q = 0.4 is removed from the infectives and isolated at a rate γ q = 1.Two simulations were performed.The first one with β = 0.5, that corresponds to R 0 = 1.5, is shown in the plot in the left of Figure 16 and has a solution that tends to the DFE.The second one, with β = 1 that gives R 0 = 3, is shown in the right plot and the solution tends to an endemic equilibrium.The figure shows that while R 0 predicts whether an epidemic starts or not, it does not indicate whether it will eventually die out or not.

Discussion
SIR-based models have been used for may infectious diseases transmitted by direct contact.In this paper, we used models with no demographics for short-lived infections such as seasonal influenza.The first reason is that the infection has a short duration and that it is not very deadly.Therefore, the number of births and deaths during the infection period is small.The second reason is that it allows us to concentrate on the effects of losing immunity.The basic reproduction number R 0 is a very important index that gives a criterion for determining whether there will be an epidemic or not.It is usually calculated as we did in this paper by assuming that at the start of the infection there is one infective and N − 1 susceptibles with N the total population.For the majority of our examples it is R 0 = βS(0)/(γN) and simplified to R 0 = β/γ.If there is also an initial time of the infection number of vaccinated individuals V(0) > 0, then R 0 = β(N − V(0))/(γN).And, if for example without vaccination R 0 = 1.5, by vaccinating more than 1/3 of the total population then R 0 < 1 and there is not an epidemic.A similar estimate can be performed if initially there is also a number of individuals S p > 0 with partial immunity and infection rate β p < β.
The models presented start with the simple SIS, SIR and SIRS, that have no immunity, total and permanent immunity, and total and temporal immunity.They were first extended to include other situations with combinations of total and partial immunity, both permanent and temporal.For diseases like seasonal influenza it is very difficult to determine if an infected and recovered individual received total or partial immunity after the infection.Or if that immunity disappears or not.The individual may not get infected again due to total immunity or partial immunity or just because the individual was not in close contact with any infective.The different models permit running scenarios with different hypotheses and comparing the results.Two population groups with different immune reactions were considered next and different scenarios were run.Finally, vaccination was introduced with the different models having different assumptions and the final model being more complete having also antivirals and isolation to reduce the infection.For the various models that had DFEs simulations were performed where the epidemic died and where it did not.The effect of introducing a delay in the time of infection was also studied.The delay was introduced in two different ways using DDEs and also adding an exposed compartment.Apart from the delayed solutions taking more time to reach the equilibrium states, the behavior is very similar to that of the corresponding ODE model.That is, they produce solutions with the same equilibrium states and the same basic reproduction number.However, adding the delay makes the model more realistic and it is expected that the delayed models will fit real data better.
The parameters related to the time between contact and infection and to the duration of the infection are well-known for seasonal influenza.The infection rates can be estimated from the number of infectives.They vary from city to city and from population group to population group.As mentioned above, the parameters related to partial immunity and vaccination effects are hard to estimate and not reported.Numerical simulations using different models and comparing with actual epidemic data can help determine the more realistic hypotheses.
Since the two delay versions of each model have different characteristic equations, their solution behaviors can be different from each other and different from the behavior of the solutions of the non-delayed model.However, for the parameter values used in the simulations that did not happen and the qualitative dynamics were all similar.Determining whether an individual has total or partial immunity or has not come into close contact is not an easy task.We proposed different models with different types of immunity.Comparisons with real data are needed to establish which one is more realistic.
Future work will consider different initial conditions, the introduction of randomness and spatial dependence.The use of other compartments such as asymptomatics is another route.

Figure 2
Figure 2 presents the plots for the SIR model, in the left for the susceptible population and for the infective in the right.Note that for I = 0, there are multiple solutions S + R = N.The simulated solutions depend on the initial values.The delay model I uses the history for both S and I, while all the other models only use S(0).

Figure 3
Figure 3 presents the plots for the susceptible (left) and infective (right) populations for the SIR model.This model has total but temporal immunity.As expected, the simulated solutions tend to the endemic equilibrium.The ODE model gives the highest infection peak, all the models have the same qualitative dynamics.

Figure 9 .
Figure 9. Simulations for the model with two classes, one with total and permanent immunity, and the other with partial and permanent immunity.(left) β = 0.5.(right) β = 4.5

Figure 10 .
Figure 10.Simulations for the model with two classes, one with total and immunity, and the other with partial immunity, both temporal.(left) β p = β/2.(right) β p = β/8.

Figure 13 .
Figure 13.Model with vaccination.Infected individuals get total immunity that wanes to partial and permanent immunity.Vaccinated individuals get total immunity that may wane or not.(left) Vaccination immunity wanes to partial permanent.(right) Vaccination immunity does not wane.

Figure 14 .
Figure 14.Model with recovered getting partial permanent immunity and vaccinated getting total immunity that in one case wanes and in the other does not.(left) Vaccination immunity wanes.(right) Vaccination immunity does not wane.

Figure 15 .
Figure 15.Model with two classes of vaccinated individuals, one with total and permanent immunity and the other with total immunity that wanes to partial and permanent immunity.(left) Waning vaccination immunity.(right) Permanent vaccination immunity.

Table 1 .
This table has the parameter values used in the simulations.