Stability analysis and optimal control of Covid-19 pandemic SEIQR fractional mathematical model with harmonic mean type incidence rate and treatment.

In the present work, we investigated the transmission dynamics of fractional order SARS-CoV-2 mathematical model with the help of Susceptible S ( t ) , Exposed E ( t ) , Infected I ( t ) , Quarantine Q ( t ) , and Recovered R ( t ) . The aims of this work is to investigate the stability and optimal control of the concerned mathematical model for both local and global stability by third additive compound matrix approach and we also obtained threshold value by the next generation approach. The author's visualized the desired results graphically. We also control each of the population of underlying model with control variables by optimal control strategies with Pontryagin's maximum Principle and obtained the desired numerical results by using the homotopy perturbation method. The proposed model is locally asymptotically unstable, while stable globally asymptotically on endemic equilibrium. We also explored the results graphically in numerical section for better understanding of transmission dynamics.

Introduction

The human society recently facing a terrible enemy is known as corona virus now a days. One of the challenging issues to the human community has to recognize the complex dynamic of Covid-19. The concern virus badly affected almost each countries and traceries around the globe. Corona virus was for first time notified in the peoples of Wuhan Chain [1]. Scientist provides different theories concerning the origination of consider virus, but it is still a mystery that from where this virus is took place. Initially, some cases of affected peoples were reported of this virus from local fish market in city of Chain Wuhan [2]. In this connection some of the scientist believed, that this virus was transformed in humans form animals. Similarly, the researches recognized that the concern virus could transmit from person to person too [3]. The latest statistics provided by World Health Organization (WHO) on January 2021, the concern virus almost affected each and every traceries and countries on the globe. According to fresh statistics of WHO [4], at 7:08pm CET, on 15th January 2021, there were 709,865 new cases reported globally, while 91,816,091 confirmed cases, including 1,986,871 deaths. The concern virus badly affected the developed countries especially UK, Spain, Italy, USA, Italy, and many more [5]. The death rates in aforementioned regions are very higher as compared to other countries. The aforesaid analysis justified the severity of proposed virus. The well known symptoms of this disease are sever coughing, regular fever and infection include respiratory issues. In addition, to this neurological sickness and gastroenteritis of contradictory strictness are also counted as symptoms of proposed virus [6].

The main source of transmission of infection is the droplets from the nose or mouth of the influence person during sneezing, coughing or speaking. Therefore, a person among the affected peoples has it high risk of been affected form the disease. As precautionary, all the traceries and countries around globe implement the police lock-down in their respective countries, in order to ensure the safety of peoples. In these circumstances paramedical staff and doctors has dedicated themselves to give health services to the affected humans. The experts of concerned area believes that root cause of consider virus, was first resulted due the bats, that is identical to SARS (Severe Acute Respiratory Syndromes), that took birth in China in the year 2003, (see [7], [8]). Some of the researchers match up the Covid-19 with SARS and MERS to categorize the family of virus from which it belonging. In [9], author’s presented that the current virus recounted to the genus of -corona virus, like SARS-Cov and MERS-Cov. For further, details see [10], [11], [12].

In this paper our aim is to study the dynamics of Covid-19 based on a epidemiological SEIQR model, where they represent respectively the susceptible, exposed, infected, quarantine, and recovered(removed) human population and with harmonic mean type incidence rate (see, [16], [17], [18]). Particularly, in this study the incidence rate is

Methodology and mathematical model formulation

Several models are introduced to simulate the dynamics of the spread and transmission of COVID-19 [19], [20], [21], [22], [23], [24]. The model, epidemiological SEIQR model is specially design for Covid-19 disease, because this model containing the quarantine compartment, which describe the proper way to decreases the concerning infection disease, that other model don’t have. The proposed model is defined in five compartment susceptible population, exposed population, infected population, quarantine population, recovered population and all of the compartments relay on time (t). In the SEIQR model have a contract between susceptible population with exposed population from which the disease is transfer from exposed into susceptible hosts . In model (1), b is the recruitment or the birth rate, is the transmission rate, susceptible population recovers at a rate of and quarantine at , while is the natural mortality rate. The parameters , and are incubation period of infected, exposed and susceptible compartment respectively, at which the population goes into the quarantine compartment . Infected and quarantine population recover at a rate of and , while is the rate of infection in exposed population and the infected population dies out at a rate of . withand

The parameters involved in SEIQR system is described in the following Table 1 :

Table 1

Description of parameter and its value.

Notation	Description of parameters
b	Fertility Rate.
β	Transmission Rate from Susceptible into Exposed Compartment.
η	Transmission Rate of Susceptible to Recovered Compartment.
μ	Natural mortality Rate.
q_1	Incubation Period of Infected Compartment.
q_2	Incubation Period of Exposed Compartment.
q_3	Incubation Period of Susceptible Compartment.
λ	Infection Rate of Exposed Compartment.
ε	Death Rate in Infected Compartment due to Disease.
γ	Recovery Rate of Infected Compartment.
τ	Recovery Rate of Quarantine Compartment

Basic reproduction number

We obtained the basic reproduction number by considering the functional equations

Thusandwhere the next generation matrix is given as

The eigenvalues of the next generation matrix arewhere , thus the basic reproduction number, is the spectral radius G is defined assuch that

Equilibrium points and their local stability

This section includes the possible fixed points of model (1). There exists two possible equilibrium points are calculated, i.e. Disease free equilibrium (DFE) and endemic equilibrium (EE). Furthermore, basic reproduction number is calculated by next generation technique and discuss the local stable analysis of these equilibrium points. We denote the infection-free equilibrium point by such that .

The infection-free equilibrium point of model (1) is locally asymptotically stable, if , otherwise unstable.

The Jacobian matrix of model (1) around the infection-free equilibrium point becomeswe easily have three eigenvalues while in matrix (10) and it takes the form

Now after some matrix operation in matrix (11), we have

We getting the rest of two eigenvalues from the matrix (13) these are: which would negative for , such that , while which is positive. Hence the model is unstable at infection-free equilibrium point for . □

The infection-endemic equilibrium point is denoted by and we have . where , and .

The infection-endemic equilibrium point of model (1) is locally asymptotically stable, if , otherwise unstable.

The Jacobian matrix of model (1) around the infection-endemic equilibrium point is given as

Clearly, we get the first two eigenvalues of matrix (17), , and which are of course negative, then the matrix (17) would becomemultiplying row second by and row third by and add row second and third, we get

Thus, the eigenvalue is . Next, we subtract row first from row second we havesubtract row second from the product of and row first, we have

Now after some matrix operation, one has

So that the eigenvalues of the matrix (22) are , and . We concluded that three of its eigenvalues are negative while the remaining is non-negative. Thus the model is unstable at disease-endemic equilibrium point, if and only if . □

Global stability analysis

In this section, we investigate the Global stability of the model (1) at disease-endemic equilibrium point.

If , then the model (1) is globally asymptotically stable at endemic equilibrium , else unstable.

Consider the non-linear equations of the model (1), and we define the functional equations for those equations,

The jacobian of system of Eqs. (23) at disease-endemic equilibrium, such that

Next, the general form of the third additive compound matrix is given by

Therefore, matrix (24) and matrix (25) yields thatwhereand

Consider , such that , and the time derivative of is defined as , thereforeand

Hence, one has , that is

Consequently, we have to find , such that from Eq. (30) to Eq. (33) implies that

Similarly for , one has

For , we have

For , we get

Let be a vector in and the Lozinski measure of B is defined as . The integration of Lozinski measure by taking limit as .

Similarly, for , we have

With same fashion, one has and similarly can be obtained as

These lead us to the following relation

The system containing only four non-linear equations of model (1) is globally asymptotically stable around its interior equilibrium . The solution of the limiting system of the remaining three equations of model (1) gives that approaches to its equilibrium point. Hence the disease endemic equilibrium point is globally asymptotically stable. □

Optimal control strategies

Consider the model (1) with modifications

Subject to initial conditions

Using the modified model (43) we define the following optimal control problem to maximize the objective functional, such thatwhere

The variable, is the control variable [14], which show the education or media campaign for the awareness of preventions (may be the lockdown) defined by the public health department;

The variable is the treatment of Susceptible individuals;

The variable is the treatment of Exposed individuals;

The variable is the treatment of Infected individuals;

The variable is the treatment of Quarantine individuals.

When they are vary from 0 to 1 then the efforts of campaign about the Covid-19 increases and also the clinical treatment also increases. Whereas, the chances of infection reduce too but when the control variables reduce towards “0” then the model (43) approaches to the original model (1) where is no campaign effort and no treatment. Where are the weight parameters for balancing. In order to measurement the control variables such that we find the control function as

Subject to the control model (43), (45), such that the control set is given by

We prove the existence of control variables to find them.

Existence of optimal control problem

Let us consider the control problem (43), (44), such that it is clear that there exists for the positive bounded solutions for the state system, if the initial conditions are non-negative and control are bounded Lebegue measurable. We established the following theorem for the existence of control problem.

There exist and optimal control

Subject to the control problem (43), (44) .

Since, the combination controls and state variables are both non-negative and non-empty, and the set of controls U is closed and convex which also satisfying the property of boundedness leads to proof of compactness required for the optimal control existence. The integrand of the objective functional (45) is also convex therefore the final condition iswhere, , while . □

Optimality conditions

We consider the optimal control problem (43)–(47) for the optimal solution and define the Lagrangian as well as Hamiltonian for the control model with initial conditions (43), (44) such that the Lagrangian and Hamiltonian are defined bywhile the Hamiltonian for the problem (43), (44)is

Equivalently, one haswhich yields that

We use the well-known Pontryagin’s Maximum Principle for the optimal solution of the proposed model (43), (44). When the optimal solution is essentially bounded for the control problem with control variables , then there exists a nontrivial vector function , the Hamiltonian system isthe optimality conditionand the adjoint equation

The optimal control variables and the solutions of the state system, then we have the adjoint variables , satisfying with transversally conditions

Furthermore, the optimal control variables can be shown as following

We establish the adjoint equations and the transversally conditions by Hamiltonian, H. Taking the derivatives of H with respect to the state variables and at endemic-equilibrium point, , such that with transversality conditions and the characteristic equations of the control space U, that is

So that the control variables are obtained as

This completes the proof. □

Numerical results

Homotopy perturbation method (HPM) (see, [13], [15], [25], [26]) used to find the semi-analytical solution of differential equations, such as ordinary differential equation(s), partial differential equation(s) and fractional order differential equation(s). In this study we also find the solution to the model (1) for our results. Thus the scheme is developed as following, . Consider the homotopy for the model (1), such that

We assume solution for the model (1) as

For the comparison of , system of Eqs. (75), (76) implies that for and similarly for , such that

First order problem

Second order problem

taking limit as in system of Eqs. (76), we have

Hence, we have semi-analytic solution for model (1) based on the system of Eqs. (86).

Solution by homotopy perturbation method

where, , and .

Discussion

In current work, we studied the transmission of SARS-CoV-2 with fractional order Susceptible, Exposed, Infected, Quarantine and Recovered Population model by means of stability, control and numerical interpretation by homotopy perturbation method. The basic reproduction number, is calculated by next generation matrix and examined the model on both disease-free and endemic equilibrium points which is locally asymptotically unstable and globally asymptotically stable. We used a third additive compound matrix for the computation of global stability and also implemented five control variables in fractional optimal control, the variable, is the control variable, which show the education or media campaign for the awareness of preventions (may be the lockdown) defined by the public health department, the variable is the treatment of Exposed individuals, the variable is the treatment of Exposed individuals, the variable is the treatment of Infected individuals, the variable is the treatment of Quarantine individuals. We have presented graphically the semi analytical results of various compartments in Fig. 1, Fig. 2, Fig. 3, Fig. 4, Fig. 5 . The respective dynamical behavior of the corresponding compartments against various fractional order has been shown. Also in Fig. 6, Fig. 7 , we have presented the behavior of basic reproductive numbers with and with contour plots against the given values of the parameters. Under the control sterility, the plot of is bonded bellow 1 which shows that stability is occurring in the dynamics with the passage of time. On the other hand in the Fig. 8, Fig. 9 , we have presented the plots of different compartments under various control parameters. We see that a decay is occurring in susceptible, exposed and infectious classes and the class of quarantined is raising due to the isolation of individuals among each other. Further the global dynamics of various compartments for different fractional order can be observed from Fig. 1, Fig. 2, Fig. 3, Fig. 4, Fig. 5. The different behavior in dynamics interpret that fractional calculus approach is an excellent way to investigate biological model of infectious disease (see Table 2, Table 3 and Fig. 10, Fig. 11, Fig. 12 ).

Fig. 1

Numerical simulation of susceptible and exposed human individuals by homotopy perturbation method for thirty days; in these plots the order of is arranged from bottom to top from  = 0.1 to 1.0.

Fig. 2

Numerical simulation of infected and quarantine human individuals by homotopy perturbation method for thirty days; in these plots the order of is arranged from bottom to top from  = 0.1 to 1.0.

Fig. 3

Numerical simulation of recovered human individuals and the behaviour of susceptible human individuals with the variations in , while ; in these plots the order of is arranged from bottom to top from  = 0.1 to 1.0.

Fig. 4

The behavior of exposed and infected human individuals with the variations in , while ; in these plots the order of is arranged from bottom to top from  = 0.01 to 0.1.

Fig. 5

The behavior of quarantine and recovered human individuals with the variations in , while ; in these plots the order of is arranged from bottom to top from  = 0.01 to 0.1.

Fig. 6

The behavior of basic reproduction number, with the variations in and .

Fig. 7

The contour plot of and to basic reproduction number, .

Fig. 8

This plots show the behavior of susceptible with the control variable the education or media campaign for the awareness of preventions, and the treatment in both Susceptible and Exposed population.

Fig. 9

This plots show the behavior of infected and quarantine population with the control variable treatment.

Table 2

Descriptions of classes.

Classes	Description of classes
S(t)	Susceptible Population at time, t.
E(t)	Exposed Population at time, t.
I(t)	Infected Population at time, t.
Q(t)	Quarantine Population at time, t.
R(t)	Recovered Population at time, t.

Table 3

Description of parameter and its value.

Notation	Description of parameters	Values
b	Fertility Rate.	3.42500
β	Transmission Rate from Susceptible into Exposed Compartment.	0.01898
η	Transmission Rate of Susceptible to Recovered Compartment.	0.1
μ	Natural mortality Rate.	0.00576
q_1	Incubation Period of Infected Compartment.	0.02200
q_2	Incubation Period of Exposed Compartment.	0.03100
q_3	Incubation Period of Susceptible Compartment.	0.03900
λ	Infection Rate of Exposed Compartment.	0.00147
ε	Death Rate in Infected Compartment due to Disease.	0.020
γ	Recovery Rate of Infected Compartment.	0.94978
τ	Recovery Rate of Quarantine Compartment	0.04990

Fig. 10

These plots show the behavior of recovered population with the control variable treatment and stability curves (Fig. 1, [27]) of susceptible population.

Fig. 11

Stability curves of exposed and infected population.

Fig. 12

Stability curves of quarantine and recovered population.

CRediT authorship contribution statement

Muhammad Sinan: Data curattion, Formal analysis, Methodology, Software, Writing-original draft, Writing - review & editing. Amjad Ali: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Visualization, Writing - original draft, Writing - review & editing. Kamal Shah: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Supervision, Validation, Visualization, Writing - original draft, Writing - review & editing. Taghreed A. Assiri: Formal analysis, Methodology, Validation, Visualization, Formal analysis, Writing - review & editing. Taher A. Nofal: Data curation, Funding acquisition, Investigation, Methodology, Resources, Software, Validation, Visualization, Writing - original draft, Writing - review & editing.

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.