Mathematical analysis of COVID-19 via new mathematical model.

We develop a new mathematical model by including the resistive class together with quarantine class and use it to investigate the transmission dynamics of the novel corona virus disease (COVID-19). Our developed model consists of four compartments, namely the susceptible class, S ( t ) , the healthy (resistive) class, H ( t ) , the infected class, I ( t ) and the quarantine class, Q ( t ) . We derive basic properties like, boundedness and positivity, of our proposed model in a biologically feasible region. To discuss the local as well as the global behaviour of the possible equilibria of the model, we compute the threshold quantity. The linearization and Lyapunov function theory are used to derive conditions for the stability analysis of the possible equilibrium states. We present numerical simulations to support our investigations. The simulations are compared with the available real data for Wuhan city in China, where the infection was initially originated.

Introduction

Communicable diseases have always been an important part of human history. Large number of outbreaks have been occurred in human history due to which millions of people lost their lives. For instance the outbreak that occurred in the previous century, known as Spanish flow during which millions of people died. Many of such diseases then became endemic like HIV/AIDS due to which thousands of people are dying each year. In the last two decades, two corona-virus epidemics have been reported [1], [2], [3], [4], [5]. Such micro organism are classes of viruses that can result infections in humans, ranging from the common cold to severe acute respiratory syndrome (SARS). The first corona-virus epidemic, known as SARS has effected more than eight thousand infections with eight hundred deaths. The second one, called MERS, spread from Saudi Arabia to number of other countries. About 25,000 individuals were infected with it among which nearly thousand people lost their life. MERS is still a root of some cases [6]. Nearly a year ago, a severe respirational infection originated in the Wuhan city of China [7]. It was reported that the main cause of this infection was the new coronavirus (COVID-9). It is reported that the disease was first transmitted from animals to humans [8]. Some researchers say that bats have the corona virus from which it was transmitted to the humans, as all the patients that were first identified working in wet market in the Wuhan city [9]. The infection was then transmitted very rapidly in the whole city of during the months of January and February 2020. Due to international travel, in March 2020 some patients were identified in USA, Thailand and Korea. Then it was also announced that the disease is contagious as it is transmitted from person to person by contact. During the mid of Mach 2020, WHO announced that the infection to be an outbreak [10]. One virologists named it ‘severe acute respiratory syndrome’ corona virus (SARS-CoV-2). There are many claims and several theories about the origin of COVID-19 in which the most probable is that it might have been originated from bat [11], or from a seafood market exposure [12]. As the disease is contagious, international travel could expedite the spread of the virus [13]. The novel corona virus outbreak currently as the most remarkable in the recent history. By April 2020, COVID-19 affected nearly the whole globe [14]. Researchers in the area of epidemiology and other branches of biology are working day and night to find treatment or vaccine for this disease. They use different tools to understand the procedure through which the disease transmits in a society and how to reduce or control it. The process of infectious diseases may be easily understood and described by using mathematical models [15], [16], [17].

Mathematical model is a powerful tool that effectively helps in investigation of real world phenomenon and processes [18], [19], [20], [21], [22]. Bernoulli was the first mathematician who gave idea about mathematical modeling of spread of an infectious disease during 1760. After that numerous researchers took interest in the said area. One can easily understand various physical and biological phenomenon and their mechanism through such models. This area has been very well extended from simple models to more complex and complicated models. With the help mathematical models, large number of infectious and other diseases have been studied e.g., see [19], [20]. Using mathematical models, researchers first try to understand the dynamics of a disease, and afterwards they develop control and curing procedures for it. For some famous study in this regard, we refer [23], [24], [25], [26].

The researchers have used the tools of nonlinear numerical analysis to establish the global, local stability for the endemic and disease free equilibrium (DFE). In same line recently researchers have greatly investigated the novel COVID-19 through mathematical models from different aspects. The concerned investigations are devoted to stability theory, numerical simulation and global local dynamics. In this regards we have referred here some good work like [27], [28], [29], [30], [31], [32], [33].

In all previous models the researchers have considered various compartments but to the best of our knowledge the resistive compartment along with quarantined class has not been reported yet. Our model makes the investigation unique by this way that we have involved the resistive class together with quarantined class in construction of our new model.

The paper is organized as follows. In Section 2 we describe mathematical model for the current COVID-19 and discuss its fundamental mathematical properties in a biologically feasible region. We perform the stability analysis in Section 3 and derive the basic reproduction number. In the same section we find conditions for the local as well as the global behaviour of the possible equilibria of the proposed model. In Section 4 we present numerical simulations to support and verify our analytical findings. These simulations are performed for biologically feasible values of the parameters of the model. Finally we conclude our work in Section 5.

Formulation of the model

Inspired from the above mentioned literature, we develop a new mathematical model based on susceptible individuals healthy or resistant individuals infected and quarantine individuals respectively. All the parameters involved in the model are assumed to be non-negative. The susceptible individuals initially move to the infectious class with a constant flow rate. The suspected or infected individuals move to the quarantine class and confirmed cases are send back to the infected compartment for further treatment. Our model under consideration is expressed in the form of the following autonomous ordinary differential equationsThe description of parameters in (1) is given in Table 1 .

Table 1

Parameters and their explanation in the model (1).

Parameters	The physical interpretation
λ	Recruitment rate susceptible
γ	Disease transmission rate
d	Natural death rate
α	Recruitment rate of healthy human
β	Transmission rate of healthy human
μ	Disease related death rate infected or suspected individuals
δ	Rate at which quarantine people get infection
θ	Cure rate of infected people in the quarantine class

A systematic diagram of the model is given in Flow Chart 2 :

Fig. 1

Flowchart of the model (1) under consideration+.

We will investigate (1) under the following biologically feasible initial conditions

Let denotes the total population at time then we haveTaking the temporal derivatives of and using (1), one observes that obeys the law of mass actionNote that (3) is an exact differential equation with solutionHence it is deduced that (4) possesses positive definite solutions for all .

The dynamical system(1)exhibits boundedness for all non-negative initial conditions which are not all identically zero in the entire region given by

Assume that be any solution set of the model (1) with some non-negative initial conditions such ascorresponding to any other non-negative initial conditions on and . Since is a positive parameter, so that from (3) one may writeSolving this equation leads towhere is the initial value of the total population of the dynamical system. Thus, for we haveHence is positive and bounded, and is the largest set for which the solutions are positive and bounded. This completes the proof of the theorem. □

Stability analysis

We will focus our attention on determining the possible stationary states of the system (1) and on deriving stability results. We consider the situation, when there is no infection of the disease in the community. The concerned state is called DFE. In the following, we denote such equilibrium by .

DFE state

Inserting in the given system and solving the autonomous differential Eq. (1) for and we obtain

To calculate the endemic equilibrium, we first need to determine the threshold quantity which plays a significant role in determining the global behavior of the given dynamical system.

Basic reproductive number

We follow the traditional technique, the next generation method [34] to derive . Let denotes the infectious class in model (1). We can writeThe Jacobians of the above matrices are respectively given byThe multiplicative inverse of the matrix is calculated as

The next generation matrix for the DFE of the proposed problem is given byThe spectral parameter of this matrix gives the threshold parameter . Thus is given byIn Fig. 2 , we plot against different parameters involved in the model under consideration.

Fig. 2

Plot of the basic reproduction number in terms of various parameters involved in the model under consideration.

Local stability at the DFE state

After determining the threshold quantity, we are now in position to derive stability conditions for our model.

The DFE of the dynamical system(1)is locally asymptotically stable equipped the threshold quantityandThe DFE is unstable forand in the contrast of.

The local stability of for the proposed model (1) can be inquired by taking the Jacobian matrix of the formExploiting elementary row operations, the last matrix yieldswhereObviously the matrix (10) has eigenvalues . The remaining two eigenvalues and are negative if and only ifHence re-arranging the inequality (11) implies that Thus condition holds along with which shows that the system (1) exhibits the local asymptotical stability. On the other hand, we see that if has non-negative real part showing instability of the given system with growing time. This completes the proof of the theorem. □

Global stability of the DFE state

In this section, we discuss the global stability of the disease free state. In the following, we derive conditions for the global stability of the disease fee equilibrium using Castillo-Chaves [35]. The global stability of the DFE depends upon the following two conditionswhere represents the susceptible and healthy classes, i.e., while denotes the infected and quarantine individuals, at the DFE point Thus the existence of the global asymptotic stability of depends upon the following conditions.

If then is globally asymptotically stable;

where for .

In the condition is an M-matrix having non-negative off diagonal entries, where represents the feasible region in biological sense. Thus the following Lemma holds.

Ifthen the equilibrium pointof the dynamical system(1)is globally asymptotically stable whenever conditionsandare satisfied.

Forthe system(1)exhibits global asymptotic stability at the DFE point

In order to prove the above theorem, we need to verify conditions and Let us use the symbols and and define where . With the use of system (1), we may writeandNow if and Thus we conclude that whenever . From which one can say that is globally asymptotically stable. Similarly for the second condition, we need to expresswhere is the Jacobian matrix of infected and quarantine classes at and Therefore, one may verify that matrix is of the formBy taking

we arrive atThe total population of dynamical system is bounded by and such that from which it follows that . Thus is positive definite. Furthermore, the off diagonal elements of the M-matrix are non-negative. Thus conditions and are satisfied. So by Lemma stated above, the DFE is globally asymptotically stable. □

Endemic equilibrium point and backward bifurcation

In this subsection we determine the endemic equilibrium point of the model (1), i.e., when the infected individuals of the system are non-zero. We take the following steps. Let represent an arbitrary endemic equilibrium point of the dynamical system (1). Solving simultaneously equations of the model(1) at steady state givesObviously for the endemic state, we have Substituting and in the third equation of the model (1) at stationary state, one obtains the following quadratic equationwhere and Surely the coefficient is always positive, as the parameters and are positive. Moreover is negative (positive) for is greater (less) than unity. Thus (19) depends upon the signs of and respectively to get its positive solutions. Now if is greater than unity, then two roots of (19) are positive real values and hence the endemic equilibrium state is unique. If is one, then is zero and hence there is either a trivial or no endemic state. The whole discussion shows that the endemic equilibrium of the model under consideration depends on . Thus there exists an interval such thatIf is non-negative and either or then there have no positive solution of (19) and have no endemic equilibrium state. Finally, we establish the following result for various range of parameters.

The model(1)has the following set of facts

Ifthen the model(1)have unique state of endemic equilibrium;

Ifandorthen the model(1)have unique state of endemic equilibrium;

Ifandthen the model(1)have two unique endemic equilibria;

Otherwise the model(1)has no endemic equilibria.

Case 3 in Theorem (3.4) illustrates the existence of phenomenon of backward bifurcation, where the local asymptotic stabilities of both DFE as well endemic equilibrium co-exist, when [36], [37]. Now to probe the backward bifurcation, we must set the discriminant to zero in order to obtain critical value of . We obtainThus backward bifurcation occurs, when or equivalently The phenomenon of backward bifurcation and its epidemiological significance needs a classical requirement of although this condition is necessary but no longer sufficient to eliminate the disease. In same scenario,the disease elimination depends upon the initial state variables or initial size of the sub population of the system. The backward bifurcation in the same model (1) suggests the feasibility of controlling the disease when threshold quantity is greater than one and would depend on the initial stature of the sub-population of the system (1).

The dynamical system(1)endures backward bifurcation when case 3 of theorem (3.4)holds along with

Atthe system(1)undergoes backward bifurcation if and only if

For the sufficient part, assume the graph of Now if implies that and hence . Consequently graph of the function passes through origin. Furthermore, has a positive root if If we exceed from zero to some value . This guarantees that there exists an open interval which contains where has two real and positive roots. In short words, we try to show if then there are two endemic equilibrium states. The necessity is conspicuous, if (19) has no real positive roots, when  □

Local stability of endemic equilibrium state

The following theorem is enough to show that the system(1) performs local stability at the endemic equilibrium point .

The unique endemic equilibrium stateof model(1)is locally asymptotically stable provided.

The Jacobian matrix of the proposed model (1) for unique endemic equilibrium state is given byTo know about the nature of the eigenvalues of the matrix (22), it is advantageous to apply the traditional row operations to obtain

where The eigenvalues of are and In addition if and only if After using simple algebra we can rewriteClearly all the coefficients of (24) are positive if and .Thus conditions of the above stated theorem are satisfied and consequently result in the local asymptotic stability of the endemic equilibrium . □

Global stability of endemic equilibrium state

In this section we explore the global stability of the endemic equilibrium point of our propose model (1)in terms of the basic reproduction number .

The endemic equilibrium stateof model(1)is globally asymptotically stable, ifotherwise unstable.

To study the global stability of the model(1) under consideration we construct the following Lyapunov function.Clearly at the endemic state the function for the prescribed values of and . Further the given function is strictly positive, whenever and . This implies that the function for the given is positive semi-definite.

The time derivative of the function results inInserting values from (1) into (26) and performing algebraic manipulation, (26) can be formatted asClearly for and This mean that the time rate derivative of the defined function is negative semi definite. From which we have that is the largest invariant sub set in support that the endemic equilibrium state of the dynamical system (1) is globally asymptotically stable under the condition . The graphical investigation of the above theorem is illustrated in Fig. 2. □

Numerical results and discussion

Here we use RK4 method to perform the numerical simulations using some real data of as in the following Table 2

Table 2

Description of the parameters used in model (1).

Parameters	The physical interpretation	Numerical value
S(t)	The susceptible population (recruit for test)	1000 thousands
H(t)	The resistant population	790 thousands
I(t)	The infected population	170 thousands
Q(t)	The quarantined population	450 thousands
λ	Recruitment rate of susceptible	0.0043217
γ	Disease transmission rate	0.125
d	Natural death rate	0.002
α	Recruitment rate of healthy human	0.535
β	Transmission rate of healthy human	0.0056
μ	Disease related death rate infected or suspected individuals	0.0008
δ	Rate at which quarantine people getting infection	0.029
θ	Cure rate of infected people in quarantine	0.35
η	Contact rate of infected and healthy people	0.025

In the first hundred days in a locality, one thousand individuals were found susceptible, among which 170 were founded infected and 790 were declared healthy. 450 individuals either died or got rid from the infection. The infected individuals were quarantined. We simulated the new constructed model under these information and found that in coming few months the infectious will grow exponentially and hence more people will be quarantined if the people and government of a locality do not properly follow SOPs. As a result the resistive population decreases and this causes increase the classes and . It can be observed from Fig. 3 that susceptibility is decreasing and as a result healthy population also declines, as shown in Fig. 4 . Hence the infected and the quarantined classes are growing up, see Figures 5 and 6 respectively. The asymptotic stability behavior is clearly observed from these plots. Next, using available data of Wuhan city (China) and using the proposed model, we compared the real and simulated data for infected class as reported in [38], [39] for initial sixty eight days from 4th January 2020 to 8th March 2020 as The concerned parametric values in the model are taken as Where for compartments, we use the initial conditions We compare the simulated and real data in Fig. 7 . We rather see a very good agreement between the simulated and real data.

Fig. 3

Dynamical behavior of the susceptible class (population).

Fig. 4

Dynamical behavior of the resistive class (population).

Fig. 5

Dynamical behavior of the infected class (population).

Fig. 6

Dynamical behavior of the quarantined class (population) (population) for the considered model (1).

Fig. 7

Comparison between simulated data and real data for the infected class using the model (1).

Conclusion

We have established a new model by considering the resistive class together with quarantine class for the transmission dynamics of COVID-19. Exploiting the Lyapunov function theory, we have developed sufficient conditions for global and local stability of the disease free and the endemic equilibria of our model in terms of the threshold quantity. We have shown the positivity and boundedness of solutions of the model in a feasible region. To support and verify our analytical work, we have performed numerical simulations. This has been achieved using the RK4 method by taking some real data of the city Wuhan in China.

It is well known fact that fractional analysis is a hot area of research in the recent time. Fractional calculus has been found very effective in modelling various real world phenomenon. Researchers have focused their attention on studying models of infectious disease like HIV, AIDS, etc under fractional order derivatives and integrals. For more detail we refer [40], [41], [42], [43], [44]. Some authors have also extended the area to classical and arbitrary order mathematical models of real world problems in physical sciences, see, for instance, [45], [46]. Our next step is to investigate the qualitative and numerical aspects of our proposed under different fractional order derivatives. Work on this is in progress and will be reported in a future publication.

CRediT authorship contribution statement

: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Software, Writing - original draft, Writing - review & editing. Saeed Ahmad: Conceptualization, Data curation, Investigation, Methodology, Software, Writing - original draft, Writing - review & editing. Saud Owyed: Data curation, Formal analysis, Methodology, Validation, Writing - review & editing. Abdel-Haleem Abdel-Aty: Data curation, Funding acquisition, Methodology, Project administration, Resources, Supervision, Validation, Writing - review & editing. Emad E. Mahmoud: Formal analysis, Funding acquisition, Resources, Software, Validation, Writing - review & editing. Kamal Shah: Conceptualization, Formal analysis, Investigation, Project administration, Supervision, Visualization, Writing - original draft, Writing - review & editing. Hussam Alrabaiah: Formal analysis, Investigation, Software, Visualization, Writing - review & editing.

Declaration of Competing Interest

I declare that I have no significant competing financial, professional, or personal interests that might have influenced the performance or presentation of the work. Also, we do not have any conflict of interest regarding this paper.