To study the transmission dynamic of SARS-CoV-2 using nonlinear saturated incidence rate.

In this work, we construct a new SARS-CoV-2 mathematical model of SQIR type. The considered model has four compartments as susceptible S , quarantine Q , infected I and recovered R . Here saturated nonlinear incidence rate is used for the transmission of the disease. We formulate our model first and then the disease-free and endemic equilibrium (EE) are calculated. Further, the basic reproduction number is computed via the next generation matrix method. Also on using the idea of Dulac function, the global stability for the proposed model is discussed. By using the Routh-Hurwitz criteria, local stability is investigated. Through nonstandard finite difference (NSFD) scheme, numerical simulations are performed. Keeping in mind the significant importance of fractional calculus in recent time, the considered model is also investigated under fractional order derivative in Caputo sense. Finally, numerical interpretation of the model by using various fractional order derivatives are provided. For fractional order model, we utilize fractional order NSFD method. Comparison with some real data is also given.

Introduction

In human history many out breaks have been reported in which millions of people died and suffered. It is a great fact that some of the worst epidemics and outbreaks in history have ruined whole civilizations and pushed dominant nations to their knees. During the mentioned outbreaks millions of people had been died. Still these dreadful outbreaks of various disease threaten humanity. Recently a deadly virus called Coronavirus (SARS-Cov-2) which has been identified in the Wuhan city of the Republic of China at the end of 2019. The said infection has been become a serious threat to human life now a days. All countries of the world are suffering with a great loss of economic as well as lives of their publics. Scientists and physicians, researchers are working continuously to make various policies, procedures for its eradication and proper cure day and night. Every state of the world has taken serious actions and made some strategies to deal with this pandemic. Health research department also is working day and night for safe vaccine for this deadly Coronavirus. For elimination of the mentioned disease different approaches have been made up to date (see [1], [2]). In start, the virus had been reported that transmitted from animal to humans. But later on researchers proved that it can also transfer from human to human (see [3], [4]). Recently scientists and researchers have gotten great success in preparing proper vaccine for the cure of said infection. Different kind of vaccines have been prepared and now available throughout the world. The vaccines have produced great impact on the reduction of the disease. In this regards some recent work published which has reported the usefulness of vaccines are refereed as [5], [6], [7]. Further a detailed survey of death ratio among vaccinee and unvaccinate people has been conducted in [8]. According to the report approximately 540 millions people have been infected and nearly 6.4 millions have been died due to infection throughout the world (see [9]).

Currently, many methods/ approaches have been adopted to study and analyze the mentioned infection. Since in last several decades advancement in epidemiology has been taken place. Therefore various tools and procedures are adopting day by day for investigating and controlling such disease in society. From research point of view mathematical modeling is one of the important tool to investigate infectious disease. For effective control of the diseases, mathematical models help to make decision for public health and also support future prediction (see few study as [10], [11], [12], [13], [14]). Mathematical models are helpful to investigate spatiotemporal patterns and dynamical behaviors of infections which spreading in communities. Having the important of mathematical models in mind, researchers have also studied SARS-Cov-2 from different point of view during the last two years. Usually mathematical models are formulated in terms of ordinary differential equations. The concerned equations are then evaluated by various tools to bring out some useful results and information. Hence researchers of this area are using different methods to ensure effective control over this disease with the help of mathematical models (see some work as [15], [16], [17], [18], [19]).

On the other hand, fractional calculus has been used as strong tools in mathematical modeling of various real world process. Greater popularity has been gotten because fractional differential operators preserve memory terms as well as past history. Further, the dynamical aspect of various real world problems like biological and physical phenomenon can be well explained via the concept of the said calculus (see some detail in [20]). For dealing the aforesaid area numerous methods and tools have been developed in previous few years (we refer few as [21], [22], [23], [24], [25]). Various methods and schemes including Adomian decomposition method, differential transform, perturbation method, etc have been used for the computation of approximate results.

For investigation of the transmission dynamics various models have been used. For instance authors [26] have used SIR type model to study the transmission dynamics of the aforesaid disease. The said model includes the classical linear incidence rate. In same line authors [27] have applied SIR type model to investigate the description of the successive waves of COVID-19. It should be kept in mind that different models of SQIR type have recently been considered under various incidence rate in literature to study COVID-19. For instance authors have considered SEIQR type models to study the situation of COVID-19 in Saudi Arabia in [28], [29]. Also SEIR type models have formulated in [30], [31] to study the said disease dynamics in Saudi Arabia. Authors [32] have studied the global dynamics of COVID-19 by using SQIR type model. Inspired from the above work, we establish a SQIR type mathematical model for SARS-Cov-2 first. Here, we consider the proposed model under the nonlinear saturated incidence rate. Then through the next generation matrix method, we compute the basic reproduction number, disease-free and endemic equilibrium points for the proposed model. Furthermore, by using the idea of Dulac function, the global stability for the proposed model is investigated. Also by the idea of the Jacobian matrix method, local stability is investigated. Since numerical simulation is an important aspect in the area of differential equations. Therefore various techniques have been developed in past for the said purposes. In this connection, one is led to perform numerical approximation using some known scheme like the Euler or Runge–Kutta 2 or Runge–Kutta 4 methods. These schemes are based on approximation theory which focus on the procedure to generate finite representation of functions involve in the considered problems. This is very tedious job mostly to get a close agreement between exact and approximate solution. This is due to the production of numerical behavior which does not present in the considered model or problem. For instance creation of ghost equilibrium points, change in the nature of stability of existing equilibrium point or destruction of domain invariance and so on. In fact these artificially produced by classical numerical schemes devoted to the non persistence of some important characteristics of the dynamics produced by the differential equations. To overcome these shortcoming, Mickens has been presented many results for it in the form of series of articles (see detail in [33], [34], [35]). Further the detail study of the scheme including stability and convergence has been given in [36], [37], [38]. The aforesaid scheme has been extensively used in various papers for numerical interpretations of various dynamical problems (see [39], [40]). Further the mentioned scheme has also applied in dealing some mathematical models (see [41]).

On the other hand, fractional-order dynamical systems have been considered for modeling majority of real world phenomenons due to their control over the memory effects. This characteristic has make fractional calculus more papular to consider in mathematical modeling of various real world problems. Keeping the importance of fractional calculus in mind, we will develop some results for numerical interpretation of the model under our consideration by using fractional order (NSFD) method. The history of fractional calculus is old like classical calculus. Different authors have investigated the said area from different aspects, we refer few as [42], [43], [44], [45], [46]. The said calculus has significant applications in various fields of science and technology (see [47], [48], [49]). The fractional order derivative has more degree of freedom due to its nonlocal nature. Therefore it is better to use fractional calculus in mathematical modeling of infectious disease as compared to classical calculus. In last few decades various authors have increasingly used the said derivative in dealing various mathematical models of physical as well biological, dynamical phenomenon, (detail can be seen in [50], [51], [52], [53], [54]). Here, we remark that COVID-19 has been investigated very well under the concept of various fractional order derivatives by applying different numerical schemes (see [55], [56], [57]).

Also for the numerical analysis of various problems, many authors have extended the NSFD method to fractional order (see [58], [59], [60], [61]). Here, we will use the mentioned method for the numerical analysis of our fractional order model also.

The mathematical model

In this part of our work, we assume to divide the entire community into four compartments. Susceptible , those people who are healthy but have the ability to gain the disease, infected , those people who are infected from the disease. In same line quarantine represents those people who have some symptoms and recovered stands for those people who have gotten ride from the disease. Here is the incidence rate of the COVID-19. Keeping these points in mind, we formulate our model as The flow chart presents the class wise dynamics of the system (1) in Fig. 1.

Fig. 1

Flow chart of the model given in (1).

The parameters involve in the above model (1) and their interpretation are given in Table 1.

Table 1

Nomenclatures and interpretation.

Symbols	Nomenclature
S	Susceptible population
I	Infected population
Q	Quarantine population
R	Recovered population
μ	mortality rate
k	Saturation constant
δ	Removal from quarantine class
α	Removal from infection
d_1	Disease death rate in infected class
d_2	Disease death rate in quarantine class
b	Individual goes to quarantine compartment
γ	Unprotected individual
β	Recruitment rate

Stability analysis

Here, we present some results regarding global and local stability of various equilibrium points. Also, we compute .

Here, we define the complete norm space by under the norm . Further the product space given by is also a Banach space with norm defined by , where .

We consider the existence of equilibrium for the system (1). The disease free equilibrium (DFE) for the Model (1) is represented by .

Also the EE is computed as In mathematical models of epidemiology, the basic reproduction number describes the transmission and predicts the control of disease. From , we can predict about the disease that weather it is going in the population or reducing. Also we make plan for the best way to protect the members of the community from this infection. The method of next generation matrix is used to compute . Let , then form the system (1), we have where and Jacobian of at the DFE is given by and for the DFE Jacobian of is given Hence We have Hence the required value of is given by From positive values of the parameters, we get . This indicated that the transmission of SARS-Cov-2 is well managed by community in their country. Here on the basis of reproduction number , we establish the following result.

There exist no positive equilibrium of the model (1) , if .

The EE of the system (1) is locally asymptotically stable if .

Local stability analysis of the model (1)

Here in this part of our paper, we discuss local stability for our proposed model (1).

DFE of the system (1) is locally asymptotically stable at with the condition and unstable if .

We calculate the variational matrix for (1) at as Jacobian matrix at equilibrium point is given by Here we have from Jacobian matrix, the following characteristic equation given by which upon calculation yields that where It is clear from the above Eq. 4 that three eigen values are negative. The 4th one is negative if or which implies that . Hence we conclude from the above result, that the equilibrium point is stable locally asymptotically if . Also from the characteristic equation, the eigen values are non-positive and unstable if .

The EE point is locally asymptotically stable, if the condition holds and vice visa.

Consider the equations We calculate the variational matrix for (1) at as Jacobian matrix at the equilibrium point is given by Let use the notions as and Then the characteristic equation from the given matrix is given by On further evaluation, one has From (4), we can see two of eigen values are negative and other two eigen values are calculated from axiality equation where Hence by using criteria of Routh–Hurwitz, the model(1) is stable locally asymptotically if and and unstable if also . So EE point is stable if Putting the value of in (5), we get that is . Hence the equilibrium point is locally asymptotically stable with the condition . □

Global stability analysis of the model (1)

In this part of our paper, we discuss global stability for our proposed system (1). First, we will reduce the system (1) and then take a Dulac function to check the global stability. The limit set of the model (1) is on the plane . Thus, we are focusing on the reduce model given by

The system (6) does not has nontrivial periodic orbits.

We take the following Dulac function for the system (6) with and as We have from (7) and Upon taking partial derivative of (8)–(10) and adding, one has This follow the conclusion. □

The EE point of the model (1) is globally asymptotically stable if the condition holds.

Take the Dulac function as From (11), (6), we have and Upon taking partial derivative of (12)–(14) and adding, one has Hence from (15), we have . □

Numerical analysis and discussion of model (1)

In this section, we present simulation for the model (1). For numerical simulations, we consider some reported data of Pakistan for the given time frame. use NSFD scheme [62], [63] to rewrite the first equation of the model (1) in the following form Assume that be the approximation of and , then we express approximation of left hand side of (16) as where ia an parameter. If we choose , then the approximation of first order derivative is given as Also the consistency condition given below should be satisfied by the function as For instance , satisfy the consistency condition as proved in [64]. Keeping the above detailed in mind, Eq. (16) can be approximated by NSFD scheme as As like (19), we approximate the whole system by using NSFD scheme as

The given NSFD scheme is convergent with order one. The proof has been given in [40].

We simulate our model (1) by using the NSFD scheme established in (20) and using the numerical values of Table 2.

Table 2

The physical interpretation of the parameters and numerical values.

Symbols	Nomenclature	Numerical value
S_0	Susceptible individual Class	217.338802 millions
I_0	Infected individual Class	1.286825 millions [65]
Q_0	Quarantine individual Class	2.56 millions
R_0	Recovered individual Class	1.274373 millions [65]
β	birth rate	0.1243
μ	Natural death rate	0.02 [65]
d_1	Death rate of corona in infected class	0.019 assumed
d_2	Death rate of corona in quarantine class	0.001 assumed
b	The rate of infection from quarantine	0.0205 assumed
κ	Infection rate	0.0032 assumed
δ	Recovery rate from quarantine	0.854302 assumed
γ	individuals rate lose immunity	0.0003 assumed
α	Rate of recovery	0.0968 [66]

Using the numerical values for initial data of compartments and parameters in Table 2, we plot the solution through NSFD method in Fig. 2, Fig. 3, Fig. 4, Fig. 5 of various compartments. Here we have taken the initial data of the Pakistan from [65], [66] and simulated the results for days.

Fig. 2

Graphical presentation of susceptible class.

Fig. 3

Graphical presentation of infected class.

Fig. 4

Graphical presentation of quarantined class.

Fig. 5

Graphical presentation of recovered class.

Here in Fig. 6, we compare our results for different compartment with that of RK4 method solution. We see that both solutions are closely agreed. As compare to RK4 method our proposed scheme NSFD is simple and does not involve complex calculation as in the RK4 method needs.

Fig. 6

Comparison of the NSFD method solution with that of the RK4 method results for the proposed model.

Numerical simulation of the system (1) under the concept of fractional calculus

This section of our paper contains numerical scheme for the proposed system (1) with the idea of fractional order derivative using NSFD method. To simulate the system (1), we take Grünwald–Letnikov approximation for Caputo derivatives as introduced in [60], [61]. Further, the numerical scheme we are using here has the ability to preserve dynamical properties. Also our proposed numerical scheme preserves the stability of the associated equilibrium point for a sufficiently small step-size . The suggested scheme has two main advantages over the traditional RK2, RK4, and Euler method which are mostly related with the time step increment. Also the said method behaves better or equivalently to RK 4 or RK2 methods. The complexity of the proposed method is comparable with Euler scheme. Hence the increase in terms of computation is significant for accuracy. Further, the computational complexity is much less than that of RKM methods. Moreover, the conservation of the dynamical constraints indicate to a scheme which yields a good dynamical behavior even for large time raise. This is a significant big application and computational advantage of the suggested method.

Let us suppose that , the Caputo fractional order derivatives is recollected as provided that integral on right is converging.

The Riemann–Liouville integral operator with fractional order is recollected as Following the procedure of [64], we have where . Here it is interesting if , we have from (22) Consider One has the following in the light of ((22)) to discretize (23) as we have where .

Keeping in mind that satisfies . The consistency rule is satisfied by some function like , etc. We present our proposed model (1) corresponding to fractional order as given by In view of (25) and by using Grünwald–Letnikov discrimination formula, we developed the required numerical scheme for the fractional order system (26) as In this part, we plot figures by using fractional order scheme by means of fractional NSFD scheme established above. Also take various fractional order and fixed values of , we simulate our results graphically. With fractional order, we get prediction for next few months from global dynamics of fractional calculus. If, we increase the order, the dynamical behavior of different classes go to the integer order.

From Fig. 7, Fig. 8, Fig. 9, Fig. 10, we see that at different fractional order the transmission dynamics behaves different. But as the fractional order goes to integer value 1, the concerned curve also tends to the curve at integer order. Further using the concept given in [67], the equilibrium point is locally asymptotically stable of the fractional order system if all eigen values of the concerned Jacobian matrix satisfies the condition given by , where . The results we established are compatible with theoretical analysis, which indicate that lower the fractional order for the system, the decay in susceptible class is faster while the growth in infection is slow. But greater the fractional order slower the decay and faster the growth. Hence choosing for suitable dynamics, we have greater degree of freedom in choosing order of the system. Therefore displaying the model for various fractional order giving accumulation which include its integer order dynamic as a special case. Further by using fractional order derivative, we remark that in physics, biology, electrochemistry and visco-elasticity, etc, the anomalous behavior of dynamical systems can be study very well. Further, to study dynamics of population densities in various phenomenon, exponential laws are traditional approach but in various systems the dynamics undergoes faster or slow. In such a circumstances where anomalous changes in the dynamics occurs, the use of fractional calculus is very best option. For instance the autonomous flow during diffusion process can be well explained via fractional order derivatives instead of classical. Due to these applications and usefulness of fractional order derivatives, we have investigated our established model under the mentioned concept. Hence we conclude that fractional calculus provides more significant results to understand the dynamics of infectious disease as compared to classical order derivative. Further here we present some real data plots for our given system to compare our numerical results with the real data reported in [68] in  Fig. 11. The confirmed cases in Pakistan per day reported from 1 March 2021 to 15th of September 2021 for 200 days are recorded as given in Table 3.

Fig. 7

Graphical presentation of susceptible class for different fractional order.

Fig. 8

Graphical presentation of infected class for different fractional order.

Fig. 9

Graphical presentation of quarantined class for different fractional order.

Fig. 10

Graphical presentation of recovered class for different fractional order.

Fig. 11

Comparison between real and simulated data at given different fractional orders for our proposed model.

Table 3

Reported cases in Pakistan from 1 March 2021 to 15th of September 2021.

Reported value	Reported value	Reported value	Reported value	Reported value	Reported value
4	3549	22037	89583	53333	9940
4	3735	23268	93233	52203	9356
5	3852	25609	97690	51057	16001
5	3902	26003	100324	50080	13706
5	4162	26230	104648	40242	13385
5	4150	27054	105087	29274	12464
5	4307	27904	106142	29626	11697
6	4362	29266	107733	27189	11542
15	4824	30503	107270	26191	10378
17	5143	31775	107023	25279	10446
18	5122	32578	107607	24983	10446
19	5660	34386	107460	24941	9940
19	6043	34642	107784	24912	9356
31	6742	36228	106023	24908	8739
51	7286	37657	106775	24935	8555
182	7703	38150	106361	24827	8585
245	8479	38900	108100	20597	8500
331	8925	39690	108466	19230	8553
439	9438	40358	103543	18253	8623
485	10103	40880	95388	17573	8633
629	10586	42687	95241	17548	8564
758	11058	44777	95219	17555	8512
856	11747	47607	94522	17588	8660
962	11996	50234	91408	17103	8883
1034	12380	53300	90358	16229	6020
1171	12900	56144	89250	16685	6234
1139	13818	59394	87345	16014	6477
1454	14498	63400	86770	16001	6545
1554	14814	57170	84234	13706	5291
1836	15716	60470	77418	13385	5546
1997	16370	75053	77360	12464	5979
2262	17574	78699	73536	11697	5786
2520	18003	83182	60234	11542	5582
2646	20267	79700	60234	10378	5525
2899	21587	84762	57668	10446
3058	22037	85321	53431	9940

Conclusion

A dynamical system addresses the transformation of SARS-Cov-2 has been considered under fractional as well as classical order derivative. We have first formulated the model. The mentioned model has been investigating for both global as well as local stability results. The Dulac function theory has been applied to prove the global stability of the proposed model. While the Routh–Hurwitz criteria has been exercised to derive results about local stability. The basic reproduction number has been calculated by the method of next generation matrix. Based on the value of , we can predict about the transmission dynamics of the disease weather it is increasing or decreasing. The value of for our model can be computed from the given values in Table 2 which is 0.127. Clearly , which indicate the future dynamics of the disease. It can be predicted from the derived value that the transmission dynamics of the disease goes out. Further on using NSFD scheme, we have simulated model by taking some real data of Pakistan. Also keeping in mind the importance of fractional calculus, we have extended our model to fractional order in Caputo sense. We have updated the classical NSFD scheme to fractional order by using Grünwald–Letnikov discretization formula for fractional order derivative. We have simulated the results against different fractional order using the aforementioned real data. Also some comparison with real data reported in Pakistan for the mentioned 200 days has been given. The continuous data graph and real data in case of infection are closely related which demonstrate the efficiency of the proposed method. Finally, we concluded that fractional order NSFD scheme worked very well for the said model in simulation. Also by using fractional order derivative to study dynamical behavior of real world problems is more realistic approach than integer order. As compared to classical simulation, we have noticed that fractional order simulations are more efficient and effective in understanding.

Funding

Currently no Funding source available.

CRediT authorship contribution statement

Kamal Shah: Established the theoretical part literature review. Thabet Abdeljawad: Supervised the project and edited the last version. Rahim Ud Din: Included table graphs and simulations.

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.