A SEIQRD EPIDEMIC MODEL TO STUDY THE DYNAMICS OF COVID-19 DISEASE

. In this paper, we propose a COVID-19 epidemic model with quarantine class. The model contains 6 sub-populations, namely the susceptible ( S ), exposed ( E ), infected ( I ), quarantined ( Q ), recovered ( R ), and death ( D ) sub-populations. For the proposed model, we show the existence, uniqueness, non-negativity, and boundedness of solution. We obtain two equilibrium points, namely the disease-free equilibrium ( DFE ) point and the endemic equilibrium ( EE ) point. Applying the next generation matrix, we get the basic reproduction number ( R 0 ) . It is found that R 0 is inversely proportional to the quarantine rate as well as to the recovery rate of infected sub-population. The DFE point always exists and if R 0 < 1 then the DFE point is asymptotically stable, both locally and globally. On the other hand, if R 0 > 1 then there exists an EE point, which is globally asymptotically stable. Here, there occurs a forward bifurcation driven by R 0 . The dynamical properties of the proposed model have been veriﬁed our numerical simulations.


INTRODUCTION
Coronavirus Disease-2019 or COVID-19, which was first found in Wuhan, China in December 2019, is an infectious disease caused by a newly discovered coronavirus. COVID-19 virus is primarily transmitted through droplets produced when an infected person coughs, sneezes, or exhales. Fever, coughing, and breathing problems are the earliest symptoms of this disease. In the following stages, the infection can cause pneumonia, serious acute respiratory syndrome, kidney damage, and even death [1]. COVID-19 is currently receiving great attention among researchers, governments and the general public due to its high rate of spread of infection and significant number of death [2]. As a result, many researchers in the Mathematics field are investigating and developing mathematical models of the spread of COVID-19 disease [3][4][5], as well as forecasting the time required for the disease to be eradicated [6][7][8].
When modelling the spread of an epidemic, population is divided into several fundamental disease states. One of simple epidemic models is the susceptible-infected-recovered (SIR) model, which consists of the susceptible class, the infected class, and the recovered class [9].
Individuals in susceptible state are at risk of getting the disease. Infected individuals are contagious and may spread the disease to others through contact with susceptible individuals. Recovered state consists of individuals who have immunity from a vaccine or acquired immunity and died people. Several investigating of COVID-19 SIR model are studied in [10][11][12]. The models are in the form of system of three nonlinear ordinary differential equations (ODEs).
It is well known that COVID-19 has latency or incubation period. Hence, some researchers have considered SEIR COVID-19 epidemic model by including an exposed (E) class, those who have been infected but have not exhibited any disease symptoms and cannot transmit the disease [13][14][15][16]. By taking into account the further characteristics of COVID-19, Peng et al. [18] extended the SEIR model by involving class of protected inviduals, quarantined individuals and closed (or death due to COVID-19 disease). In this case, Peng et al. [18] divides the population into 7 sub-population classes, namely Susceptible (S), Insusceptible, i.e. susceptible but protected because they always apply health protocols in an orderly manner (P), exposed (infected but not yet be infectious) (E), infected population (I), infected and quarantined privately or in hospital (Q), recovered (R) and closed cases (or dead) from COVID-19 (D). The quarantined class can't spread the disease to others because they have lost contact with susceptible people [17]. It is expected that the COVID-19 outbreak will be controlled by the presence of a quarantined sub-population [5]. Peng et al. [18] obtained a SEIQRD epidemic model. The SEIQRD model in [18] does not consider demographic parameters such as birth and death rates. Therefore López and Rodó [19] modified the model (2.6) by including those two demographic parameters. However, they ignored the deaths of infected individuals caused by the COVID-19 disease and assumed that the infected individuals cannot recover from the disease unless they quarantine first. Zeb et al. [2] have also proposed a COVID-19 epidemic model which includes S, E, I, Q and R classes, but without considering the closed class (D). Zeb et al. [2] assumes that disease transmission occurs through contact between susceptible individuals and exposed individuals or infected individuals where the rate of transmission follows the bilinear incidence rate.
According to Postavaru et al. [20], the exposed individuals are infected individuals but not yet infective. Thus, in this article, we propose a SEIQRD COVID-19 epidemic model by assuming that COVID-19 transmission occurs only when there is contact between susceptible individuals and infected individuals, with the transmission rate being the standard incident rate. We also consider deaths due to COVID-19 disease and recovery in both infected and quarantined subpopulations into the model, which are not considered by Zeb et al [2] nor by López and Rodó [19].
The outline of this study is given as follows. In Section 2, we present the SEIQRD model with some assumptions. We provide the basic properties of proposed model in Section 3 and 4. In Section 5 and 6, we study the equilibrium points and basic reproduction number. The stability of equilibrium points, both locally and globally, are provided in Section 7. The numerical simulation was carried out in Section 8. Finally, some conclusion will be given in Section 9.

MODEL DEVELOPMENT
Based on the assumptions described previously, we formulate a SEIQRD COVID-19 model by considering the model of Zeb et al. [2] and the model of López and Rodó [19], but without including the confined (C) compartment. The compartment diagram of the proposed model can be seen in Figure 1.
Based on the diagram, the developed model is a system of (1) where S represents the size of the susceptible sub-population, i.e. all people who are at risk of contracting COVID-19; E represents the size of the exposed sub-population, i.e. people who have been infected and have not shown symptoms of disease, I represents the size of the infected sub-population, i.e. people who have shown symptoms of the disease, Q represents the size of the quarantined sub-population, i.e. people who are infected and self-isolated or hospitalized for treatment, R represents the size of the sub-population who have fully recovered from COVID-19, and D represents the size of death sub-population caused by COVID-19. The parameters of system (1) are described in the Table 1.

THE NON-NEGATIVITY AND BOUNDEDNESS OF SOLUTION
In this section, we prove the non-negativity and boundedness of solution of model (1) to show that the model is epidemiologically meaningful. Proof. We first show the non-negativity of S. By assuming that this not the case, the intermediate value theorem in [21] guarantees the existence of τ > 0 such that This means that S > 0 in (τ, τ + ε) for arbitrary small positive constant ε. This leads to a contradiction. Hence, S ≥ 0 for all t > 0. The non-negativity of E, I, Q, and R can be proven analogously. Therefore, all solutions of model (1) are non-negative.
Total population is generally defined as the number of living humans, so that the total population (N(t)) is obtained by adding up all sub-populations in the model (1) except the sub- If all equations in the system (1) are summed up, then the following differential equation is obtained.
The solution of the equation (2) satisfies Hence, we can see that the feasible region of model (1) is which is positively invariant region.

THE EXISTENCE AND UNIQUENESS OF SOLUTION
The existence and uniqueness of solution of model (1) can be proven by Derrick and Groosman theorem in [22], which states that if Lipchitz's condition as in Definition 1 is satisfied, then the solution of model (1) exists and is unique.
The following theorem guarantee the existence and uniqueness of solution of model (1).
Theorem 2. The model (1) subject to non-negative initial values has a unique solution in Ω for all t ≥ 0.
Proof. The right side of model (1) can be written as follows.

EQUILIBRIUM POINTS
Let X = (S, E, I, Q, R) T . By setting we get E = (σ +θ +δ +µ) Then, by substituting E to the second equation of (3), we obtain Theorem 3. The model (1) has two equilibrium points as follows.

BASIC REPRODUCTION NUMBER
One of important epidemiologic metric is the basic reproduction number (R 0 ), which measures the contagiousness or transmissibility of infectious agents. The basic reproduction number can be determined by the next generation matrix method. For that aim, we consider Z = (E, I, Q) T . Then, we have The Jacobian matrices of F and V evaluated at ε 0 are respectively given by F and V as follows Hence, the next generation matrix is given by The basic reproduction number (R 0 ) is a spectral radius of the next generation matrix, that is, .
We observe that R 0 is inversely proportional to the the recovery rate of the infected subpopulation (θ ) as well as to the quarantine rate σ . The larger the value of θ or σ , the smaller the R 0 value. We also notice that I * in endemic point can be written in term of R 0 as follows Thus, the endemic equilibrium point exists if R 0 > 1.

THE STABILITY OF EQUILIBRIUM POINTS
In this section, we first provide the local asymptotic stability of disease-free equilibrium point ε 0 as follow.
Theorem 4. The disease-free point ε 0 is locally asymptotically stable in domain Ω if R 0 < 1.
In the following theorems, we show the global asymptotic stability of the disease-free equilibrium point ε 0 and the endemic equilibrium point ε * .
Theorem 5. The disease-free equilibrium point ε 0 is globally asymptotically stable in the do- Proof. To prove this theorem, we follow the method of Castilo-Chavez et al. [24]. First, we rewrite the model (1) as where the elements of Y = (S, R) ∈ R 2 + describe the number of non-infected individuals while the components of Z = (E, I, Q) ∈ R 3 + indicate the number of non-infected individuals. Let Based on the theorem in [24], the ε 0 is globally asymptotically stable if R 0 < 1 and the following conditions hold: (H.1) Y 0 is globally asymptotically stable for system dY dt = F 1 (Y, 0). (H.2) F 2 (Y, 0) = 0 and F 2 (Y, Z) = CZ −F 2 (Y, Z) whereF 2 (Y, Z) ≥ 0 for any (Y, Z) ∈ Ω and C is the Jacobian matrix ∂ F 2 ∂ Z evaluated at ε 0 .
We notice that It is clearly seen that the elements ofF 2 (Y, Z) are non-negative, and therefore condition (H.2) is satisfied.
We next consider that from which we get It is observed that S(t) → Λ µ and R(t) → 0 as t → ∞, showing that Y 0 is globally asymptotically stable. Since condition (H.1) is also satisfied, the disease-free point ε 0 is globally asymptotically stable in domain Ω.
Theorem 6. Assume that R 0 > 1, and thus the endemic equilibrium ε * exists. Then, the point ε * is globally asymptotically stable in the domain Ω.
Proof. Consider a Lyapunov function W , which is defined by The Lyapunov function W is a positive definite function in region Ω. The derivative of W with respect to t is Since the geometric mean is less than or equal to the arithmetic mean [25], we have It is clear that dW dt ≤ 0. Furthermore, dW dt = 0 is satisfied if and only if S = S * , E = E * , and I = I * . The LaSalle's Invariance Principle in [25] guarantees that S → S * , E → E * , and I → I * as t → ∞. Hence, Q → Q * and R → R * as t → ∞. In other words, ε * is globally asymptotically stable.
From discussion above, we know that the disease-free equilibrium ε 0 always exists. If R 0 < 1 then the endemic equilibrium point ε * does not exist and ε 0 is globally asymtotically stable.
However, if R 0 > 1 then ε * co-exists and is globally asymtotically stable. Hence, based on Martcheva [25], there occurs a forward bifurcation driven by R 0 . Because R 0 is directly proportional to the infection rate β , the forward bifurcation may also caused by β . For example, we perform simulation using parameters value as in Table 2 and varying the value of β . Based on this simulation, we plot the bifurcation diagram as shown in Figure 2. It can be seen that the model (1) exhibits a forward bifurcation, where the bifurcation occurs when R 0 = 1 or β = 0.2314.

NUMERICAL SIMULATIONS
To confirm our previous analytical results, we solve the SEIQRD epidemic model (1) using the fourth-order Runge-Kutta method. If not stated otherwise, our simulation uses parameter values as in Table 2   is globally asymptotically stable. This stability properties is confirmed by our simulation depicted in Figure 3, that is the numerical solution is convergent to the endemic equilibrium point  Table 2 and β = 1.7.
To see the effect of the recovery rate of the infected class (θ ), we perform simulation using parameters value as for the first simulation but with θ = 1.2. The basic reproduction number in this case is R 0 = 0.8233 < 1 and hence we can expect that the disease-free equilibrium point ε 0 = (35200, 0, 0, 0, 0) is globally asymptotically stable. This is confirmed by our numerical result shown in Figure 5.   Table 2.  Table 2.

CONCLUSIONS
A SEIQRD model to describe the spread of COVID-19 disease using the standard incidence rate has been developed in this paper. The model consists of susceptible (S), exposed (E), infected (I), quarantined (Q), recovered (R), and death caused by the COVID-19 disease(D) subpopulations. The existence, uniqueness, non-negativity, and boundedness of solution have been proven, showing that the proposed model is biologically feasible. The model has two equilibrium points, namely the disease-free equilibrium point and the endemic equilibrium point.
Using the next generation matrix method, we have determined the basic reproduction number.
The disease-free equilibrium point always exists and it is locally and globally asymptotically stable if the basic reproduction number is less than unity. If the endemic equilibrium point exists, i.e. when the basic reproduction number is greater than unity, then it is always globally asymptotically stable. The proposed model exhibits a forward bifurcation, where the basic reproduction number acts as a bifurcation parameter. The results of our dynamical analysis are confirmed by our numerical simulations. Furthermore, from the basic reproduction number formula and our numerical simulation results, the basic reproduction number can be reduced by increasing the rate of recovery or quarantine of the infected sub-population. This shows that COVID-19 disease can be controlled by treating infected individuals or by quarantining them.