
==== Front
Sci Rep
Sci Rep
Scientific Reports
2045-2322
Nature Publishing Group UK London

38589475
58192
10.1038/s41598-024-58192-7
Article
Fractional epidemic model of coronavirus disease with vaccination and crowding effects
Saleem Suhail 1
Rafiq Muhammad 27
Ahmed Nauman 37
Arif Muhammad Shoaib 1
Raza Ali 48
Iqbal Zafar 3
Niazai Shafiullah shafiullahniazai@lu.edu.af

5
Khan Ilyas i.said@mu.edu.sa

69
1 https://ror.org/03yfe9v83 grid.444783.8 0000 0004 0607 2515 Department of Mathematics, Air University, PAF Complex E-9, Islamabad, 44000 Pakistan
2 https://ror.org/04g0mqe67 grid.444936.8 0000 0004 0608 9608 Department of Mathematics, Faculty of Science and Technology, University of Central Punjab, Lahore, Pakistan
3 https://ror.org/051jrjw38 grid.440564.7 0000 0001 0415 4232 Department of Mathematics and Statistics, The University of Lahore, Lahore, Pakistan
4 https://ror.org/01b009v28 Department of Mathematics, University of Chanab, Gujrat, Pakistan
5 Department of Mathematics, Education Faculty, Laghman University, Mehtarlam City, 2701 Laghman Afghanistan
6 https://ror.org/01mcrnj60 grid.449051.d 0000 0004 0441 5633 Department of Mathematics, College of Science Al-Zulfi Majmaah University, 11952 Al-Majmaah, Saudi Arabia
7 https://ror.org/00hqkan37 grid.411323.6 0000 0001 2324 5973 Department of Computer Science and Mathematics, Lebanese American University, Beirut, 1102-2801 Lebanon
8 Department of Mathematics, Mathematics Research Center, Near East University, Near East Boulevard, 99138 Nicosia/Mersin 10, Turkey
9 grid.412431.1 0000 0004 0444 045X Department of Mathematics, Saveetha School of Engineering, SIMATS, Chennai, Tamil Nadu India
8 4 2024
8 4 2024
2024
14 815713 10 2023
26 3 2024
© The Author(s) 2024
2024
https://creativecommons.org/licenses/by/4.0/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Most of the countries in the world are affected by the coronavirus epidemic that put people in danger, with many infected cases and deaths. The crowding factor plays a significant role in the transmission of coronavirus disease. On the other hand, the vaccines of the covid-19 played a decisive role in the control of coronavirus infection. In this paper, a fractional order epidemic model (SIVR) of coronavirus disease is proposed by considering the effects of crowding and vaccination because the transmission of this infection is highly influenced by these two factors. The nonlinear incidence rate with the inclusion of these effects is a better approach to understand and analyse the dynamics of the model. The positivity and boundedness of the fractional order model is ensured by applying some standard results of Mittag Leffler function and Laplace transformation. The equilibrium points are described analytically. The existence and uniqueness of the non-integer order model is also confirmed by using results of the fixed-point theory. Stability analysis is carried out for the system at both the steady states by using Jacobian matrix theory, Routh–Hurwitz criterion and Volterra-type Lyapunov functions. Basic reproductive number is calculated by using next generation matrix. It is verified that disease-free equilibrium is locally asymptotically stable if R0<1 and endemic equilibrium is locally asymptotically stable if R0>1. Moreover, the disease-free equilibrium is globally asymptotically stable if R0<1 and endemic equilibrium is globally asymptotically stable if R0>1. The non-standard finite difference (NSFD) scheme is developed to approximate the solutions of the system. The simulated graphs are presented to show the key features of the NSFD approach. It is proved that non-standard finite difference approach preserves the positivity and boundedness properties of model. The simulated graphs show that the implementation of control strategies reduced the infected population and increase the recovered population. The impact of fractional order parameter α is described by the graphical templates. The future trends of the virus transmission are predicted under some control measures. The current work will be a value addition in the literature. The article is closed by some useful concluding remarks.

Keywords

SIVR model
Infection dynamics
Precautionary measures
Non-linear differential equation
Existence and uniqueness
Basic reproductive number
Stability
Non-standard finite difference approach
Subject terms

Biological models
Software
Mathematics and computing
issue-copyright-statement© Springer Nature Limited 2024
==== Body
pmcIntroduction

The severe acute respiratory syndrome (SARS) was initially detected in Asia in 2003. Then, this infection spread to the other parts of the world. Main symptoms of the SARS virus contain high body temperature, dysentery, discomfort, dry cough and respiratory problems etc.1. It is an air-borne disease, that spreads by the tiny droplets produced by coughing, sneezing, laughing or talking infected person. These symptoms vary from mild to severe depending upon the immunity, health condition and phase of the virus. The SARS cov2 or Covid-19 was started in December 2019 from Wuhan, China. This virus was similar to that of the SARS virus. It is considered that Covid-19 was spread in humans via the bats. The initial symptoms of this virus was similar to that of the pneumonia virus2,3. This virus propagated to the other parts of the world, eventually it reached in the United States on January 20, 2020. It was declared as pandemic on March 11, 2020 by the WHO. In the beginning, the value of R0 was 2.5, apporoximately. Because of its crown like structure, it is called as Corona virus. The spike-protein attacks on the cells of the respiratory system. Gradually, Covid-19 introduces its contaminated RNA into the human bodies. Then the virus replicates and multiplies its production in a short span of time. Everyone is susceptible to the Corona virus, but the aged and sick people are the easy target of this disease. The people with low immunity and respiratory issues, may acquire it fastly. This disease may reoccur in an individual which increases its rate of fatality. The USA acquired the highest rate of infection in the world, while Brazil was recorded with the highest mortality rate and Vietnam with the lowest in the world. Covid-19 is a different virus than influenza-A or Influenza-B viruses. The strain of it is similar to SARS Cov2 and its signs and symptoms appear with in 2 to 14 days after the vulnerability. It has been observed that, COVID-19 is more infectious than the other viruses of the family. One may lose sense of taste or smell after getting the virus. Some other issues like fever and respiratory problems may also develop in the patients. So, the importance of preventive measures for these types of diseases has been increased. Social distancing, wearing face masks and hand sanitizing can control the disease dynamics, considerably4. Vaccination imparts an important role in controlling the disease. But it is not sure that a vaccinated person cannot receive the infection. It is investigated that the people with strong immunity can restrict the propagation of the disease5. Ahmed et al. considered the SEQIIR model in 2021 to study the dynamics of the corona virus. They used ODEs and FDEs models for the deeper insight of the disease phenomena6. Likewise, In 2021, Hassan et al. suggested the SIIR compartmental model to examine the waves of the disease in Texas, USA7. Alqarni et al. described the DSIARB epidemic model to notice the dynamics of Covid-19 in the Kingdom of Saudi Arabia8. Similarly, in Brazil Savi and his co-authors studied an SEIRDC epidemic model to check the propagation of the virus9. Tiwari et al. explored the effect of quarantine on the dynamics of the disease in India. They tested the disease dispersion by analyzing the SEIRD model10. Warbhe et al. in 2021, inspected a SIRM model to study the economic losses due to Covid-1911. In 2021, Daniel studied the compartmental (SEIQCRW) model to explore the dynamics of the infection Covid-19 with diffusion process12. Prathumwan et al. considered a SLIQHR compartmental model to examine the hazards of the corona virus by adopting the safety measures13. The Covid-19 influenced the world economy badly. Balike explored the economic effects of the virus in the model in Congo by studying the SEIHQR model14. Likewise, in Egypt, Raslan discussed SEHQIR compartmental model to explore the virus propagation of the new invariant of the COVID-1915. Chen et al. proposed a new compartment model (BHRP) to simulate the disease dynamics16. In Indonesia, Sinaga et al. contemplated a model (SEIR) for examining the structural properties of the coronavirus17. James et al. highlighted the pros and cons of the mathematical modeling in designing the health policies18. Ameen et al. proposed a fractional order mathematical model (SSLLIPD) for investigating the communication of the virus19. Similarly, different compartmental models were suggested to verify the propagation of Covid-19 in the society20,21. Kahn et al. studied the model disease dynamics by including the vaccination strategy in the model22. Nave and co-authors designed a model to study the stability of the virus by adopting the fast-slow classification strategy23. Kim et al. explored that social distancing, isolation and early case detection plays a vital key role in controlling the disease dynamics of Covid-1924. Das studied the stability analysis of fractional order Covid-19 model25. Machado et al. studied the rare and extreme events during the pandemic of the Covid-1926. On the same lines, many researchers proposed different types of models with the consideration of various important parameters to examine the spread of the virus27. Quaranta and co-authors analyzed various multi-scale territorial models for the deeper understanding of the virus propagation28. Many research studies conducted to investigate the propagation of corona-like diseases29–33. Different epidemic models and strategies are designed to illustrate the transmission dynamics of the Covid-1934–37.

Mathematical models don't offer a definitive cure for infectious diseases, they only serve to simulate various scenarios and dynamics, resilience assessment, control and effective detection strategies. Timely and appropriate interventions are crucial for mitigating the social impact of diseases by controlling their spread. Many mathematical models in literature have been proposed to understand and predict the spread of infectious diseases, with the objective to flatten the infection curve and reduce mortality rates.

Fractional calculus, a field extending classical derivatives and integrals to fractional orders, has a rich history dating back to Leibniz's inquiries in 1695. This branch of applied mathematics deals with real-world phenomena modeled by non-integer-order derivatives and increasing attention globally. Various fractional derivatives, with or without singular kernels, have been developed and extensively employed to model diverse real-life problems. These derivatives include the Caputo, Riemann–Liouville, Katugampola, Caputo–Fabrizio, and Atangana–Baleanu derivatives. The applications of fractional calculus may be seen in50–54.

In38, Ali and co-authors presented SIVR (Susceptible, Infected, Vaccinated, Recovered) epidemic ODE model for analyzing and predicting the COVID-19 but as we know that the fractional order derivative is much more accurate, factual and empirical when compared with the integer order cases.

Researchers have utilized fractional calculus to enhance epidemic models, incorporating memory effects and studied the disease dynamics more accurately. The efficacy of fractional operators has led to their widespread application across various disciplines such as Finance, Engineering, Biology, and Medicine. In the context of COVID-19, where uncertainties abound, fractional calculus offers advantages over integer-order derivatives. Fractional operators, characterized by their non-local nature, are better suited to capturing complex and unpredictable systems. The memory and hereditary properties inherent in fractional order models enable more realistic representations of infectious disease dynamics, incorporating past information for better predictions. Recent developments in fractional calculus may be observed in55–59.

To fill these gaps and to address the challenges arise in the COVID-19 modeling, the present study focuses on analyzing an SIVR epidemic model using the Caputo fractional order operator. This choice is justified because the Caputo derivative has the ability to accommodate local initial conditions and it is compatible with biological and physical principles.

The primary objective of this research is to assess the impact of vaccination and crowding effect on COVID-19 dynamics in the perspective of fractional calculus. The insights gained from this study could assist in strategic planning by governments and public health authorities to bridge immunization gaps and prevent future outbreaks. Furthermore, this fractional order model will contribute to the ongoing research in COVID-19 mathematical modeling and will also assist the researchers to stimulate interest in fractional calculus modeling and mathematical epidemiology.

The design of our paper is as follows. In “Description of the model” section, we have discussed the mathematical model and performed its analysis. Then, in the sub-sections, positivity, boundedness, existence, and uniqueness are studied. In “Qualitative analysis of the proposed model” section, qualitative analysis of the proposed model (local and global stability of the model) is presented. In “Numerical simulations” section, numerical simulations are presented to analyze the dynamics of the virus, graphically. In “Numerical method” section, numerical method (Grunwald–Letnikov non-standard finite difference method) of the proposed model is presented. Then, in the sub-sections, positivity and boundedness of the numerical method are discussed. Finally, the concluding remarks are presented in the closing section. Now, we quote some basic definitions of fractional calculus.

Definition 1

(Caputo Fractional Derivative) The Caputo derivative of fractional order α of function f(t) is defined as.acDtαft=1Γ(n-α)∫atfn(x)t-xα+1-ndx,wheren-1<α<n∈N.

Laplace transform of Caputo fractional derivative

The Laplace transformation of Caputo fractional differential operator of order α is given by:L0cDtαft=sαFs-∑k=0n-1sk-n-1fk0,wheren-1<α<n∈N.

Definition 2

(The Mittag–Leffler Functions) Two parametric Mittag–Leffler function is represented by the series.Eα,βz=∑k=0∞zkΓ(αk+β),α,β>0,α,β∈R,z∈C.

Laplace transformation of the Mittag–Leffler Functions

Laplace transformation of the Mittag–Leffler function is given by.Ltβ-1Eα,βMtα=sα-βsα-M.

For further definitions and properties of fractional calculus, see39–45.

Description of the model

In this segment, we present the fractional epidemic system of coronavirus disease.

In the proposed model, the whole human population N(t) is classified into four subclasses as S(t) (Susceptible class), I(t) (Infected class), V(t) (Vaccinated class), and R(t) (Recovered or immune class). The principle of mass action is taken into account for the infection dynamics in the society. The transmission map of the model is shown in Fig. 1.Figure 1 Flow map of coronavirus model.

Parameters of the model are described as follows ΛN (the recruitment rate of the population),βαI (the force of infection of virus), 11+α1αI (the crowding effect of population on the virus), μα (the rate of mortality due to virus or natural of each subpopulation), δ1α (the rate at which infected population got vaccination during the period of quarantine, or isolation, etc.), δ2α (the rate at which susceptible population got vaccination under the program launched by World Health Organization (WHO)), γα (the rate at which infected population may recover due to its internal immunity and natural circumstances), and σα (the rate of doses in the population who recovered or got immune after vaccination). Deterministic model is based on the following assumptions. Susceptible population becomes immune against the disease after vaccination. Recovered population do not become infected. Only direct contact of infected individuals and susceptible population is considered in this model, all other types of interactions are neglected.

Model equations

The system of equations obtained from the transmission map of the virus is as follows:1 0cDtαSt=Λ-βαSI1+α1αI-δ2α+μαSfort≥0,0cDtαIt=βαSI1+α1αI-γα+δ1α+μαIfort≥0,0cDtαVt=δ2αS+δ1αI-σα+μαVfort≥0,0cDtαRt=γαI+σαV-μαRfort≥0.

Invariant region

The total dynamics of the system (1) is obtained by adding the four equations as follows0cDtαSt+0cDtαIt+0cDtαVt+0cDtαRt=Λ-μαN.

where Nt=St+It+Vt+R(t).

Finally, we have0cDtαNt=Λ-μαN.

Hence, Nt≤M, whenever t→∞.

The feasible region for the system (1) is defined byΩ=St,It,Vt,Rt∈R+4:Nt≤M.

Properties

In this section, we shall study the basic properties of our proposed model (1). The model will be biological meaningful if all the variables are non-negative for t≥0. It means that solution with non-negative initial conditions will remain non-negative for all time. We ensure this result by Theorem 1.

Theorem 1

(Positivity) For any initial data S0,I0,V0,R0∈R4+, the solution St,It,Vt,Rt for the system (1) is positive invariant in R4+.

Proof

We define the norm as f∞=supt∈Dfft.

Firstly, consider the class St,0cDtαSt=Λ-βαSI1+α1αI-δ2α+μαS,

0cDtαSt≥-δ2α+μα+βαI1+α1αIS,∀t≥0,

0cDtαSt≥-δ2α+μα+βαsupt∈DλI1+α1αsupt∈DλIS,∀t≥0,

0cDtαSt≥-δ2α+μα+βαI∞1+α1αI∞S,∀t≥0.

Let M1=δ2α+μα+βαI∞1+α1αI∞,

then0cDtαSt≥-M1St,

0cDtαSt+M1St≥0.

By taking Laplace transformation on both sides of the above expression, we getL{0cDtαSt}+M1LSt≥0,

sαLSt-sα-1S0+M1LSt≥0,

sα+M1LSt≥sα-1S0,LSt≥S0sα-1sα+M1.

By taking inverse Laplace transformation on both sides and by using the following result L-1sα-βsα-M=tβ-1Eα,βMtα, we reach at the expression given belowSt≥S0L-1sα-1sα+M1,

St≥S0Eα,1-M1tα,

St≥S0Eα-M1tα,

St≥0,∀t≥0.

Next, for the function I(t), we proceed as explained below0cDtαIt=βαSI1+α1αI-γα+δ1α+μαI,

0cDtαIt≥-γα+δ1α+μαI.

Let M2=γα+δ1α+μα, then0cDtαIt≥-M2I,

0cDtαIt+M2It≥0.

By taking Laplace transformation on both sides, we derive the following resultL{0cDtαIt}+M2L{I(t)}≥0,

sαLIt-sα-1I0+M2LIt≥0,

(sα+M2)LIt≥sα-1I0,

LIt≥I(0)sα-1sα+M2.

By taking inverse Laplace transformation on both sides and by using the following resultL-1sα-βsα-M=tβ-1Eα,βMtα,

I(t)≥I(0)L-1sα-1sα+M2,

It≥I0Eα,1-M2tα,

It≥I0Eα-M2tα,

It≥0,∀t≥0.

Similarly, for the class V(t) we have0cDtαVt=δ2αS+δ1αI-(σα+μα)V,

0cDtαVt≥-(σα+μα)V,

Let M3=σα+μα.Then, we conclude that0cDtαVt≥-M3V(t).

By solving above inequality, we getVt≥V0Eα,1-M3tα,

Vt≥V0Eα-M3tα,

Vt≥0,∀t≥0.

Finally, for the function Rt0cDtαRt=γαI+σαV-μαR,

0cDtαRt≥-μαR,

By solving above inequality, we getRt≥R0Eα,1-μαtα,

Rt≥R0Eα-μαtα,

Rt≥0,∀t≥0.

As desired.

Now, we present the second property of model (1) i.e. boundedness of the solution.

Theorem 2

(Boundedness) For any time t, the system (1) is bounded and lies in the feasible region Ω.

Proof

By letting the population function N(t)=S(t)+I(t)+V(t)+R(t), we advance as.0cDtαNt=0cDtαSt+0cDtαIt+0cDtαVt+0cDtαRt,

0cDtαNt=Λ-μαN,

0cDtαNt+μαN(t)=Λ.

By taking Laplace transform on both sides of the above expressionL0cDtαNt+μαLN(t)=ΛL1,

sαLNt-sα-1N0+μαLNt=Λs,sα+μαLN(t)=sα-1N0+Λs,

LN(t)=N(0)sα-1sα+μα+Λssα+μα,

LN(t)=N0sα-1sα+μα+Λs-1sα+μα,

LN(t)=N0sα-1sα+μα+Λsα-(1+α)sα+μα.

By taking inverse Laplace transformation on both sides and by using the following formula L-1sα-βsα-M=tβ-1Eα,βMtα, we getNt=N0L-1sα-1sα+μα+ΛL-1sα-(1+α)sα+μα,

Nt=N0Eα,1-μαtα+ΛtαEα,α+1-μαtα,

Nt=N0Eα,1-μαtα+μα.ΛμαtαEα,α+1-μαtα,

Let M=maxN0,Λμα,Nt≤MEα,1-μαtα+μαtαEα,α+1-μαtα.

By using, Eα,βz=zEα,α+βz+1Γ(β), we can write above inequality in the formNt≤M-μαtαEα,α+1-μαtα+1Γ(1)+μαtαEα,α+1-μαtα,

N(t)≤M.

Thus N(t) is bounded uniformly and hence the solution Xt=S,I,V,R of the system (1) is uniformly bounded in Ω=S,I,V,R:S+I+V+R=N≤M for t∈[0,∞).

Therefore, the model (1) is well-posed, biologically and mathematically in the invariant set Ω, as desired.

Existence and uniqueness

This subsection presents the existence and uniqueness of solution of the proposed model using the technique of fixed-point theory, for this we will apply the following Lemma.

Lemma 1

(39,46).Consider the system t0CDtaxt=gt,x,t0>0, with initial condition xt0=xt0, where, α∈0,1,g:t0,∞×Ω→R,Ω⊆C1[t0,∞), if local Lipschitz condition is satisfied by g(t, x) with respect to x, then there exists a unique solution on [t0,∞)×Ω.

Theorem 3

(Existence and Uniqueness) For any time t, the solution of system (1) will exist and the solution will be unique47.

Proof

To study the existence and uniqueness of system (1), let us consider the region Σ×[t0,γ], where.Σ=S,I,V,R∈R4,S,I,V,R∈C1t0,∞∧S,I,V,R≤Mandγ<+∞.

Let KS=Λ-βαSI1+α1αI-δ2α+μαS,KS1-KS2=Λ-βαS1I1+α1αI-δ2α+μαS1-Λ+βαS2I1+α1αI+δ2α+μαS2,=βαI1+α1αIS2-S1+δ2α+μαS2-S1,≤βα11+α1αIIS2-S1+δ2α+μαS2-S1,≤βαMS2-S1+δ2α+μαS2-S1,whereI1+α1αI≤1,≤βαM+δ2α+μαS2-S1.

Therefore, ‖KS1-KS2‖≤βαM+δ2α+μα‖S1-S2‖.

Therefore, K(S) satisfies Lipchitz condition.

For contraction mapping. βαM+δ2α+μα<1,

Let LI=βαSI1+α1αI-γα+δ1α+μαI,LI1-LI2=βαSI11+α1αI1-γα+δ1α+μαI1-βαSI21+α1αI2+γα+δ1α+μαI2,=βαSI11+α1αI1-I21+α1αI2+γα+δ1α+μαI2-I1,≤βαSI1+α1αI1I2-I2-α1αI1I21+α1αI11+α1αI2+γα+δ1α+μαI2-I1,≤βαM·11+α1αI1·11+α1αI2·I1-I2+γα+δ1α+μαI2-I1,≤βαMI1-I2+γα+δ1α+μαI2-I1,whereI1+α1αI≤1,I¯1+α1αI¯≤1,≤βαM+γα+δ1α+μαI1-I2.

Therefore, ‖LI1-LI2‖≤βαM+γα+δ1α+μα‖I1-I2‖.

Consequently, L(I) satisfies Lipchitz condition.

For contraction mapping, βαM+γα+δ1α+μα<1.

Let NV=δ2αS+δ1αI-σα+μαV,NV1-NV2=δ2αS+δ1αI-σα+μαV1-δ2αS-δ1αI+σα+μαV2,=σα+μαV2-V1,=σα+μαV2-V1,≤σα+μαV1-V2,≤σα+μαV1-V2.

Therefore, ‖NV1-NV2‖<σα+μα‖V1-V2‖.

So, N(V) satisfies Lipchitz condition.

For contraction mapping, σα+μα<1.

Let PR=γαI+σαV-μαR.PR1-PR2=γαI+σαV-μαR1-γαI-σαV+μαR2,=μαR2-R1,=μαR1-R2.

‖PR1-PR2‖<‖R1-R2‖. Where μα<1.

Therefore, P(R) satisfies Lipchitz condition.

Let, F1=βαM+δ2α+μα, F2=βαM+γα+δ1α+μα, F3=σα+μα , F4=μα.

Also let F=maxF1,F2,F3,F4.

Therefore,‖KS1-KS2‖≤F‖S1-S2‖,

‖LI1-LI2‖≤F‖I1-I2‖,

‖NV1-NV2‖≤F‖V1-V2‖,

‖PR1-PR2‖≤F‖R1-R2‖.

For F<1, K(S),L(I),N(V),P(R) are contraction mappings.

Therefore, K(S),L(I),N(V) and P(R) satisfies Lipshitz conditions and are contraction mappings. Therefore, according to Banach fixed point theorem, solution of the proposed model (1) exists and is unique, as desired.

Qualitative analysis of the proposed model

This section presents the equilibrium points, basic reproductive number and stability analysis, analytically. The order of the results for the stability is described as:Asymptotically local stability at corona virus-free equilibrium is established by Theorem 4.

Asymptotically local stability at corona existing equilibrium is ensured by Theorem 5.

Asymptotically global stability at corona virus-free equilibrium is guaranteed by Theorem 6.

Asymptotically global stability at corona existing equilibrium is confirmed by Theorem 7.

Equilibria

We determine the equilibria of the system (1) by assuming the state variables are constant and by putting the right side equal to zero. Equation (1) admits two types of equilibria as follows:(i) corona virus-free equilibrium = C0=S0,I0,V0,R0=Λδ2α+μα,0,0,0,

(ii) corona existing equilibrium = C1=S∗,I∗,V∗,R∗,whereS∗=γα+δ1α+μα1+ααI∗βα,I∗=βαΛ-δ2α+μαδ1α+μαβαδ1α+μα+ααδ2α+μαδ1α+μα,V∗=δ2αS∗+δ1αI∗σα+μα,R∗=γαI∗+σαV∗μα.

Reproduction number

We determine the reproduction number of the system (1) by using the well-known results like the next-generation matrix method after substituting the value of coronavirus free equilibrium. We get:I′V′R′=βαΛδ2α+μα00000000IVR-γα+δ1α+μα00-δ1α(σα+μα)00-γα-σαμαIVR

where A=βαΛδ2α+μα00000000,B=γα+δ1α+μα00-δ1α(σα+μα)00-γα-σαμα,AB-1=βαΛδ2α+μαγα+δ1α+μα00000000,

R0 is the spectral radius of AB-1.

Mathematically R0 is expressed by the following parametric relationR0=βαΛδ2α+μαγα+δ1α+μα.

Stability analysis

In this section, we test the local and global stability of the system (1), considering the two defined equilibria.

Theorem 4

(Local stability at C0) The system (1) at C0=S0,I0,V0,R0=Λδ2α+μα,0,0,0 is locally asymptotically stable if R0<1. Otherwise, unstable when R0>1.

Proof

The Jacobian matrix is obtained for the system (1) as follows.JS,I,V,R=-βαI1+α1αI-δ2α-μα-βαS1+α1αI200βαI1+α1αIβαS1+α1αI2-γα-δ1α-μα00δ2αδ1α-σα-μα00γασα-μα.

The Jacobian matrix at C0 is computed as followsJΛδ2α+μα,0,0,0=-δ2α-μα-βαΛδ2α+μα00βαI1+α1αIβαΛδ2α+μα-γα-δ1α-μα00δ2αδ1α-σα-μα00γασα-μα.

Consider J-λI=0, and hence-δ2α-μα-λ-βαΛδ2α+μα00βαI1+α1αIIβαΛδ2α+μα-γα-δ1α-μα-λ00δ2αδ1α-σα-μα-λ00γασα-μα-λ=0.

After solving above determinant, we get λ1=-μα<0,λ2=-σα+μα<0,λ3=-δ2α+μα<0λ4=βαΛδ2α+μα-γα-δ1α-μα<0.

Consequently, we draw the following resultβαΛδ2α+μαγα+δ1α+μα<1,R0<1.

It is clear that, system (1) is locally asymptotically stable at C0 when R0<1, as desired.

Theorem 5

(Local stability at C1) The system (1) at C1=S∗,I∗,V∗,R∗ is locally asymptotically stable if R0>1.

Proof

The Jacobian matrix at C1 of the system (1) is as follows.JS∗,I∗,V∗,R∗=-βαI∗1+α1αI∗-δ2α-μα-βαS∗1+α1αI∗200βαI∗1+α1αI∗βαS∗1+α1αI∗2-γα-δ1α-μα00δ2αδ1α-σα-μα00γασα-μα.

Consider J-λI=0, and hence-βαI∗1+α1αI∗-δ2α-μα-λ-βαS∗1+α1αI∗200βαI∗1+α1αI∗βαS∗1+α1αI∗2-γα-δ1α-μα-λ00δ2αδ1α-σα-μα-λ00γασα-μα-λ=0.

After solving above determinant, λ1=-μα<0,λ2=-σα+μα<0,λ2+βαI∗1+α1αI∗-βαS∗1+α1αI∗2+γα+δ1α+δ2α+2μαλ+γα+δ1α+μαβαI∗1+α1αI∗+δ2α+μα-δ2α+μαβαS∗1+α1αI∗2=0,

λ2+A1λ+A0=0.

where A1=βαI∗1+α1αI∗-βαS∗1+α1αI∗2+γα+δ1α+δ2α+2μα,A0=γα+δ1α+μαβαI∗1+α1αI∗+δ2α+μα-δ2α+μαβαS∗1+α1αI∗2.

Since, A1,A0 both are positive when R0>1, by the Routh-Hurwitz criterion for the second order polynomial, the system (1) is locally asymptotically stable at C1 when R0>1, as desired.

The following lemma is provided to improve the global stability analysis of the system.

Lemma 2

(25) Let x:[0,∞)→R+ be a continuous function and let t0≥0. Then, for any time t≥t0,α∈0,1andx∗∈R+, the following inequality holds.0cDtαxt-x∗-x∗lnx(t)x∗≤1-x∗(t)x0cDtαx(t).

We tackle now the global asymptotic stability of the system at the equilibrium points46.

Theorem 6

(Global stability at C0) The system (1) at C0=S0,I0,V0,R0=Λδ2α+μα,0,0,0 is globally asymptotically stable if R0<1.

Proof

Firstly, we define the Lyapunov function as.L=S+(I+V+R)-S0-S0logSS0,

L=S-S0-S0logSS0+(I+V+R).

Apply Caputo fractional derivative on both the sides,0cDtαL=0cDtαS-S0-S0logSS0+0cDtαI+0cDtαV+0cDtαR.

By using Lemma 1, following result is obtained0cDtαL≤1-S0S0cDtαS+0cDtαI+0cDtαV+0cDtαR,

0cDtαL≤S-S0SΛ-βαI1+α1αI+δ2α+μαS+βαSI1+α1αI-γα+δ1α+μαI+δ2αS+δ1αI-σα+μαV+γαI+σαV-μαR,

0cDtαL≤S-S0SΛ-βαI1+α1αI+δ2α+μαS+γα+δ1α+μαβαS1+α1αIγα+δ1α+μα-1I

-σα+μαV-δ2αS+δ1αIσα+μα-μαR-γαI+σαVμα,

0cDtαL≤S-S0SΛ-ΛSS0+γα+δ1α+μαβαΛα2α+μαγα+δ1α+μα-1I

-σα+μαV-δ2αS+δ1αIσα+μα-μαR-γαI+σαVμα,

0cDtαL≤-ΛS-S02SS0+γα+δ1α+μαR0-1I-σα+μαV-δ2αS+δ1αIσα+μα

-μαR-γαI+σαVμα.

Observe that 0cDtαL<0 when R0<1.

Hence, the system is globally asymptotically stable at the disease-free equilibrium C0.

Theorem 7

(Global stability at C1) The system (1) at C1=S∗,I∗,V∗,R∗ is globally asymptotically stable if R0>1.

Proof

The proof of this theorem is similar to the previous theorem. To prove this theorem, we construct the Lyapunov functional at C1 as:L=S-S∗-S∗logSS∗+I-I∗-I∗logII∗+V-V∗-V∗logVV∗+R-R∗-R∗logSR∗,

0cDtαL=0cDtαS-S∗-S∗logSS∗+0cDtαI-I∗-I∗logII∗+0cDtαV-V∗-V∗logVV∗+0cDtαR-R∗-R∗logSR∗.

By using Lemma 1:0cDtαL≤S-S∗S0cDtαS+I-I∗I0cDtαI+V-V∗V0cDtαV+R-R∗R0cDtαR.

After some simplifications, we get:0cDtαL≤-ΛS-S∗2SS∗-α1αβαSI-I∗21+α1αI(1+α1αI∗)-δ2αS+δ1αIV-V∗2VV∗-γαI+σαVR-R∗2RR∗.

Observe that 0cDtαL≤0 when R0>1. Moreover, 0cDtαL=0 if S=S∗,I=I∗,V=V∗,R=R∗.

Therefore, the system is globally asymptotically stable at the endemic equilibrium C1 , as desired.

Numerical simulations

In this section, parametric values for the simulations are described.

This section is devoted for investigating the key properties of the simulated graphs against the set of parametric values as mentioned in Table 1. Further, these graphs are plotted when the disease is prevailing in the human population and attains the endemic steady state with in due course of time. To analyze the dynamics of the susceptible individuals at endemic equilibrium, four appropriate values of α are selected. We have practically noticed in the outburst of COVID-19 that it propagated with different rates in the different countries of the world. These rates are biologically meaning full as every individual or group of individuals have different physical environments, immunity levels, health conditions and leaving standards etc. In this regard, fractional parameter α handles the rate of spread for the corona virus. Figure 1 shows the flow map of coronavirus model. All the graphs in Fig. 2 depicts the dynamical behavior of the susceptible population. Every trajectory of the curved graph reflects that how the susceptible compartment attains the endemic steady state for different values of α, with a particular set of parametric values including the crowding effect. It can be observed that each graph indicates the convergence towards the exact fixed point. Also, each graph possesses a specific trajectory and rate of convergence according to the fractional order parameter α.Table 1 Set of parametric values per day.

Symbols	Value/per day	Source	
Λ	μ×N(0)	Estimated	
μ	167.7×365	23	
γ	0.5000	Fitted	
β	0.4000	Fitted	
α1	0.5465	Fitted	
δ	0.5000	Fitted	
δ1	0.1000	Fitted	
δ2	5.32978	Fitted	
σ	≥0	Fitted	

Figure 2 Graphical behavior of susceptible population (disease present in the population) for various values of fractional order α.

Similarly, the curved graph in Fig. 3 shows the convergence of infected compartment, towards the exact fixed point. Infected graphs capture the dynamical evolution in the presence of viral load in the infected populace for the highlighted values of parameter. It is worth mentioning that each parameter has its own biological meaning and importance as described in the description of the model in section "Description of the model". All the graphs show the different phases of the disease dynamics according to the value of α. If value of α is less then the rate of convergence is slow and vice versa. Moreover, the graphical behavior is in line with the real behavior of the continues system.Figure 3 Graphical behavior of infected population (disease present in the population) for various values of fractional order α.

The set of the graphs in Fig. 4 demonstrate the vital behavior of vaccinated the populace in Covid-19 model. These graphs describe the evolutionary behavior of the vaccinated individuals during the propagation of Covid-19. Every graph has the specific trajectory and a different rate of convergence to attain the required steady state. Thus, fractional order parameter α plays an important role in describing the disease dynamics of every compartment.Figure 4 Graphical behavior of vaccinated population (disease present in the population) for various values of fractional order α.

Finally, the graphs in Fig. 5 shows the dynamical behavior of the recovered population, when the Covid-19 disease persists in the community. The graphs bring an important fact into the lime light that the fast recovery rate can be captured by a large value of α. Similarly, the slower rate of recovery may be represented by a smaller value of α.Figure 5 Graphical behavior of recovered population (disease present in the population) for various values of fractional order α.

Effect of vaccination

Vaccine of any disease plays a vital role in controlling the disease. Likewise, vaccination is an important strategy to control the Covid-19. The graphs in Fig. 6 depict the role of vaccinating for stopping the spread of the virus. All the graphs are drawn by choosing the α=0.8 and values of vaccination parameter δ1 are shown in the values of vaccination. When the value of δ1 is less i.e. 0.1, the number of infected individuals are greater as compared to the number of individuals against the δ1=0.2,0.3 and 0.4.Figure 6 Effect of vaccination on active cases when fractional order value α=0.8.

Consequently, the number of infected individuals It are inversely proportional to the vaccination parameters δ1. Hence, the parameter δ1 imparts a significant role to slow down the disease dynamics and the transmission of the virus can be controlled effectively by vaccinating a large part of population.

The main advantage of the theory, applied in this work is that it may be applied to examine the positivity, boundedness, stability analysis and uniquely existence of the solution for the fractional order epidemic models. On the other hand, if we look on the drawbacks or disadvantages of this study, one can notice that any epidemic model cannot be captured all the factors related to disease dynamics. So, there is a chance of error and slight deviation from the actual behavior of the disease.

Numerical method

In this section, we construct nonstandard finite difference (NSFD) scheme using Grunwald–Letnikov approximation to numerically approximate Caputo fractional derivative47–49.

Grunwald–Letnikov non-standard finite difference method

Let us now apply the Grunwald–Letnikov definition to the following fractional differential equation using the Caputo operator:0cDtαyt=fy(t),yT=y0(0<α<1),

And assume there exists a unique solution y=yT in the interval [0,T] and let yk denote the approximation of the true solution y(Tk). Then the explicit or implicit Grunwald–Letnikov method for an equidistant grid is given by:

yn+1-∑ν=1n+1Cναyn+1-ν-rn+1αy0=hαf(yn) or hαf(yn+1).

Where Cνα=-1ν-1αν,rn+1α=hαr0αTn+1=γ0,-1αn+1-α,

and the coefficients γ0,-1α=Γ(μα+1)Γ(kα+1),μ,k∈N0∪-1.

Moreover, the coefficients Cνα and rνα satisfy the relations given in the next lemma.

Lemma 3

(48) Assume that 0<α<1, then all the coefficients Cνα are positive and show the behaviour Cνα=O(1ν1+α) as ν→∞. Further, the coefficients Cνα and rνα satisfy for ν>1 the properties:0<Cν+1α<Cνα<⋯<C1α=α<1and0<rν+1α<rνα<⋯<r1α=1Γ(1-α).

Proof

It is clear that for all α∈0,1,limν→∞Cνα=0,limν→∞rνα=0.

Moreover, for α∈0,1, is α+1Γ(1-α)>1.

It also shows that ∑ν=1∞Cνα=1.

Which is the result that will used later.

Now, we provide a non-standard finite-difference discretization of the mathematical model (1).

To this, we divide the interval [0,L] into M ∈ N subintervals, respectively, and step sizes h=LM. The approximate solutions S,I,VandR of (1) will be denoted as Sn,In,Vn,Rn, respectively, for each n=0,1,⋯,N. Under the rules, the discretization of a system (1) is presented in43.

First equation of system (1):0cDtαSt=Λ-βαSI1+α1αI-(δ2α+μα)S.

Use Grunwald–Letnikov approximation 0cDtαytn+1=1ψhαyn+1-∑ν=1n+1Cναyn+1-ν-rn+1αy0 to left hand side of above equation, we have:1ψhαSn+1-∑ν=1n+1CναSn+1-ν-rn+1αS0=Λ-βαSn+1In1+α1αIn-(δ2α+μα)Sn+1,

Sn+1-∑ν=1n+1CναSn+1-ν-rn+1αS0=Λψhα-ψhαβαSn+1In1+α1αIn-ψhα(δ2α+μα)Sn+1.

By shifting Sn+1 terms to L.H.S.Sn+1+ψhαβαSn+1In1+α1αIn+ψhαδ2α+μαSn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν,

Sn+11+ψhαβαIn1+α1αIn+ψhαδ2α+μα=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν,

Sn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν1+ψhαβαIn1+α1αIn+ψhαδ2α+μα.

Second equation of system (1): 0cDtαIt=βαSI1+α1αI-(γα+δ1α+μα)I.

Use Grunwald–Letnikov approximation, 0cDtαytn+1=1ψhαyn+1-∑ν=1n+1Cναyn+1-ν-rn+1αy0 to left hand side of above equation, we have:1ψhαIn+1-∑ν=1n+1CναIn+1-ν-rn+1αI0=βαSn+1In1+α1αIn-(γα+δ1α+μα)In+1,

In+1-∑ν=1n+1CναIn+1-ν-rn+1αI0=ψhαβαSn+1In1+α1αIn-ψhα(γα+δ1α+μα)In+1.

By shifting In+1 terms to L.H.S.In+1+ψhαγα+δ1α+μαIn+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν,

In+11++ψhαγα+δ1α+μα=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν,

In+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν1+ψhαγα+δ1α+μα.

Third equation of system (1): 0cDtαVt=δ2αS+δ1αI-(σα+μα)V.

Use Grunwald–Letnikov approximation 0cDtαytn+1=1ψhαyn+1-∑ν=1n+1Cναyn+1-ν-rn+1αy0 to left hand side of above equation, we have:1ψhαVn+1-∑ν=1n+1CναVn+1-ν-rn+1αV0=δ2αSn+1+δ1αIn+1-σα+μαVn+1,

Vn+1-∑ν=1n+1CναVn+1-ν-rn+1αV0=ψhαδ2αSn+1+ψhαδ1αIn+1-ψhασα+μαVn+1,

Vn+1+ψhασα+μαVn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV0,

Vn+11+ψhασα+μα=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV0,

Vn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV01+ψhασα+μα.

Fourth equation of system (1):0cDtαRt=γαI+σαV-μαR.

Use Grunwald–Letnikov approximation 0cDtαytn+1=1ψhαyn+1-∑ν=1n+1Cναyn+1-ν-rn+1αy0 to left hand side of above equation, we have:1ψhαRn+1-∑ν=1n+1CναRn+1-ν-rn+1αR0=γαIn+1+σαVn+1-μαRn+1,

Rn+1-∑ν=1n+1CναRn+1-ν-rn+1αR0=ψhαγαIn+1+ψhασαVn+1-ψhαμαRn+1,

Rn+1+ψhαμαRn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+1,

Rn+11+ψhαμα=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+1,

Rn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+11+ψhαμα.

Now we have a system of equations:Sn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν1+ψhαβαIn1+α1αIn+ψhαδ2α+μα,

In+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν1+ψhαγα+δ1α+μα,

Vn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV01+ψhασα+μα,

Rn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+11+ψhαμα.

Next, we establish the most important properties of NSFD method48.

Properties

Theorem 8

(Positivity) The deterministic form of SIVR model system preserves the non-negativity of the solution.

Proof

We will prove this by induction method.

Suppose that S0≥0,I0≥0,V0≥0,R0≥0, then Sn≥0,In≥0,Vn≥0,Rn≥0 is satisfied for n=1,2,⋯

By induction, for n=0, above system can be written as:S1=Λψhα+r1αS0+C1αS01+ψhαβαI01+α1αI0+ψhαδ2α+μα≥0,

I1=ψhαβαS1I01+α1αI0+r1αI0+C1αI01+ψhαγα+δ1α+μα≥0,

V1=C1αV0+ψhαδ2αS1+ψhαδ1αI1+r1αV01+ψhασα+μα≥0,

R1=C1αR0+r1αR0+ψhαγαI1+ψhασαV11+ψhαμα≥0.

As all the parameters are positive, so that, S1≥0,I1≥0,V1≥0,R1≥0.

Now, we will suppose that for positive integers 1,2,⋯,n-1,Sn≥0,In≥0,Vn≥0,Rn≥0, i.e. S2,S3,⋯,Sn≥0,I2,I3,⋯.,In≥0,V1,V2,⋯.,Vn≥0,R1,R2,⋯,Rn≥0.Sn=Λψhα+rnαS0+∑ν=1nCναSn-ν1+ψhαβαIn-11+α1αIn-1+ψhαδ2α+μα≥0,∀n∈1,2,3,⋯n-1,

In=ψhαβαSnIn-11+α1αIn-1+rnαI0+∑ν=1nCναIn-ν1+ψhαγα+δ1α+μα≥0,∀n∈1,2,3,⋯n-1,

Vn=∑ν=1nCναVn-ν+ψhαδ2αSn+ψhαδ1αIn+rnαV01+ψhασα+μα≥0,∀n∈1,2,3,⋯n-1,

Rn=∑ν=1nCναRn-ν+rnαR0+ψhαγαIn+ψhασαVn1+ψhαμα≥0,∀n∈1,2,3,⋯n-1.

By using above results, for positive integer n i.e. ∈Z+ , we conclude that:Sn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν1+ψhαβαIn1+α1αIn+ψhαδ2α+μα≥0,

In+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν1+ψhαγα+δ1α+μα≥0,

Vn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV01+ψhασα+μα≥0,

Rn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+11+ψhαμα≥0.

Hence all the state variables in the discretized model guarantee the positivity of the solution, as desired.

Theorem 9

(Boundedness) Suppose that α1α>0,δ1α>0,δ2α>0,μα>0,γα>0,σα>0,ψhα>0 For all α∈(0,1), then there is a constant:NN0+1,α=N(N0,α)+1Γ(1-α)+ψhαΛ1+ψhαμα.

Such that Sn+1,In+1,In+1,Rn+1≤NN0+1,α for n=0,1,2,⋯,N0.

Proof

Sn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν1+ψhαβαIn1+α1αIn+ψhαδ2α+μα, Sn+11+ψhαβαIn1+α1αIn+ψhαδ2α+μα=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν,

Sn+11+ψhαβαIn1+α1αIn+ψhαδ2α+μα=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν,

Sn+1+ψhαβαIn1+α1αInSn+1+ψhαδ2αSn+1++ψhαμαSn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν.

In+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν1+ψhαγα+δ1α+μα,

In+11+ψhαγα+δ1α+μα=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν,

In+1+ψhαγαIn+1+ψhαδ1αIn+1+ψhαμαIn+1=ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν,=Vn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV01+ψhασα+μα,

Vn+11+ψhασα+μα=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV0,

Vn+1+ψhασαVn+1+ψhαμαVn+1=∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV0.

Rn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+11+ψhαμα,

Rn+11+ψhαμα=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+1,

Rn+1+ψhαμαRn+1=∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+1.

By adding the above equations, we have:Sn+1+ψhαβαIn1+α1αInSn+1+ψhαδ2αSn+1++ψhαμαSn+1+In+1+ψhαγαIn+1+ψhαδ1αIn+1+ψhαμαIn+1Vn+1+ψhασαVn+1+ψhαμαVn+1+Rn+1+ψhαμαRn+1=Λψhα+rn+1αS0+∑ν=1n+1CναSn+1-ν+ψhαβαSn+1In1+α1αIn+rn+1αI0+∑ν=1n+1CναIn+1-ν+∑ν=1n+1CναVn+1-ν+ψhαδ2αSn+1+ψhαδ1αIn+1+rn+1αV0+∑ν=1n+1CναRn+1-ν+rn+1αR0+ψhαγαIn+1+ψhασαVn+1.

Sn+11+ψhαμα+In+11+ψhαμα+Vn+1+Rn+11+ψhαμα=Λψhα+∑ν=1n+1CναSn+1-ν+In+1-ν+Vn+1-ν+Rn+1-ν+rn+1αS0+I0+V0+R0,

1+ψhαμαSn+1+In+1+Vn+1+Rn+1=Λψhα+∑ν=1n+1CναSn+1-ν+In+1-ν+Vn+1-ν+Rn+1-ν+rn+1αS0+I0+V0+R0,

Sn+1+In+1+Vn+1+Rn+1=Λψhα+∑ν=1n+1CναSn+1-ν+In+1-ν+Vn+1-ν+Rn+1-ν+rn+1αS0+I0+V0+R01+ψhαμα,

Sn+1+In+1+Vn+1+Rn+1=Λψhα+∑ν=1n+1CναSn+1-ν+In+1-ν+Vn+1-ν+Rn+1-ν+rn+1α1+ψhαμα,

For, n=0,S1+I1+V1+R1=Λψhα+C1αS0+I0+V0+R0+r1α1+ψhαμα,

S1+I1+V1+R1=Λψhα+C1α+r1α1+ψhαμα,

S1+I1+V1+R1=Λψhα+α+1Γ(1-α)1+ψhαμα,S1+I1+V1+R1=N(1,α),

Now for n=1 and α∈(0,1),S2+I2+V2+R2=Λψhα+∑ν=12CναS2-ν+I2-ν+V2-ν+R2-ν+r2α1+ψhαμα,

S2+I2+V2+R2=Λψhα+C1αS1+I1+V1+R1+C2α+r2α1+ψhαμα,

S2+I2+V2+R2<Λψhα+C1αN(1,α)+C2α+r1α1+ψhαμα, from lemma 3, r2α<r1α,S2+I2+V2+R2<Λψhα+N1,αC1α+C2α+1Γ(1-α)1+ψhαμα,N1,α>1,

S2+I2+V2+R2<Λψhα+N1,α∑ν=1∞Cνα+1Γ(1-α)1+ψhαμα,

S2+I2+V2+R2<Λψhα+N1,α+1Γ(1-α)1+ψhαμα, from lemma 3, ∑ν=1∞Cνα=1,S2+I2+V2+R2<N(2,α).

where N2,α=Λψhα+N1,α+1Γ(1-α)1+ψhαμα.

Similarly, for n=0, we arrive at;S3+I3+V3+R3=Λψhα+C1αS2+I2+V2+R2+C2αS1+I1+V1+R1+C3α+r3α1+ψhαμα,

S3+I3+V3+R3<Λψhα+C1αN2,α+C2αN(1,α)+C3α+r1α1+ψhαμα,

S3+I3+V3+R3<Λψhα+N2,αC1α+C2α+C3α+1Γ(1-α)1+ψhαμα,

S3+I3+V3+R3<Λψhα+N2,α∑ν=1∞Cνα+1Γ(1-α)1+ψhαμα=N(3,α),

S3+I3+V3+R3<N(3,α).

Now, we suppose that for n∈3,4,⋯,N0-1 and for α∈0,1,SN0+IN0+VN0+RN0<NN0-1,α+1Γ(1-α)+Λψhα1+ψhαμα=NN0,α.

Thus, for =N0, we obtain:Sn+1+In+1+Vn+1+Rn+1=Λψhα+∑ν=1n+1CναSn+1-ν+In+1-ν+Vn+1-ν+Rn+1-ν+rn+1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1=Λψhα+∑ν=1N0+1CναSN0+1-ν+IN0+1-ν+VN0+1-ν+RN0+1-ν+rN0+1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1=Λψhα+C1αSN0+IN0+VN0+RN0+C2αSN0-1+IN0-1+VN0-1+RN0-1+…+CN0+1αS0+I0+V0+R0+rN0+1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1<Λψhα+C1αNN0,α+C2αNN0-1,α+C3αNN0-2,α+…+CN0αN1,α+CN0+1α+r1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1<Λψhα+C1αNN0,α+C2αNN0,α+C3αNN0,α+…+CN0αNN0,α+CN0+1αNN0,α+r1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1<Λψhα+N(N0,α)∑ν=1N0+1Cνα+r1α1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1<NN0,α∑ν=1∞Cνα+r1α+Λψhα1+ψhαμα,

SN0+1+IN0+1+VN0+1+RN0+1<NN0,α+1Γ(1-α)+Λψhα1+ψhαμα=NN0+1,α.

Equivalently, it can be described as:

Sn+1+In+1+Vn+1+Rn+1<NN0+1,α ∀n=1,2,3,⋯,N0, where 0<α<1.

Finally, Sn+1,In+1,In+1,Rn+1≤NN0+1,α for =0,1,2,⋯,N0 , as desired.

Remark

(1) NN0+1,α>NN0,α>NN0-1,α>⋯>N1,α>1.

(2) NN0+1,α→1,asα→1-.

Discussion

In our current world, interdisciplinary research stands as a vital tool for ensuring the well-being and standard of living for all humanity. Over time, mathematical prediction modelling has become indispensable for understanding the behavior of epidemics, aiding policymakers in making crucial decisions and preparing for future challenges. Our study aimed to develop a mathematical model to predict the dynamics of COVID-19, employing the SIVR modelling framework. Scientists globally are striving to identify effective control measures to curb the spread of the COVID-19 virus. Strategies such as social distancing, quarantine for exposed individuals, isolation of infected persons, and mask-wearing are widely adopted. The primary objective in containing the virus is to minimize contact between susceptible and infected individuals, thereby reducing the transmission rate from the susceptible (S) to the infected (I) class.

Graphs plotted at the endemic equilibrium point illustrate the dynamical behavior of the susceptible population for varying fractional order parameter α, as depicted in Fig. 2. Each curve's trajectory demonstrates the susceptible compartment's convergence to the endemic steady state for different α values and specific parameter sets. Notably, the graphs display distinctive trajectories and convergence rates corresponding to α values, influencing the disease dynamics. Similarly, Fig. 3 exhibits curved graphs demonstrating the convergence of the infected compartment towards the fixed point, showcasing different phases of disease dynamics based on α values. Lower α values indicate slower convergence rates and vice versa, aligning with real-world system behaviors.

Figure 4 presents a set of graphs illustrating the critical role of vaccination in the COVID-19 model. These graphs depict the evolutionary behavior of vaccinated individuals during the disease propagation, with each graph exhibiting unique trajectories and convergence rates. The fractional order parameter α significantly influences the disease dynamics of each compartment.

Furthermore, Fig. 5 displays the dynamic behavior of the recovered population in persistent COVID-19 scenarios. The graphs highlight that a higher α value captures faster recovery rates, while lower α values represent slower recovery rates. Vaccination plays a pivotal role in disease control, particularly in managing COVID-19. Figure 5 illustrates the impact of vaccination on halting virus spread, with graphs drawn at α = 0.8 and varying vaccination parameter values δ1. Lower δ1 values result in a greater number of infected individuals compared to higher values of δ1 (e.g., δ1 = 0.2, 0.3, and 0.4). Thus, the number of infected individuals (I(t)) is inversely proportional to vaccination parameters δ1, underscoring its significant role in mitigating disease dynamics and effectively controlling virus transmission through widespread vaccination efforts.

Concluding remarks

The fractional order COVID-19 model is developed to study the virus propagation in the population. To this end, analytical and numerical results are established. These results ascertained the unique positive solution of the underlying model. For numerical solutions, the GL-NSFD scheme is formulated to study the numerical aspects of the proposed scheme. This scheme has provided the positive and bounded numerical solutions, which are the salient features of the compartmental models. It is also identified that the disease model has two equilibrium state i.e. virus free and virus existing state. The simulated graphs and rendered that the endemic equilibrium may be obtained physically. Moreover, the disease-free state may also be attained in real scenario. Stability analysis is arranged and results for local and global stability are established for different values of R0 i.e. R0<1 or R0>1. The results have shown that the fractional COVID-19 system preserves the stability at both the steady states depending upon the values of R0. The crowding effect of the corona virus on the disease dynamics is presented. The role of vaccination for controlling the virus spread is also studied graphically. It is observed that the virus may be controlled, significantly if maximum portion of the population is vaccinated. Moreover, the SOP’s planned by the world health organization can reduce the rate of virus propagation in the society. Consequently, our proposed fractional epidemic model for the COVID-19 may be considered for studying the disease dynamics in the community. In addition, crowding effect and vaccination strategy can play a vital role devising the health policies and pre cautionary measures. Moreover, the fractional order epidemic model will be helpful for studying the dynamics of various disease in future.

In reading through the referee reports, we've noted that a number of additional specific references have been suggested. We would like to emphasize that we do not expect a research article to be comprehensive in its referencing of the background literature and that these additions would be at your discretion. Not including references which are not directly relevant to your work would not affect the final decision. However, please do ensure that directly relevant previous work is adequately acknowledged in your paper and that the original contribution of your study to furthering the current understanding in the field is sufficiently emphasized.

Acknowledgements

The authors would like to thank the Deanship of Scientific Research at Majmaah University for supporting this work under Project R-2022-XX.

Author contributions

S.S.: Designed, and analyzed the results; R.F.: prepared figures and discussed; N.A.: software and coding; M.S.A and S.N..: analyzed the data, funding, A.R.: discussed the results; Z.I. software and coding; I.K: wrote the manuscript; supervised the research.

Data availability

The datasets used and analyzed during the current study available from the corresponding author on reasonable request.

Competing interests

The authors declare no competing interests.

Publisher's note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
==== Refs
References

1. https://www.who.int/Csr/Sars/WHOconsensus.Pdf.
2. https://www.who.int/Docs/Defaultsource/Coronaviruse/Who-China-Jointmission-on-Covid-19-Fnal-Report.Pd.
3. https://www.who.int/Emergencies/Diseases/Novel-Coronavirus-2019.
4. Baud D Qi X Nielsen-Saines K Musso D Pomar L Favre G Real estimates of mortality following Covid-19 infection Lancet Infect Dis 2020 20 7 773 10.1016/S1473-3099(20)30195-X 32171390
Baud, D. et al. Real estimates of mortality following Covid-19 infection. Lancet Infect Dis 20(7), 773 (2020).32171390 10.1016/S1473-3099(20)30195-X
5. https://www.who.int/Health-Topics/Coronavirus.
6. Ahmed I Modu GU Yusuf A Kumam P Yusuf I A mathematical model of coronavirus disease (COVID-19) containing asymptomatic and symptomatic classes Results Phys. 2021 10.1016/j.rinp.2020.103776 34367890
Ahmed, I., Modu, G. U., Yusuf, A., Kumam, P. & Yusuf, I. A mathematical model of coronavirus disease (COVID-19) containing asymptomatic and symptomatic classes. Results Phys.10.1016/j.rinp.2020.103776 (2021).34367890 10.1016/j.rinp.2020.103776
7. Hassan MN Mahmud MS Nipa KF Kamrujjaman M Mathematical modeling and Covid-19 forecast in Texas, USA: A prediction model analysis and the probability of disease outbreak Disaster Med. Public Health Prep. 2021 10.1017/dmp.2021.151 34496993
Hassan, M. N., Mahmud, M. S., Nipa, K. F. & Kamrujjaman, M. Mathematical modeling and Covid-19 forecast in Texas, USA: A prediction model analysis and the probability of disease outbreak. Disaster Med. Public Health Prep.10.1017/dmp.2021.151 (2021).34496993 10.1017/dmp.2021.151
8. Alqarni MS Alghamdi M Muhammad T Alshomrani AS Khan MA Mathematical modeling for novel coronavirus (COVID-19) and control Numer. Methods Partial Differ. Equ. 2022 38 760 776 10.1002/num.22695 33362341
Alqarni, M. S., Alghamdi, M., Muhammad, T., Alshomrani, A. S. & Khan, M. A. Mathematical modeling for novel coronavirus (COVID-19) and control. Numer. Methods Partial Differ. Equ. 38, 760–776. 10.1002/num.22695 (2022).33362341 10.1002/num.22695
9. Savi PV Savi MA Borges B A mathematical description of the dynamics of coronavirus disease 2019 (COVID-19): A case study of Brazil Comput. Math. Methods Med. 2020 10.1155/2020/9017157 33029196
Savi, P. V., Savi, M. A. & Borges, B. A mathematical description of the dynamics of coronavirus disease 2019 (COVID-19): A case study of Brazil. Comput. Math. Methods Med.10.1155/2020/9017157 (2020).33029196 10.1155/2020/9017157
10. Tiwari V Deyal N Bisht NS Mathematical modeling based study and prediction of COVID-19 epidemic dissemination under the impact of lockdown in India Front. Phys. 2020 10.3389/fphy.2020.586899
Tiwari, V., Deyal, N. & Bisht, N. S. Mathematical modeling based study and prediction of COVID-19 epidemic dissemination under the impact of lockdown in India. Front. Phys.10.3389/fphy.2020.586899 (2020).10.3389/fphy.2020.586899
11. Warbhe SD Lamba NK Deshmukh KC Impact of COVID-19: A mathematical model J. Interdiscip. Math. 2021 24 77 87 10.1080/09720502.2020.1833444
Warbhe, S. D., Lamba, N. K. & Deshmukh, K. C. Impact of COVID-19: A mathematical model. J. Interdiscip. Math. 24, 77–87. 10.1080/09720502.2020.1833444 (2021).10.1080/09720502.2020.1833444
12. Daniel Deborah O Mathematical model for the transmission of Covid-19 with nonlinear forces of infection and the need for prevention measure in Nigeria J. Infect. Dis. Epidemiol. 2020 10.23937/2474-3658/1510158
Daniel Deborah, O. Mathematical model for the transmission of Covid-19 with nonlinear forces of infection and the need for prevention measure in Nigeria. J. Infect. Dis. Epidemiol.10.23937/2474-3658/1510158 (2020).10.23937/2474-3658/1510158
13. Prathumwan D Trachoo K Chaiya I Mathematical modeling for prediction dynamics of the coronavirus disease 2019 (COVID-19) pandemic, quarantine control measures Symmetry 2020 10.3390/SYM12091404
Prathumwan, D., Trachoo, K. & Chaiya, I. Mathematical modeling for prediction dynamics of the coronavirus disease 2019 (COVID-19) pandemic, quarantine control measures. Symmetry10.3390/SYM12091404 (2020).10.3390/SYM12091404
14. Balike Dieudonné Z Mathematical model for the mitigation of the economic effects of the Covid-19 in the Democratic Republic of the Congo PLoS ONE 2021 16 e0250775 10.1371/journal.pone.0250775 33939724
Balike Dieudonné, Z. Mathematical model for the mitigation of the economic effects of the Covid-19 in the Democratic Republic of the Congo. PLoS ONE 16, e0250775. 10.1371/journal.pone.0250775 (2021).33939724 10.1371/journal.pone.0250775
15. Raslan WE Fractional mathematical modeling for epidemic prediction of COVID-19 in Egypt Ain Shams Eng. J. 2021 12 3057 3062 10.1016/j.asej.2020.10.027
Raslan, W. E. Fractional mathematical modeling for epidemic prediction of COVID-19 in Egypt. Ain Shams Eng. J. 12, 3057–3062. 10.1016/j.asej.2020.10.027 (2021).10.1016/j.asej.2020.10.027
16. Chen TM Rui J Wang QP Zhao ZY Cui JA Yin L A mathematical model for simulating the phase-based transmissibility of a novel coronavirus Infect. Dis. Poverty 2020 10.1186/s40249-020-00640-3 33380335
Chen, T. M. et al. A mathematical model for simulating the phase-based transmissibility of a novel coronavirus. Infect. Dis. Poverty10.1186/s40249-020-00640-3 (2020).33380335 10.1186/s40249-020-00640-3
17. Sinaga LP Nasution H Kartika D Stability analysis of the corona virus (Covid-19) dynamics SEIR model in Indonesia J. Phys. Conf. Ser. 2021 10.1088/1742-6596/1819/1/012043
Sinaga, L. P., Nasution, H. & Kartika, D. Stability analysis of the corona virus (Covid-19) dynamics SEIR model in Indonesia. J. Phys. Conf. Ser.10.1088/1742-6596/1819/1/012043 (2021).10.1088/1742-6596/1819/1/012043
18. James LP Salomon JA Buckee CO Menzies NA The use and misuse of mathematical modeling for infectious disease policymaking: Lessons for the COVID-19 pandemic Med. Decis. Mak. 2021 41 379 385 10.1177/0272989X21990391
James, L. P., Salomon, J. A., Buckee, C. O. & Menzies, N. A. The use and misuse of mathematical modeling for infectious disease policymaking: Lessons for the COVID-19 pandemic. Med. Decis. Mak. 41, 379–385. 10.1177/0272989X21990391 (2021).10.1177/0272989X21990391
19. Ameen IG Ali HM Alharthi MR Abdel-Aty AH Elshehabey HM Investigation of the dynamics of COVID-19 with a fractional mathematical model: A comparative study with actual data Results Phys. 2021 10.1016/j.rinp.2021.103976 33623732
Ameen, I. G., Ali, H. M., Alharthi, M. R., Abdel-Aty, A. H. & Elshehabey, H. M. Investigation of the dynamics of COVID-19 with a fractional mathematical model: A comparative study with actual data. Results Phys.10.1016/j.rinp.2021.103976 (2021).33623732 10.1016/j.rinp.2021.103976
20. Jiang S Li Q Li C Liu S He X Wang T Li H Corpe C Zhang X Xu J Mathematical models for devising the optimal SARS-CoV-2 strategy for eradication in China, South Korea, and Italy J. Transl. Med. 2020 10.1186/s12967-020-02513-7 33208161
Jiang, S. et al. Mathematical models for devising the optimal SARS-CoV-2 strategy for eradication in China, South Korea, and Italy. J. Transl. Med.10.1186/s12967-020-02513-7 (2020).33208161 10.1186/s12967-020-02513-7
21. Uddin, M. S., Nasseef, M. T., Mahmud, M., AlArjani, A. Mathematical modelling in prediction of novel CoronaVirus (COVID-19) transmission dynamics (2020).
22. Kahn R Holmdahl I Reddy S Jernigan J Mina MJ Slayton RB Mathematical modeling to inform vaccination strategies and testing approaches for coronavirus disease 2019 (COVID-19) in nursing homes Clin. Infect. Dis. 2022 74 597 603 10.1093/cid/ciab517 34086877
Kahn, R. et al. Mathematical modeling to inform vaccination strategies and testing approaches for coronavirus disease 2019 (COVID-19) in nursing homes. Clin. Infect. Dis. 74, 597–603. 10.1093/cid/ciab517 (2022).34086877 10.1093/cid/ciab517
23. Nave O Hartuv I Shemesh U Θ-SEIHRD mathematical model of Covid19-stability analysis using fast–slow decomposition PeerJ 2020 10.7717/peerj.10019 33005495
Nave, O., Hartuv, I. & Shemesh, U. Θ-SEIHRD mathematical model of Covid19-stability analysis using fast–slow decomposition. PeerJ10.7717/peerj.10019 (2020).33005495 10.7717/peerj.10019
24. Kim BN Kim E Lee S Oh C Mathematical model of Covid-19 transmission dynamics in South Korea: The impacts of travel restrictions, social distancing, and early detection Processes 2020 8 1 18 10.3390/pr8101304
Kim, B. N., Kim, E., Lee, S. & Oh, C. Mathematical model of Covid-19 transmission dynamics in South Korea: The impacts of travel restrictions, social distancing, and early detection. Processes 8, 1–18. 10.3390/pr8101304 (2020).10.3390/pr8101304
25. Das M Samanta G Stability analysis of a fractional ordered COVID-19 model Comput. Math. Biophys. 2021 9 22 45 10.1515/cmb-2020-0116
Das, M. & Samanta, G. Stability analysis of a fractional ordered COVID-19 model. Comput. Math. Biophys. 9, 22–45. 10.1515/cmb-2020-0116 (2021).10.1515/cmb-2020-0116
26. Machado JAT Lopes AM Rare and extreme events: The case of COVID-19 pandemic Nonlinear Dyn. 2020 100 2953 2972 10.1007/s11071-020-05680-w 32427206
Machado, J. A. T. & Lopes, A. M. Rare and extreme events: The case of COVID-19 pandemic. Nonlinear Dyn. 100, 2953–2972. 10.1007/s11071-020-05680-w (2020).32427206 10.1007/s11071-020-05680-w
27. Rajagopal K Hasanzadeh N Parastesh F Hamarash II Jafari S Hussain I A fractional-order model for the novel coronavirus (COVID-19) outbreak Nonlinear Dyn. 2020 101 711 718 10.1007/s11071-020-05757-6 32836806
Rajagopal, K. et al. A fractional-order model for the novel coronavirus (COVID-19) outbreak. Nonlinear Dyn. 101, 711–718. 10.1007/s11071-020-05757-6 (2020).32836806 10.1007/s11071-020-05757-6
28. Quaranta G Formica G Machado JT Lacarbonara W Masri SF Understanding COVID-19 nonlinear multi-scale dynamic spreading in Italy Nonlinear Dyn. 2020 101 1583 1619 10.1007/s11071-020-05902-1 32904911
Quaranta, G., Formica, G., Machado, J. T., Lacarbonara, W. & Masri, S. F. Understanding COVID-19 nonlinear multi-scale dynamic spreading in Italy. Nonlinear Dyn. 101, 1583–1619. 10.1007/s11071-020-05902-1 (2020).32904911 10.1007/s11071-020-05902-1
29. Aslam Noor M Raza A Arif MS Rafiq M Sooppy Nisar K Khan I Abdelwahab SF Non-standard computational analysis of the stochastic COVID-19 pandemic model: An application of computational biology Alex. Eng. J. 2022 61 619 630 10.1016/j.aej.2021.06.039
Aslam Noor, M. et al. Non-standard computational analysis of the stochastic COVID-19 pandemic model: An application of computational biology. Alex. Eng. J. 61, 619–630. 10.1016/j.aej.2021.06.039 (2022).10.1016/j.aej.2021.06.039
30. Macías-Díaz JE Raza A Ahmed N Rafiq M Analysis of a nonstandard computer method to simulate a nonlinear stochastic epidemiological model of coronavirus-like diseases Comput. Methods Programs Biomed. 2021 10.1016/j.cmpb.2021.106054 34715516
Macías-Díaz, J. E., Raza, A., Ahmed, N. & Rafiq, M. Analysis of a nonstandard computer method to simulate a nonlinear stochastic epidemiological model of coronavirus-like diseases. Comput. Methods Programs Biomed.10.1016/j.cmpb.2021.106054 (2021).34715516 10.1016/j.cmpb.2021.106054
31. Shahid N Baleanu D Ahmed N Shaikh TS Raza A Iqbal MS Rafiq M Aziz-Ur Rehman M Optimality of solution with numerical investigation for coronavirus epidemic model Comput. Mater. Contin. 2021 67 1713 1728 10.32604/cmc.2021.014191
Shahid, N. et al. Optimality of solution with numerical investigation for coronavirus epidemic model. Comput. Mater. Contin. 67, 1713–1728. 10.32604/cmc.2021.014191 (2021).10.32604/cmc.2021.014191
32. Shatanawi W Raza A Arif MS Abodayeh K Rafiq M Bibi M An effective numerical method for the solution of a stochastic coronavirus (2019-NCovid) pandemic model Comput. Mater. Contin. 2020 66 1121 1137 10.32604/cmc.2020.012070
Shatanawi, W. et al. An effective numerical method for the solution of a stochastic coronavirus (2019-NCovid) pandemic model. Comput. Mater. Contin. 66, 1121–1137. 10.32604/cmc.2020.012070 (2020).10.32604/cmc.2020.012070
33. Naveed M Rafiq M Raza A Ahmed N Khan I Nisar KS Soori AH Mathematical analysis of novel coronavirus (2019-NCov) delay pandemic model Comput. Mater. Contin. 2020 64 1401 1414 10.32604/cmc.2020.011314
Naveed, M. et al. Mathematical analysis of novel coronavirus (2019-NCov) delay pandemic model. Comput. Mater. Contin. 64, 1401–1414. 10.32604/cmc.2020.011314 (2020).10.32604/cmc.2020.011314
34. Ghosh S Samanta G Nieto JJ Application of non-parametric models for analyzing survival data of COVID-19 patients J. Infect. Public Health 2021 14 1328 1333 10.1016/j.jiph.2021.08.025 34479820
Ghosh, S., Samanta, G. & Nieto, J. J. Application of non-parametric models for analyzing survival data of COVID-19 patients. J. Infect. Public Health 14, 1328–1333. 10.1016/j.jiph.2021.08.025 (2021).34479820 10.1016/j.jiph.2021.08.025
35. Das M Samanta GP Optimal control of fractional order COVID-19 epidemic spreading in Japan and India 2020 Biophys. Rev. Lett. 2020 15 207 236 10.1142/s179304802050006x
Das, M. & Samanta, G. P. Optimal control of fractional order COVID-19 epidemic spreading in Japan and India 2020. Biophys. Rev. Lett. 15, 207–236. 10.1142/s179304802050006x (2020).10.1142/s179304802050006x
36. Saha S Samanta GP Nieto JJ Epidemic model of COVID-19 outbreak by inducing behavioural response in population Nonlinear Dyn. 2020 102 455 487 10.1007/s11071-020-05896-w 32863581
Saha, S., Samanta, G. P. & Nieto, J. J. Epidemic model of COVID-19 outbreak by inducing behavioural response in population. Nonlinear Dyn. 102, 455–487. 10.1007/s11071-020-05896-w (2020).32863581 10.1007/s11071-020-05896-w
37. Saha S Samanta GP Modelling the role of optimal social distancing on disease prevalence of COVID-19 epidemic Int. J. Dyn. Control 2021 9 1053 1077 10.1007/s40435-020-00721-z 33194535
Saha, S. & Samanta, G. P. Modelling the role of optimal social distancing on disease prevalence of COVID-19 epidemic. Int. J. Dyn. Control 9, 1053–1077. 10.1007/s40435-020-00721-z (2021).33194535 10.1007/s40435-020-00721-z
38. Raza A Rafiq M Awrejcewicz J Ahmed N Mohsin M Dynamical analysis of coronavirus disease with crowding effect, and vaccination: A study of third strain Nonlinear Dyn. 2022 107 3963 3982 10.1007/s11071-021-07108-5 35002076
Raza, A., Rafiq, M., Awrejcewicz, J., Ahmed, N. & Mohsin, M. Dynamical analysis of coronavirus disease with crowding effect, and vaccination: A study of third strain. Nonlinear Dyn. 107, 3963–3982. 10.1007/s11071-021-07108-5 (2022).35002076 10.1007/s11071-021-07108-5
39. Sontakke BR Shaikh AS Properties of Caputo operator and its applications to linear fractional differential equations Int. J. Eng. Res. Appl. 2015 5 22 27
Sontakke, B. R. & Shaikh, A. S. Properties of Caputo operator and its applications to linear fractional differential equations. Int. J. Eng. Res. Appl. 5, 22–27 (2015).
40. Albadarneh RB Batiha I Alomari AK Tahat N Numerical approach for approximating the Caputo fractional-order derivative operator AIMS Math. 2021 6 12743 12756 10.3934/math.2021735
Albadarneh, R. B., Batiha, I., Alomari, A. K. & Tahat, N. Numerical approach for approximating the Caputo fractional-order derivative operator. AIMS Math. 6, 12743–12756 (2021).10.3934/math.2021735
41. Ray SS Atangana A Noutchie SCO Kurulay M Bildik N Kilicman A Fractional calculus and its applications in applied mathematics and other sciences Math. Probl. Eng. 2014 2014 2 4 10.1155/2014/849395
Ray, S. S. et al. Fractional calculus and its applications in applied mathematics and other sciences. Math. Probl. Eng. 2014, 2–4 (2014).10.1155/2014/849395
42. Garrappa R Kaslik E Evaluation of fractional integrals and derivatives of elementary functions: Overview and tutorial Mathematics 2019 2 1 21 10.3390/math7050407
Garrappa, R. & Kaslik, E. Evaluation of fractional integrals and derivatives of elementary functions: Overview and tutorial. Mathematics 2, 1–21. 10.3390/math7050407 (2019).10.3390/math7050407
43. Calculus, F. Fractional Operators and Their (2018). (ISBN 9780128096703).
44. Atangana A Secer A A note on fractional order derivatives and table of fractional derivatives of some special functions Abstract Appl. Anal. 2013 2013 1 9
Atangana, A. & Secer, A. A note on fractional order derivatives and table of fractional derivatives of some special functions. Abstract Appl. Anal. 2013, 1–9 (2013).
45. Li Y Chen Y Podlubny I Stability of fractional-order nonlinear dynamic systems: Lyapunov direct method and generalized Mittag–Leffler stability Comput. Math. Appl. 2010 59 1810 1821 10.1016/j.camwa.2009.08.019
Li, Y., Chen, Y. & Podlubny, I. Stability of fractional-order nonlinear dynamic systems: Lyapunov direct method and generalized Mittag–Leffler stability. Comput. Math. Appl. 59, 1810–1821. 10.1016/j.camwa.2009.08.019 (2010).10.1016/j.camwa.2009.08.019
46. Baba IA Nasidi BA Fractional order epidemic model for the dynamics of novel COVID-19 Alex. Eng. J. 2021 60 537 548 10.1016/j.aej.2020.09.029
Baba, I. A. & Nasidi, B. A. Fractional order epidemic model for the dynamics of novel COVID-19. Alex. Eng. J. 60, 537–548. 10.1016/j.aej.2020.09.029 (2021).10.1016/j.aej.2020.09.029
47. Sweilam N Nagy AM Nonstandard finite difference scheme for the fractional order Salmonella transmission model J. Fract. Calculus Appl. 2019 10 1 197
Sweilam, N. & Nagy, A. M. Nonstandard finite difference scheme for the fractional order Salmonella transmission model. J. Fract. Calculus Appl. 10(1), 197 (2019).
48. Arenas AJ Gonz G Benito M Arenas AJ Gonz G Construction of nonstandard finite difference schemes for the SI and SIR epidemic models of fractional order Math. Comput. Simul. 2015 10.1016/j.matcom.2015.09.001
Arenas, A. J., Gonz, G., Benito, M., Arenas, A. J. & Gonz, G. Construction of nonstandard finite difference schemes for the SI and SIR epidemic models of fractional order. Math. Comput. Simul.10.1016/j.matcom.2015.09.001 (2015).10.1016/j.matcom.2015.09.001
49. Scherer R Kalla SL Tang Y Huang J The Grünwald–Letnikov method for fractional differential equations Comput. Math. Appl. 2011 62 902 917 10.1016/j.camwa.2011.03.054
Scherer, R., Kalla, S. L., Tang, Y. & Huang, J. The Grünwald–Letnikov method for fractional differential equations. Comput. Math. Appl. 62, 902–917. 10.1016/j.camwa.2011.03.054 (2011).10.1016/j.camwa.2011.03.054
50. Uçar E Özdemir N A fractional model of cancer-immune system with Caputo and Caputo–Fabrizio derivatives Eur. Phys. J. Plus 2021 136 1 17 10.1140/epjp/s13360-020-00966-9
Uçar, E. & Özdemir, N. A fractional model of cancer-immune system with Caputo and Caputo–Fabrizio derivatives. Eur. Phys. J. Plus 136, 1–17 (2021).10.1140/epjp/s13360-020-00966-9
51. Ucar E Özdemir N Altun E Fractional order model of immune cells influenced by cancer cells Math. Model. Natl. Phenom. 2019 14 3 308 10.1051/mmnp/2019002
Ucar, E., Özdemir, N. & Altun, E. Fractional order model of immune cells influenced by cancer cells. Math. Model. Natl. Phenom. 14(3), 308 (2019).10.1051/mmnp/2019002
52. Uçar S Analysis of hepatitis B disease with fractal–fractional Caputo derivative using real data from Turkey J. Comput. Appl. Math. 2023 419 114692 10.1016/j.cam.2022.114692
Uçar, S. Analysis of hepatitis B disease with fractal–fractional Caputo derivative using real data from Turkey. J. Comput. Appl. Math. 419, 114692 (2023).10.1016/j.cam.2022.114692
53. Özdemir N Uçar S Billur İskender Eroğlu B Dynamical analysis of fractional order model for computer virus propagation with kill signals Int. J. Nonlinear Sci. Numer. Simul. 2020 21 3–4 239 247 10.1515/ijnsns-2019-0063
Özdemir, N., Uçar, S. & Billur İskender Eroğlu, B. Dynamical analysis of fractional order model for computer virus propagation with kill signals. Int. J. Nonlinear Sci. Numer. Simul. 21(3–4), 239–247 (2020).10.1515/ijnsns-2019-0063
54. Esmehan UÇAR Examining of a tumor system with Caputo derivative Balıkesir Üniv. Fen Bilimleri Enstitüsü Dergisi 2023 25 1 37 48 10.25092/baunfbed.1113646
Esmehan, U. Ç. A. R. Examining of a tumor system with Caputo derivative. Balıkesir Üniv. Fen Bilimleri Enstitüsü Dergisi 25(1), 37–48 (2023).10.25092/baunfbed.1113646
55. Ahmad A Farman M Naik PA Zafar N Akgul A Saleem MU Modeling and numerical investigation of fractional-order bovine babesiosis disease Numer. Methods Partial Differ. Equ. 2021 37 3 1946 1964 10.1002/num.22632
Ahmad, A. et al. Modeling and numerical investigation of fractional-order bovine babesiosis disease. Numer. Methods Partial Differ. Equ. 37(3), 1946–1964 (2021).10.1002/num.22632
56. Ghori MB Naik PA Zu J Eskandari Z Naik MUD Global dynamics and bifurcation analysis of a fractional-order SEIR epidemic model with saturation incidence rate Math. Methods Appl. Sci. 2022 45 7 3665 3688 10.1002/mma.8010
Ghori, M. B., Naik, P. A., Zu, J., Eskandari, Z. & Naik, M. U. D. Global dynamics and bifurcation analysis of a fractional-order SEIR epidemic model with saturation incidence rate. Math. Methods Appl. Sci. 45(7), 3665–3688 (2022).10.1002/mma.8010
57. Naik PA Zehra A Farman M Shehzad A Shahzeen S Huang Z Forecasting and dynamical modeling of reversible enzymatic reactions with hybrid proportional fractional derivative Front. Phys. 2024 11 1307307 10.3389/fphy.2023.1307307
Naik, P. A. et al. Forecasting and dynamical modeling of reversible enzymatic reactions with hybrid proportional fractional derivative. Front. Phys. 11, 1307307 (2024).10.3389/fphy.2023.1307307
58. Naik PA Yavuz M Qureshi S Zu J Townley S Modeling and analysis of COVID-19 epidemics with treatment in fractional derivatives using real data from Pakistan Eur. Phys. J. Plus 2020 135 1 42 10.1140/epjp/s13360-020-00819-5
Naik, P. A., Yavuz, M., Qureshi, S., Zu, J. & Townley, S. Modeling and analysis of COVID-19 epidemics with treatment in fractional derivatives using real data from Pakistan. Eur. Phys. J. Plus 135, 1–42 (2020).10.1140/epjp/s13360-020-00819-5
59. Naik PA Global dynamics of a fractional-order SIR epidemic model with memory Int. J. Biomath. 2020 13 08 2050071 10.1142/S1793524520500710
Naik, P. A. Global dynamics of a fractional-order SIR epidemic model with memory. Int. J. Biomath. 13(08), 2050071 (2020).10.1142/S1793524520500710
