STABILITY ANALYSIS OF A NONLINEAR MATHEMATICAL MODEL FOR COVID-19 TRANSMISSION DYNAMICS

. The whole world had been plagued by the COVID-19 pandemic. It was ﬁrst detected in the Wuhan city of China in December 2019


INTRODUCTION
One of the most disastrous pandemics the world has witnessed throughout its history is the COVID-19 epidemic. On March 11, 2020, the World Health Organization (WHO) declared COVID-19 a pandemic. The vast majority of (COVID- 19) cases are transmitted through humanto-human contact. The virus is spread through direct contact with an infected individual. There have been more than 600 million cases of COVID-19 infection and more than 6 million deaths worldwide [1,2,3]. It is crucial to note that the number of cases is underestimated for a variety of factors, including a lack of diagnostics and asymptomatic cases [4,5,6,7,8,9,10,11,12].
The first case of COVID-19 infection was reported in the city of Wuhan, Hubei province of China [13,14,15,16,17]. Thereafter, the virus spread worldwide, affecting nearly each and every country.
The contribution of epidemiological modelling to the study of infectious disease transmission dynamics has been significant [18,19,20,21,22]. Modeling the dynamics of smallpox in 1760 marked the beginning of the study of epidemic dynamics, which has since become an essential tool for studying the spread and control of infectious diseases [23]. Kermack and McKendrick presented the the well-known Susceptible-Infected-Removed (SIR) model of ordinary differential equations in their ground-breaking paper [24]. In the 1950s and early 1960s, Bailey [25] and Bartlett [26] developed stochastic theories of disease dynamics. Bartlett also pioneered the application of Monte Carlo simulations to the study of epidemics [28,27].
Since the emergence of the COVID-19, several writers have contributed to the literature of epidemiological modelling, e.g. see [55,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54]. The classical SIR compartmental model splits the whole population into three compartments: S(t), I(t), and R(t), representing the proportion of the population susceptible to infection, infected individuals, and removed people (recovered or dead), respectively [56]. However, for the majority of infectious diseases, a latent phase precedes the transition from the diseased to the infectious state. This necessitates an additional compartment in the model, namely the exposed population E(t), making it a system of four ODEs. [57]. Effective use of SEIR has been made to comprehend the early dynamics of the COVID-19 outbreak and to assess the efficacy of various actions since the outbreak [58,59,60,61,62]. In [63], the classical SEIR model was expanded to include delays in order to incorporate the incubation period into the COVID-19 dynamics.
In this work, we extend the SEIR model by considering the class of quarantined individuals Q(t), and those individuals who are under treatment T (t). These classes are not infectious since they are under proper isolation. This is important since deterministic models with fewer compartments such as SIR or SEIR, fail to adequately characterize the evolution of COVID-19, as proved in [64]. Using dynamical system theory, this study aims to analyze the dynamics of the COVID-19 epidemic using the SEIQTR model. The discrete model is equally accurate, but we rely on the continuous deterministic model because it is easier to process. We will first evaluate the local stability of the model under both disease-free and endemic equilibrium conditions, and then we will test the model's global stability.
The rest of this paper is organized as follows. In Section 2, we formulate the mathematical model for the spread of COVID-19 infection, in Section 3, we present some basic properties such as positivity and boundedness of the model. The next section is committed to the derivation of the basic reproduction number. Section 5 is devoted to the stability analysis of the model including both local as well as global stability. Section 6 and 7 contains bifurcation analysis and sensitivity analysis of the basic reproduction number, respectively. In Section 8, we study some numerical simulations of our model substantiating results obtained in the previous sections. The last section is dedicated to discussion and conclusion.

FORMULATION OF THE MATHEMATICAL MODEL
We construct a compartmental model [65] for the spread of COVID-19 based on a continuous time deterministic nonlinear system of the differential equations. For the purpose of the model, we divide the whole population into six compartments, viz., the compartment of people susceptible to the disease S(t), people exposed to the disease who are not yet infectious E(t), the compartment of infected and infectious people I(t), the quarantined class Q(t), those are under treatment T (t), and the class R(t) of people recovered from COVID-19. Each compartment is assumed to be homogeneous. So all individuals within the same compartment are subject to the same hazards. We assume that once someone is recovered from the disease, she or he becomes immune. All newborns are assumed to be susceptible. At any time t, the total population is given by N(t) = S(t) + E(t) + I(t) + Q(t) + T (t) + R(t).
We consider the following system of non-linear differential equations for our model: with initial conditions Here, Λ is the birth rate, d is the natural death rate, β represents the effective contact rate of susceptible population with the infected population, α represents the proportion of exposed class who become infective, λ is the rate of quarantine, h is the rate at which people are being treated, γ is the recovery rate of a quarantined individual, δ represents the disease-induced death rate, and ρ represents the recovery rate due to treatment.
The diagrammatic representation of the mathematical model (1) is shown in Figure 1.

BASIC PROPERTIES OF THE MODEL
By the fundamental theory of differential equations [66,67], the solution of the system (2.1) with the initial conditions (2.2) exists for all t ≥ 0 and it is unique. For the model to be realistic, the solutions must be non-negative and uniformly bounded. We check for that in this section. Proof. From the first equation of the system (2.1), we have

Positivity of solutions.
Multiplying equation ( This implies that Integrating equation (3.2) from 0 to t, we get As every parameter is positive, equation (3.4) implies ∀t ≥ 0, Similarly, from the other equations of the system (2.1), we find Thus, all the solutions S(t), E(t), I(t), Q(t), T (t), and R(t) of the system (2.1) are positive for all t ≥ 0. This completes the proof of Theorem 1. Proof. Adding all the equations of the model (2.1), we obtain

Boundedness of solutions.
Since Using Gronwall's inequality [73] for (3.12) we get for t ≥ 0, N(0) being the total initial population. Hence, is positively invariant and the solutions of the model (2.1) remain bounded. Further, equation From inequality (3.15), we have N(t) approaches Λ d asymptotically for N(0) > Λ d . Thus, all feasible solutions of the system (2.1) are uniformly bounded on R 6 + .
By virtue of Theorem 1 and Theorem 2, for the purpose of analyzing the model (2.1), it suffices to consider the region Ω as given by (3.14).

THE BASIC REPRODUCTION NUMBER R 0
The basic reproduction number R 0 , is defined as the expected number of secondary cases that would arise from a typical primary case in a susceptible population [68,69,70,71,72].
We calculate the basic reproduction number from the first principle [70,71]. On an average, an infected person remains infectious for a period of 1/(d + λ ), the probability of the index case becoming infective rather than dying while in the class E is α/(α + d) which will transmit infection at a rate β × Λ d , so the expected number of secondary cases caused by the index case in a completely susceptible population is: probability of making it through the latent stage without dying × rate of transmission while infectious × average infectious period. Thus, The effective reproduction number, given by R 0 × S(t) N(t) , is the expected number of secondary infection at time t [71].

STABILITY ANALYSIS OF THE MODEL
In this section, we find the steady-state solutions of the model and study their stability. We will see that the stability of the equilibria is intricately related to the value of R 0 .

Equilibrium Points.
We find the equilibrium points of the model by setting all the time derivatives in (2.1) to be zero. We see that there are two possible equilibria-one corresponding to having disease and the other one to no disease at all.

Local stability.
We now proceed to study the local stability behavior of the COVID-19 disease-free equilibrium DFE and endemic equilibrium DPE derived in (5.1) and (5.2), respectively.

5.2.1.
Local stability of the disease-free equilibrium.

Global stability.
We now move on to examine the global stability of the disease-free equilibrium DFE and endemic equilibrium DPE given by (5.1) and (5.2), respectively. In order to show that the system (2.1) is globally asymptotically stable, we shall make use of the Lyapunov second method and LaSalle's principle.
Theorem 5.3. The disease-free equilibrium DFE as given by (5.1) of system (2.1) is globally asymptotically stable in Ω if R 0 ≤ 1.

5.3.2.
Global stability of the disease present equilibrium. For the global stability of the endemic equilibrium, we use the following type of function: which is non-negative for x > 0, and f (x) = 0 iff x = 1. Besides, We now present the final result on global stability.
Theorem 5.4. The endemic equilibrium point DPE of the system (2.1) given by (5.2) is globally asymptotically stable in Ω if R 0 > 1.

BIFURCATION ANALYSIS
We observe that there is an exchange of stability between the disease free equilibrium DFE (5.1) and the disease present equilibrium DPE (5.2) when R 0 crosses the threshold value 1. So, the system (2.1) passes through a bifurcation at R 0 = 1. This is the subject of our ensuing theorem.
Proof. From stability analysis in Section 5, we have found that if R 0 < 1,only the disease-free equilibrium DFE(5.1) exists, and is stable both locally and globally. Further, if R 0 > 1, the disease prsesnt equilibrium DPE(5.2) arises and it is locally as well as globally asymptotically stable. Moreover, even though the disease-free equilibrium DFE(5.1) continues to exist for R 0 > 1, it is unstable. Thus, at the threshold R 0 = 1, there is a change in both feasibility and stability of system (2.1). Following the study in [78,77,79,80,81], we therefore come to the conclusion that at R 0 = 1, the system (2.1) undergoes a transcritical bifurcation.
We graphically depict the bifurcation diagram of the system (2.1) in Figure 2.

SENSITIVITY ANALYSIS OF R 0
Sensitivity analysis measures the relative change in a state variable with the change of a parameter. The normalized forward sensitivity index of a variable to a parameter is the ratio of the relative change in the variable to the relative change in the parameter.
Following Chitnis et al. [82], we calculate the normalized forward sensitivity indices of R 0 .
This will allow us to look at the relative impact of the parameters on the reproduction number Let denotes the sensitivity index of R 0 with respect to the parameter p.
From (4.1) we have Therefore, From above, we see that R 0 , the basic reproduction number has a positive correlation with Λ, β , and α. This means that an increase or decrease in either of them will bring about a corresponding increase or decrease in the value of R 0 . This is sensible since an increase in the birth rate or effective contact rate will increase the expected number of infected cases. Moreover, the sensitivity index of R 0 with respect to Λ and β is 1; which means that the change in the value of R 0 will be directly proportional to that of Λ and of β . Since people in the Q, T, or R compartment do not spread the disease, we observe that R 0 has no sensitivity towards quarantine individual's recovery rate γ, treatment rate h, and recovery due to treatment ρ. Similarly, λ and d have negative correlation with R 0 . Hence, an increase in any one of them will bring about a decrease in R 0 . This makes sense since, for example, an increase in the rate of quarantine will decrease the number of individuals spreading disease, and hence the expected number of infected cases. Thus, the sensitivity analysis of R 0 suggests that increasing the isolation or quarantine rate λ is a good way to reduce the disease spread.
We have shown the PRCC diagram of our model (2.1) in Figure 3.

NUMERICAL SIMULATIONS
We illustrate our theoretical results via numerical simulations in this section. We will use the Runge Kutta RK4 method to solve the initial value problem. We shall modify the initial values for the populations in various compartments of the model in order to corroborate our results on local and global stability. In each case, the total initial population is calculated by N(0) = S(0) + E(0) + I(0) + Q(0) + T (0) + R(0). To substantiate the dependence of stability on the value of R 0 , we use three different sets of parameters value so that R 0 < 1 , R 0 = 1, and R 0 > 1, respectively. the disease-free equilibrium is DFE = (500, 0, 0, 0, 0, 0). Since R 0 = 1, by Theorem 5, DFE is globally asymptotically stable. In fact, we see from Figure 5 that the final susceptible population approaches the value 500 asymptotically for any initial sub-population, while the other classes viz., exposed, infectious, quarantined, recovered, and under-treatment class converge to zero.

DISCUSSIONS AND CONCLUSIONS
In this paper, we have formulated a continuous-time mathematical model for the COVID-19 disease spread, described by a system of ordinary differential equations. We have studied the basic properties of the model (2.1) viz., positivity and boundedness. The basic reproduction number has been derived to be R 0 = αβ Λ d(d+α)(d+h+γ) , and its relation to the stability of the model dynamics has been established. We calculated the equilibria of the model and found that the disease-free equilibrium is locally stable for R 0 ≤ 1, while the endemic equilibrium is locally stable for R 0 > 1. We also established that there is a unique globally stable disease-free equilibrium when R 0 ≤ 1, and a unique globally stable disease-present equilibrium when R 0 > 1. The global stability of the two equilibria has been proved by considering a suitable Lyapunov functional for the respective cases and invoking the LaSalle invariant principle. In order to see the impact of parameters on the spread dynamics of the disease, we performed sensitivity analysis of the basic reproduction number R 0 . Finally, we demonstrated our theoretical findings with the help of some numerical simulations carried out under different initial conditions. From the overall study, it is found that the transmission of COVID-19 can be controlled by maintaining the value of basic reproduction number less than or equal to one, which can be done by manipulating different model parameters. In this regard, the sensitivity analysis indicates that the transmission of COVID-19 can be lowered by increasing the rate of quarantine.
Thus the quarantining of infected population plays a crucial rule to control the transmission of COVID-19. However, there are some drawbacks in this model, for example: the immunity is assumed for lifelong once recovered, which needs further research. But, it is expected that the qualitative behavior of the model will be unchanged by this. In addition, vaccination has not been considered. We plan to incorporate these in our future work.

CONFLICT OF INTERESTS
The author(s) declare that there is no conflict of interests.