A reliable and competitive mathematical analysis of Ebola epidemic model

The purpose of this article is to discuss the dynamics of the spread of Ebola virus disease (EVD), a kind of fever commonly known as Ebola hemorrhagic fever. It is rare but severe and is considered to be extremely dangerous. Ebola virus transmits to people through domestic and wild animals, called transmitting agents, and then spreads into the human population through close and direct contact among individuals. To study the dynamics and to illustrate the stability pattern of Ebola virus in human population, we have developed an SEIR type model consisting of coupled nonlinear differential equations. These equations provide a good tool to discuss the mode of impact of Ebola virus on the human population through domestic and wild animals. We first formulate the proposed model and obtain the value of threshold parameterR0 for the model. We then determine both the disease-free equilibrium (DFE) and endemic equilibrium (EE) and discuss the stability of the model. We show that both the equilibrium states are locally asymptotically stable. Employing Lyapunov functions theory, global stabilities at both the levels are carried out. We use the Runge–Kutta method of order 4 (RK4) and a non-standard finite difference (NSFD) scheme for the susceptible–exposed–infected–recovered (SEIR) model. In contrast to RK4, which fails for large time step size, it is found that the NSFD scheme preserves the dynamics of the proposed model for any step size used. Numerical results along with the comparison, using different values of step size h, are provided.

1976, as the virus was detected in a village, in the vicinity of the Ebola River of the Democratic Republic of Congo (DRC). The first outbreak was one of the most deadly outbreaks in history having case fatality rate of 88% with 318 exposed cases and 280 deaths. This year, the second outbreak occurred approximately 850 kilometers away in South Sudan, having case fatality rate of 53% with 284 exposed cases and 151 confirmed deaths. It is remarkable that almost 25 more Ebola outbreaks have occurred after 1976 in a multitude of countries around the world including the Democratic Republic of Congo, where it originated, Liberia, South Africa, Sierra Leone, Gabon, Uganda, Sudan, Guinea, Spain, and the United States of America and decimated many people [1][2][3][4][5][6].
It is detected that a significant source of virus transmission is a contact with infected animals like fruit bats, porcupines, and non-human primates such as apes and monkeys, considered to be a natural lake for the Ebola virus. Ebola virus initially transmits to people from domestic and wild animals, potentially affecting an oversized range of individuals, and then spreads into the population through transmission among individuals via close physical contact with an infected human or with their bodily fluids like blood, saliva, mucus, tears (eye fluid), bile, secretions, breast milk, spinal column fluid, urine, and semen. It can also be transmitted to others while handling contaminated materials used by the patient like bedding and cloth. Moreover, some women may get infected in the process of breastfeeding, and therefore the virus might stay in breast milk. Aid employees and healthcare staff can also get infected during the treatment of infected patients in health centers [2,4]. Therefore, it is inevitable to take a great precaution.
Usually, symptoms of EVD include fever, muscle pain, loss of appetite, stomach pain, severe headache, sore throat, fatigue, etc., followed by vomiting, rash, bloody diarrhoea, impaired liver and kidney function, both internal and external uncontrollable bleeding, low white blood cell, increase in platelet counts, and elevated liver enzymes. Normally after an average of three days, the patient faces a rapid progression to death. Symptoms may appear anywhere, and the time interval from infection with the Ebola virus to the appearance of symptoms is called the incubation period, which normally is from 2 to 21 days [1][2][3][4]. An infected human may not lead to the spread in the host population until he shows the said symptoms of the viral disease (EVD).
The most severe, deadliest, and the largest outbreak of EVD was reported in March 2014 in Guinea which then moved to Liberia, Mali, Senegal, the United States of America, Nigeria, Spain, and across land borders. On December 2014, World Health Organization has reported overall 17,942 probable and confirmed cases in Africa which included 6388 deaths with a case fatality rate of 36%. A large number of human deaths were caused by Ebola virus in 2014-2016, where most of the infected people were foreign travelers as they traveled to the affected regions, got exposed to the virus, and showed symptoms of Ebola virus fever after they returned back to their homeland. Only in West Africa, after the first case was discovered in 2014, the EVD outbreak ended by 2016 with 113,10 confirmed deaths and nearly 286,16 suspected deaths with fatality rate of 39% (2014-2016 Ebola Outbreak in West Africa reported by CDC) officially recorded in June 2016 as discussed in [7][8][9][10][11][12] and [2][3][4].
Another outbreak of Ebola virus disease in North Kivu Province began on 1 August 2018 reported by the Ministry of Health of the DRC. The outbreak continues with moderate intensity in the eastern DRC. Almost twenty-four health zones of Ituri and North Kivu provinces of the DRC have confirmed many probable cases. On 11 June 2019, the Ugandan Ministry of Health confirmed their first and afterwards two additional imported cases from the DRC into Uganda, and the number of infectious kept quickly increasing day by day [13][14][15]. Behind the West Africa outbreak of 2014-2016, overall world's second largest outbreak of Ebola virus ever recorded was that of 2018-2019. It is the 10th and the largest Ebola outbreak ever recorded in the DRC. Guinea having about 2500 deaths due to EVD by May 2018 is one among the three hardest hit countries of West Africa. There were 2763 cases with 1841 confirmed deaths recorded from North Ituri and Kivu provinces by the DRC Ministry of Health on 4 August 2019 [16].
Clinically, it may be more troublesome to handle Ebola epidemic malady as compared to other infectious diseases in the world. Recently, a study carried out in 2016 determined the immunogen of VSV-EBOV that is 70-100% effective to shield the patient against Ebola infection. It is thought to be the first primary vaccine for the patient to fight against Ebola virus sickness. However, the United States FDA has not approved any vaccine to be used in humans. During the 2018-2019 Ebola eruption in the DRC, the first-ever multi-drug irregular management trial was conducted under an ethical framework developed in consultation with specialists within the field and by the DRC, simply to judge the safety and effectiveness of drugs employed in the treatment of Ebola patients. Today, CDC is assisting the Uganda government and the DRC, the eruption infected bordering countries, native and international partners like WHO to coordinate the activities and supply technical steering associated with surveillance, contact tracing, infection control management, risk communication, laboratory testing, vaccination, data management, health screening at the border, and health education. A key to successfully controlling outbreaks of EVD is community engagement, control practices, and protective measures. To scale back the possibilities of transmission, we should introduce a vaccination or other antiviral drugs that people will take, and also raise awareness to reduce the risk factors for Ebola infection due to human transmission. For example, we should try to avoid contact with infected human beings and infected wild animals such as infected fruit bats, porcupines, monkeys, and apes. Infected people should be treated in a very protective way especially when handling body fluids. It is very important to wash our hands after caring after patients at home or on visiting the hospital. We should treat animals only with gloves or other proper clothing on and carefully treat their blood and meat. In the affected countries, infection of Ebola virus can occur through touching the dead bodies of infected humans. To control the propagation of EVD, dead bodies of patients should be handled and kept for a minimum time and burial should be done by a group of trained people, those who have better information to conduct safe and dignified burial.
In this study, we concentrate on Ebola virus to assess the dynamics of EVD for the direct and indirect environmental transmission. For this purpose, we present an SEIR mathematical model in which the primary emanator of Ebola infection is the infected wild and domestic animals, from which the virus is transmitted slowly to household and subsequently reaches the rehabilitation places such as health care centers and hospitals, thereby taking on the medical staff. Henceforth, if not capped, the virus may transform into a wide scale epidemic, threatening the human beings at large. The important fact is that proper precautions must be used if one needs to treat infected patients or animals. The study, moreover, features a general assessment of, and different considerations with regards to, the multi-pronged and multi-faceted rampage mechanisms of Ebola virus disease pertaining to the environment and the linkage of the same to the disease broadly. Some compre-hensive mathematical techniques are used to analyze the proposed SEIR epidemic model mathematically in a reliable way.
It is always necessary to discretize the continuous model for practical purposes. The obtained discrete SEIR model should possibly contain all important dynamical properties of the corresponding continuous model. For example, most of the standard finite difference schemes, such as Runge-Kutta and Euler method and many other standard methods, when implemented to a dynamical system, can lead to major issues such as negative solutions, converging to wrong equilibrium point or wrong periodic cycle, and sometimes give numerical instabilities [2,[17][18][19] for the proposed model by increasing the time step size. In this situation, these numerical schemes become unbounded and divergent. We propose a numerical scheme to solve the proposed model by implementing a non-standard finite difference (NSFD) scheme [20]. The presented NSFD scheme is used to perform a reliable mathematical analysis of the SEIR model [2] which is available in the existing literature. Therefore, we perform a numerical analysis of the continuous SEIR model along with some comparisons with RK4 scheme for different parameter values involved to investigate the dynamic consistency of the developed NSFD scheme.
The study has been divided into various segments. Section 2 lays down a brief description of the mathematical SEIR model of Ebola virus which consists of coupled nonlinear differential equations. Equilibrium points and the reproductive number R 0 are determined in Sects. 3 and 4, respectively. Section 5 portrays the elements of scope and prospects such as positivity and boundedness of the prescribed solutions of the model. Complementarily, we shed light on local and global behaviors of SEIR model at equilibrium points in Sect. 6. Section 7 is devoted to numerical analysis and discussions. We conclude the research article in Sect. 8.

Model description and formulation
Recently, Tahir et al. [2] considered a compartmental SEIR mathematical model of Ebola virus with a closed population to describe the epidemiology and natural history of Ebola. On a careful reading of [2], we have observed some flaws and gaps in the mathematical analysis dealing with SEIR. In this paper, we not only address the problems in [2], but also include other aspects of qualitative analysis of the model. We summarize our findings as follows: 1. It is worth mentioning that in the formulation of the model in [2] the term β 3 IR appearing in the model has no explicit interpretation. Indeed, infected individuals can move to the recovered class in the presence of different measures which are missing here. Hence we have modified this term by replacing β 3 IR with β 3 I. 2. They have not discussed the positivity of variables of the epidemic problem, which is an essential property of population dynamics. 3. Disease-free and endemic equilibrium points do not satisfy the steady state epidemic model qualitatively and quantitatively. Moreover, it is observed that if the parameters used in [2] are considered, then both I and R in epidemic equilibria EE become negative, which is contrary to the positivity of solutions. 4. As there are some calculation problems in the evaluation of equilibrium points mentioned before, eigenvalues of the Jacobian at DFE and EE are misleading in connection with the stability of the model.

5.
Global stability of the model at DFE is wrongly settled because without the involvement of reproductive number R 0 , we reduce the beauty of the semi-negative property of the Lyapunov function. 6. The analytic study of equilibrium points is not compatible with graphic representation. 7. We use the NSFD scheme instead of RK4. It is worth mentioning that RK4 is not always convergent and hence not appropriate to analyze the quantitative behavior of an SEIR epidemic model [2]. However, our scheme is more reliable and hence the title of our paper.
We will rectify all of the above mistakes and give a true analysis of the epidemic model [2]. The infection of Ebola virus might be transmitted from both domestic and wild animals. The modified SEIR model divides the total population into four epidemiological classes. Susceptible humans at time instant t, that is, humans which are not yet infected but can get infected by Ebola virus, are placed in class I (say). The number of elements in class I is denoted by S. Susceptible individuals that may become exposed after an effective contact with any of Ebola infected human are placed in class II (say) which consists of exposed humans at time instant t, that is, those who still demonstrate no side effects of Ebola infection. The number of elements in class II is denoted by E. The infected humans at time instant t are included in epidemic class III. The number of elements in class III is denoted by I. Ebola infected humans that finally recovered, acquired long immunity in life, and may solely die naturally are placed in the fourth class of Ebola mathematical model. The number of elements in this class is represented by R. Thus S, E, I, and R are the variables for epidemic model (1)- (5). Compartment to compartment transmission flow of Ebola virus in SEIR model (1)-(5) is shown in Fig. 1.

Equilibrium points
SEIR model (1)-(5) admits two equilibrium points in the feasible region (7). A disease-free equilibrium (DFE) of the proposed model (1)-(5) will occur if R 0 < 1. At this stage, there is no infection in the entire population, that is, I = 0. This implies that all the Ebola infected classes will diminish, and finally human population incorporates Ebola free/susceptible humans only. However, system (1)-(5) has a unique endemic equilibrium (EE) β 2 (π + μ) , If R 0 > 1, then we have a stage when the disease will spread in the host population.
4 The basic reproduction number R 0 A threshold parameter R 0 is very important as it is used to assess the prospects or dynamics of any disease. Epidemic will exist when an infected appears into a completely susceptible host population. The basic reproductive number R 0 therefore gives the average measure of new infections produced by a primary infection. The fate of EVD and the dynamical behavior of model (1)-(5) are determined and controlled by R 0 . There will be no epidemic in the human population if R 0 < 1, and it occurs when R 0 > 1. Thus, we require a pack of control strategies if the disease becomes epidemic. Several research articles [21][22][23][24][25][26] are devoted to calculating the basic reproductive number R 0 for different epidemic models. In this section, we determine R 0 for the proposed model and calculate it by the next-generation matrix approach [2,21,22]. For the proposed model (1)-(5), the next-generation approach is performed as follows: The Jacobian of F and V is given bȳ The dominant eigenvalue of the product matrixFV -1 denoted by is the required value of threshold parameter for the proposed model (1)-(5).

Positivity and boundedness of solutions
As we are dealing with human populations, all associated parameters used in model (1)

Positivity of solutions
The solutions S, E, I, R of system (1)-(5) are positive for all t ≥ 0 with nonnegative initial conditions, when they exist. Proof Let Then clearlyt > 0. Suppose that S(0) ≥ 0. Equation (1) leads to Put f (t) = (β 1 + β 4 + β 6 )E -(β 5 + β 7 )I. Multiplying both sides by exp(μt + t 0 f (t) dt) > 0, the above equation becomes Integrating both sides from t = 0 to t =t, we obtain We multiply both sides by exp(-μtt 0 f (t) dt) > 0 to get As S(0) ≥ 0, the sum of the positive terms S is positive. Similarly, we can prove that the quantities S, E, I are positive for all t > 0. Moreover, any solution (S(t), E(t), I(t), R(t)) of model (1)-(5) satisfies the implication which completes the proof. Proof The total population is represented by Z and is defined as

Boundedness of the solutions
Differentiating the above equation with respect to t, we obtain that Using Eqs. (1)-(4), we get that Suppose, for any initial condition, It follows from Eq. (6) that By Gronwall's inequality, we have and hence Z(t) ≤ λ μ for all t ≥ 0 whenever Z(0) ≤ λ μ . Clearly, This shows that Z(t) and all other variables S, E, I, R of model (1)-(5) are bounded. Therefore, SEIR epidemic model (1)-(5) will be analyzed in a biologically feasible region The differential equation (6) for Z shows that the solution of Eqs. (1)-(4) exists in the positive orthan R 4 + , eventually enters and remains in the attracting subset B (since the set B attracts all solutions in R 4 + ). Thus the set B contains a local as well as global attractor of dynamical system (1)- (5). Moreover, the set B ⊂ R 4 + is compact and positively invariant with respect to model (1)-(5) with nonnegative initial conditions in R 4 + .

Behavior of the model at disease-free equilibrium in local sense
The local stability analysis of system (1)-(5) at point F 0 = (S 0 , E 0 , I 0 , R 0 ) = (λ/μ, 0, 0, 0) is discussed in this subsection. The Jacobian matrix for system (1)-(5) at F 0 is obtained as follows: We now prove the following important result for local stability analysis of the epidemic model at disease-free equilibrium.
For local stability analysis of the given model at F 1 , we will prove the following well-known result.

Theorem 4
If R 0 > 1, the proposed system (1)-(5) on set B is locally asymptotically stable (LAS) at endemic equilibrium F 1 ,whereas the system will be unstable if R 0 < 1.
Since In this case, each solution trajectory which starts in the feasible region B with some initial condition approaches F 0 as t → +∞. As a result, Ebola virus disease (EVD) eventually disappears from the host population. Thus by La Salle's invariance principle [57], we conclude that F 0 is GAS on B.
Theorem 6 For R 0 > 1, the endemic equilibrium point F 1 of system (1)-(5) is GAS on B if S = S 1 , E = E 1 , I = I 1 , and for R 0 < 1, the system is unstable.
Proof The following Lyapunov function V : B → R [34] is considered as a candidate to show the global stability of system (1)-(5) at endemic equilibrium F 1 , defined by the relation where K 1 , K 2 , and K 3 are positive constants to be chosen latter.
This means that each solution trajectory, which starts in the feasible region B, approaches F 1 as t → +∞ implies that Ebola virus disease (EVD) spreads in the host population. Thus by La Salle's invariance principle [57], it is concluded that F 1 is GAS on B.

Numerical analysis
This section is devoted to the numerical interpretation of SEIR model (1)-(5) using RK4 and NSFD method coded with Matlab. Different parameters and their numerical values have been taken from [2,4] as given in Table 1. First, we develop both the numerical schemes for the epidemic model, then numerical simulations are provided by the graphs to observe the dynamical behavior of EVD over time t. We also discuss the numerical results.
Proof The proof of this theorem follows similar lines as in [34].

Numerical results
This section includes numerical interpretation of system (1)-(5) using RK4 (16)  shows that the RK4 scheme fails to converge, whereas the NSFD scheme improves the result obtained by RK4 and it is seen to be convergent dynamically to the correct endemic equilibrium F 1 As the last trial, if we compare both schemes for the step size h = 2.0 and h = 2.5, RK4 exhibits negative solutions and does not converge to both F 0 and F 1 , respectively, gives negative solutions again, but NSFD preserves positivity and gives convergence of the solutions, as shown in Figs. 6-9.
The above discussion shows that the RK4 scheme is not always convergent rather conditionally convergent and depends upon the value of step size h, and fails for large step size, whereas Figs. 2-9 illustrate the power of an unconditionally convergent NSFD scheme to produce the converged and positive solutions of model (1)-(5) for any value of the step size h. Moreover, the proposed NSFD scheme is numerically stable (by Theorems 7 and 8) and easy to implement. Using different values of parameter h, numerical experiments are performed, and then the obtained results are compared for both RK4 and NSFD schemes. For both numerical schemes, the effect of different time step h is shown in Table 2.

Conclusions
In this paper, we have considered an SEIR epidemic model of Ebola virus affected by wild and domestic animals, which spread the infection within the human population at any time t. We have studied the dynamical behavior of the proposed model and the dynamics is determined by the basic reproduction number R 0 that acts virtually in controlling the infection of Ebola virus. We have proved the boundedness and nonnegativity of solutions and then well-posedness of the model. Both the disease-free and endemic equilibrium points for the epidemic model were presented and further analyzed for stability. It was proved that the disease-free equilibrium is locally and globally asymptotically stable for system (1)-(5) when R 0 < 1, which shows that EVD will die out at time instant t. Moreover, the endemic equilibrium is stable locally and globally when R 0 > 1 by using the theory of Lyapunov functions, which implies that Ebola virus will persist in the host population and will eventually lead to epidemic. Finally, we have developed the numerical schemes of RK4 and NSFD methods for the proposed model (1)-(5) to acquire numerical solutions of the SEIR model. It was observed that the NSFD numerical scheme is more reliable than RK4. The RK4 scheme is a non-preserving numerical scheme, gives negative solutions, whereas the NSFD method preserves the nonnegativity and boundedness of all solutions for different values of step size h. Graphs of the state variables against time are presented for numerical analysis of the disease.