Next Article in Journal
Study on the Nonlinear Dynamics of the (3+1)-Dimensional Jimbo-Miwa Equation in Plasma Physics
Previous Article in Journal
Knowable Moments in Stochastics: Knowing Their Advantages
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Dynamics of a Fractional-Order COVID-19 Epidemic Model with Quarantine and Standard Incidence Rate

1
Department of Mathematics, Faculty of Mathematics and Natural Sciences, University of Brawijaya, Jl. Veteran, Malang 65145, Indonesia
2
Department of Mathematics Education, Faculty of Pedagogy and Psychology, PGRI Wiranegara University, Pasuruan 67118, Indonesia
*
Author to whom correspondence should be addressed.
Axioms 2023, 12(6), 591; https://doi.org/10.3390/axioms12060591
Submission received: 8 April 2023 / Revised: 31 May 2023 / Accepted: 8 June 2023 / Published: 15 June 2023
(This article belongs to the Topic Advances in Nonlinear Dynamics: Methods and Applications)

Abstract

:
In this paper, we propose a fractional-order COVID-19 epidemic model with a quarantine and standard incidence rate using the Caputo fractional-order derivative. The model consists of six classes: susceptible (S), exposed (E), infected (I), quarantined (Q), recovered (R), and deceased (M). In our proposed model, we simultaneously consider the recovery rate and quarantine rate of infected individuals, which has not been considered in other fractional-order COVID-19 epidemic models. Furthermore, we consider the standard incidence rate in the model. For our proposed model, we prove the existence, uniqueness, non-negativity, and boundedness of the solution. The model has two equilibrium points: disease-free equilibrium and endemic equilibrium. Implementing the spectral radius of the next-generation matrix, we obtain the basic reproduction number ( R 0 ). The disease-free equilibrium always exists and is locally and globally asymptotically stable only if R 0 < 1 . On the other hand, endemic equilibrium exists and is globally asymptotically stable if R 0 > 1 . Our numerical simulation confirms the stability properties of the equilibrium. The smaller the order of the derivative, the slower the convergence of the solution of the model. Both the recovery rate and quarantine rate of the infected class are important parameters determining the stability of the equilibrium point. Based on parameter estimation from COVID-19 data in Indonesia, the fractional-order model has better performance than the first-order model for both the calibration and 20-day forecasting of confirmed daily active cases of COVID-19.

1. Introduction

Coronavirus Disease 2019, or COVID-19, is an acute respiratory disease that was first indicated in Wuhan Province, China, in 2019 [1]. COVID-19 is an infectious disease caused by infection with the SARS-CoV-2 virus [2]. Symptoms of this disease are fever, flu, cough, pneumonia, and other symptoms that cause complications such as asthma [1,3]. COVID-19 can spread quickly, and the total number of cases reached 624 million cases in the world [4]. The number of cases is also influenced by the variant of the COVID-19 virus. In general, there are five variants of COVID-19 virus: Alpha, Beta, Gamma, Delta, and Omicron [5,6]. Various studies on how COVID-19 is transmitted have been described in [7,8,9]. The mechanism of transmission most often occurs through droplets, directly or indirectly [2,10,11]. To prevent the spread of the virus, quarantine is an effective intervention [11,12,13]. Quarantine causes infected individuals to lose contact with other individuals, which can wipe out the transmission of the disease [14]. Thus, the mechanism of transmission with quarantine interventions is important to understand and study further.
Mathematical modeling is an approximation used to understand the mechanism of transmission and current state of COVID-19 [15]. If COVID-19 data are sufficient, then the forecasting of COVID-19 via parameter estimation is very helpful for policy makers, as was conducted in [16,17,18]. Otherwise, the dynamics of an epidemic model can provide conditions that lead to an endemic disease. Some researchers have studied the COVID-19 pandemic via the dynamics of mathematical models, such as in [15,19,20]. The mathematical model for a disease epidemic is basically represented by first-order differential equations.
Some researchers have developed COVID-19 epidemic models. Zeb et al. [21] proposed a COVID-19 epidemic model with a quarantined class. The paper of Zeb et al. [21] focuses on the dynamics of the COVID-19 epidemic model. In the model, there is an assumption that exposed individuals can transmit the COVID-19 virus. The current literature shows that an exposed individual is an infected individual but does not always transmit the virus [22,23]. Furthermore, the mortality of COVID-19 is not considered in the model of Zeb et al. [21]. Another COVID-19 epidemic model was proposed by López and Rodó [24] with a quarantined class and protected class. They focus on the calibration of COVID-19 data in Italy and Spain. In their model, there is an assumption that the death class only comes from the quarantined class. Musafir et al. [15] stated that infected individuals can also cause death from COVID-19. In both models, the incidence rate used is bilinear, with a maximum proportion of the total of living classes. A recent model has been proposed to resolve the weak assumption for these models. The proposed model by Darti et al. [25] assumes that the exposed class is not contagious and the infected class may induce death from COVID-19. These assumptions improve the proposed models by Zeb et al. [21] and López and Rodó [24]. All three models are represented by nonlinear first-order differential equations.
In recent decades, fractional calculus has been developed to provide some integral and derivative options [26,27]. A fractional derivative is a non-local operator, which means the rate of the function depends on previous histories from the initial value, which is called the memory effect [26]. In epidemiology, the future spread of a disease depends on recent and previous conditions because of the impact of both heredity and history on individual habits [28]. Hence, the fractional derivative is more realistic than the integer derivative. Moreover, some researchers, such as [2,20,29,30], state that fractional derivatives provide better accuracy and flexibility than integer derivatives. There are two types of fractional operators obtained by generalizing the integer-order derivative, namely, the Riemann–Liouville operator and the Caputo operator [31]. Li et al. [32] show that the Caputo operator is better than the Riemann–Liouville operator in terms of applied mathematics.
Recently, many researchers have proposed a fractional-order COVID-19 epidemic model. Yadav et al. [33] introduced a fractional-order SEIR-type model by dividing the infected class into two subclasses: the infected subclass without symptoms and the infected subclass with symptoms. The proposed model was analyzed dynamically and solved using q-homotopy analysis with the Sumudu transform method. With the increasing incidences of COVID-19, various measures have been taken by many countries to control the spread of COVID-19. One of the well-known preventive measures is a lockdown. The impact of lockdowns has been modeled by [34], i.e., by considering a fractional-order SI-type epidemic model where the population is split into the subpopulation not under lockdown and the subpopulation under lockdown. In [35], the authors considered a fractional SEIR model incorporating public risk awareness and the spread of disease caused by both human-to-human and zoonotic transmissions. Then, in [36], the authors extended the model of [33] by considering the classes of quarantine, isolation, and environmental impacts.
We notice that although the fractional-order COVID-19 epidemic models above were derived based on the Caputo fractional derivative operator, they did not pay attention to the inconsistency in the time dimension. Such dimensional inconsistency is caused by replacing the first-order derivative with a Caputo fractional-order derivative. Moreover, none of those models have considered the effects of fractional-order derivatives. In this study, we propose a fractional-order COVID-19 epidemic model by including the dimensional issue. We implement the proposed model using real COVID-19 data. Since the available data only consist of the confirmed active COVID-19 cases (the infected subclass), without noting whether this is with or without symptoms, we consider a simpler model than proposed by [36]. For this aim, we modify the COVID-19 epidemic model proposed by Darti et al. [25] to a fractional-order model with a Caputo fractional operator. The proposed COVID-19 epidemic model by Darti et al. [25] is given as follows, with the parameter definitions listed in Table 1:
d S ( t ) d t = Λ ^ β ^ S ( t ) I ( t ) N ( t ) μ ^ S ( t ) , d E ( t ) d t = β ^ S ( t ) I ( t ) N ( t ) ( γ ^ + μ ^ ) E ( t ) , d I ( t ) d t = γ ^ E ( t ) ( σ ^ + η ^ + δ ^ + μ ^ ) I ( t ) , d Q ( t ) d t = σ ^ I ( t ) ( ϑ ^ + κ ^ + μ ^ ) Q ( t ) , d R ( t ) d t = η ^ I ( t ) + ϑ ^ Q ( t ) μ ^ R ( t ) , d M ( t ) d t = δ ^ I ( t ) + κ ^ Q ( t ) .
In this study, we first propose the fractional-order COVID-19 epidemic model in Section 2. The basic properties of the proposed model are presented in Section 3 and Section 4. The equilibrium points of the model and their stability are given in Section 5 and Section 6, respectively. To measure how sensitive the model parameters are to COVID-19 transmission, in Section 7, we provide a sensitivity analysis. In Section 8, we numerically simulate the properties of the dynamic analysis of our proposed model. The numerical solutions are obtained by adopting the predictor–corrector scheme proposed by [37]. In Section 9, we implement our proposed model using COVID-19 data in Indonesia. Parameter estimation is performed using the nonlinear least square method. This can be achieved using MATLAB’s built-in function lsqcurvefit. Here, we minimize the squared error between numerical solutions for infected subpopulations and confirmed active daily data. Finally, we remark upon the conclusions of this study in Section 10.

2. Model Formulation

We modify the model proposed by Darti et al. [25] to a fractional-order model with the Caputo operator. The Caputo fractional operator is given as follows.
Definition 1 
([38]). Let t 0 and α ( 0 , 1 ) . The Riemann–Liouville integral of order α is defined by
R L I t α x ( t ) = R L D t α x ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 x ( τ ) d τ
where x C n [ 0 , ) , R and Γ ( · ) is a Gamma function.
Definition 2 
([38]). Let t 0 and α ( 0 , 1 ) . The Caputo derivative of order α is defined by
C D t α x ( t ) = 1 Γ ( 1 α ) 0 t ( t τ ) α x ( τ ) d τ ,
where x C n [ 0 , ) , R .
Our proposed model has five population classes: the susceptible class (S), the exposed class (E), the infected class (I), the quarantined class (Q), the recovered class (R), and the death class (M). We define variable N as the total of living individuals, which is denoted as N = S + E + I + Q + R . The individuals of the susceptible class (S) who are in close contact with infected individuals have the potential to be exposed individuals. The exposed class is infected but not yet infectious. After the incubation period of the virus, exposed individuals move to the infected class (I), which is composed of the infected individuals. The infected individuals can be either quarantined, recovered, or dead from illness. The quarantined individuals (Q) can be either recovered or dead. The interaction of all classes is depicted in Figure 1. The proposed fractional-order SEIQRM model with the Caputo derivative operator is derived by replacing the first-order derivatives in the SEIQRM model (1) with the Caputo fractional derivatives. However, if we perform this step, then there is a time-dimensional inconsistency on both sides of all the equations because the fractional derivatives of order α have the dimension of ( t i m e ) α . To handle the dimension inconsistency, we rescale the parameters and obtain the following fractional-order SEIQRM model:
C D t α S ( t ) = Λ β S ( t ) I ( t ) N ( t ) μ S ( t ) , C D t α E ( t ) = β S ( t ) I ( t ) N ( t ) ( γ + μ ) E ( t ) , C D t α I ( t ) = γ E ( t ) ( σ + η + δ + μ ) I ( t ) , C D t α Q ( t ) = σ I ( t ) ( ϑ + κ + μ ) Q ( t ) , C D t α R ( t ) = η I ( t ) + ϑ Q ( t ) μ R ( t ) , C D t α M ( t ) = δ I ( t ) + κ Q ( t ) ,
where Λ = Λ ^ α , β = β ^ α , μ = μ ^ α , γ = γ ^ α , σ = σ ^ α , η = η ^ α , δ = δ ^ α , ϑ = ϑ ^ α , and κ = κ ^ α .

3. Existence and Uniqueness of Solution

In this section, we explicitly prove the existence and uniqueness of the solution of model (2).
Theorem 1. 
Let the initial values of model (2) be non-negative. The model (2) has a unique solution in R + 5 for all t 0 .
Proof. 
If the right-hand side of model (2) is written by
f 1 = Λ β S ( t ) I ( t ) N ( t ) μ S ( t ) , f 2 = β S ( t ) I ( t ) N ( t ) ( γ + μ ) E ( t ) , f 3 = γ E ( t ) ( σ + η + δ + μ ) I ( t ) , f 4 = σ I ( t ) ( ϑ + κ + μ ) Q ( t ) , f 5 = η I ( t ) + ϑ Q ( t ) μ R ( t ) ,
since I N 1 , we find for every S , S ¯ R + 5 that
f 1 ( S , t ) f 1 ( S ¯ , t ) = β S I N μ S + β S ¯ I N + μ S ¯ p 1 S S ¯ ,
with p 1 = ( β + μ ) . In the same way, we can show that
f 2 ( E , t ) f 2 ( E ¯ , t ) p 2 E E ¯ , f 3 ( I , t ) f 3 ( I ¯ , t ) p 3 I I ¯ , f 4 ( Q , t ) f 4 ( Q ¯ , t ) p 4 Q Q ¯ , f 5 ( R , t ) f 5 ( R ¯ , t ) p 5 R R ¯ ,
where p 2 = γ + μ , p 3 = σ + η + δ + μ , p 4 = ϑ + κ + μ , p 5 = μ . It is clear that if 0 < p i < 1 , i = 1 , 2 , , 5 , then all kernels of f i contract and satisfy the Lipschitz condition.
If both sides of model (2) are integrated by Definition 1, then we obtain
S ( t ) S ( 0 ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 S , τ d τ , E ( t ) E ( 0 ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 2 E , τ d τ , I ( t ) I ( 0 ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 3 I , τ d τ , Q ( t ) Q ( 0 ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 4 Q , τ d τ , R ( t ) R ( 0 ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 5 R , τ d τ .
Hence, we obtain the iterative scheme
S n ( t ) = S ( 0 ) + 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 S n 1 , τ d τ , E n ( t ) = E ( 0 ) + 1 Γ ( α ) 0 t ( t τ ) α 1 f 2 E n 1 , τ d τ , I n ( t ) = I ( 0 ) + 1 Γ ( α ) 0 t ( t τ ) α 1 f 3 I n 1 , τ d τ , Q n ( t ) = Q ( 0 ) + 1 Γ ( α ) 0 t ( t τ ) α 1 f 4 Q n 1 , τ d τ , R n ( t ) = R ( 0 ) + 1 Γ ( α ) 0 t ( t τ ) α 1 f 5 R n 1 , τ d τ .
Next, we construct a recursive formulation based on the Equation (4) as follows:
ψ 1 n ( t ) = S n ( t ) S n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 ( S n 1 , τ ) f 1 ( S n 2 , τ ) d τ , ψ 2 n ( t ) = E n ( t ) E n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 2 ( E n 1 , τ ) f 2 ( E n 2 , τ ) d τ , ψ 3 n ( t ) = I n ( t ) I n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 3 ( I n 1 , τ ) f 3 ( I n 2 , τ ) d τ , ψ 4 n ( t ) = Q n ( t ) Q n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 4 ( Q n 1 , τ ) f 4 ( Q n 2 , τ ) d τ , ψ 5 n ( t ) = R n ( t ) R n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 5 ( R n 1 , τ ) f 5 ( R n 2 , τ ) d τ .
The norm of the first equation of (5) can be expressed as
ψ 1 n ( t ) = S n ( t ) S n 1 ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 ( S n 1 , τ ) f 1 ( S n 2 , τ ) d τ 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 ( S n 1 , τ ) f 1 ( S n 2 , τ ) d τ p 1 Γ ( α ) 0 t ( t τ ) α 1 S n 1 ( τ ) S n 2 ( τ ) d τ = p 1 Γ ( α ) 0 t ( t τ ) α 1 ψ 1 ( n 1 ) ( τ ) d τ = p 1 t Γ ( α + 1 ) n S ( 0 ) .
Similarly, we have
ψ j n ( t ) p 1 t Γ ( α + 1 ) n X j ( 0 ) , f o r   a l l   j = 2 , 3 , 4 , 5 ,
where ( X 1 , X 2 , X 3 , X 4 , X 5 ) = ( S , E , I , Q , R ) . Assume that there is t = T such that
p j T Γ ( α + 1 ) < 1 , f o r   a l l   j = 1 , 2 , , 5 .
It is clear that for every j = 1 , 2 , , 5 , ψ j n ( t ) 0 as n . Therefore, the solution of model (2) exists; that is, it is constructed from the Equation (3).
Next, we show that every solution with the same initial value is unique. Suppose the model (2) has two different solutions with the same initial value, where the other solutions are denoted by S ^ , E ^ , I ^ , Q ^ , and R ^ . We obtain
S ( t ) S ^ ( t ) = 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 ( S , τ ) f 1 ( S ^ , τ ) d τ .
The norm of the above equation is
S ( t ) S ^ ( t ) 1 Γ ( α ) 0 t ( t τ ) α 1 f 1 ( S , τ ) f 1 ( S ^ , τ ) ,
so that by using the Lipschitz condition, we obtain
S ( t ) S ^ ( t ) p 1 T Γ ( α + 1 ) S ( t ) S ^ ( t ) 1 p j T Γ ( α + 1 ) S ( t ) S ^ ( t ) 0 .
The inequality (6) is satisfied if and only if S ( t ) S ^ ( t ) = 0 , which means S ( t ) = S ^ ( t ) . This leads to a contradiction, so that S ( t ) is unique. In a similar way, we can show that E ( t ) , I ( t ) , Q ( t ) , and R ( t ) are unique. □

4. Non-Negativity and Boundedness of the Solution

In this section, we show the non-negativity and boundedness of the solution of model (2). We also obtain the feasible region of model (2). The following theorem guarantees the non-negativity and boundedness of the solution.
Theorem 2. 
All solutions of model (2) subject to non-negative initial values are non-negative and ultimately bounded. The feasible region of model (2) is given by
Ω = ( S , E , I , Q , R ) R + 5 0 : N Λ μ .
Proof. 
We first show that S is non-negative. Assume S < 0 for some t > 0 . Then, there is at least one τ > 0 such that
S ( τ ) 0 , S ( τ ) = 0 , S ( τ + ) < 0 .
From model (2), we obtain
C D t α S ( t ) | t = τ = Λ > 0
According to the corollary of the non-decreasing function in [39], S > 0 in ( τ , τ + ε ) for arbitrary small positive constant ε . This leads to a contradiction. Hence, S 0 for all t > 0 . Similarly, we can prove E, I, Q, and R are non-negative. So, all solutions of model (2) are non-negative.
We next show the boundedness of the solution. From model (2), we have
C D t α N ( t ) = Λ μ N ( t ) δ I ( t ) κ Q ( t ) Λ μ N ( t ) .
According to comparison theorem in [39], we obtain
N ( t ) Λ μ + N ( 0 ) Λ μ E α ( μ t α ) .
Since the Mittag–Leffler expression E α ( μ t α ) 0 as t , it is clear that lim t N ( t ) Λ μ , and thus, N ( t ) is bounded with N ( t ) Λ μ . Hence, all solutions of model (2) are ultimately bounded. Then, we obtain the feasible region of model (2) as follows:
Ω = ( S , E , I , Q , R ) R + 5 0 : N Λ μ .

5. The Equilibrium Points and Basic Reproduction Number

To obtain the equilibrium point of model (2), we set
C D t α X = 0 ,
where X = ( S ,   E ,   I ,   Q ,   R ) T . Based on the third equation to fifth equation of (7), we obtain E = ( σ + η + δ + μ ) γ I ,   Q = σ ( ϑ + κ + μ ) I ,   R = η I + ϑ Q μ , and N = Λ δ I κ Q μ = Λ δ + κ σ η + κ + μ I μ . If E is substituted to the second equation of (7), then we have either
I = 0   o r   S = ( γ + μ ) ( σ + η + δ + μ ) β γ N .
For I = 0 , we have E = 0 , Q = 0 , R = 0 , N = Λ μ , and S = Λ μ , from which we obtain disease-free equilibrium point X 0 = ( Λ μ , 0 , 0 , 0 , 0 ) . From X 0 , we establish the basic reproduction number. The basic reproduction number, denoted R 0 , is the average number of secondary individuals infected by one primary infected individual in the susceptible class [40]. The basic reproduction number can be determined by the spectral radius of the next-generation matrix.
Let Z = ( E , I , Q ) T . From model (2), we have
C D t α Z = F ( Z ) V ( Z ) ,
where
F ( Z ) = β S I N 0 0 a n d V ( Z ) = ( γ + μ ) E γ E + ( σ + η + δ + μ ) I σ I + ( ϑ + κ + μ ) Q .
By calculating the spectral radius of next generation matrix, we obtain
R 0 = β γ ( γ + μ ) ( σ + η + δ + μ ) ;
see [25] for more details. Notice that R 0 is inversely proportional to the recovery rate of infected class η and quarantine rate σ . Therefore, if θ or σ increase, then R 0 decreases.
We next determine the second equilibrium point. For S = ( γ + μ ) ( σ + η + δ + μ ) β γ N = N R 0 , based on the first equation of (7), we obtain
I = Λ ( 1 R 0 ) δ + κ σ η + κ + μ β ,
from which we obtain the endemic equilibrium point ( S , E , I , Q , R ) with
I = Λ ( 1 R 0 ) δ + κ σ η + κ + μ β , S = Λ δ + κ σ η + κ + μ I μ R 0 , E = ( σ + η + δ + μ ) γ I , Q = σ ϑ + κ μ I , R = η I + ϑ Q μ ,
which exists if R 0 > 1 .
Theorem 3. 
The equilibrium points of model (2) are
  • X 0 = ( S 0 , 0 , 0 , 0 , 0 ) as the disease-free equilibrium with S 0 = Λ μ , which always exists;
  • X = ( S , E , I , Q , R ) as the endemic equilibrium, with S , E , I , Q , R given by Equations (8), which exists if R 0 > 1 .

6. The Stability of Equilibrium Points

In this section, we first study the local asymptotic stability of the equilibrium points of model (2). The following theorem provides the local stability of the equilibrium point.
Theorem 4. 
The disease-free equilibrium point X 0 is locally asymptotically stable if R 0 < 1 and is unstable if R 0 > 1 .
Proof. 
The Jacobian matrix of model (2) evaluated at X 0 is
J ( X 0 ) = μ 0 β 0 0 0 ( γ + μ ) β 0 0 0 γ ( σ + η + δ + μ ) 0 0 0 0 σ ( ϑ + κ + μ ) 0 0 0 η ϑ μ .
There are five eigenvalues ξ j ; that is, ξ 1 = ξ 2 = μ < 0 , ξ 3 = ( ϑ + κ + μ ) < 0 , and ξ 4 , ξ 5 are eigenvalues of matrix J S ( X 0 ) as follows:
J S ( X 0 ) = ( γ + μ ) β γ ( σ + η + δ + μ ) .
Notice that T r a c e ( J S ) = ( σ + η + δ + γ + 2 μ ) < 0 and d e t ( J S ) = ( γ + μ ) ( σ + η + δ + μ ) β γ . It is clear that d e t ( J S ) > 0 if R 0 < 1 . Thus, the real parts of both ξ 4 and ξ 5 are negative if R 0 < 1 . Moreover, if R 0 > 1 , then d e t ( J S ) < 0 . Thus, ξ 4 is positive and ξ 5 is negative. Hence, | arg ( ξ 4 ) | < α π 2 for all α ( 0 , 1 ) , which is the unstable condition of the equilibrium point. Therefore, X 0 is locally asymptotically stable if R 0 < 1 and is unstable if R 0 > 1 . □
We next study the global asymptotic stability of equilibrium points. The following theorems show the condition of the global stability of the disease-free equilibrium X 0 and the endemic equilibrium X .
Theorem 5. 
The disease-free equilibrium point X 0 is globally asymptotically stable if R 0 < 1 .
Proof. 
Assume that R 0 < 1 . Consider a Lyapunov function L, which is defined by
L ( E , I ) = E + γ + μ γ I
If we apply fractional derivative to L with respect to t [41], we obtain
C D t α L = C D t α E + γ + μ γ C D t α I = β S I N ( γ + μ ) E + γ + μ γ γ E ( σ + η + δ + μ ) I = β S I N ( γ + μ ) ( σ + η + δ + μ ) γ I β I ( γ + μ ) ( σ + η + δ + μ ) γ I = R 0 1 R 0 β I .
If R 0 < 1 , then C D t α L 0 . Moreover, the condition C D t α L = 0 is satisfied if and only if E = 0 and I = 0 . The Lemma of the LaSalle Invariance Principle in [42] guarantees E 0 and I 0 as t . Hence, Q 0 and R 0 as t , and thus, S Λ μ = S 0 . Therefore, X 0 is globally asymptotically stable. □
Theorem 6. 
Let the endemic equilibrium X exists, that is, R 0 > 1 . X is globally asymptotically stable in the domain Ω.
The proof of Theorem 6 uses same Lyapunov function W as in [25] by considering the inequality of the fractional derivative of the Lyapunov function [41]:
C D t α W 1 S S C D t α S + 1 E E C D t α E + γ + μ γ 1 I I C D t α I ,
where we simply change the equal sign in [25] to an inequality sign.

7. Sensitivity Analysis

In this section, we perform sensitivity analysis of basic reproduction number R 0 . This analysis aims to identify the role of each parameter of model (2) to R 0 , which is represented by the sensitivity index. The normalized forward sensitivity index of R 0 for parameter p ^ is formulated by
I p ^ R 0 = p ^ R 0 × R 0 p ^ .
The greater the absolute value of the index, the more influential the parameter is to R 0 . The sign of the index shows the enlargement and reduction of R 0 , that is, a positive sign enlarges R 0 while a negative sign reduces R 0 . We analytically calculate the sensitivity index for each parameter as follows:
I Λ R 0 = 0 , I β R 0 = 1 , I μ R 0 = μ ( σ + η + δ + μ ) + ( γ + μ ) ( γ + μ ) ( σ + η + δ + μ ) , I γ R 0 = μ ( γ + μ ) , I σ R 0 = σ ( σ + η + δ + μ ) , I η R 0 = η ( σ + η + δ + μ ) , I δ R 0 = δ ( σ + η + δ + μ ) , I ϑ R 0 = 0 , I κ R 0 = 0 .
To display the sensitivity index of R 0 for each parameter, we substitute the value of the parameter given in Table 2 into the index formula for each parameter. The sensitivity indices of R 0 are depicted in Figure 2, and the values are given in Table 3. It can be seen that parameter β , σ , and η have important role in reducing R 0 . However, the sensitivity of β is clear because the definition of R 0 is very close to disease transmission β . Thus, we discuss the sensitivity indices of both σ and η , i.e., 0.4585 and 0.4458 , respectively. Since both these indices are negative, the greater the values of σ or η , the smaller R 0 is. The results show both parameters σ and η have an important role in reducing the basic reproduction number and thus have an impact on the stability of the equilibrium point.

8. Numerical Simulation

In this section, we confirm our analytical results by setting the parameter values and implementing a predictor–corrector scheme in [37] to obtain the numerical solution of model (2). The values of the parameters in the first simulation are given by Table 2. Note that the parameter units shown in Table 2 are due to parameter scaling. Therefore, to maintain consistency with the physical dimensions, the units of all model parameters (2) are adjusted to the dimensions d a y α .
Since the quarantined class has better prevention than the infected one [45], we set the recovery rate to ϑ = 0.1 > η and the death rate to κ = 0.005 . Here, we use (38,000, 2000, 1000, 3000, 50, and 3000) as initial values in the model (2). For these parameters value, we obtain R 0 = 1.1881 > 1 , and thus two equilibrium points exist, namely, the unstable disease-free equilibrium X 0 = ( 220,410,536.4 ,   0 ,   0 ,   0 ,   0 ) and the globally asymptotically stable endemic equilibrium X = ( 184,244,619.17 ,   3,856,537.79 ,   2,198,361.42 ,   1,834,462.79 ,   26,769,323.88 ) . These stability properties are verified via our simulation depicted in Figure 3. In Figure 3, there is different behavior for different orders of fractional derivative ( α ) . It is clear that the smaller the order α , the slower the convergence of the solutions is. The peaks are also lower when α is smaller. Furthermore, for the numerical solution of R at α = 0.8 , there is no peak in the time interval [ 0 , 3500 ] . However, the curves of the solution for each α have the same behavior and preserve their convergence.
For the second simulation, we execute the model using the parameter values in Table 2, except the quarantine rate is σ = 0.2 . For the value of the parameters, we have R 0 = 0.81456 < 1 , and therefore, only one equilibrium exists, that is, a globally asymptotically stable disease-free equilibrium X 0 = ( 220,410,536.4 ,   0 ,   0 ,   0 , 0 ) . This stability equilibrium is confirmed by our numerical simulation, depicted in Figure 4; that is, the disease-free equilibrium X 0 is globally asymptotically stable for all orders of fractional derivative ( α ) . Again, the smaller the order of α , the slower the convergence of the solution. Therefore, endemic stability from the first simulation can be prevented by increasing the quarantine rate, which is actually supported by education and lockdown.
For the third simulation, the parameter values listed in Table 2 are used, except the recovery rate of infected individuals is η = 0.2 . We have R 0 = 0.8075 < 1 , and therefore, only one equilibrium point exists, that is, a globally asymptotically stable disease-free equilibrium X 0 = ( 220,410,536.4 ,   0 ,   0 ,   0 ,   0 ) . This global asymptotic stability of the disease-free equilibrium is confirmed by our numerical simulation, depicted in Figure 5. We observe that all numerical solutions are convergent to X 0 . Again, the order α affects the convergence rate of the solution. This simulation is the difference between our paper and those of Zeb et al. [21] and López and Rodó [24], both of whom did not consider the recovery rate of infected class η . We conclude that the recovery rate of the infected class is important to consider in the COVID-19 epidemic model because it can change the stability of the equilibrium point. Epidemiologically, the recovery rate can be increased by treatment based on Molnupiravir-assisted pharmacotherapy [46].

9. Parameter Estimation

In this section, we estimate the parameter values of model (2) using COVID-19 data from Indonesia. The COVID-19 data are adopted from the daily distribution map of confirmed active COVID-19 cases from the sixth wave, i.e., from 20 October 2022 to 10 January 2023, which is available on the official website of the Indonesian government at Distribution Map of COVID-19. Available online: https://covid19.go.id/peta-sebaran (accessed on 15 January 2023) [47]. The number of confirmed active daily cases corresponds to the infected subpopulations of model (2). The model’s parameter values are then estimated using the least square method, namely by minimizing the squared error between the numerical solution of the infected subpopulations and the data of confirmed active daily cases. To perform this process, we employ a built-in MATLAB function, that is, lsqcurvefit, which is based on the non-linear least square method. The initial value of model (2) is S ( 0 ) = 270 × 10 6 , E ( 0 ) = 50,000, I ( 0 ) = 18,919, Q ( 0 ) = 17,000, R ( 0 ) = 5 × 10 6 , M ( 0 ) = 150,000, where I ( 0 ) corresponds to first-day case of the data used. We also set the initial estimates for each parameter, that is, Λ = 270 × 10 6 70 × 365 , β = 10 , μ = 0.00001 , γ = 0.2 , σ = 0.4 , η = 0.9 , δ = 0.001 , ϑ = 0.9 , and κ = 0.001 . Each parameter value is estimated in the range from 0 to 1, except for recruitment rate Λ and disease transmission β . Ideally, recruitment rate Λ is around Λ = 270 × 10 6 70 × 365 = 10,567 , where the total population of Indonesia is around 270 × 10 6 and the average age of the Indonesian population is around 70 years. Hence, we set Λ in the range of 0 to 20,000 . We assume the estimated value of disease transmission β is in the range 0 to 50.
We perform parameter estimation for fractional-order α = 1 ,   0.97 ,   0.94 ,   0.9 ,   0.87 ,   0.84 ,   0.8 . Notice that estimated values of the parameter for α = 1 correspond to the parameter of the first-order model proposed by Darti et al. [25]. For each value of α , model calibration is performed for 63 days of data, i.e., from 20 October 2022, to 21 December 2022, whereas the additional data from 22 December 2022, to 10 January 2023, is used for 20-day forecasting. We compare the calibration and forecasting of both the fractional-order model and the first-order model. To perform this comparison, we determine the root mean square error (RMSE), expressed by
R M S E = k = 0 N I ( t k ) I d a t a ( t k ) 2 ,
where N is the number of data points, I ( t k ) is the numerical solution of the infected compartment of model (2) at t k , and I d a t a ( t k ) is the data of confirmed daily active cases at t k .
The estimated values of the parameters for different α are shown in Table 4. Notice that the estimated values of the parameter for each α are different. These cannot be compared between α because different alpha values have different units. The fitted curve for α = 1 , 0.97 , 0.93 , 0.9 , 0.87 , 0.83 , 0.8 is depicted in Figure 6. The confirmed daily active data used for calibration are marked by red dots, and the data for the 20-day forecasting are marked by blue dots. We discuss the results of calibration and forecasting from the viewpoint of the fitted curve and the obtained RMSE.
For the fitted curve of a first-order model, i.e., when α = 1 , it can be seen that there are large deviations almost everywhere except for around the inflection point and peak point. Fitted curves for the fractional-order model typically have some deviation at the beginning and end of the period. Except at the end of the period, the fitted curves for α = 0.97 and α = 0.94 depicted in Figure 6 are good. It should be noted that not all orders α are good for data fitting. For example, the fitted curve for α = 0.9 appears to have a large deviation from the data used. Thus, the value of α also influences the data fitting performance. Nonetheless, based on the fitted curves, a fractional-order curve is typically better than a first-order curve for calibration and 20-day forecasting. These results show that a fractional-order model has greater flexibility for data fitting, and thus the best-fit curve for the data can be provided.
To confirm the observation of the fitted curves for different values of α , we also determine the RMSE for the calibration and 20-day forecasting of confirmed daily active cases given by Table 5 and Table 6, respectively. It can be seen that the RMSE of calibration for the fractional-order model is smaller than the first-order model, except for α = 0.9 . The smallest RMSE of calibration occurs when fractional order α = 0.97 , where, in this case, the fractional-order model has a RMSE that is 2.7 times smaller than the first-order model. These results verify that a fractional-order model generally performs better than a first-order model in calibrating the COVID-19 data. Similarly, the RMSE of the 20-day forecasting for the fractional-order model is also smaller than the first-order model, except for α = 0.9 . Again, a fractional-order model generally provides better 20-day forecasts than a first-order model, especially the fractional-order model with α = 0.97 , which has an RMSE that is 2.3 times smaller than the first-order model. Finally, we can say a fractional-order model has better accuracy and flexibility than a first-order model when implemented using confirmed daily active cases of COVID-19 in Indonesia. This also confirms the statement of Alrabaiah et al. [29] that a fractional-order model is more flexible than a first-order model.

10. Conclusions

We construct a fractional-order COVID-19 epidemic model with quarantine and standard incidence rate using the Caputo fractional-order derivative. We modify earlier models by simultaneously considering the recovery rate and quarantine rate of infected individuals. We also consider the standard incidence rate in the model. We have proven the existence, uniqueness, non-negativity, and boundedness of the solution. The model has two equilibrium points, i.e., a disease-free equilibrium and an endemic equilibrium. We obtain the basic reproduction number ( R 0 ) by implementing the spectral radius of the next-generation matrix. The disease-free equilibrium always exists and is locally and globally asymptotically stable only if R 0 < 1 . In contrast, the endemic equilibrium exists and is globally asymptotically stable if R 0 > 1 . The stability properties of the equilibrium points are confirmed by numerical simulation. From the numerical simulation, we can see that the smaller the order of the derivative, the slower the convergence of the solution of the model and vice versa. In addition, we conclude that both the recovery rate and quarantine rate of the infected class have an important role in reducing R 0 , as well as for the stability of the equilibrium point. Based on parameter estimation results, a fractional-order model typically performs better than a first-order model for both calibrating and forecasting confirmed cases of COVID-19 in Indonesia. In our case, the fractional-order model with α = 0.97 has an RMSE that is twice as small as the first-order model for both the calibration and forecasting of confirmed COVID-19 cases.

Author Contributions

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

Funding

This research is funded by the Faculty of Mathematics and Natural Sciences (FMIPA) through Public Funds DPA (Dokumen Pelaksanaan Anggaran) Perguruan Tinggi Berbadan Hukum (PTNBH), University of Brawijaya, and based on FMIPA Professor Grant Contract No. 3084.19/UN10.F09/PN/2022.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data used in this study are available on the official website of the Indonesian government at https://covid19.go.id/peta-sebaran, accessed on 15 January 2023.

Acknowledgments

The authors thank the anonymous reviewers for their valuable suggestions and comments on our paper.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Riyapan, P.; Shuaib, S.E.; Intarasit, A. A mathematical model of COVID-19 pandemic: A case study of Bangkok, Thailand. Comput. Math. Methods Med. 2021, 2021, 6664483. [Google Scholar] [CrossRef] [PubMed]
  2. Ahmad, S.; Ullah, A.; Al-Mdallal, Q.M.; Khan, H.; Shah, K.; Khan, A. Fractional order mathematical modeling of COVID-19 transmission. Chaos Solitons Fractals 2020, 139, 110256. [Google Scholar] [CrossRef] [PubMed]
  3. Ahmad, S.W.; Sarwar, M.; Shah, K.; Ahmadian, A.; Salahshour, S. Fractional order mathematical modeling of novel corona virus (COVID-19). Math. Methods Appl. Sci. 2021, 46, 7847–7860. [Google Scholar] [CrossRef] [PubMed]
  4. World Health Organization. COVID-19 Weekly Epidemiological Update, Edition 115, 26 October 2022. 2022. Available online: https://apps.who.int/iris/bitstream/handle/10665/363853/nCoV-weekly-sitrep26Oct22-eng.pdf?sequence=2 (accessed on 5 September 2022).
  5. Burki, T.K. Omicron variant and booster COVID-19 vaccines. Lancet Respir. Med. 2022, 10, e17. [Google Scholar] [CrossRef]
  6. Araf, Y.; Akter, F.; Tang, Y.D.; Fatemi, R.; Parvez, M.S.A.; Zheng, C.; Hossain, M.G. Omicron variant of SARS-CoV-2: Genomics, transmissibility, and responses to current COVID-19 vaccines. J. Med. Virol. 2022, 94, 1825–1832. [Google Scholar] [CrossRef]
  7. Lotfi, M.; Hamblin, M.R.; Rezaei, N. COVID-19: Transmission, prevention, and potential therapeutic opportunities. Clin. Chim. Acta 2020, 508, 254–266. [Google Scholar] [CrossRef]
  8. Salian, V.S.; Wright, J.A.; Vedell, P.T.; Nair, S.; Li, C.; Kandimalla, M.; Tang, X.; Carmona Porquera, E.M.; Kalari, K.R.; Kandimalla, K.K. COVID-19 transmission, current treatment, and future therapeutic strategies. Mol. Pharm. 2021, 18, 754–771. [Google Scholar] [CrossRef]
  9. Tellier, R. COVID-19: The case for aerosol transmission. Interface Focus 2022, 12, 20210072. [Google Scholar] [CrossRef]
  10. Akindeinde, S.O.; Okyere, E.; Adewumi, A.O.; Lebelo, R.S.; Fabelurin, O.O.; Moore, S.E. Caputo fractional-order SEIRP model for COVID-19 Pandemic. Alex. Eng. J. 2022, 61, 829–845. [Google Scholar] [CrossRef]
  11. Tuan, N.H.; Mohammadi, H.; Rezapour, S. A mathematical model for COVID-19 transmission by using the Caputo fractional derivative. Chaos Solitons Fractals 2020, 140, 110107. [Google Scholar] [CrossRef]
  12. Bushnaq, S.; Khan, S.A.; Shah, K.; Zaman, G. Mathematical analysis of HIV/AIDS infection model with Caputo-Fabrizio fractional derivative. Cogent Math. Stat. 2018, 5, 1432521. [Google Scholar] [CrossRef]
  13. Khan, A.A.; Ullah, S.; Amin, R. Optimal control analysis of COVID-19 vaccine epidemic model: A case study. Eur. Phys. J. Plus 2022, 137, 1–25. [Google Scholar] [CrossRef] [PubMed]
  14. Khajji, B.; Kouidere, A.; Elhia, M.; Balatif, O.; Rachik, M. Fractional optimal control problem for an age-structured model of COVID-19 transmission. Chaos Solitons Fractals 2021, 143, 110625. [Google Scholar] [CrossRef] [PubMed]
  15. Musafir, R.R.; Suryanto, A.; Darti, I. Dynamics of COVID-19 epidemic model with asymptomatic infection, quarantine, protection and vaccination. Commun. Biomath. Sci. 2021, 4, 106–124. [Google Scholar] [CrossRef]
  16. Rayungsari, M.; Aufin, M.; Imamah, N. Parameters estimation of generalized Richards model for COVID-19 cases in indonesia using genetic algorithm. Jambura J. Biomath. (JJBM) 2020, 1, 25–30. [Google Scholar] [CrossRef]
  17. Darti, I.; Suryanto, A.; Panigoro, H.S.; Susanto, H. Forecasting COVID-19 epidemic in Spain and Italy using a generalized Richards model with quantified uncertainty. Commun. Biomath. Sci. 2020, 3, 90–100. [Google Scholar] [CrossRef]
  18. Musafir, R.R.; Anam, S. Parameter Estimation of COVID-19 Compartment Model in Indonesia Using Particle Swarm Optimization. J. Berk. Epidemiol. 2022, 10, 283–292. [Google Scholar] [CrossRef]
  19. Chu, Y.M.; Ali, A.; Khan, M.A.; Islam, S.; Ullah, S. Dynamics of fractional order COVID-19 model with a case study of Saudi Arabia. Results Phys. 2021, 21, 103787. [Google Scholar] [CrossRef]
  20. Paul, S.; Mahata, A.; Mukherjee, S.; Roy, B. Dynamics of SIQR epidemic model with fractional order derivative. Partial Differ. Equ. Appl. Math. 2022, 5, 100216. [Google Scholar] [CrossRef]
  21. Zeb, A.; Alzahrani, E.; Erturk, V.S.; Zaman, G. Mathematical model for coronavirus disease 2019 (COVID-19) containing isolation class. BioMed Res. Int. 2020, 2020, 3452402. [Google Scholar] [CrossRef]
  22. Postavaru, O.; Anton, S.; Toma, A. COVID-19 pandemic and chaos theory. Math. Comput. Simul. 2021, 181, 138–149. [Google Scholar] [CrossRef]
  23. Rangasamy, M.; Chesneau, C.; Martin-Barreiro, C.; Leiva, V. On a novel dynamics of SEIR epidemic models with a potential application to COVID-19. Symmetry 2022, 14, 1436. [Google Scholar] [CrossRef]
  24. López, L.; Rodó, X. A modified SEIR model to predict the COVID-19 outbreak in Spain and Italy: Simulating control scenarios and multi-scale epidemics. Results Phys. 2021, 21, 103746. [Google Scholar] [CrossRef]
  25. Darti, I.; Trisilowati; Rayungsari, M.; Musafir, R.R.; Suryanto, A. A SEIQRD Epidemic Model to Study the Dynamics of COVID-19 Disease. Commun. Math. Biol. Neurosci. 2023, 2023, 5. [Google Scholar] [CrossRef]
  26. Baleanu, D.; Abadi, M.H.; Jajarmi, A.; Vahid, K.Z.; Nieto, J. A new comparative study on the general fractional model of COVID-19 with isolation and quarantine effects. Alex. Eng. J. 2022, 61, 4779–4791. [Google Scholar] [CrossRef]
  27. Diethelm, K. The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Operators of Caputo Type; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2010. [Google Scholar]
  28. Farman, M.; Aslam, M.; Akgül, A.; Ahmad, A. Modeling of fractional-order COVID-19 epidemic model with quarantine and social distancing. Math. Methods Appl. Sci. 2021, 44, 9334–9350. [Google Scholar] [CrossRef] [PubMed]
  29. Alrabaiah, H.; Arfan, M.; Shah, K.; Mahariq, I.; Ullah, A. A comparative study of spreading of novel corona virus disease by ussing fractional order modified SEIR model. Alex. Eng. J. 2021, 60, 573–585. [Google Scholar] [CrossRef]
  30. Basti, B.; Hammami, N.; Berrabah, I.; Nouioua, F.; Djemiat, R.; Benhamidouche, N. Stability analysis and existence of solutions for a modified SIRD model of COVID-19 with fractional derivatives. Symmetry 2021, 13, 1431. [Google Scholar] [CrossRef]
  31. Podlubny, I. Fractional Differential Equations: Mathematics in Science and Engineering; Academic Press: Cambridge, MA, USA, 1999. [Google Scholar]
  32. Li, C.; Qian, D.; Chen, Y. On Riemann-Liouville and caputo derivatives. Discret. Dyn. Nat. Soc. 2011, 2011, 562494. [Google Scholar] [CrossRef] [Green Version]
  33. Yadan, S.; Kumar, D.; Singh, J.; Baleanu, D. Analysis and dynamics of fractional order COVID-19 model with memory effect. Results Phys. 2021, 24, 104017. [Google Scholar]
  34. Denu, D.; Kermausuor, S. Analysis of a Fractional-Order COVID-19 Epidemic Model with Lockdown. Vaccines 2022, 10, 1773. [Google Scholar] [CrossRef] [PubMed]
  35. Danane, J.; Hammouch, Z.; Allali, K.; Rashid, S.; Singh, J. A fractional-order model of coronavirus disease 2019 (COVID-19) with governmental action and individual reaction. Math. Methods Appl. Sci. 2023, 46, 8275–8288. [Google Scholar] [CrossRef]
  36. Oud, M.A.; Ali, A.; Alrabaiah, H.; Ullah, S.; Khan, M.A.; Islam, S. A fractional order mathematical model for COVID-19 dynamics with quarantine, isolation, and environmental viral load. Adv. Differ. Equ. 2021, 2021, 106. [Google Scholar] [CrossRef] [PubMed]
  37. Diethelm, K.; Ford, N.J.; Freed, A.D. A predictor-corrector approach for the numerical solution of fractional differential equations. Nonlinear Dyn. 2002, 29, 3–22. [Google Scholar] [CrossRef]
  38. Petráš, I. Fractional-Order Nonlinear Systems: Modeling, Analysis and Simulation; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2011. [Google Scholar]
  39. Li, H.L.; Zhang, L.; Hu, C.; Jiang, Y.L.; Teng, Z. Dynamical analysis of a fractional-order predator-prey model incorporating a prey refuge. J. Appl. Math. Comput. 2017, 54, 435–449. [Google Scholar] [CrossRef]
  40. Murray, J.D. Mathematical Biology I. An Introduction, 2nd ed.; Interdisciplinary Applied Mathematics; Springer Science & Business Media: Berlin/Heidelberg, Germany, 2002; Volume 17. [Google Scholar] [CrossRef]
  41. Vargas-De-León, C. Volterra-type Lyapunov functions for fractional-order epidemic systems. Commun. Nonlinear Sci. Numer. Simul. 2015, 24, 75–85. [Google Scholar] [CrossRef]
  42. Huo, J.; Zhao, H.; Zhu, L. The effect of vaccines on backward bifurcation in a fractional order HIV model. Nonlinear Anal. Real World Appl. 2015, 26, 289–305. [Google Scholar] [CrossRef]
  43. Ahmad, Z.; Arif, M.; Ali, F.; Khan, I.; Nisar, K.S. A report on COVID-19 epidemic in Pakistan using SEIR fractional model. Sci. Rep. 2020, 10, 22268. [Google Scholar] [CrossRef]
  44. Carcione, J.M.; Santos, J.E.; Bagaini, C.; Ba, J. A simulation of a COVID-19 epidemic based on a deterministic SEIR model. Front. Public Health 2020, 8, 230. [Google Scholar] [CrossRef]
  45. Marshall, P. The impact of quarantine on COVID-19 infections. Epidemiol. Methods 2021, 10, 1–8. [Google Scholar] [CrossRef]
  46. Pourkarim, F.; Pourtaghi-Anvarian, S.; Rezaee, H. Molnupiravir: A new candidate for COVID-19 treatment. Pharmacol. Res. Perspect. 2022, 10, e00909. [Google Scholar] [CrossRef] [PubMed]
  47. Distribution Map of COVID-19. Available online: https://covid19.go.id/peta-sebaran (accessed on 15 January 2023).
Figure 1. Compartment diagram of the proposed model.
Figure 1. Compartment diagram of the proposed model.
Axioms 12 00591 g001
Figure 2. Bar plot of sensitivity indices of R 0 for each parameter.
Figure 2. Bar plot of sensitivity indices of R 0 for each parameter.
Axioms 12 00591 g002
Figure 3. Numerical solution of model (2) for different fractional derivative orders ( α ) using the parameter values in Table 2.
Figure 3. Numerical solution of model (2) for different fractional derivative orders ( α ) using the parameter values in Table 2.
Axioms 12 00591 g003
Figure 4. Numerical solution of model (2) for different fractional derivative order ( α ) with σ = 0.2 and other parameter values in Table 2.
Figure 4. Numerical solution of model (2) for different fractional derivative order ( α ) with σ = 0.2 and other parameter values in Table 2.
Axioms 12 00591 g004
Figure 5. Numerical solution of model (2) for different fractional derivative orders ( α ) with η = 0.2 and other parameter values from Table 2.
Figure 5. Numerical solution of model (2) for different fractional derivative orders ( α ) with η = 0.2 and other parameter values from Table 2.
Axioms 12 00591 g005
Figure 6. Fitted curve of numerical solution of model (2) to confirmed COVID-19 data in Indonesia for different values of α . Part (a) shows the fitted curves for derivative order α = 1 , 0.97 , 0.94 , 0.9 . Part (b) shows the fitted curves for derivative order α = 1 , 0.87 , 0.84 , 0.8 . Notice that the fitted curve for α = 1 corresponds to a first-order model.
Figure 6. Fitted curve of numerical solution of model (2) to confirmed COVID-19 data in Indonesia for different values of α . Part (a) shows the fitted curves for derivative order α = 1 , 0.97 , 0.94 , 0.9 . Part (b) shows the fitted curves for derivative order α = 1 , 0.87 , 0.84 , 0.8 . Notice that the fitted curve for α = 1 corresponds to a first-order model.
Axioms 12 00591 g006
Table 1. Description of variables and parameters.
Table 1. Description of variables and parameters.
ParameterDescription
S ( t ) Susceptible class
E ( t ) Exposed class
I ( t ) Infected class
Q ( t ) Quarantined class
R ( t ) Recovered class
M ( t ) Death class
N ( t ) Total of living class
Λ ^ Recruitment rate
β ^ Rate of disease transmission
μ ^ Natural death rate
γ ^ Incubation rate
σ ^ Quarantine rate
η ^ Recovery rate of infected class
δ ^ Mortality rate of infected class
ϑ ^ Recovery rate of quarantined class
κ ^ Mortality rate of quarantined class
Table 2. The parameter values of the sensitivity index.
Table 2. The parameter values of the sensitivity index.
ParameterUnitValueSource
Λ i n d i v i d u a l d a y α 3,270,186 [43]
β 1 d a y α 0.29 [43]
μ 1 d a y α 1 67.4 [43]
γ 1 d a y α 0.1243 [43]
σ 1 d a y α 0.1 Assumed
η 1 d a y α 0.09722 [43]
δ 1 d a y α 0.006 [44]
ϑ 1 d a y α 0.1 Assumed
κ 1 d a y α 0.005 Assumed
Table 3. The value of sensitivity index of R 0 with respect to parameter values in Table 2.
Table 3. The value of sensitivity index of R 0 with respect to parameter values in Table 2.
ParameterValue of Sensitivity Index
Λ 0
β 1
μ 0.1746
γ 0.1066
σ 0.4585
η 0.4458
δ 0.0275
ϑ 0
κ 0
Table 4. Estimated parameter value for different values of derivative order.
Table 4. Estimated parameter value for different values of derivative order.
ParameterDerivative Order
α = 1 α = 0 . 97 α = 0 . 94 α = 0 . 90 α = 0 . 87 α = 0 . 84 α = 0 . 80
Λ 20 , 000 1 , 853.4473 0.1936 5 , 259.6532 12.6850 446.9792 11.1194
β 3.2961 1.8559 2.9466 9.5676 4.0059 4.9821 5.7791
μ 0.0412 0.2193 0.3010 0.3459 0.2764 0.3126 0.3188
γ 1 0.1081 0.1247 0.0943 0.3604 0.3966 0.5526
σ 1 0.2452 × 10 4 0.0184 0.1745 0.5946 0.7644 0.9999
η 1 0.0930 0.1304 0.9224 0.5185 0.5876 0.7543
δ 0.8812 0.1312 0.1891 0.0042 0.4550 0.6031 0.9764
ϑ 0.7624 0.0372 0.1012 0.8045 0.0655 0.0005 0.0004
κ 2.751 × 10 14 0.898 × 10 6 0.3903 0.1947 0.9111 0.9749 0.9998
Table 5. Root mean square error of the calibration of confirmed COVID-19 cases.
Table 5. Root mean square error of the calibration of confirmed COVID-19 cases.
Derivative Order
α = 1 α = 0.97 α = 0.94 α = 0.90 α = 0.87 α = 0.84 α = 0.80
2.5880 × 10 4 0.9701 × 10 4 1.5819 × 10 4 3.4915 × 10 4 1.3084 × 10 4 1.4788 × 10 4 1.6319 × 10 4
Table 6. Root mean square error of 20-day forecasting of confirmed COVID-19 cases.
Table 6. Root mean square error of 20-day forecasting of confirmed COVID-19 cases.
Derivative Order
α = 1 α = 0.97 α = 0.94 α = 0.90 α = 0.87 α = 0.84 α = 0.80
4.1023 × 10 4 1.7944 × 10 4 2.7407 × 10 4 7.4194 × 10 4 1.8617 × 10 4 2.2247 × 10 4 2.7200 × 10 4
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

Trisilowati; Darti, I.; Musafir, R.R.; Rayungsari, M.; Suryanto, A. Dynamics of a Fractional-Order COVID-19 Epidemic Model with Quarantine and Standard Incidence Rate. Axioms 2023, 12, 591. https://doi.org/10.3390/axioms12060591

AMA Style

Trisilowati, Darti I, Musafir RR, Rayungsari M, Suryanto A. Dynamics of a Fractional-Order COVID-19 Epidemic Model with Quarantine and Standard Incidence Rate. Axioms. 2023; 12(6):591. https://doi.org/10.3390/axioms12060591

Chicago/Turabian Style

Trisilowati, Isnani Darti, Raqqasyi Rahmatullah Musafir, Maya Rayungsari, and Agus Suryanto. 2023. "Dynamics of a Fractional-Order COVID-19 Epidemic Model with Quarantine and Standard Incidence Rate" Axioms 12, no. 6: 591. https://doi.org/10.3390/axioms12060591

APA Style

Trisilowati, Darti, I., Musafir, R. R., Rayungsari, M., & Suryanto, A. (2023). Dynamics of a Fractional-Order COVID-19 Epidemic Model with Quarantine and Standard Incidence Rate. Axioms, 12(6), 591. https://doi.org/10.3390/axioms12060591

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