A fractional order mathematical model for COVID-19 dynamics with quarantine, isolation, and environmental viral load

COVID-19 or coronavirus is a newly emerged infectious disease that started in Wuhan, China, in December 2019 and spread worldwide very quickly. Although the recovery rate is greater than the death rate, the COVID-19 infection is becoming very harmful for the human community and causing financial loses to their economy. No proper vaccine for this infection has been introduced in the market in order to treat the infected people. Various approaches have been implemented recently to study the dynamics of this novel infection. Mathematical models are one of the effective tools in this regard to understand the transmission patterns of COVID-19. In the present paper, we formulate a fractional epidemic model in the Caputo sense with the consideration of quarantine, isolation, and environmental impacts to examine the dynamics of the COVID-19 outbreak. The fractional models are quite useful for understanding better the disease epidemics as well as capture the memory and nonlocality effects. First, we construct the model in ordinary differential equations and further consider the Caputo operator to formulate its fractional derivative. We present some of the necessary mathematical analysis for the fractional model. Furthermore, the model is fitted to the reported cases in Pakistan, one of the epicenters of COVID-19 in Asia. The estimated value of the important threshold parameter of the model, known as the basic reproduction number, is evaluated theoretically and numerically. Based on the real fitted parameters, we obtained \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathcal{R}_{0} \approx 1.50$\end{document}R0≈1.50. Finally, an efficient numerical scheme of Adams–Moulton type is used in order to simulate the fractional model. The impact of some of the key model parameters on the disease dynamics and its elimination are shown graphically for various values of noninteger order of the Caputo derivative. We conclude that the use of fractional epidemic model provides a better understanding and biologically more insights about the disease dynamics.


Introduction
The novel virus (2019-nCoV) that is highly transmissible and pathogenic was first identified from a single individual in Wuhan city in China. This novel infection causes a severe acute respiratory syndrome and it has spread across the world. The reported COVID-19 confirmed cases are over 10.27 million, and there have been more than 0.5 million deaths till 30 June 2020 globally so far [1]. The worst affected regions due to coronavirus are America, Europe, Africa, South-East Asia, Western Pacific, Eastern, and Mediterranean. The initial symptoms of a COVID-19 infection include dry cough, fever, fatigue, and breath shortening that appear in 2-10 days and further cause pneumonia, SARS, kidney failure, and even death [2]. The pandemic has continuously spread due to absence of vaccine and antiviral treatments. Thus WHO announced it a global issue. The policy makers have implemented the non-pharmaceutical intervention like social distancing, self-quarantine, isolation of infected, wearing mask, protective kits for medical personnel, and travel restrictions to minimize the disease incidence. It is also a challenging problem for scientists and virologist evaluating potential treatments based on ongoing clinical trials.
Researchers suggested many mathematical models to analyze the dynamical behavior and spread of the novel virus which can help to predict the future situation and even control of the COVID-19 pandemic [3]. In the analysis of mathematical models of coronavirus, the reproductive number has a significant role in describing the nonlinear dynamics of physical and biological engineering problems. The reproduction number indicates that COVID-19 has been continuously increasing or has been controlled. In Pakistan, 209,337 confirmed infected cases have been reported and about 4304 have lost their lives out of over 220 million population to date [4,5]. The first case was reported in Karachi on 26 February 2020, and day by day situation is getting worse and virus is spreading quickly due to limiting testing. The government is unable to maintain strict lockdown and has imposed a smart lockdown by easing restrictions due to severe economic hardships, especially for labor community who earns for living to survive every day. To study the dynamics of COVID-19 transmission pattern, many mathematical models provide more insight on how to control the disease spread to health authorities [6][7][8]. Fanelli and Piazza [9] studied a novel compartmental model describing the transmission patterns of COVID-19 in three highly infected countries. The dynamics of COVID-19 with an impact of non-pharmaceutical interventions was studied by Ullah and Khan [10] on Pakistani data. The fractional mathematical models rendering the natural fact in a systematic way as in [11,12] and [13,14] are used to simulate the transmission of coronavirus. Different mathematical models with an effect of nonlocality and fading memory process by using differential operators have been presented [14,15]. The fractional order epidemic models are more helpful and reliable in analyzing the dynamics of an infectious disease than the classical integer order models [16,17]. The fractional order models for different diseases show cooperatively better fit to the real data. In [18,19] a different fractional operator is suggested, and applications of these fractional operators are found in [20,21]. Recently a Caputo fractional order COVID-19 model has been studied in [22]. Some other fractional mathematical models for the investigations of infectious diseases have been studied in [15,[23][24][25][26]. For example, the fractional diffusion equations and their analysis are studied in [15]. A new scheme to solve numerically the fractional order differential equations is utilized in [23]. New numerical investigations for the fluid in nonconventional media are suggested in [24]. Coronavirus modelings, simulations, and their possible control through a mathematical model are studied by the author in [25]. The spread of coronavirus in South Africa and Turkey with detailed statistical and mathemati-cal results is studied in [26]. Recently, the authors have studied the analysis of coronavirus model in fractional derivative [27]. The use of quarantine and isolations in the modeling of coronavirus is investigated in [28]. A mathematical model for the dynamical analysis of coronavirus and its control analysis is studied by the authors in [10]. The notified cases of coronavirus in Saudi Arabia through a mathematical model are considered in [29], where the authors provide suggestions on possible controls based on the parameters.
Environmental viral load plays an essential role in the disease incidence and is considered to be one of the main transmission routes of COVID-19. In this study, we reformulate the model [28] with the impact of quarantine, isolation, and environmental effects on the transmission dynamics of coronavirus with the application of Caputo derivative. The parameter values are estimated from the cumulative COVID-19 cases reported in Pakistan. The fractional order models provide better understanding and give more insights about the pandemic. The rest of the work is arranged as follows: In Sect. 2 basics preliminaries are presented, while the model formulation for integer case with parameter estimation and curve fitting is presented in Sect. 3. Model derivation and basics properties are presented in Sects. 4 and 5, respectively. In Sect. 6 we present the analysis of the model, while the numerical simulations are depicted in Sect. 7. Brief concluding remarks are presented in Sect. 8.

Preliminaries on fractional derivative
In order to proceed, first we recall some basic definitions regarding fractional calculus.

Definition 1
The fractional order derivative in the Caputo case with order α for a function g ∈ C n is defined as follows [19]: Clearly, C D α t (g(t)) tends to g (t) as α → 1.

Definition 5
The constant point z * is an equilibrium point of the Caputo-fractional model, then [30] if and only if g(t, z * ) = 0.

The classical integer order model formulation
In this section, we briefly discuss the integer order model of the dynamics of COVID-19 with quarantine, isolation, and environmental load which is mainly studied for the Chinese population data in the case of fractal-fractional Atangana-Baleanu derivative [28]. For the mathematical model formulation, the net population N(t) at time t is further divided into mutually-exclusive sub-compartments as susceptible or healthy S(t), exposed E(t), infected with clinical symptoms I(t), asymptomatically infected A(t), quarantined class Q(t), hospitalized H(t), and the recovered R(t) individuals. The class M(t) denotes the environmental viral load due to the infected people. It is to be that the asymptomatic infected individuals are also capable of transmitting the infection. The nonlinear ordinary differential equations governed by these assumptions are described as follows: where subject to the following nonnegative initial conditions: The biological description of the parameters involved in COVID-19 model (6) is given in Table 1.

Parameter estimation
In this section, we estimate the model parameters with the help of the well-known statistical technique known as nonlinear least square curve fitting approach. The confirmed  infected cases from 1 March to 30 June 2020 in Pakistan are used to estimate the model parameter values [31]. The recruitment rate and the natural mortality rate are estimated from literature, and other parameters are fitted from real data. The updated reproduction number evaluated by using the estimated and fitted parameters is R 0 ≈ 1.50. The predicted curve having a better agreement to the actual reported cases is depicted in Fig. 1. The model parameter values are shown in Table 1. The initial conditions used in data fitting as well as in numerical results are given as S(0) = 220,870,336, E(0) = 20,000,

Model derivation in Caputo operator
We reformulate the COVID-19 fractional order (6) model using a Caputo fractional derivative in order to observe the memory effects and gain more insights about the pandemic. System (6) in terms of integral form is followed by substituting the value of kernel as a power-law correlation function. After applying the Caputo fractional derivative of order α -1, we get Both are the inverse operators, we get the COVID-19 model in the Caputo operator form The initial conditions are as follows: 5 Basic properties

Invariant region
The dynamics of Caputo fractional model (10) is explored in a feasible region ⊂ R 8

Lemma 1
The region ⊂ R 7 + is positively invariant with nonnegative initial conditions for model (10) Proof After summing the components of human population in model (10), we obtain a total population as follows: and then we have After applying the Laplace transform, we have By considering the inverse Laplace, we arrive at where the Mittag-Leffler function is shown by , and the Laplace transform is Thus, the solution of the model with the nonnegative conditions in remains in . So, the region is positively invariant and attracts all the solutions in R 7 + .
Now, for the positivity of the system solution, let

and y(t) = S(t), E(t), I(t), A(t), Q(t), H(t), R(t) T .
Corollary 1 [32] Suppose that g(t) ∈ C[m, n] and C D α t g(t) ∈ (m, n],where α ∈ (0, 1]. Then if (i). C D α t g(t) ≥ 0, ∀y ∈ (m, n), then g(t) is nondecreasing; (ii). C D α t g(t) ≤ 0, ∀y ∈ (m, n), then g(t) is nonincreasing. Proof To show that the solution of the model is nonnegative, it is required to show that on each hyper-plane bounding the positive orthant there is a the vector field point R 7 + . From system (10), we have

Positivity and boundedness
The solution is Also, the Mittag-Leffler function is bounded ∀ t > 0. Therefore, we have The solution of the system will remain in R 7 + by using Corollary 1, and hence a biologically feasible region is constructed as follows: Since all the terms are positive, the solution of model (10) is bounded.
6 Analysis of the model

Disease-free equilibrium (DFE)
For the equilibrium points of fractional model (10), we have The disease-free equilibrium is denoted by E * 0 = (S 0 , E 0 , I 0 , A 0 , Q 0 , H 0 , R 0 , M 0 ) and is given by

The reproduction number
The next generation technique is used to derive the reproduction number R 0 for the dynamics of disease. Let x = (E, I, A, Q, H, M) T , then we have where The next generation matrix is of the form The spectral radius of the next generation matrix is shown by
In a similar fashion, we also find that the arguments of the equation (λ 4q + a 1 λ 3q + a 2 λ 2q + a 3 λ q + a 4 ) = 0 are all greater than π 2M if R 0 < 1, having an argument less than π 2M for R 0 > 1. The DFE is locally asymptotically stable for R 0 < 1.

Global stability of DFE
The Lyapunov function approach is used to proceed to the result for the GAS of the proposed model at the disease-free state. For this, we have the following theorem.

Theorem 2 The DFE of the Caputo COVID-19 model is GAS if
Proof Consider the following appropriate Lyapunov function: where the coefficients Yj, for j = 1, 2, 3, 4, are unknown positive constants, and they will be chosen later. The Caputo-fractional derivative of F(t), along model (10), yields Let us choose Hence it follows that C D α t F ≤ 0 for R 0 ≤ 1, all the parameters and variables are nonnegative with C D α t F = 0 iff E = I = A = M = 0. Thus (E, I, A, M) → (0, 0, 0, 0) as t → ∞. By using E = I = A = M = 0 in the first and second last three equations of (10) implies that S → δ , and Q, H, R → 0 as t → ∞. Thus, using Lyapunov stability theorems for the fractional case developed in [30], the solution of Caputo model (10) with nonnegative initial conditions approaches to E * 0 as t → ∞ in a feasible region. Thus, it follows that the DFE of model (10) is GAS.

Theorem 3
The proposed model has: and the discriminant is greater than zero, then there exist two endemic equilibria. (iv) Otherwise, no endemic equilibrium exists.
Thus from case (i) of Theorem 3, it is clear that model (10) has a unique positive endemic equilibrium when R 0 > 1, and point (iii) indicates the possibility of backward bifurcation when R 0 < 1.

Numerical results
This section presents simulation and discussion for the Caputo COVID-19 model (10). The proposed fractional model is solved numerically using a generalized predictorcorrector of the Adams-Bashforth-Moulton method [33,34]. The biological parameters estimated from the actual data reported in Pakistan from 1 March to 30 June 2020 and tabulated in Table 1 are utilized to obtain the simulation results. We varied different model parameters and the order α of the Caputo operator in order to explain the role of various parameters and memory index on the disease transmission patterns and control. The impact of contact rates ζ 1 , ζ 2 , ψ (a measure of social distancing effects), quarantine or contact-tracing policy ε 1 , and hospitalization or self-isolation rate τ is depicted graphically in Figs. 2-6. The impact of environmental viral load due to the symptomatic (m 1 ) and asymptomatic infected individuals (m 2 ) is depicted in Figs. 7 and 8, respectively. Finally, the role of removal rate of virus from the environment (via disinfection spray etc.) is shown in Fig. 9.

Effect of reduction in contact rates (a measure of social distancing)
First of all, we simulate COVID-19 model (10) for different values of effective contacts ζ 1 in order to evaluate the impact of social distancing on the dynamics of cumulative symptomatic, asymptomatic, and hospitalized infected individuals. The simulation results for different four values of fractional order are shown in Fig. 2(a)-(d). It can be observed from these interpretations that the reduction in the effective contact rates ζ 1 up to 50 percent to its baseline values reduces the peak of cumulative COVID-19 infected cases very well. The behavior is the same, seen for all four values of fractional order, although slightly faster decrease in pandemic peak is observed for smaller values of α. The impact of disease transmission rate due to the environmental viral load ζ 2 is demonstrated in Fig. 3(a)-(d).
With the reduction in this rate, a rapid decrease is observed in the total of symptomatic, asymptomatic, and hospitalized infected individuals and even becomes more biologically significant for smaller values of fractional order. This interpretation shows that the environmental transmission plays a significant role in disease prevalence. Finally, the influence  Fig. 4(a)-(d) for four different values of α. A reasonable decrease in the pandemic peak is observed with a reduction of contact rate ψ. Overall, from these interpretations, we conclude that maintaining a strict social distancing and avoiding public gathering are needed to control the transmission of infection in the future.

Effect of contact-tracing and hospitalization/self-isolation
The quarantine contact-tracing policy of exposed individuals and hospitalization/isolation of confirmed infected cases is another effective strategy to control the COVID-19 incidence. Here, we simulate the model with different levels of quarantine and hospitalization rates. The resulting graphical interpretation is shown in Figs. 5(a)-(d) and 6(a)-(d), respectively. A dramatic decrease in pandemic peak is observed with an increase (up to 50 percent) in the contact-tracing rate as shown in Fig. 5(a)-(d). The same behavior is obtained for all values of α. The influence of the increase in hospitalization of confirmed infected individuals is analyzed in Fig. 6(a)-(d). These results suggest that in the absence of treatment or vaccine, a proper and effective contact-tracing and self-isolation policy should be adopted until the elimination of pandemic.

Effect of environmental viral load
The environmental viral load is one of the main routes of the COVID-19 infection transmission. In this section, we simulate model (10) with the reduction in the viral contribution to the environment by symptomatic and asymptomatic COVID-19 infective individuals at different rates. Figure 7(a)-(d) demonstrates the impact of the reduction in m 1 with different levels. It is observed that a reduction in m 1 (up to 50 percent) reduces the infected curve with a reasonable rate as seen in Fig. 7(a)-(d). The decrease in pandemic peaks is relatively more significant for smaller values of the fractional order of Caputo operator. The role of m 2 (viral contribution of asymptomatic individuals to the environment) is analyzed in Fig. 8(a)-(d). It can be seen that a reduction in m 3 at 50 percent to its baseline values given in Table 1 reduces the total number of symptomatic, asymptomatic, and hospitalized infected individuals very well. The similar behavior is obtained for all values of fractional order α. Finally, the influence of the removal rate of virus from the environment is depicted in Fig. 9(a)-(d). It is observed that enhancement in m 3 at 40 percent to its baseline value dramatically reduces the infection curves as shown in Fig. 9(a)-(d). This interpretation suggests that proper disinfection spray is necessary to remove the viruses from the environment in order to reduce the disease prevalence.

Conclusion
The environmental viral load plays an essential role in the spread of the current COVID-19 pandemic. Limited research has been done on the environmental impacts on disease dynamics. In this paper, we studied the transmission dynamics of COVID-19 pandemic with quarantine, hospitalization (or self-isolation), and the environmental viral load impacts through a Caputo fractional model. Initially, the model is presented with classical integer order differential equations and then reformulated using fractional order operator with the power-law kernel. The basic mathematical analysis of the Caputo COVID-19 model is performed and the local and global stability of the disease-free case is proven. Utilizing the nonlinear least square approach, the model parameters are estimated and fitted for the COVID-19 actual cases recorded in Pakistan. The most important threshold quantity, known as the basic reproduction number, is presented both theoretically and numerically. The predictor-corrector iterative scheme of Adams-Moulton type is applied in order to solve the model numerically. The real estimated and fitted parameters are used to simulate the Caputo COVID-19 model. We depicted and discussed the simulation results for three different sets of model parameters showing the impact on the dynamics cumulative symptomatic, asymptomatic, and hospitalized infected cases. Firstly, we analyzed the influence of reduction in contact rates (a measure of social distancing) on the infected curves. The graphical interpretation of this case revealed that with a reduction of effective contacts with infected individuals, a significant decrease in the cumulative infected cases is seen. It is further observed that the reduction in the environmental viral load dramatically reduces the pandemic peak. Secondly, we explored the impact of the most commonly used interventions i.e., quarantine (or contact-tracing policy) and hospitalization (or self-isolation) on cumulative infected population. These graphical interpretations demonstrated that contact-tracing of exposed people is more effective in the elimination of infection. It can be seen that with an increase in the quarantine rate up to 50 percent to the base line value reduces the infected cases significantly. Finally, we have shown the impact of variation in the environmental viral load due to symptomatic and asymptomatic COVID-19 infected individuals. It is found that reducing the rate of viral release into the environment by asymptomatic infected individuals decreases the disease burden very well. It is because the asymptomatic infected individuals are unaware of the infection and freely continue the daily routine without any necessary safety precautions like the use of face mask and sanitizer. Moreover, we have also depicted the impact of variation in the removal rate of virus from the environment (or surfaces) on the disease prevalence. A dramatic reduction in the infected cases is found with an enhancement in this rate, which demonstrates that disinfection spray is very helpful in the disease elimination. We believe that the investigations in the paper will be beneficial for the health and decision-making authorities to combat the disease in the community. In the future, the present model can be reformulated by incorporating some suitable time-dependent control interventions using optimal control theory.