Produce forecasts for the following series using whichever of NAIVE(y), SNAIVE(y) or RW(y ~ drift()) is more appropriate in each case:
It seems if seasonality is stronger than trend then Seasonal Naïve works best and if trend is stronger than seasonality then Drift works best as a baseline comparison forecast.
Australian Population (global_economy)
From the graph below, Drift works best. Notably there is no seasonality to population.
### Exercise 5.1 Part a
# Check out packages
library(fpp3)
# Produce forecast
train <- global_economy |>
filter(Country == "Australia") |>
filter_index("1960" ~ "2017")
fit <- train |>
model(
`Drift` = NAIVE(Population ~ drift()),
`Naïve` = NAIVE(Population),
`Seasonal naïve` = SNAIVE(Population)
)
fc <- fit |> forecast(h = 14)
fc |>
autoplot(train, level = NULL) +
autolayer(
filter_index(train, "2007 Q1" ~ .),
colour = "black"
) +
labs(
y = "Total Population",
title = "Forecasts for Australia's Total Population"
) +
guides(colour = guide_legend(title = "Forecast"))
Bricks (aus_production)
From the graph below, Seasonal naïve works best.
### Exercise 5.1 Part b
# Produce forecast
train <- aus_production |>
filter_index("1992 Q1" ~ "2005 Q2")
fit <- train |>
model(
`Drift` = NAIVE(Bricks ~ drift()),
`Naïve` = NAIVE(Bricks),
`Seasonal naïve` = SNAIVE(Bricks)
)
fc <- fit |> forecast(h = 14)
fc |>
autoplot(train, level = NULL) +
autolayer(
filter_index(train, "2007 Q1" ~ .),
colour = "black"
) +
labs(
y = "Millions",
title = "Forecasts for quarterly brick production"
) +
guides(colour = guide_legend(title = "Forecast"))
NSW Lambs (aus_livestock)
From the graph below, Seasonal naïve works best as a baseline comparison forecast.
### Exercise 5.1 Part c
# Produce forecast
train <- aus_livestock |>
filter(Animal == "Lambs") |>
filter(State == "New South Wales") |>
filter_index("2010 Jul" ~ "2018 Nov")
fit <- train |>
model(
`Drift` = NAIVE(Count ~ drift()),
`Naïve` = NAIVE(Count),
`Seasonal naïve` = SNAIVE(Count)
)
fc <- fit |> forecast(h = 14)
fc |>
autoplot(train, level = NULL) +
autolayer(
filter_index(train, "2007 Q1" ~ .),
colour = "black"
) +
labs(
y = "Lambs slaughtered",
title = "Forecasts for monthly lamb slaughter in New South Wales"
) +
guides(colour = guide_legend(title = "Forecast"))
Household wealth (hh_budget).
From the graph below, Drift works best as a baseline comparison forecast.
### Exercise 5.1 Part d
# Produce forecast
train <- hh_budget |>
filter_index("1995" ~ "2016")
fit <- train |>
model(
`Drift` = NAIVE(Wealth ~ drift()),
`Naïve` = NAIVE(Wealth),
`Seasonal naïve` = SNAIVE(Wealth)
)
fc <- fit |> forecast(h = 14)
fc |>
autoplot(train, level = NULL) +
autolayer(
filter_index(train, "2007 Q1" ~ .),
colour = "black"
) +
labs(
y = "Wealth as a percentage of net disposable income",
title = "Forecasts for household wealth by year"
) +
guides(colour = guide_legend(title = "Forecast"))
Australian takeaway food turnover (aus_retail).
Drift seems to be better than Seasonal naïve as our baseline forecast comparison because this doesn’t look like true seasonality and so trend is more important.
### Exercise 5.1 Part e
# Produce forecast
train <- aus_retail |>
filter_index("2010 Jan" ~ "2018 Dec") |>
filter(Industry == "Takeaway food services") |>
filter(State == "Australian Capital Territory")
fit <- train |>
model(
`Drift` = NAIVE(Turnover ~ drift()),
`Naïve` = NAIVE(Turnover),
`Seasonal naïve` = SNAIVE(Turnover)
)
fc <- fit |> forecast(h = 14)
fc |>
autoplot(train, level = NULL) +
autolayer(
filter_index(train, "2007 Q1" ~ .),
colour = "black"
) +
labs(
y = "$Million AUD",
x = "",
title = "Australian takeaway food turnover"
) +
guides(colour = guide_legend(title = "Forecast"))
Use the Facebook stock price (data set gafa_stock) to do the following:
Produce a time plot of the series.
# Exercise 5.2 Part a
# Plot series
gafa_stock |>
filter(Symbol == "FB") |>
autoplot(Adj_Close)
Produce forecasts using the drift method and plot them.
### Exercise 5.2 Part b
# Re-index based on trading days
facebook_stock <- gafa_stock |>
filter(Symbol == "FB") |>
mutate(day = row_number()) |>
update_tsibble(index = day, regular = TRUE)
# Fit the drift model
facebook_fit <- facebook_stock |>
model(
Drift = NAIVE(Close ~ drift())
)
# Produce forecasts
facebook_fc <- facebook_fit |>
forecast(h = 300)
# Plot the forecasts
facebook_fc |>
autoplot(facebook_stock, level = NULL) +
labs(y = "$US",
title = "Facebook closing with 300 day forecast",
subtitle = "(Jan 2014 - Dec 2018)") +
guides(colour = guide_legend(title = "Forecast"))
Show that the forecasts are identical to extending the line drawn between the first and last observations.
I’ve drawn a line between the first and last observations and it’s close to showing that the forecast is an extension of the slope from the start to the end of the historical time series.
### Exercise 5.2 Part c
# Plot the forecasts
facebook_fc |>
autoplot(facebook_stock, level = NULL) +
labs(y = "$US",
title = "Facebook closing: Forecast as slope of Beginning to End",
subtitle = "(Jan 2014 - Dec 2018)") +
guides(colour = guide_legend(title = "Forecast")) +
geom_segment(aes(x = 1, y = 55, xend = 1253, yend = 131, colour = "segment"))
Try using some of the other benchmark functions to forecast the same data set. Which do you think is best? Why?
The Naive is best because this is a random walk so there is no trend to continue. Additionally, while the mean is close it’s an accident of the numbers and we need the naive because it starts from the last historical date.
### Exercise 5.2 Part d
# Re-index based on trading days
facebook_stock <- gafa_stock |>
filter(Symbol == "FB") |>
mutate(day = row_number()) |>
update_tsibble(index = day, regular = TRUE)
# Fit the models
facebook_fit <- facebook_stock |>
model(
Mean = MEAN(Close),
Naïve = NAIVE(Close),
Drift = NAIVE(Close ~ drift())
)
# Produce forecasts
facebook_fc <- facebook_fit |>
forecast(h = 300)
# Plot the forecasts
facebook_fc |>
autoplot(facebook_stock, level = NULL) +
labs(y = "$US",
title = "Facebook closing with 300 day forecast",
subtitle = "(Jan 2014 - Dec 2018)") +
guides(colour = guide_legend(title = "Forecast"))
Apply a seasonal naïve method to the quarterly Australian beer production data from 1992. Check if the residuals look like white noise, and plot the forecasts. The following code will help. What do you conclude?
The residuals are stationary so it looks like there is no trend to extract from the residuals, however in the ACF plot there is a significant spike at four which represents to me that there is seasonal quarterly autocorrelation left in the residuals that could be explained by a better model. Also, while the histogram looks like a normal distribution, it could also be two normal distributions next to each other like a bactrian camel.
Conclusion? We can do better.
### Exercise 5.3
# Extract data of interest
recent_production <- aus_production |>
filter(year(Quarter) >= 1992)
# Define and estimate a model
fit <- recent_production |> model(SNAIVE(Beer))
# Look at the residuals
fit |> gg_tsresiduals()
# Look at some forecasts
fit |> forecast() |> autoplot(recent_production)
Repeat the previous exercise using the Australian Exports series from global_economy and the Bricks series from aus_production. Use whichever of NAIVE() or SNAIVE() is more appropriate in each case.
Australian Exports series from
global_economy
In the case of the Australian exports series from
global_economy
the NAIVE() model is more appropriate.
Additionally the residuals look like normally distributed white noise. One of the critical values in the ACF plot is outside of the threshold but just by a little bit and it’s one out of 17 which is reasonable for 95%. If we wanted to investigate further we could run a Ljung-box test.
### Exercise 5.4 Part a
# Extract data of interest
australian_exports <- global_economy |>
filter(Country == "Australia")
# Define and estimate a model
fit <- australian_exports |>
model(NAIVE(Exports))
# Look at the residuals
fit |> gg_tsresiduals()
# Look at some forecasts
fit |> forecast() |> autoplot(australian_exports)
In the case of the Bricks series from aus_production
the
SNAIVE() model is more appropriate, however when we look at the
residuals we can see that there is autocorrelation outside of the
threshold in the ACF plot which could be extracted in a better model. It
looks like there’s a second seasonality with a period of six or seven
quarters that maybe could be captured with a higher order ARIMA
model.every six quarters there’s a mesocycle. Also note the histogram is
skewed and not normal.
### Exercise 5.4 Part b
# Extract data of interest
recent_production <- aus_production |>
filter_index("1992 Q1" ~ "2005 Q2")
# Define and estimate a model
fit <- recent_production |>
model(SNAIVE(Bricks))
# Look at the residuals
fit |> gg_tsresiduals()
# Look at some forecasts
fit |> forecast(h=10) |> autoplot(recent_production)
For your retail time series (from Exercise 7 in Section 2.10):
Create a training dataset consisting of observations before 2011 using
-Done, see code below with Show
button on right-
### Exercise 5.7 Part a
# Create training dataset
set.seed(20240915)
myseries <- aus_retail |>
filter(`Series ID` == sample(aus_retail$`Series ID`,1))
myseries_train <- myseries |>
filter(year(Month) < 2011)
Check that your data have been split appropriately by producing the following plot.
### Exercise 5.7 Part b
# Plot split
autoplot(myseries, Turnover) +
autolayer(myseries_train, Turnover, colour = "red")
Fit a seasonal naïve model using SNAIVE() applied to your training data (myseries_train).
-Done, see code below with Show
button on right-
### Exercise 5.7 Part c
# Fit seasonal naive model
fit <- myseries_train |>
model(SNAIVE(Turnover))
Check the residuals. Do the residuals appear to be uncorrelated and normally distributed?
The residuals show a significant number of critical values outside of the threshold in the ACF plot so there is significant correlation left to be captured in the model. Also the residuals do not look normally distributed as they are centered on a positive, non-zero value. Which means that there is trend left in the residuals that could be captured in a better model.
### Exercise 5.7 Part d
# Check residuals
fit |> gg_tsresiduals()
Produce forecasts for the test data
### Exercise 5.7 Part e
# Produce forecast
fc <- fit |>
forecast(new_data = anti_join(myseries, myseries_train))
fc |> autoplot(myseries)
Compare the accuracy of your forecasts against the actual values.
The RMSE for the fitted model is ~6.29. The RMSE for the forecasted values is ~12.24. This suggests that the model is less accurate when predicting future values (the forecasted values) compared to how well it performed on historical data (the fitted model).
We would need additional context to evaluate the model however since we’ve previously established there is trend and autocorrelation in the residuals, we know we can do better with another type of model.
### Exercise 5.7 Part f
# Compare forecast accuracy
fit |> accuracy()
## # A tibble: 1 × 12
## State Industry .model .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
## <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 South… Departm… SNAIV… Trai… 2.64 6.30 5.03 2.91 5.47 1 1 1.19e-4
fc |> accuracy(myseries)
## # A tibble: 1 × 12
## .model State Industry .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
## <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 SNAIVE(T… Sout… Departm… Test -10.5 12.2 10.7 -9.77 9.92 2.13 1.95 0.230
How sensitive are the accuracy measures to the amount of training data used?
If there’s a small training data set it’s not possible to obtain a reliable forecast. What we’re using seems to be sufficient training data for modeling, and I do not believe results would be significantly different if we increased or decreased the amount of training data but rather we should pursue better types of models to reasonably forecast this phenomenon. Normally more training data should help but SNAIVE() isn’t sophisticated enough of a model to take advantage of more training data then we already have.
These exercises come from ‘Forecasting: Principles and Practice’ by
Rob J Hyndman and George Athanasopoulos, 3rd Ed
https://otexts.com/fpp3/