DATA624 - Homework 5
Glen Dale Davis
2023-10-02
library(fpp3)
library(RColorBrewer)
library(knitr)
library(pracma)palette <- brewer.pal(n = 11, name = "BrBG")
#function taken from here: https://bookdown.org/yihui/rmarkdown-cookbook/font-color.html
colorize <- function(x, color){
if (knitr::is_latex_output()) {
sprintf("\\textcolor{%s}{%s}", color, x)
}
else if (knitr::is_html_output()){
sprintf("<span style='color: %s;'>%s</span>", color, x)
}
else x
}Exercise 8.1:
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 \(\alpha\) and \(\ell_0\), and generate forecasts for the next four months.
data(aus_livestock)
pigs <- aus_livestock |>
filter(Animal == "Pigs" & State == "Victoria") |>
mutate(Count_Thousands = Count / 1000)
fit <- pigs |>
model(SimExpSm = ETS(Count_Thousands ~ error("A") + trend("N") + season("N")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| Animal | State | .model | term | estimate |
|---|---|---|---|---|
| Pigs | Victoria | SimExpSm | alpha | 0.32 |
| Pigs | Victoria | SimExpSm | l[0] | 100.49 |
pigs_since_2010 <- pigs |>
filter(year(Month) >= 2010)
fit_augment <- augment(fit)
fit_augment_since_2010 <- fit_augment |>
filter(year(Month) >= 2010)
fit |>
forecast(h = 4) |>
autoplot(pigs_since_2010, level = 95, color = "#80CDC1") +
geom_line(aes(y = .fitted), color = "#DFC27D",
data = fit_augment_since_2010) +
labs(title = "Pigs Slaughtered", y = "Number in Thousands")
In the above chart, only data from 2010 onward are displayed so that the Fitted Values and the Forecasts from Simple Exponential Smoothing can be seen more clearly.
- Compute a \(95\%\) prediction interval for the first forecast using \(\hat{y} \pm 1.96s\) where \(s\) is the standard deviation of the residuals. Compare your interval with the interval produced by R.
acc <- accuracy(fit)
RMSE <- round(as.numeric(acc$RMSE), 2)
kable(acc, format = "simple")| Animal | State | .model | .type | ME | RMSE | MAE | MPE | MAPE | MASE | RMSSE | ACF1 |
|---|---|---|---|---|---|---|---|---|---|---|---|
| Pigs | Victoria | SimExpSm | Training | -0.0295501 | 9.336336 | 7.190129 | -1.056173 | 8.346368 | 0.7761893 | 0.7508488 | 0.0105559 |
\(s = 9.34\)
fc <- fit |>
forecast(h = 1)
y_hat <- round(as.numeric(fc$.mean), 2)\(\hat{y} = 95.19\)
l <- y_hat - (1.96 * RMSE)
u <- y_hat + (1.96 * RMSE)
pi <- as.data.frame(t(c(l, u)))
colnames(pi) <- c("Lower", "Upper")
rownames(pi) <- "**Calculated Prediction Interval**"
kable(pi, format = "simple")| Lower | Upper | |
|---|---|---|
| Calculated Prediction Interval | 76.8836 | 113.4964 |
The \(95\%\) calculated prediction interval aligns well with the interval produced by setting the level argument of the autoplot function to \(95\).
Exercise 8.5:
Data set global_economy contains the annual Exports from many countries. Select one country to analyse.
data(global_economy)
keep <- c("Country", "Year", "Exports")
canada <- global_economy |>
select(all_of(keep)) |>
filter(Country == "Canada")- Plot the Exports series and discuss the main features of the data.
canada |>
autoplot(Exports) +
labs(title = "Canadian Exports")
In Canada, there was a general upward trend for Exports (as percentage of GDP) from 1960 to 1990, after which the rate of the upward trend increased until 2000. Then there was a steep decline in Exports (as percentage of GDP) until 2009, after which the decline has leveled off. We could be seeing the previous general upward trend recommence or just a pause before further decline.
- Use an ETS(A,N,N) model to forecast the series, and plot the forecasts.
fit <- canada |>
filter(Year <= 2012) |>
model(SimExpSm = ETS(Exports ~ error("A") + trend("N") + season("N")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| Country | .model | term | estimate |
|---|---|---|---|
| Canada | SimExpSm | alpha | 1.00 |
| Canada | SimExpSm | l[0] | 16.97 |
fit_augment <- augment(fit)
fit |>
forecast(h = 10) |>
autoplot(canada, level = 95, color = "#80CDC1") +
geom_line(aes(y = .fitted), color = "#DFC27D",
data = fit_augment) +
labs(title = "Canadian Exports",
subtitle = "Simple Exponential Smoothing Forecasts",
y = "Percentage of GDP")
- Compute the RMSE values for the training data.
fc <- fit |>
forecast(h = 1)
y_hat1 <- round(as.numeric(fc$.mean), 2)
acc <- accuracy(fit)
RMSE1 <- round(as.numeric(acc$RMSE), 2)
kable(acc, format = "simple")| Country | .model | .type | ME | RMSE | MAE | MPE | MAPE | MASE | RMSSE | ACF1 |
|---|---|---|---|---|---|---|---|---|---|---|
| Canada | SimExpSm | Training | 0.2499818 | 1.679944 | 1.266428 | 0.9283237 | 4.417591 | 0.9811533 | 0.9905539 | 0.3144439 |
\(RMSE = 1.68\)
- Compare 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.
fit <- canada |>
filter(Year <= 2012) |>
model(HoltLinTrend = ETS(Exports ~ error("A") + trend("A") + season("N")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| Country | .model | term | estimate |
|---|---|---|---|
| Canada | HoltLinTrend | alpha | 1.00 |
| Canada | HoltLinTrend | beta | 0.23 |
| Canada | HoltLinTrend | l[0] | 16.46 |
| Canada | HoltLinTrend | b[0] | 0.08 |
fit_augment <- augment(fit)
fit |>
forecast(h = 10) |>
autoplot(canada, level = 95, color = "#80CDC1") +
geom_line(aes(y = .fitted), color = "#DFC27D",
data = fit_augment) +
labs(title = "Canadian Exports",
subtitle = "Holt's Linear Trend Method Forecasts",
y = "Percentage of GDP")
fc <- fit |>
forecast(h = 1)
y_hat2 <- round(as.numeric(fc$.mean), 2)
acc <- accuracy(fit)
RMSE2 <- round(as.numeric(acc$RMSE), 2)
kable(acc, format = "simple")| Country | .model | .type | ME | RMSE | MAE | MPE | MAPE | MASE | RMSSE | ACF1 |
|---|---|---|---|---|---|---|---|---|---|---|
| Canada | HoltLinTrend | Training | -0.0517386 | 1.678591 | 1.245532 | -0.0624154 | 4.43095 | 0.9649638 | 0.9897558 | 0.1679992 |
\(RMSE = 1.68\)
Both models use an alpha coefficient close to 1, so more weight is given to more recent observations. The RMSE is the same for both models, but they have different numbers of parameters, so neither is more accurate by that standard. The fitted values from the ETS(A,A,N) model, Holt’s Linear Trend Method, demonstrate higher peaks and lower valleys than the ETS(A,N,N) model, Simple Exponential Smoothing. Holt’s Linear Trend Method has a beta coefficient of 0.23, indicating somewhat frequent changes in trend.
The merits of using a model that puts more weight on more recent observations and allows for a linear trend in the forecasts vs. one that doesn’t is: the analyst can see how things might go assuming the recent conditions continue. If there were reason to believe recent conditions would not continue, Simple Exponential Smoothing might be preferable.
- Compare the forecasts from both methods. Which do you think is best?
The ETS(A,A,N) model is probably overestimating a downward trend for the point forecasts, and the prediction interval it provides gets very wide. The ETS(A,N,N) model forecasts are closer to the test data, and its prediction interval still allows for (less extreme) uncertainty. The ETS(A,N,N) model seems best here.
- Calculate a \(95\%\) prediction interval for the first forecast for each model, using the RMSE values and assuming normal errors. Compare your intervals with those produced using R.
l <- y_hat1 - (1.96 * RMSE1)
u <- y_hat1 + (1.96 * RMSE1)
pi <- as.data.frame(t(c(l, u)))
colnames(pi) <- c("Lower", "Upper")
rownames(pi) <- "**Calculated Prediction Interval: Simple Exponential Smoothing**"
kable(pi, format = "simple")| Lower | Upper | |
|---|---|---|
| Calculated Prediction Interval: Simple Exponential Smoothing | 26.9172 | 33.5028 |
l <- y_hat2 - (1.96 * RMSE2)
u <- y_hat2 + (1.96 * RMSE2)
pi <- as.data.frame(t(c(l, u)))
colnames(pi) <- c("Lower", "Upper")
rownames(pi) <- "**Calculated Prediction Interval: Holt's Linear Trend Method**"
kable(pi, format = "simple")| Lower | Upper | |
|---|---|---|
| Calculated Prediction Interval: Holt’s Linear Trend Method | 26.3872 | 32.9728 |
Again, the \(95\%\) calculated prediction intervals align well with the intervals produced by setting the level argument of the autoplot functions to \(95\).
Exercise 8.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 |>
select(all_of(keep)) |>
filter(Country == "China")
china |>
autoplot(Exports)
lam <- china |>
features(Exports, features = guerrero) |>
pull(lambda_guerrero)
fit <- china |>
filter(Year <= 2012) |>
model(HoltLinTrend = ETS(Exports ~ error("A") + trend("A") + season("N")),
DampedHolt = ETS(Exports ~ error("A") + trend("Ad") + season("N")),
HoltLinTrend_Log = ETS(log(Exports) ~
error("A") + trend("A") + season("N")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| Country | .model | term | estimate |
|---|---|---|---|
| China | HoltLinTrend | alpha | 1.00 |
| China | HoltLinTrend | beta | 0.00 |
| China | HoltLinTrend | l[0] | 4.07 |
| China | HoltLinTrend | b[0] | 0.40 |
| China | DampedHolt | alpha | 1.00 |
| China | DampedHolt | beta | 0.22 |
| China | DampedHolt | phi | 0.80 |
| China | DampedHolt | l[0] | 4.15 |
| China | DampedHolt | b[0] | -0.11 |
| China | HoltLinTrend_Log | alpha | 1.00 |
| China | HoltLinTrend_Log | beta | 0.00 |
| China | HoltLinTrend_Log | l[0] | 1.52 |
| China | HoltLinTrend_Log | b[0] | 0.03 |
fit |>
forecast(h = 25) |>
autoplot(china, level = NULL) +
labs(title = "Chinese Exports",
subtitle = "Forecasts",
y = "Percentage of GDP")
The Damped Holt’s Linear Trend Method is influenced more by the recent downward trend than Holt’s Linear Trend Method, and it continues it, but with the smallest leveling off phi coefficient allowed of 0.8. Forecasts from Holt’s Linear Trend Method are more influenced by the general upward trend from the data prior to the recent decline. So they continue an upward trend. But they are not as optimistic as the forecasts from log transforming the Exports variable, then using the same method. These forecasts assume the rate of the trend will increase exponentially.
Here, the forecasts from the Damped Holt’s Linear Trend Method appear most appropriate, as they fit the test data best.
Exercise 8.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?
data(aus_production)
keep <- c("Quarter", "Gas")
gas <- aus_production |>
select(all_of(keep))
fit <- gas |>
filter(year(Quarter) <= 2007) |>
model(HoltWintAdd = ETS(Gas ~ error("A") + trend("A") + season("A")),
HoltWintMult = ETS(Gas ~ error("M") + trend("A") + season("M")),
DampHoltWintMult = ETS(Gas ~ error("M") + trend("Ad") + season("M")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| .model | term | estimate |
|---|---|---|
| HoltWintAdd | alpha | 0.50 |
| HoltWintAdd | beta | 0.00 |
| HoltWintAdd | gamma | 0.39 |
| HoltWintAdd | l[0] | 12.17 |
| HoltWintAdd | b[0] | 1.16 |
| HoltWintAdd | s[0] | -9.02 |
| HoltWintAdd | s[-1] | 5.29 |
| HoltWintAdd | s[-2] | 8.28 |
| HoltWintAdd | s[-3] | -4.56 |
| HoltWintMult | alpha | 0.66 |
| HoltWintMult | beta | 0.15 |
| HoltWintMult | gamma | 0.09 |
| HoltWintMult | l[0] | 5.88 |
| HoltWintMult | b[0] | 0.03 |
| HoltWintMult | s[0] | 0.93 |
| HoltWintMult | s[-1] | 1.18 |
| HoltWintMult | s[-2] | 1.07 |
| HoltWintMult | s[-3] | 0.82 |
| DampHoltWintMult | alpha | 0.65 |
| DampHoltWintMult | beta | 0.15 |
| DampHoltWintMult | gamma | 0.08 |
| DampHoltWintMult | phi | 0.98 |
| DampHoltWintMult | l[0] | 5.88 |
| DampHoltWintMult | b[0] | 0.07 |
| DampHoltWintMult | s[0] | 0.93 |
| DampHoltWintMult | s[-1] | 1.18 |
| DampHoltWintMult | s[-2] | 1.07 |
| DampHoltWintMult | s[-3] | 0.82 |
fit |>
forecast(h = 24) |>
autoplot(gas, level = NULL) +
labs(title = "Gas Production",
subtitle = "Forecasts",
y = "Petajoules")
We plot the full data here to answer the first question. Multiplicative seasonality is required because the seasonal variations are not constant throughout the series. They change proportional to the level of the series, as can clearly be seen when looking at the difference in seasonal variation between the beginning, middle, and end of the time series. Seasonal variation goes from pretty small to pretty large.
Next we zoom in on the most recent data to compare forecasts.
gas_since_2005 <- gas |>
filter(year(Quarter) >= 2005)
fit |>
forecast(h = 24) |>
autoplot(gas_since_2005, level = NULL) +
labs(title = "Gas Production",
subtitle = "Forecasts",
y = "Petajoules")
Even though multiplicative seasonality is necessary here, using additive seasonality actually comes closer to fitting the test data. Damping the trend does improve the multiplicative seasonality forecasts slightly, but the damped multiplicative forecasts are still not as close to the test data as the additive forecasts. We aren’t suggesting we should use additive seasonality when we know it’s not appropriate, however.
Exercise 8.8:
Recall your retail time series data (from Exercise 7 in Section 2.10).
data(aus_retail)
set.seed(1221)
my_series <- aus_retail |>
filter(`Series ID` == sample(aus_retail$`Series ID`, 1))
my_series |>
autoplot(Turnover)
- Why is multiplicative seasonality necessary for this series?
Multiplicative seasonality is again necessary here because the seasonal variation is not constant throughout the time series. The seasonal variation at the beginning levels of the time series is smaller than at the middle/end levels of the time series.
- Apply Holt-Winters’ multiplicative method to the data. Experiment with making the trend damped.
fit <- my_series |>
filter(year(Month) < 2011) |>
model(HoltWintMult = ETS(Turnover ~ error("M") + trend("A") + season("M")),
DampHoltWintMult = ETS(Turnover ~
error("M") + trend("Ad") + season("M")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| State | Industry | .model | term | estimate |
|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | alpha | 0.60 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | beta | 0.00 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | gamma | 0.13 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | l[0] | 2.30 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | b[0] | 0.04 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[0] | 1.13 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-1] | 0.95 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-2] | 0.87 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-3] | 1.43 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-4] | 0.99 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-5] | 0.92 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-6] | 0.94 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-7] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-8] | 0.95 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-9] | 0.90 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-10] | 0.98 |
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | s[-11] | 0.93 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | alpha | 0.58 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | beta | 0.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | gamma | 0.13 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | phi | 0.98 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | l[0] | 2.45 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | b[0] | 0.05 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[0] | 1.16 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-1] | 0.94 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-2] | 0.87 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-3] | 1.42 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-4] | 0.98 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-5] | 0.93 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-6] | 0.93 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-7] | 1.01 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-8] | 0.96 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-9] | 0.90 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-10] | 0.96 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-11] | 0.93 |
fit |>
forecast(h = 168) |>
autoplot(my_series, level = NULL) +
labs(title = "Australian Turnover (My Series)",
subtitle = "Forecasts",
y = "Turnover")
- Compare the RMSE of the one-step forecasts from the two methods. Which do you prefer?
acc <- accuracy(fit)
kable(acc, format = "simple")| State | Industry | .model | .type | ME | RMSE | MAE | MPE | MAPE | MASE | RMSSE | ACF1 |
|---|---|---|---|---|---|---|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | HoltWintMult | Training | -0.0373398 | 0.6889633 | 0.5020404 | -0.9923539 | 7.025620 | 0.4540452 | 0.4442017 | 0.0351070 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | Training | 0.0115891 | 0.6893085 | 0.4991133 | -0.3248674 | 6.957762 | 0.4513980 | 0.4444243 | 0.0476988 |
The RMSE is practically the same for the two methods. The Damped method provides forecasts that are closer to the test data, and it would be preferable here.
- Check that the residuals from the best method look like white noise.
fit_dhwm_only <- my_series |>
filter(year(Month) < 2011) |>
model(DampHoltWintMult = ETS(Turnover ~
error("M") + trend("Ad") + season("M")))
fit_dhwm_only |>
gg_tsresiduals()
There are a few spikes outside the blue bounds of the ACF plot. None of them are very large, but the presence of multiple spikes indicates the data are not white noise.
aug <- fit_dhwm_only |>
augment() |>
features(.innov, ljung_box, lag = 24)
kable(aug, format="simple")| State | Industry | .model | lb_stat | lb_pvalue |
|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | 48.86096 | 0.0019718 |
The Ljung-Box test also confirms that the residuals are distinguishable from white noise.
- Now find the test set RMSE, while training the model to the end of 2010. Can you beat the seasonal naive approach from Exercise 7 in Section 5.11?
fit <- my_series |>
filter(year(Month) < 2011) |>
model(DampHoltWintMult = ETS(Turnover ~
error("M") + trend("Ad") + season("M")),
SeasonalNaive = SNAIVE(Turnover))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| State | Industry | .model | term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | alpha | 0.58 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | beta | 0.00 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | gamma | 0.13 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | phi | 0.98 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | l[0] | 2.45 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | b[0] | 0.05 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[0] | 1.16 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-1] | 0.94 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-2] | 0.87 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-3] | 1.42 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-4] | 0.98 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-5] | 0.93 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-6] | 0.93 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-7] | 1.01 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-8] | 0.96 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-9] | 0.90 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-10] | 0.96 | NA | NA | NA |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-11] | 0.93 | NA | NA | NA |
fit |>
forecast(h = 168) |>
autoplot(my_series, level = NULL) +
labs(title = "Australian Turnover (My Series)",
subtitle = "Forecasts",
y = "Turnover")
acc <- accuracy(fit)
kable(acc, format = "simple")| State | Industry | .model | .type | ME | RMSE | MAE | MPE | MAPE | MASE | RMSSE | ACF1 |
|---|---|---|---|---|---|---|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | Training | 0.0115891 | 0.6893085 | 0.4991133 | -0.3248674 | 6.957762 | 0.451398 | 0.4444243 | 0.0476988 |
| Australian Capital Territory | Newspaper and book retailing | SeasonalNaive | Training | 0.2258258 | 1.5510144 | 1.1057057 | 2.1931710 | 14.952835 | 1.000000 | 1.0000000 | 0.8094843 |
my_series_since_2011 <- my_series |>
filter(year(Month) >= 2011)
fit |>
forecast(h = 168) |>
autoplot(my_series_since_2011, level = NULL) +
labs(title = "Australian Turnover (My Series)",
subtitle = "Forecasts",
y = "Turnover")
Yes, the Damped Holt’s Linear Trend Method beats the Seasonal Naive method.
Exercise 8.9:
For the same retail data, try an STL decomposition applied to the Box-Cox transformed series, followed by ETS on the seasonally adjusted data. How does that compare with your best previous forecasts on the test set?
lam <- my_series |>
features(Turnover, features = guerrero) |>
pull(lambda_guerrero)
my_series_transformed <- my_series |>
mutate(Turnover_Transformed = nthroot(Turnover, n = 4)) |>
model(STL = STL(Turnover_Transformed)) |>
components() |>
select(-.model)
fit <- my_series_transformed |>
filter(year(Month) < 2011) |>
model(DampHoltWintMult = ETS(season_adjust^4 ~
error("M") + trend("Ad") + season("M")))
fit_print <- tidy(fit)
fit_print$estimate <- format(round(fit_print$estimate, 2), scientific = FALSE)
kable(fit_print, format = "simple")| State | Industry | .model | term | estimate |
|---|---|---|---|---|
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | alpha | 0.67 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | beta | 0.04 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | gamma | 0.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | phi | 0.91 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | l[0] | 2.37 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | b[0] | 0.09 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[0] | 1.01 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-1] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-2] | 1.01 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-3] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-4] | 0.98 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-5] | 0.99 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-6] | 0.99 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-7] | 0.99 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-8] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-9] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-10] | 1.00 |
| Australian Capital Territory | Newspaper and book retailing | DampHoltWintMult | s[-11] | 1.01 |
my_series_transformed_since_2011 <- my_series_transformed |>
filter(year(Month) >= 2011)
fit |>
forecast(h = 168) |>
autoplot(my_series_transformed_since_2011, level = NULL) +
labs(title = "Australian Turnover (My Series)",
subtitle = "Forecasts",
y = "Turnover")
These forecasts are better than any previous forecasts on the test set.