library(fpp3)

1 Exercise 8.1 (aus_livestock: pigs slaughtered in Victoria)

pigs_vic <- aus_livestock |>
  filter(Animal == "Pigs", State == "Victoria") |>
  select(Month, Count)

head(pigs_vic)
## # A tsibble: 6 x 2 [1M]
##      Month  Count
##      <mth>  <dbl>
## 1 1972 Jul  94200
## 2 1972 Aug 105700
## 3 1972 Sep  96500
## 4 1972 Oct 117100
## 5 1972 Nov 104600
## 6 1972 Dec 100500

1.1 (a) ETS equivalent to simple exponential smoothing; α, ℓ0; 4-month forecasts

fit_ses <- pigs_vic |>
  model(SES = ETS(Count ~ error("A") + trend("N") + season("N")))

report(fit_ses)
## 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.07
fc_ses <- fit_ses |> forecast(h = "4 months")
fc_ses
## # A fable: 4 x 4 [1M]
## # Key:     .model [1]
##   .model    Month
##   <chr>     <mth>
## 1 SES    2019 Jan
## 2 SES    2019 Feb
## 3 SES    2019 Mar
## 4 SES    2019 Apr
## # ℹ 2 more variables: Count <dist>, .mean <dbl>
autoplot(fc_ses, pigs_vic) +
  labs(title = "Pigs slaughtered in Victoria: SES via ETS(A,N,N)", x = "Month", y = "Count")

Write-up (8.1a):
I fit an ETS(A,N,N) model, which is equivalent to simple exponential smoothing. I recorded the estimated alpha and initial level l[0] from report(), then generated forecasts for the next four months.

1.2 (b) Manual 95% PI for first forecast: ŷ ± 1.96 s; compare to R

aug <- augment(fit_ses)
s <- sd(aug$.resid, na.rm = TRUE)

first_mean <- fc_ses |>
  as_tibble() |>
  slice(1) |>
  pull(.mean)

manual_lo <- first_mean - 1.96 * s
manual_hi <- first_mean + 1.96 * s

c(s = s, first_mean = first_mean, manual_lo = manual_lo, manual_hi = manual_hi)
##          s first_mean  manual_lo  manual_hi 
##   9344.666  95186.557  76871.012 113502.102
r_int <- fc_ses |> hilo(95) |> as_tibble() |> slice(1)
r_int |> select(.mean, `95%`)
## # A tibble: 1 × 2
##    .mean                  `95%`
##    <dbl>                 <hilo>
## 1 95187. [76854.79, 113518.3]95

Write-up (8.1b):
I computed the first 95% prediction interval using ŷ ± 1.96·s, where s is the residual standard deviation from the fitted model. The result is very close to the 95% interval produced by hilo(95).

2 Exercise 8.5 (global_economy: Exports for one country)

country_choice <- "United States"

exports_cty <- global_economy |>
  filter(Country == country_choice) |>
  select(Year, Exports)

glimpse(exports_cty)
## Rows: 58
## Columns: 2
## $ Year    <dbl> 1960, 1961, 1962, 1963, 1964, 1965, 1966, 1967, 1968, 1969, 19…
## $ Exports <dbl> 4.969630, 4.899698, 4.809122, 4.870028, 5.103529, 4.988571, 5.…

2.1 (a) Plot exports & discuss features

autoplot(exports_cty, Exports) +
  labs(title = paste("Exports:", country_choice), x = "Year", y = "Exports")

Write-up (8.5a):
The exports series shows long-run growth with periods of faster/slower change. There are also noticeable dips around major economic disruptions.

2.2 (b) ETS(A,N,N) forecast + plot

fit_ann <- exports_cty |>
  model(ANN = ETS(Exports ~ error("A") + trend("N") + season("N")))

fc_ann <- fit_ann |> forecast(h = 10)

autoplot(fc_ann, exports_cty) +
  labs(title = paste("ETS(A,N,N) Forecast:", country_choice), x = "Year", y = "Exports")

2.3 (c) Training RMSE

accuracy(fit_ann) |> select(any_of(c(".model","RMSE","MAE","MAPE")))
## # A tibble: 1 × 4
##   .model  RMSE   MAE  MAPE
##   <chr>  <dbl> <dbl> <dbl>
## 1 ANN    0.627 0.465  5.10

2.4 (d) Compare to ETS(A,A,N); discuss merits

fit_compare <- exports_cty |>
  model(
    ANN = ETS(Exports ~ error("A") + trend("N") + season("N")),
    AAN = ETS(Exports ~ error("A") + trend("A") + season("N"))
  )

accuracy(fit_compare) |> select(any_of(c(".model","RMSE","AICc","AIC")))
## # A tibble: 2 × 2
##   .model  RMSE
##   <chr>  <dbl>
## 1 ANN    0.627
## 2 AAN    0.615
fc_compare <- fit_compare |> forecast(h = 10)

autoplot(fc_compare, exports_cty) +
  labs(title = paste("ANN vs AAN ETS Forecasts:", country_choice), x = "Year", y = "Exports")

Write-up (8.5d):
ANN is the simpler model, while AAN adds a trend component (one more parameter). I compared the training RMSE and the forecast shapes. If the series shows a consistent trend, AAN can fit better; otherwise ANN may be enough.

2.5 (e) Compare forecasts; which is best?

Write-up (8.5e):
I prefer the model that balances fit and simplicity. If RMSE improves only slightly with AAN, I would stick with ANN; if AAN clearly improves accuracy and produces reasonable forecasts, I would choose AAN.

2.6 (f) Manual 95% PI for first forecast for each model using RMSE; compare to R

acc_tbl <- accuracy(fit_compare) |> select(any_of(c(".model","RMSE")))
first_fc <- fc_compare |> as_tibble() |> group_by(.model) |> slice(1) |> ungroup()

manual_pi <- first_fc |>
  left_join(acc_tbl, by = ".model") |>
  mutate(
    manual_lo = .mean - 1.96 * RMSE,
    manual_hi = .mean + 1.96 * RMSE
  ) |>
  select(.model, .mean, RMSE, manual_lo, manual_hi)

manual_pi
## # A tibble: 2 × 5
##   .model .mean  RMSE manual_lo manual_hi
##   <chr>  <dbl> <dbl>     <dbl>     <dbl>
## 1 AAN     12.1 0.615      10.9      13.3
## 2 ANN     11.9 0.627      10.7      13.1
r_pi <- fc_compare |> hilo(95) |> as_tibble() |> group_by(.model) |> slice(1) |> ungroup()
r_pi |> select(.model, .mean, `95%`)
## # A tibble: 2 × 3
##   .model .mean                  `95%`
##   <chr>  <dbl>                 <hilo>
## 1 AAN     12.1 [10.88421, 13.36085]95
## 2 ANN     11.9 [10.65073, 13.13064]95

Write-up (8.5f):
I built an approximate 95% interval using ±1.96·RMSE and compared it to R’s interval. The RMSE-based interval is a rough approximation, so small differences from R’s intervals are expected.

3 Exercise 8.6 (global_economy: Chinese GDP, ETS with options)

china_gdp <- global_economy |>
  filter(Country == "China") |>
  select(Year, GDP)

glimpse(china_gdp)
## Rows: 58
## Columns: 2
## $ Year <dbl> 1960, 1961, 1962, 1963, 1964, 1965, 1966, 1967, 1968, 1969, 1970,…
## $ GDP  <dbl> 59716467625, 50056868958, 47209359006, 50706799903, 59708343489, …
h_years <- 15  # keep it lighter for Posit

fit_china <- china_gdp |>
  model(
    ETS_auto   = ETS(GDP),
    ETS_damped = ETS(GDP ~ error("A") + trend("Ad") + season("N")),
    ETS_boxcox = ETS(box_cox(GDP, lambda = 0) ~ error("A") + trend("A") + season("N"))
  )

fc_china <- fit_china |> forecast(h = h_years)

autoplot(fc_china, china_gdp) +
  labs(title = "China GDP: ETS options comparison", x = "Year", y = "GDP")

Write-up (8.6):
The damped-trend model produces more conservative long-run forecasts (trend growth tapers off). The Box–Cox option changes the scale and can reduce the impact of very large values, which can change the curvature of the forecast path.

4 Exercise 8.7 (aus_production: Gas; multiplicative seasonality and damped trend)

gas <- aus_production |> select(Quarter, Gas)
head(gas)
## # A tsibble: 6 x 2 [1Q]
##   Quarter   Gas
##     <qtr> <dbl>
## 1 1956 Q1     5
## 2 1956 Q2     6
## 3 1956 Q3     7
## 4 1956 Q4     6
## 5 1957 Q1     5
## 6 1957 Q2     7
fit_gas <- gas |>
  model(
    ETS_auto   = ETS(Gas),
    ETS_damped = ETS(Gas ~ trend("Ad"))
  )

accuracy(fit_gas) |> select(any_of(c(".model","RMSE","AICc","AIC")))
## # A tibble: 2 × 2
##   .model      RMSE
##   <chr>      <dbl>
## 1 ETS_auto    4.60
## 2 ETS_damped  4.59
fc_gas <- fit_gas |> forecast(h = "3 years")

autoplot(fc_gas, gas) +
  labs(title = "aus_production Gas: ETS forecasts", x = "Quarter", y = "Gas")

Write-up (8.7):
Multiplicative seasonality is appropriate when seasonal swings increase as the series level increases. A damped trend can help avoid unrealistically fast long-run growth if the trend is not expected to continue indefinitely.

5 Exercise 8.8 (Retail series + Holt-Winters multiplicative)

# Single series to keep memory low
retail <- aus_retail |>
  filter(State == "New South Wales", Industry == "Department stores") |>
  select(Month, Turnover)

glimpse(retail)
## Rows: 441
## Columns: 2
## $ Month    <mth> 1982 Apr, 1982 May, 1982 Jun, 1982 Jul, 1982 Aug, 1982 Sep, 1…
## $ Turnover <dbl> 178.3, 202.8, 176.3, 172.6, 169.6, 181.4, 173.9, 206.6, 346.6…

5.1 (a) Why multiplicative seasonality?

Write-up (8.8a):
Multiplicative seasonality makes sense here because the seasonal changes are larger when the overall level of turnover is higher (the peaks and troughs scale with the level).

5.2 (b–e) Holt-Winters multiplicative; damped trend; RMSE; residual checks; test RMSE (train through 2010)

train <- retail |> filter(year(Month) <= 2010)
test  <- retail |> filter(year(Month) >= 2011)

fit_snaive <- train |> model(SNAIVE = SNAIVE(Turnover))

fit_hw <- train |>
  model(
    HW_mult        = ETS(Turnover ~ error("M") + trend("A")  + season("M")),
    HW_mult_damped = ETS(Turnover ~ error("M") + trend("Ad") + season("M"))
  )

# One-step (training) accuracy
accuracy(fit_hw) |> select(any_of(c(".model","RMSE","AICc","AIC")))
## # A tibble: 2 × 2
##   .model          RMSE
##   <chr>          <dbl>
## 1 HW_mult         19.2
## 2 HW_mult_damped  19.2
# Residual diagnostics (run for the model you prefer)
fit_hw |> select(HW_mult) |> gg_tsresiduals()

# Test RMSE vs seasonal naive
fc_hw <- fit_hw |> forecast(new_data = test)
fc_snaive <- fit_snaive |> forecast(new_data = test)

bind_rows(
  accuracy(fc_hw, test),
  accuracy(fc_snaive, test)
) |> select(.model, RMSE)
## # A tibble: 3 × 2
##   .model          RMSE
##   <chr>          <dbl>
## 1 HW_mult         24.5
## 2 HW_mult_damped  36.1
## 3 SNAIVE          25.5

Write-up (8.8b–e):
I compared Holt–Winters multiplicative with and without a damped trend using training RMSE, then checked residual diagnostics for the preferred model. Finally, I compared test RMSE (post-2010) against the seasonal naïve benchmark to see whether the ETS/HW approach improves accuracy.

6 Exercise 8.9 (STL + Box-Cox + ETS on seasonally adjusted; compare test RMSE)

lambda <- train |> features(Turnover, features = guerrero) |> pull(lambda_guerrero)

fit_stl_ets <- train |>
  model(
    STL_ETS = decomposition_model(
      STL(box_cox(Turnover, lambda) ~ season(window = "periodic")),
      ETS(season_adjust)
    )
  )

fc_stl_ets <- fit_stl_ets |> forecast(new_data = test)

bind_rows(
  accuracy(fc_stl_ets, test),
  accuracy(fc_hw, test)
) |> select(.model, RMSE)
## # A tibble: 3 × 2
##   .model          RMSE
##   <chr>          <dbl>
## 1 STL_ETS         29.6
## 2 HW_mult         24.5
## 3 HW_mult_damped  36.1

Write-up (8.9):
I compared the STL+ETS test RMSE to the best Holt–Winters/ETS model from 8.8. I would choose the method with lower test RMSE unless the difference is tiny, in which case I would prefer the simpler approach.

7 Session Info

sessionInfo()
## R version 4.5.2 (2025-10-31)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.4 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=C.UTF-8       LC_NUMERIC=C           LC_TIME=C.UTF-8       
##  [4] LC_COLLATE=C.UTF-8     LC_MONETARY=C.UTF-8    LC_MESSAGES=C.UTF-8   
##  [7] LC_PAPER=C.UTF-8       LC_NAME=C              LC_ADDRESS=C          
## [10] LC_TELEPHONE=C         LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C   
## 
## time zone: UTC
## tzcode source: system (glibc)
## 
## 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.6.1  ggtime_0.2.0     
##  [5] tsibbledata_0.4.1 tsibble_1.2.0     ggplot2_4.0.2     lubridate_1.9.5  
##  [9] tidyr_1.3.2       dplyr_1.2.0       tibble_3.3.1      fpp3_1.0.3       
## 
## loaded via a namespace (and not attached):
##  [1] ggdist_3.3.3         rappdirs_0.3.4       sass_0.4.10         
##  [4] utf8_1.2.6           generics_0.1.4       anytime_0.3.12      
##  [7] digest_0.6.39        magrittr_2.0.4       evaluate_1.0.5      
## [10] grid_4.5.2           timechange_0.4.0     RColorBrewer_1.1-3  
## [13] fastmap_1.2.0        jsonlite_2.0.0       purrr_1.2.1         
## [16] scales_1.4.0         numDeriv_2016.8-1.1  jquerylib_0.1.4     
## [19] cli_3.6.5            rlang_1.1.7          crayon_1.5.3        
## [22] withr_3.0.2          cachem_1.1.0         yaml_2.3.12         
## [25] tools_4.5.2          vctrs_0.7.1          R6_2.6.1            
## [28] lifecycle_1.0.5      pkgconfig_2.0.3      progressr_0.18.0    
## [31] pillar_1.11.1        bslib_0.10.0         gtable_0.3.6        
## [34] glue_1.8.0           Rcpp_1.1.1           xfun_0.56           
## [37] tidyselect_1.2.1     rstudioapi_0.18.0    knitr_1.51          
## [40] farver_2.1.2         htmltools_0.5.9      rmarkdown_2.30      
## [43] labeling_0.4.3       compiler_4.5.2       S7_0.2.1            
## [46] distributional_0.6.0