6. Forecast the Chinese GDP from the
global_economy
data set using an ETS model. Experiment with
the various options in the ETS()
function to see how much
the forecasts change with damped trend, or with a Box-Cox
transformation. Try to develop an intuition of what each is doing to the
forecasts.
[Hint: use a relatively large value of h when forecasting, so you can clearly see the differences between the various options when plotting the forecasts.]
china <- global_economy |>
filter(Country == "China")
china |> autoplot(GDP)
It clearly needs a relatively strong transformation due to the increasing variance.
china |> autoplot(box_cox(GDP, 0.2))
china |> features(GDP, guerrero)
Country <fctr> | lambda_guerrero <dbl> | |||
---|---|---|---|---|
China | -0.03446284 |
Making λ=0.2 looks ok.
The Guerrero method suggests an even stronger transformation. Let’s also try a log.
fit <- china |>
model(
ets = ETS(GDP),
ets_damped = ETS(GDP ~ trend("Ad")),
ets_bc = ETS(box_cox(GDP, 0.2)),
ets_log = ETS(log(GDP))
)
fit
Country <fctr> | ets <S3: lst_mdl> | ets_damped <S3: lst_mdl> | ets_bc <S3: lst_mdl> | ets_log <S3: lst_mdl> |
---|---|---|---|---|
China | <ETS(M,A,N)> | <ETS(M,Ad,N)> | <ETS(A,A,N)> | <ETS(A,A,N)> |
augment(fit)
Country <fctr> | .model <chr> | Year <dbl> | GDP <dbl> | .fitted <dbl> | .resid <dbl> | .innov <dbl> |
---|---|---|---|---|---|---|
China | ets | 1960 | 5.971647e+10 | 4.900169e+10 | 1.071478e+10 | 2.186614e-01 |
China | ets | 1961 | 5.005687e+10 | 6.634664e+10 | -1.628977e+10 | -2.455252e-01 |
China | ets | 1962 | 4.720936e+10 | 5.160737e+10 | -4.398009e+09 | -8.522057e-02 |
China | ets | 1963 | 5.070680e+10 | 4.738649e+10 | 3.320305e+09 | 7.006860e-02 |
China | ets | 1964 | 5.970834e+10 | 5.191909e+10 | 7.789252e+09 | 1.500267e-01 |
China | ets | 1965 | 7.043627e+10 | 6.335042e+10 | 7.085845e+09 | 1.118516e-01 |
China | ets | 1966 | 7.672029e+10 | 7.628919e+10 | 4.310994e+08 | 5.650858e-03 |
China | ets | 1967 | 7.288163e+10 | 8.270838e+10 | -9.826744e+09 | -1.188120e-01 |
China | ets | 1968 | 7.084654e+10 | 7.580482e+10 | -4.958286e+09 | -6.540858e-02 |
China | ets | 1969 | 7.970591e+10 | 7.222226e+10 | 7.483647e+09 | 1.036197e-01 |
fit |>
forecast(h = "20 years") |>
autoplot(china, level = NULL)
The transformations have a big effect, with small lambda values creating big increases in the forecasts. The damping has relatively a small effect.
7.Find an ETS model for the Gas data from
aus_production
and forecast the next few years. Why is
multiplicative seasonality necessary here? Experiment with making the
trend damped. Does it improve the forecasts?
aus_production |> autoplot(Gas)
There is a huge increase in variance as the series increases in level which makes it necessary to use multiplicative seasonality.
fit <- aus_production |>
model(
hw = ETS(Gas ~ error("M") + trend("A") + season("M")),
hwdamped = ETS(Gas ~ error("M") + trend("Ad") + season("M")),
)
fit |> glance()
.model <chr> | sigma2 <dbl> | log_lik <dbl> | AIC <dbl> | AICc <dbl> | BIC <dbl> | MSE <dbl> | AMSE <dbl> | |
---|---|---|---|---|---|---|---|---|
hw | 0.003244894 | -831.4643 | 1680.929 | 1681.794 | 1711.389 | 21.11507 | 32.16296 | |
hwdamped | 0.003285280 | -832.0139 | 1684.028 | 1685.091 | 1717.873 | 21.08500 | 32.04319 |
The non-damped model seems to be doing slightly better here, probably because the trend is very strong over most of the historical data.
fit |>
select(hw) |>
gg_tsresiduals()
fit |> tidy()
.model <chr> | term <chr> | estimate <dbl> | ||
---|---|---|---|---|
hw | alpha | 0.65285449 | ||
hw | beta | 0.14416745 | ||
hw | gamma | 0.09784922 | ||
hw | l[0] | 5.94559214 | ||
hw | b[0] | 0.07062881 | ||
hw | s[0] | 0.93092356 | ||
hw | s[-1] | 1.17788274 | ||
hw | s[-2] | 1.07485100 | ||
hw | s[-3] | 0.81634270 | ||
hwdamped | alpha | 0.64890445 |
# A tibble: 19 × 3
fit |>
augment() |>
filter(.model == "hw") |>
features(.innov, ljung_box, lag = 24)
.model <chr> | lb_stat <dbl> | lb_pvalue <dbl> | ||
---|---|---|---|---|
hw | 57.1136 | 0.0001613395 |
There is still some small correlations left in the residuals, showing the model has not fully captured the available information. There also appears to be some heteroskedasticity in the residuals with larger variance in the first half the series.
fit |>
forecast(h = 36) |>
filter(.model == "hw") |>
autoplot(aus_production)
While the point forecasts look ok, the intervals are excessively wide.
8.Recall your retail time series data (from Exercise 7 in Section 2.10).
a.Why is multiplicative seasonality necessary for this series?
set.seed(12345678)
myseries <- aus_retail |>
filter(`Series ID` == sample(aus_retail$`Series ID`, 1))
myseries |> autoplot(Turnover)
The variation in the seasonal pattern increases as the level of the series rises. (This may not be true for every series, but is true for almost all of them.)
b. Apply Holt-Winters’ multiplicative method to the data. Experiment with making the trend damped.
fit <- myseries |>
model(
hw = ETS(Turnover ~ error("M") + trend("A") + season("M")),
hwdamped = ETS(Turnover ~ error("M") + trend("Ad") + season("M"))
)
fc <- fit |> forecast(h = 36)
fc |> autoplot(myseries)
c.Compare the RMSE of the one-step forecasts from the two methods. Which do you prefer?
accuracy(fit)
State <chr> | Industry <chr> | .model <chr> | |
---|---|---|---|
Northern Territory | Clothing, footwear and personal accessory retailing | hw | |
Northern Territory | Clothing, footwear and personal accessory retailing | hwdamped |
The non-damped method is doing slightly better (on RMSE), but the damped method is doing better on most other scores. I’d be inclined to use the damped method here as the trend at the end of the series seems to be flattening.
d.Check that the residuals from the best method look like white noise.
fit |>
select("hwdamped") |>
gg_tsresiduals()
There are significant spikes at lags 8 and 18 in the ACF, but they are relatively small and probably of no consequence.
augment(fit) |>
filter(.model == "hwdamped") |>
features(.innov, ljung_box, lag = 36)
State <chr> | Industry <chr> | .model <chr> | |
---|---|---|---|
Northern Territory | Clothing, footwear and personal accessory retailing | hwdamped |
Overall, the autocorrelation in the residuals is not significant.
e.Now find the test set RMSE, while training the model to the end of 2010. Can you beat the seasonal naïve approach from Exercise 7 in Section 5.11?
myseries |>
filter(year(Month) < 2011) |>
model(
snaive = SNAIVE(Turnover),
hw = ETS(Turnover ~ error("M") + trend("A") + season("M")),
hwdamped = ETS(Turnover ~ error("M") + trend("Ad") + season("M"))
) |>
forecast(h = "7 years") |>
accuracy(myseries)
.model <chr> | State <chr> | Industry <chr> | |
---|---|---|---|
hw | Northern Territory | Clothing, footwear and personal accessory retailing | |
hwdamped | Northern Territory | Clothing, footwear and personal accessory retailing | |
snaive | Northern Territory | Clothing, footwear and personal accessory retailing |
The SNAIVE model is doing much better than the HW model for this data set.