A novel discrete-time COVID-19 epidemic model including the compartment of vaccinated individuals

1.   Introduction

Epidemic models have received considerable attention over the last years ^[1]. They can be represented by either dynamic continuous-time or discrete-time systems, i.e., they are described by either differential or difference equations ^[1]. Recently, several epidemic models have been proposed with the aim to understand, describe and control the spread of the COVID-19 pandemic. These models analyse the evolution of the disease over time by dividing the communities into some compartments, which mainly include the susceptible class (S), the exposed class (E), the infected class (I) and the removed class (R). For example, in ^[2] a SIR model is utilized for estimating the infectivity and recovery rates from COVID-19 real data. Then, the estimated rates have been exploited to analyse the evolution of the pandemic over time. In ^[3] another SIR model is presented, with the aim to describe the spread of the COVID-19 epidemic in Wuhan. The model is updated with real-time input data in order to derive clinical parameters that could support public officials in decision-makings. In ^[4] a dynamic continuous-time model of the COVID-19 pandemic is introduced. By using phase portraits and time-series plots, the authors of reference [4] highlight the chaotic dynamics of their epidemic model. In ^[5] a nonlinear dynamic SIR model, which includes the effect of social distancing, is proposed to carefully analyse the spread of the COVID-19 pandemic. In ^[6] an accurate model that includes eight compartments (i.e., susceptible (S), infected (I), diagnosed (D), ailing (A), recognized (R), threatened (T), healed (H) and extinct (E)) is presented. The model is called SIDARTHE and introduces the distinction between diagnosed/non-diagnosed individuals, being the former isolated and, consequently, less likely to spread the infection.In ^[7] a COVID-19 model is illustrated, which includes the effects of media on public awareness as well as a removal rate based on the hospital-bed population ratio. The approaches to modelling the spread of the COVID-19 pandemic have become to change, at the beginning of the year 2021, by virtue of the introduction of the vaccines. Namely, a number of papers have been published, which take into account (in different ways) the role of vaccination for controlling the pandemic. Most of these models are described by differential equations, whereas very few models involve difference equations. Referring to continuous-time models, in [8] a novel SEIRS system that includes the vaccine rate is illustrated. The model clearly shows that, when the vaccine rate increases, the infection decreases and the recovered population increases over time. In ^[9] a continuous-time model based on eight state variables is presented, where the vaccination is considered as preventive action, rather than as a system variable. In particular, the approach in ^[9] highlights that vaccination is a key solution to eliminate the coronavirus among healthcare workers. In ^[10] an extended Kalman filter is developed to estimate the state variables of the model. In particular, the spread of COVID-19 is modelled by a set of ordinary differential equations where vaccinations represent the control input signals. In ^[11] a nonlinear robust control policy is exploited to analyse the role of vaccinations on the spread of the COVID-19 pandemic. The proposed continuous-time dynamic model includes eight state variables, i.e., susceptible, exposed, infected, quarantined, hospitalized, recovered, deceased and insusceptible populations. In ^[12] an innovative SEIR model that takes into account the effect of vaccination on the spread of COVID-19 is illustrated. Some properties of the model, including output-reach ability and non-negativity of the solution, are carefully analysed. In ^[13] a continuous-time SEIR model, which includes both a feedback vaccination law and an antiviral treatment control law, is illustrated. In particular, the existence of a unique disease-free equilibrium point is proved, along with the attainability of the endemic equilibrium point. In ^[14] a SEIR-type epidemic model with time delay and vaccination control is developed. The vaccination strategy is expressed as a state delayed feedback, which is related to the current and previous state of the epidemic model.Reference [15] is devoted to modelling and predicting COVID-19 confirmed cases through a multiple linear regression. In particular, a novel trial is illustrated, which combines both growth rates and vaccination rates in modelling the spread of the COVID-19 pandemic.

It is worth noting that discrete-time epidemic models present the advantage, over the continuous-time model, to allows a direct acquisition of recorded data, given that the sampling period (i.e., the unity), may correspond to one day or one week, depending on the available data ^[16]. In spite of these considerations, just a few papers have been published so far (to the best of the authors' knowledge) regarding discrete-time models that take into account the effect of vaccination on the COVID-19 pandemic. Namely, in ^[17] a novel discrete SEIR COVID-19 model is developed. In particular, a feedback vaccination control law is introduced to control the \textquotedblleft susceptible class\textquotedblright, in ^[18] the authors present a separate new epidemiological model for susceptible to have recovered infection (SEIR) subject to a two-dose feedback vaccination efforts, in ^[19] the authors developed a susceptible- infected- recovered- vaccine-1- vaccine-2- death model to analyze the behavior of the epidemic in Japan. Besides the fact that the topic of discrete-time models including vaccination is almost unexplored, it should be also noted that most the COVID-19 models (both continuous-time and discrete-time) published to date consider the vaccination as a control law or as a preventive action. Based on these considerations, this paper makes a contribution to the topic of discrete mathematical modelling of epidemics by introducing a new discrete-time COVID-19 model, which includes the number of vaccinated individuals as an additional state variable describing the system dynamics and it has the same characteristics as the Susceptible class. The paper shows the existence of two fixed points at most in the proposed compartment model described by difference equations, i.e., the disease-free fixed point and the epidemic fixed point. By considering both the forward difference system and the backward difference system, some stability analyses of the disease-free fixed point are carried out. In particular, for the backward difference system a theorem is proved, which gives a condition for the disappearance of the pandemic when an inequality involving some epidemic parameters is satisfied. The paper is organized as follows. In Section 2 a new discrete-time compartment model for describing the spread of the COVID-19 pandemic is illustrated. The discrete system involves five state variables, i.e., the Susceptible class S, the Recovered class R, the Infection class I, the Infection dangerous class Id and the Vaccinated class V.In Section 3it is shown that the system possesses two fixed points, i.e., the disease-free fixed point and the epidemic fixed point. Moreover, the basic reproduction number is computed. In Section 4 a stability analysis for the forward difference system is conducted, whereas in Section 5 a similar analysis is carried out for the backward difference system. In particular, in Section 5 a novel theorem is proved, which assures the global stability of the disease-free fixed point, indicating that the pandemic disappears provided that an inequality involving some epidemic parameters is satisfied. Finally, in Section 6 simulation results of the conceived discrete COVID-19 compartment model are carried out, along with comparisons regarding the performances of both the forward difference system and the backward difference system.

2.   A novel discrete-time compartment model including the vaccinated class

For the proposed model, the study population (N) is divided into two main classes: the class of individuals exposed to infection and the class of individuals infected. Each of the previous classes is also divided into secondary class, so the class of people who are exposed to infection is divided into three sub-class: people who are exposed to infection and were not previously infected and did not receive vaccination (S), people who were previously infected and recovered from the disease and are at risk of being infected again (R), and the class of vaccinated people against the epidemic (V). As for the class of infected persons, it is divided into two secondary classes: The people with good immunity and for whom infection does not pose a great risk and suppose that their ratio in the society is $\lambda, (\lambda\leq1)$, so the class of infected people from this group are going to I. The infected persons for whom the infection is dangerous and consists of the elderly, pregnant women and people with chronic diseases (whose their ratio in the society is $(1-\lambda)) I_d$. The transition between classes is described in the following diagram: So the mathematical model that we will be interested in, which is given in detail in ^[20] as follow:

The proposed model's flowchart and parameters descriptors are well explained in (2.2).

Where $r'_{i} = \frac{p_{i}k}{N}, i = 1, 2, 3$, $k$ is the average numbers of contacts per capita (per unit of time), $p_{i}$ is the probabilities of contagion ($p_{1} > p_{2} > p_{3}$) and $N$ is the total population (It can be considered as the maximum value of the population)

Adding up the equations given in (2.1), we find

In addition, the following initial conditions take into consideration:

In real life, the statistics are discrete, so the discrete system is closer to modeling the spread of Corona virus, and in light of this we will use the following approximation:

System 2.1 becomes:

We multiply both sides by $h$, and put

and

where $n\in \mathbb{N}, X(t) = \left(S\left(t\right), R\left(t\right), V\left(t\right), I\left(t\right), I_{d}\left(t\right) \right) ^{t}.$ And we also put

We finally get the following system:

where

and where $\Delta$ is the forward difference operator ($\Delta X\left(n\right) = X\left(n+1\right) -X\left(n\right)$).

We can also use the following approximation:

and by following similar steps to the above, we get:

where $\nabla$ is the backward difference operator ($\nabla X\left(n\right) = X\left(n\right) -X\left(n-1\right)$).

3.   Fixed points and basic reproduction number

3.1.   Fixed points

To study the dynamics of the systems 2.5 and 2.7, one must first find the fixed points, and to find the fixed points must be solved the following nonlinear system:

The previous equation has the point $E_{0} = \left(\frac{\Omega }{\left(\mu +\upsilon \right) }, 0, \frac{\upsilon \Omega }{\mu \left(\mu +\upsilon \right) }, 0, 0\right)$ as a solution. It can be seen that the disease at this point is non-exist and therefore it is called the disease-free fixed point which we will focus on studying later.

If we suppose that $\left(I^{\ast }+I_{d}^{\ast }\right) \neq 0,$ So we will get:

This system is complex and difficult to solve in the abstract case. Even if we can solve it, studying the stability of this fixed point contains many obstacles. In general, this point is called the endemic equilibrium point $E^{\ast } = \left(S^{\ast }, R^{\ast }, V^{\ast }, I^{\ast }, I_{d}^{\ast }\right)$ and its existence and stability can be studied numerically.

3.2.   The basic reproduction number

The basic reproduction number is a very important number for studying the behavior of epidemic systems. According to ^[21] the basic reproduction number is defined as the spectral radius of the next generation matrix $FV^{-1}$, where $F$ and $V$ are the Jacobian matrices for $\mathcal{F}$ and $\mathcal{V}$ respectively at the disease-free fixed point, and where $\mathcal{F}$ the rate of appearance of new infections:

and $\mathcal{V}$ the disappearance rate of infections:

The Jacobite matrices are:

and

So the next generation matrix is:

According to it the basic reproductive number is given by:

4.   Stability analysis for the forward difference system

4.1.   Existence and uniqueness

System 2.5 can be written as:

Thus, it is defined by a regression relationship, we can notice that knowing each term allows us to know the term after it, so existence and uniqueness are trivial in this case. It cannot be shown that the solution is positive for system 2.5. In fact, the solution to system 2.5 is not always positive even when the initial conditions are positive.

4.2.   Stability analysis

We are now studying the stability of the disease-free fixed point, we mention the following theorem that we need in the study of stability:

Theorem 1. ^[22] System 4.1 is asymptotically stable if each eigenvalue $\lambda_*$ of the Jacobite matrix $J_*$ of $X\left(n\right) +F\left(X\left(n\right) \right)$ at the fixed point is satisfied

Using the above theorem we result:

Theorem 2. Suppose that $R_{0} < 1,$ if

where

Then the disease-free fixed point $E_{0}$ of system 2.5 is locallyasymptotically stable.

Proof. System 2.5 is equivalent to system 4.1, and the Jacobian matrix of $X\left(n\right) +F\left(X\left(n\right) \right)$ at $E_{0}$ is:

where $I$ is the identity matrix and J is the Jacobian matrix of $F(X(n))$ at $E_0$:

The characteristic polynomial of $J_*$ is:

So the matrix $J_*$ has $\lambda_1 = 1-\mu$ as normal eigenvalue and $\lambda_2 = 1-\left(\mu +\upsilon \right)$ as a double eigenvalue, so according to conditions 4.3 and because all parameters are positives ($\mu +\upsilon < 2$ and $\mu, \upsilon > 0$), we get:

so $\lambda_1$ and $\lambda_2$ satisfy the condition of Theorem 1. Stayed the roots of the polynomial:

If $R_{0} < 1$ then $A, B > 0$, we have from the conditions 4.3 ($A^{2} < 4B$):

So $B-A+1 > \frac{\left(A-2\right)^{2}}{4} > 0$, and the roots are complex and conjugated, we not them $\lambda$ and $\bar{\lambda}$ and we have

(because $\frac{\Omega \delta \left(\upsilon r_{3}+\lambda \mu r_{1}\right) }{\mu \left(\mu +\rho \right) \left(\mu +\upsilon \right)} > 0$). So all the eigenvalues of the characteristic equation of $J_*$ satisfy the condition of Theorem 1. Accordingly the disease-free fixed point $E_{0}$ of system 2.5 is locally asymptotically stable.

The structure of this system does not allow us to study the global stability, therefore we will move to the backward difference system.

5.   Stability analysis for the backward difference system

5.1.   Existence and uniqueness

System 2.7 can be written as:

So, to calculate each term, the term before it must belong to the image of the function $x\rightarrow x-F(x)$, and this makes the existence not trivial like the previous case, for this purpose, we give the following result:

Theorem 3. For any $n^{\ast }\in \mathbb{N}$, there exists $\delta ^{\ast } > 0$, so that the solution of system 2.7 (for $n\leq n^{\ast }$) is exists when $\left\Vert x\left(0\right) -E_{0}\right\Vert < \delta ^{\ast }.$ Also, this solution is positive when the initialcondition is positive.

Proof. System 2.7 is equivalent to system 5.1, and from the system 4.5:

so the solution is the fixed point of $T\left(X\left(n\right) \right)$, where:

in the space $\ell ^{\infty }\left(\left[ 0, n^{\ast }\right] \cap \mathbb{N} \right)$, with its norm $\left(\left\Vert {.}\right\Vert _{\infty }\right)$, so for $x_{1}, x_{2}\in \ell ^{\infty }\left(\left[ 0, n^{\ast }\right] \cap \mathbb{N} \right)$, we have:

since $F$ continues in the neighborhood of $E_{0},$ we have $\forall \epsilon > 0, \exists \delta > 0:x_{1}\left(i\right), x_{2}\left(i\right) \in B\left(E_{0}, \delta \right)$ for $i\leq n\Longrightarrow \left\Vert F\left(x_{1}\left(i\right) \right) -F\left(x_{2}\left(i\right) \right) \right\Vert < \epsilon \left\Vert x_{1}\left(i\right) -x_{2}\left(i\right) \right\Vert.$ Take $\epsilon ^{\ast } > 0$ check that $n\epsilon ^{\ast } < 1$ for $n\in \left[ 0, n^{\ast }\right] \cap \mathbb{N}$, then $\exists \delta ^{\ast } > 0,$ check that:

So $T$ is contraction map in $\ell ^{\infty }\left(\left[ 0, n^{\ast }\right] \cap \mathbb{N} \right)$, and from it, by Banach fixed point theorem there exists one solution for the system 5.1 for $n\leq n^{\ast }$.

To prove the positivity of the solution we assume the opposite. Let $S$ be the first component to become negative at $n_{1}\in \mathbb{N},$ (i.e $R\left(n_{1}\right), V\left(n_{1}\right), I\left(n_{1}\right), I_{d}\left(n_{1}\right) \geq 0$ and $S\left(n_{1}\right) \leq 0$), then:

this means

a contradiction so $S\left(n\right) \geq 0$ for $n\leq n^{\ast }.$ If $R$ the first component to become negative at $n_{2}\in \mathbb{N},$ (i.e $S\left(n_{2}\right), V\left(n_{2}\right), I\left(n_{2}\right), I_{d}\left(n_{2}\right) \geq 0$ and $R\left(n_{2}\right) \leq 0$), then:

this means

a contradiction so $R\left(n\right) \geq 0$ for $n\leq n^{\ast }.$ If $V$ the first component to become negative at $n_{3}\in \mathbb{N},$ (i.e $S\left(n_{3}\right), R\left(n_{3}\right), I\left(n_{3}\right), I_{d}\left(n_{3}\right) \geq 0$ and $V\left(n_{3}\right) \leq 0$) then:

this means

a contradiction so $V\left(n\right) \geq 0$ for $n\leq n^{\ast }.$ If $I$ the first component to become negative at $n_{4}\in \mathbb{N},$ (i.e $S\left(n_{4}\right), R\left(n_{4}\right), V\left(n_{4}\right), I_{d}\left(n_{4}\right) \geq 0$ and $I\left(n_{4}\right) \leq 0$) then:

so

because we know that:

a contradiction, so $I\left(n\right) \geq 0$ for $n\leq n^{\ast }.$ If $I_{d}$ the first component to become negative at $n_{5}\in \mathbb{N},$ (i.e $S\left(n_{5}\right), R\left(n_{5}\right), V\left(n_{5}\right), I\left(n_{5}\right) \geq 0$ and $I_{d}\left(n_{5}\right) \leq 0$) then:

so

because we know that:

a contradiction, so $I\left(n\right) \geq 0$ for $n\leq n^{\ast }.$In the end, we conclude that:

5.2.   Invariant region

Theorem 4. The System 2.7 have

as invariant region, where $\mathbb{R} _{+}^{5} = \left\{ \left(x_{1}, x_{2}, x_{3}, x_{4}, x_{5}\right) \in \mathbb{R} ^{5}\; {\mathit{and}}\; x_{i}\geq 0\; {\mathit{for}}\; i = 1..5\right\}.$

Proof. Adding the equations of system 2.7, we get:

Since $I_{d}$ is positive, we get:

where

and from it we find:

Let $N(0)\leq \frac{\Omega }{\mu },$ and suppose that $N(n)\leq \frac{\Omega }{\mu }$ for a natural number $n,$ we get:

then by induction for all $n$ when the solution exists:

Therefore, the solution belongs to the invariant region:

5.3.   Local stability

We are now studying the stability of the disease-free fixed point, so we present the following result:

Theorem 5. Suppose that $R_{0} < 1$. Then the disease-free fixed point $E_{0}$ of thesystem 2.7 is locally asymptotically stable.

Proof. Locally, system 2.7 behaves the same as that of:

$\Leftrightarrow$

since

when $R_{0} < 1,$ so

According to Theorem 1, the system 5.3 is stable when each eigenvalue $\lambda_{*}$ of $\left(I-J\right)^{-1}$ is satisfied

its meaning that each eigenvalue $\lambda^{-1}_{*}$ of $(I-J)$ is satisfied

and whereas

where $\lambda$ an eigenvalue for $J$, the condition of stability becomes:

The characteristic polynomial of $J$ is :

So the matrix $J$ has $\lambda_1 = -\mu$, as normal eigenvalue and $\lambda_2 = -\left(\mu +\upsilon \right)$, as a double eigenvalue, and from it:

and

Stayed the roots of the polynomial:

if $R_{0} < 1,$ then $A, B > 0$, so according to Routh-Hurwitz criterion, the roots of this polynomial called $\lambda_3$ and $\lambda_3$ have a negative real part, so:

and

where $Re(\lambda)$ is the real part of $\lambda$. Therefore, all eigenvalues of matrix $J$ satisfying the condition 5.6, so the disease-free fixed point is locally asymptotically stable.

5.4.   A condition for the disappearance of the disease

The main purpose of studying epidemic models is to find applicable conditions in the real life and contribute to the disappearance of the epidemic, so we will now give a condition related to reducing friction to get rid of the epidemic.

Theorem 6. If

then the disease will disappear.

Proof. Adding the last two equations in the system 2.7 together, we find that the infected class is described by the following equation:

Since $I$ is positive and $r_{i} = \frac{p_{i}k}{N}, i = 1, 2, 3, \left(p_{3}\leq p_{2}\leq p_{1}\right) :$

$\Leftrightarrow$

$\Leftrightarrow$

$\Leftrightarrow$

We note that if $\left(1-p_{1}k+\left(\mu +\rho \right) \right) ^{-1} < 1,$ then $\left(I+I_{d}\right) (t)\rightarrow 0$ when $t\rightarrow \infty.$

6.   Numerical simulations and comparisons

In this section we will compare systems 2.5 and 2.7, note that the solution of the system 2.5 has always existed, but it is not always positive and also it cannot be proven that it is bounded. In system 2.7 the solution does not always exist, but if it is exists by taking a positive initial condition the solution is positive, and it also belongs to an invariant region. As for stability, the system 2.5 and even when $R_{0} < 1$, stability is not guaranteed except by adding other conditions, contrary to system 2.7 the only condition for local stability is $R_{0} < 1$ and therefore it is the one who maintains the characteristics of the continuous system 2.1. In the system 2.5 We couldn't find a condition for global stability, which was so easily formulated in the system 2.7.

In the following, we will do some numerical simulations with the same initial value and the same parameters, and note the essential differences between the two systems 2.5 and 2.7. Therefore, we choose a divided population as follows:

And we take the parameters as:

We get the simulation as:

We notice that there is no difference in the previous example. In fact, essential differences sometimes do not appear in short periods of time. We now take the same system as before and make

we get

There is no difference at the beginning, but when we take long periods of time, we notice the differences. We have increased the rate of spread of the disease. We note that the disease is spreading quickly. In this case, we notice an increase in $I$ and a decrease in $S$ for the forward difference system, but the opposite happens to the backward difference system.

7.   Conclusion

Referring to the epidemic models that take into account the effect of vaccination on the COVID-19 pandemic, all the models (both continuous-time and discrete-time) published to date consider the vaccination as a control law or as a preventive action. No paper includes the vaccination as an additional compartment, to be added to the remaining compartments when describing the system dynamics. This manuscript has presented a novel discrete-time COVID-19 model, which includes the number of vaccinated individuals as an additional state variable in the system equations. The paper has shown that the proposed compartment model, described by difference equations, has a disease-free fixed point and an epidemic fixed point. By considering both the forward difference system and the backward difference system, some stability analyses of the disease-free fixed point have been carried out. In particular, for the backward difference system a novel theorem has been proved, which has provided a condition for the disappearance of the pandemic when an inequality involving some epidemic parameters is satisfied.This result represents a remarkable finding of the proposed approach, which may help decision-makers to better understand the epidemiological behaviour of the COVID-19 over time.Finally, numerical simulations of the proposed discrete model have been carried out, along with comparisons regarding the performances of both the forward difference system and the backward difference system.

Acknowledgments

This research has been funded by Scientific Research Deanship at University of Ha'il - Saudi Arabia through project number RG-21067.

Conflict of interest

The authors declare no conflict of interest.