Discussion 2

Pick a time series of interest to you. Build an auto.arima model as well as an ETS model. Which performed better? Now hold out 6 months of data for a test set and try to forecast using the ETS and the auto.arima. Which performs better on the hold-out set? # Collect data from FRED

Seasonally un-adjusted quarterly PCE

library(tidyverse)
library(fpp3)
library(fredr)
library(ggthemes)
nonfarm_payroll <- fredr(
  series_id = "PAYNSA",
  observation_start = as.Date("1950-01-01"),
  observation_end = as.Date("2023-10-01"),
  frequency = "m",
  units = "lin"
) |> 
  select(date, value) |> 
  mutate(date = tsibble::yearmonth(date)) |> 
  as_tsibble()

Plot

Change points -

• Decline during GFC
• Increase post-COVID pandemic recession

nonfarm_payroll |> 
  ggplot(aes(date, value, color = "temp")) +
  geom_line(linewidth = 1.25, show.legend = FALSE) +
  scale_x_yearquarter(
    breaks = scales::breaks_width("10 years")
  ) +
  scale_y_continuous(
    breaks = scales::breaks_width(25000),
    labels = scales::label_comma()
  ) +
  scale_color_economist() +
  labs(
    title = "Total nonfarm payrolls",
    subtitle = "Not seasonally adjusted, in thousands"
  ) +
  theme_minimal() +
  theme(axis.title = element_blank())

STL decomposition

decomp <- nonfarm_payroll |> 
  model(stl = STL(value))
decomp |> 
  components() |> 
  autoplot()

nonfarm_payroll |> 
  autoplot(value) +
  autolayer(components(decomp), season_adjust, color = "darkred") +
  scale_y_continuous(
    breaks = scales::breaks_width(25000),
    labels = scales::label_comma()
  ) +
  theme_minimal() +
  theme(axis.title = element_blank())

Train/test splits

train <- nonfarm_payroll |> 
  filter_index("1950 Q1" ~ "2000 Q1")

test <- nonfarm_payroll |> 
  filter_index("2000 Q2" ~ .)

fit <- train |> 
  model(
    multiplicative = ETS(value ~ error("M") + trend("M") + season("M")),
    auto_arima = ARIMA(value, stepwise = FALSE, approximation = FALSE)
  )

fit |> report()

# A tibble: 2 × 11
  .model   sigma2 log_lik    AIC   AICc    BIC     MSE    AMSE      MAE ar_roots
  <chr>     <dbl>   <dbl>  <dbl>  <dbl>  <dbl>   <dbl>   <dbl>    <dbl> <list>  
1 multip… 7.04e-5  -5803. 11641. 11642. 11715. 226780. 393133.  0.00514 <NULL>  
2 auto_a… 4.07e+4  -3956.  7922.  7922.  7943.     NA      NA  NA       <cpl>   
# ℹ 1 more variable: ma_roots <list>

Forecast on full data set with additive model

auto_arima <- train |> 
  model(auto_arima = ARIMA(value, stepwise = FALSE, approximation = FALSE))

auto_arima |> 
  forecast(h = 36) |> 
  autoplot(
    linewidth = 1.25
  ) +
  autolayer(
    nonfarm_payroll |> 
      filter_index("1996 Q1" ~ "2020 Q1"),
    level = NULL,
    linewidth = 1.25,
    show_gap = FALSE
  ) +
  labs(
    title = "Forecasted nonfarm payrolls 12 quarters ahead using auto ARIMA"
  ) +
  theme_minimal()

Collect fitted values

fitted_arima <- fitted(auto_arima, 36)

Check against test set

train |> 
  autoplot(value) +
  autolayer(fitted_arima, .fitted, color = "darkred", linewidth = 1.25) +
  theme_minimal() +
  theme(legend.text = element_text()) +
  labs(
    title = "Fitted forecast vs. actual nonfarm payroll"
  )