{'fraction_non_alphanumeric': 0.03208444737664502, 'fraction_numerical': 0.02361193096982534, 'mean_word_length': 4.620608899297424, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 6, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Mathematical modelling and analysis of COVID-19 epidemic and predicting its future situation in Ethiopia

The epidemic of the coronavirus disease 2019 (COVID-19) has been rising rapidly and life-threatening worldwide since its inception. The lack of an established vaccine for this disease has caused millions of illnesses and hundreds of thousands of deaths globally. Mathematical models have become crucial tools in determining the potential and seriousness of the disease and in helping the types of strategic intervention measures to be taken to prevent and control the intensity of the spread of the disease. In this study, a compartmental epidemic model of COVID-19 is proposed and analyzed to predict the transmission dynamics of the disease in Ethiopia. Analytically, the basic reproduction number is determined. To observe the dynamics of the system, a detailed stability analysis of the disease-free equilibrium (DFE) of the proposed model is carried out. Our result shows that the DFE is stable if the basic reproduction number is less than unity and unstable otherwise. Also, the parameters of the assumed model are estimated using the actual data of COVID-19 from Ethiopia reported for three months between March and June 2020. Furthermore, we performed a sensitivity analysis of the basic reproductive number and found that reducing the rate of transmission is the most important factor in achieving disease control. Numerical simulations demonstrate the suitability of the proposed model for the actual COVID-19 data in Ethiopia. In particular, the numerical simulation shows an increase in the rate of transmission leads to a significant increase in the infected individuals. Thus, results of the numerical simulations are in agreement with the sensitivity results of the system. The possible implication of this is that declining the rate of transmission to the desired level could enable us to combat the disease. Numerical simulations are also performed to forecast the disease prevalence in the community.


Introduction
The outbreak of coronavirus 2019 (COVID-19) pandemic emerged in December 2019 in Wuhan city, China, and was named a novel COVID- 19. The disease was declared as a global pandemic by the World Health Organization (WHO) on 11 March 2020 [1]. According to the WHO daily situational report [1], the outbreak of the virus remains rapidly increasing and life-threatening worldwide since its starting time. It is claimed in suffering more than 200 countries and territories worldwide [2]. Researches indicated that the main transmission route of COVID-19 is through people to people contact [3][4][5]. Its transmission shows almost an exponential growth in the early phase of the epidemic as daily global reported data indicated. It results in mild to high severe clinical outcomes around 20% of all the confirmed cases [6][7][8]. As of 12 June 2020, the COVID-19 causes for more than 7.4 million cases and more than 481 thousand deaths globally; and Africa shares 155,762 cases and 3,700 deaths. On 14 March 2020, the first 5 confirmed cases were detected in the Addis Ababa city, Ethiopia. In Ethiopia, the epidemic grows slowly near to mid of June 2020 since its onset on 14 March 2020 compared to the other countries. As of 12 June 2020, Ethiopia reported 2,915 cases and 47 deaths [1].
Since no medication for COVID-19 has yet been found, the WHO along with the concerned bodies continues to announce and release nonmedical potential intervention strategic measures to combat COVID-19 all the time intensively in addition to the vaccine discovery trials. Social distancing, contact tracing, quarantine, isolation, cleaning hands with soaps often, staying home, wearing a mask, and so on are some of the strategic intervention measures being practiced to mitigate and suppress the outbreak of the COVID-19 pandemic worldwide. Ethiopia adopts some of the intervention measures such as closing schools, lockdown international travels, restricted social gathering, and so on. However, administrative problems are observed in the strict implementation of the measures, and it causes a rapid increase in the infection.
Mathematical models have the potential to understand the dynamic transmission of communicable and non-communicable diseases and forecasting their long-run transmission behaviors [9,3,[10][11][12]. Mathematical models of disease dynamics could address the following questions: when will the epidemic be occurring? How many individuals will be infected next time? At what time will the epidemic peak be attained? How many individuals will get infected at the peak period? When will the epidemic last? Several scholars have developed mathematical models of the novel COVID-19 to study the dynamic transmission in Wuhan, China, and the other countries in the world [9,3,[10][11][12][13][14][15][16][17][18][19][20].
Recently, the authors [21] analyzed a Susceptible-Infected-Recovered-Dead (SIRD) COVID-19 using a discrete-time model assumption to forecast its outbreak. In this study, we propose a continuous nonlinear deterministic ordinary differential equation Susceptible-Infected-Recovered-Dead (SIRD) system to study the outbreak of the COVID-19 pandemic and further investigate its future trends in Ethiopia using COVID-19 daily reported data. In the process, we investigate the global stability of the disease-free equilibrium, perform the parameter estimation of the proposed model, and the sensitivity analysis of the basic reproduction number to the parameters. The rest part of this paper is presented as follows. Section 2 presents the formulation of the model and further studies the qualitative analysis. In Section 3, we study the basic reproduction number, the existence of a disease-free equilibrium, its local and global stability. Section 4 presents the estimation of the parameters. The numerical simulation and discussions are presented in Section 5. In the last section, a brief conclusion is given.

Model description
The total population size N at each time t with respect to their disease status is divided into four mutually-exclusive compartments of susceptible population (S), infected population (I), recovered population (R), and dead population (D). In the formulation of the method, we make the following assumptions: (i) The total population is closed. (ii) The population is homogeneously mixed. iii Recruitment and natural deaths in the population are not taken into consideration. (iv) The exposed compartment is not taken into account. v The recovered population develop permanent immunity. (vi) All variables and parameters of the model are non-negative.
Based on the assumptions and population scheme in Fig. 1, the differential equations that describe the transmission of COVID-19 become dS dt subject to the initial conditions where the total population of the model is given by The susceptible population is declined following the infection at a rate μ (the effective transmission rate of COVID-19 resulted from interactions of infected and susceptible human individuals per unit time) while the infected population is generated by same rate μ. The infected compartment is reduced by the recovery rate and disease-induced death rate, (β + γ), where, 1 β+γ is the average infectious period of the disease. Further, the population of recovery from COVID-19 infection is increased proportional to the infected individuals (I) at a rate β. The dead population is recruited at a rate γ. Descriptions of the state variables and parameters are presented in Table 1, and the flowchart diagram for the system is illustrated in Fig. 1.

Basic properties of the model
In this section, we show that all solutions of the proposed model with positive initial values remain positive.

Positivity of the solutions
For the system Eq. 1 to be biologically feasible, it is important to prove that all its state solutions are positive at all time t > 0.
Hence, T > 0. From the first equation of the system Eq. 1, we have Since T is the maximum of all time t, it can be written as It follows that

S(t)
Susceptible human individuals at a time t.

I(t)
Infected human individuals at a time t.

R(t)
Recovered human individuals at a time t.

D(t)
Dead human individuals at a time t.

N
The total size of human population. μ The rate of transmission from susceptible to infectious individuals. β The rate of recovery from infection. γ The disease induced death rate.
Using a similar method, it can shown that all other state variables of the model remain positive for all time t > 0.

Basic reproduction number
The basic reproduction number, R 0 , is the most important term in disease dynamics. It is used to quantify the transmission dynamics of an infectious disease. It is defined as the number of the secondary infected people produced by one infected individual introduced into a completely susceptible population. R 0 is used to decide if a disease persists in a community in a given population. If R 0 is less than one, each infected individual infects less than one individual on average and therefore the disease will cease to persist in the community. However, if R 0 is greater than one the disease persists and can spread [22]. The basic reproduction number, R 0 , is determined by the next generation matrix method. It can be determined from the general classical epidemiological models using [23] approach and summarized as follows: In addition, x, E, I denote the number of non-infected, latent, and infected individuals, respectively. Let the fixed point , 0) and it can be also further written as A = M − D, where M and D are the M-matrix (M ij ≥ 0,i ∕ = j) and diagonal matrix (D > 0), respectively. Hence, the basic reproduction number, R 0 , is given as the spectral radius of MD − 1 . That is, Using [23] approach, we can easily determine R 0 of system Eq. 1. That is, x = (S, R), I = (I, D). Furthermore, we have

We then get
Hence, using Eq. 3, one can easily find R 0 as R 0 := μS0 β+γ = μN β+γ , where S 0 ≊N is the initial population during the onset of infection.

Local stability of the disease-free equilibrium
In epidemiological modelling, we admit two equilibrium points in general: Disease-free equilibrium (DFE) and endemic equilibrium (EE). Disease-free equilibrium is in the absence of infected individuals in the population. On the other hand, endemic equilibrium is the equilibrium solution in the presence of infected individuals. The model presented in Eq. 1 has a DFE, computed by equating the right-hand side of system Eq. 1 to zero, given by E 0 = (N, 0, 0, 0). We have the following theorem.
Proof. The local stability of the system Eq. 1 is investigated using the jacobian matrix of the model Eq. 1 computed at E 0 , given by Following Eq. 4, the eigenvalues are λ 1 = 0, which is simple, and λ 2 = (β + γ)(R 0 − 1). This turns out that the disease-free equilibrium, E 0 , is locally asymptotically stable if R 0 < 1, and unstable if R 0 > 1. This completes the proof.

Global stability of disease-free equilibrium
In this section, we consider two conditions, (H 1 ) and (H 2 ), that would guarantee the global asymptotically stable of the DFE stated in [24]. Now, we write the system Eq. 1 in the form: where, x ∈ R m denotes the components of uninfected individuals and I ∈ R n denotes the components of infected individuals. Let U 0 = (x * , 0) be denoted the disease-free equilibrium of the system Eq. 1. Moreover, the conditions (H 1 ) and (H 2 ) below are supposed to be held true.
where, A = D I G(x * , 0) is an M-matrix and Ω is the region where the model makes biological sense.
Proof. To investigate the global stability of the DFE, we rewrite the system Eq. 1 in the form of Eq. 6. Then, x = (S, R) and I = (I, D). Moreover, we then have Since N⩾S⩾0, it is clear that Ĝ (x, I)⩾0. Thus, the disease-free equilibrium is globally stable.

Parameter estimation
Now we return to the estimation of the constant parameters given in the system Eq. 1. We now formulate the system as Here, y is the state variables and θ is the constant parameter values to be determined. To measure our fit to the real data, we define a least squares objective function y(i) represents the real data and y(i) is the solution to Eq. 7. To obtain the optimum parameter values we minimize the objective function .
The algorithm is presented below: 1. Guess initial parameter values p 0 . Set p = p 0 . 2. Using MATLAB version 2013a ode45 routine, solve Eq. 7 using p to find the solution y(i). 3. Evaluate error using Eq. 8. 4. Use p to minimize Eq. 9 using an optimization algorithm nlinfit to find the parameters with 95% confidence interval p and the covariance matrix. Update p =p. 5. Check for the convergence. If the convergence is not satisfied go to 2. 6. On convergence, set p =p.
Here we have used the inverse problem which requires an experimental data and an iterative solution. This inverse problem is an iterative procedure for solving the forward problem and hence, the ode45 solver will be nested within nlinfit.
Using the above algorithm the estimated parameter values are given in Table 2.
Using the estimated parameter values presented in Table 2, the basic reproduction number is R 0 = 8.45. This indicates every infected individual can probably infect on average 8 individuals per unit time.
Similar to [25], we also calculate the mean absolute error (MAE) and the root mean square error (RMSE). The MAE and RMSE are used to determine the accuracy of the predictions which are also called the performance metrics. The computed values of MAE and RMSE are 2.464511452058917e − 10 and 12.77, respectively. The values for MAE and RMSE show noticeably small compared to the results observed in [21,25]. This indicates that the model simulation fitted very well to the real data. The MAE is given as: and the RMSE: where N p represents the data sample size.

Model validation
Here, we use the graphical approach to validate the pattern. The basic principle is that if the residuals appear to be behaving arbitrarily, it is sound to conclude that the model matches the actual data. However, if the non-random structure is noticeable in the residuals, the model presumed is not good enough to represent the data. Or the data are not suitable for the model. Fig. 2 shows the fitting of the susceptible and infected population. Fig. 2a depicts the fitting of the susceptible model output to the real susceptible population and the COVID-19 model simulation to the COVID-19 observed data is presented Fig. 2b. Furthermore, the fitting of the recovery confirmed population to the recovery model simulation is seen in Fig. 2c. In all cases, the continuous line depicts the output of the model simulation, and the dotted points represent the actual observed data taken from Ethiopia. Even though some data points do not fit well, overall the model simulation and the actual data fit well.
Residuals plots for the model and the real observed data are presented in Fig. 2d. As given in Fig. 2d, the residuals seem to be random and the errors are reasonably small for the model. The possible implication of this is the model fits very well to the field data of Ethiopia. In general, the outcome of the model validation shows that the estimated parameters are reasonable and hence the assumed model fits well to the real data.

Uncertainties analysis
It is natural that some of the parameters in the disease dynamics maybe related and others be unrelated. An increase in one parameter may cause in an increase in the other parameter. During the parameter  estimation of the parameters some of these could covary. Covariance is the measure of how the parameter p i is related with p j . The covariance matrix approximates the uncertainties in the model parameters. The interrelations between the parameters is based on the signs. If the sign of the covariance is positive, then an increase in one parameter causes in an increase in the other parameter. If the covariance between any two parameters is 0, then parameters are not related. However, if the covariance between two parameters is negative, then an increase in one parameter will cause in a decrease in the other parameter. The diagonal elements of the covariance matrix is the variances of the parameters and the off-diagonal elements are the covariances between all potential pairs of parameter. The matrix Π is positive semidefinite and symmetric. The covariance matrix Π for the model is: The covariance between the transmission coefficient rate (μ) and the recovery rate (β), the transmission coefficient rate (μ) and the rate of death due to infection (γ) is negative. The possible implication is that an increase in one parameter value may result in a decrease of the other. On the other hand, the covariance between the recovery rate (β) and the rate of disease induced death rate (γ) is zero. This indicates that an increase/or decrease in one parameter does not cause in an increase or decrease of the other parameter. The summary of the covariance of the model is depicted in Table 4.

Sensitivity analysis
Sensitivity analysis is the vital tool to investigate the importance of each parameter for the system Eq. 1 in prevalence and transmission of the disease. This tool is usually applied to determine the robustness of model predictions to parameters as errors occurred in collected data and presumed parameter values [26]. In particular, sensitivity analysis is used to identify parameters of the system that would have great influence to the basic reproduction number, described in Section 3.2, and should be targeted by intervention strategies through the stakeholders. To this end, we employ the normalized forward sensitivity index of a variable for the system relative to each parameter, that could be defined  Table 3 and the residuals plot.

Table 4
Covariance of the parameters.

Cμ,γ Negative
Cβ,γ Zero as the ratio of the relative change in the variable to the relative change in the parameter [26,27].
Definition 4.1. [26][27][28]. The normalized forward sensitivity index of the basic reproduction number, R 0 , which is differentiable with respect to each parameter of the system is defined by Following Eq. 13, we calculate the normalized forward sensitivity index of R 0 with respect to each parameter as follows: The estimated values of model parameters and sensitivity indices of R 0 are presented in Table 5 and Table 6, respectively.
Note that one should carefully estimate the sensitive parameter, since a small change in the parameter causes a big quantitative changes. In contrast, it could be given less emphasizes in estimation of a parameter with small value for the sensitive index, because a small variation in that parameter could not bring a significant quantitative changes.
The sensitivity index of R 0 remains independent of any parameter with respect to the parameter μ. On the other hand, it depends on the parameters β and γ with respect to the parameters β and γ, as it can be seen in Eq. 14. Hence, Eq. 6 shows that R 0 is highly sensitive to the parameters μ and β. It is clear that R 0 attains its maximum value with the parameter, μ. In particularly, an increase (decrease) the value of μ by 100% will increase (decrease) R 0 by same value 100%, and an increase (decrease) the value of β by 100% will increase (decrease) the value of  Table 6 Sensitivity indices of R 0 calculated at the estimated parameter values presented in Table 5.  Table 2.
R 0 by 91.42%. On the other hand, an increase (decrease) the value of γ by 100% will increase (decrease) R 0 by 8.58%, respectively.

Numerical simulation and discussion
Numerical results of system Eq. 1 are discussed in this section. Based on the real data collected from 14 March to 12 June 2020 in Ethiopia presented in Table 3, we estimated the parameters' values of the system Eq. 1, and perform numerical simulations to investigate the impact of the most and least sensitive parameters at the estimated values of parameters. We also carry out numerical simulations and predictions using the estimated parameters' values. For simulation purpose, we set the following values as initial conditions: S 0 = N − 1, I 0 = 1, R 0 = 0, and D 0 = 0, where N is the total size of Ethiopian population given by N = 115 millions [2].

The effect of the most sensitive parameters
This subsection is devoted to show the numerical simulation results are in accordance with the sensitivity analysis results in Section 4.3. As can be observed, for example, from Fig. 3b a small increase in the transmission parameter results in a great effect in the transmission of the disease in the community which is in agreement with the result in Table 6. Fig. 3 shows the solution trajectories of the system by varying the rate of the transmission of the disease, μ, and the recovery rate, β. An increase in the rate of transmission μ causes a decrease in the susceptible population (Fig. 3a). However, an increase in the rate of this parameter increases in the dynamics of the disease (Fig. 3b). More precisely, increasing μ from 0.000000000376963 to 0.000000000476963 and 0.000000000576963 resulted in significant decline in the susceptible individuals. In contrast, increasing the parameter μ from 0.000000000376963 to 0.000000000476963 and 0.000000000576963 causes in significant increase in the infected population. On the other hand, Fig. 3c reveals that decreasing the rate of recovery β from 0.007176440779023 to 0.006176440779023 and 0.005176440779023 results in significant decrease in the number of susceptible individuals. Conversely, the infected population increases while the rate of recovery decreases, as depicted in Fig. 3d. Thus, our results are in agreement with results reported in the sensitivity analysis (Table 6). Therefore, to minimize or eradicate the disease from the population policymakers are recommended to minimize the rate of transmission of the disease.

The effect of the least sensitive parameter
Based on the results in the sensitivity analysis, the disease-induced death rate of COVID-19 is the least sensitive parameter. The numerical simulation of the system by varying this parameter, γ, is presented in Fig. 4. It is evident from Fig. 4a and Fig. 4c that the susceptible and recovered populations associated with decreasing of the disease induced death rate is not significant. Precisely, decreasing γ from

Numerical simulations using the estimated parameters
The behavior of state variables in the proposed model using the estimated parameters is described in Fig. 5. Fig. 5a represents the human susceptible population to CODID-19. Fig. 5b depicts the population infected with COVID-19. It shows that the infected population with COVID-19 is small at the beginning and are monotonically increasing eventually. Fig. 5c describes the humans size recovered from COVID-19. The number of recovered individuals continues to be high in number as the number of days increases. Moreover, Fig. 5d exhibits that the number of dead individuals become monotonically increasing.

COVID-19 forecasts for Ethiopia
Ethiopia declared its first confirmed COVID-19 case on 14 March 2020 [1], and continues reporting daily new cases. Prediction of the pandemic is most essential to take effective and efficient strategic intervention measures. Hence, this analysis aims to investigate the trend of the pandemic COVID-19 in Ethiopia using the value of the estimated parameters in Table 2. It attempts to carry out a numerical simulation of the system Eq. 1 for the prediction of COVID-19 for 60 days. To achieve this, we use the final values of the state variables as an initial condition of the system Eq. 1. The number of state variables for the next 60 days is depicted in Fig. 6. It can be observed from this that the number of susceptible individuals decreases monotonically (Fig. 6a). And the number of infected, recovered, and dead individuals are increasing rapidly (Fig. 6b, Fig. 6c, and Fig. 6d). According to this numerical simulation, policymakers need to enforce preventive measures such as wearing a face mask and disinfecting and decontamination of contaminated areas to flatten the curve. The analysis from the prediction of our system reveals high variations in the cumulative number of infected, recovered, and died individuals relative to the official data published on August 11, 2020 [2]. On the 11th of August 2020, the system predicted 113,877 cumulative number of infections, while the reported cumulative number of infected individuals was 24,175. On 11 August 2020, the projected cumulative numbers of recoveries and deaths were 13,951 and 1,310, respectively, while the confirmed cumulative numbers recorded were 10,696 for the recovered and 440 for the dead individuals. The variation between the real reported cumulative data and the forecasted cumulative value could be attributed to the lack of inclusion of protection measures in the proposed system.

Conclusion
In this paper, we proposed a SIRD simple mathematical model to address the current COVID-19 transmission dynamics and predicting its future trends in Ethiopia. We carried out the local and global stability analysis of the disease free-equilibrium point E 0 in case of R 0 < 1 based on the geometric method presented in [29]. Based on the reported data of COVID-19 in Ethiopia, we also estimated the value of the parameters transmission rate, the recovery rate, and the disease-induced death rate.  Table 2. A graphical method is applied to validate the model. Our result shows that the proposed model fits well with the real data reported in the country for 91 days. It is found out from the sensitivity analysis that the rate of transmission parameter is the most sensitive parameter. The epidemiological implication of this is that an increase in this parameter could aggravate the disease in the community. On the other hand, the disease-induced death rate is the least sensitive parameter. Meaning an increase or decrease in this parameter does not result in a substantial spread of the disease. Furthermore, the reproduction number estimated for the case of Ethiopia is found out to be R 0 = 8.45. The biological significance of this finding is that the elimination of COVID-19 may not be achieved in the population in a short period. Unlike the authors [21] who have considered a discrete differential equation with the frequencydependent force of infection, we proposed an ordinary differential equation with a linear force of infection. The numerical experiment, based on the estimated parameter values, of the system exhibited a small number of infected individuals at the beginning and are increasing eventually. Numerical results of the model in the prediction also indicated a significant variation in the cumulative number of infected variables relative to the cumulative number of confirmed reported cases. The variation may be resulted from the lack of incorporating the intervention strategic measures to the system. In addition, our findings on the short-term prediction of the disease show that the disease will aggravate in the community. This result is in agreement with results obtained in [30].
The study has limitations. First, the values of the parameters defining the system are determined by applying field data to COVID-19.
Estimating the parameters is a difficult task due to the absence of a significant part of the infectious cycle. The most important limitation is the use of ode45 and the optimization method nlinfit to find the optimum parameters. nlinfit is an optimization strategy that can obtain both positive and negative optimum parameters, but the required parameter values are non-negative. In this case, we have to simulate the system using different initial values of parameters until we get positive optimum parameter values. In addition, this paper considers a deterministic model with a homogeneous population mixing assumption, ignoring the stochastic model assumption and the spatial distribution of the population. It has also failed to include intervention measures like vaccination, quarantine, and social distancing. Although our assumed system answers relevant questions, extending the system into an agestructured and model with implementing vaccination may result in a wide range of global dynamics [31][32][33]. Therefore, incorporating these into the system would be an exciting area of future research to foresee their impacts on the transmission and mitigation of the epidemic in Ethiopia.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.  Table 2 when the initial population is (S 0 ,I 0 ,R 0 ,D 0 ) = (114999781, 193, 24, 2).