Q1

Produce forecasts for the following series using whichever of NAIVE(y), SNAIVE(y) or RW(y ~ drift()) is more appropriate in each case:

  • Australian Population (global_economy)

From the graph below we can see that we have a strong upward trend with no seasonal trends. Because of the most appropriate model to choose is the Drift model.

global_economy %>% 
  filter(Country == "Australia") %>% 
  select(Year,Population) -> aus_pop
aus_pop %>% autoplot(Population)

From the graph below note that the model captures the strong upward trend.

aus_fit <- aus_pop %>% 
  model(Drift = RW(Population ~ drift()))

aus_fit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.var = "Population", level = NULL)+
  autolayer(filter_index(aus_pop, "1960" ~ .))
## Plot variable not specified, automatically selected `.vars = Population`

  • Bricks (aus_production)

As we can see from the graph below we want to capture seasonal trend so we will use the seasonal naive model.

Bricks <- aus_production %>% select(Bricks,Quarter)
Bricks %>% autoplot(Bricks)

As we can see from the graph below the Seasonal Naive model looks like it did a good job at continuing the seasonal pattern.

Bricks_fit <- aus_production %>% 
  filter(!is.na(Bricks)) %>% 
  model(SNAIVE(Bricks ~ lag("year"))) 

Bricks_fit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.var = "Bricks", level = NULL)+
  autolayer(filter_index(Bricks, "1960" ~ .))
## Plot variable not specified, automatically selected `.vars = Bricks`

  • NSW Lambs (aus_livestock)

From the plot below there doesn’t seem to be an upward or downward trend or a seasonal trend so we will use the Naive model.

aus_lambs <- aus_livestock %>% 
  filter(Animal == "Lambs",
         State == "New South Wales")
aus_lambs %>% autoplot(Count)

# lambs_fit <- aus_lambs %>% 
#   model(NAIVE(Count))
lambs_fit <- aus_lambs %>% 
  model(SNAIVE(Count ~ lag("year")))

lambs_fit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.vars = "Count", level = NULL)+
  autolayer(filter_index(aus_lambs, "1960" ~ .))
## Plot variable not specified, automatically selected `.vars = Count`

  • Household wealth (hh_budget)

From the graph below we can see that there seem to be a business cycle but no seasonal trend. There is also a strong upward trend for the last 7-8 years so we will perform a drift model.

HW <- hh_budget %>% 
  index_by(Year) %>% 
  summarise(Wealth = sum(Wealth)/n()) #average household wealth between countries
HW %>% autoplot(Wealth)

Looking at our model below it doesn’t look too bad.

HW_fit <- HW %>% 
  model(Drift = RW(Wealth ~ drift()))

HW_fit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.vars = "Wealth", level = NULL)+
  autolayer(filter_index(HW, "1960" ~ .))
## Plot variable not specified, automatically selected `.vars = Wealth`

  • Australian takeaway food turnover (aus_retail)

From the graph below we can see we have a general upward trend. A seasonal trend is not clear just looking at the graph (appears to fluctate randomly).

aus_tft <-
  aus_retail %>% 
  filter(stringr::str_detect(State,"Australian") &
           stringr::str_detect(Industry,"takeaway food")) %>% 
  select(c(Month,Turnover))
aus_tft %>% autoplot(Turnover)

Below we can see our drift forecast.

aus_tfit <- aus_tft %>% 
  model(Drift = RW(Turnover ~ drift()))

aus_tfit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.vars = "Turnover", level = NULL)+
  autolayer(filter_index(aus_tft, "1983" ~ .))
## Plot variable not specified, automatically selected `.vars = Turnover`

Q2

Apply a seasonal naive 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.

The residuals look random and there mean seems to be close to 0 so it looks like our seasonal naive model accounts for all the available information.

# 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()

Then looking at our forecast below we can see that our the seasonal naive method works well in this case.

fit %>% 
  forecast(h = "5 years") %>% 
  autoplot(.vars = "Beer", level = NULL)+
  autolayer(filter_index(recent_production, "1992" ~ .))
## Plot variable not specified, automatically selected `.vars = Beer`

Q3

Are the following statements true or false? Explain your answer.

  1. Good forecast methods should have normally distributed residuals.

True, Our forecast should account for all available information. That is, the remaining error (residuals) should not be correlated and thus making them normally distributed.

  1. A model with small residuals will give good forecasts.

True (with condition). Given that the model is not overfitted this is true. The best forecast will come from a model that is not overfitted or under fitted.

  1. The best measure of forecast accuracy is MAPE.

False, All the different measures of error are good for somethings. Some are better in certain scenarios.

  1. If your model doesn’t forecast well, you should make it more complicated.

False. You don’t necessarily want an overly complicated model. The best thing to do is to find a simple model that works and try to account for more variability if you can.

  1. Always choose the model with the best forecast accuracy as measured on the test set.

False. If you know you model is overfitted (which will produce high accuracy) then you should not predict with that model.

Q4

For your retail time series (from Exercise 8 in Section 2.10):

set.seed(12345678)
myseries <- aus_retail %>%
  filter(`Series ID` == sample(aus_retail$`Series ID`, 1))
  1. Create a training dataset consisting of observations before 2011.
myseries_train <- myseries %>%
  filter(year(Month) < 2011)
myseries_train %>% tail()
## # A tsibble: 6 x 5 [1M]
## # Key:       State, Industry [1]
##   State              Industry                      `Series ID`    Month Turnover
##   <chr>              <chr>                         <chr>          <mth>    <dbl>
## 1 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Jul     16.1
## 2 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Aug     13.8
## 3 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Sep     13.6
## 4 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Oct     12.3
## 5 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Nov     11.7
## 6 Northern Territory Clothing, footwear and perso~ A3349767W   2010 Dec     17.9
  1. Check that your data have been split appropriately by producing the following plot.

The plot below colors the training data set red. As we can see our data has been split appropriately.

autoplot(myseries, Turnover) +
  autolayer(myseries_train, Turnover, colour = "red")

  1. Calculate seasonal naive forecasts using SNAIVE() applied to your training data (myseries_train).
fit <- myseries_train %>%
  model(SNAIVE(Turnover))
  1. Check the residuals. Do the residuals appear to be uncorrelated and normally distributed?

The acf plot shows some unaccounted for trend. The distributions of the residuals are normally distributed but it appears to not be cantered at 1.

fit %>% gg_tsresiduals()

  1. Produce forecasts for the test data.

As we can see below the CI are much larger than what actually happened.

fc <- fit %>%
  forecast(new_data = anti_join(myseries, myseries_train))
## Joining, by = c("State", "Industry", "Series ID", "Month", "Turnover")
fc %>% autoplot(myseries)

  1. Compare the accuracy of your forecasts against the actual values.

One thing to note is that we were predicting 8 years out. As we can see from the error measurements below our trainging data set (for the most part) had a lower error then the test data. This is to be expected since it was used to train the model. However, Some of the differences are large so we might suspect that the model was under fitted or over fitted.

bind_rows(
  accuracy(fit),
  accuracy(fc, myseries)
) %>%
  select(-State, -Industry, -.model)
## # A tibble: 2 x 9
##   .type       ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>    <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Training 0.439  1.21 0.915  5.23 12.4   1     1    0.768
## 2 Test     0.836  1.55 1.24   5.94  9.06  1.36  1.28 0.601

Q5

  1. Create a training set for Australian takeaway food turnover (aus_retail) by withholding the last four years as a test set.
set.seed(3.14159)

tft <- 
  aus_retail  %>%
  filter(`Series ID` == sample(aus_retail$`Series ID`, 1))
tft %>% autoplot(Turnover)

tft_train <- 
  tft %>% 
  filter(year(Month) < 2015) #with hold last 4 years
tft_train %>% autoplot(Turnover)

  1. Fit all the appropriate benchmark methods to the training set and forecast the periods covered by the test set.

We can see from above that we might want to transform our data to have equal variance throughout the time. Below we used a box-cox transformation to forecast. The only benchmark method that makes sense to look at is seasonal naive since there is no strong trend but the graph has a lot of variablility.

Lambda <- features(tft_train, Turnover, features = guerrero)

tft_fit <- 
  tft_train %>% 
  model(SNAIVE(box_cox(Turnover,lambda = Lambda$lambda_guerrero)))

# tft_fit %>% 
#   forecast(h = "4 years") %>% 
#   autoplot(.vars = "Turnover_box_cox") +
#   autolayer(filter_index(tft))
# code above produces same graph as below
tft_fc <- tft_fit %>% 
  forecast(new_data = anti_join(tft, tft_train))
## Joining, by = c("State", "Industry", "Series ID", "Month", "Turnover")
tft_fc %>% autoplot(tft)

Our residuals are close to being normally distributed but they it seems to be a little off. Then our acf plot shows that there is still a lot of patterns that are unaccounted for in this simple model.

tft_fit %>% gg_tsresiduals()

  1. Compute the accuracy of your forecasts. Which method does best?

We only did one type of forecast mentioned above. This forecast has a nice RMSE since the training is close to the Test.

bind_rows(
  accuracy(tft_fit),
  accuracy(tft_fc, tft)
) %>%
  select(-State, -Industry, -.model)
## # A tibble: 2 x 9
##   .type       ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>    <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Training  2.12  7.16  5.51  4.47  8.80 1      1    0.615
## 2 Test     -2.93  7.37  5.46 -3.96  7.11 0.991  1.03 0.154
  1. Do the residuals from the best method resemble white noise?

As mentioned in part b they don’t form white noise so we can try to combine some of the simple models to get a more complicated model and see how it does.