Coronavirus disease (COVID-19) is a severe acute respiratory syndrome that occurred for the first time in Wuhan, the capital city of Hubei Province in China, in December 2019. In late February, the virus was detected even in Europe, and from that moment on, it has continuously spread across the world. As of 31 May, according to Worldometer's COVID-19 data, the pandemic affected 213 countries and territories, and two conveyances, with more than 6,2 million confirmed cases and a death toll higher than 370,000 people. The main epicenters of the pandemic are mostly advanced countries, i.e. Brazil, France, Germany, Italy, Russia, Spain, the UK, and the USA. However, in the last few weeks, other relatively poorer countries, such as India, Peru, and Turkey, are progressively climbing the world rankings. When an epidemic occurs, one crucial issue is determining its evolution and inflection point. Therefore, the aim of this paper is to provide a short to medium-term forecast of the spread of the COVID-19 disease, and its inflection point, in three of the most affected countries worldwide, Italy, Russia, and the USA. The prediction will be estimated by using an autoregressive integrated moving average (ARIMA) model, applied to the number of new daily confirmed positive cases of COVID-19. The specific time frame for each country will be chosen by considering, like departure point, the moment when daily cases began to show signs of stabilization. For the remainder, this paper is organized as follows. In section 2, I will introduce ARIMA model and related literature. In section 3, I will present the data used to forecast and discuss the empirical strategy. In section 4, I will discuss the results. In section 5, I will provide some concluding considerations and possible policy implications.
In the last few months, an increasing body of literature has attempted to forecast the incidence and prevalence of the COVID-2019 pandemic by using different approaches, such as the ARIMA (Benvenuto et al., 2020; Chakraborty and Ghosh, 2020; 2 Ribeiro et al., 2020 3 , Singh et al., 2020 ), the exposed-identified-recovered (EIR) model (Xiong and Yan, 2020) , the susceptible-exposed-infected-recovered (SEIR) model (Wu et al., 2020) , the segmented Poisson model (Zhang et al., 2020) , the SIDARTHE 4 model (Giordano et al., 2020) , the susceptible-infected-recovered (SIR) model (Nesteruk, 2020) , the SIR/death model (Anastassopoulou et al., 2020; Fanelli and Piazza, 2020) , and mathematical methods based on travel volume (Tuite et al., 2020) . The ARIMA model is clearly one of the most preferred because of its good properties. In fact, it is quite easy to fit, manage, and its mathematical interpretation is easy and immediate also for non-academics. It was introduced for the first time, as the Box-Jenkins approach, by statisticians Box and Jenkins in a highly influential seminal work published in 1970. It could be considered one of the most used prediction models for epidemic time series (Rios et al., 2000; Li et al., 2012; Zhang et al., 2014; Zheng et al., 2015) , and it is frequently used with non-stationary time series to capture the linear trend of an epidemic or disease. It allows the prediction of a given time series by considering its own lags, i.e. the previous values of the time series, and the lagged forecast errors. Moreover, it is very flexible and can be easily adapted to any kind of data, considering trend, cyclicity, seasonality, calendar variation, randomness disturbances like other diseases, external or exogenous interventions, outliers, and other relevant real aspects of time series (Pack, 1990; Barnett and Dobson, 2010) . Therefore, the ARIMA models allow in a simple way to investigate COVID-2019 trends, which are currently of huge economic and social impact, by helping the health authority to continuously monitor the epidemic and to better allocate the available resources. In Table 1 , I selected 19 studies on the prediction of the spread of various diseases by using an ARIMA framework. 
The analysis is carried on the number of new daily confirmed positive cases of COVID-19 in Italy, Russia, and the USA, collected by the data aggregator website named Worldometer. 5 As suggested by several authors (Box and Tiao, 1975; McCleary et al., 1980; Box et al., 1994) , a reasonable ARIMA model requires at least 40-50 observations. Therefore, to meet this criterion, I considered a minimum of 53 and a maximum of 69 data points. The specific investigated time window for each country was chosen by identifying the period when daily cases began to stabilize, which can be approximately seen as the peak of the epidemic (Figures 1, 2 In particular, to forecast the COVID-19 outbreak spread in Italy, Russia, and the USA, I used a nonseasonal ARIMA model. In fact, the autocorrelation function (ACF) and the partial autocorrelation function (PACF) of the raw time series showed a lack of any seasonal significative pattern. 6 The nonseasonal ARIMA models are generally classified as "AR-I-MA (p, d, q)", where: i) p is the order of autoregressive terms (AR); ii) d is the order of nonseasonal differences (I) required to make the time series stationary; and iii) q is the order of moving average (MA). The optimal parameters for ARIMA models were chosen by considering the following four sequential steps: 7 i. First, I applied the auto.arima() function, included in the "forecast package" of the widely known R Project for Statistical Computing, developed by Hyndman and Khandakar (2008) . This algorithm allows to identify the best order of an ARIMA process by considering a unit root test to identify the appropriate degree of differencing, 8 and the minimization of the corrected Akaike's Information Criterion (AICc) 9 and the log-likelihood function (MLE) to identify the AR and MA parameters; ii. Second, I evaluated the forecast accuracy by implementing the following four accuracy measures: the mean absolute error (MAE), mean absolute percentage error (MAPE), mean absolute scaled error (MASE), and root mean squared error (RMSE). 10 These accuracy measures allow the validation of the models selected through the auto.arima() algorithm. iii. Third, I examined the autocorrelation function (ACF) and the partial autocorrelation function (PACF) of the residuals from ARIMA models, to detect the presence of significant spikes. iv.
Finally, I verified the fundamental statistical assumptions, homoscedasticity, and serial correlation. 11 The first was controlled by using Engle's Lagrange Multiplier (1982) test for autoregressive conditional heteroskedasticity (ARCH), and the former was detected by implementing the Ljung- Box's (1978) test for serial correlation. The basic ARIMA forecasting estimated equation is the following [1]:
Where Δ means the degree of nonseasonal differences, ŷ is the time series that will be predicted at time t, p is the lag order of AR, λ is the coefficient for each parameter p, q is the lag order of MA, ψ is the coefficient of each parameter q, and means the residuals of errors in time t. 10 The mathematical formulas for each measure of forecast accuracy are provided in Table A1 in Appendix A. 11 The normality assumption of the residuals is not generally necessary (Hyndman and Athanasopoulos, 2018, 3 .3). 
The minimization of AICc and forecast error measurements suggests the following parameters (Table 2) : Italy (4, 2, 4), Russia (1, 2, 1), and USA (6, 2, 3). 12 In the case of Italy, ARIMA (4, 2, 4) is the absolute best according to all criteria. For Russia, MAE, MASE, and RMSE suggest ARIMA (1, 2, 2); however, the AICc and MAPE are definitely lower for the ARIMA (1, 2, 1) . For the USA, ARIMA (6, 2, 3) performs far better in all measures of forecast accuracy, but the AICc, which is slightly higher than the ARIMA (2, 2, 4). I chose the former, also based on the residuals auto-correlation structure.
Moreover, MAPE and MASE may give other useful information on the fitting accuracy.
In particular, MAPEwhich is the most widely used error measure (Goodwin and Lawton, 1999; Ren and Glasure, 2009; Moreno et al., 2013; Kim H. and Kim S., 2016 )indicates an overall forecast accuracy (100-MAPE) of 86.96% for Italy, 88.61% for Russia, and 90.41% for USA. 13 According to Lewis (1982, p. 40) , since the values of MAPE are close to 10, the forecasting is definitely good, especially in the case of the USA that shows a value lower than this "limit" (9.59). 14 MASE, which was more recently proposed by Hyndman and Koehler (2008) , seems to be more versatile and reliable. Since it is always lower than 1 for all three countries, it indicates that the computed forecast performs definitely better than in-sample one-step forecasts from the naïve method. The ACF and PACF correlogram plots of the residuals for Italy ( Figure C1 ) and the USA ( Figure C3 ) show a white noise process; in fact, any of the associated spikes go beyond the 95% confidence intervals. Meanwhile, the ACF and PACF correlogram plots for Russia exhibit a significant spike at lag 10, which reaches the 95% confidence limits ( Figure C2 ). Even if it should not be a particular matter of concern, I further examine the ARIMA model process by applying the Ljung-Box's test for autocorrelation and Engle's LM test for the Arch effect. The Ljung-Box's test shows that the null hypothesis of serial independence of the residuals can be accepted at each distinct lag for all three countries (Table B2 ). In the same way, Engle's LM test for the ARCH effect shows that the null hypothesis of no autoregressive conditional heteroscedasticity can always be accepted (Table B3) . Therefore, the three ARIMA models show a very good fit. In Figures 4, 5, and 6 , I present the forecast results for 30 days in each country. From April 15 to May 14, Italy's coronavirus infectionthat is the "benchmark" of the analysis showed a clear sinusoidal and declining trend (Figure 4) . The new daily COVID-19 12 The estimated parameters are reported in Table B1 in Appendix B. 13 It is necessary to stress that the MAPE has also been affected by some criticalities. In fact, it put greater penalty on negative errors than on positive errors and may be very problematic with time points close to zero (Hyndman and Athanasopoulos, 2018, 3.4) . However, if the first is mostly out of control, the second aspect is not worry about because all time points considered in the analysis are far from zero. 14 Specifically, according to the interpretation of Lewis, the models for Italy and Russia are good forecasting, and that for the USA can be considered a highly accurate forecasting.
cases should drop near zero in mid-June 2020, with an estimated epidemic final size of about 246,627 people. The outcomes for Russia ( Figure 5 ) and the USA ( Figure 6 ) are also characterized by a declining trend, indicating that they have probably already reached the peak of the COVID-19 outbreak. The major difference between the Italian and Russia/USA cases is that the slowdown in new COVID-19 cases seems to require more time for the latter. Specifically, Russia will reach zero local COVID-19 cases in late August 2020, while in the USA, the new COVID-19 cases will drop to zero between the end of September 2020 and the start of October 2020. This may confirm the importance of the rapid implementation of strict lockdown measures to contain the spread of the epidemic. In fact, Italy has been under one of the world's strictest lockdown for about 2 months, while Russia and USA have applied easier lockdown restrictions To control the reliability of my approach, in Tables 3 and 4 , I compared the predicted and latest actual daily cases in Italy and the USA, respectively. The results show that in Italy, the overall predicted new daily cases are overestimated, and this deviation grows over time. In the first 10 days of the forecast, the percentage deviation was about 7.27, and it increased until the 30 th day of forecast, when it stabilized around 22%. On May 26, i.e. after 42 days, the overall deviation between the predicted and actual values was equal to 15,326 cases, with a change of 22.52 percentage point. After 10 days, the MAPE was equal to 13.75% and it reached 33.32% at the end of the forecast (after 42 days). Therefore, for Italy, ARIMA performs better at predicting the final size rather than the daily cases. Regarding the USA, Table 4 shows that the overall predicted cases are a little bit underestimated. In particular, after 5 days, the percentage deviation was about -5.86% and it decreased distinctively after 19 days, reaching 4.11%. After 19 days, the overall deviation between predicted and actual values was equal to 17,128 cases. MAPE, which was equal to 5.73% after 5 days, slightly increased to 5.91% on June 4. Therefore, ARIMA models seem to perform better for the USA than Italy.
About Russia, Table 4 shows that, after 5 days, the percentage deviation was about -2.82%, and it increased to -5.82% after 11 days. MAPE was equal 3.51% after 5 days, and 6.06% after 11 days. These values can be considered reasonable errors, all things considered. 15 In fact, if we look at Figures 7 and 8 , we see that in Italy and the USA, the predicted and actual daily cases follow a very similar and comparable trend. About Russia, the trends seem to differ a little after 5 days (Figure 9 ). This could indicate that the epidemic is still in the stabilisation phase. Even if the analytical approach seems to work well enough, these results are a rough guide only and should be additionally verified and updated. 
To the best of my knowledge, this is the first study to implement an ARIMA model to predict the incidence of COVID-19 in the short/mid-term in Italy, Russia, and the USA. ARIMA models can be viewed as an easy and immediate tool to program health monitoring systems and to better allocate the available resources. In particular, it can help public health decision makers to plan the number of beds needed for moderate and critical care, and to better allocate and manage medical staff, healthcare devices, and healthcare facilities.
Results suggest that Italy, the USA, and Russia reached the peak of COVID-19 infections in mid-April, mid-May, and late May, respectively. In particular, the USA and Russia will require a considerable length of time to drop near zero daily cases, if compared to Italy. This could be because Italy imposed stricter nationwide lockdown measures such as severe traffic and travel restrictions, bans on gatherings, and closure of commercial activitiesto mitigate the spread of the outbreak. Moreover, the comparison between the predicted and more recent actual values showed that forecasts are reliable enough, even in the mid-term, especially for Italy and the USA. Therefore, ARIMA may be considered a good model for short to medium-term forecasting, but the results should be interpreted with thriftiness. Finally, further useful and more precise forecasting may be provided by continuously updating these data, adding interventions and other real aspects, and applying the model to other countries and/or regions.
APPENDIX A Figure C1 . ACF and PACF plots of residuals for Italy. Figure C2 . ACF and PACF plots of residuals for Russia. Figure C3 . ACF and PACF plots of residuals for the USA.