Dynamics of fractional order COVID-19 model with a case study of Saudi Arabia

The novel coronavirus disease or COVID-19 is still posing an alarming situation around the globe. The whole world is facing the second wave of this novel pandemic. Recently, the researchers are focused to study the complex dynamics and possible control of this global infection. Mathematical modeling is a useful tool and gains much interest in this regard. In this paper, a fractional-order transmission model is considered to study its dynamical behavior using the real cases reported in Saudia Arabia. The classical Caputo type derivative of fractional order is used in order to formulate the model. The transmission of the infection through the environment is taken into consideration. The documented data since March 02, 2020 up to July 31, 2020 are considered for estimation of parameters of system. We have the estimated basic reproduction number (R0) for the data is 1.2937. The Banach fixed point analysis has been used for the existence and uniqueness of the solution. The stability analysis at infection free equilibrium and at the endemic state are presented in details via a nonlinear Lyapunov function in conjunction with LaSalle Invariance Principle. An efficient numerical scheme of Adams-Molten type is implemented for the iterative solution of the model, which plays an important role in determining the impact of control measures and also sensitive parameters that can reduce the infection in the general public and thereby reduce the spread of pandemic as shown graphically. We present some graphical results for the model and the effect of the important sensitive parameters for possible infection minimization in the population.


Introduction
The fatal disease called a coronavirus  has grown beyond the expectation of everyone and seems to stay around for months or even years. The rapid spread of the virus with a high level of severity has become an unprecedented threat for the public health. To understand the mechanism of this virus and its possible eliminations, the researchers and scientists are looking to discover some treatments and vaccinations for its eliminations. The source of transmission is mainly through the droplets and person to person contact. To reduce the infection spread, the wearing of mask, stay at home, social distancing and wash the hands regularly are the highly recommended. Still the coronavirus infection spreading and according to the recent updates, over 20.2 million are infected with 787,909 deaths globally so far [20,13].
Like other countries, the rapid spreading of COVID-19 in Saudi Arabia; in which by the 30 of July 2020, an accounted 288609 confirmed cases of infected people, with 252035 recovered, and 3167 lost their lives as well as the severe economic hardship. The first case reported in Saudi Arabia was on second March 2020 while the second case has been reported on March 14 by the entry of the person from Iran via Bahrain [11]. After these cases, the government make plans and other preventive measure for the infection further spread. In this regard, a strict lock-down and curfew were announced, as a results the educations sector and the business markets were closed, flights suspended, Ummrah and Hajj are restricted to a certain number of people. The development of the vaccine and the possible control of the infections, the researchers and the biologists around the globe started thinking to have some vaccine for the disease elimination. The virologists are focusing their attention on vaccine ongoing trials, many to developing the vaccine and also a specific treatment mechanism to prevent the spread of this deadly virus. The researchers around the world focus to study the coronavirus infection in different scenarios. Some were studying the statistical tools to understand the infected cases and to find some relationships that can be used further to help the infection minimization. Some studying the virus biologically and to think for some vaccine developments. Besides this, the others who developing mathematical models in terms of mathematics and to predict the infection eradication stage. The mathematical models are were considered in fast for many infectious diseases and its outbreak. Like corona infection, the authors used the theory of differential equations to formulate mathematical models with different characteristics and provided useful results for its eliminations.
Mathematical modeling is a valuable tool to study very effectively the disease spread and its control. In this regard the researchers around the world presented many useful mathematical models in the last few decades to study the dynamics of infectious diseases and to determine some useful strategies for the better elimination of infection. The compartmental models along with real cases are more helpful to provide useful information about a particular disease outbreak. Several mathematical models are considered in literature in order to investigate and analyze the complex transmission pattern of the novel ongoing COVID-19 pandemic, see [18,19,28], using ordinary, stochastic and delay differential equations. More briefly, the COVID-19 dynamics through a mathematical with Pakistani data and with the impact of public interventions are studied in [23]. The most three high infected countries due to corona virus and its dynamics are briefly studied in [15]. Recently, the authors in [3] investigated a new COVID-19 compartmental model involving the environmental viral concentration in order to study the corona virus infection and its dynamics. The authors obtained the model parameters in [3] for the coronavirus data of Saudi Arabia. Further, they constructed the transmission model under the assumption that the people exposed to the infection can also transmit the virus to the community.
The fractional mathematical models that are regarding the generalized model and considered useful for modeling purposes in epidemiology. Various benefits can be obtained from a fractional order system in the sense of best data fitting, information about its memory and to identify the best possible value of the fractional order than can best describe the model for the real cases. In addition, the hereditary properties, it makes the models constructed in fractional derivatives more strong and useful for describing the real phenomenon. It is useful than the ordinary order and with the benefit of cross over behaviors. Different mathematical models with the effect of non-locality and with other interesting results are proposed in [16,2] and the references therein. Mathematical models where the authors proposed different fractional derivatives, see [5,21,17], and its applications, see [9,6,25]. Among the existing fractional order derivatives, the Caputo derivative gaining much interest and is widely used in problems arising in science and engineering. A vast literature is available on the application of Caputo type derivative. For instance, in [24], the Caputo fractional derivative is utilized in order to analyze the tuberculoses infection. The application of Caputo fractional derivative to a non-linear model of biochemical see [1] and the Hepatitis B virus with the use of hospitalization class, see [25]. In epidemiology, the Caputo fractional derivative is applied as a useful tool to provides a better understanding of disease dynamics. The dynamics of tuberculous (TB) in Pakistan was studied using the Caputo mathematical model in [24]. The authors in [24] showed that the Caputo model provides a better approximation to the reported infected cases. Recently, a number of Caputo-fractional order model are formulated to explore the transmission dynamics of novel coronavirus infection [22,7,8]. Some other mathematical models that addressed the coronavirus disease are studied in [4]. For example, the dynamical analysis of coronavirus model using lockdown effect has been analyzed in [4]. Using the real cases reported in Pakistan have been analyzed through a fractional mathematical model given in [7].
Motivated by the above discussion, the goal of the present work is to adopt the use of Caputo fractional-order time derivative to analyze the dynamics of COVID-19 in Saudi Arabia within the given period March 02, 2020 till July 31, 2020. We will study the model qualitatively to show the disease eradication using the important sensitive parameters. We extend the integer-order mathematical model recently studied in [3] to fractional order using Caputo derivative. The reproduction number shall be evaluated for the updated estimated and fitted model parameters. We can show more interesting results and its theoretically and its justifications numerically in the form of graphs. The effect of fractional order values on model as a graphically will be shown. This paper can be categorized with the sections is as follows: In Section "Basics of fractional calculus and model descriptions", we give the basics related to the fractional calculus and the model descriptions. The results of the system, existence and uniqueness shown in Section "Analysis of the fractional model" while the stability analysis of the model are discussed in Section "Iterative solution and stability analysis". We presented the graphical results with discussion in Section "Results and discussion" and concluding the findings in Section "Conclusion".

Basics of fractional calculus and model descriptions
We give here in brief the essential definitions regarding fractional calculus and the model description in fractional derivative. The following important definitions are considered that will be utilized as a application to our proposed study.

Definition 2.1.
Consider y ∈ C n be function, then the well known classical Caputo derivative having fractional order α in (n − 1, n) where n ∈ N [21] is defined as: clearly C D α t (y(t)) tends to y ′ (t) as α→1.
Definition 2.2. The Corresponding integral with fractional order α > 0 of the function y : R + →R is described as follows:

Definition 2.5.
Let y * denotes the equilibrium point of the Caputo fractional model then:

Model descriptions
We consider here to generalize the COVID-19 model with the description of memory effects by using the Caputo derivative. For this purpose, we divide the total populations of humans at any time t in five different sub-classes. These sub-classes include, the susceptible S(t), exposed E(t), asymptomatic infected (no clinical symptoms but infect healthy people) A(t), having disease symptoms and infected people or symptomatic infected I(t) and the recovered individuals R(t). The environmental class given by B(t) shows the concentration of virus in the environment. Then, the following is the representations of the corona virus model in terms of fractional differential equations: The COVID-19 model in Caputo operator form is; with the corresponding initial conditions (ICs) defined as: The birth rate of humans population and the natural mortality rate is considered as Λ and ν respectively. The parameters β i for i = 1, 2, 3 are used to describe the disease transmission rates through the direct transmission from the exposed, symptomatic and asymptomatic infected individual respectively. Whereas, β 4 shows the virus transmission rate from the environment or contaminated surfaces. The symptoms develop in exposed individuals at the rate (1 − ω)μ they become infected and join I(t) class while the remaining do not show any (or having mild) symptoms enter to asymptomatic class. The recovery rate in infected and asymptomatic class is ϕ 1 and ϕ 2 respectively. The people die at the rate contributed by E(t), I(t) and asymptomatic infected people to the environment and removal rate from the environment is ψ. In system (2), the net human population at time t is described by the term N(t). Moreover, By adding them, we observe we obtain Therefore, we have The feasible region Ω, shown by

Parameter estimation
In order to parameterize the biological parameters of the proposed model by considering the reported infected cases of corona virus in Saudi Arabia for the given period starting from March 02 till 31 July 2020 (the peak time of COVID-19), we make use of the well-known least square fitting technique. The birth rate denoted by Λ and the natural mortality rate ν are estimated from literature [12] as can be found in the table given below. The total population of KSA is N(0) = 34813871 and the average life span is 1/74.87 years. The other parameters of the model under consideration are estimated from the real data of the aforementioned period in KSA. Consequently, using the estimated and fitted parameters the reproduction number is evaluated R 0 ≈ 1.2937. The predicted curve is depicted in Fig. 1, which provides a reasonable fitting curve to the actual cumulative infected cases. The Table 1 describes the parameters values estimated from the real data.

Analysis of the fractional model
In this section, we investigate some basic and necessary mathematical features of the proposed model (2). To present the existence as well as the uniqueness (EU) of the problem solutions, we proceed as follow: Proof. System (2) gets the following for after utilizing the Caputo integral

Existence and uniqueness (EU) of the model solution
where the kernels are, and B(t) having an upper bound. Let the two function S(t) and S * (t) into consideration, and in similar manner for other functions. we have continuing in the same way for the remaining equations, we get where l 1 , l 2 , l 3 , l 4 , l 5 and l 6 denote the respective Lipschitz constants to the functions K 1 , K 2 , K 3 , K 4 , K 5 and K 6 . Hence, the Lipschitz condition are satisfied. The equations in (4) can be shown recursively as: The corresponding difference among the successive terms along with the initial conditions of model will takes the form; ∫ t 0 (t− ς) α− 1 (K 6 (ς,B n− 1 (ς))− K 6 (ς,B n− 2 (ς)))dς.
So by hypothesis, 1 Γ(α) l j m < 1, and S n , E n , I n , A n , R n , B n is the cauchy sequence. Hence, complete the desired result for the fractional-order model (2) is reached. □

Iterative solution and stability analysis
To establish the results in details, we provide the following results: . Then, the iterative approach, y n+1 = h(G * , y n ), G * is stable if lim n→∞ c n = 0, that is lim n→∞ c * n = p * for z n+1 = G * where, n is taken as the Picard's iteration then G * iteration is stable. The theorem can be summarized as below: (B , ‖.‖) defines a Banach space and G * be a self-map upon B , then for all x, y ∈ B , we have Let P defines a self map, then we have Proof. Since the map P is a fixed point, therefore, we evaluate the following equation Taking the norm of both sides of above equation After simplification, we have Now it assumed that, after substituting the above relation, we have Because the sequences S n (t), E m (t), I m (t), A m (t), and B m (t) are convergent and bounded, their exist five different constant S 1 >0,S 2 >0, S 3 >0,S 4 >0 and S 5 >0 for all t. Hence we have ‖S n (t)‖<|S 1 ‖, ‖E m (t)‖<‖ S 2 ‖, ‖I m (t)‖<‖S 3 ‖, ‖A m (t)‖<‖S 4 ‖, ‖B m (t)‖<‖S 5 ‖, (m, n)∈N× N.
By the relation we can attain Hence, the proof is complete. □

Disease free equilibrium (DFE) and the basic reproduction number
We obtain the basic reproduction number of the system (2) and its equilibrium point at the disease free case. Let the disease free equilibrium DFE of the system (2) denoted by E ** 0 = (S 0 , E 0 , I 0 , A 0 , R 0 , B 0 ), and can be obtained by setting, Solving these equations at the disease free state, we get In order to get the expression for the basic reproduction number R 0 of the system (2), we apply the technique known as the next generation method [26]. The well-known next-generation procedure is utilized in order to derive the basic reproductive number R 0 for the model (2). Let x = (E, I, A, B) T , then the necessary matrices are Now, for linearization the Jacobian of above matrices at disease freestate is: We have Where, for simplicity we take Thus, the expression of R 0 in term of model parameters is obtained as: where

Local stability of DFE
In order to study the local asymptotical stability of the system (2), we give the following result: Proof. Now the Jacobian of system (2) at E ** 0 is: Expanding the matrix (9) gives, The repeated two eigenvalues are − ν, − ν, contains negative real part. We used the polynomial in (10) and estimating the remaining eigenvalues by obtaining the coefficients given below: Clearly, a i are all positive for i = 1, .., 4 if R 0 < 1. The Routh-Hurtwiz criteria can be used further to satisfy the conditions and make it sure the asymptotical stability of the system (2). This computation can be done easily through some computational software. The argument of the roots of equation (λ m + ν) 2 = 0 are , where k = 0, 1, 2, …., (m − 1).
Similar approach can be adopted to show that the argument of the roots of Eq. (10) are all greater than π 2M if R 0 < 1. Moreover, we can find an argument less than π 2M for R 0 > 1. Thus disease-free equilibrium is LAS for R 0 < 1. □ Lemma 4.1. The system (2) is unstable at the disease-free equilibrium if R 0 > 1.

Global stability of DFE
The Lyapunov function approach will be used to proceed the results for GAS of proposed model at disease free-state. For this we state and prove theorem: Proof. We present the Lyapunov function defined below by establishing the above result: where the coefficients f j > 0, for j = 1, …, 4, denoting the arbitrary positive constants, and later in the proof we can assign appropriate values. Applying the Caputo-fractional derivative on G(t), together with the use of the model (2), we have, . Thus (E, I, A, B)→(0, 0, 0, 0) as t→∞ and model (2) implies that S→ Λ ν , and R→0 as t→∞. It can be followed from the results given in [27], the solution of system (2) with non-negative ICs tends to E ** 0 whenever

The endemic equilibrium
In order to find the endemic equilibrium (EE) of the model (2), we denote it by where, N ** = (S ** + E ** + I ** + A ** + R ** ). The system (2) gives the following results at a steady state: Further, at endemic steady-state, Substituting (14) in (15), and after simplification, we get where, If R 0 > 1, then a unique EE exist.
Proof. Defining a non-linear Lyapunov function for the system (16), and then by taking its time fractional derivative, we have the following result, From (19), we have . . .
after substituting, we get We have the following fact, Hence, C D α t G (t) ≤ 0 for R 01 . Thus, by Lasalle invariance principle, So, every solution associated to the model (16) approaches to its unique endemic equilibrium. So, whenever R 0 > 1, the system (16) at (E ** 1 ) is globally asymptotically stable. □ Conjecture 1. The unique EE of the system described in (2) is GAS in Ω, if R 0 > 1, in Ω.

Results and discussion
We present the simulation and discussion for the Caputo COVID-19 model (2) in the given section. We consider the numerical results given in Table 1 for the numerical results. The proposed fractional model is solved numerically using the method described in details in [10,14]. The biological parameters estimated from COVID-19 confirmed cases and tabulated in Table 1 to assess the dynamics of novel COVID-19 pandemic in KSA. In order to analyze the role of different parameters and memory index in the dynamical behavior of the disease incidence, we varied various values of the key parameters of the model and the α.  Fig. 3 and 4 respectively. This behavior is interpreted for two values of the order of Caputo operator i.e., α = 1, 0.90. It is found that by decreasing the contact rates to its estimated baseline for two distinct values of α, we can see a significant decrease. Although, a comparative faster decrease is observed in the infected population for α = 0.90. But, the peaks of infection curves occur at the longer period of time in this case. The role of viral concentration on the COVID-19 incidence is depicted in Fig. 5. By decreasing the virus concentration rate, one can       observe the decrease in the infected compartments as shown in Fig. 5. Fig. 6 shows that decay in the infected population due to decreasing the removal rate ψ. The behavior is even more prominent for a smaller value of fractional order α as can be seen in 6(b-d).
Moreover, we also investigated the impact of some parameters on the dynamical behavior of cumulative symptomatic and asymptomatic infected cases. The impact of contact rates β 1 and β 4 on the cumulative symptomatic and asymptomatic people is shown in Fig. 7 and Fig. 8 respectively. The reduction in these parameters can reduce the cumulative symptomatic and asymptomatic people significantly as can be seen in Fig. 7 and 8. This reveals that the disease incidence can be reduced following protective measures i.e., using a facemask, sanitizer, glues, etc. through disinfection spry in order to reduce the environmental viral load. Furthermore, the decrease in the cumulative symptomatic and asymptomatic population on the reduction in the virus concentration rate φ 1 due to exposed people is depicted in 9. This behavior becomes more prominent for smaller values of α as shown in 9 (b-c). From Fig. 10 the role of removal rate of the virus from the environmental (via disinfection spray etc.), is observed and seems to be more feasible for small values of the fractional operator α.

Conclusion
We established a mathematical model for the understanding the corona virus infection in the Kingdom of Saudi Arabia through a fractional mathematical model in Caputo sense. We considered the real cases reported in the Kingdom have been analyzed and obtained a reasonable fit to the data. The environmental impact on the COVID-19 incidence is taken into consideration. Initially, the model in integer order were considered and then the fractional order operator with the power law kernel were applied for its generalization. We provided the related properties for the fractional model. The threshold quantity R 0 and parameters were estimated from the reported COVID-19 cases in Saudi Arabia with the help of nonlinear least square procedure and we found that R 0 = 1.2937. The stability analysis at the DFE and EEP are explored in detailed. We used a numerical scheme for the solution of the fractional order model and presented the graphical results. The predictor corrector scheme of Adams Molten type is utilized and studied the influence of various parameters by varying to its baseline value for distinct values of α. The dramatic reduction in the disease burden was observed with an enhancement in contact tracing policy. The impact of variation in the environmental viral load due to symptomatic and asymptomatic COVID-19 infected individuals is analyzed graphically. It is observed that by reducing the viral contribution in the environment by asymptomatic infected individuals (i.e.,φ 3 ) decrease the disease burden significantly. We also depicted the impact of variation in the removal rate of virus from the environment (or surfaces) graphically on the disease prevalence. In overall simulation results, that the sensitive parameters decrease the infection in the population very fast. Such important parameters that can be considered as a control for the infection eradication in the population. The present work can be extended by using the cases in the state of Saudi Arabia with the second wave in more recent introduced fractal-fractional operators.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.