Forecasting: Principles and Practice Chapter 8 Exercises
Neil Martin
2024-10-02
Consider the the number of pigs slaughtered in Victoria, available in the aus_livestock dataset.
Use the ETS() function to estimate the equivalent model for simple exponential smoothing. Find the optimal values of α and ℓ0, and generate forecasts for the next four months.
vic_pigs <- aus_livestock |> filter(Animal == "Pigs", State == "Victoria")
vic_pigs |> autoplot()## Plot variable not specified, automatically selected `.vars = Count`
vicpig_mod <- vic_pigs |>
model(ETS(Count ~ error("A") + trend("N") + season("N")))
report(vicpig_mod)## Series: Count
## Model: ETS(A,N,N)
## Smoothing parameters:
## alpha = 0.3221247
##
## Initial states:
## l[0]
## 100646.6
##
## sigma^2: 87480760
##
## AIC AICc BIC
## 13737.10 13737.14 13750.07vicpig_mod## # A mable: 1 x 3
## # Key: Animal, State [1]
## Animal State `ETS(Count ~ error("A") + trend("N") + season("N"))`
## <fct> <fct> <model>
## 1 Pigs Victoria <ETS(A,N,N)>vicpig_fc <- vicpig_mod |>
forecast(h = 4) |>
autoplot(vic_pigs)
vicpig_fc
#### Data set global_economy contains the annual Exports from many
countries. Select one country to analyse.
india <- global_economy |>
filter(Country == "India")
india |> autoplot(Exports)
From this plot we can see that the German exports have been following a significant linear trend from the early 1990s. There is a slight dip around the 2008 financial crash but there seems to be an uptick in exports from Germany after this.
Use an ETS(A,N,N) model to forecast the series, and plot the forecasts.
india_mod <- india |>
model(ETS(Exports ~ error("A") + trend("N") + season("N")))
report(india_mod)## Series: Exports
## Model: ETS(A,N,N)
## Smoothing parameters:
## alpha = 0.9998965
##
## Initial states:
## l[0]
## 4.507901
##
## sigma^2: 1.4709
##
## AIC AICc BIC
## 261.8526 262.2971 268.0339india_mod## # A mable: 1 x 2
## # Key: Country [1]
## Country `ETS(Exports ~ error("A") + trend("N") + season("N"))`
## <fct> <model>
## 1 India <ETS(A,N,N)>india_fc <- india_mod |>
forecast(h = 4) |>
autoplot(india)
india_fc
fitted_india <- fitted(india_mod)
residuals <- india$Exports - fitted_india$.fitted
rmse <- sqrt(mean(residuals^2))
print(rmse)## [1] 1.191728Compare the results to those from an ETS(A,A,N) model. (Remember that the trended model is using one more parameter than the simpler model.) Discuss the merits of the two forecasting methods for this data set.
Compare the forecasts from both methods. Which do you think is best?
india_mod <- india |>
model(ETS(Exports ~ error("A") + trend("A") + season("N")))
report(india_mod)## Series: Exports
## Model: ETS(A,A,N)
## Smoothing parameters:
## alpha = 0.8289068
## beta = 0.2270934
##
## Initial states:
## l[0] b[0]
## 4.598138 -0.1008885
##
## sigma^2: 1.4854
##
## AIC AICc BIC
## 264.3120 265.4659 274.6143india_mod## # A mable: 1 x 2
## # Key: Country [1]
## Country `ETS(Exports ~ error("A") + trend("A") + season("N"))`
## <fct> <model>
## 1 India <ETS(A,A,N)>india_fc <- india_mod |>
forecast(h = 4) |>
autoplot(india)
india_fc
fitted_india <- fitted(india_mod)
residuals <- india$Exports - fitted_india$.fitted
rmse <- sqrt(mean(residuals^2))
print(rmse)## [1] 1.176006The RMSE shows that the AAN model is roughly 0.02 better than the ANN model.
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.
ch_gdp <- global_economy |>
filter(Country == "China") |>
mutate(GDP_Billions = GDP / 1e9)
ch_gdp |> autoplot(GDP_Billions)
lambda <- ch_gdp |>
features(GDP_Billions, features = guerrero) |>
pull(lambda_guerrero)
ch_gdp_mod <- ch_gdp |>
model(`Simple` = ETS(GDP ~ error("A") + trend("N") + season("N")),
`Holt's method` = ETS(GDP ~ error("A") + trend("A") + season("N")),
`Damped Holt's method` = ETS(GDP ~ error("A") + trend("Ad", phi = 0.8) + season("N")),
`Box-Cox` = ETS(box_cox(GDP,lambda) ~ error("A") + trend("A") + season("N")),
`Box-Cox Damped` = ETS(box_cox(GDP,lambda) ~ error("A") + trend("Ad", phi = 0.8) + season("N")),
`Log` = ETS(log(GDP) ~ error("A") + trend("A") + season("N")),
`Log Damped` = ETS(log(GDP) ~ error("A") + trend("Ad", phi = 0.8) + season("N"))
)
ch_gdp_mod## # A mable: 1 x 8
## # Key: Country [1]
## Country Simple `Holt's method` `Damped Holt's method` `Box-Cox`
## <fct> <model> <model> <model> <model>
## 1 China <ETS(A,N,N)> <ETS(A,A,N)> <ETS(A,Ad,N)> <ETS(A,A,N)>
## # ℹ 3 more variables: `Box-Cox Damped` <model>, Log <model>,
## # `Log Damped` <model>ch_gdp_fc <- ch_gdp_mod |>
forecast(h = 20) |>
autoplot(ch_gdp) + labs(title = "Chinese 2040 GDP Forecast")
ch_gdp_fc
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_gas <- aus_production |>
select(Gas)
aus_gas |> autoplot()## Plot variable not specified, automatically selected `.vars = Gas`
aus_gas_mod <- aus_gas|>
model(additive = ETS(Gas ~ error("A") + trend("A") + season("A")),
multiplicative = ETS(Gas ~ error("M") + trend("A") + season("M")),
`damped multiplicative` = ETS(Gas ~ error("M") + trend("Ad", phi = 0.9) + season("M")))
aus_gas_mod## # A mable: 1 x 3
## additive multiplicative `damped multiplicative`
## <model> <model> <model>
## 1 <ETS(A,A,A)> <ETS(M,A,M)> <ETS(M,Ad,M)>report(aus_gas_mod)## Warning in report.mdl_df(aus_gas_mod): Model reporting is only supported for
## individual models, so a glance will be shown. To see the report for a specific
## model, use `select()` and `filter()` to identify a single model.## # A tibble: 3 × 9
## .model sigma2 log_lik AIC AICc BIC MSE AMSE MAE
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 additive 23.6 -927. 1872. 1873. 1903. 22.7 29.7 3.35
## 2 multiplicative 0.00324 -831. 1681. 1682. 1711. 21.1 32.2 0.0413
## 3 damped multiplicative 0.00340 -835. 1688. 1689. 1719. 21.0 32.4 0.0424aus_gas_fc <- aus_gas_mod |>
forecast(h = 10) |>
autoplot(aus_gas) + labs(title = "Australian Gas Production")
aus_gas_fc