A review of commonly used

compartmental models in epidemiology

Mohamed Mehdaoui 111E-mail<EMAIL_ADDRESS>
_Faculty of Sciences and Technics of Settat, B.P. 577, Settat 26000, Morocco_

Abstract

In order to model an epidemic, different approaches can be adopted. Mainly,
the deterministic approach and the stochastic one. Recently, a large amount of
literature has been published using the two approaches. The aim of this paper
is to illustrate the usual framework used for commonly adopted compartmental
models in epidemiology and introduce variant analytic and numerical tools that
interfere on each one of those models, as well as the general related types of
existing, ongoing and future possible contributions.

_Keywords_ : Dynamical systems; Epidemiological modeling; Differential
equations; Optimal control; Stability theory.

##  Introduction

Epidemiology is a scientific discipline that studies the distribution (who,
when, and where) patterns and determinants of health and disease conditions in
defined populations. It is a cornerstone of public health, and shapes policy
decisions and evidence-based practice by identifying risk factors for disease
and targets for preventive healthcare. The mathematical modeling of an
epidemic allows to describe the spread of an infectious disease and predicts
its evolution. This knowledge is crucial to determine politics and measures
that can control the fast spread of the disease. Mathematical models allow to
estimate the effects of the epidemic, for instance the number of total people
that are infected and the duration of the epidemic, moreover, the effects of
prevention measures such as social distancing, vaccination or confinement. The
study of the dynamics of an epidemic requires the knowledge of many structural
variables: the absolute time, the age of the hosts, the duration of the
infection, the immunity status. etc. Time is essential because it describes
the continuous variation of a type of systems called dynamical systems. In the
simple case, the system is formalized in terms of ordinary differential
equations (ODEs). In some cases, to illustrate the ”memory” effect of the
dynamics, the system is formalized by fractional differential equations
(FDEs), while delay differential equations (DDEs) interfere to capture the
incubation period for instance. As for the partial differential equations
(PDEs), they interfere when the spatial factor is included into the model.
Deterministic dynamical systems, whether they’re based on differential
equations or partial differential equations, are easy to simulate. Their
smoothness allows to explore a diversity of scenarios. Moreover, many
theoretical and numerical tools allow to exploit these systems in a way to
obtain a formula that expresses the reproduction number in terms of the
model’s parameters. However, because the spread of a disease is a random
process, sometimes It’s pertinent to write the model in a probabilistic point
of view. Moreover, the deterministic epidemiological models are formed under
the assumption of a large-sized population. When it’s not the case, the
interactions between individuals are not uniform but have a random effect. In
this case, the model is governed by stochastic differential equations (SDEs).

Before choosing one of the two approaches of epidemiological modeling, the
restriction to compartmental models is usually done. The logic behind this
type of models is to express the dynamics of the population by dividing this
latter into as many compartments as the clinical states. Usually, these
compartments correspond to the population of susceptible (S), those who are
infected (I) and those who don’t participate to the epidemic anymore (R for
recovered or removed individuals who are either dead or immune for a period of
time). Variant considerations can be added to make the model more realistic
and complex (exposed individuals but not yet infectious, infectious
individuals but not yet detected, etc.). The transition from a clinical state
to another is described by an incidence function.

Fig. 1: The relation between different approaches of compartmental modeling

##  Common biological/mathematical background

The compartmental mathematical models reviewed in this paper are expressed in
either of the following forms
$\begin{cases}\dfrac{du(t)}{dt}=f(u(t))\;t>0,\\\ \text{+ positive initial
conditions}.\end{cases}$ (1)
$\begin{cases}du(t)=f(u(t))dt+A(u(t))\;dW(t)\;t>0,\\\ \text{+ positive initial
conditions}.\end{cases}$ (2)
$\begin{cases}D^{\alpha}u(t)=f(u(t))\;t>0,\\\ \text{+ positive initial
conditions}.\end{cases}\;\;\;\;\;$ (3)
$\begin{cases}\dfrac{du(t)}{dt}=f(u_{1}(t-\tau_{1}),\cdots,u_{m}(t-\tau_{m}))\;\;t>0\\\
\text{+ positive initial functions}.\end{cases}\;\;\;\;\;$ (4)

$\begin{cases}-\underbrace{\mathcal{D}\Delta u(x,t)}_{\text{pointwise
product}}+\dfrac{\partial u(x,t)}{\partial t}=f(u(x,t))\;x\in\Omega,t>0,\\\
\partial_{\eta}u_{1}(x,t)=\cdots=\partial_{\eta}u_{m}(x,t)=0\;x\in\partial\Omega,t>0,\\\
\text{+ positive initial conditions}.\end{cases}$ (5)

such that

  * •

$u:=(u_{1},\cdots,u_{m})$ with $m\in\mathbb{N}^{*}$ representing the number of
compartments.

  * •

$f:(0,\infty)^{m}\longrightarrow\mathbb{R}^{m}$ is a continuously
differentiable function.

  * •

$\dfrac{d}{dt}$ is the classical time derivative which will be denoted by ′
for $\eqref{eq:ode}$ and $\eqref{eq:dde1}.$

  * •

$\dfrac{\partial}{\partial t}$ is the partial derivative with respect to time.

  * •

$\Omega\subset\mathbb{R}^{p}$ ($p\in\\{1,2\\}$) is an open bounded set,
$\partial\Omega$ its boundary and $\partial_{\eta}$ is the outward normal
derivative.

  * •

$\Delta u:=(\Delta u_{1},\cdots,\Delta u_{m})$ where
$\Delta:=\dfrac{\partial^{2}}{\partial
x_{1}^{2}}+\dfrac{\partial^{2}}{\partial x_{2}^{2}}$ is the Laplace operator.

  * •

$\mathcal{D}:=(d_{1},\cdots,d_{m})$ is a vector of strictly positive diffusion
coefficients.

  * •

$W_{t}$ is an $n$-dimensional Wiener process (will also be denoted by
$B_{t}$).

  * •

$A(u(t))\in\mathcal{M}_{m\times n}(\mathbb{R})$ is a $m\times n$ matrix with
locally Lipschitz coefficients.

  * •

$D^{\alpha}$ is some fractional time derivative of order $\alpha\in(0,1)$ (see
de Oliveira & Tenreiro Machado (2014)).

  * •

$\tau_{i}\geq 0\;\forall i\in\\{1,\cdots,m\\}$ are the discrete time delays.

For more terminology on each of the previous differential equations, we refer
to Hsieh (2013); Oksendal (2013); Hale (1977); Evans (2016).
All the considered compartmental models share a common biological/mathematical
framework as well as some steps followed in their mathematical and qualitative
analysis such as:

#### The basic reproduction number $\mathcal{R}_{0}$:

One quantity that plays a crucial role on all of the approaches of
epidemiological modeling illustrated in Figure 1, is the basic reproduction
number, which is biologically speaking, defined as the average number of
secondary infected caused by a single infected individual on the period of his
contagiousness. One question that intuitively arises is: How this quantity can
be expressed mathematically based on the chosen compartmental model? The
answer was given by Driesschea & Watmough (2002). The problem was reduced to a
computation of the spectral radius of the product of two matrices $F$ and
$V^{-1}$, such that $F$ represents the rate of appearance of new infections
and $V$ is the rate of transfer of individuals into infected compartments by
all other means.

#### The study of mathematical and biological well-posedness:

The mathematical well-posedeness of a compartmental model relies on proving
the existence and uniqueness of the solution. Depending on the type of the
differential equation that governs the model, there exist many theorems that
can be applied to achieve this purpose. For instance, the Cauchy-Lipschitz
theorem for (1) (see Hsieh (2013) pp. 1-8). For (2), this is established by
means of the Itô existence theorem (see Oksendal (2013) pp. 65-84). For (4),
(see Hale (1977) pp. 13-23). Finally for (5), (see Vrabie (2003) pp. 249-261).
Once the mathematical well-posedeness is achieved, the biological well-
posedness relies on proving the boundedness of the unique solution as well as
its positiveness. Given a positive initial condition, for
$\eqref{eq:ode},\eqref{eq:sde},\eqref{eq:fde}$ and $\eqref{eq:dde1},$ this is
usually proven by a contradiction argument. For $\eqref{eq:pde}$, the maximum
principle approach is adopted (see Protter & Weinberger (1984)).

#### The study of different disease equilibria:

In dynamical systems, it is of great importance to study the solutions that
don’t change in time. Such solutions are called equilibria. For compartmental
models, two equilibria are crucial: the Free Disease Equilibrium (FDE) and the
Endemic Equilibrium (EE). The first equilibrium, is obtained by assigning the
value zero to all infected classes, and then solving the equation $f(u)=0,$ to
deduce the rest of the terms of the FDE. On the other hand, for the EE, under
the hypothesis $\mathcal{R}_{0}>1.$ Solving the system $f(u)=0$ yields the
other corresponding equilibrium. We emphasize that for
$\eqref{eq:ode},\eqref{eq:pde},\eqref{eq:fde}$ and $\eqref{eq:sde}$ the EE can
be obtained by simple or long algebraic manipulations depending on the number
of compartments, whereas for $\eqref{eq:dde1},$ even for a small number of
compartments, further analysis is required. Once the existence of the FDE and
EE equilibria is shown, their local and global stability must be studied. For
most of the previous models, the corresponding system is linearized around the
corresponding equilibrium. For the FDE (resp EE), under the hypothesis
$\mathcal{R}_{0}<1,$ (resp $\mathcal{R}_{0}>1$) it can be shown that the
characteristic equation associated to the linearized system has roots with
strictly negative real part which yields the local asymptotic stability of the
FDE (resp EE). One should know that the conditions under which the roots have
strictly negative real parts are not always easy to obtain, especially when
the number of compartments increases. But, In virtue of the Routh-Hurwitz
criterion (see Gantmakher (2005)), assumptions on the parameters can be
concluded. On the other hand, the global stability can be very challenging for
all the previous models, given that the latter requires a suitable choice of
the Lyapunov function and unfortunately, no mathematical method can be used to
get the suitable choice, but several authors have tackled the problem and
proposed a general form that can work as Lyapunov functions for some
particular epidemiological models. For instance, Boukhouima et al. (2020);
Syafruddin & Noorani (2013); Hattaf & Yousfi (2020). If the Lyapunov function
is found, under the same assumptions on $\mathcal{R}_{0},$ and sometimes
further assumptions depending on the model, the global stability of the
equilibria can be proved using the classical Lyapunov stability theory (see
Bellman (2008)).

#### Numerical simulation of the model:

Once the mathematical analysis is concluded, to support the theoretical
results, a numerical simulation is in order. Before the simulation, the
question of the values assigned to the parameters arises. One approach to
solve this problem is to consider an estimation of the parameters, but such an
approach requires availability of the data up to the present within
governments’ institutions. As an example, in the case of models taking the
form (1), we refer to Marinov et al. (2014), where the authors used the method
of variational imbedding (MVI) to identify the transmission rate $\beta>0$ as
well as the recovery rate $\gamma>0$ of the SIR model. We highly recommend
Tarantola (1987) to the reader to get familiar with the approach used in the
cited paper. A different approach was used by Medvedeva et al. (2020), the
authors developed a numerical formula estimating the parameters of the SIR
model, based on the approximation of the classical derivative (first order and
second order) of the model and the minimization of a least square sum. Once
the parameters are estimated, the numerical simulation is carried on and a
variety of numerical methods can be used for this purpose. For instance, for
(1), the fourth order Runge-Kutta method is preferred due to its high
accuracy. for (3), the Adams-Bashforth method or its generalized version can
be used. For (5), one can think of the Finite Difference Method or the
Operator Splitting Method. For (4), a modified version of Runge-Kutta methods
can be adopted. Finally, for (2), one can use the Euler–Maruyama method or the
Milstein’s higher order method. To get a better understanding of these methods
we refer to Bellen & Zennaro (2013); Griffiths & Higham (2011); Owolabi
(2020); Diethelm et al. (2002); Kloeden & Platen (1999); Higham (2001);
Workman & Lenhart (2007).

#### The study of the optimal control problem:

When the model is simulated numerically, particularly, for the case
$\mathcal{R}_{0}>1,$ the problem of finding control strategies that stop the
epidemic from spreading further is studied. Usually this control measures
reside on introducing treatment and vaccination controls to the studied model.
The study of the optimal control problem is carried on theoretically at first
by proving the existence and uniqueness of such an optimal control, and then
by using direct or indirect numerical methods to compute it. For (1), (2) (3)
and (4), usually the indirect method is preferred due the reduction of the
problem to a system of initial value problems in virtue of the Pontryagin
maximum principle and some of its variations, this system is then solved by
means of the Forward-Backward Sweep Method. We refer to Tröltzsch (2010) to
gain more knowledge on both theoretical and numerical aspects of optimal
control theory and Lenhart & Workman (2007) for its application to biological
models.

##  Review of some of the main contributions in compartmental models

###  The foundation

Strictly speaking, the birth of compartmental epidemiological models was
attributed Kermack & McKendrick (1927). They formalized the concept of
compartmental models by using a set of ordinary differential equations to
describe the behavior of an epidemic. Their method is still considered valid,
and is largely used in recent research. The main objective of the two
researchers was to understand the reasons why the pandemic of Spanish Flu
didn’t infect the whole population. The model is composed into three
compartments as the following figure shows

Fig. 2: Diagram of the SIR model.

The three compartments are considered as follows

$\bullet$

$S(t):$ the class of those individuals who are capable of contracting the
disease and becoming infected.

$\bullet$

$I(t):$ the class of infected individuals, those who capable of transmitting
the disease to others.

$\bullet$

$R(t):$ the class of those individuals who have had the disease and are dead,
or have recovered and are permanently immune or are isolated.

This model has two parameters defined by

$\bullet$

$\beta>0:$ the infection/transmission rate.

$\bullet$

$\gamma>0:$ the recovery/removal rate.

The differential system governing the original SIR model is:

$\begin{cases}S^{\prime}(t)=-\beta I(t)S(t)\;\;t>0,\\\ I^{\prime}(t)=\beta
I(t)S(t)-\gamma I(t)\;\;t>0,\\\ R^{\prime}(t)=\gamma I(t)\;\;t>0,\\\ \text{+
Initial positive conditions.}\end{cases}$ (6)

such that $N=S(t)+I(t)+R(t)\;\;\;\forall t\geq 0$ is the total population.

It is important to state that the original SIR model has some limitations. For
instance, the vital dynamics (natural birth/death rate) are not captured by
the model, also, the model assumes that all the infected people are infectious
and are spreading the disease among the susceptible population, moreover, it
doesn’t consider the incubation period. Consequently, the diseases for which
this model can be applied are limited and further contributions to make it
more practicable are required.

Notation: All the reviewed models’ parameters will be renamed if necessary to
unify the notation. Unless stated otherwise,
$\Lambda,\beta,\gamma,\mu,d,\mathcal{R}_{0}$ and the letters
S,I,R,E,V,D,A,H,Q,C,T and M hold the same definitions that will be stated
thereafter.

###  Some of the main contributions to (1)

The contributions of this form take into account the following

  * •

Extension of the SIR model: This is usually done by considering more clinical
states resulting in the addition of more compartments to the model. Examples
of such added compartments are the exposed (E), the vaccinated (V), the
deceased (D), the asymptomatic (A), the hospitalized (H), the quarantined (Q),
the cross-immune (C), the treated (T), the maternally-derived immune (M), etc.
To each interaction of a compartment with another, a positive parameter is
added. Furthermore, to make the modified models more realistic, further vital
parameters (natural birth rate $\lambda>0$, natural mortality rate $\mu>0$,
death rate caused by the disease $d>0$) are taken into account. Another type
of extension concerns diseases in which different species may intervene (ex:
Humans and Mosquitoes for Dengue fever), or while considering heterogeneous
host populations. In this case, multi-group SIR and their extensions are
considered.

  * •

Modification/Generalization of the bilinear incidence rate: It is known that
the bilinear incidence rate $"\beta SI"$ of the original SIR model is based on
the law of the chemical mass action stating that the variation of the rate of
a chemical reaction is directly proportional to the product of the activities
or concentrations of the reactants. This incident rate is classic, and does
not take into account the saturation phenomena for a large number of infected
individuals. One way to remedy this problem is to consider other rates defined
by $g(I)S,$ $g(I)h(S),$ or in a more general way, one of the form $f(S,I),$
where $f,h$ and $g$ are functions satisfying adequate assumptions Capasso &
Serio (1978); Korobeinikov & Maini (2005). Some of the commonly well-known
incidence rates since the $20^{th}$ century are:

    * $\bullet$

$\dfrac{\beta SI}{N}:$ Standard Incidence rate Heesterbeek et al. (1996).

    * $\bullet$

$\dfrac{\beta SI}{1+aI}:$ Holling type incidence rate Anderson & May (1992).

    * $\bullet$

$\dfrac{\beta SI}{1+aS+bI}:$ Beddington-DeAngelis incidence rate Beddington
(1975).

    * $\bullet$

$\dfrac{\beta SI}{(1+aS)(1+\gamma I)}:$ Crowley-Martin incidence rate Crowley
& Martin (1989).

Such that $a$ and $b$ are strictly positive constants, $N$ is the sum of all
the considered compartments, while $\beta$ and $\gamma$ are respectively the
transmission and the recovery rates. Note that for a fixed value of $S,$ and
for a very large value of the infected class, all the above incidence rates
become bounded and the crowding effect is taken into account which is
practically acceptable. Not the same can be said in the case for the classical
bilinear rate. Below, we present a figure showing the four incidence rates
along with the classical incident rate for the SIR model.

Fig. 3: Some of the early proposed incidence rates compared to the bilinear
incidence rate

The work on new proposed incidence rates and extended SIR models continued in
the $21^{st}$ century. Indeed, Xiao & Ruan (2007) worked with the incidence
rate given by $\dfrac{\beta IS}{1+\alpha I^{2}}$ with $\alpha>0.$ Their aim
was to take into account the change on the behavior of the susceptible
population when they become aware of the high number of infected individuals.
Their new proposed SIRS model was expressed in the following system of odes

$\begin{cases}\begin{aligned}
&S^{\prime}(t)=\Lambda-\mathrm{\mu}S(t)-\frac{\beta S(t)I(t)}{1+\alpha
I^{2}(t)}+\kappa R(t),\\\ &I^{\prime}(t)=\frac{\beta S(t)I(t)}{1+\alpha
I^{2}(t)}-(\mu+\kappa)I(t),\\\ &R^{\prime}(t)=\gamma I(t)-(\mu+\kappa)R(t),\\\
&\text{+ positive initial conditions.}\end{aligned}\end{cases}$ (7)

with one added parameter $\kappa>0,$ being the rate in which the recovered
lose immunity and return to the susceptible class.

In the same year, the research related to incidence rates begun to steer to
the generalization of the previous proposed ones. In fact, Korobeinikov (2006)
proposed the following SIR model with a general incidence rate $f(S,I)$

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\mu-f(S(t),I(t))-\mu S(t),\\\
&I^{\prime}(t)=f(S(t),I(t))-(\mu+\gamma)I(t),\\\ &R^{\prime}(t)=\gamma
I(t)-\mu R(t),\\\ &\text{+ positive initial
conditions.}\end{aligned}\end{cases}$ (8)

where the birth/death rates are taken to be equal to $\mu>0$ and such that $f$
is is a continuously differentiable function satisfying the following
assumptions

    1. 1.

$f:(S,I)\in Q\rightarrow f(S,I)\in\mathbb{R}^{+}$ is a monotonically
increasing function $\forall S,I>0.$

    2. 2.

$f(S,0)=f(0,I)=0\;\forall S,I>0.$

    3. 3.

$f$ is concave with respect to I.

where $Q:=[0,1]\times[0,1]\subset\mathbb{R}_{+}^{2}.$

######  Remarks 3.1.

It should be noted that while the assumptions are biologically intuitive, each
one serves a mathematical purpose. Assumption 1 is for the mathematical well
posedness, while 2 is for the existence of the FDE and 3 is for the global
stability analysis of the EE.

Korobeinikov (2006) proved the global assymptotic stability of the FDE (resp
EE) under the condition $\mathcal{R}_{0}:=\dfrac{1}{\mu+\gamma}\dfrac{\partial
f(S,I)}{\partial I}\bigg{|}_{FDE}\leq 1$ (resp $\mathcal{R}_{0}>1$) by means
of the Lyapunov function. Later, Buonomo & Rionero (2010) considered the same
SIRS model with general incidence rate, and they proved the local stability of
the FDE and the EE under weaker assumptions on the general incidence rate
$f(S,I).$
The work on the development of adequate incidence rates for the SIRS model
continued. For instance, Li et al. (2017) considered a general incidence rate
of the form $f(I)S$ with $f$ verifying some adjusted assumptions.

  * $\bullet$

Introduction of immigration: Parallel progress in ameliorating the SIR model
from a different perspective begun when Sigdel & McCluskey (2014) proposed an
SEI model with immigration. Their aim was to take into account individuals
that enter a population and their contribution into maintaining the presence
of tuberculosis disease. The associated model was expressed in the following
system of odes

$\begin{cases}\begin{aligned}
&S^{\prime}(t)=(1-p-q)\Lambda-c\beta\frac{S(t)I(t)}{N(t)}-\mu
S(t)+r_{1}E(t)+r_{2}I(t),\\\
&E^{\prime}(t)=p\Lambda+c\beta\frac{S(t)I(t)}{N(t)}-\left(k+\mu+r_{1}\right)E(t),\\\
&I^{\prime}(t)=q\Lambda+kE(t)-\left(\mu+d+r_{2}\right)I(t),\\\ &\text{+
positive initial conditions.}\end{aligned}\end{cases}$ (9)

with the following assumptions

    * $\bullet$

$c>0$, $\Lambda$ is the total rate in which new individuals enter the whole
population while $(1-p-q)\Lambda$ in the rate in which they enter the
susceptible class, $p\Lambda$ the one in which they enter the exposed class
and finally $q\Lambda,$ the one is which they enter the infected class. With
the assumption $p+q\in[0,1].$

    * $\bullet$

$r_{1}$ (resp $r_{2}$) is a strictly positive rate in which exposed (resp
infected) are treated .

In the same paper, Sigdel & McCluskey (2014) tried to generalize the
recruitment rate of immigrants by considering the recruitment function $B(N),$
where $N:=S+E+I.$ The following assumptions were made

    * $\bullet$

(H1) There exists a unique $N_{0}>0$ such that $B\left(N_{0}\right)-\mu
N_{0}=0$.

    * $\bullet$

(H2) $B^{\prime}(N)<\frac{B(N)}{N}$ for $N\in\left(0,N_{0}\right]$.

    * $\bullet$

(H3) $b_{0}=\lim_{N\rightarrow 0+}\frac{B(N)}{N}>\mu$, allowing that $b_{0}$
may be infinite.

$(H1)$ and $(H3)$ are biologically related assumptions while $(H2)$ is to
ensure the existence and uniqueness of the endemic equilibrium.
Later on, several authors’ aim was to extend the previous models into ones
that take into account the concept of immigration.

For instance, Cui et al. (2020) considered an SIRS model with a constant
immigration and an incidence rate of the form $g(I)S$ with $g$ verifying two
assumptions.

Driven by the same motivation, Khan, Alaoui, Zeb, Tilioua & Djilali (2021)
considered an SEI model with immigration and a general incidence rate of the
form $f(S,I)$.

  * $\bullet$

Introduction of the treatment: One last angle of contribution that we discuss
in this section is the one related to the finding of a suitable treatment
function. This begun when Wang & Ruan (2004) remarked that the removal rate
being linear in terms of the infected in not practicable because of the
limited capacity of treatment in each country. This resulted in the proposal
of the following treatment function

$h(I)=\begin{cases}r&\text{ for }I>0,\\\ 0&\text{ for }I=0.\end{cases}$

where $r>0$ is a constant representing the capacity of treatment for
infectives.
The model is formulated in the following system of odes

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\Lambda-\mu S(t)-\beta
S(t)I(t),\\\ &I^{\prime}(t)=\beta S(t)I(t)-(\mu+\gamma)I(t)-h(I(t)),\\\
&R^{\prime}(t)=\gamma I(t)+h(I(t))-\mu R(t),\\\ &\text{+ positive initial
conditions.}\end{aligned}\end{cases}$ (10)

One limitation of the proposed treatment function is that a constant treatment
is always considered even when the treatment capacity is not reached.
Wang (2006) solved this limitation by considering the following treatment
function

$T(I)=\begin{cases}rI&\text{ if }0\leqslant I\leqslant I_{0},\\\ k&\text{ if
}I>I_{0}.\end{cases}$

where $k=rI_{0}.$ The modified model was given by

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\Lambda-\mu S(t)-\beta
S(t)I(t),\\\ &I^{\prime}(t)=\beta S(t)I(t)-(\mu+\gamma+d)I(t)-T(I(t)),\\\
&R^{\prime}(t)=\gamma I(t)+T(I(t))-\mu R(t),\\\ &\text{+ positive initial
conditions.}\end{aligned}\end{cases}$

For an outbreak disease such as SARS, Zhonghua & Yaohong (2009) thought
replacing the treatment rate function $T$ with a saturated rate of Holling
type would be more pertinent. The considered treatment function was

$h(I)=\frac{\beta I}{1+\alpha I},\;\;\;\alpha>0.$

The related contributions done after were mostly proposing new models obtained
by taking combinations of saturated rates for the treatment function and some
of the recent studied incidence rates introduced in this paper.

We give a figure illustrating the reviewed treatment rates as well as some
that have been inspired from the well-known incidence rates.

Fig. 4: Some of the commonly used treatment rates

Below, we present a table summarizing the previous discussed contributions
falling into $\eqref{eq:ode}$

Table 1: Some of the main contributions falling into category $\eqref{eq:ode}$ Ref. | Model | Incidence |  Treatment |  Contributions
---|---|---|---|---
Xiao & Ruan (2007) | SIR | $\dfrac{\beta IS}{1+\alpha I^{2}},\;\;\alpha>0$ |  None |  Proposing a new incidence rate taking into account the behavioral change towards the infected.
Korobeinikov (2006) | SIR | $f(S,I)$ |  None |  Proposing a general incidence rate for the SIR model, global stability of the DFE & EE for the model under some assumptions on $f.$
Buonomo & Rionero (2010) | SIRS | $f(S,I)$ |  None |  Local stability of the SIRS model with general incidence rate under weakened conditions on $f.$
Table 2: Some of the main contributions falling into category $\eqref{eq:ode}$ continued Ref. | Model | Incidence |  Treatment |  Contributions
---|---|---|---|---
Sigdel & McCluskey (2014) | SEI | Bilinear |  Constant treatment rate, General treatment rate |  Proposition of a new SEI model taking into account the affect of immigration, generalization of the model by considering a general immigration rate, global stability with constant/general immigration rate.
Khan, Alaoui, Zeb, Tilioua & Djilali (2021) | SEI | $f(S,I)$ |  None |  Proposing a SEI model with a general incidence rate and constant immigration rates within the the infected and the exposed, global stability of the FDE and EE, Numerical simulation for the particular bilinear incidence rate.
Wang & Ruan (2004) | SIR | Bilinear |  treatment rate with limited capacity |  Proposing an SIR model with a more realistic treatment rate taking into account the limited resources, bifurcation analysis, introducing conditions under which the disease is instinct and limited treatment is sufficient.
Wang (2006) | SIR | Bilinear |  Infected-depending treatment rate with limited capacity |  Tackling the limitation of the treatment rate in Wang & Ruan (2004) by introducing a rate proportional to the infected class as long as the capacity is not reached, study of existence and uniqueness of the DFE & EE, bifurcation analysis.

It should be noted that the efficiency of contributions falling into category
(1) depends on the better understanding of the modeled disease. For a simple
one, the limitation to a reduced number of compartments and a bilinear
incidence rate may suffice Khan et al. (2019). But for complicated diseases,
such as Covid-19 or seasonal diseases, a large number of compartments is
required for the first and an adequate choice of the incidence rate is needed
for the second to capture the maximum of possible intervening factors in each
one of the modeled diseases Ndaïrou et al. (2020); Xing & Li (2021).

As far as the contributions of category $\eqref{eq:ode}$ can go, this type of
modeling will always have several limitations, some of which are:

  1. $\bullet$

Limitation 1: They don’t capture the uncertainty and variability that is
inherent in real-life epidemics due to factors such as the unpredictability of
person-to-person contact.

  2. $\bullet$

Limitation 2: They don’t take into account the memory effects on their
dynamics (dependency on the past). As a matter of fact, given the knowledge of
the history of a disease, people can use different precautions that affect the
change in the dynamics.

  3. $\bullet$

Limitation 3: They’re based on the assumption of instant cause and effect
which is not correct in real life. For instance, an infected individual needs
some time-period before becoming infectious and starts spreading the disease.
Furthermore, when treatment is applied, a time-period is also needed in order
for the effect of the treatment to start showing.

  4. $\bullet$

Limitation 4: They don’t take into account the spatial factor. For infectious
diseases, the movement of individuals is a major factor of the disease spread,
thus the densities of different compartments of the model shouldn’t be just in
terms of time but also space.

###  Some of the main contributions to (2)

The essential aim of contributions falling into category $\eqref{eq:sde}$ is
to overcome ”Limitation 1” by adapting all the previous contributions of
category $\eqref{eq:ode}$ to a stochastic point of view.

Tornatore et al. (2005) adjusted the deterministic SIR model to the following
stochastic one

$\begin{cases}\begin{aligned} &\mathrm{d}S(t)=(-\beta S(t)I(t)-\mu
S(t)+\mu)\mathrm{d}t-\sigma S(t)I(t)\mathrm{d}W(t),\\\ &\mathrm{d}I(t)=(\beta
S(t)I(t)-(\gamma+\mu)I(t))\mathrm{d}t+\sigma S(t)I(t)\mathrm{d}W(t),\\\
&\mathrm{d}R(t)=(\gamma I(t)-\mu R(t))\mathrm{d}t,\\\ &\text{+ positive
initial conditions.}\end{aligned}\end{cases}$ (11)

with the natural birth and death rates being equal to $\mu>0.$ The authors
were able to prove that the global stability of the DFE (resp EE) holds under
the condition
$0<\beta<\min\left\\{\gamma+\mu-\frac{\sigma^{2}}{2},2\mu\right\\}.$ (resp
$\beta>\lambda+\mu+\frac{\sigma^{2}}{2}$). One can remark that for $\sigma=0,$
this condition takes into account that
$\mathcal{R}_{0}:=\dfrac{\beta}{\gamma+\mu}<1$ (resp $\mathcal{R}_{0}>1$)
which is exactly the same condition that assures the asymptotic global
stability of the DFE (resp EE) associated the deterministic counterpart of the
considered stochastic SIR model. Numerical simulations showed important
oscillations of the obtained solutions, which captures the random behavior of
disease dynamics in reality.

Lu (2009) extended (11) by taking into account the loss of immunity to the
disease. The author proposed the following SIRS model

$\begin{cases}\begin{aligned} &\mathrm{d}S(t)=(-\beta S(t)I(t)-\mu S(t)+\gamma
R(t)+\mu)\mathrm{d}t-\sigma S(t)I(t)\mathrm{d}W(t),\\\
&\mathrm{\leavevmode\nobreak\ d}I(t)=(\beta
S(t)I(t)-(\gamma+\mu)I(t))\mathrm{d}t+\sigma S(t)I(t)\mathrm{d}W(t),\\\
&\mathrm{\leavevmode\nobreak\ d}R(t)=(\gamma
I(t)-(\mu+\lambda)R(t))\mathrm{d}t,\\\ &\text{+ positive initial
conditions.}\end{aligned}\end{cases}$ (12)

Under the assumption $\beta<\lambda+\mu-\frac{\sigma^{2}}{2}$, stochastic
asymptotic stability of the FDE has been proved based on the Lyapunov method.

Jiang et al. (2011) proposed the following SIR model with a different added
noise term

$\begin{cases}\begin{aligned} &\mathrm{d}S(t)=(\Lambda-\beta S(t)I(t)-\mu
S(t))\mathrm{d}t+\sigma_{1}S(t)\mathrm{d}B_{1}(t),\\\ &\mathrm{d}I(t)=(\beta
S(t)I(t)-(\mu+d+\gamma)I(t))\mathrm{d}t+\sigma_{2}I(t)\mathrm{d}B_{2}(t),\\\
&\mathrm{d}R(t)=(\gamma I(t)-\mu
R(t))\mathrm{d}t+\sigma_{3}R(t)\mathrm{d}B_{3}(t),\\\ &\text{+ positive
initial conditions}.\end{aligned}.\end{cases}$ (13)

with natural birth rate (resp death rate) $\Lambda>0$ (resp $\mu>0$) and a
death rate caused by the disease $d>0$ as well as different stochastic
intensities $\sigma_{i}>0$ $i\in\\{1,2,3\\}$ and different Wiener processes
$B_{i}$ $i\in\\{1,2,3\\}.$ The authors proved the existence and uniqueness of
a global positive solution to the model. Moreover, it was shown that if
$\sigma_{1}^{2}<\mu$ and $\sigma_{2}^{2}<2(\mu+d+\gamma)$ then the boundedness
of solutions is assured, they also studied the asymptotic behavior of the
solution around the deterministic FDE and EE equilibria associated the
deterministic counterpart of the proposed model. Precisely, if the boundedness
assumptions are verified and
$\mathcal{R}_{0}:=\dfrac{\beta\Lambda}{\mu(\mu+d+\gamma)}<1$ then the
oscillations of solutions around the DFE associated to the deterministic
counterpart of the model increase in terms of $\sigma_{1}$ and $\sigma_{2}.$
For $\mathcal{R}_{0},$ the oscillations also occur around the state
$\left(\dfrac{2\mu}{2\mu-\sigma_{1}^{2}}S^{*},\dfrac{2\mu\left(\mu+d+\gamma-p\gamma^{2}\right)}{2\mu(\mu+d+\gamma)-p\gamma^{2}-\mu\sigma_{2}^{2}}I^{*},\dfrac{\mu}{\mu-\sigma_{3}^{2}}R^{*}\right)$
provided that
$\sigma_{1}^{2}<2\mu,\quad\sigma_{2}^{2}<2(\mu+d+\gamma),\quad\sigma_{3}^{2}<\mu$.
Such that the EE is the endemic equilibrium of the deterministic counterpart
of the model.

The numerical simulations were done using the Milstein’s Higher Order Method
Higham (2001) to support the obtained theoretical results and showed that for
a small value of the stochastic intensities, the boundedness of the solution
is assured and for a large time, the different stochastic equilibria approach
their deterministic counterparts.

Rao (2014) modified $\eqref{eq:lu}$ by taking a Holling-type incidence rate,
the noise was added to the three compartments and was supposed proportional to
the distances of the states from the steady ones corresponding to the endemic
equilibrium of the deterministic counterpart of the model. The proposed model
was as follows

$\begin{cases}\begin{aligned} &\mathrm{d}S(t)=\left(\Lambda-dS(t)-\frac{\beta
S(t)I(t)}{1+aI(t)}+\lambda
R(t)\right)\mathrm{d}t+\sigma_{1}\left(S(t)-S^{*}\right)\mathrm{d}B(t),\\\
&\mathrm{d}I(t)=\left(\frac{\beta S(t)I(t)}{1+\alpha
I(t)}-(\mu+\gamma)I(t)\right)\mathrm{d}t+\sigma_{2}\left(I(t)-I^{*}\right)\mathrm{d}B(t),\\\
&\mathrm{d}R(t)=(\gamma
I(t)-(\mu+\lambda)R(t))\mathrm{d}t+\sigma_{3}\left(R(t)-R^{*}\right)\mathrm{d}B(t),\\\
&\text{+ positive initial conditions.}\end{aligned}\end{cases}$ (14)

With a natural birth rate $\Lambda>0,$ a natural death rate $\mu>0$ and
different stochastic intensities $\sigma_{1},\sigma_{2},\sigma_{3}>0,$ and
such that $(S^{*},I^{*},R^{*})$ is the EE of the deterministic counterpart of
the model obtained by taking $\sigma_{1}=\sigma_{2}=\sigma_{3}=0.$
The author established the mathematical well posedness of the model as well as
the biological one. The stochastic permanence of the disease was proven to
occur if
$\max\left\\{\sigma_{1}^{2},\sigma_{2}^{2},\sigma_{3}^{2}\right\\}<2(\Lambda-\mu)$.
Numerical simulations using the Milstein’s Higher Order Method Higham (2001)
have been done for different values of the intensities to show the
contribution of the additional noise to the model.

El Koufi et al. (2019) tried to adapt a deterministic SIR model with
vaccination and vertical transmission to a stochastic point of view with a
more generalized incidence function. Their proposed model was

$\begin{cases}\begin{aligned} &dS(t)=\bigg{(}-\frac{\beta
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)+\alpha_{3}S(t)I(t)}-bS(t)+(1-m)pdI(t)\\\
&\hskip 22.76228pt+b(1-m)(S(t)+R(t))\bigg{)}dt-\frac{\sigma
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)+\alpha_{3}S(t)I(t)}dB(t),\\\
&dI(t)=\left(\frac{\beta
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)+\alpha_{3}S(t)I(t)}-(pd+\gamma)I(t)\right)dt\\\
&\hskip 19.91684pt+\frac{\sigma
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)+\alpha_{3}S(t)I(t)}dB(t),\\\
&dR(t)=(\gamma I(t)-bR(t)+dmpI(t)+mb(S(t)+R(t)))dt,\\\ +&\text{ positive
initial conditions.}\end{aligned}\end{cases}$ (15)

With the additional parameters $b>0$ (resp $d>0$) is the mortality rate in the
susceptible and the recovered (resp infective) individuals, $p$ is the
proportion of the offspring of the infective and $q>0$ is the proportion of
the rest that are born infected, with $p+q=1$ and
$\alpha_{1},\alpha_{2},\alpha_{3}>0.$ Finally $m\in[0,1]$ is the successful
vaccination proportion to the newborn from the susceptible and the recovered.
The value of the work presented by the authors is that it implicitly deals
with all other known transmition rates (Beddington-DeAngellis, Crowley Martin,
Holling-type) discussed earlier in this paper.
The authors proved the existence and uniqueness of a global positive solution
to the model, they also proved that if either
$\sigma^{2}>\dfrac{\beta^{2}}{(2(pd+r))}$ or ( $R_{s}<1$ and
$\sigma^{2}<\beta$ ) the disease dies out, while if $R^{*}_{s}>1,$ then the
disease persists in the mean. Where $R_{s}$ (resp $R^{*}_{s}$) is the
extinction (resp persistence) threshold given by

$\begin{cases}R_{s}=\dfrac{R_{0}}{(1-m)}-\dfrac{\sigma^{2}}{2\left(1+\alpha_{1}(1-m)\right)^{2}(pd+r)}.\\\
R_{s}^{*}=(1-m)R_{0}-\dfrac{\sigma^{2}}{2\left(1+\alpha_{1}(1-m)\right)^{2}(1-m)(pd+r)}.\end{cases}$

such that $\mathcal{R}_{0}$ is the basic reproduction number associated to the
deterministic counterpart of the proposed model (i.e for $\sigma=0$). The
numerical simulation of the model confirmed the theoretical results, and a
comparison of solutions with the deterministic case $\sigma=0$ was
illustrated.

###  Some of the main contributions to (3)

The purpose of this type of contributions is to deal with Limitation 2 by
including the ”memory effect” to the dynamics of the and adapting the
previously introduced models and the contributions made as well. This is done
by replacing the classical derivative with a well chosen fractional one.
Ozalp & Demi̇rci̇ (2011) proposed the following SEIR model

$\begin{cases}\begin{aligned} &D^{\alpha}S(t)=\Lambda-\frac{p\Lambda
E(t)}{N(t)}-\frac{q\Lambda I(t)}{N(t)}-r\frac{S(t)I(t)}{N(t)}-\mu S(t),\\\
&D^{\alpha}E(t)=\frac{p\Lambda E(t)}{N(t)}+\frac{q\Lambda
I(t)}{N(t)}+r\frac{S(t)I(t)}{N}-\mu E(t)-\beta E(t),\\\ &D^{\alpha}I(t)=\beta
E(t)-\mu I(t)-dI(t)-\gamma I(t),\\\ &D^{\alpha}R(t)=\gamma I(t)-\mu R(t),\\\
&\text{+positive initial conditions.}\end{aligned}\end{cases}$ (16)

such that $r>0$ is the horizontal transmission rate between the susceptible
and the exposed. $p\in[0,1]$ resp ($q\in[0,1]$) is the probability of the
offspring of exposed (resp infected) being born in the exposed class,
$\beta>0$ is the rate of exposed individuals becoming infectious. and
$D^{\alpha}$ is the Caputo fractionnal order derivative of order
$\alpha\in(0,1)$ de Oliveira & Tenreiro Machado (2014).

The authors proved the mathematical well posedness of the model as well as the
positiveness and boundedness of the unique solution. Under the assumption
$\mathcal{R}_{0}:=\dfrac{(q\mu+r)\beta}{(\mu+d+\gamma)(\mu-p\mu+\beta)}<1$,
the DFE was proven to be asymptotically stable. For the endemic equilibrium,
the application of the fractional Routh-Hurwitz criterion (see Ahmed et al.
(2006)) permitted the authors to obtain assumptions on the model’s parameters
under which the EE is asymptotically stable. The numerical simulation was
carried on by means of the generalized Adam-Bashford method and an adequate
choice of the parameters illustrated the theoretical stability results.

In the same year, to model Influenza, El-Shahed & Alsaedi (2011) replaced the
classical derivative of the model proposed by Casagrandi et al. (2006) with
the Caputo fractional one to get the following modified SIRC model

$\begin{cases}\begin{aligned} &D^{\alpha}S(t)=\mu(1-S(t))-\beta
S(t)I(t)+\gamma C(t),\\\ &D^{\alpha}I(t)=\beta S(t)I(t)+\sigma\beta
C(t)I(t)-(\mu+\theta)I(t),\\\ &D^{\alpha}R(t)=(1-\sigma)\beta C(t)I(t)+\theta
I(t)-(\mu+\delta)R(t),\\\ &D^{\alpha}C(t)=\delta R(t)-\beta
C(t)I(t)-(\mu+\gamma)C(t),\\\ &\text{+ positive initial
conditions.}\end{aligned}\end{cases}$ (17)

Such that the birth and death dates are equal to $\mu>0$, $\gamma^{-1}$ is the
cross-immune period $\theta^{-1}$ is the infectious period, $\delta^{-1}$ is
the total immune period and finally $\sigma$ is the fraction of the exposed
cross-immune recruited into the infective subpopulation.

All the usual steps of study were followed, the existence, uniqueness,
boundedness and positiveness of the solution were proven. The local asymptotic
stability of the FDE is assured under the condition
$R_{0}:=\dfrac{\beta}{(\mu+\theta)}<1$. For the endemic equilibrium, the
authors uded the fractional Routh-Hurwitz criterion (see Ahmed et al. (2006))
which allowed them to obtain conditions under which the EE is locally
asymptotically stable.

Mouaouine et al. (2018) considered the following SIR model with fractional
Caputo derivative.

$\begin{cases}\begin{aligned} &D^{\alpha}S(t)=\Lambda-\mu S(t)-\frac{\beta
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I+\alpha_{3}S(t)I(t)},\\\
&D^{\alpha}I(t)=\frac{\beta
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)+\alpha_{3}S(t)I(t)}-(\mu+d+\gamma)I(t),\\\
&D^{\alpha}R(t)=\gamma I(t)-\mu R(t),\\\ &\text{+ positive initial
conditions.}\end{aligned}.\end{cases}$ (18)

With $\alpha_{1},\alpha_{2},\alpha_{3}>0.$ Note that this model holds the
value of studying several models that can be considered with each incidence
rate by studying only one with a generalized incidence rate.
As always, the mathematical and biological well posedness were established.
The local asymptotic stability of the DFE (resp EE) was proven under the
condition
$\mathcal{R}_{0}:=\dfrac{\beta\Lambda}{\left(\mu+\alpha_{1}\lambda\right)(\mu+d+r)}\leq
1$ (resp $\mathcal{R}_{0}>1$), while using the Routh-Hurwitz criterion (see
Ahmed et al. (2006)) for the EE case. The global stability of the FDE and EE
was proven under the same previous conditions by means of the Lyapunov method.
Numerical simulations were carried out by using an algorithm based on the
fractional Euler’s method (see Odibat & Momani (2008)) and supported the
theoretical results obtained.

In the aim of adapting the contribution of the generalization of the incidence
rate to the fractional case, Altaf Khan et al. (2020) proposed the two
following SIR models

$\begin{cases}\begin{aligned} &D_{C}^{\alpha}S(t)=\Lambda-\mu
S(t)-S(t)f(I(t))+\mu_{1}I(t)+\theta R(t)\\\
&D_{C}^{\alpha}I(t)=S(t)f(I(t))-\left(\mu+\mu_{1}+\mu_{2}+d\right)I(t)\\\
&D_{C}^{\alpha}R(t)=\mu_{2}I(t)-(\mu+\theta)R(t),\\\ &\text{+ positive initial
conditions}.\end{aligned}\end{cases}$ (19) $\begin{cases}\begin{aligned}
&D_{ABC}^{\alpha}S(t)=\Lambda-\mu S(t)-S(t)f(I(t))+\mu_{1}I(t)+\theta R(t)\\\
&D_{ABC}^{\alpha}I(t)=S(t)f(I(t))-\left(\mu+\mu_{1}+\mu_{2}+d\right)I(t)\\\
&D_{ABC}^{\alpha}R(t)=\mu_{2}I(t)-(\mu+\theta)R(t),\\\ &\text{+ positive
initial conditions}.\end{aligned}\end{cases}$ (20)

Such that $D^{\alpha}_{C}$ (resp $D^{\alpha}_{ABC}$) is the Caputo (resp
Atangana-Baleanu-Caputo) fractional derivative, $\mu_{1}>0$ resp ($\mu_{2}>0$)
is the rate of movement of the infected to the susceptible (resp recovered)
class and $\theta>0$ is the rate in which the recovered population becomes
susceptible. $f$ was supposed to be non-negative and continuously
differentiable in $\mathbb{R}^{*}_{+},$ locally lipschitz on $\mathbb{R}_{+}$
and satisfying

  1. $\bullet$

$f(0)=0$ and $f(I)>0\;\;\forall I>0.$

  2. $\bullet$

$\dfrac{f(I)}{I}$ is continuous and monotonically non-increasing $\forall I>0$
and
$\lim_{I\rightarrow 0^{+}}\dfrac{f(I)}{I}$ exists.

######  Remarks 3.2.

Although this model was considered by the authors to be an SIR model, the
susceptibility of the recovered population in real life is generally due to
immunity-loss to the disease. Taken that into account, the model is rather an
SIRS model.

For both models, The authors proved the mathematical and biological well
posedness of the two models. The DFE was proven to be asymptotically locally
and globally stable under the condition
$\mathcal{R}_{0}:=\dfrac{\Lambda}{\mu\left(\mu+\mu_{1}+\mu_{2}+d\right)}\dfrac{\partial
f(0)}{\partial I}<1$, whereas for the EE, it was proven to be asymptotically
globally stable under the condition $\mathcal{R}_{0}>1.$ To show the
importance of the choice of the incidence functions, numerical simulations
were done for all three of them (Bilinear, Holling type and Beddington-
DeAngellis). Interpretations in terms of the fractional order $\alpha$ were
given, precisely, the decreasing (resp increasing) of the susceptible and
recovered (resp the infected) for increased values of $\alpha.$ The numerical
simulation also proved the feasibility of the Atangana-Baleanu-Caputo
fractional derivative over the Caputo one.

###  Some of the contributions to (4)

For this kind of contributions, Limitation 3 is dealt with by adding of the
notion of delay to the model to get a model falling into category (4).

Kaddar (2009) considered the following SIR model with a Beddington-DeAngelis
incidence rate

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\Lambda-\mu S(t)-\frac{\beta
S(t-\tau)I(t-\tau)}{1+\alpha_{1}S(t-\tau)+\alpha_{2}I(t-\tau)},\\\
&I^{\prime}(t)=\frac{\beta
S(t)I(t)}{1+\alpha_{1}S(t)+\alpha_{2}I(t)}-(\mu+d+\gamma)I(t),\\\
&R^{\prime}(t)=\gamma I(t)-\mu R(t),\\\ &\text{+ positive initial
functions.}\end{aligned}\end{cases}$ (21)

With $\alpha_{1},\alpha_{2}>0$ and $\tau>0$ representing the incubation
period.

The author proved that the FDE is locally asymptotically stable under the
condition
$\mathcal{R}_{0}:=\dfrac{\Lambda\left(\beta-\alpha_{1}(\mu+\alpha+\gamma)\right)}{\mu(\mu+\alpha+\gamma)}<1.$
One can remark that the incubation period does not affect the stability of the
FDE. For the EE, under the assumptions $\mathcal{R}_{0}<1$ and
$\alpha_{2}\mu<\beta-\alpha_{1}(\mu+\alpha+\gamma)$, the author established
the existence of a critical incubation period $\tau_{0}$ such that if
$\tau\in[0,\tau_{0})$ the EE is locally asymptotically stable and unstable for
$\tau>\tau_{0}.$ For the critical case $\tau=\tau_{0}$ bifurcation from the EE
occurs. Numerical simulation was done separately for the three possible values
of $\tau$ and confirmed the obtained theoretical results.

Global stability of an SIR model with Holling type incidence rate was studied
by McCluskey (2010). The considered model was

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\Lambda-\mu_{1}S(t)-\frac{\beta
S(t-\tau)I(t-\tau)}{1+aI(t-\tau)},\\\ &I^{\prime}(t)=\frac{\beta
S(t-\tau)I(t-\tau)}{1+aI(t)}-(\mu_{2}+\gamma)I(t),\\\ &R^{\prime}(t)=\gamma
I(t)-\mu_{3}R(t),\\\ &\text{+ positive initial
functions.}\end{aligned}\end{cases}$ (22)

Such that $\tau>0$ is the latent period and $\mu_{1},\mu_{2},\mu_{3}$ are
respectively the natural death rates of the susceptible, infected and
recovered population. For biological reasons, it was supposed that the removal
of infectives is at least as fast as the removal of susceptibles the
assumption. Mathematically speaking, $\mu_{1}\leq\mu_{2}+\gamma$.

It should be known that the model was studied earlier by Xu & Ma (2009) in
which they proved that the FDE is globally asymptotically stable under the
condition
$\mathcal{R}_{0}:=\dfrac{\Lambda\beta}{\mu_{1}\left(\mu_{2}+\gamma\right)}<1$
. Whereas, if $\mathcal{R}_{0}>1$ and $\underbrace{\Lambda
a>\mu_{2}+\gamma}_{(H)}$, the EE is locally stable. McCluskey’s contribution
was discarding $(H)$ by considering the function $g(x):=x-1-ln(x)$ and using
it to define the following carefully chosen Lyapunov function

$V(t):=\frac{1}{\beta
f\left(I^{*}\right)}g\left(\frac{S(t)}{S^{*}}\right)+\frac{I^{*}}{\beta
S^{*}f\left(I^{*}\right)}g\left(\frac{I(t)}{I^{*}}\right)+\int_{0}^{\tau}g\left(\frac{I(t-s)}{I^{*}}\right)\mathrm{d}s$

where

$f(x)=\dfrac{x}{1+\alpha x},\quad\alpha>0.$

In the aim of adapting the generalization of the incidence rate to this type
of models, Li & Liu (2014) proposed in 2014 the following model with a general
incidence and recovery rates.

$\begin{cases}\begin{aligned} &S^{\prime}(t)=\Lambda-\mu S(t)-F(S(t),I(t))\\\
&I^{\prime}(t)=e^{-\mu\tau}F(S(t-\tau),I(t-\tau))-(\mu+d)I(t)-G(I(t)),\\\
&R^{\prime}(t)=G(I(t))-\mu R(t),\\\ &\text{+ positive initial
functions.}\end{aligned}\end{cases}$

Such that $\tau\geq 0$ is the latent period and the term $e^{-\mu\tau}$ is
used to model the survival rate of the population. All along the authors’
study the total assumptions fixed on the functions $F$ and $G$ are

  1. 1.

$F(S,I)>0,\dfrac{\partial F(S,I)}{\partial I}>0,\dfrac{\partial
F(S,I)}{\partial S}>0$ $\forall$ $S,I>0$.

  2. 2.

$F(S,0)=F(0,I)=0,\dfrac{\partial F(S,I)}{\partial
S}\bigg{|}_{(S,0)}=0,\dfrac{\partial F(S,I)}{\partial I}\bigg{|}_{(S,0)}>0$
$\forall S,I>0$.

  3. 3.

$G(0)=0,\;G^{\prime}(I)>0$ $\forall I\geq 0.$

  4. 4.

$F^{\prime}(S,0)$ is increasing $\forall S>0$.

  5. 5.

$F(S,I)\leq I\dfrac{\partial F(S,I)}{\partial I}\bigg{|}_{(S,0)}\forall I>0.$

  6. 6.

$G^{\prime}(0)\leq\dfrac{G(I)}{I}\;\forall I>0.$

  7. 7.

$\dfrac{I}{I^{*}}\leq\dfrac{F(S,I)}{F\left(S,I^{*}\right)}$ $\forall
I\in\left(0,I^{*}\right),\dfrac{F(S,I)}{F\left(S,I^{*}\right)}\leq\dfrac{I}{I^{*}}$
$\forall I\geq I^{*}.$

  8. 8.

$\dfrac{G(I)}{G\left(I^{*}\right)}\leq\dfrac{I}{I^{*}}$ $\forall
I\in\left(0,I^{*}\right)$ and
$\dfrac{I}{I^{*}}\leq\dfrac{G(I)}{G\left(I^{*}\right)}$ $\forall I\geq I^{*}$.

The first condition of assumption 2 allowed the authors to prove the existence
of the FDE. $\forall\tau>0,$ the assumptions $1,2$ and $3$ along with
$\mathcal{R}_{0}:=\dfrac{e^{-\mu\tau}\dfrac{\partial F(S,I)}{\partial
I}\bigg{|}_{FDE}}{\mu+\alpha+G^{\prime}(0)}<1$ allowed them to prove the
existence of the EE. For $\mathcal{R}_{0}\leq 1,$ the local asymptotic
stability of the FDE was established under the assumptions $1,2,3$ while it
was disproven if $\mathcal{R}_{0}>1.$ For $\mathcal{R}_{0}\leq 1,$ the global
stability of the FDE was proven and required the assumptions $1-6.$ Whereas
for $\mathcal{R}_{0}>1$ the assumptions $1,2$ and $3$ (resp $1,2,3$ and $7$)
were used to prove the existence (resp global asymptotic stability) of the EE.

######  Remarks 3.3.

While the term $e^{-\mu\tau}$ only shows in $\mathcal{R}_{0}.$ Comparing these
results to the ones obtained by Kaddar in the model introduced in the
beginning of the section by taking $F(S,I)=\dfrac{\beta
SI}{1+\alpha_{1}S+\alpha_{2}I}$ and $G(I)=\gamma I,$ we deduce that no
critical latent period showed in the analysis of Li & Liu (2014) and the FDE
was proven to be globally asymptotically stable only under the assumption
$\mathcal{R}_{0}>1$. Implying that the survival of the population plays an
important role in the dynamics of the disease.

###  Some of the main contributions to (5)

For this category of contribution, the aim is to deal with Limitation 4 by
introducing the spacial factor to the considered models.
Chinviriyasit & Chinviriyasit (2010) proposed the following SIR model

$\begin{cases}\begin{aligned} &\dfrac{\partial S(x,t)}{\partial
t}-\alpha\Delta S(x,t)=\mu N(x,t)-\mu S(x,t)-\beta
S(x,t)I(x,t),&x\in\Omega,t>0,\\\ &\dfrac{\partial I(x,t)}{\partial
t}-\alpha\Delta I(x,t)=-(\mu+\gamma)I(x,t)+\beta
S(x,t)I(x,t),&x\in\Omega,t>0,\\\ &\dfrac{\partial R(x,t)}{\partial
t}-\alpha\Delta R(x,t)=\gamma I+\mu R(x,t),&x\in\Omega,t>0,\\\
&\partial_{\eta}S(x,t)=\partial_{\eta}I(x,t)=\partial_{\eta}R(x,t)=0,&x\in\partial\Omega,t>0,\\\
&\text{+ positive initial conditions.}\end{aligned}\end{cases}$ (23)

With birth and death rates equal to $\mu>0,$ and $\alpha>0$ being the rate of
propagation of the individuals.
The authors proved the existence and uniqueness of the solution for a maximum
time of existence $T_{max}.$ From the boundedness of the solution it was
deduced that $T_{max}=\infty.$ The positiveness of the solution was an
immediate consequence of the maximum principle. The local asymptotic stability
of the DFE was proven (resp disproven) under the assumption
$\mathcal{R}_{0}:=\dfrac{N\beta}{\mu+\gamma}<1$ (resp $\mathcal{R}_{0}>1$),
the global asymptotic stability was carried out by the Lyapunov method. The EE
was proven to be locally (resp globally) asymptotically stable under the
condition $\mathcal{R}_{0}>1.$ The numerical simulation was restricted to the
one dimensional case in space. The authors proposed a numerical scheme based
on the forward approximation on time and a $\theta$ approximation in space.
The numerical results supported the theoretical ones, and the infected
population was proven to spread more in space in terms of increased values of
the diffusion rate $\alpha$.

######  Remarks 3.4.

The conditions
$\partial_{\eta}S(x,t)=\partial_{\eta}I(x,t)=\partial_{\eta}R(x,t)=0$ model
the confinement of the population, on the other hand, while the diffusion
rates are taken in this model to be constant, for contagious diseases, they
can depend on the infected population. When the infected increase, individuals
tend not to propagate in space to not get infected. It would be interesting to
restudy the model by considering infected depending diffusion rates. One last
remark concerns the biological feasibility of the proposed scheme, although
not done by the authors, the numerical scheme proposed is proven to be
positive invariant, that is, the iterates of the approximate solution are
always positive.

In the aim of contribution to the generalization of the incidence rate to
models falling into category $\eqref{eq:pde},$ Lotfi et al. (2014) proposed
the following generalized diffusive SIR model

$\begin{cases}\begin{aligned} &\dfrac{\partial S(x,t)}{\partial t}=d_{S}\Delta
S(x,t)+\Lambda-\mu S(x,t)\\\ &\hskip 42.67912pt-\dfrac{\beta
S(x,t)I(x,t)}{1+\alpha_{1}S(x,t)+\alpha_{2}I(x,t)+\alpha_{3}S(x,t)I(x,t)}\;x\in\Omega,t>0,\\\
&\dfrac{\partial I(x,t)}{\partial t}=d_{I}\Delta
I(x,t)-(\mu+d+\gamma)I(x,t)\\\ &\hskip 42.67912pt+\frac{\beta
S(x,t)I(x,t)}{1+\alpha_{1}S(x,t)+\alpha_{2}I(x,t)+\alpha_{3}S(x,t)I(x,t)}\;x\in\Omega,t>0,\\\
&\dfrac{\partial R(x,t)}{\partial t}=d_{R}\Delta R(x,t)+\gamma I(x,t)-\mu
R(x,t)\;x\in\Omega,t>0,\\\
&\partial_{\eta}S(x,t)=\partial_{\eta}I(x,t)=\partial_{\eta}R(x,t)=0,\;x\in\partial\Omega,t>0,\\\
&\text{+ positive initial conditions.}\end{aligned}\end{cases}$ (24)

with $d_{S},d_{I}>0$ and $d_{R}>0$ being the diffusion rates associated to the
susceptible, infected and recovered population respectively while
$\alpha_{1},\alpha_{2},\alpha_{3}>0.$

By the same analogy of the previous model (23), the authors proved the
existence and uniqueness of a global, positive and bounded solution to the
model. Knowing that
$\mathcal{R}_{0}:=\dfrac{\beta\Lambda}{\left(\mu+\alpha_{1}\Lambda\right)(\mu+d+\gamma)},$
a linearization of the system around the FDE (resp EE) proved the latter to be
asymptotically locally stable under the condition $\mathcal{R}_{0}<1$ (resp
$\mathcal{R}_{0}>1$) while it was disproven for the FDE in the case
$\mathcal{R}_{0}>1.$ Inspired by the Lyapunov function given by Hattaf et al.
(2013), the authors were able to prove the global asymptotic stability of the
FDE (resp EE) under the condition $\mathcal{R}_{0}\leq 1$ (resp
$\mathcal{R}_{0}>1$). The numerical simulation supported the theoretical
stability results, moreover, as in (23), for the case $\mathcal{R}_{0}>1$ and
a high diffusion rate $d_{I},$ the population of the infected spreads in
space.

Knowing that the diffusion plays an important role, the question of
controlling the disease in the case of fast spread arises. Taking that into
account, Adnaoui et al. (2020) considered the same previous model but were
interested into studying the existence of a spatiotemporal optimal control
representing the vaccination. The modified model was expressed by

$\begin{cases}\begin{aligned} &\dfrac{\partial S(x,t)}{\partial t}=d_{S}\Delta
S(x,t)+\Lambda-(\mu+v(x,t))S(x,t)\\\ &\hskip 42.67912pt-\dfrac{\beta
S(x,t)I(x,t)}{1+\alpha_{1}S(x,t)+\alpha_{2}I(x,t)+\alpha_{3}S(x,t)I(x,t)}\;x\in\Omega,t>0,\\\
&\dfrac{\partial I(x,t)}{\partial t}=d_{I}\Delta
I(x,t)-(\mu+d+\gamma)I(x,t)\\\ &\hskip 42.67912pt+\frac{\beta
S(x,t)I(x,t)}{1+\alpha_{1}S(x,t)+\alpha_{2}I(x,t)+\alpha_{3}S(x,t)I(x,t)}\;x\in\Omega,t>0,\\\
&\dfrac{\partial R(x,t)}{\partial t}=d_{R}\Delta R(x,t)+\gamma I(x,t)-\mu
R(x,t)+v(x,t)S(x,t)\;x\in\Omega,t>0,\\\
&\partial_{\eta}S=\partial_{\eta}I=\partial_{\eta}R=0,\;x\in\partial\Omega,t>0,\\\
&\text{+ positive initial conditions.}\end{aligned}\end{cases}$ (25)

The authors proved the mathematical and biological well posedness of the model
with vaccination which allowed to prove the existence of an optimal solution
to the optimal control problem in which the aim is the minimization (over a
biologically feasible set) of the following objective functional

$J(v)=\int_{0}^{T}\int_{\Omega}\left(\rho_{1}S(x,t)+\rho_{2}I(x,t)\right)dxdt+\frac{\alpha}{2}\|v\|_{L^{2}(\Omega\times[0,T])}^{2}.$
(26)

where $\alpha,\rho_{1},\rho_{2}>0$ are chosen weighted constants. The
characterization of the optimal control $v^{*}$ was given in terms of the
optimal state by the following formula
$v^{*}=\min\left(v^{\max},\max\left(0,\dfrac{y_{1}^{*}p_{1}-y_{1}^{*}p_{3}}{\alpha}\right)\right),$
where $p$ is the solution to the adjoint problem. The numerical simulation was
carried on by the forward-backward sweep method in which the direct problem is
solved forward in time, and the adjoint problem backwards. A choice of
parameters corresponding to a situation in which the fast spread of the
disease occurs was illustrated by simulating the model without control. In the
case with control, interesting results have been deduced. Primarily, the major
role of the immediate vaccination (1 day after the beginning of the disease)
of the population in controlling the spread of the disease over the delayed
one (20 days after).

##  Conclusion

The primal goal of this paper was to familiarize future researchers with the
existing types of compartmental models in epidemiology as well as the
contributions’ scope of each type. The paper takes into account some of the
main contributions in each category, but by no means covers them all, since
that would require a whole book. The reader should know that the ongoing and
future works are in the aim of developing more complex models combining the
reviewed types $\eqref{eq:sde},\eqref{eq:fde},\eqref{eq:dde1}$ and
$\eqref{eq:pde}.$ We refer for instance to the following recent works Khan,
Ikram, Din, Humphries & Akgul (2021); Macías‐Díaz et al. (2020); Chinnathambi
et al. (2019); Ge & Chen (2021).

Knowing all these types, it can be concluded that, while modeling a disease,
the pertinence of the chosen model is indeed based on its capacity of
maximizing the interpretation of reality but also on its ability to minimize
the difficulty of the required mathematical and numerical tools for its
analysis. And the level that one can sacrifice on each side is based on the
better understanding of the modeled disease, such understanding can be gained
in terms of contributions with experts in the modeled field, in this case,
biologists.

## References

  * (1)
  * Adnaoui et al. (2020) Adnaoui, K., Elberrai, I., Laaroussi, A. E. A. & Hattaf, K. (2020), ‘A spatiotemporal sir epidemic model two-dimensional with problem of optimal control’, Bol. Soc. Paran. Mat. 2020.
  * Ahmed et al. (2006) Ahmed, E., El-Sayed, A. & El-Saka, H. A. (2006), ‘On some routh–hurwitz conditions for fractional order differential equations and their applications in lorenz, rössler, chua and chen systems’, Phys. Lett. A 358, 1–4.
  * Altaf Khan et al. (2020) Altaf Khan, M., Ismail, M., Ullah, S. & Farhan, M. (2020), ‘Fractional order sir model with generalized incidence rate’, AIMS Math. 5, 1856–1880.
  * Anderson & May (1992) Anderson, R. M. & May, R. M. (1992), Infectious diseases of humans : dynamics and control, Oxford University Press.
  * Beddington (1975) Beddington, J. R. (1975), ‘Mutual interference between parasites or predators and its effect on searching efficiency’, The Journal of Animal Ecology 44, 331.
  * Bellen & Zennaro (2013) Bellen, A. & Zennaro, M. (2013), Numerical methods for delay differential equations, Oxford University Press.
  * Bellman (2008) Bellman, R. (2008), Stability theory of differential equations, Dover Publications.
  * Boukhouima et al. (2020) Boukhouima, A., Hattaf, K., Lotfi, E. M., Mahrouf, M., Torres, D. F. & Yousfi, N. (2020), ‘Lyapunov functions for fractional-order systems in biology: Methods and applications’, Chaos Solitons Fractals 140, 110224.
  * Buonomo & Rionero (2010) Buonomo, B. & Rionero, S. (2010), ‘On the lyapunov stability for sirs epidemic models with general nonlinear incidence rate’, Appl. Math. Comput 217, 4010–4016.
  * Capasso & Serio (1978) Capasso, V. & Serio, G. (1978), ‘A generalization of the kermack-mckendrick deterministic epidemic model’, Math. Biosci. Eng. 42, 43–61.
  * Casagrandi et al. (2006) Casagrandi, R., Bolzoni, L., Levin, S. A. & Andreasen, V. (2006), ‘The sirc model and influenza a’, Math. Biosci. 200, 152–169.
  * Chinnathambi et al. (2019) Chinnathambi, R., Rihan, F. A. & Alsakaji, H. J. (2019), ‘A fractional‐order model with time delay for tuberculosis with endogenous reactivation and exogenous reinfections’, Math. Methods Appl. Sci. 44, 8011–8025.
  * Chinviriyasit & Chinviriyasit (2010) Chinviriyasit, S. & Chinviriyasit, W. (2010), ‘Numerical modelling of an sir epidemic model with diffusion’, Appl. Math. Comput. 216, 395–409.
  * Crowley & Martin (1989) Crowley, P. H. & Martin, E. K. (1989), ‘Functional responses and interference within and between year classes of a dragonfly population’, Journal of the North American Benthological Society 8, 211–221.
  * Cui et al. (2020) Cui, Q., Du, Q. & Wang, L. (2020), ‘Global dynamics of a generalized sirs epidemic model with constant immigration’, Math. Probl. Eng. 2020, 1–9.
  * de Oliveira & Tenreiro Machado (2014) de Oliveira, E. C. & Tenreiro Machado, J. A. (2014), ‘A review of definitions for fractional derivatives and integral’, Math. Probl. Eng. 2014, 1–6.
  * Diethelm et al. (2002) Diethelm, K., Ford, N. & Freed, A. (2002), ‘A predictor–corrector approach for the numerical solution of fractional differential equations’, Nonlinear Dynam. 29, 3–22.
  * Driesschea & Watmough (2002) Driesschea, P. & Watmough, J. (2002), ‘Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission’, Math. Biosci. 180, 29–48.
  * El Koufi et al. (2019) El Koufi, A., Adnani, J., Bennar, A. & Yousfi, N. (2019), ‘Analysis of a stochastic sir model with vaccination and nonlinear incidence rate’, Int. J. Differ. Equ. 2019, 1–9.
  * El-Shahed & Alsaedi (2011) El-Shahed, M. & Alsaedi, A. (2011), ‘The fractional sirc model and influenza a’, Math. Probl. Eng. 2011, 1–9.
  * Evans (2016) Evans, L. C. (2016), Partial differential equations, Providence, Ri American Math. Soc. Himayatnagar, Hyderabad (Telangana) Universities Press (India.
  * Gantmakher (2005) Gantmakher, F. (2005), Applications of the theory of matrices, Dover Publications.
  * Ge & Chen (2021) Ge, F. & Chen, Y. (2021), ‘Optimal vaccination and treatment policies for regional approximate controllability of the time-fractional reaction–diffusion sir epidemic systems’, ISA Trans. 115, 143–152.
  * Griffiths & Higham (2011) Griffiths, D. F. & Higham, D. J. (2011), Numerical methods for ordinary differential equations : initial value problems, Springer.
  * Hale (1977) Hale, J. K. (1977), Theory of functional differential equations, Springer-Verlag.
  * Hattaf et al. (2013) Hattaf, K., Lashari, A., Louartassi, Y. & Yousfi, N. (2013), ‘A delayed sir epidemic model with a general incidence rate’, Electron. J. Qual. Theory Differ. Equ. pp. 1–9.
  * Hattaf & Yousfi (2020) Hattaf, K. & Yousfi, N. (2020), ‘Global stability for fractional diffusion equations in biological systems’, Complexity 2020, 1–6.
  * Heesterbeek et al. (1996) Heesterbeek, J., Mart De Yong, C. & Diekmann, O. (1996), ‘How does transmission of infection depend on population size?’, Pub. of the Newton Ins. pp. 84–94.
  * Higham (2001) Higham, D. J. (2001), ‘An algorithmic introduction to numerical simulation of stochastic differential equations’, SIAM Rev. 43, 525–546.
  * Hsieh (2013) Hsieh, P.-F. (2013), Basic theory of ordinary differential equations., Springer.
  * Jiang et al. (2011) Jiang, D., Yu, J., Ji, C. & Shi, N. (2011), ‘Asymptotic behavior of global positive solution to a stochastic sir model’, Math. Comput. Modelling 54, 221–232.
  * Kaddar (2009) Kaddar, A. (2009), ‘On the dynamics of a delayed sir epidemic model with a modified saturated incidence rate’, J. Differential Equations 133, 1–7.
  * Kermack & McKendrick (1927) Kermack, W. O. & McKendrick, A. G. (1927), ‘A contribution to the mathematical theory of epidemics’.
  * Khan, Ikram, Din, Humphries & Akgul (2021) Khan, A., Ikram, R., Din, A., Humphries, U. W. & Akgul, A. (2021), ‘Stochastic covid-19 seiq epidemic model with time-delay’, Results Phys. 30, 104775.
  * Khan et al. (2019) Khan, T., Ullah, Z., Ali, N. & Zaman, G. (2019), ‘Modeling and control of the hepatitis b virus spreading using an epidemic model’, Chaos Solitons Fractals 124, 1–9.
  * Khan, Alaoui, Zeb, Tilioua & Djilali (2021) Khan, Z. A., Alaoui, A. L., Zeb, A., Tilioua, M. & Djilali, S. (2021), ‘Global dynamics of a sei epidemic model with immigration and generalized nonlinear incidence functional’, Results Phys. 27, 104477.
  * Kloeden & Platen (1999) Kloeden, P. E. & Platen, E. (1999), Numerical solution of stochastic differential equations, Springer.
  * Korobeinikov (2006) Korobeinikov, A. (2006), ‘Lyapunov functions and global stability for sir and sirs epidemiological models with non-linear transmission’, Bull. Math. Biol. 68, 615–626.
  * Korobeinikov & Maini (2005) Korobeinikov, A. & Maini, P. K. (2005), ‘Non-linear incidence and stability of infectious disease models’, Math. Med. Biol. 22, 113–128.
  * Lenhart & Workman (2007) Lenhart, S. & Workman, J. T. (2007), Optimal control applied to biological models, Chapman and Hall/Crc.
  * Li & Liu (2014) Li, M. & Liu, X. (2014), ‘An sir epidemic model with time delay and general nonlinear incidence rate’, Abstr. Appl. Anal. 2014, 1–7.
  * Li et al. (2017) Li, T., Zhang, F., Liu, H. & Chen, Y. (2017), ‘Threshold dynamics of an sirs model with nonlinear incidence rate and transfer from infectious to susceptible’, Appl. Math. Lett. 70, 52–57.
  * Lotfi et al. (2014) Lotfi, E. M., Maziane, M., Hattaf, K. & Yousfi, N. (2014), ‘Partial differential equations of an epidemic model with spatial diffusion’, Int. J. Partial Differ. Equ. 2014, 1–6.
  * Lu (2009) Lu, Q. (2009), ‘Stability of sirs system with random perturbations’, Phys. A: Stat. Mech. Appl. 388, 3677–3686.
  * Macías‐Díaz et al. (2020) Macías‐Díaz, J. E., Ahmed, N., Jawaz, M., Rafiq, M. & Aziz ur Rehman, M. (2020), ‘Design and analysis of a discrete method for a time‐delayed reaction–diffusion epidemic model’, Math. Methods Appl. Sci. 44, 5110–5122.
  * Marinov et al. (2014) Marinov, T. T., Marinova, R. S., Omojola, J. & Jackson, M. (2014), ‘Inverse problem for coefficient identification in sir epidemic models’, Comput. Math. with Appl. 67, 2218–2227.
  * McCluskey (2010) McCluskey, C. C. (2010), ‘Global stability for an sir epidemic model with delay and nonlinear incidence’, Nonlinear Anal. Real World Appl. 11, 3106–3109.
  * Medvedeva et al. (2020) Medvedeva, M., Simos, T. E., Tsitouras, C. & Katsikis, V. (2020), ‘Direct estimation of sir model parameters through second‐order finite differences’, Math. Methods Appl. Sci. 44, 3819–3826.
  * Mouaouine et al. (2018) Mouaouine, A., Boukhouima, A., Hattaf, K. & Yousfi, N. (2018), ‘A fractional order sir epidemic model with nonlinear incidence rate’, Adv. Differ. Equ. 2018.
  * Ndaïrou et al. (2020) Ndaïrou, F., Area, I., Nieto, J. J. & Torres, D. F. (2020), ‘Mathematical modeling of covid-19 transmission dynamics with a case study of wuhan’, Chaos Solitons Fractals 17, 2708–2724.
  * Odibat & Momani (2008) Odibat, Z. & Momani, M. (2008), ‘An algorithm for the numerical solution of differential equations of fractional order’, J. Appl. Math. Inform. 26, 15–27.
  * Oksendal (2013) Oksendal, B. K. (2013), Stochastic differential equations : an introduction with applications, Springer.
  * Owolabi (2020) Owolabi, K. M. (2020), Numerical Methods For Fractional Differentiation, Springer.
  * Ozalp & Demi̇rci̇ (2011) Ozalp, N. & Demi̇rci̇, E. (2011), ‘A fractional order seir model with vertical transmission’, Math. comput. model. 54, 1–6.
  * Protter & Weinberger (1984) Protter, M. H. & Weinberger, H. F. (1984), Maximum principles in differential equations., Springer.
  * Rao (2014) Rao, F. (2014), ‘Dynamics analysis of a stochastic sir epidemic model’, Abstr. Appl. Anal. 2014, 1–9.
  * Sigdel & McCluskey (2014) Sigdel, R. P. & McCluskey, C. C. (2014), ‘Global stability for an sei model of infectious disease with immigration’, Appl. Math. Comput. 243, 684–689.
  * Syafruddin & Noorani (2013) Syafruddin, S. & Noorani, M. S. M. (2013), ‘Lyapunov function of sir and seir model for transmission of dengue fever disease’, Int. J. Simul. Model. 8, 177.
  * Tarantola (1987) Tarantola, A. (1987), Inverse problem theory: Methods for data fitting and model parameter estimation, Elsevier Science Pub. Co. Inc.,New York, Ny.
  * Tornatore et al. (2005) Tornatore, E., Maria Buccellato, S. & Vetro, P. (2005), ‘Stability of a stochastic sir system’, Phys. A: Stat. Mech. Appl. 354, 111–126.
  * Tröltzsch (2010) Tröltzsch, F. (2010), Optimal control of partial differential equations : theory, methods, and applications, American Mathematical Society.
  * Vrabie (2003) Vrabie, I. I. (2003), Co-semigroups and applications, Elsevier Science.
  * Wang (2006) Wang, W. (2006), ‘Backward bifurcation of an epidemic model with treatment’, Math. Biosci. 201, 58–71.
  * Wang & Ruan (2004) Wang, W. & Ruan, S. (2004), ‘Bifurcations in an epidemic model with constant removal rate of the infectives’, J. Math. Anal. 291, 775–793.
  * Workman & Lenhart (2007) Workman, J. T. & Lenhart, S. (2007), Numerical Solution of Stochastic Differential Equations, Taylor and Francis Group.
  * Xiao & Ruan (2007) Xiao, D. & Ruan, S. (2007), ‘Global analysis of an epidemic model with nonmonotone incidence rate’, Math. Biosci. 208, 419–429.
  * Xing & Li (2021) Xing, Y. & Li, H.-X. (2021), ‘Almost periodic solutions for a svir epidemic model with relapse’, Math. Biosci. Eng. 18, 7191–7217.
  * Xu & Ma (2009) Xu, R. & Ma, Z. (2009), ‘Global stability of a sir epidemic model with nonlinear incidence rate and time delay’, Nonlinear Anal. Real World Appl. 10, 3175–3189.
  * Zhonghua & Yaohong (2009) Zhonghua, Z. & Yaohong, S. (2009), ‘Qualitative analysis of a sir epidemic model with saturated treatment rate’, J. Appl. Math. Comput. 34, 177–194.