key: cord-0783099-3hfxv6h8
authors: Nave, OPhir; Hartuv, Israel; Shemesh, Uziel
title: Θ-SEIHRD mathematical model of Covid19-stability analysis using fast-slow decomposition
date: 2020-09-21
journal: PeerJ
DOI: 10.7717/peerj.10019
sha: 0fc56685fa8429341608111613127475d05d050e
doc_id: 783099
cord_uid: 3hfxv6h8

In general, a mathematical model that contains many linear/nonlinear differential equations, describing a phenomenon, does not have an explicit hierarchy of system variables. That is, the identification of the fast variables and the slow variables of the system is not explicitly clear. The decomposition of a system into fast and slow subsystems is usually based on intuitive ideas and knowledge of the mathematical model being investigated. In this study, we apply the singular perturbed vector field (SPVF) method to the COVID-19 mathematical model of to expose the hierarchy of the model. This decomposition enables us to rewrite the model in new coordinates in the form of fast and slow subsystems and, hence, to investigate only the fast subsystem with different asymptotic methods. In addition, this decomposition enables us to investigate the stability analysis of the model, which is important in case of COVID-19. We found the stable equilibrium points of the mathematical model and compared the results of the model with those reported by the Chinese authorities and found a fit of approximately 96 percent.

The coronavirus belongs to the severe acute respiratory syndrome (SARS) family. It usually does not cause disease in humans, but infects animals-mammals and poultry. If humans are infected with the virus, the disease usually causes mild cooling, and the infection passes without any treatment. Conversely, if the infected person has a weak immune system for some reason, the coronavirus can be fatal. It is a common virus that attacks almost every person at least once in his/her lifetime-especially during childhoodand, as is usually the case, it is not dangerous. However, the current epidemic is a recently mutated virus that has become fatal (which is why the World Health Organization initially called it "the new coronavirus" 2019). It should be noted that, in 2002-2003 , there was an outbreak of a similar virus from the SARS family, causing the death of more than 800 people (General Health Fund in Israel, https://www.clalit.co.il/he/yourhealth/family/ The rest of the paper is organized as follows. In "The θ-SEIHRD Model of COVID-19", we describe the mathematical model of the coronavirus called θ-SEIHRD model. In "Results and Analysis", we present the results of our comparative analysis of the cases reported by the Chinese authorities. In addition, we apply the SPVF model to the θ-SEIHRD model and find the stability of the equilibrium points. "Discussion" presents our discussion. Finally, "Conclusions" presents the conclusion of this paper.

In this section, we introduce the mathematical model of the θ-SEIHRD model as presented in (Ivorra et al., 2020) . This model is based on the Be−CoDiS model presented in (Ivorra, Ramos & Ngom, 2015) . It is important at this point to note that this model is not the SIR-SEIR standard model. It considers the known special characteristics of the considered disease, such as the existence of infectious undetected cases and the sanitary conditions and infectiousness of hospitalized patients. The assumptions of the model can be found in (Ivorra et al., 2020) . Here, we present the mathematical formulation of the COVID-19 spread. The system of equations includes nine nonlinear ordinary differential equations and has the following forms:

whereṼ ¼ ðS; E; I; I u ; H R ; H D ; R d ; R u ; DÞ; g E ; g I ; g H R and g I u are the transition rates of E, I, H R , and I u , respectively  1 day à ; β(·) is the disease contact rate  1 day à ; m is the fatality rate  1 day à ; and τ indicates the infected people, who move from one territory to another per day. The initial conditions of the system at time t 0 (for each country, however, the initial time can be changed. For example, in Wuhan, China, t 0 = 10 A. CET (technical adjustment) because it was on 7 December 2019 when the first case of patients with symptoms was confirmed) are as follows:

As can be seen above, the presentation of the above mathematical model has no explicit hierarchy and, hence, it is impossible to apply different asymptotic methods. In the next section, we apply the SPVF methods to expose the hierarchy of the considered model.

In this section, we apply the SPVF method to the system of Eqs. (1)-(9), as presented in (Nave, 2017) . While we apply the SPVF method, we expose the hierarchy of the model. This procedure enables us to rewrite the model under consideration in new coordinates and, hence, the "new" model can be decomposed into the so-called fast and slow subsystems.

The set of parameters and functions used in our calculations is as follows:

Parameters: 

Functions:

b I u ¼ 0:375 þ 0:375 ð1 À vðtÞÞ ð1 À hÞ; vðtÞ ¼ 0:003 Á m I þ 0:997 Á ð1 À m I Þ;

Here, we present only the absolute values of the eigenvalues of the SPVF method (and not the eigenvectors) to avoid overwhelming the readers with too much information: 

As we observed from the above results, and according to the SPVF method, the maximal gap of the eigenvalues is between λ 1 and λ 2 . This gap implies that the first dynamic variable of the system, which is written in the new coordinates, is fast compared with the rest of the variables. The results of the SPVF algorithm are presents in Figs. 1 and 2. We compute the value of |λ i + 1 |/|λ i | (i = 1,…,9), for 365 days and plot only the maximum values of this quotient for every day. According to Fig. 1 we can see that for every day the gap indeed exists (for every day the graph is not zero). This means that the SPVF algorithm exposes the hierarchy of the system at the new coordinates. But it is not clear from these results what is the (1)-(6) and the slow subsystem are Eqs. (7)- (9). The results are present in Fig. 2 . We can see that at some days the fast subsystem includes only the first equation and for some other days the fast subsystem includes the first two equations of the model in the new coordinates. One can analyze the results of the SPVF algorithm from the other direction. That is, if we look at Fig. 2 on the 100th day, for example, then we see that only the first equation is the fast subsystem. But if we look at the 250th day for example, then the first two equations are the fast subsystem. The SPVF algorithm aims to find a frame of coordinates for which the system has an SPS form. Hence, we have a degree of freedom to choose any frame of coordinates for this purpose. At our analysis, we have chosen the frame of coordinates, that is, eigenvectors, that belong to the 49th day where on this day we received the maximal gap of all the gaps that we obtained from the SPVF algorithm. On this day the fast subsystem is the first equation of the model.

We The new model is rewritten in the new coordinates using the eigenvectorsũ 1 ;ũ 2 ; …;ũ 9 , which correspond to the eigenvalues λ 1 , λ 2 ,…,λ 9 . This means that the new dynamic variables of the model are linear combinations of its old variables when the coefficients of the linear combinations are taken from the eigenvectors, that is, 

where t refers to the transpose operation. In matrix form, system (15) can be written as

Here, we denote the matrix B as a matrix whose columns are the eigenvectors of the SPVF method. The right-hand side (RHD) of this system is a function of the variablesṼ, whereas its left-hand side is a function of the new variablesṼ. To write the system of Eq. (16), with the same variablesṼ, we should differentiate this system with respect to time, that is,

is the vector field of systems (1)-(9). To rewrite both sides of the model using the new variables, we should express the old variables as a function of the new variables using the inverse matrix of the eigenvectors, that is,

Substituting Eqs. (19) into (17) 

is the vector field of the new model (the RHD of the model is written in the new coordinates).

After we transformed and presented the model in the new coordinates using the eigenvectors of the SPVF method, the model can be decomposed into the fast and slow subsystems based on the gap of the eigenvalues. The first eigenvalue, λ 1 , indicate that the first variable of the system, sayS, is the fast variable of the system and thatẼ,Ĩ,Ĩ u ,H R ,H D , R d ,R u andD are the slow ones. Hence, the model presented in the system of Eq. (20) can be decomposed as follows:

According to the SPVF method, the stability analysis can be investigated on the fast subsystem. To find the equilibrium points of the model, we should solve the following linear system:

This implies that we look for the equilibrium points of the fast variable (denoted by an asterisk), that is,S Ã , while the rest of the variables are frozen, and we take their values as the initial conditions of the system. System (22) is a system of one equation with one unknown (fast) variables:S Ã . We solve this system and substitute the valuesS Ã in the slow subsystem and solve the system for the slow variables, that is,Ẽ Ã ,Ĩ Ã ,Ĩ Ã u ,H Ã R ,H Ã D ,R Ã d ,R Ã u and D Ã . Now, we have eight equations with eight unknown (slow) variables.

To check the stability of these equilibrium points, we compute the Jacobian matrix calculated at the equilibrium point. It is important at this point to note that the equilibrium points in the new coordinates have no biological meaning. Moreover, to provide a biological interpretation to the equilibrium points that we received, we have to transfer them to the old coordinates for which we need to calculate the inverse matrix of the matrix with the eigenvectors, that is,

Here, we inverse transform only the stable equilibrium points. We received only one stable point in the new coordinates, and, therefore, we should transform only one stable equilibrium point.

According to our numerical simulations, if the initial condition is on 12 January 2020, then the equilibrium point is reached after 49 days. According to official reports (from the Wuhan County authorities), the closure of Wuhan ended after 62 days. This means that the situation stabilized after 62 days. As we can see from the results of the mathematical model and what actually happened, there is a difference of 13 days. Medically, another 13 days of lockdown would not have affected the civilians. However, apparently, the economic considerations outweighed the medical considerations and the authorities decided to lift the lockdown earlier than expected. In Fig. 3 , we present the evolution of the number of cases in China for different values of the parameters. As we can see from these results, the equilibrium point is attained at approximately 83K. The results are summarized in Table 1 . In this table, we present only the important variables that are stable. 

As we have shown in the previous section, we obtain the stable equilibrium points of the mathematical model owing to the application of the SPVF method. The application of this method enables us to decompose the system into a fast and a slow subsystems. After the decomposition, we explored only the fast subsystem. In general, instead of decomposing a given system into fast and slow subsystems, one can solve the original system and determine the stable equilibrium points numerically. However, the big problem with the numerical method is that the equilibrium points are represented by their values and not as points, and they depend on the original system's parameters as can be obtained by the considered decomposition. That is, if we want to change the system parameters and determine new equilibrium points, we have to resolve the mathematical model numerically each time, which is time consuming. Moreover, the SPVF method allows us to first find all the equilibrium points of the original model analytically. In addition, the equilibrium points depend on the parameters of the original system. Therefore, if we want to change the model parameters, we do not have to resolve the mathematical model again but only need to change the parameters in the stable equilibrium points that have been determined. As can be seen from the results obtained from the mathematical model and from the results reported by the authorities, the relative error of the model can be calculated for which, we define the relative error of the equilibrium points for each dynamic variable of the system as

where U model is the dynamic variable of the system, that is, the stable equilibrium point of the original model (after the inverse transform of the stable equilibrium points) and U real denotes the values reported by the authorities. The results are as follows:

As can be seen from the results of the relative error, the reported cases are indeed small in number, as was expected from the decomposition method.

In this paper, we investigated the stability of the θ-SEIHRD mathematical model of COVID-19. The mathematical model of coronavirus is presented with a hidden hierarchy. This implies we cannot know which of the dynamic variables of the model is progressing fast and which one is progressing slowly. The hidden hierarchy of the model neither allows for an asymptotic analysis in general, nor an analytical investigation in particular, but only the running of numerical simulations. In particular, the equilibrium points of the system and their stability cannot be investigated and obtained from the present model. Therefore, in this study, we implemented the SPVF method, which explicitly exposes the system hierarchy. This method transforms the model into new coordinates using the eigenvectors of the given vector field. In the new coordinates, the model is presented in the form of the SPS, that is, with an explicit hierarchy. The hierarchy of the new mathematical model enables us to decompose the model and divide it into subsystems. In our case, we only divided it into two subsystems: a fast and a slow subsystem. According to the SPVF theory, the fast subsystem can be investigated while the slow system is frozen. We found the equilibrium points the fast system in the new coordinates analytically. We analyzed the stability of the equilibrium points and obtained a single equilibrium point for the new model. To obtain a biological interpretation of the equilibrium point, we inverse transformed the stable equilibrium point into the old coordinates. According to our analytical results and the official reports of the Chinese authorities, we found that the stable equilibrium points obtained from the mathematical model are very close to those of the official reports.

Susceptible: The person is not infected by the disease pathogen E Exposed: The person is in the incubation period after being infected by the disease pathogen and has no visible clinical signs. The individual could infect other people but with a lower probability than that by the people in the infectious compartments. After the incubation period, the person moves to the infectious compartment I I Infectious: After the incubation period, the first compartment of the infectious period in which nobody is expected to be detected yet. The person has completed the incubation period, may infect other people, and start developing clinical signs. After this period, people in this compartment can be either taken in charge by sanitary authorities (and we classify such people as hospitalized) or may not be detected by the said authorities and continue to be infectious (but in another compartment, I u ) I u Infectious but undetected: After being in compartment I, the person can still infect other people and have clinical signs; however, they have not yet been detected and reported by the authorities. We assume that only people with low or medium symptoms can be included in this compartment, and not the people who die. After this period, people in this compartment move to the recovered compartment R u .

Hospitalized: The person is in the hospital (or under quarantine at home) and can still infect other people. At the end of this stage, a person moves on to the Recovered compartment R d H D Hospitalized who will die: The person is hospitalized and can still infect other people. At the end of this state, the person is transferred to the Dead compartment R d

Recovered after being previously detected as infectious: The person was previously detected as infectious, survived the disease, is no longer infectious, and has developed a natural immunity to the virus. When a person enters this

Singularly perturbed vector fields

Mathematical modeling of the spread of the coronavirus disease 2019 (COVID-19) taking into account the undetected infections. The case of China

Mathematical formulation and validation of the Be-FAST model for classical swine fever virus spread between and within farms

Validation of the forecasts for the international spread of the coronavirus disease 2019 (COVID-19) done with the Be-CoDis mathematical model

Application of the Be-CoDis mathematical model to forecast the international spread of the 2019 Wuhan coronavirus outbreak

Be-CoDiS: a mathematical model to predict the risk of human diseases spread between countries-validation and application to the 2014 Ebola virus disease epidemic

2020. Coronavirus COVID-19 Global Cases by the Center for Systems Science and Engineering (CSSE)

Early dynamics of transmission and control of COVID-19: a mathematical modelling study

Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (SARS-CoV2)

Transmission dynamics of 2019 novel coronavirus

Singularly perturbed vector field method (SPVF) applied to combustion of monodisperse fuel spray

The global dynamics for an age-structured tuberculosis transmission model with the exponential progression rate

compartment, he/she remains in the hospital for a convalescence period of C 0 days (average) R u Recovered after being previously infectious but undetected: The person was not previously detected as infectious, survived the disease, is no longer infectious, and has developed a natural immunity to the virus D Dead by COVID-19: The person did not survive the disease t duration

The authors received no funding for this work.

The authors declare that they have no competing interests.

OPhir Nave conceived and designed the experiments, performed the experiments, analyzed the data, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft. Israel Hartuv conceived and designed the experiments, performed the experiments, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft. Uziel Shemesh conceived and designed the experiments, performed the experiments, analyzed the data, prepared figures and/or tables, authored or reviewed drafts of the paper, and approved the final draft.

The following information was supplied regarding data availability: Data is available at WorldMeters.info: https://www.worldometers.info/coronavirus/ country/china/