A shiny app for modeling the lifetime in primary breast cancer patients through phase-type distributions

1.   Introduction

Mathematical modeling becomes an essential tool in the description and analysis of diseases dynamics ^[1,2]. In the last decades, many models have been used to understand the dynamics of infectious diseases with the purpose of controlling them ^[3,4]. In this context, the choose of the variables that characterize population dynamics, as the epidemiological state of individual, the tracking of the relevant biological processes, for example, how infection is transmitted, and the decision about the type of model to be used among the available options such as discrete models ^[5,6], ordinary differential systems ^[7,8,9,10], age-structured PDE ^{[10,11,12,10,13,14]}, delay differential systems ^{[15,16,17,18]}, stochastic models ^[19,20,21], and etcetera, depend on the main question addressed.

Vaccination is the one of the most efficient way to halt disease transmission through promoting population immunity ^[22,23]. The designed of vaccination strategies are based on the type of the infectious agent (viruses, bacteria, fungi, protozoa, or worms) and always search for risk groups ^[24], thresholds such as the proportion of the population to vaccinate ^[24,2], and the optimum age for vaccination ^[14,25], with the aim of optimize disease control. The duration of immunity promoted by vaccination and its efficacy determine the number of doses and the interval among them to ensure that the individuals are protected. For example, recently, human papillomavirus (HPV) vaccination moves from a 3–dose schedule to a 2–dose schedule ^[26]. The change in vaccination schedule was due to evidences that the antibody response generated after 2 doses is enough to prevent virus infection ^[27]; besides vaccination cost is significantly reduced in this 2 dose–scheme. HPV vaccination target-individuals between 9 to 14 years age because exposure to infection is higher at younger ages with a peak after the debut of sexual activity ^[28]. To prolong the immunity conferred by certain vaccines, it is sometimes necessary to update them. It is easier to focus on the individuals that are already vaccinated to incite them to update their vaccine. Indeed, these individuals are already known and easier to encourage to be vaccinate again.

For kids, the immunizations schedule depends on where they live, the child's health, the type of vaccine, and also which vaccines are available. For example, the DTaP vaccine (diphtheria, tetanus and acellular pertussis) is applied at 2, 4, 6, 15 months, and 4 years old ^[29]. For influenza, starting after six months old, the flu vaccine is recommended every year as the main virus that are circulating change and evolve ^[29]. According to WHO (World Health Organization), a collaborative global vaccination program was able to eradicate smallpox in 1980. But currently, diseases like whooping cough, polio, measles, and rubella that were controlled or almost eradicated are appearing again because of vaccine refusals, under-vaccination, waning immunity, less effective immunizations, and imported cases ^[30,31].

The duration of protection provided by any mechanism (vaccination is one of them) plays an important role on the evolution and control of epidemics. Although a lot of models in the literature address one or several questions related to vaccination, few of them considered the lost of immunity ^[23,32,33], and to our knowledge, no one considered the already vaccinated individuals that need to update their vaccine. Mathematical models of vaccination have been studied since 1760 ^{[34,35,36,37,38,39]}. To our knowledge, D. Bernoulli ^[34] proposed the first mathematical model of vaccination. He studied the impact of smallpox vaccination on the life expectancy of the immunized population. The works of ^{[35,36,37,39]} deal with age structured models of vaccination. They considered vaccine-induced temporary or permanent immunity. They studied the asymptotic behavior of the steady states. Recently, ^[38] have considered a model of vaccination with temporary immunity described by a delay differential system. The delay represents the length of immunity period. However, the population of individuals that update their vaccine at the end of their period of protection has never been explicitly incorporated in these models. It is sometimes difficult to reach a reasonable percentage of people to vaccinate in the total population to halt the disease transmission. Therefore, It would be interesting to combine vaccination of a part of total population with a proportion of individuals that were previously vaccinated.

With this in mind, we propose a new mathematical model that take into account the temporary protection and the specific individuals that are at the end of their previous period of protection. The model is an extension of the classical Kermack and McKendrick model ^[9] which includes a compartment of individuals with a temporary protection. Although, the model is formulated as a direct disease transmission model, it can be adapted to take into account the saturation in the transmission rate observed in vector-borne diseases. The epidemic model is present in section 2. Thus, we reduce the model (section 3) by using the method of characteristics to a delay differential-difference system. We present also some results about the existence, uniqueness, positivity and uniform boundedness of the solutions. Section 4 is devoted to the study of the existence of the steady states: disease-free and endemic. Section 5 concerns the computation of the basic reproduction number $\mathcal{R}_0$ of the model and its comparison with the $\mathcal{R}_0$ of the classical Kermack-McKendrick model. In Section 6, the local asymptotic stability of the steady states is established. In section 7, we establish the uniform persistence of the infected individuals. In section 8, we show that if $\mathcal{R}_0 > 1$, then the endemic steady state is globally asymptotically stable, and if $\mathcal{R}_0 < 1$, then the disease-free steady state is globally asymptotically stable. Finally, section 9 summarizes the main results and conclusion.

2.   Mathematical model

The epidemiological model splits the total population $N$ in four classes: susceptible ($S$), infected ($I$), recovered ($R$) and protected individuals ($P$). It is an extension of the classical SIR Kermack and McKendrick model ^[9] that includes a compartment of protected individuals with limited duration $\tau$. Let $p: = p(t, a)$ be the age distribution of the population of protected individuals (the age in this model is the time since an individual is temporarily protected). So, the total number of protected individuals at time $t$ is

The model is given by

The evolution of the density of the protected individuals is given by

The boundary condition is

The system is combined with nonnegative initial conditions

All the parameters of the model are nonnegative constants, and they are described in Table 1. The parameter $\alpha \in (0, 1)$ represents a specific protection rate which corresponds to the individuals that get protected again at the end of the previous period of protection. These individuals may represent a population of volunteers or a specific age group, for example children, that are always vaccinated. Figure 1 provides a schematic representation of the epidemiological model.

3.   Reduction to a differential-difference system

Using the characteristics method, see for instance ^[40], we obtain, for $t>0$ and $a \in [0, \tau]$,

We put

Then the expression of $p(t, a)$ becomes, for $t>\tau$ and $a\in [0, \tau]$,

Consequently,

Finally, we obtain the following system

By doing a translation in time ($t\mapsto t-\tau$), we can consider the initial conditions

By adding the equations of $S$, $I$, $R$ and $P$, we show that the total population $N = S + I + R + P$ satisfies

where

It follows that

$P$ and $R$ depend on $S$, $I$ and $u$. However, the equations of $S$, $I$ and $u$ are independent on $P$ and $R$. Then, we can focus on the reduced system

which is completed by the initial conditions

System (3.2)–(3.3) is a coupled system of differential and difference equations with discrete delay.

Remark 1. The particular case when $\tau = 0$ reduces the system (3.2) to the classical Kermak and McKendrick model ^[9] given by

Let us introduce $C: = C([-\tau, 0], \mathbb{R})$, the space of continuous functions on $[-\tau, 0]$ and $C^+: = C([-\tau, 0], \mathbb{R}^+)$, the space of nonnegative continuous functions on $[-\tau, 0]$. Throughout this paper, we assume $\varLambda$, $\gamma_S$, $h$, $\gamma_p$, $\tau$, $\mu$, $\gamma_I$, $\beta\geq0$, $\alpha \in (0, 1)$, $S_0\geq0$, $I_0\geq0$, $\phi\in C^+$. The existence and uniqueness of nonnegative solutions of (3.2)–(3.3) can be obtained as in ^[41]. Also, we observe that by the method of steps we can solve the system (3.2)–(3.3) in each interval $[k\tau, (k+1)\tau]$, for $k = 0, 1, \dots$.

Consider the auxiliary linear homogeneous difference equation

where the function $u_t \in C$ is defined, for $t\geq 0$ and $u \in C([-\tau, +\infty), \mathbb{R})$, by $u_t(\theta) = u(t+\theta)$ for $\theta \in [-\tau, 0]$, and the operator $\mathcal{D} \colon C \to\mathbb{R}$ is given, for $\psi \in C$, by

Remark that

with $\|\psi\| = \sup_{\theta\in[-\tau, 0]}\lvert\psi(\theta)\rvert$. The condition (3.5) says that the zero solution of the linear difference equation (3.4) is globally asymptotically stable ^[42].

Now, we deal with the nonnegativity of the solutions of the system (3.2).

Proposition 1. All the solutions $(S, I, u)$ of the system (3.2) with nonnegative initial conditions are nonnegative. Furthermore, $(S, I)$ has a continuous first derivative for all $t > 0$ and $u$ is continuous for all $t\geq-\tau$ if and only if the initial condition $(S_0, I_0, \phi)$ satisfies the compatibility condition

Proof. Let $(S, I, u)$ be a solution of (3.2) associated to the initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+ \times \mathbb{R}^+\times C^+$. We first prove the nonnegativity on the interval $[0, \tau]$, and we apply the same reasoning by steps on each interval $[k\tau, (k+1)\tau]$, for $k = 1, 2 \dots$. For $t\in[0, \tau]$, we have $t-\tau\in[-\tau, 0]$. Then, the system (3.2) becomes

The idea is to extend the analogous result known for ODE to our system as established in (^[43], Theorem 3.4). We have the following implications

and

This implies $S(t)\geq 0$ and $I(t)\geq0$ for $t\in[0, \tau]$. Then, $u(t) = h S(t)+\alpha e^{-\gamma_p \tau}\phi(t-\tau)\geq 0$, for $t \in [0, \tau]$. Hence, one just repeats the same argument by steps. We conclude that $S$, $I$ and $u$ are nonnegative on $[0, +\infty)$. We obtain a nonnegative piecewise solution $(S, I, u)$ of (3.2). We can easily prove that $(S, I)$ has a continuous first derivative for all $t > 0$ and $u$ is continuous for all $t\geq-\tau$ if and only if $\phi(0) = h S_0+\alpha e^{-\gamma_p \tau}\phi(-\tau)$.

Next, we investigate the boundedness of the solutions of (3.2). We propose to prove the following result.

Proposition 2. The solutions of the system (3.2) are uniformly bounded.

Proof. Let $(S, I, u)$ be the solution of (3.2) associated to the initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+\times \mathbb{R}^+ \times C^+$. First, one can observe from (3.1) that

Then, $S$ and $I$ are uniformly bounded. Moreover, for $t>0$, we have (see ^[44] or Lemma 3.5 in ^[45])

with $\nu >0$, $C>0$. This implies that $u$ is bounded. On the other hand, we have from the equation of $u$

Then,

This complete the proof.

4.   Existence of steady states

In this section, we establish the existence of the steady states of the system (3.2). Let $(S^*, I^*, u^*)$ be a steady state of (3.2). Then,

The third equation of (4.1) implies that

From the second equation of (4.1), we have for $\mu+\gamma_I>0$ and $\beta>0$

Suppose that $I^* = 0$. Then, $S^*$ satisfies the following equation

Consequently, we obtain for $\gamma_S+h-(\alpha\gamma_S+h) e^{-\gamma_p \tau} >0$,

Remember that $\alpha \in (0, 1)$. Then by assuming $\gamma_S+h>0$, we have $\gamma_S+h-(\alpha\gamma_S+h) e^{-\gamma_p \tau} >0$. We put, for $\gamma_S+h>0$,

So, under the condition $\gamma_S+h>0$, $(S^0, 0, u^0)$ is always a steady state of (3.2). It describes the disappearance of the epidemic. We will refer to this steady state as the disease-free steady state.

Suppose now that $I^*> 0$. Then, $S^* = (\mu+\gamma_I)/\beta$ with $\mu+\gamma_I>0$ and $\beta>0$. We have

Moreover, the first equation of (4.1) implies that

In fact, $I^*$ is given by the following expression

Then, the existence of a positive steady state is equivalent to

We set

We will refer to this steady state as the endemic steady state.

We summarize the existence of steady states in the following result.

Theorem 4.1. Assume that (4.3) holds. Then, the system (3.2) has two distinct steady states: a disease-free steady state $(S^0, 0, u^0)$, which is given by (4.2), and an endemic steady state$(\overline{S}, \overline{I}, \overline{u})$, which is given by (4.4). If (4.3) does not hold, then $(S^0, 0, u^0)$ is the only steady state.

In the next section, we derive the basic reproduction number $\mathcal{R}_0$. We study also the influence of some parameters on this threshold.

5.   The basic reproduction number $\mathcal{R}_0$

The number $\mathcal{R}_0$ is defined as the average number of secondary infections that occur when one infective individual is introduced into a completely susceptible population. By dividing the equation of $I(t)$, in the system (3.2), by $(\mu+\gamma_I)I$ we get

The fraction $\beta/(\mu+\gamma_I)$ can be interpreted as the number of contacts per infected individuals during their infectious period that lead to the transmission of the disease. If

the disease persist, otherwise, it disappears. Then, the basic reproduction number of the disease is defined by

Remark that the condition (4.3) is equivalent to $\mathcal{R}_0>1$. Then, the existence of the endemic stable state is guaranteed by the condition $\mathcal{R}_0>1$. In comparison with the classical Kermack-McKendrick SIR epidemic model, the new parameters in our model are $h$, $\alpha$ and $\tau$. In the next proposition we study the behavior of $\mathcal{R}_0$ in terms of these new parameters.

Proposition 3. $\mathcal{R}_0$ is a decreasing function with respect to $h$, $\alpha$ and $\tau$. Furthermore,

where $m: = h \in J: = [0, +\infty) \mapsto \mathcal{R}_0(h)$, $m: = \alpha \in J: = [0, 1] \mapsto \mathcal{R}_0(\alpha)$ or $m: = \tau \in J: = [0, +\infty) \mapsto \mathcal{R}_0(\tau)$, and

We remark that ${\inf_{m\in J} \mathcal{R}_0(m)}$ is in fact ${\min_{m\in J} \mathcal{R}_0(m)}$ for $h$ and $\alpha$. The proof of this proposition is easy, so we drop the details. The behavior of $\mathcal{R}_0$ as a function of $\tau, \alpha$ and $h$ is shown in Figures 2 and 3.

Remark 2. Observe that:

● If $\tau = 0$ or $h = 0$, then $\mathcal{R}_0$ becomes

which is the basic reproduction number associated to the classical Kermack-McKendrick model (see Remark 1). If we increase $\tau$ or $h$, $\mathcal{R}_0$ decreases until the value ${\min\limits_{m\in J} \mathcal{R}_0(m)}$, $m = \tau$ or $h$, which is given in Proposition 3 (see Figures 2 and 3).

● If $\alpha = 0$, then $\mathcal{R}_0$ becomes

We can decrease $\mathcal{R}_0$ by increasing $\alpha$ until the threshold given by ${\min\limits_{\alpha\in J} \mathcal{R}_0(\alpha)}$ (see Figures 2 and 3).

6.   Local asymptotic stability

The purpose of this section is to study the local asymptotic stability of each steady state in the cases $\mathcal{R}_0<1$ and $\mathcal{R}_0>1$. We use the results in ^[44] to linearize the differential-difference system (3.2) about the steady states and to derive the characteristic equations. As we will often vary the delay $\tau$ to study the stability of the steady states, we consider when it is necessary the dependence of $\mathcal{R}_0$ in terms of the parameter $\tau$, $\mathcal{R}_0: = \mathcal{R}_0(\tau)$.

The linearized system of (3.2) (see ^[44]) about any steady state $(S^*, I^*, u^*)$ is given by

We can check (see ^[44]) that the characteristic equation of this system is

Remark 3. Let $\tau>0$. The characteristic equation (6.1) has the form

Then, any sequence of distinct roots $\{\lambda_n\}$ of this last equation satisfies $\lim\limits_{n \to +\infty}\left|\lambda_n \right| = +\infty$. Furthermore, this sequence approaches the roots of the equation

which are given by

As $\alpha \in (0, 1)$, we conclude that all the branches of eigenvalues that appear from infinity have negative real parts.

6.1.   Local asymptotic stability of the disease-free steady state

The linearized system of (3.2) about the equilibrium $(S^0, 0, u^0)$ is

and the characteristic equation is given by

The following proposition deals with the instability of disease-free steady state.

Proposition 4. Assume that $\mathcal{R}_0>1$. Then, there exists a positive real root of (6.2), and the steady state $(S^0, 0, u^0)$ is unstable.

Proof. From the characteristic equation (6.2), we have the following eigenvalue

Clearly, this eigenvalue is real and positive when $\mathcal{R}_0>1$. The proof is complete.

Suppose now that $\mathcal{R}_0<1$. Then,

Thus, the local stability of $(S^0, 0, u^0)$ is determined by the sign of the real part of $\lambda \in \mathbb{C}$ satisfying

We have the following theorem.

Theorem 6.1. Assume that $\mathcal{R}_0<1$. Then, all roots of the characteristic equation (6.2) have negative real parts, and the steady state $(S^0, 0, u^0)$ is locally asymptotically stable.

Proof. Our approach is to see the stability of $(S^0, 0, u^0)$ when the delay is equal to zero and, by using the continuity and Remark 3, we check if the stability can be lost by the appearance of a pure imaginary roots. We consider in this proof $\mathcal{R}_0$ as a function of $\tau \in [0, +\infty) \mapsto \mathcal{R}_0(\tau)$. Setting $\tau = 0$ with $\mathcal{R}_0(0)<1$, there exists only one root of (6.3) given by $\lambda = -\gamma_S$ which is negative. We conclude that $(S^0, 0, u^0)$ is locally asymptotically stable when $\tau = 0$.

Hence, we look for purely imaginary roots $\pm i\omega$, $\omega \in \mathbb{R}$. Remark that if $\lambda$ is a root of (6.3) then its conjugate $\overline{\lambda}$ is also a root of (6.3). Then, we can look for purely imaginary roots $i\omega$ with $\omega>0$. We put

Then, by separating real and imaginary parts in (6.3), we obtain

This last system is equivalent to

It follows, by taking $\cos^2(\omega\tau)+\sin^2(\omega\tau) = 1$, that

We can observe that $\rho-(\gamma_S+h)<0$ and $1-\eta>0$, which is absurd. Then, no $\pm i\omega$ satisfying (6.3) exist. Hence, when $\mathcal{R}_0(\tau)<1$ all roots of (6.2) have negative real parts. Then, $(S^0, 0, u^0)$ is locally asymptotically stable.

6.2.   Local asymptotic stability of the endemic steady state

In this section, we show the local asymptotic stability of the endemic steady state. We assume that

The linearized system of (3.2) about $(\overline{S}, \overline{I}, \overline{u})$ is given by

The characteristic equation of this system is given by

To study the local asymptotic stability of the endemic steady state, we use the same technique as for the disease-free steady state. We consider again $\mathcal{R}_0$ as a function of $\tau \in [0, +\infty) \mapsto \mathcal{R}_0(\tau)$. We have the following lemma.

Lemma 6.2. For $\tau = 0$ and under the condition $\mathcal{R}_0(0)>1$, the characteristic equation (6.1) has only roots with negative real parts.

Proof. As $\tau \in [0, +\infty) \mapsto \mathcal{R}_0(\tau)$ is a decreasing function, the assumption $\mathcal{R}_0(0)>1$ implies that $\mathcal{R}_0(\tau)>1$, for all $\tau \geq 0$. For $\tau = 0$, the characteristic equation (6.1) becomes

It is clear that all coefficients of the above equation are positive. Then, the Routh-Hurwitz criterion implies that all the roots have negative real parts. This is corresponding, when $\mathcal{R}_0(0)>1$, to the local asymptotic stability of $(\overline{S}, \overline{I}, \overline{u})$. The previous lemma states that $(\overline{S}, \overline{I}, \overline{u})$ is locally asymptotically stable for $\tau = 0$. Using this assertion and Remark 3, we show in the following theorem that the eigenvalues of (6.4) stay in the left half plane for all $\tau>0$ with $\mathcal{R}_0(\tau)>1$.

Theorem 6.3. Assume that $\mathcal{R}_0(\tau)>1$. Then, all roots of (6.1) have negative real parts, and the steady state $(\overline{S}, \overline{I}, \overline{u})$ is locally asymptotically stable.

Proof. As in the proof of Theorem 6.1, we show that no purely imaginary roots $\lambda = \pm i\omega$, $\omega \in \mathbb{R}$ exist. It is sufficient to look for purely imaginary roots $i\omega$ with $\omega>0$. Then, by separating real and imaginary parts in (6.4), we obtain

where

It follows, from the system (6.5), that $\omega$ satisfies

This implies that

Remember that $1-(\alpha c)^2>0$. We put

Then, $x = \omega^2>0$ satisfies

If $D$ is nonnegative, then the Routh-Hurwitz criterion implies that no positive real roots of (6.6) exist.

Suppose that $D$ is negative. The discriminant of (6.6) is given by

It is clear that $D-2b <0$. On the other hand, we have

and

We conclude that $D+2b>0$, and then $\Delta_x<0$. Consequently, there is no real root of (6.6). In all cases, no $x: = \omega^2$ exits.Consequently, the equation (6.4) has no imaginary root $i\omega$. This means that no stability switches occurs. Hence, all roots of (6.4) have negative real parts. Then, the steady state $(\overline{S}, \overline{I}, \overline{u})$ is locally asymptotically stable.

7.   Uniform persistence

In all this section, we assume that

This condition implies that the steady state $(S^0, 0, u^0)$ is unstable. The instability of the disease-free steady state is not in contradiction with the existence of a nonnegative initial condition for which

We therefore need to prove the persistence of the component $I$ of infected individuals ^[46] (we will prove in fact the uniform persistence), that ensures survival of the infected individuals. First, we prove the following uniform weak persistence result.

Lemma 7.1. Assume that $\mathcal{R}_0>1$. Then, there exists a constant $\epsilon > 0$ such that for any initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+ \times \mathbb{R}^{+*} \times C^+$

Proof. Assume that

where

We can choose $\epsilon > 0$ sufficiently small such that

where

Remark that $S^0 > S^0_{\epsilon} >0$, for all $\epsilon >0$. With the choice of $\epsilon >0$ satisfying (7.2), we are going to show that (7.1) holds true. On the contrary, suppose that $\limsup_{t \to +\infty} I(t) \leq \epsilon$. Then, there exists a sufficiently large $T_\epsilon>0$ such that $I(t) \leq \epsilon$ for all $t \geq T_\epsilon$. Then, we have for all $t \geq T_\epsilon$

We will use a comparison principle. Then, we consider the following problem

The system (7.4) has a unique steady state $(S^0_{\epsilon}, u^0_{\epsilon})$, with $S^0_{\epsilon}$ given by (7.3) and

The corresponding $\mathcal{R}_0$ is given by (7.2). We will use the following result.

Lemma 7.2. Under the assumption (7.2) and for all initial condition $(S_0, \phi) \in \mathbb{R}^+ \times C^+$, we have $(S_\epsilon^+(t), u_\epsilon^+(t)) \to (S^0_\epsilon, u^0_\epsilon)$ as $t \to +\infty$.

Proof. See the proof of Theorem 8.1.

On the other hand, we can choose a sufficiently small $\tilde{\epsilon}>0$ such that

In the same time, we can choose a large $T'_{\tilde{\epsilon}}>T_{\epsilon}$ such that $S^+_\epsilon(t) > S^0_\epsilon-\tilde{\epsilon}$, for all $t \geq T'_{\tilde{\epsilon}}$. Furthermore, by comparison principle we have $S(t) \geq S^+_\epsilon(t) > S^0_\epsilon-\tilde{\epsilon}$, for all $t \geq T'_{\tilde{\epsilon}}$. Let $\xi>0$ be sufficiently small such that

Multiplying the equation of $I(t)$ by $e^{-\xi t}$ and integrating from $T'_{\tilde{\epsilon}}$ to $+\infty$, we obtain

This is a contradiction with (7.6). Hence (7.1) holds true.

In the next result, we use Lemma 7.1 to prove the following stronger result about the uniform strong persistence of infected individuals. The idea is based on the paper ^[47].

Theorem 7.3. Assume that $\mathcal{R}_0>1$. Then, there exists a constant $0<\epsilon' \leq \epsilon$, where $\epsilon$ is given by Lemma 7.1, such that for any initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+ \times \mathbb{R}^{+*} \times C^+$

Proof. From Lemma 7.1, we have $\limsup_{t \to +\infty} I(t) > \epsilon$. Then, there exists an increasing positive sequence $\{ \eta_k \}_{k = 0}^{+\infty}$, $\eta_k \to + \infty$ such that $I(\eta_k) > \epsilon$.

We prove Theorem 7.3 by contradiction. Suppose that for all $\epsilon' \in (0, \epsilon]$ there exists an initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+ \times \mathbb{R}^{+*} \times C^+$, such that

Then, there exist a positive increasing sequence $\{ t_k \}_{k = 0}^{+\infty}$ and a positive decreasing sequence $\{ \mu_k \}_{k = 0}^{+\infty}$ such that $t_k > \eta_k$, $\lim_{k\to+\infty} \mu_k = 0$ and

Then, $I(t_k) < \mu_k < \epsilon$. By the continuity of $I$, there exists a sequence $\{ \nu_k \}_{k = 0}^{+\infty}$, $\nu_k \in (\eta_k, t_k)$ such that

Let $\{ I_k \}_{k = 0}^{+\infty}$ and $\{ S_k \}_{k = 0}^{+\infty}$ be the sequences such that $I_k : = I(\nu_k) = \epsilon$ and $S_k : = S(\nu_k) \in \mathbb{R}^+$. The sequence $\{ I_k \}_{k = 0}^{+\infty}$ is constant and since the sequence $\{ S_k \}_{k = 0}^{+\infty}$ is uniformly bounded, it follows that there exist a convergent subsequences of $\{ I_k \}_{k = 0}^{+\infty}$ and $\{ S_k \}_{k = 0}^{+\infty}$ (denoted again $\{ I_k \}_{k = 0}^{+\infty}$ and $\{ S_k \}_{k = 0}^{+\infty}$) and $\rho \in \mathbb{R}^+$ such that $I_k = \epsilon$ and $\lim_{k \to +\infty} S_k = \rho$. Let consider the following problem

For each initial condition $\phi \in C: = C([-\tau, 0], \mathbb{R})$, the difference equation (7.9) has a unique solution $w$, which is continuous on $(0, +\infty)$. Let $\left\{ u_k \right\}_{k = 0}^{+\infty}$ be the functional sequence in $C$ defined by $u_k(\theta) : = w(\nu_k + \theta)$, $\theta\in [-\tau, 0]$, with $\nu_k > \tau$, for $k$ large enough (we make a translation of $k$ to have $\nu_k > \tau$ for all $k \in \mathbb{N}$). Then,

From the Proposition 2, we have the uniform boundedness of the sequence $\{ u_k \}_{k = 0}^{+\infty}$. For $\theta, \theta' \in [-\tau, 0]$, we have

where $n_k : = \lfloor \nu_k /\tau \rfloor$. Since $\phi$ is uniformly continuous on $[-\tau, 0]$, the inequality (7.10) implies that the sequence $\{ u_k \}_{k = 0}^{+\infty}$ is equicontinuous. Hence, it follows from the Ascoli-Arzela theorem that there exists $u^* \in C^+$ such that $\lim_{k \to +\infty} u_k = u^*$ (otherwise, we can choose a convergent subsequence).

Let consider now the solution of (3.2) with the initial condition $S_0 = \rho$, $I_0 = \epsilon$ and $\phi = u^* \in C^+$. We denote this solution by $(S^\infty, I^\infty, u^\infty)$. From Lemma 7.1, there exists $\sigma > 0$ and we can find $0<m<\epsilon$ (because $I^\infty$ is positive), such that

Next, we prove that we obtain a contradiction. For each $k \in \mathbb{N}$, we put $\tilde{I}_k(t) : = I(\nu_k+t)$, $t>0$. From (7.11), the continuity and the fact that

we have (recall that $\lim_{k\to+\infty} \mu_k = 0$), for $k$ large enough

On the other hand, for $\tilde{t}_k : = t_k - \nu_k$, we have from (7.7) and (7.8) that

We distinguish three cases, if $\sigma < \tilde{t}_k$, then the second inequality in (7.13) gives $\tilde{I}_k (\sigma) < \epsilon$ which contradicts the first inequality in (7.12). If $\sigma = \tilde{t}_k$, then the first inequality in (7.12) contradicts the first inequality in (7.13). If $\tilde{t}_k < \sigma$, then the second inequality in (7.12) gives $\tilde{I}_k (\tilde{t}_k)>\mu_k$ which contradicts the first inequality in (7.13). Consequently, there exists $\epsilon' \in (0, \epsilon]$ such that for any initial condition $(S_0, I_0, \phi) \in \mathbb{R}^+ \times \mathbb{R}^{+*} \times C^+$,

This completes the proof.

8.   Lyapunov functional and global asymptotic stability

In this section, we construct Lyapunov functionals to prove that the disease-free steady state is globally asymptotically stable when the basic reproduction number $\mathcal{R}_0<1$ and that the unique endemic steady state is globally asymptotically stable when $\mathcal{R}_0>1$.

8.1.   The global asymptotic stability of the disease-free steady state

In this part, we assume that $\mathcal{R}_0<1$ and we prove the global asymptotic stability of the disease-free steady state $(S^0, 0, u^0)$ of the system (3.2):

The solutions of this system satisfy, for all $t > 0$,

By the comparison principle, we have $S(t) \leq S^+(t)$ and $u(t) \leq u^+(t)$ for all $t>0$, where $(S^+, u^+)$ is the solution of the following problem

The system (8.1) has a unique steady state $(S^0, u^0)$, where $S^0$ and $u^0$ are the first and third components of the disease-free steady state of the system (3.2). In the next result, we show that the steady state $(S^0, u^0)$ of (8.1) is globally asymptotically stable.

Theorem 8.1. The unique steady state $(S^0, u^0)$ of (8.1) is globally asymptotically stable.

Proof. We put, for $t>0$,

Then, we get the linear differential-difference system

Let's consider the following Lyapunov functional

defined by

This functional satisfies, for $\nu_1(s) = s^2/2$ and $\nu_2(s) = ((1/2)+\tau\vartheta)s^2$, the inequalities

Moreover, the system (8.2) is input-to-state stable (see ^[44,48], for the definition and some properties of the notion of input-to-state stable). More precisely, there exist constants $C>0$ and $\sigma>0$ such that the solution $(\hat{S}, \hat{u})$ of (8.2) satisfies

The above estimation is an immediate consequence of (^[42], Theorem 3.5, page 275). By differentiating the function $t \mapsto V (\hat{S}(t), \hat{u}_t)$ along the solution $(\hat{S}, \hat{u})$ of the system (8.2), we obtain, for $t>0$

We want to find $\epsilon >0$ such that

We consider $\dfrac{d}{dt}V(\hat{S}, \hat{u}_t)+\epsilon \hat{S}^2(t)$ as a second order polynomial function of $\hat{S}(t)$. We compute the discriminant

and the expression

It is clear that we can choose $\epsilon>0$ small enough such that

We conclude that

Hence $(0, 0)$ is a globally asymptotically stable steady state of (8.2) (see, ^[44,48]). This completes the proof of Theorem 8.1.

Let $\epsilon >0$ and consider the set

Lemma 8.2. For any sufficiently small $\epsilon>0$, the subset $\Omega_{\epsilon}$ of $\mathbb{R}^+ \times \mathbb{R}^+ \times C^+$ is a global attractor for the system (3.2).

Proof. The solutions of (3.2) satisfy, for all $t > 0$,

By the comparison principle, we have $S(t) \leq S^+(t)$ and $u(t) \leq u^+(t)$ for all $t>0$, where $(S^+, u^+)$ is the solution of the system (8.1). Theorem 8.1 shows that $S^+(t) \to S^0$ and $u^+(t) \to u^0$ as $t \to +\infty$. This convergence implies that $\Omega_\epsilon$ is a global attractor for the system (3.2) in $\mathbb{R}^+ \times \mathbb{R}^+ \times C^+$. This completes the proof.

Thanks to Lemma 8.2, we can restrict the global stability analysis of the disease-free steady state of (3.2) to the set $\Omega_\epsilon$.

Theorem 8.3. Assume that $\mathcal{R}_0<1$. Then, the disease-free steady state $(S^0, 0, u^0)$ of (3.2) is globally asymptotically stable.

Proof. It suffices to consider the solutions in $\Omega_\epsilon$ for any sufficiently small $\epsilon>0$. We then have, for $t>0$,

Since ${R}_0 < 1$, we can choose $\epsilon>0$ such that the right-hand side of the above inequality is negative. This implies that $\lim_{t \to +\infty} I(t) = 0$.

From the above result, we see that for any $\epsilon > 0$, there exists a $T_\epsilon>0$ such that $I(t) \leq \epsilon$ for all $t \geq T_\epsilon$. We then have, for $t > T_\epsilon$

Then, we have $S(t) \geq S_\epsilon (t)$ and $u(t) \geq u_\epsilon (t)$ for all $t \geq T_\epsilon$, where $(S_\epsilon, u_\epsilon)$ is the solution of the following problem

As in the proof of Lemma 7.1 and Theorem 8.1, we can show that $S_\epsilon (t) \to S^0_\epsilon$ and $u_\epsilon (t) \to u^0_\epsilon$ as $t \to +\infty$, where $(S^0_\epsilon, u^0_\epsilon)$ is the steady state of (8.3), given by (7.3) and (7.5). Then, there exists a $\tilde{T}_\epsilon>T_\epsilon>0$ such that, for $t \geq \tilde{T}_\epsilon$,

Since $\epsilon > 0$ is arbitrary, $S^0_\epsilon \to S^0$ and $u^0_\epsilon \to u^0$ as $\epsilon \to 0$, we have that $\lim_{t \to +\infty} S(t) = S^0$ and $\lim_{t \to +\infty} u(t) = u^0$. Recalling from Theorem 6.1 that $(S^0, 0, u^0)$ is locally asymptotically stable. Then, it is globally asymptotically stable. This completes the proof.

8.2.   The global asymptotic stability of the endemic steady state

In this section, we assume that

Let $(\overline{S}, \overline{I}, \overline{u})$ be the endemic steady state. This equilibrium satisfies $\overline{S}>0$, $\overline{I}>0$ and $\overline{u}>0$. Let $\tilde{S}(t) : = S(t)-\overline{S}$ and $\tilde{u}(t) : = u(t) - \overline{u}$. Then, the system (3.2) transforms to

Note that $\beta \overline{S} = \mu + \gamma_I$.

Theorem 8.4. Assume that $\mathcal{R}_0>1$. Then, the steady state $(\overline{S}, \overline{I}, \overline{u})$ is globally asymptotically stable.

Proof. Let consider the following Lyapunov function

defined by

where

We point that the function

satisfies $G(I_0) \geq 0$ for all $I_0>0$ and $G(I_0) = 0$ if and only if $I_0 = \overline{I}$. This observation means that $W(S_0, I_0, u_0) = 0$ if and only if $(S_0, I_0, u_0, ) = (0, \overline{I}, 0)$.

We set

Then, the derivative of $t \mapsto W (\tilde{S}(t), I(t), \tilde{u}_t)$ along the solution trajectory is given by

Since $c>0$, we obtain

with

Then, the function $t \mapsto W(\tilde{S}(t), I(t), \tilde{u}_t)$ is nonincreasing and we have

Furthermore, by integration (8.5), we get

The both sides of the inequality (8.6) are nondecreasing functions. Then, the limits exist and satisfy

As the function $\tilde{S}'(t)$ is uniformly bounded $\tilde{S}(t)$ is uniformly continuous. Then, the Barbalat's Lemma ^[49] applied to the function $t \mapsto \int_0^t \tilde{S}^2(s)ds$, shows that

Using ^[45], Lemma 3.5, we obtain

Then, the expression of the function $W$ implies that

Furthermore, the function $\tilde{S}(t)$ is bounded and differentiable, then the fluctuations Lemma implies that there exists a sequence $t_k \to +\infty$ such that $\lim\limits_{k \to +\infty}\tilde{S}'(t_k) = 0$. Then, the first equation of (8.4), implies that $\lim\limits_{k \to +\infty}I(t_k) = \overline{I}$. The continuity of the function $G$ gives $\lim\limits_{k \to +\infty}G(I(t_k)) = G(\overline{I}) = 0$. Then, $W^* = 0$. From the properties of the function $G$, we conclude that

This prove the global asymptotic stability of $(\overline{S}, \overline{I}, \overline{u})$.

9.   Discussion

In this work, an epidemiological SIR model with a class of age structured temporary protected individuals is presented. A coupled system of differential-difference equations with delay is derived from this model by using the method of the characteristics. We followed the same idea as in our recent work ^[44]. The model presents two steady states: disease-free and endemic. The condition for the existence and the asymptotic behavior of these steady states (local and global asymptotic stability) are discussed. The global asymptotic stability was proved by using Lyapunov functions. In summary, if $\mathcal{R}_0<1$ the disease-free steady state is globally asymptotic stable, otherwise if $\mathcal{R}_0>1$, which ensures also its existence, the endemic steady state is globally asymptotic stable. The threshold $\mathcal{R}_0$ is the well-known "basic reproduction number of the disease". Comparing the new $\mathcal{R}_0$ with the one obtained by the classical Kermack-McKendrick model, the new parameters $h$, $\alpha$ and $\tau$ decrease the $\mathcal{R}_0$ value, which means that the disease is easier to control when protection, through vaccination or drugs, is take into account (Figures 2 and 3). This emphasize the importance of the compliance with the adopted control strategies. The parameters $h$, $\tau$ and $\alpha$ are, respectively, the protection rate, the duration of the temporary protection phase, and the specific protection rate at the end of the previous period of protection. Keeping one of the three parameters fixed and varying the others, we can observe that if one parameter decrease, the other has to increase to keep the transmission of the disease under control. In all cases, control is achieved if $\mathcal{R}_0(m)<1$, with $m = h$, $\alpha$ or $\tau$, (Figures 2 and 3). Finally, we showed that when short immunity is considered, the fraction of individuals that update their vaccine at the end of their period of protection has an important impact on disease dynamics.

Acknowledgments

CPF acknowledges Grant 18/24058-1 from São Paulo Research Foundation (FAPESP). Part of this research was developed during the visit of Dr Adimy to IBB/UNESP by PROPG 02/2018. We also want to thank the anonymous referees for their careful reading that helped us to improve the manuscript.

Conflict of interest

The authors declare no conflicts of interest.