DATA 624 — Assignment 5 (Forecasting, FPP3 Chapter 8)
Cailinn
2026-10-05
Exercises 8.1, 8.5, 8.6, 8.7, 8.8 and 8.9 from Hyndman & Athanasopoulos, Forecasting: Principles and Practice (3rd ed.), Chapter 8 (Exponential smoothing).
Exercise 8.1 — Pigs slaughtered in Victoria
(a) Simple exponential smoothing is
ETS(A,N,N):
pigs <- aus_livestock |>
filter(Animal == "Pigs", State == "Victoria")
fit <- pigs |> model(ETS(Count ~ error("A") + trend("N") + season("N")))
report(fit)## 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.07fc <- fit |> forecast(h = 4)
fc## # A fable: 4 x 6 [1M]
## # Key: Animal, State, .model [1]
## Animal State .model Month
## <fct> <fct> <chr> <mth>
## 1 Pigs Victoria "ETS(Count ~ error(\"A\") + trend(\"N\") + season(\"… 2019 Jan
## 2 Pigs Victoria "ETS(Count ~ error(\"A\") + trend(\"N\") + season(\"… 2019 Feb
## 3 Pigs Victoria "ETS(Count ~ error(\"A\") + trend(\"N\") + season(\"… 2019 Mar
## 4 Pigs Victoria "ETS(Count ~ error(\"A\") + trend(\"N\") + season(\"… 2019 Apr
## # ℹ 2 more variables: Count <dist>, .mean <dbl>The optimal smoothing parameter is alpha ~ 0.32 and the initial level
l0 ~ 100,647 (see report() above). The forecasts are flat
at ~95,187 for all four months, as expected for SES with no trend.
(b) Manual 95% prediction interval for the first forecast:
s <- sd(residuals(fit)$.resid, na.rm = TRUE)
yhat1 <- fc$.mean[1]
c(lower = yhat1 - 1.96 * s, upper = yhat1 + 1.96 * s)## lower upper
## 76871.01 113502.10hilo(fc, 95)[1, ]## # A tsibble: 1 x 7 [1M]
## # Key: Animal, State, .model [1]
## Animal State .model Month
## <fct> <fct> <chr> <mth>
## 1 Pigs Victoria "ETS(Count ~ error(\"A\") + trend(\"N\") + season(\"… 2019 Jan
## # ℹ 3 more variables: Count <dist>, .mean <dbl>, `95%` <hilo>The manually computed interval (about 76,871 to 113,502) is very close to the interval R computes, though not identical — R’s interval accounts for the estimation uncertainty in alpha and l0, not just residual variance, so it is very slightly different (and based on a t/normal approximation that can differ in its exact bounds).
Exercise 8.5 — Australian exports
aus_exports <- global_economy |> filter(Country == "Australia")
aus_exports |>
autoplot(Exports) +
labs(title = "Australian exports (% of GDP)")
(a) Exports as a share of GDP trend gently upward over the full series, with a dip in the 1980s-90s recessions and a clear rise after 2000 (the mining boom); there is no seasonality since the data are annual.
(b) ETS(A,N,N) forecast:
fit_ann <- aus_exports |> model(ANN = ETS(Exports ~ error("A") + trend("N") + season("N")))
fit_ann |> forecast(h = 5) |> autoplot(aus_exports) +
labs(title = "Australian exports: ETS(A,N,N) forecast")
(c) Training RMSE:
accuracy(fit_ann) |> select(.model, RMSE)## # A tibble: 1 × 2
## .model RMSE
## <chr> <dbl>
## 1 ANN 1.15(d) Compare to ETS(A,A,N):
fit_cmp <- aus_exports |> model(
ANN = ETS(Exports ~ error("A") + trend("N") + season("N")),
AAN = ETS(Exports ~ error("A") + trend("A") + season("N"))
)
accuracy(fit_cmp) |> select(.model, RMSE, MASE)## # A tibble: 2 × 3
## .model RMSE MASE
## <chr> <dbl> <dbl>
## 1 ANN 1.15 0.928
## 2 AAN 1.12 0.907The trended AAN model has a slightly lower training RMSE
(1.12 vs 1.15), at the cost of one extra parameter (the trend slope,
b). The improvement is small, since the data has only a
weak long-run trend relative to its variability.
(e) Compare forecasts:
fit_cmp |> forecast(h = 5) |> autoplot(aus_exports, level = NULL) +
labs(title = "ANN vs AAN forecasts")
ANN forecasts a flat line, while AAN
extrapolates the recent upward trend. Given the modest trend visible in
the data and its small RMSE improvement, AAN is
mildly preferred, but the difference is small enough that
ANN is a defensible, simpler choice.
(f) Manual 95% prediction intervals (using training RMSE as sigma, normal errors) for the first forecast of each model, compared to R’s intervals:
rmse_vals <- accuracy(fit_cmp)$RMSE
fc_cmp <- fit_cmp |> forecast(h = 5)
yhat_ann <- fc_cmp |> filter(.model == "ANN") |> pull(.mean) |> head(1)
yhat_aan <- fc_cmp |> filter(.model == "AAN") |> pull(.mean) |> head(1)
tibble(
model = c("ANN", "AAN"),
manual_lower = c(yhat_ann, yhat_aan) - 1.96 * rmse_vals,
manual_upper = c(yhat_ann, yhat_aan) + 1.96 * rmse_vals
)## # A tibble: 2 × 3
## model manual_lower manual_upper
## <chr> <dbl> <dbl>
## 1 ANN 18.4 22.9
## 2 AAN 18.6 23.0fc_cmp |> group_by(.model) |> filter(row_number() == 1) |> hilo(95) |>
as_tibble() |> select(.model, `95%`)## # A tibble: 2 × 2
## # Groups: .model [2]
## .model `95%`
## <chr> <hilo>
## 1 ANN [18.31970, 22.89462]95
## 2 AAN [18.57028, 23.10700]95The manual intervals are close to R’s but not identical, for the same reason as 8.1(b): R’s intervals incorporate parameter estimation uncertainty, not just the training residual spread.
Exercise 8.6 — Chinese GDP
china <- global_economy |> filter(Country == "China")
china_fit <- china |> model(
ets = ETS(GDP),
damped = ETS(GDP ~ trend("Ad")),
bc_ets = ETS(box_cox(GDP, 0.2)),
bc_damped = ETS(box_cox(GDP, 0.2) ~ trend("Ad"))
)
china_fit## # A mable: 1 x 5
## # Key: Country [1]
## Country ets damped bc_ets bc_damped
## <fct> <model> <model> <model> <model>
## 1 China <ETS(M,A,N)> <ETS(M,Ad,N)> <ETS(A,A,N)> <ETS(A,Ad,N)>china_fit |> forecast(h = 50) |>
autoplot(china, level = NULL) +
labs(title = "Chinese GDP: ETS options", y = "US$")
Letting ETS() choose automatically selects
ETS(M,A,N) — multiplicative errors with a
linear trend — because GDP growth is roughly proportional to its level,
not additive. Over a long horizon this forecasts unbounded
exponential-looking growth. Adding a damped
trend flattens the long-run forecast to a more conservative,
levelling-off growth rate — more plausible 50 years out, since no
economy grows at a constant percentage forever. The Box-Cox
transformation (lambda = 0.2) stabilizes the variance so
ETS() instead chooses additive error/trend on the
transformed scale; back-transformed, this also curbs some of the
explosive growth compared to the untransformed model, though less
aggressively than damping. Combining Box-Cox with
damping produces the most conservative, realistic-looking
long-run forecast of the four.
Exercise 8.7 — Australian gas production
fit_gas <- aus_production |> model(
ets = ETS(Gas),
damped = ETS(Gas ~ trend("Ad"))
)
fit_gas## # A mable: 1 x 2
## ets damped
## <model> <model>
## 1 <ETS(M,A,M)> <ETS(M,Ad,M)>accuracy(fit_gas) |> select(.model, RMSE, MASE)## # A tibble: 2 × 3
## .model RMSE MASE
## <chr> <dbl> <dbl>
## 1 ets 4.60 0.542
## 2 damped 4.59 0.544fit_gas |> select(ets) |> forecast(h = "3 years") |>
autoplot(aus_production) +
labs(title = "Australian gas production: ETS forecast")
Letting ETS() choose automatically picks
multiplicative seasonality (ETS(M,A,M)).
This is necessary because the size of the seasonal swings grows
along with the level of the series (gas production has
increased roughly ten-fold since the 1950s, and the seasonal
peaks/troughs have grown proportionally) — additive seasonality assumes
constant-size seasonal swings, which would badly under-forecast the
seasonal amplitude in recent years.
fit_gas |> forecast(h = "3 years") |> autoplot(aus_production, level = NULL) +
labs(title = "ETS vs damped-trend ETS for gas production")
Adding a damped trend barely changes training RMSE (4.60 vs 4.59) and the forecasts are nearly indistinguishable — the gas trend is close enough to linear over the forecast horizon that damping makes little practical difference here.
Exercise 8.8 — Holt-Winters on retail data
set.seed(12345678)
myseries <- aus_retail |>
filter(`Series ID` == sample(aus_retail$`Series ID`, 1))(a) This series (Northern Territory clothing/footwear/accessory retail) needs multiplicative seasonality for the same reason as Gas: the size of the seasonal (Christmas) spike grows as the overall level of turnover grows, so a fixed additive seasonal amount would not scale with the trend.
(b) Holt-Winters’ multiplicative method, with and without a damped trend:
fit_hw <- myseries |> model(
hw = ETS(Turnover ~ error("M") + trend("A") + season("M")),
hw_damped = ETS(Turnover ~ error("M") + trend("Ad") + season("M"))
)
fit_hw |> forecast(h = "3 years") |>
autoplot(myseries, level = NULL) +
labs(title = "Holt-Winters multiplicative: damped vs not")
(c) Compare one-step training RMSE:
accuracy(fit_hw) |> select(.model, RMSE, MASE)## # A tibble: 2 × 3
## .model RMSE MASE
## <chr> <dbl> <dbl>
## 1 hw 0.613 0.513
## 2 hw_damped 0.616 0.507RMSE is nearly identical (0.613 vs 0.616), with the damped version very slightly better on MASE. The damped trend is mildly preferred — it is a safer long-run assumption with essentially no cost in fit quality.
(d) Residual diagnostics for the preferred (damped) model:
fit_hw |> select(hw_damped) |> gg_tsresiduals()
The residuals are centred near zero with only minor ACF spikes, close enough to white noise for a reasonable, though not perfect, model.
(e) Test-set RMSE (train through 2010) vs seasonal naive from Exercise 5.7:
train <- myseries |> filter(year(Month) < 2011)
fit_test <- train |> model(
hw_damped = ETS(Turnover ~ error("M") + trend("Ad") + season("M")),
snaive = SNAIVE(Turnover)
)
fc_test <- fit_test |> forecast(new_data = anti_join(myseries, train))
fc_test |> accuracy(myseries) |> select(.model, RMSE, MASE)## # A tibble: 2 × 3
## .model RMSE MASE
## <chr> <dbl> <dbl>
## 1 hw_damped 1.15 0.960
## 2 snaive 1.55 1.36Holt-Winters damped (RMSE 1.15, MASE 0.96) is substantially more accurate than the seasonal naive benchmark (RMSE 1.55, MASE 1.36) from Exercise 5.7 — as expected, since ETS captures both the trend and the growing seasonal pattern that seasonal naive ignores.
Exercise 8.9 — STL + ETS vs Holt-Winters
fit_stl <- train |> model(
stlf = decomposition_model(
STL(box_cox(Turnover, 0.1) ~ season(window = "periodic")),
ETS(season_adjust)
)
)
fc_stl <- fit_stl |> forecast(new_data = anti_join(myseries, train))
fc_stl |> accuracy(myseries) |> select(.model, RMSE, MASE)## # A tibble: 1 × 3
## .model RMSE MASE
## <chr> <dbl> <dbl>
## 1 stlf 1.14 0.968The STL + ETS approach (RMSE 1.14, MASE 0.97) performs about as well as the damped Holt-Winters model from Exercise 8.8 (RMSE 1.15, MASE 0.96) — a very slightly lower RMSE but a very slightly higher MASE, so the two methods are essentially tied on this series. Both comfortably beat the seasonal naive benchmark.
Session info
sessionInfo()## R version 4.6.1 (2026-06-24)
## Platform: aarch64-apple-darwin25.4.0
## Running under: macOS Tahoe 26.7
##
## Matrix products: default
## BLAS: /opt/homebrew/Cellar/openblas/0.3.34/lib/libopenblasp-r0.3.34.dylib
## LAPACK: /opt/homebrew/Cellar/r/4.6.1/lib/R/lib/libRlapack.dylib; LAPACK version 3.12.1
##
## locale:
## [1] C.UTF-8/C.UTF-8/C.UTF-8/C/C.UTF-8/C.UTF-8
##
## time zone: America/New_York
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] fable_0.5.0 feasts_0.5.0 fabletools_0.8.0 ggtime_1.0.0
## [5] tsibbledata_0.4.1 tsibble_1.2.0 ggplot2_4.0.3 lubridate_1.9.5
## [9] tidyr_1.3.2 dplyr_1.2.1 tibble_3.3.1 fpp3_1.0.3
##
## loaded via a namespace (and not attached):
## [1] ggdist_3.3.3 utf8_1.2.6 rappdirs_0.3.4
## [4] sass_0.4.10 generics_0.1.4 anytime_0.3.13
## [7] digest_0.6.39 magrittr_2.0.5 evaluate_1.0.5
## [10] grid_4.6.1 timechange_0.4.0 RColorBrewer_1.1-3
## [13] mixtime_0.3.0 fastmap_1.2.0 jsonlite_2.0.0
## [16] purrr_1.2.2 scales_1.4.0 numDeriv_2016.8-1.1
## [19] jquerylib_0.1.4 cli_3.6.6 rlang_1.3.0
## [22] crayon_1.5.3 vecvec_1.3.0 withr_3.0.3
## [25] cachem_1.1.0 yaml_2.3.12 tools_4.6.1
## [28] tzdb_0.5.0 vctrs_0.7.3 R6_2.6.1
## [31] lifecycle_1.0.5 pkgconfig_2.0.3 progressr_1.0.0
## [34] pillar_1.11.1 bslib_0.12.0 gtable_0.3.6
## [37] glue_1.8.1 Rcpp_1.1.2 xfun_0.60
## [40] tidyselect_1.2.1 rstudioapi_0.19.0 knitr_1.52
## [43] farver_2.1.2 htmltools_0.5.9 labeling_0.4.3
## [46] rmarkdown_2.32 compiler_4.6.1 S7_0.2.2
## [49] distributional_0.9.0