# How to improve forecast

I have the number of hospitalizations due to a particular disease and I'm trying to forecast the number of hospitalizations for the next 3 months.

I've been reading about ARIMA and time series but I'm mainly learning on my own and I have some questions and I would also like to check what I've done is correct. First I remove the known last 4 values to compare the forecasted values with the observed values.

The plot of observed values:

I checked the ACF and PACF plots and since the data wasn't stationary I transformed the data with diff.

Then I proceeded to the model. I fitted an ARIMA with seasonal pattern and differencing equal to 1.

> fit_sarima <- auto.arima(ts, d = 1, seasonal = T, stepwise = F, approximation = F)
> summary(fit_sarima)

Series: ts
ARIMA(0,1,0)(1,1,0)[12]

Coefficients:
sar1
-0.5835
s.e.   0.1413

sigma^2 estimated as 31949:  log likelihood=-219.95
AIC=443.9   AICc=444.3   BIC=446.9

Training set error measures:
ME     RMSE      MAE       MPE     MAPE
Training set -24.67867 149.0825 97.59303 -2.522906 8.812427
MASE        ACF1
Training set 0.5977595 -0.05551552

> checkresiduals(fit_sarima)

Ljung-Box test

data:  Residuals from ARIMA(0,1,0)(1,1,0)[12]
Q* = 8.9865, df = 8, p-value = 0.3434

Model df: 1.   Total lags used: 9


The model seems to fit the data. I forecasted the next 4 months and compared them with the observed values.

> forecast(fit_sarima, h=4, level=95)

Point Forecast    Lo 95    Hi 95
Nov 2019       823.9865 473.6554 1174.318
Dec 2019      1092.5124 597.0694 1587.955
Jan 2020      1578.2382 971.4469 2185.029
Feb 2020      1478.9154 778.2532 2179.578

> obs
time      hosp
1: nov 2019    863
2: dez 2019   1029
3: jan 2020   1186
4: fev 2020   1038


However, when I plot the fitted values with the observed values I get some weird values for the forecast. Maybe it's expected, I'm not expecting a perfect fit, but I was expecting a somewhat more smoothed curve. Should I have more data available to forecast? (I don't have more but I'm assuming this is a limitation for forecast).

I also used a ets model to fit the data and although the model seemed less appropriate since the Ljung-Box is significant and the MPE and AIC are higher the curve is smoother compared with the ARIMA, and the confidence intervals smaller.

> fit_exp <- ets(ts)
> summary(fit_exp)

ETS(M,N,M)

Call:
ets(y = ts)

Smoothing parameters:
alpha = 0.5077
gamma = 3e-04

Initial states:
l = 1156.097
s = 1.1451 0.8585 0.832 0.7256 0.8352 0.8213
0.848 0.9036 1.0145 1.1833 1.3147 1.5183

sigma:  0.1088

AIC     AICc      BIC
629.7517 645.7517 657.1813

Training set error measures:
ME     RMSE      MAE       MPE     MAPE
Training set -8.479541 127.0509 84.58677 -1.216796 7.032338
MASE      ACF1
Training set 0.5180958 0.1998706

> checkresiduals(fit_exp)

Ljung-Box test

data:  Residuals from ETS(M,N,M)
Q* = 27.678, df = 3, p-value = 4.244e-06

Model df: 14.   Total lags used: 17

> forecast(fit_exp, h=4, level=95)

Point Forecast     Lo 95    Hi 95
Nov 2019       874.5813  688.0925 1061.070
Dec 2019      1166.6136  887.2940 1445.933
Jan 2020      1546.8097 1139.8502 1953.769
Feb 2020      1339.2916  957.7705 1720.813

> obs
time      hosp
1: nov 2019    863
2: dez 2019   1029
3: jan 2020   1186
4: fev 2020   1038


Again, looking at the fitted values the ARIMA (red) looks a bit strange. Is there any way to improve this model? Or should I go with the exponential smoothing model instead (blue)?

**I ended up comparing the accuracy of both models considering the 3 extra months the test dataset and the remaining the training dataset. The accuracy was better for the exponential smoothing and I ended up selecting that model instead, despite the fit of the ARIMA being better.

• Thank you for your comment. I'm using auto.arima. Regarding the data, I imagined the historical data is too short but I don't have more data available. I only need to forecast 4 or so months. I tested the data with difference equal to 1 and it was okay. Hence using d = 1, but I only used the adf.test. My understanding is that auto.arima also doesn't work for d higher than 2. Did I misunderstood the Ljung box test? No significant is white noise so the model should be fine. No correlated residuals,right? Aug 10 '20 at 21:27
• The exponential smoothing was significant for the Ljung box test, but accuracy and the fitted data looks better with the exponential smoothing. Hence my confusion. ARIMA has these weird peaks in the fitted data (last plot, red line). Also, I'm not sure if d equal to 1 is enough to create stationarity and I thought ets also worked with non-stationary data. Thank you! I will check kpss. Aug 11 '20 at 4:59
• I tried kpss and the data is stationary. Even if I remove d = 1 auto.arima adds the differencing so I guess it's appropriate. I also used your suggestion and did the training with three years and the test with one year and accuracy is much better for arima now. I guess between the two it would be best to use arima. I need to look into combine models. How would you combine the CI? Just average as well, or take the lower and upper? ARIMA gives me wide CIs but it's probably to be expected since the historical data is only 4 years. Aug 11 '20 at 8:34