Fractional-order SEIQRDP model for simulating the dynamics of COVID-19 epidemic
Mohamed Bahloul, Abderrazak Chahid, Taous Meriem Laleg-Kirati
FF RACTIONAL - ORDER
SEIQRDP
MODEL FOR SIMULATING THEDYNAMICS OF
COVID-19
EPIDEMIC
Mohamed A. Bahloul
Computer, Electrical and Mathematical Sciences and Engineering Division (CEMSE)King Abdullah University of Science and Technology (KAUST)Thuwal 23955-6900, Makkah Province, Saudi Arabia [email protected]
Abderrazak Chahid
Computer, Electrical and Mathematical Sciences and Engineering Division (CEMSE)King Abdullah University of Science and Technology (KAUST)Thuwal 23955-6900, Makkah Province, Saudi Arabia [email protected]
Taous-Meriem Laleg-Kirati
Computer, Electrical and Mathematical Sciences and Engineering Division (CEMSE)King Abdullah University of Science and Technology (KAUST)Thuwal 23955-6900, Makkah Province, Saudi Arabia [email protected]
June 9, 2020 A BSTRACT
The novel corona-virus disease (COVID-19), known as the causative virus of outbreak pneumoniainitially recognized in the mainland of China, late December 2019. COVID-19 reaches out tomany countries in the world, and the number of daily cases continues to increase rapidly. Inorder to simulate, track, and forecast the trend of the virus spread, several mathematical andstatistical models have been developed.
Susceptible-Exposed-Infected-Quarantined-Recovered-Death-Insusceptible (SEIQRDP) model is one of the most promising dynamic systems thathas been proposed for estimating the transmissibility of the COVID-19. In the present study,we propose a Fractional-order SEIQRDP model to analyze the COVID-19 epidemic. TheFractional-order paradigm offers a flexible, appropriate, and reliable framework for pandemicgrowth characterization. In fact, fractional-order operator is not local and consider the memory ofthe variables. Hence, it takes into account the sub-diffusion process of confirmed and recoveredcases growth. The results of the validation of the model using real COVID-19 data are presented,and the pertinence of the proposed model to analyze, understand and predict the epidemic is discussed.
COVID-19 is a respiratory disease caused by the new coronavirus that was first identified in Wuhan, Hubei province,China, late December 2019 [1]. This novel virus soon began to spread out around the world, and on 30 January, WHOdeclared the outbreak a Public Health Emergency of International Concern (PHEIC). On 11 March, WHO Director-General marked COVID-19 as a pandemic [2, 3]. The transmission of COVID-19 is primarily through respiratorydroplets and contact routes [4]. Accordingly, ideal interventions to control the spread include: quarantine, isolation,increase home confinement, promoting the wearing of face masks, travel restrictions, the closing of public space, and a r X i v : . [ q - b i o . P E ] J un PREPRINT - J
UNE
9, 2020cancellation of events. The number of cases increased rapidly to more than 3.25 million cases, including around 231,000deaths worldwide as of April, , 2020.The novel corona-virus caused-pneumonia has attracted the attention of scientists with different backgrounds rangedfrom epidemiology to data science, mathematics, and statistics. Accordingly, several studies based on either statistics ormathematical modeling have been proposed for better analysis and a deep understanding of the evolution of this epidemic.Besides, enormous efforts have been devoted to predicting the inflection point and ending time of this epidemic in orderto help make decisions concerning the different measures that have been taken by different governments.Mathematical models have always played a crucial role in the understanding of the spread of the virus and in providingrelevant guidelines for controlling the pandemic. Basically, mathematical paradigms are considered to be very useful inthis context, as they provide detailed mechanisms for the epidemic dynamics. Among the most widely investigatedmodel for the characterization of the COVID-19 outbreak in the world is the classic Susceptible-Exposed-Infectious-Recovered (SEIR) model [1]. Since the outbreak of the virus, the SEIR model has been intensively utilized to evaluatethe effectiveness of multiple measures, which seems to be a challenging task for general other estimation methods [5–7].For instance, it has been employed to evaluate the effects of lock-down on the transmission dynamics between provincesin China, such as the effect of the lock-down in Hubei province on the transmission in Wuhan and Beijing [7]. Inaddition, a cascading scheme of the SEIR model has been studied to emulate the process of transmission from infectionsources to humans. This approach was efficient to reach useful conclusions on the outbreak dynamics [8]. The workof [9] presented a generalization of the classical SEIR known as
Susceptible-Exposed-Infected-Quarantined-Recovered-Death-Insusceptible (SEIQRDP) model, for the epidemic analysis of COVID-19 in China. This generalization is basedon the introduction of a new quarantined state and of the effect of preventive actions, which are considered as crucialepidemic parameters for COVID-19, such as the latent and quarantine time. As a result, considering this generalizationof SEIR, the estimation of the inflection point, ending time, and total infected cases in widely affected regions wasaccurately determined and verified.Over the last decade, the fractional-order derivative (FD), defined as a generalization of the conventional integerderivative to a non-integer order (arbitrary order) operator, has been used to simulate many phenomena involvingmemory and delays including epidemic behavior [10, 11]. FD models offer a promising tool for the description ofcomplex systems, in addition to their potential to incorporate accurately the memory and delay involved in the systems,it also provides more flexibility than classical integer-order models in fitting the data accurately [12].In this paper, weinvestigate a fractional version of the (
SEIQRDP ) model in modeling the COVID-19 epidemic’s dynamics. We will usereal data to fit the model and analyze the results, and we will provide some insights on the interpretation and role of thefractional derivatives [11].This paper is organized as follows. In Section II, we will recall some basic concepts from the
SEIQRDP epidemicmodel and the fractional-order derivatives. Section III is devoted to the presentation of fractional-order (
SEIQRDP )epidemic model (
F-SEIQRDP ). In Section IV, we present the materials and methods. Section v present the estimationresults. The last section discusses the obtained results and provides some future directions on the use of the model foranalyzing and controlling the COVID-19 epidemic.
In this section, we recall some basic concepts from the
SEIQRDP epidemic model and fractional-order derivativestheory.
As described in [9] the
SEIQRDP model consists of seven sub-populations (states), i.e { S ( t ) , E ( t ) , I ( t ) , Q ( t ) , R ( t ) , D ( t ) , P ( t ) } denoting at time t the followings: • S ( t ) : Susceptible cases. • E ( t ) : Exposed cases, which are infected but not yet be infectious, in a latent period. • I ( t ) : Infectious cases, which have infectious capacity but not yet be quarantined. • Q ( t ) : Quarantined cases, which are confirmed and infected. • R ( t ) : Recovered cases. • D ( t ) : Dead cases. • P ( t ) : Insusceptible cases. 2 PREPRINT - J
UNE
9, 2020It contains also six parameters defined as follows: • α : the protection rate. • β : the infection rate. • γ : the inverse of the average latent time. • δ : the rate at which infectious people enter in quarantine. • λ ( t ) : a time-dependant coefficient used in the description of the cure rate. It is expressed as: λ ( t ) = λ (cid:2) − e − λ t (cid:3) , (1)where λ and λ are empirical coefficients. • κ ( t ) : time-dependant coefficient used in the description of the mortality rate. It is expressed as: κ ( t ) = κ e − κ t , (2)where κ and κ are empirical coefficients.It is worth to note that the time-dependent expressions of the cure rate, λ ( t ) , and the mortality rate, κ ( t ) , are assumed tobe in the above forms based on the analysis of real data collected in some provinces in China, in January 2020, [13].The plot and analysis of this data showed the gradually increase of the cure rate and the quick decrease of the mortalityrate. Furthermore, this assumptions are very reasonable by nature as the function of death rate in such pandemic alwaysconverges to zero while the cure rate continues increasing toward a consistent level. The other parameters are assumedto be constant as they are not fluctuating over time.The above parameters are controlled by the application of the preventive interventions as well as the effectiveness of thehealth systems in the investigated region. Figure 1 illustrates the relations between all the states. The dynamic of eachstate is mathematically characterized by ordinary differential equations (ODE) as follows:Figure 1: Flowchart illustrates the COVID-19 based SEIQRDP epidemic model, adopted from Fig.1 in [9].3 PREPRINT - J
UNE
9, 2020 d S ( t )d t = − αS ( t ) − β S ( t ) I ( t ) N d E ( t )d t = − γE ( t ) + β S ( t ) I ( t ) N d I ( t )d t = γE ( t ) − δI ( t )d Q ( t )d t = δI ( t ) − λ ( t ) Q ( t ) − κ ( t ) Q ( t )d R ( t )d t = λ ( t ) Q ( t )d D ( t )d t = κ ( t ) Q ( t )d P ( t )d t = αS ( t ) (3)where N represents the total population in the studied region expressed as N = S + E + I + Q + R + D + P .Comparing to the classical SEIR model,
SEIQRDP is augmented by three new states, {
Q(t) , D(t) & P(t) }. This newquarantined state
Q(t) and the recovery state
R(t) constitute, originally, the recovery state of the classical
SEIR model.
In the past few decades, the theory of fractional calculus (FC) has gained significant research attention in several fieldssuch as biology and epidemic modeling [14–16]. This is originated from the interdisciplinary nature of this field as wellas the flexibility and effectiveness of FC in describing complex physical systems. For example, the characterizationof bio-impedance, modeling of the viscoelasticity and biological cells, and representing the mechanical propertiesof the arterial system, as well as respiratory systems, have been investigated extensively through the exploring ofFC [12, 17–19]. The concept of FC is not new dating from the pioneer conversation between
L’Hopital and
Leibniz in that yielded to the generalization of the conventional integer derivative to a non-integer order operator [20], asfollows: D qt = d q d t α if q > , if q = 0 , (cid:82) t ( df ) − q if q < (4)where q ∈ R is the order of the operator known as the fractional-order, and df is the derivative function.There are several fractional-order derivative definitions. In this work, we introduce the three most frequently usedones in the sense of the Riemann–Liouville , Caputo and
Grünwald–Letnikov
FD-based definitions [21–23]. The
Grünwald–Letnikov scheme based on finite differences has been adopted in the numerical implementation of theproposed
F-GESIR .For a function g ( t ) that satisfies some smoothness conditions then: • The
Riemann–Liouville definition is given as:
RLa D qt g ( t ) = 1Γ ( n − q ) d n d t n (cid:90) t (1 − τ ) − q − n g ( τ ) dτ. (5) • The
Caputo definition for FD is expressed as follows: Ca D qt g ( t ) = 1Γ ( n − q ) (cid:90) t (1 − τ ) − q − n d n d t n g ( τ ) dτ, (6)where Γ is the Euler gamma function and ( n − < q < n ).4 PREPRINT - J
UNE
9, 2020 • The GL definition is given as:
GLa D qt g ( t ) = lim h → h − q [ t − ah ] (cid:88) j =0 ( − j (cid:18) qj (cid:19) g ( t − jh ) , (7)where a is the terminal point and [.] means the integer part. Similar to the
SEIQRDP , [13], the
F-SEIQRDP epidemiological model considers that the total population ( N ) isdivided into seven sub-populations i.e { S ( t ) , E ( t ) , I ( t ) , Q ( t ) , R ( t ) , D ( t ) , P ( t ) } . As the fractional-order derivativetakes into account the history of the state, we believe that this operator is more suitable to describe the dynamics ofthe epidemic COVID-19. Using the definition of the fractional-order derivative operator, we consider that each statefollows a fractional-order behavior.Considering the nonlinear FODEs in this matrix form: D qt X ( t ) = AX ( t ) + L ( X ) , (8)where, D q = [ D q S , D q E , D q I , D q Q , D q R , D q D , D q P ] T is the fractional-order derivative operator for all the states and X = [ S, E, I, Q, R, D, P ] T represents the state vector. A = − α q S − γ q E γ q I − δ q I δ q Q − κ q Q ( t ) − λ q Q ( t ) 0 0 00 0 0 λ q R ( t ) 0 0 00 0 0 κ q D ( t ) 0 0 0 α q P represent the parameters. L ( X ) depicts the nonlinear term that is function of the susceptible and L = S ( t ) I ( t ) − β q S N − β q E N The updated epidemic data of different countries around the world is collected from authoritative and known sources asfollows: • France : The data is gathered from three main sources: "Agence Regionale de Sante", "Santé Publique France"and "Geodes". This data is publicly available. https://github.com/cedricguadalupe/FRANCE-COVID-19 PREPRINT - J
UNE
9, 2020 • Italy : This data is provided by the Italian government and it is publicly available. • Other countries : The data is gathered from different official sources: World Health Organization (WHO),Center of Disease Control and Prevention (CDC), the COVID Tracking Project (testing and hospitalizations),etc. The data repository is operated by the Johns Hopkins University Center for Systems Science andEngineering (JHU CSSE) and Supported by ESRI Living Atlas Team and the Johns Hopkins UniversityApplied Physics Lab (JHU APL). The repository is publicly available. The parameters of the proposed fractional-order model were estimated by a non-linear least square minimization routine,making use of the well-known
MATLAB − R2019b , function lsqnonlin . This function is based on the trust-regionreflective method [24]. The steps used to obtain the optimal estimates are outlined in Algorithm 1.
Algorithm 2
Parameter estimation of epidemic data
Input : t : Time in days R : Recovered cases I : Confirmed cases D : Dead cases guess : The initial guess of the parameters f unmodel : The model to be fitted Output: param : Fitted parameter of the funmodel - Set the initial conditions
E = I; (cid:46)
Unknown but unlikely to be zero.Q = I-R-D;input= [E; I; Q; R; D] - Run the fitting optimization param=lsqcurvefit (t, funmodel , guess , input)The fitting performances are evaluated using the followings metrics: RM SE = (cid:114) l Σ li =1 (cid:16) y ( i ) − ˆ y ( i ) (cid:17) , (9)and ReM SE = l Σ li =1 | y ( i ) − ˆ y ( i ) | max ( y ) , (10)where y and ˆ y are the real and fitted data, respectively. l is the length of the data. The fitting performance of predicting the dynamics of Q ( T ) , R ( t ) , and D ( t ) populations using SEIQRDP and theproposed
F-SEIQRDP are presented in Table I. It is worth to note that
SEIQRDP epidemiological model can beconsidered as a special case of the
F-SEIQRDP , where all the fractional differentiation orders are equal to 1. From thereported results, it is clear that for all the studied cities the
ReRMSE as well as the
RMSE based on
F-SEIQRDP modelare less than the ones reported using
SEIQRDP . These results show the usefulness of the fractional-order derivativeoperator in fitting real data of the pandemic. In addition, it demonstrates the potential of the fractional-order frameworkin estimating the size and the key milestones of the spread of the epidemic-COVID-19. The appropriateness of thefractional-order paradigm can be visualized from the fact that: FD operator is not local and depends on the strength ofthe memory that is controlled by the fractional differentiation order. On the other hand, the epidemiological dynamicalprocess is involving the memory effect within the sub-diffusion process of confirmed and recovered cases growth.The parameter estimates for all the studied populations in different cities are reported in Table II. In this study, we choosedifferent cities that present different circumstances in terms of the total number of population, the number of infected https://github.com/pcm-dpc/COVID-19 https://github.com/CSSEGISandData/COVID-19 PREPRINT - J
UNE
9, 2020Table 1: Estimation error comparison between the GESIR model and the proposed F-GESIR model for differentcountries.
Model EstimationError China Italy France
Beijing Guangdong Henan Hubei Lombardia Veneto Emilia-Romagna Piemonte Nouvelle-Aquitaine
F-SEIQRDP
ReMSE
RMSE
SEIQRDP
ReMSE
RMSE cases, and the lock-down schemes. Figures 2,3 and 4 show examples of the predicted dynamic of the quarantined,recovered and death sub-populations in {Beijing, Guangdong, Henan, and Hubei} cities in China, {Lombardia, Veneto,Figure 2: Predictions of the proposed fractional model using data from China.7
PREPRINT - J
UNE
9, 2020Figure 3: Predictions of the proposed fractional model using data from Italy.Figure 4: Predictions of the proposed fractional model using data from France.8
PREPRINT - J
UNE
9, 2020Emilia-Romagna, and Piemonte} cities in Italy, and {Nouvelle Aquitaine} region in France, respectively. It is apparentin all the figures that we present the trend of the epidemiological dynamic till May th , which is the future concerningthe date of simulation April, th . This shows the potential of the model in predicting the trend of the pandemicdynamic in the future. Mathematical models are considered vital in every stage of the epidemic evolution. In fact, the simulation of theepidemic dynamics helps to track and monitor the spread of the virus. In addition, models are very effective in estimatingthe size of the pandemics and hence they might assist the specialized health actors to make the right decisions andTable 2: Epidemic spread parameters estimation using the proposed F-SEIQRDP model for different countries.
Country China Italy FranceCity
Beijing Guangdong Henan Hubei Lombardia Veneto Emilia-Romagna Piemonte NouvelleAquitaine
Population rate(million) 21.54 113.46 94 58.50 10.04 4.90 4.45 4.37 5.98Protection rate ( α ) ( β ) ( γ − ) ( δ − ) λ λ κ κ q S q E q I q Q q R q D q P PREPRINT - J
UNE
9, 2020minimize the losses. This paper proposes a new general fractional-order model for the evolution of the COVID-19pandemic. The first validation results show the accurate model fitting using real COVID-19 data from different countries.The fractional-order derivatives provide new and pertinent parameters for the control of the epidemic. However, themodel has some limitations that we can list as follows. • We observe that for countries with less data (reported cases in less than 30 days), the estimation is less accuratebecause the trend does not appear yet (more than 28 days). • From country to country, the initial guess of parameters should be chosen carefully to guarantee the bestpossible fitting of the used optimization solver. • The time-variable parameters need to be considered carefully because the countries do not have similar medicalfacilities and expertise, and they perform different numbers of tests per day. Besides, some countries haveadopted some precautions strategies earlier than others such as quarantine, lock-down...etc.
Acknowledgment
Research reported in this publication was supported by King Abdullah University of Science and Technology (KAUST).
Funding
This research project has been funded by King Abdullah University of Science and Technology (KAUST) Base ResearchFund (BAS/1/1627-01-01).
References [1] Chaolin Huang, Yeming Wang, Xingwang Li, Lili Ren, Jianping Zhao, Yi Hu, Li Zhang, Guohui Fan, JiuyangXu, Xiaoying Gu, et al. Clinical features of patients infected with 2019 novel coronavirus in wuhan, china.
TheLancet , 395(10223):497–506, 2020.[2] Juliet Bedford, Delia Enria, Johan Giesecke, David L Heymann, Chikwe Ihekweazu, Gary Kobinger, H CliffordLane, Ziad Memish, Myoung-don Oh, Anne Schuchat, et al. Covid-19: towards controlling of a pandemic.
TheLancet , 395(10229):1015–1018, 2020.[3] Yan-Rong Guo, Qing-Dong Cao, Zhong-Si Hong, Yuan-Yang Tan, Shou-Deng Chen, Hong-Jun Jin, Kai-SenTan, De-Yun Wang, and Yan Yan. The origin, transmission and clinical therapies on coronavirus disease 2019(covid-19) outbreak–an update on the status.
Military Medical Research , 7(1):1–10, 2020.[4] J Liu, X Liao, S Qian, J Yuan, F Wang, Y Liu, Z Wang, FS Wang, L Liu, and Z Zhang. Community transmissionof severe acute respiratory syndrome coronavirus 2, shenzhen, china, 2020.
Emerging Infectious Diseases , 26(6),2020.[5] Biao Tang, Nicola Luigi Bragazzi, Qian Li, Sanyi Tang, Yanni Xiao, and Jianhong Wu. An updated estimation ofthe risk of transmission of the novel coronavirus (2019-ncov).
Infectious disease modelling , 5:248–255, 2020.[6] Biao Tang, Xia Wang, Qian Li, Nicola Luigi Bragazzi, Sanyi Tang, Yanni Xiao, and Jianhong Wu. Estimationof the transmission risk of the 2019-ncov and its implication for public health interventions.
Journal of ClinicalMedicine , 9(2):462, 2020.[7] Mingwang Shen, Zhihang Peng, Yuming Guo, Yanni Xiao, and Lei Zhang. Lockdown may partially halt thespread of 2019 novel coronavirus in hubei province, china. medRxiv , 2020.[8] Tian-Mu Chen, Jia Rui, Qiu-Peng Wang, Ze-Yu Zhao, Jing-An Cui, and Ling Yin. A mathematical model forsimulating the phase-based transmissibility of a novel coronavirus.
Infectious diseases of poverty , 9(1):1–8, 2020.[9] Liangrong Peng, Wuyue Yang, Dongyan Zhang, Changjing Zhuge, and Liu Hong. Epidemic analysis of covid-19in china by dynamical modeling. arXiv preprint arXiv:2002.06563 , 2020.[10] Gilberto González-Parra, Abraham J Arenas, and Benito M Chen-Charpentier. A fractional order epidemicmodel for the simulation of outbreaks of influenza a (h1n1).
Mathematical methods in the Applied Sciences ,37(15):2218–2226, 2014.[11] Amjad S Shaikh, Iqbal N Shaikh, and Kottakkaran Sooppy Nisar. A mathematical model of covid-19 usingfractional derivative: Outbreak in india with dynamics of transmission and control. 2020.[12] Richard L Magin.
Fractional calculus in bioengineering , volume 2. Begell House Redding, 2006.10
PREPRINT - J
UNE
Retrieved April , 6:2020, 2020.[14] E Ahmed and AS Elgazzar. On fractional order differential equations model for nonlocal epidemics.
Physica A:Statistical Mechanics and its Applications , 379(2):607–614, 2007.[15] Ivan Area, Hanan Batarfi, Jorge Losada, Juan J Nieto, Wafa Shammakh, and Ángela Torres. On a fractional orderebola epidemic model.
Advances in Difference Equations , 2015(1):278, 2015.[16] Md Rafiul Islam, Angela Peace, Daniel Medina, and Tamer Oraby. Integer versus fractional order seir deterministicand stochastic models of measles.
International Journal of Environmental Research and Public Health , 17(6):2014,2020.[17] Mohamed A Bahloul and Taous Meriem Laleg-Kirati. Arterial viscoelastic model using lumped parameter circuitwith fractional-order capacitor. In , pages 53–56. IEEE, 2018.[18] Rudolf Hilfer et al.
Applications of fractional calculus in physics , volume 35. World scientific Singapore, 2000.[19] Mohamed A Bahloul and Taous-Meriem Laleg-Kirati. Assessment of fractional-order arterial windkessel as amodel of aortic input impedance.
IEEE Open Journal of Engineering in Medicine and Biology , 2020.[20] Igor Podlubny.
Fractional differential equations: an introduction to fractional derivatives, fractional differentialequations, to methods of their solution and some of their applications , volume 198. Elsevier, 1998.[21] Ivo Petrás. Fractional derivatives, fractional integrals, and fractional differential equations in matlab. In
Engineer-ing education and research using MATLAB . IntechOpen, 2011.[22] Igor Podlubny. An introduction to fractional derivatives, fractional differential equations, to methods of theirsolution and some of their applications.
Mathematics in science and engineering , 198:xxiv+–340, 1999.[23] KS Miller and B Ross. An introduction to the fractional calculus and fractional differential equations: John willeyand sons inc.
New York, NY , 1993.[24] Thomas F Coleman and Yuying Li. An interior trust region approach for nonlinear minimization subject to bounds.