We used the Box-Jenkins approach to ARIMA (p, d, q) modeling of time series [20] . This model-building process is designed to take advantage of associations in the sequentially lagged relationships that usually exist in periodically collected data [21] . The following were the parameters selected when fitting the ARIMA model: p, the order of autoregression; d, the degree of difference; q, the order of moving average. The annual data used in this study did not show seasonal pattern, so the series was differenced at the nonseasonal level to induce stationarity. Autocorrelation function (ACF) graph and Partial autocorrelation function (PACF) graph were used to identify the order of moving average (MA) and autoregressive (AR) terms included in the ARIMA model. Estimates of the model's parameters were obtained by the conditional least squares method. Diagnostic checking including residual analysis and the Akaike Information Criterion (AIC) was used to compare the goodness-of-fit among ARIMA models. The Ljung-Box test was used to measure the ACF of the residuals. In addition, we used the mean absolute percentage error (MAPE) and fitting effect diagram to assess forecast accuracy. 