key: cord-0754007-pqkw9g1s authors: Abuhasel, Khaled Ali; Khadr, Mosaad; Alquraish, Mohammed M. title: Analyzing and forecasting COVID‐19 pandemic in the Kingdom of Saudi Arabia using ARIMA and SIR models date: 2020-10-05 journal: Comput Intell DOI: 10.1111/coin.12407 sha: 7117b45c004ac818a537dad852dd46b1d6fdf05b doc_id: 754007 cord_uid: pqkw9g1s The novel coronavirus COVID‐19 is spreading all across the globe. By June 29, 2020, the World Health Organization announced that the number of cases worldwide had reached 9 994 206 and resulted in more than 499 024 deaths. The earliest case of COVID‐19 in the Kingdom of Saudi Arabia (KSA) was registered on March 2 in 2020. Since then, the number of infections as per the outcome of the tests increased gradually on a daily basis. The KSA has 182 493 cases, with 124 755 recoveries and 1551 deaths on June 29, 2020. There have been significant efforts to develop models that forecast the risks, parameters, and impacts of this epidemic. These models can aid in controlling and preventing the outbreak of these infections. In this regard, this article details the extent to which the infection cases, prevalence, and recovery rate of this pandemic are in the country and the predictions that can be made using the past and current data. The well‐known classical SIR model was applied to predict the highest number of cases that may be realized and the flattening of the curve afterward. On the other hand, the ARIMA model was used to predict the prevalence cases. Results of the SIR model indicate that the repatriation plan reduced the estimated reproduction number. The results further affirm that the containment technique used by Saudi Arabia to curb the spread of the disease was efficient. Moreover, using the results, close interaction between people, despite the current measures remains a great risk factor to the spread of the disease. This may force the government to take even more stringent measures. By validating the performance of the applied models, ARIMA proved to be a good forecasting method from current data. The past data and the forecasted data, as per the ARIMA model provided high correlation, showing that there were minimum errors. reproduction number. The results further affirm that the containment technique used by Saudi Arabia to curb the spread of the disease was efficient. Moreover, using the results, close interaction between people, despite the current measures remains a great risk factor to the spread of the disease. This may force the government to take even more stringent measures. By validating the performance of the applied models, ARIMA proved to be a good forecasting method from current data. The past data and the forecasted data, as per the ARIMA model provided high correlation, showing that there were minimum errors. ARIMA model, COVID-19, forecasting, Saudi Arabia, SIR model INTRODUCTION The first case of the ongoing epidemic of coronavirus (COVID- 19) disease was first reported in Wuhan City, in China in early December 2019. 1 In the subsequent weeks, this virus spread inside China and proceeded to spread widely in other countries and infected more than 1 996 681 people and taken more than 127 590 lives as on the second week of April, 2020. 2 Estimating trends in the spread of any infection over time may provide information and insight into the epidemiological situation and ascertain if outbreak control strategies have a significant effect. 3, 4 Such approaches may provide decision-makers with the scenarios of the expected potential future progress and help to quantify risks and direct mitigation strategies. Models of epidemiologic research assume an exponential increase, 5 however, some other methods for estimating and predicting the rate of infections are based on mechanistic-statistical approaches. 6, 7 Moreover, mathematical and statistical modeling can also be beneficial in producing accurate, short-term forecasts of confirmed cases. 8 Several research articles have been released in recent weeks on epidemic prediction of COVID-19 [9] [10] [11] [12] [13] [14] and most of the researchers concentrated on developing a new framework based on artificial intelligence tools. [14] [15] [16] [17] Time-series forecasting models, based on statistics and stochastic processes, forecast epidemic outbreaks, mainly involving autoregressive integrated moving average, exponential smoothing method, Gray model, and Markov chain model. 12, 18, 19 Since the initial spread of the COVID-19 epidemic, many approaches have been developed and applied to assess and forecast potential numbers of the confirmed cases. These attempts are based either on data mining mostly or on health models to investigate outbreak data and to predict the range of new confirmed cases. The forecasting should provide users with the most likely outcome and probabilities that different scenarios could occur with an objective level of confidence. Models that simulate and forecast the spread of COVID-19 are essential tool that can improve the informed decision taking. These models should provide users with the most possible outcome and probabilities with a reasonable degree of confidence that various scenarios will occur. Various studies have been on prediction of the COVID-19 pandemic have been brought forward. One of these studies is DeFelice et al 20 provided a prediction model that uses the spillover transmission risk and also the West Nile virus (WNV) case to predict the COVID-19 in the form of a compartmental model. In this model, the predicted trend was based on data that had been previously collected on the trend of WNV, specifically in Long Island and New York between the years 2001 and 2014. Moreover, 21 to ensure that there was consistency in the model, there were parallel comparisons between the model and other models such as time series outcomes on Hepatitis A virus. Furthermore, there was a 13-year data from the Turkey infection that had been modeled using a multilayer perceptron (MLP), which was also used to compare the reliability of the model. Other models, such as time-delay neural networks, radial basis function, and ARIMA were also part of the prediction process. In the comparison process, the MLP proved to more accurate than the other models. To further ensure that the comparison was robust, the MLP was tested for short-term data, such as 6-years data on influenza ad also weekly data on the same influenza. 22, 23 Various researchers have proposed different trend on the COVID-19 cases and prevalence rates since the first case was reported. Such predictions have been made public by being published in platforms such as websites, newspapers, and also be propagated by media companies. Although some of the researchers have employed almost similar models, the prediction proves a wide range of differences. [24] [25] [26] [27] [28] [29] [30] [31] There has been a rich history as far as research on epidemics in the field research of physics. 32 Moreover, other mathematical models have also been employed in fitting the data on the spread of epidemics. The balancing of these models provides a balanced system that can be used to provide a reliable prediction as far epidemic spread based on past data is concerned. Due to the usability of these methods, the models have been proposed as appropriate models for the COVID-19 pandemic, which is currently the novel issue that need to be predicted to put down the best protective measures. When the Kingdom of Saudi Arabia (KSA) reported the first case of the COVID-19 pandemic, the country undertook several containment strategies that have been part of the parameters in the prediction method due to their ability to lessen the rate of infection. The main measure was to ensure that the movement of people is halted by putting people under lockdown to reduce the frequency of interactions. There are people of interest that the government is mostly concerned with, including the elderly and people with underlying problems that may adversely affect their immune system. Using all the variabilities that come with this pandemic, the main point of focus in this research is to use the data from already infected persons and use such data to create short-term prediction on how the pandemic is shaping itself in the country. These models will be fed and utilized using the ARIMA model due to its reliability compared with the other models that have been used in the past to establish future trends of epidemics and pandemics. This work use the classical SIR model (Figure 1 ), popularly described as the Kermack-Mckendrick to describe the transmission of the COVID-19 virus pandemic. 33 The SIR model is a compartmental model for modeling how a disease spreads through a population. It is an acronym for susceptible, infected, and recovered. The population under study is presumed to be invariant in which three groups are derived from it. In this model, natural death and birth rate are not considered because the outstanding period of the disease is much shorter than the lifetime of the human. The model that will be employed is as follows: . (1) In the model above, the response variables are the number of people in the population that is at risk in a given time t, the number of those people in the population that have the capability of spreading the virus, and the number of individuals that have recovered in a given time t. Moreover, the response variables must satisfy the condition that their sum is equal to the population. On the other hand, all the response variables assume nonnegativity as the initial conditions. The contact rate among the members of the population is represented by , while the mean infections period is represented by 1/ . For the population at risk, model (S(t)) narrows to>. while the population that recovers model (R(t)) is an increasing function represented by: If R 0 > N/S 0 , then number of infectious people I(t) rise to a maximum value. If R 0 < N = S 0 , then I(t) decreases to 0. For more details about SIR, see Reference 34. During the initial period, before the reporting of the first case, the number of people at risk can be assumed to be equal to the whole population under study so that the ratio become 1 for the initial condition. The SIR model cannot be used in isolation since it does not take into account the number of deaths due to the pandemic although such occurrences are inevitable. Therefore, the model is being combined with appropriate models, such as the death equation so that the number of death can also be accounted for. When the SIR model is combined with the death model, the following equation are produced. where, the number of deaths in the death equation is represented by D. From the combined equation, R(t) + D(t) represents the number of individuals removed from the population. A description of the probabilistic structure in a sequence is known as the stochastic process. 35 There have been various studies by different researches that aim at describe the probabilistic structure of sequences. One of the main model is the ARIMA which was developed by Box and Jenkins. 36 A process whose statistical properties do not change with time, that is, process with constant mean and constant variance, known as a stationary process, is a crucial collection of stochastic processes. Mathematically, the joint distribution of Y (t 1 ), ...,Y (t k ) and Y(t 1 + ),......,Y(t k + ) is the same for all t 1 , t 2 , … , t k of a stationary process. Simply put, shifting the origin of time by a quantity does not change the statistical properties of the process. Usually, dealing with real-time data, most time series does not exhibit stationarity in nature as they have no fixed mean. The features of the vital collection of that shows stationarity in the differences between time series models are generally referred to as ARIMA models. Autoregressive in the models are the time lags while the error, which are associated with these models are known as moving averages. For the model to be modeled into a time series model, a term by term differentiation is executed to make the model stationary. 37 The ARIMA may be described by specific parameters depending on whether or not it is seasonal. In this case ARIMA (p, d, q). 38, 39 The following equation are associated with the ARIMA model above: In the equations above, Y t represents the observed series, while B acts as the operator for backshift process. Moreover, t is a sequence of errors that are normally and independently distributed N(0, 2 ). The parameters for whether or not the autoregressive is seasonal is represented by while the same parameter for moving averages is represented by . Furthermore, parameters p and q are the powers of and , respectively. In addition, d denotes nonseasonal and differences. In real-time data, taking the first difference (d = 1) is usually found to be sufficient and occasionally second difference (d = 2) would be enough to achieve stationarity. While developing the Time series model, there are three stages, which include identification, approximation, and debugging. 40 This is the first step in developing the model, which entails analysis of the data to check for its stationarity and normality in its distribution. 41 If the model does not prove to be stationary, it can be made stationary by data transformation, such as using the log transform. After the model has been examined for both stationarity and normality, the second step involves checking for the existence of correlation and autocorrelation between the data as well as partial autocorrelation. After that, the parameter for nonseasonality and their orders are found so that the model can be fully developed. The second stage, which involves approximation, several models are that seem efficient and a clear representation of the developed model are chosen to represent the model. The parameters of the approximating models are developed form the main model through least squares. Using different combinations of AR and MA models, the result can be a series of different models that can be used for the forecasting process. 42 However, an ideal model from the combinations must satisfy certain conditions. To check the usability of the model, the following diagnosis are executed: (i) The efficient model should have the minimum Akaike information criteria (AIC) or Bayesian information criteria or the Schwarz-Bayesian information criteria (SBC). AIC = (−2 log L + 2 m). In which, the value of m is gotten by summing P, Q p, and q, while L is the chance function. Alternatively, AIC is given by {n (1 + log 2 ) + n log 2 + 2 m}. The constant is eliminated and therefore the model lacks the first term, which happens to the constant. When AIC is omitted, SBC can be considered, which is represented as: After the ARIMA model, which is considered appropriate among the alternatives, is put in place, it can be tested for a goodness fit, which entails testing its efficiency. The model is assumed to be a good fit if the residuals are approximately equal to the white noise. The most essential tools are the plots of ACF and PACF. These plots can also be used to test whether or not the model is appropriate by checking whether the coefficients are tending to zero. (ii) Tests normality of the residuals A plot that can be used to show how the data is distributed can also be used to investigate the distribution of errors. The graphs that can be used for this purpose include histograms and the normal probability plots. If the residuals are normally distributed, the normal probability plot shows the data in a straight line. However, a few residuals may lie outside the straight line, which may be seen as the line of best fit. A single variable time-series dataset was sampled to analyze a real-time forecasting of the novel COVID-19 cases in the KSA. The set of data include the numbers of daily confirmed diagnosis, recoveries and fatalities cases from March 2, 2020 to June 30, 2020 as reported by The Saudi Ministry of Health 43 ( Figure 2 ). The models presented in this article were performed using Matlab. Initially, the total number of infections after time t are represented by: In estimation, the function used will be: In which the constant a, b, and c are the factors used for estimation in the least square method. Using the information above, and with TNI(t) we obtain: From Equation (12), TNI(t) = b e ct − a, in which a = 13 370; b = 6852; and c = 0.062 (Figure 3 ). This function is fitting the observed data with R 2 = 0.9985. When = 1/14, we have that t 0 = 23; Figure 4 illustrate the results from the SIR models with predefined cases in which the predictions of the infection in terms of the rates and peaks. In Figure 4 , an illustration of how the SIR is fitted to the model is shown. Figure 4 shows that the number of infections are increasing between March 2, 2020 to June 30, 2020, while the number of those who recovered is at the rate of 13%. The data is in a span of 120 days due to the date that the first case was reported and the time the model is being developed. The prediction model was efficient in producing data that can be used as the expectation of the number of infection in the country in the days to come. According to the results, the peak number of infections seem to be up to 4 × 10 6 . However, the results considers the measures that the government has put in place to curb the spread of the pandemic. From the model as of 30.06, the estimated number of infections was 4842, while the number of infection were 3943. The estimation on the stochastic processes were reproduced using a reproduction number of R o = 2.4, as shown in Figure 5 . When the reproduction number increases, the models shows that the number of infections are bound to increase up to 75% with the main reason being that there is a population of asymptomatic people, which can be seen as a career and are capable of spreading the disease to members of the society that are at high risk of being defeated by the virus. Figure 6 shows that the number infected people from the prediction model can be up to 7 × 10 6 where an R o = 2.4 is used. In addition, the rate of infection and that of recovery 0.171 and 0.173, respectively. When the prediction model is compared with the actual values, the number of infections seem to have reduced by the same rate that the model provided as an increase level, while the actual peak showed a level up to 40% of the whole population. 44 It is advisable to plot the associated data and pay attention to the unique features exhibited by a time series. It gives guide for choosing an appropriate modeling approach that directly captures identified features. Before starting the procedure, there is a need to make the time series stationary. To stabilize the variance, we used square root transformation on the time series of the infected number of cases per day. For investigating the stationarity of time series, we take the support of the KPSS and ADF test, and results are shown in Table 1 . The first difference of series, that is, d = 1, is optimum, as mentioned earlier, to make series reasonably stationary. Based on a 5% significance level both the tests, the Dickey-Fuller, which is an augmented model, the Schmidt, Kwiatkowski, Shin test (KPSS), and Phillips reject the hypothesis of stationarity of time series without making any difference. Afterward taking the first difference, both the criteria agree on the stationarity of time series. Furthermore, the ACF and PACF series' first difference and square root transformation are used to estimate another two parameters of the candidate model. From Figure 6A ,B, the ACF display one spike, and the PACF also displays one significant spike. As seen on of the series is showing a decaying trend in a sine-wave sequence. The series is the ACF. Considering the PACF, the starting four values are essential. This shows that the general model can be expressed as a combination of the two models, TA B L E 1 T A B L E 2 AIC vs SBC for some of the qualified models that is, MA and AR. Visual inspection of the ACF and PACF can be quite inconclusive and misleading. This can be seen by simulating ARMA models and observing the unexpected inconclusive appearances of their ACFs and PACFs. Therefore, to select the most appropriate model, we investigated combinations of models with p = 1:6, and q = 1:6 to select the most appropriate model among the qualified models. Results from the investigated models together with corresponding AIC values are in Table 2 . In this model, ARIMA (1,1,1), which showed the minimum AIC, meaning that its residuals of satisfy the condition to a large extent compared with the other models is used. We examine that all the residuals are scattered around zero mean with constant variance. After identifying the models and estimating the parameters, the next stage of verification follows. This stage serves to check the patter followed by the residuals and whether or not they form a specified trend. From the model above, the residuals seem to be within the required limits in which they are independently and normally distributed ( Figure 7) . Moreover, the residuals are equal to the models white noise. Normality is shown since the plot shows a typical straight line through most of the residuals with a few of them lying outside the line (Figure 8 ). From the results, it is sufficient to accept the null hypothesis that the series of residuals displays no relations among them. When developing ARIMA models, the sole reason is to predict a variable that correspond to the given data. In this case, the model take two approaches, namely, the sample period forecast and the postsample period forecasting. In the model discussed, the prediction was for a period of 1 day lead time as shown in Figure 9 . The predicted values and the actual cases seem to have a close relation. Furthermore, the root mean square error, the mean absolute deviation, R-squared (R 2 ), and Nash-Sutcliffe coefficient (E) were used to evaluate the results of the developed models. The forecast is done with 1-day to 10-days lead-time to investigate the forecast accuracy for longer lead-time. Results of model evaluation are shown in Table 3 where it can be observed that the forecast accuracy decreases with longer lead-time. Table 3 clearly shows that the ARIMA models performance is more satisfactory as it results in minimum values of the residuals, therefore providing almost the same values as the observed ones. Two important models have been considered in the prediction of the infections and recovery rates of the novel COVID-19 pandemic. First, the classical deterministic SIR model, deals with fitting and the analysis of the COVID-19 spread in Saudi Arabia. The second model is the ARIMA model, which is stochastic process used to forecast whether or not the model is linear or stationary as exhibited in the one-variable time series data and is used to show the exact behavior of actual numbers associated with the daily COVID-19 cases. The two models were performed using the daily data from March 3 in 2020 until June 30 of the same year. Using this model, the reproduction number was 1.86. On the other hand, the measures taken by the government to curb the spread of the pandemic was tested for effectiveness using SIR model. Moreover, the repatriation plan as of the current data was also assessed showing a higher reproduction number than that of the measures taken by the government. Therefore, using the results from the data, intervention measures were more effective in controlling the spread of the disease in the country. However, if the rate of interaction between members of the population remain uncontrolled, the number of infection would hike. This will be facilitated by the existence of asymptomatic persons that will transmit the disease to other in the population. The model was chosen since it is much accurate than other candidate models. In this model, the predicted values had a good correlation with the observed values showing a decrease level of errors while using it. Due to its effectiveness, the forecasting model based on the stochastic process, such as ARIMA proved to be a good predictor of the infection in the pandemic and capable of being used in the medical field to facilitate the planning of how the number may trend. Multiple-input deep convolutional neural network model for COVID-19 forecasting in China, medRxiv The novel zoonotic COVID-19 pandemic: an expected global health concern Early dynamics of transmission and control of COVID-19: a mathematical modelling study A cascade of causes that led to the COVID-19 tragedy in Italy and in other European Union countries Fractal kinetics of COVID-19 pandemics (with update 3/1/20), medRxiv Novel coronavirus 2019-nCoV: early estimation of epidemiological parameters and epidemic predictions, medRxiv Estimating the number of infections and the impact of non-pharmaceutical interventions on COVID-19 Real-time forecasts of the COVID-19 epidemic in China from Neural network based country wise risk prediction of COVID-19 Epidemic analysis of COVID-19 in China by dynamical modeling Data visualization analysis and simulation prediction for COVID-19 Prediction and analysis of coronavirus disease Predictions of 2019-ncov transmission ending via comprehensive methods Probabilistic prediction of COVID-19 infections for China and Italy, using an ensemble of stochastically-perturbed logistic curves Finding an accurate early forecasting model from small dataset: a case of 2019-ncov novel coronavirus outbreak AI-driven tools for coronavirus outbreak: need of active learning and cross-population train/test models on Multitudinal/multimodal data Optimization method for forecasting confirmed cases of COVID-19 in China Time series analysis and forecasting of chlamydia trachomatis incidence using surveillance data from Epidemiological features and time-series analysis of influenza incidence in urban and rural areas of Shenyang, China Ensemble forecast of human West Nile virus cases and mosquito infection rates Comparison of four different time series methods to forecast hepatitis a virus infection Forecasting seasonal outbreaks of influenza Real-time influenza forecasts during the 2012-2013 season Short-term forecasts of the COVID-19 epidemic in Guangdong and Zhejiang Forecasting the novel coronavirus COVID-19 Estimation of the reproductive number of novel coronavirus (COVID-19) and the probable outbreak size on the Diamond princess cruise ship: a data-driven analysis Predicting the cumulative number of cases for the COVID-19 epidemic in Real-time forecasts and risk assessment of novel coronavirus (COVID-19) cases: a data-driven analysis Modeling and forecasting the early evolution of the COVID-19 pandemic in Brazil COVID-19 disease outbreak forecasting of registered and recovered cases after sixty day lockdown in Italy: a data driven model approach Forecasting of COVID-19 confirmed cases in different countries with Arima models, medRxiv Epidemic processes in complex networks Analysis of the COVID-19 pandemic by SIR model and machine learning technics for forecasting The SIR model and the foundations of public health Time Series Analysis: Forecasting and Control Streamflow drought time series forecasting Analysis of the trade pattern of Pakistan: past trends and future prospects Time Series Analysis and its Applications Meteorological drought forecasting using stochastic models Time Series Analysis: Forecasting and Control SAS for Forecasting Time Series Antieigenvalues and antieigenvectors in statistics The Saudi Ministry of Health Epidemic situation and forecasting of COVID-19 in Saudi Arabia using the SIR model, medRxiv; 2020. How to cite this article The authors extend their appreciation to the Deanship of Science and Research at University of Bisha Saudi Arabia for funding this work through COVID-19 Initiative Project under Grant Number (UB -COVID-07-1441). The authors declare no conflict of interest. The data that support the findings of this study are openly available at https://www.moh.gov.sa/. Khaled Ali Abuhasel https://orcid.org/0000-0003-4419-3708