Next Article in Journal
Pathway to Fractional Integrals, Fractional Differential Equations, and Role of the H-Function
Next Article in Special Issue
Necessary and Sufficient Criteria for a Four-Weight Weak-Type Maximal Inequality in the Orlicz Class
Previous Article in Journal
One Turán Type Problem on Uniform Hypergraphs
Previous Article in Special Issue
Concerning Transformations of Bases Associated with Unimodular diag(1, −1, −1)-Matrices
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Review

Exploring the Landscape of Fractional-Order Models in Epidemiology: A Comparative Simulation Study

by
Ritu Agarwal
1,†,
Pooja Airan
1,† and
Ravi P. Agarwal
2,*
1
Department of Mathematics, Malaviya National Institute of Technology, Jaipur 302017, India
2
Department of Mathematics and Systems Engineering, Florida Institute of Technology, Melbourne, FL 32901, USA
*
Author to whom correspondence should be addressed.
These authors contributed equally to this work.
Axioms 2024, 13(8), 545; https://doi.org/10.3390/axioms13080545
Submission received: 21 June 2024 / Revised: 31 July 2024 / Accepted: 9 August 2024 / Published: 11 August 2024
(This article belongs to the Special Issue Theory of Functions and Applications II)

Abstract

:
Mathematical models play a crucial role in evaluating real-life processes qualitatively and quantitatively. They have been extensively employed to study the spread of diseases such as hepatitis B, COVID-19, influenza, and other epidemics. Many researchers have discussed various types of epidemiological models, including deterministic, stochastic, and fractional order models, for this purpose. This article presents a comprehensive review and comparative study of the transmission dynamics of fractional order in epidemiological modeling. A significant portion of the paper is dedicated to the graphical simulation of these models, providing a visual representation of their behavior and characteristics. The article further embarks on a comparative analysis of fractional-order models with their integer-order counterparts. This comparison sheds light on the nuances and subtleties that differentiate these models, thereby offering valuable insights into their respective strengths and limitations. The paper also explores time delay models, non-linear incidence rate models, and stochastic models, explaining their use and significance in epidemiology. It includes studies and models that focus on the transmission dynamics of diseases using fractional order models, as well as comparisons with integer-order models. The findings from this study contribute to the broader understanding of epidemiological modeling, paving the way for more accurate and effective strategies in disease control and prevention.

Graphical Abstract

1. Introduction

The word epidemiology comes from the Greek words epi, meaning on or upon, demos, meaning people, and logos, meaning the study of. Hence, ‘Epidemiology’ is defined as the study of the distribution and determinants of health-related states or events in specified populations, and the application of this study to the control of health problems [1]. Epidemiology is a distinct branch of research that identifies the distribution of diseases, their causes, and ways of controlling them. To conduct this, epidemiologists must understand how political, social, and scientific factors interact to increase disease risk. It is a very complex science since it must take into account a variety of factors that affect human diseases, including microorganisms, human social or travel dynamics, vaccination, quarantine period, immunity, climate, and others. Epidemiological research will undoubtedly continue to serve as the basis for formulating public health strategies in the future [2].
Reviewing the SIR (Susceptible-Infected-Recovered) model in the context of existing knowledge is important due to several key reasons. Since its development in the 1927s, the SIR model has served as a fundamental epidemiological tool, influencing our understanding of the dynamics of infectious diseases. Validating its accuracy with real-world data from past outbreaks like influenza and COVID-19 demonstrates its reliability and relevance. To understand the model’s limitations and the specific situations in which it works best, it is essential to examine its underlying assumptions, which include a closed population and homogeneous mixing. The SIR model has evolved into extensions like SIRS and SEIR (Susceptible-Exposed-Infected-Recovered), which address its limitations and adapt to different disease scenarios. Reviewing these modifications provides insights into its adaptability and ongoing development. Reviewing the fractional SIR model, which incorporates fractional calculus into the traditional SIR model, allows a more nuanced representation of disease dynamics by incorporating memory effects. The review of the fractional models will help us to better understand the many operators that may be applied to generalize these models as well as the differences between the integer and fractional order models.
The objective of this review is to evaluate the effectiveness and accuracy of the fractional SIR model in modeling infectious disease dynamics. Its specific objectives are to evaluate the applicability of the fractional model to different infectious diseases, compare the performance of the fractional SIR model with the classical SIR model using available data, and investigate the theoretical developments and real-world applications of using fractional calculus in epidemiological modeling.

1.1. Epidemic Modeling

The term ‘epidemic modeling’ refers to a collection of methods for analyzing the transmission of communicable infections among host populations using mathematical, statistical, and computational tools. It does so by utilizing information and theories that outline the demographic changes, environmental factors, possibilities for disease transmission, and health effects of various conditions [3]. Some of the notable reasons to model the dynamics of the epidemics are as follows:
  • Direct observation can be used to evaluate the impact of interventions.
  • The dynamics of epidemics are non-linear, which presents a second challenge [4] and the significance of its stochastic elements (beginning and end of a wave, super-propagation) [5] simple proportionality criteria are unable to forecast the trajectory of an epidemic; instead, more computation carried out by models are necessary.
  • A third difficulty is that accurate descriptions of an epidemic would require complete knowledge of the entire population.
Modeling is useful for assessing the current state, predicting the future, and fully understanding the past course of an epidemic. Modeling in infectious disease epidemiology today depends on a variety of methods based on decades of developments at the interface of demography, medicine, and biology. The vast majority of models adopt a compartmental description of the disease as a sequence of different stages encountered upon infection to cure or death [3]. An example of this model is SIR.
The SIR model is helpful in analyzing how a virus spreads among diverse populations, especially during times of increased infection rates. The model can provide information about how the virus spreads over time that is not available from the data alone. The SIR model can estimate how long the virus will persist within a population and predict future infection and death rates. Moreover, it provides information about reducing the virus’s impact, which is difficult to determine from recorded data alone.
The SIR model can be used to predict an outbreak of an infectious illness like measles, COVID-19, Hepatitis B, etc. It is especially useful for tracking epidemic trends because it can be adjusted to take into account rises and is based on recorded data.
One reason to examine infectious diseases is to refine control and then completely demolish the infection from the population. The dynamics of infectious virus infections are better understood through mathematical modeling, which also provides control measures [6]. It makes predictions about how a disease will spread over time by taking into consideration the key variables that influence how a disease develops, such as transmission and recovery rates. Mathematical models are now important tools for understanding how infectious diseases spread and are controlled. According to epidemiological models, there are three categories of people in the population: those who are susceptible, those who are infected, and those who have recovered. For some diseases, the principle of mass action [7] governs how an epidemic spreads.
The Hepatitis B virus is a contagious illness that affects people all over the world and causes a large number of infections. Along with other illnesses, the virus also results in liver damage. In humans, there are two possible transmission routes: horizontal and vertical. Many HBV liver infection and death cases were reported in 2015, and many nations throughout the world, including China, are still facing this issue with numerous cases and fatalities.
Various scientists have created mathematical models based on the properties of the Hepatitis B virus. For instance, HBV with control strategies [8], information related vaccination [9], effects of vaccination on HBV infected peoples [10], HBV with age structure analysis [11], HBV virus as a multi-group model [12], the mathematical model for HBV with cost-effective analysis [13], etc.
A novel strain of Coronavirus (SARS-CoV-2) was recognized in Wuhan, Hubei Province, China in December 2019 causing a severe and potentially fatal respiratory syndrome, i.e., COVID-19. Since then, the World Health Organization (WHO) declared it a pandemic on March 11, which has spread around the world [14,15,16]. Coronavirus is considered to be the most challenging virus in the community of flu virus. Millions of individuals were forced by national governments to live in isolation. The sickness is spreading fast around the world. The greatest methods for controlling the pandemic in the absence of an effective vaccine include wearing a face mask, isolating oneself, and avoiding close contact with others.
Mathematical models are essential for independently determining key parameters like disease transmission, recovery, and mortality for different populations. To stop the disease from spreading, various countries have taken a variety of actions. The construction of mathematical models that involve parameters and various compartments gives us information about the transmission dynamics of infectious diseases.

1.2. Historical Review of the Epidemic Model

The modeling of infectious diseases is a tool that has been used to study the mechanisms by which diseases spread, to guess the future course of an outbreak, and to evaluate strategies to control an epidemic [17].
The outbreak and spread of diseases have been questioned and studied for many years. In fact, the study of infectious disease data started with the work of John Graunt (1620–1674) in 1662 [18]. He looked at the causes of death and provided a method for calculating the relative risks of dying from various illnesses. He also provided the first theoretical framework for competing hazards.
In 1760, Daniel Bernoulli (1700–1782) [19] claimed to develop the first mathematical epidemiological model on smallpox vaccination. In the eighteenth century, smallpox was a common disease. To produce lifelong immunity against smallpox, inoculation with a mild strain was introduced but with a small risk of infection and death. However, his study has gained more recognition in the actuarial field than the epidemiological literature, while receiving a largely positive response. However, Dietz et al. [20] more recently generalized his strategy.
Then, in 1855, John Snow was able to notice the temporal and spatial-temporal pattern of cholera cases epidemics in London [21]. The source of infection was the street water pump. In 1873, a similar understanding of the spread of typhoid was achieved by William Budd [22]. In 1840, William Farr [23] analyzed statistical data to uncover the underlying principles governing the rise and decline of epidemics.
In 1906, Hamer proposed that the transmission of infection is affected by the number of susceptible populations and the number of infectious populations [24]. He proposed a mass action law (the number of new infections when susceptible individuals (S) come in contact with the infectious individuals (I) is proportional to SI) for the rate of new infections.
A very general model of infectious disease was defined using the first stage of compartmental models in epidemiology. In the compartmental model, the population is divided into different compartments representing their status of diseases, demographics, and other factors, and where mathematical equations are used to model transitions between different compartments. The compartmental model is a very common modeling technique, often applied to the mathematical modeling of infectious diseases, chemical kinetics, and others. In the SIR model, the population is divided into compartments with labels, for example, S, I, or R, (Susceptible, Infectious, or Recovered). Individuals may move from one compartment to another. The sequence of the labels indicates the flow between compartments; for instance, SEIS represents the progression from susceptible to exposed to infectious, and then back to susceptible.
In 1902, Dr. Ross won the second Nobel Prize in Medicine for his work illustrating the dynamics of malaria transmission between mosquitoes and humans. However, Ross provided a simple compartmental model [25] including mosquitoes and humans, showing that it would be enough to lower the mosquito population below a critical threshold. Since that time, the concept of basic reproduction number was first introduced in mathematical epidemiology.
In the early 20th century, Ross [25] in 1916 and Hudson [26] in 1917, Kermack and Mckendrick [27] in 1927 and George Kendall [28] in 1956 proposed such models. The research Ross did on malaria is one example of such a model.
The mathematical model to describe an influenza epidemic was first developed by Kermack and Mckendrick, which is popularly known as the Susceptible-Infectious-Recovered (or SIR) model. Kermack and McKendrick’s theoretical papers on infectious disease models published between 1927 and 1933 had a significant impact on the development of mathematical epidemiology models. The SIR model effectively predicted outbreak behaviors that closely matched those observed in many recorded epidemics [29].
The two primary methods of disease transmission are the horizontal and vertical transmissions:
  • Horizontal transmission: Horizontal transmission occurs when a parasite moves from an infected to an uninfected individual, whether by direct contact or an infectious particle [30]. When illness is spread via direct contact, this implies that the infected party spreads the disease to another through close contact. This can include touching and activities that involve the exchange of bodily fluids (such as kissing or sexual activities). Indirect contact, such as breathing in the same atmospheric air as the infected person, can also result in the spread of an illness. When the sick person is coughing and sneezing, viruses are forcibly coming out of their lungs, ready to be transmitted to, and infect, the next susceptible individual in close proximity.
  • Vertical transmission: Vertical transmission occurs when an infected individual reproduces (either sexually or asexually), giving rise to progeny that also harbors the infection [30]. An example of a virus that can be spread via vertical transmission is HIV (human immunodeficiency virus). This virus is capable of crossing the trans-placental barrier, meaning that if the mother is infected, she may spread this virus to her unborn child. HIV may also be spread to a newborn via breastfeeding, as lactating mothers can spread the pathogen to their children via breast milk containing viral particles. In addition to vertical transmission, HIV may also be spread via horizontal transmission. This virus is spread by bodily fluids, such as blood, semen, and vaginal secretions.

1.3. Types of Epidemic Models

Epidemic models are divided into two types: deterministic and stochastic (random) models.
  • Deterministic Models Deterministic or Compartmental mathematical models are often used in dealing with diseases that spread in large populations such as tuberculosis. In this model, individuals in the population are allotted to different subgroups or compartments, each representing a specific stage of the epidemic [31]. The transition rates between classes are mathematically represented as derivatives, leading to the formulation of the model using differential equations. When developing such models, it is assumed that the population size within each compartment changes continuously over time, and the epidemic process follows deterministic rules. Changes in the population of a compartment can be calculated using the history only that was used to develop the model [29].
  • Stochastic Models Stochastic means having a random variable. A stochastic model uses random variables in one or more inputs over time to estimate probability distributions of potential outcomes. These models rely on random fluctuations in exposure risk, disease dynamics, and other factors related to illness. Stochastic methods can be used to assess disease spread at the individual level within small or large populations [32].

1.4. Fractionalization of the Epidemiological Models

A number of fractional order epidemiological models have been studied by the researchers. Creating a comprehensive list of such contributions is a big task. For the sake of brevity, we have made efforts to include a list of important contributions, which may be referred to for a more detailed study.
Rafei et al. [33] formulated a SIR model which was solved by the homotopy perturbation method in 2007. Doungmo Goufo [34] derived the SIR model in epidemiology with non-linear incidence and fractional derivative order using Caputo derivative in 2014. Okyere et al. [35] studied the fractional order SIR model with constant population and the mathematical formulation of the non-integer order initial value problem based on the famous fractional order Caputo derivative in 2016. A fractional order SIS pandemic model was created by Hassouna et al. [36] in 2018. In this paper, mathematical formulation is based on the fractional Caputo derivative.
Rostamy and Mottaghi [37] formulated the SIRV model in 2016. They extended the SIR model with vaccination into a fractional-order model by using a system of fractional ordinary differential equations using the Caputo derivative of order ρ ( 0 , 1 ] . Oke [38] formulated a SIRV epidemic model with non-linear Force of Infection and Treatment in 2019 and discussed the stability analysis. Qureshi and Jan [39] formulated the SIRV model of the measles epidemic with optimized fractional order under the Caputo differential operator in 2021.
Koca [40] formulated the SIRD model for the spread of Ebola virus and used the method Atangana–Baleanu derivative in Caputo sense (ABC) in 2018. Dokuyucu and Dutta [41] give a fractional order SIRD model for Ebola Virus with the new Caputo fractional derivative without a singular kernel in 2020. Singh [42] developed the SIRD epidemic model for the disease Ebola in 2020. In this paper, the Grunwald–Letnikov (GL) fractional derivative has been used for the fractional Ebola virus model.
Ozalp and Demirci [43] gave a fractional order SEIR model with vertical transmission within a nonconstant population in 2011. Area [44] formulated an SEIR fractional order model for the disease Ebola in 2015. They used the fractional Riemann–Liouville derivative of order ρ . Almeida [45], formulated the SEIR model using the Caputo fractional derivative in 2018.
Casagrandi et al. [46] developed an SIRC model for the influenza A disease in 2006. El and Alsaedi [47] described a SIRC model for the influenza A disease using the Caputo fractional derivative in 2011. Agarwal et al. [48] formulated an epidemic model on the Hepatitis B Virus Infection in 2020. Arfan et al. [49] formulated a fractional-type modified SEIR model under the Atangana–Baleanu–Caputo (ABC) derivative in 2021. Agarwal et al. [50] formulated a fractional-order mathematical model for the cell cycle of a tumor cell in 2020.
Asfour and Ibrahim [51] formulated the MSEIR epidemic model using the differential Riemann–Liouville and the Caputo fractional derivative in the epidemic in 2015. Almeida et al. [52] formulated an epidemiological MSEIR model using Caputo fractional derivative in 2019. Qureshi and Jan [39] formulated an MSEIR model of the chickenpox disease with fractional derivatives in 2021. In this paper, they used a fractional order model under the Caputo fractional, under the Caputo–Fabrizio and under the Atangana–Baleanu–Caputo. Khan et al. [53] gave Hepatitis B Virus transmission with fractional analysis in 2022 and used the Caputo–Fabrizio (C–F) operator.
Qureshi and Atangana [54] formulated a SITR model in 2020 for diarrhea disease using fractal fractional differentiation under the Riemann–Liouville operator [55]. Atangana [56] defined the model for COVID-19 using the fractal-fractional operator in 2020. Baleanu et al. [57] formulated a fractional differential equation model for the COVID-19 transmission by using the Caputo–Fabrizio derivative in 2020. Islam et al. [58] discussed the SEIR stochastic model for fractional order of the disease measles in 2020. Pandey et al. [59] described a model for crowding effects on the COVID-19 virus using the Atangana–Baleanu (AB) fractional operator in 2022. Owolabi and Pindza [60] studied a fractional order non-linear epidemic model in 2022 under the Caputo operator and fixed point theory to describe the evolution of tuberculosis in the population. The applications of fractional order models are illustrated in the study of bone mineralization in [61,62,63].
The article attempts to give a thorough introduction to the fractional SIR model and its extensions, emphasizing certain important aspects that improve our understanding of the dynamics of infectious diseases. The focus is on determining the equilibrium points and calculating the reproduction number, both of which are crucial for analyzing the stability and dynamics of disease spread. We proceed by presenting the theorems related to the stability of the system, along with those addressing the existence and uniqueness of solutions. Furthermore, detailed numerical and graphical simulations employing the Lagrange two-step method are conducted to illustrate the models’ behavior under various conditions. By applying the same parametric values across different scenarios, we validate the models and ensure consistency in their behavior. This method shows how the dynamics of the fractional-order and integer-order models differ from one another and offers a thorough comparison.
This article is organized into various sections, as follows: Section 2 introduces the mathematical preliminaries in which we discuss the fractional operators and numerical schemes for solving the Caputo fractional differential equations. Section 3 describes the mathematical models for epidemics which include qualitative analysis of the model. Section 4 focuses on the Kermack–McKendrick SIR model, detailing its formulation and significance; Section 5 explores extensions of the SIR model, such as SEIR, SIRD, SVIR, and MSIR, etc., along with equilibrium points, reproduction number, numerical and graphical simulations. In Section 6, the SIR model with delay is discussed. Section 7 introduces stochastic SIR models to account for randomness and variability in disease spread. Section 8 describes the methodology adopted for the review and flow diagrams to illustrate the research process. Finally, Section 9 concludes the paper by summarizing key findings and suggesting future research directions.

2. Mathematical Preliminaries

2.1. Fractional Operators

In the last thirty years, the study of fractional calculus, particularly the calculus of derivatives and integrals of arbitrary order (both real and complex), has expanded due to its obvious applications in a number of scientific and engineering domains. Without a doubt, this field of mathematics has offered a variety of tools that could be beneficial in solving differential and integral equations, as well as a variety of other issues involving special functions of mathematical physics and all of their expansions and generalizations.
It is important to note that traditional mathematical models with integer order derivatives frequently fall short. The beauty of this topic is in its fractional derivatives and integrals, which are not just a local property but also take into account history and locally dispersed effects. That is why this subject better translates nature’s reality. Non-integer order fractional derivatives and integrals can be utilized to explain processes with memory. If at each time t, the output of a system depends only upon the input at time t, such a system is called a memoryless system while to find the current value of the output of the system needs to remember previous values of the input, non-memoryless or memory systems.
We shall require the following definitions in the sequel:
Definition 1 
(Riemann–Liouville fractional integral, page 45 of [64]). Let R e ( ρ ) > 0 , and let f be a piecewise continuous function on J = ( 0 , ) and integrable on any finite sub-interval of J = [ 0 , ) . Then, for x > 0 , we call
  0 R L I x ρ f ( x ) =   0 R L D x ρ f ( x ) = 1 Γ ( ρ ) 0 x ( x t ) ρ 1 f ( t ) d t ,
the Riemann–Liouville fractional integral of f of order ρ.
Equivalently, in the convolution form,
  0 R L D x ρ f ( x ) = 1 Γ ( ρ ) f ( x ) x ρ 1 , x > 0 , R e ( ρ ) > 0 .
Laplace transform (page 69 of [64]) of the Riemann–Liouville fractional integral is
L {   0 R L D x ρ f ( x ) } = 1 Γ ( ρ ) L { t ρ 1 } L { f ( x ) } , ( using convolution theorem ) = s ρ F ( s ) , ρ > 0 ,
where F ( s ) is the Laplace transform of the function f ( x ) .
Further, the Riemann–Liouville fractional derivative associated with the Riemann–Liouville fractional integral (1) is defined as,
Definition 2 
(Riemann–Liouville fractional derivative). Let f be a piecewise continuous function on J = ( 0 , ) . Then the Riemann–Liouville fractional derivative of order ρ, R e ( ρ ) > 0 is,
  0 R L D x ρ f ( x ) = 1 Γ ( n ρ ) d d x n 0 x ( x t ) n ρ 1 f ( t ) d t , n 1 ρ < n , n N , x > 0 = d n d x n   0 R L I x n ρ f ( x ) .
Equivalently, in the convolution form,
  0 R L D x ρ f ( x ) = 1 Γ ( n ρ ) d d x n ( f ( x ) x n ρ 1 ) , x > 0 , ρ > 0 .
Laplace transform of the Riemann–Liouville fractional derivative is
L {   0 R L D x ρ f ( x ) } = s ρ L { f ( x ) } k = 0 n 1 s k f ( ρ k 1 ) ( 0 ) .
If n = 1 , that is 0 ρ < 1 , then
  0 R L D x ρ f ( x ) = 1 Γ ( 1 ρ ) d d x 0 x ( x t ) ρ f ( t ) d t ,
and
L {   0 R L D x ρ f ( x ) } = s ρ F ( s ) f ( ρ 1 ) ( 0 ) ,
where F ( s ) is the Laplace transform of the function f ( x ) .
While studying the earthquake problem in 1967, Caputo [65] defined and gave the definition which is now known as the Caputo fractional derivative. This operator resolves the limitation of the Riemann–Liouville fractional derivative and also it permits the utilization of classical initial conditions associated with differential equations of integer order.
In the following term, Caputo introduced the concept of the fractional derivative:
Definition 3 
(Caputo Fractional Derivative). Suppose that ρ > 0 , x > 0 , and the fractional derivative of real-valued function f ( x ) AC n [ a , b ]
  0 C D x ρ f ( x ) = 1 Γ ( n ρ ) 0 x ( x t ) n ρ 1 f ( n ) ( t ) d t , n 1 ρ < n d n f d x n , ρ = n
is called the Caputo fractional derivative of order ρ.
The Caputo fractional differential operator can also be written as
  0 C D x ρ f ( x ) =   R L   0 I x n ρ f ( n ) ( x ) .
Equivalently, in the convolution form,
  0 C D x ρ f ( x ) = 1 Γ ( n ρ ) f ( n ) ( x ) x n ρ 1 , x > 0 .
The Laplace transform [66] of the Caputo fractional derivative is
L {   0 C D x ρ f ( x ) } = s ρ L { f ( x ) } k = 0 n 1 s ρ k 1 f ( k ) ( 0 ) .
Following relations between the Caputo derivative and the R-L integral are important for solving the fractional order differential equations (page 95–96 of [66]):
(   0 C D x ρ   0 R L I x ρ f ) ( x ) = f ( x ) , and (   0 R L I x ρ   0 C D x ρ f ) ( x ) = f ( x ) k = 0 n 1 f ( k ) ( 0 ) k ! ( x ) k
under the suitable conditions of existence.
The application of fractional calculus to solve problems in the real world is becoming more and more important. The memory effect and inherited qualities can be better understood using fractional-order models than integer-order models. Fractional calculus has been increasingly significant in a number of disciplines recently, including signal and image processing, chemistry, biology, mechanics, and electricity. A table of correspondence between the type of phenomenon and fractional operators is given in [67].

2.2. Numerical Schemes for Solving Caputo Fractional Differential Equation

When the number of variables is less, the system of non-linear first-order differential equations can be solved using analytical, semi-analytical and numerical methods. But as the number of variables increases, the complexity of the problem increases, and handling the solution with analytical, and semi-analytical methods becomes difficult. In this paper, we have used two numerical techniques: Euler’s method and Lagrange’s two-point method to solve the various epidemic models.
Euler’s Method
Let us consider the initial value problem involving the Caputo fractional derivative of order ρ
  0 C D x ρ y ( x ) = H ( x , y ( x ) ) y ( 0 ) = y 0 , where 0 < ρ 1 , x > 0 .
Let [ 0 , a ] be the interval over which we find the solution of the Equation (14). Subdivide the interval [ 0 , a ] into n sub-intervals [ x n , x n + 1 ] of equal length h = a n by using the nodes x n = n h , for n = 0 , 1 , 2 , , n .
Assuming y ( x ) , D ρ y ( x ) , D 2 ρ y ( x ) are continuous on the interval [ 0 , a ] and using the generalized Taylor’s formula to expand y ( x ) about x 0 , we obtain
y ( x ) = y ( x 0 ) + ( x x 0 ) ρ Γ ( ρ + 1 ) ( D ρ y ( x 0 ) ) + ( x x 0 ) 2 ρ Γ ( 2 ρ + 1 ) ( D 2 ρ y ( x 0 ) ) .
When ( D ρ y ( x 0 ) ) = H ( x 0 , y ( x 0 ) ) and x x 0 = h are substituted then for y ( x 1 ) , we have
y ( x 1 ) = y ( x 0 ) + h ρ Γ ( ρ + 1 ) ( D ρ y ( x 0 ) ) + h 2 ρ Γ ( 2 ρ + 1 ) ( D 2 ρ y ( x 0 ) ) .
If the step size h is very small then we neglect the second-order term involving h 2 ρ and we obtain
y ( x 1 ) = y ( x 0 ) + h ρ Γ ( ρ + 1 ) ( D ρ y ( x 0 ) ) .
The process is repeated and generates a sequence of points that approximates the solution y ( x ) . The general formula for the fractional Euler’s method [68] over the interval [ 0 , a ] is
y ( x n + 1 ) = y ( x n ) + h ρ Γ ( ρ + 1 ) H ( x n , y ( x n ) ) ,
where x n + 1 = x n + h , n = 0 , 1 , , n 1 .
For the value ρ = 1 , it reduces to the classical Euler’s method.
Lagrange’s two-step Method
Integrating (14), we obtain
y ( x ) y ( 0 ) = 1 Γ ( ρ ) 0 x ( x t ) ρ 1 H ( t , y ( t ) ) d t .
At the point x n + 1 = ( n + 1 ) h , n = 0 , 1 , 2 ,
y ( x n + 1 ) y ( 0 ) = 1 Γ ( ρ ) 0 x n + 1 ( x n + 1 t ) ρ 1 H ( t , y ( t ) ) d t
and we have
y ( x n + 1 ) = y ( 0 ) + 1 Γ ( ρ ) m = 0 n x m x m + 1 ( x n + 1 t ) ρ 1 H ( t , y ( t ) ) d t .
To simplify the integral on the right-hand side of Equation (21), we replace the function H ( t , y ( t ) ) with two-point Lagrange polynomial, and we obtain
y ( x n + 1 ) = y ( 0 ) + 1 Γ ( ρ ) m = 1 n x m x m + 1 ( x n + 1 t ) ρ 1 1 h H ( x m , y m ) ( t x m 1 ) H ( x m , y m 1 ) ( t x m ) d t y ( x n + 1 ) = y ( 0 ) + 1 Γ ( ρ ) m = 1 n H ( x m , y m ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m 1 ) d t 1 Γ ( ρ ) m = 1 n H ( x m 1 , y m 1 ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m ) d t
Now, take
x m x m + 1 ( t x m 1 ) ( x n + 1 t ) ρ 1 d t = ( t x m 1 ) ( x n + 1 t ) ρ ρ ( x n + 1 t ) ρ + 1 ρ ( ρ + 1 ) | x m x m + 1 = ( x n + 1 x m + 1 ) ρ ρ ( x m 1 x m + 1 ) ( x n + 1 x m + 1 ) ( ρ + 1 ) + ( x n + 1 x m ) ρ ρ ( x m x m 1 ) + ( x n + 1 x m ) ( ρ + 1 ) = ( n m ) ρ h ρ + 1 ρ ( ρ + 1 ) 2 ρ 2 n + m + ( n m + 1 ) ρ h ρ + 1 ρ ( ρ + 1 ) ρ 2 + n m = h ρ + 1 ρ ( ρ + 1 ) ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) .
Similarly,
x m x m + 1 ( t x m ) ( x n + 1 t ) ρ 1 d t = h ρ + 1 ρ ( ρ + 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) .
After substituting the above values in Equation (22), we obtain the Lagrangian approximation [69] of y, i.e., the solution of the Equation (14),
y n + 1 = y 0 + h ρ Γ ( ρ + 2 ) m = 1 n H ( x m , y m ) ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ Γ ( ρ + 2 ) m = 1 n H ( x m 1 , y m 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ )
The error for the Lagrange’s two-step method is obtained as
| E n ρ | h ρ + 2 2 Γ ( ρ + 2 ) sup t [ 0 , x n + 1 ] | 2 t 2 ( H ( t , y ( t ) ) | × [ ( n + 1 ) ρ ρ h ρ ] × n ( n + 4 + 2 ρ ) 2 .
Euler–Maruyama scheme for Ito stochastic differential equations
Let us consider the Ito Stochastic differential equation
d y ( x ) = f ( x , y ( x ) ) d x + σ g ( x , y ( x ) ) d B ( x ) , a x b , y ( 0 ) = y 0 ,
with y ( x ) R for each x, σ is the intensity of the white noise, d B is the white noise, f : R R and g : R × [ a , b ] R .
By using the definition of stochastic differential, the Equation (27) is equivalent to the following stochastic integral equation
y ( x ) = y 0 + x 0 x f ( t , y ( t ) ) d t + x 0 x σ g ( t , y ( t ) ) d B ( t ) , a x b .
For a discretization a = x 0 < x 1 < x 2 < < x n = b of the interval [ a , b ] , a Euler approximation is a continuous time stochastic process y = { y ( x ) , a x b } satisfying the iterative scheme
y n + 1 = y n + f ( x n , y n ) ( x n + 1 x n ) + g ( x n , y n ) ( B x n + 1 B x n ) ,
where y n = y ( x n ) for the values of the approximation at the discretization x n = x 0 + n h with step size h = b a N = Δ x for some integer N (page 110 of [70]). The sequence { y n , n = 1 , , N } of values of Euler’s approximation (29) at the points x n , n = 0 , 1 , , N of the time discretization can be computed similar to those of the deterministic case. The main difference is to generate random increments
Δ B n = B x n + 1 B x n ,
for n = 0 , 1 , , N 1 of the Brownian process B = { B x , x 0 } . The increments are independent Gaussian random variables with mean E ( Δ B n ) = 0 and variance E ( ( Δ B n ) 2 ) = Δ x . The increments (30) of the Brownian process can be generated by some random number generator [71,72]. So that
Δ B n η Δ x ,
η N o r m a l ( 0 , 1 ) is the standard normal random number.
Lagrange’s scheme for Caputo fractional stochastic differential equation
We here develop the two-point Lagrange’s scheme for solving the stochastic Caputo fractional differential equations based on the Euler–Maruyama scheme. Consider the Ito stochastic Caputo fractional differential equation
  0 C D x ρ y ( x ) = f ( x , y ( x ) ) d x + σ g ( x , y ( x ) ) d B ( x ) , y ( x 0 ) = y 0 .
Integrating (32), we obtain
y ( x ) y ( 0 ) = 1 Γ ( ρ ) 0 x ( x t ) ρ 1 f ( t , y ( t ) ) d t + σ Γ ( ρ ) 0 x ( x t ) ρ 1 g ( t , y ( t ) ) d B ( t ) .
At the point x n + 1 = ( n + 1 ) h , n = 0 , 1 , 2 ,
y ( x n + 1 ) y ( 0 ) = 1 Γ ( ρ ) 0 x n + 1 ( x n + 1 t ) ρ 1 f ( t , y ( t ) ) d t + σ Γ ( ρ ) 0 x n + 1 ( x n + 1 t ) ρ 1 g ( t , y ( t ) ) d B ( t ) .
and we have
y ( x n + 1 ) y ( 0 ) = 1 Γ ( ρ ) m = 0 n x m x m + 1 ( x n + 1 t ) ρ 1 f ( t , y ( t ) ) d t + σ x m x n + 1 ( x n + 1 t ) ρ 1 g ( t , y ( t ) ) d B ( t ) = 1 Γ ( ρ ) m = 0 n x m x m + 1 ( x n + 1 t ) ρ 1 f ( t , y ( t ) ) d t + σ x m x n + 1 ( x n + 1 t ) ρ 1 g ( t , y ( t ) ) B ˙ ( t ) d t .
To simplify the integral on the right-hand side of Equation (35), we replace the functions f ( t , y ( t ) ) and g ( t , y ( t ) ) with the two-point Lagrange polynomial, and we obtain
y ( x n + 1 ) = y ( 0 ) + 1 Γ ( ρ ) m = 1 n x m x m + 1 ( x n + 1 t ) ρ 1 1 h f ( x m , y m ) ( t x m 1 ) f ( x m , y m 1 ) ( t x m ) d t + σ Γ ( ρ ) m = 1 n x m x m + 1 ( x n + 1 t ) ρ 1 1 h Φ ( x m , y m ) ( t x m 1 ) Φ ( x m , y m 1 ) ( t x m ) d t y ( x n + 1 ) = y ( 0 ) + 1 Γ ( ρ ) m = 1 n f ( x m , y m ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m 1 ) d t 1 Γ ( ρ ) m = 1 n f ( x m 1 , y m 1 ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m ) d t + 1 Γ ( ρ ) m = 1 n Φ ( x m , y m ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m 1 ) d t 1 Γ ( ρ ) m = 1 n Φ ( x m 1 , y m 1 ) h x m x m + 1 ( x n + 1 t ) ρ 1 ( t x m ) d t
where Φ ( t , y ( t ) ) = g ( t , y ( t ) ) B ˙ ( t ) . On solving, we obtain
y n + 1 = y 0 + h ρ Γ ( ρ + 2 ) m = 1 n f ( x m , y m ) ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ Γ ( ρ + 2 ) m = 1 n f ( x m 1 , y m 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) + h ρ Γ ( ρ + 2 ) m = 1 n Φ ( x m , y m ) ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ Γ ( ρ + 2 ) m = 1 n Φ ( x m 1 , y m 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ )
The Lagrangian approximation [73,74] of y, representing the solution to Equation (27), is provided by
y n + 1 = y 0 + h ρ Γ ( ρ + 2 ) m = 1 n f ( x m , y m ) ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ Γ ( ρ + 2 ) m = 1 n f ( x m 1 , y m 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) + h ρ 1 2 Γ ( ρ + 2 ) m = 1 n g ( x m , y m ) ( B ( x n ) B ( x n 1 ) ) × ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ 1 2 Γ ( ρ + 2 ) m = 1 n g ( x m 1 , y m 1 ) ( B ( x n 1 ) B ( x n 2 ) ) × ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) .

3. Mathematical Models for Epidemics

Deterministic models are a very popular modeling technique. They are often used in the mathematical modeling of infectious diseases. When determining disease control strategies in animal health, there remain gaps in understanding how the interaction between the environment and the host affects infection transmission and disease progression. Epidemiological disease models provide a way to address these uncertainties by integrating data from field and experimental research with expert insights into infection dynamics and disease management strategies [75]. Epidemiological models divide the complete population into three classes which change with time t: susceptible, infective, and recovered individuals and the spread of an epidemic is governed by the principle of mass action [7], for a certain disease. In the SIS model recovery does not lead to immunity. In this model, individuals move from the susceptible class to the infective class and back again. But in the SIR model individuals recover with permanent immunity. If the immunity is temporary in the recovered class, it is an SIRS model. If there is no recovery, an SI model is used. In SI models, people never leave the infectious state and have lifelong infections. For example, herpes is a disease with lifelong infectiousness. SIR models are effective for viral agents such as Influenza and SIS models are useful for bacterial agents such as plague and meningitis. Hethcote [76] discussed the qualitative analysis of communicable disease models. He determined the asymptotic stability regions for the equilibrium points for models involving temporary immunity, disease-related fatalities, carriers, migration, dissimilar interacting groups, and transmission by vectors.
Models are used to predict the progress of epidemics and to examine the impact of various model parameters. By analyzing these parameters, researchers can understand how changes might impact various control measures, such as vaccination, treatment, and other preventive strategies [77]. Mathematical models are essential for analyzing, explaining, and forecasting the dynamics of disease transmission. For globally significant diseases, these models also help in controlling outbreaks and formulating public health strategies as preventive measures [78]. These models are useful in comparing the effects of prevention or control procedures.
Definitions and Notation
The various variables, their notations and parameters used in this manuscript are described here.
  • Susceptible s : refers to the group of individuals who have not been infected yet but are susceptible to becoming infected. These individuals may remain susceptible.
    The fraction of susceptible individuals, denoted by S, is defined as the ratio of the number of susceptible individuals s to the total population N, i.e., S = s N .
  • Infected i : represents the group of people who are infected with the virus and have the ability to spread it to other people who are susceptible. An infected person may remain infectious for a very long time and may be removed from the infected population to improve or disappear completely.
    The fraction of infected individuals, denoted by I, is defined as the ratio of the number of infected individuals i to the total population N, i.e., I = i N .
  • Removed r : is the class of the individuals who have recovered from the virus and are assumed to be immune, r ( t ) or have died, d ( t ) .
    The fraction of recovered individuals, denoted by R, is defined as the ratio of the number of susceptible individuals r to the total population N, i.e., R = r N .
  • Exposed E : is the fraction of the individuals whose body is a host for infection but are not yet able to transmit the disease.
  • Vaccinate V : is the fraction of the individuals who are vaccinated and are removed from the susceptible compartment.
  • Quarantine Q : represents those who have been placed in isolation in order to avoid the disease from spreading to other people.
  • Maternally Derived Immunity M : is the fraction of individuals with passive immunity, protected by maternal antibodies.
All the variables are the function of time t. The parametric values used in the graphical simulation are given in the Table 1.

3.1. Qualitative Analysis of the Model

The objective of this section is to perform a qualitative analysis of the dynamical systems discussed in the manuscript. We begin by identifying the equilibrium points of the dynamical system and subsequently calculate the basic reproduction number at the disease-free equilibrium point. Following this, we present the theorems that define the stability of the system. Additionally, we discussed the theorem required for establishing the existence and uniqueness of the solution.

3.1.1. Equilibrium Points and Stability of the System

For every dynamical system, the study of equilibrium points and the stability of the system is very important. An equilibrium (or equilibrium point) of a dynamical system is a solution that does not change with time. In mathematical terms, for the system differential equations
d y 1 d t = f 1 ( y 1 , y 2 , , y n ) , d y 2 d t = f 2 ( y 1 , y 2 , , y n ) , d y n d t = f n ( y 1 , y 2 , , y n ) ,
the equilibrium points are evaluated by putting d y k d t = 0 , k = 1 , 2 , . . . , n .
For the epidemiological models, there are two types of equilibrium points: disease-free and endemic.
Disease-Free Equilibrium Point When there are no infected people in the population, the equilibrium point is disease-free. At this point, the number of infected individuals I ( t ) is equal to zero, indicating that the entire population remains susceptible to the disease.
Endemic Equilibrium Point When the number of infected individuals in the model is not equal to zero, i.e., I 0 , the model attains its equilibrium in the state of endemic.
The reproduction number found using the next-generation matrix is evaluated at the disease-free equilibrium point.

3.1.2. Basic Reproduction Number ( R 0 )

For the study of an epidemic, the study of its reproduction number is highly significant. The basic reproduction number ( R 0 ) [79] is the expected number of secondary cases produced in a completely susceptible population, by a typical infected individual during its entire period of infectiousness. In general, it is the measure of how rapidly an infectious disease spreads through a population. If the transmission rate increases, the number of infected persons increases and hence the value of R 0 increases. If the recovery rate increases, the number of infected persons decreases and hence the value of R 0 decreases. R 0 is dimensionless because it is a ratio of two rates.
There are three cases for the value of R 0 , which are
  • when R 0 < 1 means the number of secondary infections will decrease over time and eventually the outbreak will end on its own,
  • when R 0 = 1 means the number of infectives is stable,
  • when R 0 > 1 means the outbreak is self-sustaining unless effective control measures are implemented.
Thus, when R 0 > 1 the disease can invade but when R 0 < 1 then the diseases die out.
MacDonald [80] was the first person who identified it and gave the name, the threshold quantity explicitly as Basic Reproduction Number (BRN) in his work on malaria.
Theorem 1 
([81,82]). The disease-free equilibrium point of the system of the differential equation is locally asymptotically stable if and only if R 0 < 1 .
Theorem 2. 
The endemic equilibrium point of the system of the differential equation is locally asymptotically stable if and only if R 0 > 1 .
Computation of R 0 by the next generation matrix method
The next-generation matrix method is used to determine the value of R 0 for an epidemic model. This was very first studied by Diekmann in 1990 [79].
Let X = ( x 1 , x 2 , , x n ) T represent the number of individuals in each compartment, where the first m < n compartments contain the infected individuals. Assuming that the disease-free equilibrium point x 0 exists and is stable. The differential equations x 1 , x 2 , , x m at the disease-free equilibrium points decouple from the other equations. A more detailed description of the assumptions may be found in [83,84]. Take these equations and write in the form d x i d t = F i ( x ) V i ( x ) for i = 1 , 2 , , m . Here F i ( x ) is the rate of appearance of new infections in compartment i, and V i ( x ) denotes the rate of other transitions between compartment i and other infected compartments. It is assumed that each function is continuously differentiable at least twice in each variable [85].
Now, we define F = F i ( x 0 ) x j m × m and V = V i ( x 0 ) x j m × m for 1 i , j m such that F is entry-wise non-negative and V is a non-singular matrix, hence V 1 is entry-wise non-negative. Let I ( 0 ) be the number of initially infected individuals, then F V 1 I ( 0 ) is an entry-wise non-negative vector giving the expected number of new infections. Matrix F V 1 has ( i , j ) entry equal to the expected number of the secondary infections in compartment i produced by an infected individual introduced in compartment j.
Thus, F V 1 is the next generation matrix and the reproduction number is given as R 0 = ρ ( F V 1 ) , where ρ is the spectral radius.

4. Kermack–McKendrick SIR Model

The SIR model is a traditional model of disease transmission within a population. The SIR models consist of the variables S, I, and R which are functions of time. The SIR model tracks the number of susceptible, infected and recovered individuals during an epidemic with the help of ordinary differential equations. It provides a conceptual framework for identifying disease transmission throughout a population. It gives us a clear situation and prediction of the spread of the virus in communities. Its purpose is to predict how many people are at risk for infection, have an active infection, or have recovered from infection at any particular time. It is the simplest compartmental model.
The SIR model assumes the homogeneous mixing of the population, meaning that all individuals in the population are assumed to have an equal probability of coming in contact with one another. It assumes a population that is closed, without migration, births, or deaths from causes other than the epidemic itself. The system is dynamic, described by three ordinary differential equations (ODEs) that depict the population’s temporal changes. This epidemiological model captures the dynamics of acute infections that confer lifelong immunity once recovered. This model can be used to treat diseases including measles, smallpox, chickenpox, mumps, typhoid fever, and diphtheria, for individuals who develop lifelong immunity. Generally, the total population size is considered to be constant, i.e., N = s + i + r .
Kermack–McKendrick proposed the following SIR model without vital dynamics. The system of fractional order differential equations (FODEs) describes the rates of change for three populations,
  0 C D t ρ s ( t ) = β s ( t ) i ( t ) N ,   0 C D t ρ i ( t ) = β s ( t ) i ( t ) N γ i ( t ) ,   0 C D t ρ r ( t ) = γ i ( t ) ,
given to the initial conditions s ( 0 ) > 0 , i ( 0 ) > 0 and r ( 0 ) 0 , β is the transmission rate which affects the transition from the susceptible compartment to the infected compartment and γ is the recovery rate which affects the transition from the infected compartment to the recovered compartment.
For the sake of simplicity, consider the proportions by redefining S = s N , I = i N , R = r N . The system (40) reduced to
  0 C D t ρ S ( t ) = β S ( t ) I ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) γ I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) ,
with the initial conditions S ( 0 ) > 0 , I ( 0 ) > 0 and R ( 0 ) 0 and S + I + R = 1 .
Hethcote [76] in 1976 defined the SIR model with vital dynamics. To incorporate demography into the SIR model, the most straightforward approach is to assume a natural lifespan for hosts, which lasts 1 d years. Then, the rate of natural mortality is given by d. It is important to observe that this factor is independent of the disease and is not intended to reflect the pathogenicity of the infectious agent. Here, b represents the birth rate. A schematic representation of the SIR model is given in Figure 1. The normalized initial value problem (IVP) for the SIR model with vital dynamics is given by,
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
where S ( 0 ) > 0 , I ( 0 ) > 0 and R ( 0 ) 0 .
In the SIR model, it is assumed that the individuals recover with permanent immunity.
Remark 1. 
In our scope of study, we have considered all the models with vital dynamics. Since the case, b = d = 0 without vital dynamics is always a particular case.
To find the equilibrium point, equating to zero all the derivatives of the system of Equation (42), we obtain
b β S I d S = 0 ,
β S I γ I d I = 0 ,
γ I d R = 0 .
From the Equation (44),
( β S γ d ) I = 0 ,
we obtain I = 0 or S = γ + d β . On substituting I = 0 in (43) and (45), we get S = b d and R = 0 . Here, the number of infected individuals are zero and whole population is susceptible. Hence, this is the disease free equilibrium point ( S 0 , I 0 , R 0 ) = ( b d , 0 , 0 ) .
Now, we will find the another equilibrium point known as endemic equilibrium point. Now, put S = γ + d β in (43), we obtain
b β γ + d β I d γ + d β = 0 .
On solving, we obtain I = b β d ( d + γ ) β ( d + γ ) . Put the value of I in (45), we obtain R = γ d b β d ( d + γ ) β ( d + γ ) .
Hence, the endemic equilibrium point for the SIR model is
( S * , I * , R * ) = γ + d β , b β d ( d + γ ) β ( d + γ ) , γ d b β d ( d + γ ) β ( d + γ ) .
The system of differential Equation (42) has two equilibrium points.
Reproduction number
We determine the reproduction number at the disease free equilibrium point through the next generation matrix approach. Here, the disease class is I only. Therefore, F = [ β S I ] and V = [ ( γ + d ) I ] . Evaluating the matrices F = β b d 1 × 1 and V = [ γ + d ] 1 × 1 at the disease free equilibrium point, we obtain
F V 1 = β b d ( γ + d ) .
The spectral radius of the matrix F V 1 is the basic reproduction number, R 0 = β b d ( γ + d ) .
Remark 2. 
In case the model is having equal birth and death rates b = d , then R 0 = β ( γ + d ) . Further, if the model is without vital dynamics b = d = 0 then R 0 = β γ .

4.1. Existence and Uniqueness of the Solution of the SIR Model

The existence and uniqueness of the solutions of a system of differential equations are very commonly available in the literature. Here, for the convenience of reader, we have included the existence and uniqueness of the Caputo fractional SIR model. In case of extensions and variations in the SIR model, the same approach can be used to prove the corresponding results.
Consider Ξ = { ( S , I , R ) R 3 + ; S + I + R b d } , we will show that the closed set Ξ is the feasible region of the system (42).
Lemma 1 
([86]). The closed set Ξ is positively invariant with respect to the fractional system (42).
Proof. 
Add all the relations given in the system (42), we obtain
  0 C D t ρ N ( t ) = b d ( S ( t ) + I ( t ) + R ( t ) ) = b d N ( t ) .
Taking the Laplace transform using the Equation (12) for 0 ρ 1 of (46), we obtain
s ρ L { N ( t ) } s ρ 1 N ( 0 ) = b d L { N ( t ) } ( s ρ + d ) L { N ( t ) } = b s + s ρ 1 N ( 0 ) L { N ( t ) } = b s ( s ρ + d ) + s ρ 1 N ( 0 ) s ρ + d N ( t ) = L 1 b s ( s ρ + d ) + s ρ 1 N ( 0 ) s ρ + d = L 1 b d s b s ρ 1 d ( s ρ + d ) + N ( 0 ) L 1 s ρ 1 s ρ + d
After solving the inverse Laplace, we get
N ( t ) = b d b d E ρ , 1 ( d t ρ ) + N ( 0 ) E ρ , 1 ( d t ρ ) ,
where E ρ , π ( x ) = n = 0 x n Γ ( ρ n + π ) , x R , ρ , π > 0 is the two parameter Mittag–Leffler function and E ρ , 1 ( x ) = E ρ ( x ) = n = 0 x n Γ ( ρ n + 1 ) .
Hence,
N ( t ) = b d + E ρ ( d t ρ ) N ( 0 ) b d .
Thus, if N ( 0 ) b d , then for t > 0 , N ( t ) b d . Hence, the closed set Ξ is positively invariant with respect to fractional model (42). □
Existence and Uniqueness of Solution of the SIR Model
C ( J ) denotes the set of real-valued functions that are continuous and defined on the interval J, which is a subset of the real numbers, R . Then, V = C ( J ) C ( J ) C ( J ) is the Banach space with the norm for ( S , I , R ) V defined as ( S , I , R ) = S + I + R , where S , I , R C ( J ) and · = sup t J | · | . We will use the fixed point theorem to demonstrate that the system of differential equations describing the SIR model has a solution.
Applying the integral operator (1) to the system of Equation (42),
  0 R L I t ρ   0 C D t ρ S ( t ) =   0 R L I t ρ ( b β S ( t ) I ( t ) d S ( t ) ) ,   0 R L I t ρ   0 C D t ρ I ( t ) =   0 R L I t ρ ( β S ( t ) I ( t ) γ I ( t ) d I ( t ) ) ,   0 R L I t ρ   0 C D t ρ R ( t ) =   0 R L I t ρ ( γ I ( t ) d R ( t ) ) ,
we obtain
S ( t ) S ( 0 ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ( b β S ( u ) I ( u ) d R ( u ) ) d u ,
I ( t ) I ( 0 ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ( β S ( u ) I ( u ) γ I ( u ) d I ( u ) ) d u ,
R ( t ) R ( 0 ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ( γ I ( u ) d S ( u ) ) d u .
Let us denote the following as kernels:
K 1 ( t , S ) = b β S I d S , K 2 ( t , I ) = β S I γ I d I , K 3 ( t , R ) = γ I d R .
The following theorem states that the kernels K i , i = 1 , 2 , 3 , satisfy certain conditions.
Theorem 3. 
Let S ( t ) k 1 , I ( t ) k 2 , R ( t ) k 3 . The Lipschitz condition and contraction would be satisfied by K 1 , K 2 , K 3 , define in (54) for the Lipschitz constants 0 l 1 = | β k 2 + d | < 1 , 0 l 2 = | β k 1 γ d | < 1 , 0 l 3 = | d | < 1 .
Proof. 
To begin, let us use K 1 . Given two functions S and S ( 1 ) that satisfy the Equation (42), then
K 1 ( t , S ) K 1 ( t , S ( 1 ) ) = b β S I d S b + β S ( 1 ) I + d S ( 1 ) = ( β I d ) ( S S ( 1 ) ) = | β I + d | ( S S ( 1 ) ) | β k 2 + d | S S ( 1 ) .
K 1 ( t , S ) K 1 ( t , S ( 1 ) ) l 1 S S ( 1 ) .
Clearly, l 1 is a fixed parameter and S is a bounded function. Hence, the Lipschitz condition is satisfied for K 1 , and it is contraction mapping when 0 l 1 = β k 2 + d < 1 .
Similarly, the other two kernels also satisfy the Lipschitz condition, i.e.,
K 2 ( t , I ) K 2 ( t , I ( 1 ) ) l 2 I I ( 1 ) , K 3 ( t , R ) K 3 ( t , R 1 ) l 3 R R ( 1 ) .
Using the aforementioned kernels from (54) in the Equations (51)–(53), we obtain
S ( t ) = S ( 0 ) + 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 1 ( u , S ( u ) ) d u , I ( t ) = I ( 0 ) + 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 2 ( u , I ( u ) ) d u , R ( t ) = R ( 0 ) + 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 3 ( u , R ( u ) ) d u .
The corresponding recursive formulas are provided by
S ( n ) ( t ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 1 ( u , S ( n 1 ) ( u ) ) d u , I ( n ) ( t ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 2 ( u , I ( n 1 ) ( u ) ) d u , R ( n ) ( t ) = 1 Γ ( ρ ) 0 t ( t u ) ρ 1 K 3 ( u , R ( n 1 ) ( u ) ) d u .
The initial conditions are S ( 0 ) = S ( 0 ) , I ( 0 ) = I ( 0 ) , R ( 0 ) = R ( 0 ) .
The difference of each term in (59) with the succeeding term is represented by the following equations, respectively,
ψ 1 n ( t ) = S ( n ) ( t ) S ( n 1 ) ( t ) = 1 Γ ( ρ ) 0 t K 1 ( u , S ( n 1 ) ( u ) ) K 1 ( u , S ( n 2 ) ( u ) ) ( t u ) ρ 1 d u ,
ψ 2 n ( t ) = I ( n ) ( t ) I ( n 1 ) ( t ) = 1 Γ ( ρ ) 0 t K 2 ( u , I ( n 1 ) ( u ) ) K 2 ( u , I ( n 2 ) ( u ) ) ( t u ) ρ 1 d u ,
ψ 3 n ( t ) = R ( n ) ( t ) R ( n 1 ) ( t ) = 1 Γ ( ρ ) 0 t K 3 ( u , R ( n 1 ) ( u ) ) K 3 ( u , R ( n 2 ) ( u ) ) ( t u ) ρ 1 d u .
Taking the norm of the above system
ψ 1 n ( t ) = S ( n ) ( t ) S ( n 1 ) ( t ) = 1 Γ ( ρ ) 0 t K 1 ( u , S ( n 1 ) ( u ) ) K 1 ( u , S ( n 2 ) ( u ) ) ( t u ) ρ 1 d u 1 Γ ( ρ ) 0 t K 1 ( u , S ( n 1 ) ( u ) ) K 1 ( u , S ( n 2 ) ( u ) ) ( t u ) ρ 1 d u .
As the kernel K 1 fulfil the Lipschitz condition, we have
ψ 1 n ( t ) l 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ψ 1 ( n 1 ) ( u ) d u .
Similarly,
ψ 2 n ( t ) l 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ψ 2 ( n 1 ) ( u ) d u ,
ψ 3 n ( t ) l 1 Γ ( ρ ) 0 t ( t u ) ρ 1 ψ 3 ( n 1 ) ( u ) d u .
Hence,
S ( n ) ( t ) = i = 1 n ψ 1 i ( t ) ,
I ( n ) ( t ) = i = 1 n ψ 2 i ( t ) ,
R ( n ) ( t ) = i = 1 n ψ 3 i ( t ) .
Theorem 4. 
The system of fractional SIR model (42) has an exact coupled solution under the condition that we can find t 1 such that
l 1 t 1 ρ Γ ( ρ + 1 ) 1 , l 2 t 1 ρ Γ ( ρ + 1 ) 1 , l 3 t 1 ρ Γ ( ρ + 1 ) 1 .
Proof. 
The functions S ( t ) , I ( t ) , R ( t ) are bounded and the kernels K i , i = 1 , 2 , 3 satisfied the Lipschitz condition,
S ( t ) S ( 0 ) = S ( n ) ( t ) H 1 ( n ) ( t ) , I ( t ) I ( 0 ) = I ( n ) ( t ) H 2 ( n ) ( t ) , R ( t ) R ( 0 ) = R ( n ) ( t ) H 3 ( n ) ( t ) .
Now,
H 1 ( n ) ( t ) = S ( n ) ( t ) S ( t ) + S ( 0 ) H 1 ( n ) ( t ) = 1 Γ ( ρ ) 0 t K 1 ( u , S ( n 1 ) ( u ) ) ( t u ) ρ 1 d u S ( t ) + S ( 0 ) = 1 Γ ( ρ ) 0 t K 1 ( u , S ( n 1 ) ( u ) ) ( t u ) ρ 1 d u S ( 0 ) 0 t K 1 ( u , S ( u ) ) ( t u ) ρ 1 d u + S ( 0 ) = 1 Γ ( ρ ) 0 t ( K 1 ( u , S ( n 1 ) ( u ) ) K 1 ( u , S ( u ) ) ) ( t u ) ρ 1 d u l 1 Γ ( ρ ) ( S ( n 1 ) S ) 0 t ( t u ) ρ 1 d u l 1 t ρ ρ Γ ( ρ ) S ( n 1 ) S .
Repeating the above process, we obtain
H 1 ( n ) ( t ) l 1 t ρ Γ ( ρ + 1 ) n + 1 λ ,
where λ = S ( n 1 ) S and ∃ t 1 such that
H 1 ( n ) ( t ) l 1 t 1 ρ Γ ( ρ + 1 ) n + 1 λ .
Taking the limit n , since 0 l 1 t 1 ρ Γ ( ρ + 1 ) < 1 ,
H 1 ( n ) ( t ) 0 . S ( t ) S ( 0 ) = lim n S ( n ) ( t ) .
Similarly,
H 2 ( n ) ( t ) l 2 t 1 ρ Γ ( ρ + 1 ) I ( n 1 ) I , H 3 ( n ) ( t ) l 3 t 1 ρ Γ ( ρ + 1 ) R ( n 1 ) R ,
and hence, we have
H 2 ( n ) ( t ) 0 . I ( t ) I ( 0 ) = lim n I ( n ) ( t ) , H 3 ( n ) ( t ) 0 . R ( t ) R ( 0 ) = lim n R ( n ) ( t ) .
Hence, the existence of the solution is proved.
Theorem 5. 
The solution of the system of fractional SIR model (42) is unique under the condition (70).
Proof. 
For the sake of convenience t 1 = t . To demonstrate the uniqueness of the solution, let us suppose that S , I , R be the other set of solutions of the system (42). Then, from (59), we have
S ( t ) S ( t ) = 1 Γ ( ρ ) 0 t ( K 1 ( u , S ) K 1 ( u , S ) ) ( t u ) ρ 1 d u S ( t ) S ( t ) l 1 t ρ Γ ( ρ + 1 ) S ( t ) S ( t ) S ( t ) S ( t ) = 0 , sin ce l 1 t ρ Γ ( ρ + 1 ) < 1 S ( t ) = S ( t ) .
Similarly, we can show for the variable I and R. Hence, the uniqueness of the solution is proved. □
In [86], the authors have proved the existence and uniqueness of the solution for the SEIR epidemic model for the spread of COVID-19 using the Caputo fractional derivative. Using the fixed-point theorem, the authors in [41] proved that there exists a unique solution for the fractional Ebola virus model. [60], the researchers established the existence and uniqueness of the solution for the non-linear tuberculosis epidemic model using the Caputo operator. The above theorem can be proved for other epidemic models in a similar way.
Graphical simulation
The graphical solution of the SIR model obtained using the Lagrange’s two point method is illustrated in Figure 2. Figure 2a represents the SIR model without vital dynamics. Here, we draw the graphs for the integer order and fractional order ( ρ = 0.9 ) models to show the comparison between them. The proportion of susceptibility decreases as time increases and becomes stable after 20 days. Initially, the proportion of infected is very small which increases gradually to attain a peak (within 10–12 days depending on the value of R 0 , here R 0 = 3.182), and then starts decreasing due to the recovering population.
When the model is considered with vital dynamics, due to natural births and deaths rates, the graph changes its shape as shown in Figure 2b. It is observed that the changes in factional order ( ρ = 0.9 ) model are comparatively at a slower pace with reference to integer order. The order of the fractional model can be chosen to fit the real data.
Tian et al. [16] discussed the SIR epidemic model in 2020 and discussed its global stability. In [87], Koziol et al. discussed the fractional-order generalization of the susceptible-infected-recovered (SIR) epidemic model for predicting the spread of the COVID-19 disease in 2020. Alshomrani et al. [88] gave a SIR model on the disease COVID-19 in the year 2021. In this paper, they analyzed the qualitative behavior of the fractional SIR model involving the Caputo derivative, uniqueness conditions under the Banach contraction principle and stability analysis with basic reproduction number.
The two models SIS and IR, which are very different from the SIR model but whose study is vital since many diseases fit into these two categories are covered in the next two subsections.

4.2. The SIS Model

If the recovery from a disease does not provide immunity, then the model is converted to SIS model, since people move from susceptible class to infected then the return back to susceptible class. A schematic representation of the SIS model is given in Figure 3.
In general, SIR models are appropriate for viral agent diseases such as measles, mumps, and smallpox, while SIS models are appropriate for some bacterial agent diseases such as meningitis, plague, and venereal diseases, and for protozoan agent diseases such as malaria and sleeping sickness [76,89].
The dynamics of the normalized SIS model is given as follows:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) + ω I ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) ω I ( t ) d I ( t )
where S ( 0 ) > 0 , I ( 0 ) > 0 , S ( t ) + I ( t ) = 1 .
The disease free equilibrium point for the SIS model is ( S 0 , I 0 ) = ( b d , 0 ) .
The disease class is I here. Hence, F = [ β S I ] and V = [ ( ω + d ) I ] . The matrix F = β b d 1 × 1 and V = [ ω + d ] 1 × 1 at the disease free equilibrium point.
F V 1 = β b d ( ω + d ) . The eigen-value of the matrix F V 1 is β b d ( ω + d ) only.
Hence, the basic reproduction number, R 0 = β b d ( ω + d ) .
Graphical simulation
In Figure 4 the graphical solution for the SIS model is given. The figure shows that when the fraction of susceptible decreases, the fraction of infective increases. Initially, there is a fast decline in the fraction of the susceptible population, which stabilizes after approximately 10 days. Initially infectivity starts increasing with lag and then increases sharply till a major portion of the population is infected. It keeps on increasing gradually due to difference in the death and birth rate. The inverse relationship between the two curves highlights the transition dynamics from susceptible to infected as the disease spreads through the population.
Liu et al. in 2019 [90] discussed the fractional order SIS model with Caputo derivative and investigated the stability of the fractional SIS model on complex networks with linear treatment function. Angran et al. [91] in 2023 described the analytical solution of a non-linear fractional order SIS epidemic model using the Laplace residual power series (LRPS) method.

4.3. The IR Model

The infected-recovered mathematical model in which individuals are being infected and then recover The population dynamics of the IR model of epidemics is given in Figure 5. The dynamics of two populations are determined by the following system of fractional order differential equations (FODEs):
  0 C D t ρ I ( t ) = γ I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) ,
where γ is recovery rate and with initial conditions I ( 0 ) = I 0 > 0 and R ( 0 ) = R 0 > 0 . The larger the value of γ the more quickly people go from I to R compartment.
For ρ = 1 , the differential Equation (80) can be easily solved analytically, hence we get an analytical solution of the IR model.
Consider,
d I d t = γ I d I I = γ d t log ( I ) = γ t + log ( I 0 ) log I I 0 = γ t I ( t ) = I 0 e γ t .
Now take,
d R d t = γ I d R d t = γ I 0 e γ t R ( t ) = γ I 0 e γ t γ + κ R ( t ) = I 0 e γ t + κ
where κ is the arbitrary constant. Using initial condition R ( 0 ) = R 0 , we get κ = R 0 + I 0 .
The system of Equation (80) can be solved as follows:
I ( t ) = I 0 e γ t R ( t ) = R 0 + I 0 ( 1 e γ t ) .
As t number of infected individuals is zero, i.e., infected individuals are removed from the population and the whole population is recovered.
Graphical simulation
The graphical representation of the solution of the IR model is given in the Figure 6. The figure shows that the infected and recovered populations are inversely proportional to each other. As the proportion of infected decreases, the fraction of recovery increases. It is observed that in the fractional order model, the infection occurs in a lesser population compared to the integer order model and similarly in the recovered.

5. Extensions of the SIR Model

The SIR model is a very simple epidemic model, but we can extend it by including other epidemiological factors such as loss of immunity, Maternal immunity, vaccination, quarantine, etc. In this section, we discuss various extensions of the fractional order SIR models.

5.1. The SIRS Model

The SIRS model is exactly the same as the Kermack–McKendrick model with the addition of δ , the rate at which recovered individuals return back to the susceptible category. In the SIR model recovery implies permanent acquired immunity, but here removed individuals return to the susceptibles compartment with the rate δ . The population dynamics of the SIRS model of epidemics are given in Figure 7. This model uses the following system of differential equations in normalized form:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) + δ R ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) δ R ( t ) d R ( t ) ,
where β is the transmission rate of infection, γ is the recovery rate and δ is the rate of loss of immunity.
The SIRS model can apply to many infectious diseases such as polio, tetanus, diphtheria, measles, hepatitis, and influenza, chickenpox, influenza, measles, mumps, rubella, AIDS and others [92,93].
Disease free equilibrium point
To find the equilibrium point put the system of Equation (82) equal to zero. As we are interested in finding the disease-free equilibrium point, so put I = 0 . Hence, we obtain the disease-free equilibrium point for the SIRS model as ( S 0 , I 0 , R 0 ) = ( b d , 0 , 0 ) .
Reproduction number
We determine the reproduction number at the disease-free equilibrium point using the next-generation matrix method. Here, the reproduction number is the same as the SIR model, i.e., R 0 = β b d ( γ + d ) because the disease class and disease-free equilibrium point in SIR and SIRS are the same.
Graphical simulation
The graphical representation of the solution for the SIRS model can be shown in Figure 8. In the SIRS model, initially the fraction of susceptible is high, then decreases rapidly as infections rise and then fluctuates as the fraction of population moves between compartments due to immunity loss. The fraction of infective rises to a peak, representing the height of the epidemic, and then decreases as the recovery increases. The fraction of the recovered population starts at zero, then increases rapidly, and then decreases and increases depending on the rate of immunity loss relative to recovery. It is observed that the fractional order ( ρ = 0.9 ) model shows a small fluctuation in the curves as compared to integer order as shown in the graph.
Korobeinikov and Wake [94] discussed the SIR, SIRS and SIS epidemiological models in 2002 and discussed its Lyapunov and global stability. Zhen et al. [95] and Nakata et al. [96] discussed the SIRS epidemic model in 2006 and discussed its global stability analysis. El et al. in 2029 [97] introduced a fractional order SIRS model with homogeneous networks. In this paper, he discussed the effect of fractional order equations with the Caputo derivative and the role of networks, especially in local and global stability.

5.2. The SEIR Model

The SIR model is extended by adding a compartment E (exposed). A schematic representation of the SEIR model is given in Figure 9. A significant amount of time passes after someone catches a major infection during which they are infected but not yet infectious. So, during that time, people move into the exposed compartment. The total number of population is 1, i.e.,
S ( t ) + E ( t ) + I ( t ) + R ( t ) = 1 .
The dynamics of the normalized SEIR model is given as follows:
  0 C D t ρ S ( t ) = b α S ( t ) I ( t ) d S ( t ) ,   0 C D t ρ E ( t ) = α S ( t ) I ( t ) ϵ E ( t ) d E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
with the initial conditions, S ( 0 ) > 0 , I ( 0 ) > 0 , E ( 0 ) 0 , R ( 0 ) 0 .
Disease free equilibrium point
To determine the equilibrium point in the SEIR model, set the right side of the system of Equation (83) to zero and for the disease-free equilibrium point put E = 0 and I = 0 . We obtain the disease-free equilibrium point for the SEIR model, i.e., ( S 0 , E 0 , I 0 , R 0 ) = ( b d , 0 , 0 , 0 ) .
Reproduction number
Using the next generation matrix method we find the reproduction number at the disease-free equilibrium point. Here, the disease classes are
  0 C D t ρ E ( t ) = α S ( t ) I ( t ) ϵ E ( t ) d E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) γ I ( t ) d I ( t ) .
Now, F = α S I 0 and V = ϵ E + d E ϵ E + γ I + d I .
At the disease-free equilibrium point, we have F = 0 α b d 0 0 and V = ϵ + d 0 ϵ γ + d .
Therefore, F V 1 = b α ϵ d ( ϵ + d ) ( γ + d ) α ( γ + d ) [ 4 p t ] 0 0 .
The dominant eigen-value of matrix F V 1 is b α ϵ d ( ϵ + d ) ( γ + d ) .
Hence, the basic reproduction number, R 0 = b α ϵ d ( ϵ + d ) ( γ + d ) .
Graphical simulation
Figure 10 illustrates the graphical solution of the SEIR model. The figure shows that when the fraction of the susceptible population becomes infectious, it starts to decrease and the fraction of exposed starts to increase, attaining a peak and then decreasing over time. The infected population simultaneously rises, attains a peak and decreases as recovery begins. From the graph, it is observed that the fractional order model has a slower pace than the integer order.
Li and Muldowney [98] discussed the SEIR Model in 1994 and discussed its stability analysis. In 1999, Michael et al. [99] talked about the global dynamics of an SEIR model with a variable total population size. Korobeinikov [100] discussed the SEIR and SEIS epidemic models in 2004 and discussed the Lyapunov functions. El-Sheikh and El-Marouf [101] discussed the SEIR epidemic model with the vertical transmission in 2004 and discussed its stability and bifurcation of solutions.
Sene et al. [102] discussed the SEIR epidemic models in 2022 using the Caputo derivative and discussed the local stability and global asymptotic stability of the equilibrium points of the SEIR model with the Matignon criterion and Lyapunov direct method.

5.3. The SEIRS Model

In the SEIR model when the immunity is lost then recovered individuals come back to the susceptible class again and the mode we obtained is the SEIRS model. A schematic representation of the SEIRS model is given in Figure 11. Diseases like some viral infections that have a time of incubation and temporary immunity are studied using the SEIRS model.
The dynamics of the normalized SEIRS model [103] is given as
  0 C D t ρ S ( t ) = b α S ( t ) I ( t ) + δ R ( t ) d S ( t ) ,   0 C D t ρ E ( t ) = α S ( t ) I ( t ) ϵ E ( t ) d E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) δ R ( t ) d R ( t ) ,
with the initial conditions S ( 0 ) > 0 , E ( 0 ) > 0 , I ( 0 ) > 0 , R ( 0 ) 0 .
Disease free equilibrium point
The disease free equilibrium point is ( S 0 , E 0 , I 0 , R 0 ) = ( b d , 0 , 0 , 0 ) .
Reproduction number
Since the equations of disease class are the same for both SEIR and SEIRS models, the basic reproduction number of the model SEIRS model is the same as the SEIR model, i.e., R 0 = b α ϵ d ( ϵ + d ) ( γ + d ) .
Graphical simulation
Figure 12 illustrates the graphical solution of the SEIRS model. The SEIRS model, which includes temporary immunity, provides a more realistic representation of certain infectious diseases than the SEIR model. Due to temporary immunity, the recovered become susceptible after a certain period of time, hence the peak of the recovered is lower compared to the SEIR model (Figure 10). The fractional order model provides slower dynamics and more gradual changes in the compartments.
Examples of diseases with a latent period (the time that passes between being exposed to something that can cause disease and having symptoms) are TB, Marek’s disease virus (MDV) poultry infection, chicken pox, Epstein Barr virus (spread through saliva), Herpesviridae (a large family of DNA viruses that cause infection and certain diseases).
Van den Driessche et al. [104] discussed the SEIRS models in epidemiology in 1999 and discussed its stability analysis. Zheng et al. [105] discussed the SEIRS models with nonlinear incidence rate function in 2017 and discussed its stability analysis.

5.4. The SVIR Model

The Susceptible-Vaccinated-Infection-Recovered (SVIR) model is an extended SIR model in which there is a compartment for the vaccination of susceptible population [106,107]. A schematic representation of the SVIR model is given in Figure 13. This model uses the following system of differential equations:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) ν S ( t ) d S ( t ) ,   0 C D t ρ V ( t ) = ν S ( t ) β 1 V ( t ) I ( t ) γ 1 V ( t ) d V ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) + β 1 V ( t ) I ( t ) γ I ( t ) d I ( t )   0 C D t ρ R ( t ) = γ 1 V ( t ) + γ I ( t ) d R ( t ) ,
where b is the birth rate, β and β 1 are the transmission rate from susceptible and vaccination, respectively, ν is the vaccination rate, γ and γ 1 are the recovery rate from vaccination and infected, respectively, and d is the death rate. The initial conditions S ( 0 ) > 0 , V ( 0 ) 0 , I ( 0 ) > 0 , R ( 0 ) 0 .
Disease free equilibrium point
To find the equilibrium point put the system of Equation (86) equal to zero, i.e., we obtain
b β S I ν S d S = 0 , ν S β 1 V I γ 1 V d V = 0 , β S I + β 1 V I γ I d I = 0 γ 1 V + γ I d R = 0 .
For the disease free equilibrium point, put I = 0 , we obtain three equations
b ν S d S = 0 ,
ν S γ 1 V d V = 0 ,
γ 1 V d R = 0 .
On solving the above equations, we obtain S = b d + ν , V = b ν ( d + ν ) ( d + γ 1 ) , R = b ν γ 1 d ( d + ν ) ( d + γ 1 ) . Hence, the disease-free equilibrium point for the SVIR model is
( S 0 , V 0 , I 0 , R 0 ) = b d + ν , b ν ( d + ν ) ( d + γ 1 ) , 0 , b ν γ 1 d ( d + ν ) ( d + γ 1 ) .
Reproduction number
Here, the disease class is
  0 C D t ρ I ( t ) = β S ( t ) I ( t ) + β 1 V ( t ) I ( t ) γ I ( t ) d I ( t ) .
Now, F = β S I + β 1 V I and V = ( γ + d ) I .
The matrices F = [ β S 0 + β 1 V 0 ] and V = [ γ + d ] at the disease free equilibrium point given in (90) The eigen-value of the matrix F V 1 is β S 0 + β 1 V 0 γ + d only.
Hence, the basic reproduction number, R 0 = β S 0 + β 1 V 0 γ + d where the value of S 0 and V 0 is given in (90).
Graphical simulation
The graphical representation of the solution of the SVIR model is presented in Figure 14. The figure shows that the fraction of the susceptible population decreases fast in the beginning due to the effects of vaccination and the onset of infections. A small proportion of the susceptible population is infected, so the peak of infectivity is small. It is observed that the recovery phase progresses rapidly as there is recovery from both infection and due to the immunity gained through vaccination. It can be seen in the graph that for the fractional order, the changes in the susceptible and recovered curves are slower than in the integer order.
Chauhan et al. [108] discussed the SIR model with vaccination in 2014 and discussed its stability analysis. Nezihal et al. [109] in 2021 introduced the fractional order SVIR model for COVID-19 under the Caputo derivative. Mahayana [110] discussed the SIRV Model for the disease COVID-19 in 2022 and performed its Lyapunov stability analysis.

5.5. SIRD Model

The only difference between the Kermack–McKendrick SIR model and the SIRD model is the addition of a new compartment D for the individuals who died due to the infection. A schematic representation of the SIRD model is given in Figure 15. This model was studied by Chatterjee [111] in 2021 on COVID-19 disease.
The rates of change in the populations are given by the following system of equations
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) ( γ + μ ) I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,   0 C D t ρ D ( t ) = μ I ( t ) ,
where β is the transmission rate of infection, γ is the recovery rate and μ is the death rate due to infection. Total population is constant, i.e., S + I + R + D = 1 .
Disease free equilibrium point The disease-free equilibrium point is ( S 0 , I 0 , R 0 , D 0 ) = ( b d , 0 , 0 , 0 ) .
Reproduction number
The next-generation matrix method is used to find the reproduction number at the disease-free equilibrium point. Here, the disease class is I only. Therefore, F = [ β S I ] and V = [ ( γ + μ + d ) I ] . Evaluating the matrices F = β b d and V = [ γ + μ + d ] at the disease free equilibrium point, we obtain
F V 1 = β b d ( γ + μ + d ) .
The spectral radius of the matrix F V 1 is the basic reproduction number, R 0 = β b d ( γ + μ + d ) .
Graphical simulation
Figure 16 provides a graphical representation of the solution of the SIRD model. In the graphical representation of the SIRD model, the susceptible curve rapidly declines as the susceptible population becomes infected and then rises due to the impact of the natural births. The infectious curve rises to a peak and then declines as infectives recover or die. The recovered curve representing the fraction of individuals who gain immunity, increases and tends to become steady. The curve representing deaths due to infection gradually rises and becomes stable as all other fractions of the population tend to stabilize.
Koca [40] formulated the SIRD Model in 2018 for the spread of Ebola virus and used the method Atangana–Baleanu derivative in Caputo (ABC). Dokuyucu and Dutta [41] gives a fractional order SIRD model in 2020 for Ebola Virus with the new Caputo fractional derivative without a singular kernel. Singh [42] developed the SIRD epidemic model for the disease Ebola in 2020. In this paper, the Grunwald–Letnikov (GL) fractional derivative has been used for the fractional Ebola virus model.

5.6. Maternal Immunity

A new compartment M for maternal immunity is added in the SIR model. Maternal immunity can be modeled in two ways as shown in Figure 17.
The two ways of modeling maternal immunity are illustrated through the two models MSIR and MSEIR, respectively.

5.6.1. The MSIR Model

The MSIR model [112] is an extension of the classic SIR model, which includes a compartment for individuals with maternal immunity (M) as illustrated in the Figure 17a. According to this model, infants with passive immunity from their mothers are born immune to infectious diseases. A schematic representation of the MSIR model is given in Figure 18.
The transmission dynamics of the MSIR model is described as
  0 C D t ρ M ( t ) = b θ M ( t ) d M ( t )   0 C D t ρ S ( t ) = θ M ( t ) β S ( t ) I ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
where b is the birth rate, θ is the loss of temporal immunity period, d represents the natural mortality rate, β is the rate of transmission of infection and γ is the recovery rate.
Disease free equilibrium point
By equating the system (93) equal to zero and for the disease-free equilibrium point and putting I = 0 , the system (93) reduces to
b θ M d M = 0 ,
θ M d S = 0 ,
d R = 0 ,
This gives R = 0 , M = b d + θ , and S = b θ d ( d + θ ) .
Hence, the disease free equilibrium point is ( M 0 , S 0 , I 0 , R 0 ) = b d + θ , b θ d ( d + θ ) , 0 , 0 .
Reproduction number
Since the disease class is the same as the SIR model, the basic reproduction number at the given equilibrium point is given as, R 0 = b β θ d ( d + θ ) ( d + γ ) .
Graphical simulation
Figure 19 provides a graphical representation of the solution of the MSIR model. The graph shows that the maternal immunity curve is increasing over time because a portion of newborns is directly added to the maternal immunity compartment. This means that each new birth contributes to the maternal immunity group, leading to a rising trend in the curve. The fractional model shows comparatively lower peaks in everything.
By using the Laplace–Adomian decomposition approach in 2014, Mafuta et al. [113] approximate the series solution of the MSIR epidemic model.

5.6.2. MSEIR Model

The MSIR Model discussed in Section 5.6.1 has additional compartments added to make the fundamental epidemiological models more realistic. This model is appropriate for illnesses like measles, varicella, rubella, or mumps, among others, that result in lifelong immunity.
Due to the fact that only susceptible mothers may give birth to susceptible offspring, the susceptible population grows at a rate of b S . The other newborns b ( 1 S ) enter the passively immune class M because all other classes were exposed to the disease and have immunity. At various rates, people can be moved from one compartment to another. A schematic representation of the MSEIR model is given in Figure 20. The total population is constant, i.e., M + S + E + I + R = 1 .
The following set of fractional order differential equations determines the flow of the infection:
  0 C D t ρ M ( t ) = b ( 1 S ( t ) ) ( θ + d ) M ( t ) ,   0 C D t ρ S ( t ) = b S ( t ) + θ M ( t ) α S ( t ) I ( t ) d S ( t ) ,   0 C D t ρ E ( t ) = α S ( t ) I ( t ) ( ϵ + d ) E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) ( γ + d ) I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) .
By substituting S = 1 M E I R , the above system (97) reduces to
  0 C D t ρ M ( t ) = b ( M + E + I + R ) ( t ) ( θ + d ) M ( t ) ,   0 C D t ρ E ( t ) = α ( 1 M ( t ) E ( t ) I ( t ) R ( t ) ) I ( t ) ( ϵ + d ) E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) ( γ + d ) I ( t )   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) .
Disease free equilibrium point
For the equilibrium point put the system (98) equal to zero, we obtain
b ( M + E + I + R ) ( θ + d ) M = 0 , α ( 1 M E I R ) I ( ϵ + d ) E = 0 , ϵ E ( γ + d ) I = 0 , γ I d R = 0 .
For the disease-free equilibrium point put I = 0 , E = 0 , we obtain
b ( M + R ) ( θ + d ) M = 0 ,
d R = 0 ,
which give R = 0 , M = 0 . As S = 1 M E I R , so the value of S = 1 . Hence, the disease-free equilibrium point is ( M 0 , S 0 , E 0 , I 0 , R 0 ) = ( 0 , 1 , 0 , 0 , 0 ) .
Reproduction number
We find the reproduction number at the disease-free equilibrium point using the next-generation matrix method. Here, the disease classes are
  0 C D t ρ E ( t ) = α ( 1 M ( t ) E ( t ) I ( t ) R ( t ) ) I ( t ) ( ϵ + d ) E ( t ) ,   0 C D t ρ I ( t ) = ϵ E ( t ) ( γ + d ) I ( t ) .
Now, F = α ( 1 M E I R ) I 0 and V = ( ϵ + d ) E ϵ E + γ I + d I .
At the disease-free equilibrium point, we have F = 0 α 0 0 and V = ϵ + d 0 ϵ γ + d .
F V 1 = α ϵ ( ϵ + d ) ( γ + d ) α ( γ + d ) 0 0 .
The dominant eigenvalue of F V 1 is α ϵ ( ϵ + d ) ( γ + d ) .
Hence, the basic reproduction number, R 0 = α ϵ ( ϵ + d ) ( γ + d ) .
Graphical simulation
The graphical solution of the MSEIR model is given in Figure 21. The graph shows that new births contribute to both maternal immunity and susceptible compartments, leading to a slower rise in the maternal immunity curve compared to the MSIR model (Figure 19). It is observed that there is a delay in the rise of the maternal immunity curve due to the slow initial accumulation of maternal immunity. In the fractional order model, this change is more gradual than in the integer order model.
El-Doma [114] discussed the MSEIR age-structured epidemic model in 2006 and discussed its stability analysis. Ibrahim et al. [115] discussed the MSEIR epidemic model in 2014 and found the solution using the Homotopy Analysis Method. Almeida et al. [52] formulated an epidemiological MSEIR Model using Caputo fractional derivative in 2019. Qureshi and Jan [39] formulated the MSEIR Model of the chickenpox disease and used fractional operators like Caputo, Caputo–Fabrizio and Atangana–Baleanu–Caputo in 2021. Khan et al. [53] discussed the model for the disease Hepatitis B Virus Transmission using the Caputo–Fabrizio (C–F) operator in 2022.

5.7. Non-Linear Incidence Rate Model

In mathematical models, the rate at which a susceptible population becomes infected, representing the transmission of a disease, is referred to as the incidence rate function. Until now we discussed all the models with linear incidence rate β I . In this section, we have illustrated the case of different incidence rates which can be non-linear. There are various types of incidence rates discussed in Table 1 of [116]. We consider a non-linear saturated incidence rate. Firstly, Capasso and Serio [117] introduced the non-linear saturated incidence rate, g ( I ) S into epidemic models after the cholera epidemic in Bari in 1973. The g ( I ) is defined as
β I 1 + υ I ,
which tends to a saturation level when I becomes large. Here, β I measures the force of infection of the disease and 1 1 + υ I ( υ shows the saturation constant), measures the inhibition effect from the behavioral change in the susceptible individuals when their number increases or from the crowding effect of the infective individuals [118].
There are the following hypotheses about the transmission process of the infectious disease and its host population.
  • The transfer rates from one compartment to another are supposed to be proportional to the compartment’s population. The number of contacts of a susceptible per unit of time cannot always increase linearly with I. To make it more realistic the interaction term is of the form g ( I ) S ( t ) , where the function g characterizes the non-linearity property and it is assumed to be at least C 3 ( 0 , 1 ] with g ( 0 ) = 0 and g ( I ) > 0 for 0 < I 1 , with N 0 = 1 the initial total population.
  • The law of mass action remains valid since the group of individuals in the host population is well-mixed and homogeneous.
  • The mode of transmission is horizontal and happens when hosts come into direct contact.
  • There is no latent period and following an infection, all infected individual hosts become infectious.
  • Since there is no immunity loss, there is no transfer from compartment R to compartment S.

5.7.1. Kermack–McKendrick Epidemic Model with Non-Linear Saturated Incidence Rate

A schematic representation of the modified SIR model is given in Figure 22.
The system of differential equation of the modified SIR model [117] is
  0 C D t ρ S ( t ) = b β g ( I ) S ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) g ( I ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
subject to the initial conditions S ( 0 ) 0 , I ( 0 ) > 0 , R ( 0 ) = 0 , ( S 0 , I 0 , R 0 ) R + 3 and g ( I ) = β I 1 + υ I .
Disease free equilibrium point
The disease-free equilibrium point of the modified SIR model is ( S 0 , I 0 , R 0 ) = ( 1 , 0 , 0 ) .
Reproduction number
The next-generation matrix method is used to find the reproduction number at the disease-free equilibrium point. Here, the disease class is
  0 C D t ρ I ( t ) = β S ( t ) I ( t ) 1 + υ I ( t ) ( γ + d ) I ( t ) .
Now, F = β S I 1 + υ I and V = [ ( γ + d ) I ] . The matrix F = [ β ] and V = [ γ + d ] at the disease free equilibrium point.
F V 1 = β γ + d . The eigenvalue of the matrix F V 1 is β γ + d only.
Hence, the basic reproduction number, R 0 = β γ + d .
Graphical simulation
The graphical solution of the SIR model with a non-linear saturated incidence rate is given in Figure 23. The graph shows how the non-linear incidence rate affects the shapes of the curves for susceptible (S), infectious (I), and recovered (R). The curves of the SIR model are influenced by the parameter υ (coincidence rate). The non-linear incidence rate results in a less pronounced peak in infections compared to the linear incidence rate model. A sharp increase in recovery is seen as the infection reaches its maximum. As the parametric value of υ is varied to increase, the peak of infected goes down. The fractional order model results in slower changes in the compartment dynamics.
Ruan [119] in 2003, introduced the global dynamics of an epidemic model with vital dynamics and the non-linear incidence rate of saturated mass action. Xu [118] in 2009 discussed the global stability of a SIR epidemic model with non-linear incidence rate and time delay. Khan et al. [120] in 2018 discussed the complex dynamics of an SEIR epidemic model with saturated incidence rate and treatment.

5.7.2. The SIQR Model

Nirwani [121] studied an SIQR model in 2015 with a saturated incident rate and discussed its stability analysis. Cao [122] defined the SIQR epidemic model in 2019 with standard incidence and discussed the existence and uniqueness of the global positive solution to the stochastic system. In the basic SIR model, there is a new compartment, quarantine Q which assumes that a portion of the infective individuals pass through the Q compartment before going to the recovered compartment. Although I and Q are both infected in this case, the separation is appropriate since it simulates the fact that infectious cases that have been proven to be contagious are forced to self-isolate from the public. The population dynamics of the SIQR model of epidemics is given in Figure 24.
Ahmed et al. [123] studied the SIQR model and this model uses the following system of differential equations:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) 1 + υ I ( t ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) 1 + υ I ( t ) ( γ + Ω + d + μ ) I ( t ) ,   0 C D t ρ Q ( t ) = Ω I ( t ) ( d + μ + φ ) Q ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) + φ Q ( t ) d R ( t ) ,
where υ represents the coincidence rate, Ω is the rate at which infected individuals are quarantined, φ is the rate at which quarantined individuals recover, γ represents the rate at which infected individuals recover, μ and μ represent the disease-related death rate in I and Q, respectively. The initial conditions are S ( 0 ) = S 0 0 , I ( 0 ) = I 0 0 , Q ( 0 ) = Q 0 0 , R ( 0 ) = R 0 0 .
The saturation incidence rate, developed by Esteva and Matias [124] in 2001, is represented as β S I 1 + υ I in this case.
Disease free equilibrium point
For the equilibrium points, put the above system of differential Equation (105) equals to zero, i.e.,   0 C D t ρ S ( t ) =   0 C D t ρ I ( t ) =   0 C D t ρ Q ( t ) =   0 C D t ρ R ( t ) = 0 . To find the disease-free equilibrium point put I = 0 , the system (105) reduces to
b d S = 0 , ( d + μ + φ ) Q = 0 , φ Q d R = 0 .
On solving, we obtain S = b d , Q = 0 , R = 0 . Hence, the disease-free equilibrium point is ( S 0 , I 0 , Q 0 , R 0 ) = ( b d , 0 , 0 , 0 ) .
Reproduction number
Using the next generation matrix method we find the reproduction number at the disease-free equilibrium point. Here, the disease classes are
  0 C D t ρ I ( t ) = β S ( t ) I ( t ) 1 + υ I ( t ) ( γ + Ω + d + μ ) I ( t ) ,   0 C D t ρ Q ( t ) = Ω I ( t ) ( d + μ + φ ) Q ( t ) .
Now, F = β S I 1 + υ I 0 and V = ( γ + Ω + d + μ ) I Ω I + ( d + μ + φ ) Q . At the disease-free equilibrium point, the matrices are
F = b β d 0 0 0 and V = γ + Ω + d + μ 0 Ω d + μ + φ . Hence
F V 1 = b β d ( γ + Ω + d + μ ) 0 0 0 .
The eigen-values of F V 1 are λ 1 = 0 and λ 2 = b β d ( γ + Ω + d + μ ) .
Since the dominant eigen-value (spectral radius) is b β d ( γ + Ω + d + μ ) . Hence, the basic reproduction number is R 0 = b β d ( γ + Ω + d + μ ) .
Graphical Simulation
The solution for the SIQR model is graphically presented in Figure 25. In this model, a fraction of the infected individuals are quarantined and recovered, which effectively reduces the peak of the infected curve. It is observed that the curve of quarantine rises as infection increases, attains the peak and then declines as the infection reduces. The non-linear incidence rate with the parameter υ significantly impacts the dynamics of the SIQR model. As the parametric value of υ increases then the peak of infected and quarantine goes down. The impact of fractional order can be seen in the susceptible and recovered curves.

6. SIR Model with Delay

In mathematics, a delay differential equation (DDE) is a type of differential equation in which the derivative of the unknown function depends not only on the current state but also on the history of the system.
Beretta and Takeuchi [125] discussed the SIR epidemic model with time delays in 1995 and discussed its global stability. Ma et al. [126] formulated the SIR epidemic model with time delay in 2004 and discussed its global stability. Wang [127] formulated the SEIRS epidemic model with time delay in 2002 and discussed its global stability. Eyaran [128] discussed the SEIR model with delay differential equation and discussed the stability of the equilibrium points in 2019. The SIRC epidemic model with a time delay for the disease COVID-19 was discussed by Rihan et al. [129] in 2020. Liu et al. [130] studied the time delay model for the diseases COVID-19 in 2020. Devipriya [131] studied the SEIR Model for the COVID-19 epidemic using a delay differential equation and discussed its stability analysis in 2021. The time-delay differential equation in its generic form is
d x ( t ) d t = f ( t , x ( t ) , x t ) ,
where x ( t ) R n and x t = { x ( τ ) : τ t ; t , τ R } R n represents the trajectory of the solution in the past.
There are two types of delay in the differential equation, i.e., discrete delay (lag) and distributed delay.
Delay modeling in epidemiology is crucial in a variety of fields, including medicine, engineering, chemistry, physics, mathematics, economics, and many more [132]. The delay is the amount of time infected people need to spread the infection. Every infectious disease has a built-in lag time. Infectious illness dynamics are slowed down by artificial delays as well. A delay is a time lag, which can never be negative in either situation. It cannot even be zero because it usually takes some time for a communicable disease to spread throughout a community. Delay is, therefore, a more feasible course of action, which is good. Delay modeling is a more general form of deterministic modeling. The delayed analysis is a more effective method to curb the spread of infectious diseases since it is closer to reality, more authentic, and more genuine.
Let us use coronavirus as an example. The only preventative measure is a delay strategy, such as wearing a face mask, public holidays, travel restrictions, and social isolation. Other than vaccinations, medications, etc., these are the most effective methods for controlling the infection.
The time delay effect is included in a more general epidemic model for both infected and recovered people. This suggests that there should be a time delay throughout the disease incubation period so that susceptible persons are not instantaneously infected whenever they come in contact with infected individuals. Individuals who have recovered from infection can also experience the time delay effect, which means that at any time, infected people do not instantly become recovered people. Maleewong [133] talked about the time delay epidemic model for COVID-19 in his article from 2020, which we can see here.
The SIR model can then be generalized with discrete time delay effects as shown by the following:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t τ ) d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) I ( t τ ) γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
where τ represents the effects of time delay for the infective individuals.
However, it may be more realistic to assume that τ is a distributed parameter [125]. The SIR model with distributed time delay effect is shown by the following differential equations:
  0 C D t ρ S ( t ) = b β S ( t ) 0 f ( τ ) I ( t τ ) d τ d S ( t ) ,   0 C D t ρ I ( t ) = β S ( t ) 0 f ( τ ) I ( t τ ) d τ γ I ( t ) d I ( t ) ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) ,
where f ( τ ) is the fraction of the vector population in which the time taken to become infectious is τ . It is assumed that f ( τ ) is non-negative, square-integrable on [ 0 , ) and satisfies
0 f ( τ ) d τ = 1 , 0 τ f ( τ ) d τ < .
Disease free equilibrium point
The equilibrium points of a system with time delays are identical to those in a system with no delays [134]. Hence, the disease-free equilibrium point for the delayed SIR model is ( S 0 , I 0 , R 0 ) = ( b d , 0 , 0 ) .
Reproduction number
The basic reproduction is the same as the SIR model as the disease class is the same as the SIR model with delay τ . Hence, the basic reproduction number is, R 0 = β b d ( γ + d ) .
Graphical simulation
The graphical solution of the SIR model with delay is illustrated in Figure 26. Figure 26a shows the dynamics of the SIR model remain similar generally, but adding a delay causes an apparent change in the infectious population’s peak. The delay affects the timing of the infection spread, causing the peak of the infected population to occur later and the peak is small compared to the model without delay. However, the fractional-order model with delay exhibits a higher infection peak compared to the model without delay as shown in Figure 26b.

7. Stochastic SIR Model

The deterministic method has certain drawbacks when it comes to mathematically simulating the spread of an infectious disease, and it can be challenging to predict the system’s dynamics with any degree of accuracy. Since stochastic differential equation models provide a degree of realism that is somewhat higher than that of their deterministic counterparts, they are widely used in many applied scientific fields, including infectious dynamics [135]. As a result of continuous environmental variations, the parameters used in the modeling of biological systems are not absolute constants; they constantly vary around an average value.
Here, we transform the deterministic SIR model into a corresponding stochastic problem by incorporating white noise for the parameter β of the system (42). Let the fluctuation be in the transmission parameter β , i.e.,
β β + σ d B ( t ) d t ,
where B ( t ) is the white noise with B ( 0 ) = 0 , and σ is the intensity of the white noise. The stochastic version corresponding to the deterministic Caputo fractional model (42) is defined as follows:
  0 C D t ρ S ( t ) = b β S ( t ) I ( t ) d S ( t ) σ S ( t ) I ( t ) d B ( t ) d t ,   0 C D t ρ I ( t ) = β S ( t ) I ( t ) γ I ( t ) d I ( t ) + σ S ( t ) I ( t ) d B ( t ) d t ,   0 C D t ρ R ( t ) = γ I ( t ) d R ( t ) .
Disease free equilibrium point
The disease-free equilibrium point for the stochastic SIR model is the same as the deterministic SIR model, i.e., ( S 0 , I 0 , R 0 ) = ( b d , 0 , 0 ) .
Reproduction number
The basic reproduction number [136] of an epidemic model R 0 is given by
R 0 = 0 b ( a ) F ( a ) d a ,
where b ( a ) is the average number of newly infected individuals by infected individuals if it is infectious during all the time between 0 and a. F ( a ) is the probability of newly infected individuals continuously infecting during the time interval between 0 and a. This is also called the underlying survival probability.
For the SIR model, we consider b ( a ) = b β d and F ( a ) = e ( d + γ ) a given in [137]. By putting these values in (112), we obtain the reproduction number,
R 0 = b β d ( d + γ ) .
Graphical simulation
Using Lagrange’s two-point numerical approximation (38) for the stochastic Caputo fractional SIR model (111),
S n + 1 = S 0 + h ρ Γ ( ρ + 2 ) m = 1 n ( b β S n I n d S n ) × ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ Γ ( ρ + 2 ) m = 1 n ( b β S n 1 I n 1 d S n 1 ) ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) + h ρ 1 2 Γ ( ρ + 2 ) m = 1 n ( σ S n I n ) ( B ( x n ) B ( x n 1 ) ) × ( n m + 1 ) ρ ( n m + 2 + ρ ) ( n m ) ρ ( n m + 2 + 2 ρ ) h ρ 1 2 Γ ( ρ + 2 ) m = 1 n ( σ S n 1 I n 1 ) ( B ( x n 1 ) B ( x n 2 ) ) × ( n m + 1 ) ( ρ + 1 ) ( n m ) ρ ( n m + 1 + ρ ) ,
and similarly for I and R.
Figure 27 displays the graphical representation of the solution for the stochastic Caputo fractional SIR model. The graphical interpretation of a stochastic SIR model provides insights into how an infectious disease might spread in a population when randomness and variability are taken into account. As we can see, deterministic models produce smooth curves while in stochastic there are noisy lines. The stochastic SIR model considers the probabilistic nature of disease transmission and recovery. There is more noise in the susceptible and infected curves than the recovered curve because the fluctuation is only in the transmission rate ( β ).

8. Methods

This work falls under the category of scoping review. The methodology adopted for the review and flow diagram are explained in this section.
  • Eligibility criteria: The review underscores the importance of these models in the realm of epidemiology, demonstrating their efficacy in analyzing and predicting the spread of contagious diseases. Studies considered for inclusion are peer-reviewed journal articles, conference proceedings, and book chapters that utilize fractional order models in the context of diseases such as hepatitis B, COVID-19, influenza, and other epidemics. The review focuses on research that reports on the transmission dynamics, accuracy, and effectiveness of these models, including graphical simulations, and those that provide a comparative analysis with deterministic, stochastic, and integer-order models. The review excludes real-world data in favor of numerical simulations to focus on the theoretical examination of the fractional order models.
  • Information Sources: Research articles and books on epidemiology specifically mathematical modeling based on the Kermack–McKendrick type models. Google Scholar serves as a tool to establish a chronological database, facilitating the search for both historical studies and anticipated future developments in a given field.
  • Search Strategy: We initiated our research by examining the oldest sources such as the works of Kermack and McKendrick, as well as papers by Hethcote. Our investigation commenced with this foundational research, which gave us an understanding of the basic ideas and uses of epidemiological modeling. This approach allowed us to build a comprehensive understanding of the evolution of these models from their development to their contemporary applications in epidemiology.
  • Selection process: Keywords for the selection process were as follows:
    • Epidemiology
    • Mathematical modeling
    • SIR model: Deterministic and Stochastic
    • Fractional epidemic models
    The research articles with these keywords are searched through various database searching websites.
  • Data collection process: Data collection was collaboratively carried out by researchers who consulted regularly to ensure methodological consistency and accuracy. We carefully gathered data from a wide range of sources to ensure a robust and comprehensive dataset. This included sourcing information from peer-reviewed journals, which provided high-quality, verified research findings. To ensure the accuracy of the data collected, the researchers cross-checked each other’s extractions. If any differences were found, we worked together to reach an agreement through careful conversation.
  • Data items:
    (a)
    Data were sought for research articles based on the mathematical model. Approximately 250 papers were reviewed, from which specific papers proved particularly valuable for the following models:
    TopicArticles Reviewed
    Epidemiology44
    SIS/IR5
    SIR/SIRS10
    SEIR/SEIRS13
    SVIR6
    SIRD4
    MSIR/MSEIR6
    Fractional Epidemic Models30
    Stochastic Model5
    Epidemic Models with Delay8
    Non-Linear Incidence Rate Models6
    (b)
    The value of the parameters b , d were taken from https://www.macrotrends.net/global-metrics/countries/IND/india/birth-rate accessed on 1 January 2024; the values of the parameters β , γ were taken from [35]; the rest of the parameters were estimated to validate the model.
    We developed the corresponding fractional model and simulated graphs for the integer and fractional order models so that a comparison can be made.
  • Synthesis methods:
    Data Extraction: Relevant data from each study were systematically extracted and recorded in a structured format. This included key model parameters, outcome measures, study design characteristics, and simulation results.
    Categorization and Classification: Extracted data were categorized based on the type of epidemiological model (e.g., SEIR, SIRD, SVIR, MSIR), and type of analysis (deterministic, stochastic, fractional order). This classification helped in organizing the data for more straightforward comparison and synthesis.
    Normalization or Standardization: Depending on the variables’ distributions and scale, normalization (scaling to a common range) or standardization (scaling to have zero mean and unit variance) may be necessary. The parameters used in the article are normalized to provide a general overview of the model.

8.1. Flow Diagram

After duplicates were removed, a total of 300 citations were identified from searches of electronic databases and review article references. Based on the title and the abstract, 40 were excluded, with 260 full-text articles to be retrieved and assessed for eligibility. Of these, 123 full-text articles were excluded because of not getting any relevant things aligned with our work. The remaining 137 studies were considered eligible for this review. The flow diagram is shown in the Figure 28.
Study characteristics: In this scoping review, the sources of evidence were evaluated based on several key characteristics to ensure a thorough and reliable analysis of the epidemiological models discussed. Here is an elaborate explanation of these characteristics:
Study Design: Type of study (e.g., theoretical, empirical, simulation-based), the model used (e.g., SEIR, SIRD, SVIR, MSIR), and order of the model (integer or fractional).
Model Parameters: Key parameters used in the models (e.g., transmission rate, recovery rate, vaccination rate, exposed rate, etc.), their values, and sources of these parameters.
Numerical Methods: Numerical methods and techniques used to solve the models (e.g., Euler’s Method, Lagrange’s two-step method, stochastic simulation algorithms).
Results of syntheses:
Graphical Representation: Graphical simulations and visualizations are employed to illustrate the behaviors and characteristics of different models. Transmission dynamics, encompassing real-life scenarios like vaccination, quarantine, exposure, non-linear incidence rates, and uncertain parameters, are elucidated through pictorial representations. Their effects are demonstrated with numerical simulations and graphs, aiding readers in systematically comprehending epidemiological models and in formulating their own.
Comparative Analysis: A comparison between integer order and fractional order models is conducted through numerical simulation. Research has demonstrated that fractional order models can yield better fits to empirical data than the traditional integer order models by incorporating memory and history effects, improving the accuracy of long-term predictions.
Fractional models incorporate memory effects, allowing the model to account for past states of the system. This is important in infectious disease modeling, where the disease’s spread can be influenced by past infections and recovery patterns. The findings of the review will assist in directing public health policy by highlighting which models are most effective in different scenarios. This can lead to better preparedness and response strategies for future epidemics and pandemics.

8.2. Discussion

This study concentrates on the types of models rather than the methods applied to solve them. We utilized a singular fractional operator, the Caputo fractional derivative, alongside Lagrange’s two-point numerical method for all deterministic models. Various other fractional derivative operators and numerical methods documented in the literature could be used to develop new models. Despite attempts to perform an exhaustive evaluation, there are multiple factors that could have influenced the quality of this work. Some of these factors are outlined below:
  • Data Availability and Quality: The effectiveness of any model heavily depends on the availability and quality of data used for calibration and validation. Limitations in data availability or reliability may affect the accuracy and generalizability of the findings presented in a review.
  • Parameter Estimation Challenges: Variability in parameter values across different populations or diseases could impact the model’s predictive capability.
  • Model Validation: Validating the performance of fractional order models against real-world datasets remains challenging.
  • Selection Bias: The research articles were leading in many directions we limited ourselves to a particular scope of study.
  • Time Constraints: Due to time constraints, the review process may not have encompassed all relevant literature up to the most current publications.
The results of this review underscore the potential of the fractional SIR model to enhance our understanding and management of infectious diseases. By integrating these findings into practice, policy, and future research, we can improve public health outcomes and better prepare for future epidemics and pandemics.

8.3. Registration

The scoping review is registered on Open Science Framework (OSF) and can be found at https://doi.org/10.17605/OSF.IO/DJ45R accessed on 31 July 2024.

9. Conclusions

Numerous epidemiological models that are used to analyze contagious diseases have been studied in the literature. However, efforts have been made here for the compilation of various compartmental models that have been defined in epidemiology. Through this study, young researchers will be able to gain a comprehensive understanding of the different models and the articles that have used these models. The review of fractional models will help them better understand the many operators that may be applied to generalize these models as well as the differences between the integer and fractional order models.
The paper underscores the importance of these models in the realm of epidemiology, demonstrating their efficacy in analyzing and predicting the spread of contagious diseases. It also opens up avenues for future research, suggesting that these models could be fitted to real data available for various diseases. This could potentially lead to more accurate predictions and better disease management strategies.

Author Contributions

Conceptualization, R.A. and R.P.A.; methodology, R.A. and P.A.; software, P.A.; validation, R.A., P.A. and R.P.A.; formal analysis, R.A. and P.A.; investigation, P.A.; resources, R.P.A.; data curation, R.A. and P.A.; writing—original draft preparation, P.A.; writing—review and editing, R.A. and R.P.A.; visualization, R.A. and P.A.; supervision, R.A.; project administration, R.P.A.; funding acquisition, R.P.A. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

No new data were created or analyzed in this study. Data sharing is not applicable to this article.

Acknowledgments

The author Pooja Airan is thankful to the Council of Scientific and Industrial Research, India for providing fellowship under research scheme No. 09/0964(16111)/2022-EMR-I.

Conflicts of Interest

The authors declare no conflicts of interest.

References

  1. Last, J.M. A Dictionary Of Epidemiology; Oxford University Press: Oxford, UK, 2001. [Google Scholar]
  2. Hotez, P.J. The rise of neglected tropical diseases in the ‘new Texas’. PLoS Neglected Trop. Dis. 2018, 12, e0005581. [Google Scholar] [CrossRef] [PubMed]
  3. Keeling, M.J.; Rohani, P. Modeling Infectious Diseases in Humans and Animals; Princeton University Press: Princeton, NJ, USA, 2011. [Google Scholar]
  4. Bénéteau, T.; Elie, B.; Sofonea, M.T.; Alizon, S. Estimating dates of origin and end of COVID-19 epidemics. Peer Community J. 2021, 1. [Google Scholar] [CrossRef]
  5. Salje, H.; Tran Kiem, C.; Lefrancq, N.; Courtejoie, N.; Bosetti, P.; Paireau, J.; Andronico, A.; Hozé, N.; Richet, J.; Dubost, C.L. Estimating the burden of SARS-CoV-2 in France. Science 2020, 369, 208–211. [Google Scholar] [CrossRef] [PubMed]
  6. Dokuyucu, M.A.; Celik, E. Analyzing a novel coronavirus model (COVID-19) in the sense of Caputo-Fabrizio fractional operator. Appl. Comput. Math. 2021, 20, 49–69. [Google Scholar]
  7. Anderson, R.M.; May, R.M. Infectious Diseases of Humans: Dynamics and Control; Oxford University Press: Oxford, UK, 1992. [Google Scholar]
  8. Thornley, S.; Bullen, C.; Roberts, M. Hepatitis B in a high prevalence New Zealand population: A mathematical model applied to infection control policy. J. Theor. Biol. 2008, 254, 599–603. [Google Scholar] [CrossRef] [PubMed]
  9. Maynard, J.E.; Kane, M.A.; Hadler, S.C. Global control of Hepatitis B through vaccination: Role of Hepatitis B vaccine in the Expanded Programme on Immunization. Clin. Infect. Dis. 1989, 11, S574–S578. [Google Scholar] [CrossRef] [PubMed]
  10. Pang, J.; Cui, J.A.; Zhou, X. Dynamical behavior of a Hepatitis B virus transmission model with vaccination. J. Theor. Biol. 2010, 265, 572–578. [Google Scholar] [CrossRef] [PubMed]
  11. Mann, J.; Roberts, M. Modelling the epidemiology of Hepatitis B in New Zealand. J. Theor. Biol. 2011, 269, 266–272. [Google Scholar] [CrossRef] [PubMed]
  12. Wang, J.; Pang, J.; Liu, X. Modelling diseases with relapse and nonlinear incidence of infection: A multi-group epidemic model. J. Biol. Dyn. 2014, 8, 99–116. [Google Scholar] [CrossRef]
  13. Wang, J.; Zhang, R.; Kuniya, T. The stability analysis of an SVEIR model with continuous age-structure in the exposed and infectious classes. J. Biol. Dyn. 2015, 9, 73–101. [Google Scholar] [CrossRef]
  14. WHO. World health organization. Responding to Community Spread of COVID-19. Reference WHO/COVID-19/Community _Transmission/2020.1; WHO: Geneva, Switzerland, 2020. [Google Scholar]
  15. Tang, B.; Wang, X.; Li, Q.; Bragazzi, N.L.; Tang, S.; Xiao, Y.; Wu, J. Estimation of the transmission risk of the 2019-nCoV and its implication for public health interventions. J. Clin. Med. 2020, 9, 462. [Google Scholar] [CrossRef] [PubMed]
  16. Tian, H.; Liu, Y.; Li, Y.; Wu, C.H.; Chen, B.; Kraemer, M.U.; Li, B.; Cai, J.; Xu, B.; Yang, Q.; et al. An investigation of transmission control measures during the first 50 days of the COVID-19 epidemic in China. Science 2020, 368, 638–642. [Google Scholar] [CrossRef] [PubMed]
  17. Daley, D.J.; Gani, J. Epidemic Modeling: An Introduction; Cambridge University Press: Cambridge, UK, 2005. [Google Scholar]
  18. Graunt, J. Natural and Political Observations Made upon the Bills of Mortality; Johns Hopkins Press: Baltimore, MD, USA, 1939; Volume 2. [Google Scholar]
  19. Bernoulli, D. Reflexions sur les avantages de l’inoculation. Mercur. Fr. 1760, 141, 173–190. [Google Scholar]
  20. Dietz, K.; Heesterbeek, J. Daniel Bernoulli’s epidemiological model revisited. Math. Biosci. 2002, 180, 1–21. [Google Scholar] [CrossRef] [PubMed]
  21. Snow, J. On the Mode of Communication of Cholera; John Churchill: London, UK, 1855. [Google Scholar]
  22. Budd, W. Typhoid Fever: Its Nature, Mode of Spreading, and Prevention; Longmans, Green: London, UK, 1873. [Google Scholar]
  23. Farr, W. Progress of epidemics. In Second Report of the Registrar General of England and Wales; His majesty’s Stationery Office: London, UK, 1840; pp. 16–20. [Google Scholar]
  24. Hamer, W.H. Epidemic Disease in England: The Evidence of Variability and of Persistency of Type; Bedford Press: London, UK, 1906. [Google Scholar]
  25. Ross, R. An application of the theory of probabilities to the study of a priori pathometry.—Part I. Proc. R. Soc. Lond. Ser. A Contain. Pap. Math. Phys. Character 1916, 92, 204–230. [Google Scholar]
  26. Ross, R.; Hudson, H.P. An application of the theory of probabilities to the study of a priori pathometry.—Part II. Proc. R. Soc. Lond. Ser. A Containing Pap. Math. Phys. Character 1917, 93, 212–225. [Google Scholar]
  27. Kermack, W.O.; McKendrick, A.G. A contribution to the mathematical theory of epidemics. Proc. R. Soc. Lond. Ser. A Containing Pap. Math. Phys. Character 1927, 115, 700–721. [Google Scholar]
  28. Kendall, D.G. Deterministic and stochastic epidemics in closed populations. In Proceedings of the Third Berkeley Symposium on Mathematical Statistics and Probability; University of California Press Berkeley: Oakland, CA, USA, 1956; Volume 4, pp. 149–165. [Google Scholar]
  29. Brauer, F.; Castillo-Chavez, C. Mathematical Models in Population Biology and Epidemiology; Springer: Berlin/Heidelberg, Germany, 2001. [Google Scholar]
  30. Turner, P.E.; Cooper, V.S.; Lenski, R.E. Tradeoff between horizontal and vertical modes of transmission in bacterial plasmids. Evolution 1998, 52, 315–329. [Google Scholar] [CrossRef]
  31. Dietz, K. Epidemics and rumours: A survey. J. R. Stat. Soc. Ser. A (Gen.) 1967, 130, 505–528. [Google Scholar] [CrossRef]
  32. Nakamura, G.M.; Monteiro, A.C.P.; Cardoso, G.C.; Martinez, A.S. Efficient method for comprehensive computation of agent-level epidemic dissemination in networks. Sci. Rep. 2017, 7, 40885. [Google Scholar] [CrossRef] [PubMed]
  33. Rafei, M.; Ganji, D.; Daniali, H. Solution of the epidemic model by homotopy perturbation method. Appl. Math. Comput. 2007, 187, 1056–1062. [Google Scholar] [CrossRef]
  34. Doungmo Goufo, E.F.; Maritz, R.; Munganga, J. Some properties of the Kermack-McKendrick epidemic model with fractional derivative and nonlinear incidence. Adv. Differ. Equ. 2014, 2014, 278. [Google Scholar] [CrossRef]
  35. Okyere, E.; Oduro, F.T.; Amponsah, S.K.; Dontwi, I.K.; Frempong, N.K. Fractional order SIR model with constant population. Br. J. Math. Comput. Sci 2016, 14, 1–12. [Google Scholar] [CrossRef] [PubMed]
  36. Hassouna, M.; Ouhadan, A.; El Kinani, E. On the solution of fractional order SIS epidemic model. Chaos Solitons Fractals 2018, 117, 168–174. [Google Scholar] [CrossRef] [PubMed]
  37. Rostamy, D.; Mottaghi, E. Stability analysis of a fractional-order epidemics model with multiple equilibriums. Adv. Differ. Equ. 2016, 2016, 170. [Google Scholar] [CrossRef]
  38. Oke, M.; Ogunmiloro, O.; Akinwumi, C.; Raji, R. Mathematical modeling and stability analysis of a SIRV epidemic model with non-linear force of infection and treatment. Commun. Math. Appl. 2019, 10, 717–731. [Google Scholar] [CrossRef]
  39. Qureshi, S.; Jan, R. Modeling of measles epidemic with optimized fractional order under Caputo differential operator. Chaos Solitons Fractals 2021, 145, 110766. [Google Scholar] [CrossRef]
  40. Koca, I. Modelling the spread of Ebola virus with Atangana-Baleanu fractional operators. Eur. Phys. J. Plus 2018, 133, 133. [Google Scholar] [CrossRef]
  41. Dokuyucu, M.A.; Dutta, H. A fractional order model for Ebola Virus with the new Caputo fractional derivative without singular kernel. Chaos Solitons Fractals 2020, 134, 109717. [Google Scholar] [CrossRef]
  42. Singh, H. Analysis for fractional dynamics of Ebola virus model. Chaos Solitons Fractals 2020, 138, 109992. [Google Scholar] [CrossRef] [PubMed]
  43. Özalp, N.; Demirci, E. A fractional order SEIR model with vertical transmission. Math. Comput. Model. 2011, 54, 1–6. [Google Scholar] [CrossRef]
  44. Area, I.; Batarfi, H.; Losada, J.; Nieto, J.J.; Shammakh, W.; Torres, Á. On a fractional order Ebola epidemic model. Adv. Differ. Equ. 2015, 2015, 278. [Google Scholar] [CrossRef]
  45. Almeida, R. Analysis of a fractional SEIR model with treatment. Appl. Math. Lett. 2018, 84, 56–62. [Google Scholar] [CrossRef]
  46. Casagrandi, R.; Bolzoni, L.; Levin, S.A.; Andreasen, V. The SIRC model and influenza A. Math. Biosci. 2006, 200, 152–169. [Google Scholar] [CrossRef]
  47. El-Shahed, M.; Alsaedi, A. The fractional SIRC model and influenza A. Math. Probl. Eng. 2011, 2011, 480378. [Google Scholar] [CrossRef]
  48. Agarwal, R.; Kritika; Purohit, S.D.; Mishra, J. A Mathematical Fractional Model to Study the Hepatitis B Virus Infection. In Mathematical Modeling and Soft Computing in Epidemiology; CRC Press: Boca Raton, FL, USA, 2020; pp. 273–290. [Google Scholar]
  49. Arfan, M.; Alrabaiah, H.; Rahman, M.U.; Sun, Y.L.; Hashim, A.S.; Pansera, B.A.; Ahmadian, A.; Salahshour, S. Investigation of fractal-fractional order model of COVID-19 in Pakistan under Atangana-Baleanu Caputo (ABC) derivative. Results Phys. 2021, 24, 104046. [Google Scholar] [CrossRef] [PubMed]
  50. Agarwal, R.; Kritika; Purohit, S.D. Fractional Order Mathematical Model for the Cell Cycle of a Tumour Cell. In Fractional Calculus in Medical and Health Science; CRC Press: Boca Raton, FL, USA, 2020; pp. 129–147. [Google Scholar]
  51. Asfour, H.A.; Ibrahim, M. On the differential fractional transformation method of MSEIR epidemic model. Int. J. Comput. Appl. 2015, 113, 10–16. [Google Scholar]
  52. Almeida, R.; Brito da Cruz, A.; Martins, N.; Monteiro, M.T.T. An epidemiological MSEIR model described by the Caputo fractional derivative. Int. J. Dyn. Control 2019, 7, 776–784. [Google Scholar] [CrossRef]
  53. Khan, M.; Khan, T.; Ahmad, I.; Shah, Z.; Khan, A. Modeling of Hepatitis B virus transmission with fractional analysis. Math. Probl. Eng. 2022, 2022, 6202049. [Google Scholar] [CrossRef]
  54. Qureshi, S.; Atangana, A. Fractal-fractional differentiation for the modeling and mathematical analysis of nonlinear diarrhea transmission dynamics under the use of real data. Chaos Solitons Fractals 2020, 136, 109812. [Google Scholar] [CrossRef]
  55. Atangana, A. Fractal-fractional differentiation and integration: Connecting fractal calculus and fractional calculus to predict complex system. Chaos Solitons Fractals 2017, 102, 396–406. [Google Scholar] [CrossRef]
  56. Atangana, A. Extension of rate of change concept: From local to nonlocal operators with applications. Results Phys. 2020, 19, 103515. [Google Scholar] [CrossRef]
  57. Baleanu, D.; Mohammadi, H.; Rezapour, S. A fractional differential equation model for the COVID-19 transmission by using the Caputo-Fabrizio derivative. Adv. Differ. Equ. 2020, 2020, 99. [Google Scholar] [CrossRef]
  58. Islam, M.R.; Peace, A.; Medina, D.; Oraby, T. Integer versus fractional order SEIR deterministic and stochastic models of Measles. Int. J. Environ. Res. Public Health 2020, 17, 2014. [Google Scholar] [CrossRef] [PubMed]
  59. Pandey, R.M.; Chandola, A.; Agarwal, R. Mathematical model and interpretation of crowding effects on SARS-CoV-2 using Atangana-Baleanu fractional operator. In Methods of Mathematical Modelling; Elsevier: Amsterdam, The Netherlands, 2022; pp. 41–58. [Google Scholar]
  60. Owolabi, K.M.; Pindza, E. A nonlinear epidemic model for tuberculosis with Caputo operator and fixed point theory. Healthc. Anal. 2022, 2, 100111. [Google Scholar] [CrossRef]
  61. Agarwal, R.; Airan, P.; Midha, C. Mathematical analysis of the non-linear dynamics of the bone mineralization. In Mathematical Methods in Medical and Biological Sciences; Singh, H., Ed.; Elsevier: Amsterdam, The Netherlands, 2023. [Google Scholar]
  62. Agarwal, R.; Midha, C. Study and mathematical analysis of the novel fractional bone mineralization model. J. Comput. Anal. Appl. 2024, 33, 289–310. [Google Scholar]
  63. Agarwal, R.; Airan, P.; Sajid, M. Numerical and graphical simulation of the non-linear fractional dynamical system of bone mineralization. Math. Biosci. Eng. 2024, 21, 5138–5163. [Google Scholar] [CrossRef]
  64. Miller, K.S.; Ross, B. An Introduction to the Fractional Calculus and Fractional Differential Equations; Wiley: Hoboken, NJ, USA, 1993. [Google Scholar]
  65. Caputo, M. Linear models of dissipation whose Q is almost frequency independent—II. Geophys. J. Int. 1967, 13, 529–539. [Google Scholar] [CrossRef]
  66. Kilbas, A.A.; Srivastava, H.M.; Trujillo, J.J. Theory and Applications of Fractional Differential Equations; Elsevier: Hoboken, NJ, USA, 2006; Volume 204. [Google Scholar]
  67. Tarasov, V.E.; Tarasova, S.S. Fractional derivatives and integrals: What are they needed for? Mathematics 2020, 8, 164. [Google Scholar] [CrossRef]
  68. Odibat, Z.M.; Momani, S. An algorithm for the numerical solution of differential equations of fractional order. J. Appl. Math. Inform. 2008, 26, 15–27. [Google Scholar]
  69. Atangana, A.; Araz, S.I. New Numerical Scheme with Newton Polynomial: Theory, Methods, and Applications; Academic Press: Cambridge, MA, USA, 2021. [Google Scholar]
  70. Kloeden, P.E.; Platen, E.; Schurz, H. Numerical solution of SDE through Computer Experiments; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2012. [Google Scholar]
  71. Allen, L.J. Stochastic population and epidemic models. Math. Biosci. Lect. Ser. Stochastics Biol. Syst. 2015, 1, 120–128. [Google Scholar]
  72. Allen, L.J. A primer on stochastic epidemic models: Formulation, numerical simulation, and analysis. Infect. Dis. Model. 2017, 2, 128–142. [Google Scholar] [CrossRef]
  73. Alkahtani, B.S.T.; Koca, I. Fractional stochastic SIR model. Results Phys. 2021, 24, 104124. [Google Scholar] [CrossRef]
  74. Atangana, A.; Araz, S.İ. Fractional Stochastic Differential Equations: Applications to COVID-19 Modeling; Springer Nature: Berlin/Heidelberg, Germany, 2022. [Google Scholar]
  75. Garner, M.; Hamilton, S. Principles of epidemiological modelling. Rev. Sci. Tech.-OIE 2011, 30, 407–416. [Google Scholar] [CrossRef]
  76. Hethcote, H.W. Qualitative analyses of communicable disease models. Math. Biosci. 1976, 28, 335–356. [Google Scholar] [CrossRef]
  77. Gran, J.M. Infectious Disease Modelling and Causal Inference; Institutt for Klinisk Medisin: Copenhagen, Denmark, 2011. [Google Scholar]
  78. Wu, J.; Dhingra, R.; Gambhir, M.; Remais, J.V. Sensitivity analysis of infectious disease models: Methods, advances and their application. J. R. Soc. Interface 2013, 10, 20121018. [Google Scholar] [CrossRef]
  79. Diekmann, O.; Heesterbeek, J.A.P.; Metz, J.A. On the definition and the computation of the basic reproduction ratio R0 in models for infectious diseases in heterogeneous populations. J. Math. Biol. 1990, 28, 365–382. [Google Scholar] [CrossRef]
  80. Macdonald, G. The Epidemiology and Control of Malaria; Oxford University Press: London, UK, 1957. [Google Scholar]
  81. Shaikh, A.S.; Shaikh, I.N.; Nisar, K.S. A mathematical model of COVID-19 using fractional derivative: Outbreak in India with dynamics of transmission and control. Adv. Differ. Equ. 2020, 2020, 373. [Google Scholar] [CrossRef]
  82. Nabti, A.; Ghanbari, B. Global stability analysis of a fractional SVEIR epidemic model. Math. Methods Appl. Sci. 2021, 44, 8577–8597. [Google Scholar] [CrossRef]
  83. Van den Driessche, P.; Watmough, J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Math. Biosci. 2002, 180, 29–48. [Google Scholar] [CrossRef] [PubMed]
  84. Van den Driessche, P.; Watmough, J. Further notes on the basic reproduction number. In Mathematical Epidemiology; Springer: Berlin/Heidelberg, Germany, 2008; pp. 159–178. [Google Scholar]
  85. Van den Driessche, P. Reproduction numbers of infectious disease models. Infect. Dis. Model. 2017, 2, 288–303. [Google Scholar] [CrossRef] [PubMed]
  86. Rezapour, S.; Mohammadi, H.; Samei, M.E. SEIR epidemic model for COVID-19 transmission by Caputo derivative of fractional order. Adv. Differ. Equ. 2020, 2020, 490. [Google Scholar] [CrossRef] [PubMed]
  87. Kozioł, K.; Stanisławski, R.; Bialic, G. Fractional-order SIR epidemic model for transmission prediction of COVID-19 disease. Appl. Sci. 2020, 10, 8316. [Google Scholar] [CrossRef]
  88. Alshomrani, A.S.; Ullah, M.Z.; Baleanu, D. Caputo SIR model for COVID-19 under optimized fractional order. Adv. Differ. Equ. 2021, 2021, 185. [Google Scholar] [CrossRef] [PubMed]
  89. Hethcote, H.W. Three basic epidemiological models. In Applied Mathematical Ecology; Springer: Berlin/Heidelberg, Germany, 1989; pp. 119–144. [Google Scholar]
  90. Liu, N.; Fang, J.; Deng, W.; Sun, J.w. Stability analysis of a fractional-order SIS model on complex networks with linear treatment function. Adv. Differ. Equ. 2019, 2019, 1–10. [Google Scholar] [CrossRef]
  91. Liu, A.; Yasin, F.; Afzal, Z.; Nazeer, W. Analytical solution of a non-linear fractional order SIS epidemic model utilizing a new technique. Alex. Eng. J. 2023, 73, 123–129. [Google Scholar] [CrossRef]
  92. Tang, Y.; Huang, D.; Ruan, S.; Zhang, W. Coexistence of limit cycles and homoclinic loops in a SIRS model with a nonlinear incidence rate. SIAM J. Appl. Math. 2008, 69, 621–639. [Google Scholar] [CrossRef]
  93. Greenhalgh, D.; Khan, Q.; Lewis, F. Hopf bifurcation in two SIRS density dependent epidemic models. Math. Comput. Model. 2004, 39, 1261–1283. [Google Scholar] [CrossRef]
  94. Korobeinikov, A.; Wake, G.C. Lyapunov functions and global stability for SIR, SIRS, and SIS epidemiological models. Appl. Math. Lett. 2002, 15, 955–960. [Google Scholar] [CrossRef]
  95. Zhen, J.; Ma, Z.; Han, M. Global stability of an SIRS epidemic model with delays. Acta Math. Sci. 2006, 26, 291–306. [Google Scholar] [CrossRef]
  96. Nakata, Y.; Enatsu, Y.; Muroya, Y. On the global stability of an SIRS epidemic model with distributed delays. In Proceedings of the Conference Publications; American Institute of Mathematical Sciences: Pasadena, CA, USA, 2011; p. 1119. [Google Scholar]
  97. El-Saka, H.; Arafa, A.; Gouda, M. Dynamical analysis of a fractional SIRS model on homogenous networks. Adv. Differ. Equ. 2019, 2019, 144. [Google Scholar] [CrossRef]
  98. Li, M.Y.; Muldowney, J.S. Global stability for the SEIR model in epidemiology. Math. Biosci. 1995, 125, 155–164. [Google Scholar] [CrossRef] [PubMed]
  99. Li, M.Y.; Graef, J.R.; Wang, L.; Karsai, J. Global dynamics of a SEIR model with varying total population size. Math. Biosci. 1999, 160, 191–213. [Google Scholar] [CrossRef] [PubMed]
  100. Korobeinikov, A. Lyapunov functions and global properties for SEIR and SEIS epidemic models. Math. Med. Biol. A J. IMA 2004, 21, 75–83. [Google Scholar] [CrossRef] [PubMed]
  101. El-Sheikh, M.; El-Marouf, S. On stability and bifurcation of solutions of an SEIR epidemic model with vertical transmission. Int. J. Math. Math. Sci. 2004, 2004, 2971–2987. [Google Scholar] [CrossRef]
  102. Sene, N. Numerical methods applied to a class of SEIR epidemic models described by the Caputo derivative. In Methods of Mathematical Modelling; Elsevier: Amsterdam, The Netherlands, 2022; pp. 23–40. [Google Scholar]
  103. Bjørnstad, O.N.; Shea, K.; Krzywinski, M.; Altman, N. The SEIRS model for infectious disease dynamics. Nat. Methods 2020, 17, 557–559. [Google Scholar] [CrossRef] [PubMed]
  104. Van den Driessche, P.; Li, M.; Muldowney, J. Global stability of SEIRS models in epidemiology. Can. Appl. Math. Q. 1999, 7, 409–425. [Google Scholar]
  105. Zheng, L.; Yang, X.; Zhang, L. On global stability analysis for SEIRS models in epidemiology with nonlinear incidence rate function. Int. J. Biomath. 2017, 10, 1750019. [Google Scholar] [CrossRef]
  106. Liu, X.; Takeuchi, Y.; Iwami, S. SVIR epidemic models with vaccination strategies. J. Theor. Biol. 2008, 253, 1–11. [Google Scholar] [CrossRef]
  107. Schlickeiser, R.; Kröger, M. Analytical modeling of the temporal evolution of epidemics outbreaks accounting for vaccinations. Physics 2021, 3, 386–426. [Google Scholar] [CrossRef]
  108. Chauhan, S.; Misra, O.P.; Dhar, J. Stability analysis of SIR model with vaccination. Am. J. Comput. Appl. Math. 2014, 4, 17–23. [Google Scholar]
  109. Gokbulut, N.; Amilo, D.; Kaymakamzade, B. Fractional SVIR model for COVID-19 under Caputo derivative. J. Biometry Stud. 2021, 1, 58–64. [Google Scholar] [CrossRef]
  110. Mahayana, D. Lyapunov Stability Analysis of COVID-19 SIRV Model. In Proceedings of the 2022 IEEE 18th International Colloquium on Signal Processing & Applications (CSPA), Selangor, Malaysia, 12 May 2022; IEEE: Piscataway, NJ, USA, 2022; pp. 287–292. [Google Scholar]
  111. Chatterjee, S.; Sarkar, A.; Chatterjee, S.; Karmakar, M.; Paul, R. Studying the progress of COVID-19 outbreak in India using SIRD model. Indian J. Phys. 2021, 95, 1941–1957. [Google Scholar] [CrossRef] [PubMed]
  112. Somma, S.A.; Akinwande, N.I.; Gana, P.; Abdulrahman, S.; Ashezua, T.T. Modified Maternally-Derived-Immunity Susceptible Infectious Recovered (MSIR) Model of Infectious Disease: Existence of Equilibrium and Basic Reproduction Number. Niger. J. Technol. Res. (NJTR) 2015, 10, 40–43. [Google Scholar] [CrossRef]
  113. Mafuta, P.; Mushanyu, J.; Nhawu, G. Extension of the analysis of the MSIR model by applying the Laplace-Adomian decomposition method. World J. Model. Simul. 2014, 10, 185–192. [Google Scholar]
  114. El-Doma, M. Stability analysis for an MSEIR age-structured epidemic model. Dyn. Contin. Discret. Impuls. Syst. Ser. A 2006, 13, 85. [Google Scholar]
  115. Ibrahim, M.; Raji, I.; Aladesuyi, A.; Nwagwo, A. On the Homotopy analysis method of MSEIR epidemic model. IOSR J. Appl. Phys. 2014, 6, 55–61. [Google Scholar] [CrossRef]
  116. Khan, M.A.; Ismail, M.; Ullah, S.; Farhan, M. Fractional order SIR model with generalized incidence rate. AIMS Math. 2020, 5, 1856–1880. [Google Scholar] [CrossRef]
  117. Capasso, V.; Serio, G. A generalization of the Kermack-McKendrick deterministic epidemic model. Math. Biosci. 1978, 42, 43–61. [Google Scholar] [CrossRef]
  118. Xu, R.; Ma, Z. Global stability of a SIR epidemic model with nonlinear incidence rate and time delay. Nonlinear Anal. Real World Appl. 2009, 10, 3175–3189. [Google Scholar] [CrossRef]
  119. Ruan, S.; Wang, W. Dynamical behavior of an epidemic model with a nonlinear incidence rate. J. Differ. Equ. 2003, 188, 135–163. [Google Scholar] [CrossRef]
  120. Khan, M.A.; Khan, Y.; Islam, S. Complex dynamics of an SEIR epidemic model with saturated incidence rate and treatment. Phys. A Stat. Mech. Its Appl. 2018, 493, 210–227. [Google Scholar] [CrossRef]
  121. Nirwani, N.; Badshah, V.; Khandelwal, R. Dynamical study of an SIQR model with saturated incidence rate. J. Comput. Math. Sci. 2015, 6, 564–570. [Google Scholar] [CrossRef]
  122. Cao, Z.; Feng, W.; Wen, X.; Zu, L.; Cheng, M. Dynamics of a stochastic SIQR epidemic model with standard incidence. Phys. A Stat. Mech. Its Appl. 2019, 527, 121180. [Google Scholar] [CrossRef]
  123. Ahmed, N.; Raza, A.; Rafiq, M.; Ahmadian, A.; Batool, N.; Salahshour, S. Numerical and bifurcation analysis of SIQR model. Chaos Solitons Fractals 2021, 150, 111133. [Google Scholar] [CrossRef]
  124. Esteva, L.; Matias, M. A model for vector transmitted diseases with saturation incidence. J. Biol. Syst. 2001, 9, 235–245. [Google Scholar] [CrossRef]
  125. Beretta, E.; Takeuchi, Y. Global stability of an SIR epidemic model with time delays. J. Math. Biol. 1995, 33, 250–260. [Google Scholar] [CrossRef]
  126. Ma, W.; Song, M.; Takeuchi, Y. Global stability of an SIR epidemicmodel with time delay. Appl. Math. Lett. 2004, 17, 1141–1145. [Google Scholar] [CrossRef]
  127. Wang, W. Global behavior of an SEIRS epidemic model with time delays. Appl. Math. Lett. 2002, 15, 423–428. [Google Scholar] [CrossRef]
  128. Eyaran, W.E.; Osman, S.; Wainaina, M. Modelling and analysis of SEIR with delay differential equation. Glob. J. Pure Appl. Math. 2019, 15, 365–382. [Google Scholar]
  129. Rihan, F.A.; Alsakaji, H.J.; Rajivganthi, C. Stochastic SIRC epidemic model with time-delay for COVID-19. Adv. Differ. Equ. 2020, 2020, 502. [Google Scholar] [CrossRef] [PubMed]
  130. Liu, X.; Zheng, X.; Balachandran, B. COVID-19: Data-driven dynamics, statistical and distributed delay models, and observations. Nonlinear Dyn. 2020, 101, 1527–1543. [Google Scholar] [CrossRef] [PubMed]
  131. Devipriya, R.; Dhamodharavadhani, S.; Selvi, S. SEIR Model for COVID-19 Epidemic Using Delay Differential Equation. In Proceedings of the Journal of Physics: Conference Series; IOP Publishing: Bristol, UK, 2021; Volume 1767, p. 012005. [Google Scholar]
  132. Naveed, M.; Baleanu, D.; Rafiq, M.; Raza, A.; Soori, A.H.; Ahmed, N. Dynamical behavior and sensitivity analysis of a delayed coronavirus epidemic model. Comput. Mater. Contin. 2020, 65, 225–241. [Google Scholar] [CrossRef]
  133. Maleewong, M. Time delay epidemic model for COVID-19. MedRxiv 2020. [Google Scholar] [CrossRef]
  134. Kumar, A.; Goel, K.; Nilam. A deterministic time-delayed SIR epidemic model: Mathematical modeling and analysis. Theory Biosci. 2020, 139, 67–76. [Google Scholar] [CrossRef] [PubMed]
  135. Ji, C.; Jiang, D.; Shi, N. Multigroup SIR epidemic model with stochastic perturbation. Phys. A Stat. Mech. Its Appl. 2011, 390, 1747–1762. [Google Scholar] [CrossRef]
  136. Heffernan, J.M.; Smith, R.J.; Wahl, L.M. Perspectives on the basic reproductive ratio. J. R. Soc. Interface 2005, 2, 281–293. [Google Scholar] [CrossRef]
  137. Ríos-Gutiérrez, A.; Torres, S.; Arunachalam, V. Studies on the basic reproduction number in stochastic epidemic models with random perturbations. Adv. Differ. Equ. 2021, 2021, 288. [Google Scholar] [CrossRef]
Figure 1. Population dynamics of the SIR model of epidemics.
Figure 1. Population dynamics of the SIR model of epidemics.
Axioms 13 00545 g001
Figure 2. Graphical representation of the SIR model over time showing solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 .
Figure 2. Graphical representation of the SIR model over time showing solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 .
Axioms 13 00545 g002
Figure 3. Population dynamics of the SIS model of epidemics.
Figure 3. Population dynamics of the SIS model of epidemics.
Axioms 13 00545 g003
Figure 4. The SIS model with vital dynamics is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , ω = 0.06 .
Figure 4. The SIS model with vital dynamics is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , ω = 0.06 .
Axioms 13 00545 g004
Figure 5. Population dynamics of the IR model of epidemics.
Figure 5. Population dynamics of the IR model of epidemics.
Axioms 13 00545 g005
Figure 6. The IR model is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with respect to time when γ = 0.33 .
Figure 6. The IR model is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with respect to time when γ = 0.33 .
Axioms 13 00545 g006
Figure 7. Population dynamics of the SIRS model of epidemics.
Figure 7. Population dynamics of the SIRS model of epidemics.
Axioms 13 00545 g007
Figure 8. The SIRS model with vital dynamics is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with regard to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , δ = 0.026 .
Figure 8. The SIRS model with vital dynamics is represented graphically for integer order (solid lines) and fractional order ρ = 0.9 (dashed lines) with regard to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , δ = 0.026 .
Axioms 13 00545 g008
Figure 9. Population dynamics of the SEIR model of epidemics.
Figure 9. Population dynamics of the SEIR model of epidemics.
Axioms 13 00545 g009
Figure 10. The behavior of the SEIR model with vital dynamics over time is graphically shown by dashed lines for fractional order with ρ = 0.9 and solid lines for integer order when b = 0.01675 , d = 0.007473 , α = 1.05 , γ = 0.33 , ϵ = 0.3 .
Figure 10. The behavior of the SEIR model with vital dynamics over time is graphically shown by dashed lines for fractional order with ρ = 0.9 and solid lines for integer order when b = 0.01675 , d = 0.007473 , α = 1.05 , γ = 0.33 , ϵ = 0.3 .
Axioms 13 00545 g010
Figure 11. Population dynamics of the SEIRS model of epidemics.
Figure 11. Population dynamics of the SEIRS model of epidemics.
Axioms 13 00545 g011
Figure 12. The behavior of the SEIRS model with vital dynamics over time is graphically shown by dashed lines for fractional order with ρ = 0.9 and solid lines for integer order when b = 0.01675 , d = 0.007473 , α = 1.05 , γ = 0.33 , ϵ = 0.3 , δ = 0.026 .
Figure 12. The behavior of the SEIRS model with vital dynamics over time is graphically shown by dashed lines for fractional order with ρ = 0.9 and solid lines for integer order when b = 0.01675 , d = 0.007473 , α = 1.05 , γ = 0.33 , ϵ = 0.3 , δ = 0.026 .
Axioms 13 00545 g012
Figure 13. Population dynamics of the SVIR model of epidemics.
Figure 13. Population dynamics of the SVIR model of epidemics.
Axioms 13 00545 g013
Figure 14. Graphical illustration of the SVIR model with vital dynamics over time, featuring solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , β 1 = 0.725 , γ = 0.33 , γ 1 = 0.66 , ν = 0.08 .
Figure 14. Graphical illustration of the SVIR model with vital dynamics over time, featuring solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , β 1 = 0.725 , γ = 0.33 , γ 1 = 0.66 , ν = 0.08 .
Axioms 13 00545 g014
Figure 15. Population dynamics of the SIRD model of epidemics.
Figure 15. Population dynamics of the SIRD model of epidemics.
Axioms 13 00545 g015
Figure 16. The SIRD model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , μ = 0.034 .
Figure 16. The SIRD model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , μ = 0.034 .
Axioms 13 00545 g016
Figure 17. Extended SIR model with maternal immunity compartment, illustrating two different approaches to incorporate newly born individuals. (a) In this type of maternal immunity model, the whole newly born individual enters into the maternal immunity compartment with the rate b because for a few months they have maternal immunity after that they come into the susceptible compartment. (b) In this type of model, not every newborn individual has maternal immunity, some may be susceptible since birth. So, a rate of b S will come in the susceptible class and with the rate b ( 1 S ) will come in the maternal immunity compartment.
Figure 17. Extended SIR model with maternal immunity compartment, illustrating two different approaches to incorporate newly born individuals. (a) In this type of maternal immunity model, the whole newly born individual enters into the maternal immunity compartment with the rate b because for a few months they have maternal immunity after that they come into the susceptible compartment. (b) In this type of model, not every newborn individual has maternal immunity, some may be susceptible since birth. So, a rate of b S will come in the susceptible class and with the rate b ( 1 S ) will come in the maternal immunity compartment.
Axioms 13 00545 g017
Figure 18. Population dynamics of the MSIR model of the epidemics.
Figure 18. Population dynamics of the MSIR model of the epidemics.
Axioms 13 00545 g018
Figure 19. The MSIR model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , θ = 0.05 .
Figure 19. The MSIR model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , θ = 0.05 .
Axioms 13 00545 g019
Figure 20. Population dynamics of the MSEIR model of epidemics.
Figure 20. Population dynamics of the MSEIR model of epidemics.
Axioms 13 00545 g020
Figure 21. The MSEIR model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , α = 1.05 , γ = 0.33 , θ = 0.05 , ϵ = 0.3 .
Figure 21. The MSEIR model with vital dynamics is represented graphically with respect to time for both fractional order ρ = 0.9 (dashed lines) and integer order (solid lines) when b = 0.01675 , d = 0.007473 , β = 1.05 , α = 1.05 , γ = 0.33 , θ = 0.05 , ϵ = 0.3 .
Axioms 13 00545 g021
Figure 22. Population Dynamics of the SIR model of epidemics.
Figure 22. Population Dynamics of the SIR model of epidemics.
Axioms 13 00545 g022
Figure 23. The graph illustrates the SIR model with vital dynamics with a non-linear saturated incidence rate over time, represented by solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , υ = 0.002 .
Figure 23. The graph illustrates the SIR model with vital dynamics with a non-linear saturated incidence rate over time, represented by solid lines for integer order and dashed lines for fractional order with ρ = 0.9 when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , υ = 0.002 .
Axioms 13 00545 g023
Figure 24. Population Dynamics of the SIQR model of epidemics.
Figure 24. Population Dynamics of the SIQR model of epidemics.
Axioms 13 00545 g024
Figure 25. An illustration of the SIQR model with vital dynamics with non-linear incidence rate for a fractional order ρ = 0.9 (dashed lines) and an integer order (solid lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , υ = 0.002 , μ = 0.034 , μ = 0.0006 , Ω = 0.1 , φ = 0.2 .
Figure 25. An illustration of the SIQR model with vital dynamics with non-linear incidence rate for a fractional order ρ = 0.9 (dashed lines) and an integer order (solid lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 , υ = 0.002 , μ = 0.034 , μ = 0.0006 , Ω = 0.1 , φ = 0.2 .
Axioms 13 00545 g025
Figure 26. Graphical representation of the SIR model without delay (dashed lines) and with delay (solid lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 .
Figure 26. Graphical representation of the SIR model without delay (dashed lines) and with delay (solid lines) with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , γ = 0.33 .
Axioms 13 00545 g026
Figure 27. Graphical representation of the stochastic fractional SIR model (noisy lines) and deterministic fractional SIR model (dashed lines) for fractional order ρ = 0.9 with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , σ = 0.01 , γ = 0.33 .
Figure 27. Graphical representation of the stochastic fractional SIR model (noisy lines) and deterministic fractional SIR model (dashed lines) for fractional order ρ = 0.9 with respect to time when b = 0.01675 , d = 0.007473 , β = 1.05 , σ = 0.01 , γ = 0.33 .
Axioms 13 00545 g027
Figure 28. A visual guide to study selection in the review.
Figure 28. A visual guide to study selection in the review.
Axioms 13 00545 g028
Table 1. Description of the model parameters and the parameter values used in the graphical simulation. See Section 8 (for the data Sources see 6. Data items (b)).
Table 1. Description of the model parameters and the parameter values used in the graphical simulation. See Section 8 (for the data Sources see 6. Data items (b)).
ParameterDescriptionEstimated Value (per Day)
bthe natural birth rate0.01675
dthe natural death rate0.007473
α the rate at which susceptible is exposed1.05
β the transmission rate of infection1.05
β 1 the rate of disease transmission between vaccinated and infected individuals0.725
ω the rate at which infected is susceptible0.06
γ the rate at which infected recovered0.33
γ 1 the rate at which a vaccinated individual obtains immunity0.66
δ the rate of loss of immunity0.026
ν the rate of vaccination0.08
ϵ the rate at which exposed is infected0.3
μ death rate due to infection0.034
Ω the rate at which infected is quarantined0.1
μ death rate due to quarantine0.0006
φ the rate at which quarantine is recovered0.2
θ rate of loss of temporal immunity0.05
υ coincidence rate0.002
σ the intensity of the white noise0.01
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Agarwal, R.; Airan, P.; Agarwal, R.P. Exploring the Landscape of Fractional-Order Models in Epidemiology: A Comparative Simulation Study. Axioms 2024, 13, 545. https://doi.org/10.3390/axioms13080545

AMA Style

Agarwal R, Airan P, Agarwal RP. Exploring the Landscape of Fractional-Order Models in Epidemiology: A Comparative Simulation Study. Axioms. 2024; 13(8):545. https://doi.org/10.3390/axioms13080545

Chicago/Turabian Style

Agarwal, Ritu, Pooja Airan, and Ravi P. Agarwal. 2024. "Exploring the Landscape of Fractional-Order Models in Epidemiology: A Comparative Simulation Study" Axioms 13, no. 8: 545. https://doi.org/10.3390/axioms13080545

APA Style

Agarwal, R., Airan, P., & Agarwal, R. P. (2024). Exploring the Landscape of Fractional-Order Models in Epidemiology: A Comparative Simulation Study. Axioms, 13(8), 545. https://doi.org/10.3390/axioms13080545

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop