Question 8.1

Consider the the number of pigs slaughtered in Victoria, available in the aus_livestock dataset.

  1. Use the \(ETS()\) function to estimate the equivalent model for simple exponential smoothing. Find the optimal values of \(α\) and \(ℓ_0\), and generate forecasts for the next four months.

  2. 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.

Question 8.1a Answer

victoria_pigs_slaughtered <- aus_livestock %>%
  filter(State == "Victoria" & Animal == "Pigs")

victoria_pigs_slaughtered_fit <- victoria_pigs_slaughtered %>%
  model(ETS(Count ~ error("A") + trend("N") + season("N")))

victoria_pigs_slaughtered_forecast <- victoria_pigs_slaughtered_fit %>%
  forecast(h = 4)

victoria_pigs_slaughtered_forecast %>%
  autoplot(victoria_pigs_slaughtered %>% filter(Month >= yearmonth("2010 Jan"))) +
  labs(y = "Count", title = "Victoria Pigs Slaughtered") +
  guides(colour = "none")

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

The output above shows us that the parameter estimate for \(\hat\alpha = 0.322\). The estimate for the initial state \(ℓ_0\) is 100646.6.

Question 8.1b Answer

resid_std <- sd(augment(victoria_pigs_slaughtered_fit)$.resid)

sprintf(
  "Calculated confidence interval: [%f, %f]", 
  victoria_pigs_slaughtered_forecast$.mean[1] - (resid_std * 1.96),
  victoria_pigs_slaughtered_forecast$.mean[1] + (resid_std * 1.96)
  )
## [1] "Calculated confidence interval: [76871.012478, 113502.102384]"

The output above shows us that the 95% prediction interval for the first forecast is between 76871.01 and 113502.1 when computed using the \(\hat{y} \pm 1.96s\) formula.

hilo(victoria_pigs_slaughtered_forecast$Count, 95)
## <hilo[4]>
## [1] [76854.79, 113518.3]95 [75927.17, 114445.9]95 [75042.22, 115330.9]95
## [4] [74194.54, 116178.6]95

The output above shows us that for that the first prediction interval using a 95% CI is from 76854.79 to 113518.3. So this interval is pretty close to the interval calculated from the \(\hat{y} \pm 1.96s\) formula.

Question 8.5

Data set global_economy contains the annual Exports from many countries. Select one country to analyse.

  1. Plot the Exports series and discuss the main features of the data.

  2. Use an ETS(A,N,N) model to forecast the series, and plot the forecasts.

  3. Compute the RMSE values for the training data.

  4. 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.

  5. Compare the forecasts from both methods. Which do you think is best?

  6. 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.

Question 8.5a Answer

Since this textbook likes Australia so much, I decided to pick Australia to answer all of the parts to this question.

aus_economy <- global_economy %>%
  filter(Country == "Australia")

aus_economy %>%
  autoplot(Exports)

So this plot shows us an increasing trend in the series. The trend looks like it has begun to taper off starting in 2010.

Question 8.5b Answer

The ETS(A, N, N) model is also known as Holt’s linear trend method. ETS(A, Ad, N) is damped Holt’s method. The forecasts resulting from these models are shown in the plot below.

aus_economy %>%
  model(
    SES = ETS(Exports ~ error("A") + trend("N") + season("N"))
  ) %>%
  forecast(h = 15) %>%
  autoplot(aus_economy) +
  labs(title = "Australian population") +
  guides(colour = guide_legend(title = "Forecast"))

Question 8.5c Answer

aus_economy %>%
  stretch_tsibble(.init = 10) %>%
  model(
    SES = ETS(Exports ~ error("A") + trend("N") + season("N"))
  ) %>%
  forecast(h = 15) %>%
  accuracy(aus_economy)
## Warning: The future dataset is incomplete, incomplete out-of-sample data will be treated as missing. 
## 15 observations are missing between 2018 and 2032
## # A tibble: 1 × 11
##   .model Country   .type    ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>  <fct>     <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 SES    Australia Test   1.55  2.38  1.90  8.43  10.4  1.93  1.93 0.667

The output above shows us that the RMSE for the training data is 2.3812.

Question 8.5d Answer

aus_economy %>%
  model(
    SES = ETS(Exports ~ error("A") + trend("N") + season("N")),
    Holt = ETS(Exports ~ error("A") + trend("A") + season("N"))
  ) %>%
  forecast(h = 15) %>%
  autoplot(aus_economy, level = NULL) +
  labs(title = "Australian population") +
  guides(colour = guide_legend(title = "Forecast"))

Simple exponential smoothing is basically the naive method. It forecasts the very last observation for all forecast periods. This would be beneficial for random walk data (stock data). However, the increasing trend of this data makes a SES model inappropriate for this type of data. Holt’s method looks like it works well with this data, since the exports show an increasing trend as the years go by.

Question 8.5e Answer

aus_economy %>%
  stretch_tsibble(.init = 10) %>%
  model(
    SES = ETS(Exports ~ error("A") + trend("N") + season("N")),
    Holt = ETS(Exports ~ error("A") + trend("A") + season("N"))
  ) %>%
  forecast(h = 15) %>%
  accuracy(aus_economy)
## Warning: The future dataset is incomplete, incomplete out-of-sample data will be treated as missing. 
## 15 observations are missing between 2018 and 2032
## # A tibble: 2 × 11
##   .model Country   .type    ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>  <fct>     <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Holt   Australia Test  0.932  2.30  1.84  5.47  10.3  1.87  1.86 0.776
## 2 SES    Australia Test  1.55   2.38  1.90  8.43  10.4  1.93  1.93 0.667
aus_economy_model <- aus_economy %>%
  model(
    SES = ETS(Exports ~ error("A") + trend("N") + season("N")),
    Holt = ETS(Exports ~ error("A") + trend("A") + season("N"))
  )

aus_economy_forecasts <- aus_economy_model %>%
  forecast(h = 15)

aus_economy_forecasts %>%
  group_by(.model) %>%
  mutate(interval = hilo(Exports, 95)) %>%
  filter(row_number() == 1)
## # A tsibble: 2 x 6 [1Y]
## # Key:       Country, .model [2]
## # Groups:    .model [2]
##   Country   .model  Year    Exports .mean               interval
##   <fct>     <chr>  <dbl>     <dist> <dbl>                 <hilo>
## 1 Australia SES     2018 N(21, 1.4)  20.6 [18.31970, 22.89462]95
## 2 Australia Holt    2018 N(21, 1.3)  20.8 [18.57028, 23.10700]95

If we go off of the RMSE, than we see that Holt’s method has a better fit to the training data than SES. If we go off of the 95% CI prediction interval for the first forecast, we can see that Holt’s method is still better (by a small margin). Therefore, I think that Holt’s method is the better one out of the two.

Question 8.5f Answer

aus_economy_augment <- aus_economy_model %>%
  augment()

resid_std_ses <- sd((aus_economy_augment %>% filter(.model == "SES"))$.resid)
resid_std_holt <- sd((aus_economy_augment %>% filter(.model == "Holt"))$.resid)

sprintf(
  "Calculated confidence interval (SES): [%f, %f]", 
  aus_economy_forecasts$.mean[1] - (resid_std_ses * 1.96),
  aus_economy_forecasts$.mean[1] + (resid_std_ses * 1.96)
  )
## [1] "Calculated confidence interval (SES): [18.386715, 22.827603]"
sprintf(
  "Calculated confidence interval (Holt's): [%f, %f]", 
  aus_economy_forecasts$.mean[1] - (resid_std_holt * 1.96),
  aus_economy_forecasts$.mean[1] + (resid_std_holt * 1.96)
  )
## [1] "Calculated confidence interval (Holt's): [18.399257, 22.815061]"

Both of the calculated confidence intervals that were calculated above are roughly the same as the calculated confidence intervals that are shown in the interval column in the tsibble in “Question 8.5e Answer”. However, the left and right bounds in the calculated SES confidence interval are closer in magnitude to the SES confidence intervals generated from R, compared to the left and right bounds to the calculated and R-computed Holt’s confidence intervals.

Question 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.]

Question 8.6 Answer

Since we are predicting based on GDP, we should probably change it to GDP per capita (3.1 - Population adjustments), since GDP is affected by population changes.

china_economy <- global_economy %>% 
  filter(Country == "China") %>%
  mutate(GDP_per_capita = GDP/Population)

china_economy %>%
  autoplot(GDP_per_capita)

china_economy %>%
  model(
    `SES` = ETS(GDP_per_capita ~ error("A") + trend("N") + season("N")),
    `Damped Holt's method (phi = 0.8)` = ETS(GDP_per_capita ~ error("A") + trend("Ad", phi = 0.8) + season("N")),
    `Damped Holt's method (phi = 0.85)` = ETS(GDP_per_capita ~ error("A") + trend("Ad", phi = 0.85) + season("N")),
    `Damped Holt's method (phi = 0.9)` = ETS(GDP_per_capita ~ error("A") + trend("Ad", phi = 0.9) + season("N")),
    `Damped Holt's method (phi = 0.98)` = ETS(GDP_per_capita ~ error("A") + trend("Ad", phi = 0.98) + season("N")),
  ) %>%
  forecast(h = 15) %>%
  autoplot(china_economy, level = NULL) +
  labs(title = "China GDP Per Capita") +
  guides(colour = guide_legend(title = "Forecast"))

As we can see here, as we decrease the damping parameter \(\phi\), the forecasts begin to level out to a constant at an earlier time. By forecasting 15 years into the future, we can see that we see significant tapering at \(\phi = 0.8\). This tapering becomes less evident as we increase \(\phi\).

lambda <- china_economy %>%
  features(GDP_per_capita, features = guerrero) |>
  pull(lambda_guerrero)

china_economy %>%
  model(
    `SES` = ETS(box_cox(GDP_per_capita, lambda) ~ error("A") + trend("N") + season("N")),
    `Damped Holt's method (phi = 0.8)` = ETS(box_cox(GDP_per_capita, lambda) ~ error("A") + trend("Ad", phi = 0.8) + season("N")),
    `Damped Holt's method (phi = 0.85)` = ETS(box_cox(GDP_per_capita, lambda) ~ error("A") + trend("Ad", phi = 0.85) + season("N")),
    `Damped Holt's method (phi = 0.9)` = ETS(box_cox(GDP_per_capita, lambda) ~ error("A") + trend("Ad", phi = 0.9) + season("N")),
    `Damped Holt's method (phi = 0.98)` = ETS(box_cox(GDP_per_capita, lambda) ~ error("A") + trend("Ad", phi = 0.98) + season("N")),
  ) %>%
  forecast(h = 15) %>%
  autoplot(china_economy, level = NULL) +
  labs(title = "China GDP Per Capita") +
  guides(colour = guide_legend(title = "Forecast"))

Applying a Box-Cox Transformation has resulted in the following forecasts. As we increase \(phi\), the slopes of the forecasts increase. With that being said, current UN projections suggest that the population of China will stagnate by 2029. Therefore, the forecasts from Box-Cox transforming the original data should definitely not be reported.

Question 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?

Question 8.7 Answer

Let’s use the methodology outlined in section 8.6 of the textbook to pick out a model through minimizing the AICc.

aus_production %>%
  autoplot(Gas)

This dataset would probably benefit from a Box-Cox transformation. But let’s create two different ETS models; one where the data was transformed with Box-Cox, and one where the data was left alone. For this first ETS model, we will leave the data alone.

aus_gas <- aus_production %>%
  select(Gas)

fit <- aus_gas %>%
  model(ETS(Gas))

report(fit)
## Series: Gas 
## Model: ETS(M,A,M) 
##   Smoothing parameters:
##     alpha = 0.6528545 
##     beta  = 0.1441675 
##     gamma = 0.09784922 
## 
##   Initial states:
##      l[0]       b[0]      s[0]    s[-1]    s[-2]     s[-3]
##  5.945592 0.07062881 0.9309236 1.177883 1.074851 0.8163427
## 
##   sigma^2:  0.0032
## 
##      AIC     AICc      BIC 
## 1680.929 1681.794 1711.389
fit %>%
  accuracy()
## # A tibble: 1 × 10
##   .model   .type        ME  RMSE   MAE   MPE  MAPE  MASE RMSSE    ACF1
##   <chr>    <chr>     <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>   <dbl>
## 1 ETS(Gas) Training -0.115  4.60  3.02 0.199  4.08 0.542 0.606 -0.0131
fit %>% 
  forecast(h = 15) %>%
  autoplot(aus_gas) +
  labs(title = "Australia Gas Production (Original Data)") +
  guides(colour = guide_legend(title = "Forecast"))

Now let’s transform the data using Box-Cox (using the methodology outlined in section 3.1 - Transformations and Adjustments) and report the model coefficients for that model.

lambda <- aus_gas |>
  features(Gas, features = guerrero) |>
  pull(lambda_guerrero)

aus_gas <- aus_gas %>%
  mutate(aus_gas_transformed = box_cox(Gas, lambda))

fit_box_cox <- aus_gas %>%
  model(ETS(aus_gas_transformed))

report(fit_box_cox)
## Series: aus_gas_transformed 
## Model: ETS(A,A,A) 
##   Smoothing parameters:
##     alpha = 0.81796 
##     beta  = 0.1216283 
##     gamma = 0.0001000577 
## 
##   Initial states:
##      l[0]       b[0]       s[0]     s[-1]     s[-2]      s[-3]
##  1.888944 0.03198818 -0.1065298 0.2508721 0.1003452 -0.2446875
## 
##   sigma^2:  0.0059
## 
##      AIC     AICc      BIC 
## 63.28945 64.15484 93.74991
fit_box_cox %>%
  accuracy()
## # A tibble: 1 × 10
##   .model         .type       ME   RMSE    MAE      MPE  MAPE  MASE RMSSE    ACF1
##   <chr>          <chr>    <dbl>  <dbl>  <dbl>    <dbl> <dbl> <dbl> <dbl>   <dbl>
## 1 ETS(aus_gas_t… Trai… -9.34e-4 0.0751 0.0570 -0.00928  1.52 0.471 0.402 -0.0283
fit_box_cox %>% 
  forecast(h = 15) %>%
  autoplot(aus_gas) +
  labs(title = "Australia Gas Production (Box Cox)") +
  guides(colour = guide_legend(title = "Forecast"))

After transforming the data using Box-Cox, we can see that we have a much smaller RMSE. With that being said, the Box-Cox transformation gives as an additive seasonality model while the original data gives us a multiplicative seasonality model. The heights of the seasonal periods change over time in the original data. This is why in the model with the original data, multiplicative seasonality is used. Since one part of the question explicitly says to explain why multiplicative seasonality is necessary, we will proceed with the model produced from the original data.

damped_fit <- aus_gas %>%
  model(fit = ETS(Gas  ~ trend('Ad', phi = 0.92)))

damped_fit %>% 
  forecast(h = 15) %>%
  autoplot(aus_gas) +
  labs(title = "Australia Gas Production (Damped Model)") +
  guides(colour = guide_legend(title = "Forecast"))

The confidence intervals are still pretty high. Also, there seems to have not been much of a change when comparing the forecasts to the damped model to the forecasts using the model constructed from the original data.

Question 8.8

Recall your retail time series data (from Exercise 8 in Section 2.10).

  1. Why is multiplicative seasonality necessary for this series?
  2. Apply Holt-Winters’ multiplicative method to the data. Experiment with making the trend damped.
  3. Compare the RMSE of the one-step forecasts from the two methods. Which do you prefer?
  4. Check that the residuals from the best method look like white noise.
  5. Now find the test set RMSE, while training the model to the end of 2010. Can you beat the seasonal naïve approach from Exercise 7 in Section 5.11?

Question 8.8a Answer

set.seed(23987)
myseries <- aus_retail %>%
  filter(`Series ID` == sample(aus_retail$`Series ID`,1))
myseries %>%
  autoplot(Turnover)

The plot above shows us that multiplicative seasonality is necessary for this dataset since the heights of the seasonal periods are changing as the years go by.

Question 8.8b Answer

fit <- myseries |>
  model(
    Multiplicative = ETS(Turnover ~ error("M") + trend("A") + season("M")),
    Damped = ETS(Turnover ~ error("M") + trend("Ad") + season("M"))
  )
fc <- fit |> forecast(h = "3 years")
fc |>
  autoplot(myseries, level = NULL) +
  labs(title="Turnover Time Series",
       y="turnover") +
  guides(colour = guide_legend(title = "Forecast"))

Question 8.8c Answer

fit %>%
  accuracy()
## # A tibble: 2 × 12
##   State   Indus…¹ .model .type    ME  RMSE   MAE   MPE  MAPE  MASE RMSSE    ACF1
##   <chr>   <chr>   <chr>  <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>   <dbl>
## 1 Wester… Other … Multi… Trai… 0.264  8.64  6.07 0.239  3.48 0.481 0.496  0.0972
## 2 Wester… Other … Damped Trai… 0.726  8.32  5.84 0.401  3.35 0.462 0.477 -0.0416
## # … with abbreviated variable name ¹​Industry

The RMSE is slightly lower for the Damped model compared to the Multiplicative model, which makes the Damped model the preferable model.

Question 8.8d Answer

fit %>%
  select(Damped) %>%
  gg_tsresiduals()

The innovation residuals plot looks like white noise. Therefore, we satisfied the nonconstant variance assumption.

Question 8.8e Answer

We can use the methodology outlined in Section 5.8 - Evaluating point forecast accuracy, to create training and testing data.

myseries_train <- myseries %>%
  filter(year(Month) < 2011)

myseries_test <- myseries %>%
  filter(year(Month) >= 2011)
myseries_fit <- myseries_train %>%
  model(
    Multiplicative = ETS(Turnover ~ error("M") + trend("A") + season("M")),
    Damped = ETS(Turnover ~ error("M") + trend("Ad", phi = 0.8) + season("M"))
  )

myseries_forecasts <- myseries_fit |>
  forecast(h = nrow(myseries_test))

myseries_forecasts |>
  autoplot(
    myseries,
    level = NULL
  ) +
  labs(
    y = "Megalitres",
    title = "Forecasts for quarterly beer production"
  ) +
  guides(colour = guide_legend(title = "Forecast"))

accuracy(myseries_forecasts, myseries_test)
## # A tibble: 2 × 12
##   .model     State Indus…¹ .type    ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>      <chr> <chr>   <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Damped     West… Other … Test  115.  120.  115.   32.5  32.5   NaN   NaN 0.773
## 2 Multiplic… West… Other … Test   59.4  62.3  59.4  16.7  16.7   NaN   NaN 0.657
## # … with abbreviated variable name ¹​Industry

It looks like the Multiplicative method has a much lower RMSE than the Damped method. However, the question refers to the model that was best in Question 8.8b, which is the damped method. In question 5.7f, the RMSE I got was 87.35882. This RMSE that was calculated from the seasonal naive approach is lower than the RMSE for the Damped model. Therefore, we can conclude that for my dataset, the seasonal naive approach is better than the Damped model, but the Multiplicative model has the overall lowest RMSE.

Question 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?

Question 8.9 Answer

First thing we have to do is transform the data using Box-Cox (section 3.1 - Transformations and adjustments). Then we can perform the STL decomposition using the methodology outlined in section 3.6 - STL decomposition. We can extract the STL components using the methodology outlined in section 3.2 - Time series components.

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

myseries_transformed <- myseries %>%
  mutate(Turnover_transformed = box_cox(Turnover, lambda))

myseries_train <- myseries_train %>%
  mutate(Turnover_transformed = box_cox(Turnover, lambda))

myseries_test <- myseries_test %>%
  mutate(Turnover_transformed = box_cox(Turnover, lambda))

myseries_stl_components <- myseries_transformed %>%
  model(
    STL(Turnover_transformed ~ trend(window = 7) +
                   season(window = "periodic"),
    robust = TRUE)) %>%
  components()

myseries_stl_components_train <- myseries_train %>%
  model(
    STL(Turnover_transformed ~ trend(window = 7) +
                   season(window = "periodic"),
    robust = TRUE)) %>%
  components()

myseries_stl_components_test <- myseries_test %>%
  model(
    STL(Turnover_transformed ~ trend(window = 7) +
                   season(window = "periodic"),
    robust = TRUE)) %>%
  components()

myseries_stl_components_train %>%
  autoplot()

head(myseries_stl_components_train)
## # A dable: 6 x 9 [1M]
## # Key:     State, Industry, .model [1]
## # :        Turnover_transformed = trend + season_year + remainder
##   State            Indus…¹ .model    Month Turno…² trend seaso…³ remai…⁴ seaso…⁵
##   <chr>            <chr>   <chr>     <mth>   <dbl> <dbl>   <dbl>   <dbl>   <dbl>
## 1 Western Austral… Other … "STL(… 1982 Apr    5.00  5.19 -0.161  -0.0241    5.16
## 2 Western Austral… Other … "STL(… 1982 May    5.12  5.16 -0.0606  0.0177    5.18
## 3 Western Austral… Other … "STL(… 1982 Jun    5.02  5.14 -0.155   0.0372    5.18
## 4 Western Austral… Other … "STL(… 1982 Jul    5.00  5.11 -0.0874 -0.0247    5.09
## 5 Western Austral… Other … "STL(… 1982 Aug    5.00  5.08 -0.0484 -0.0307    5.05
## 6 Western Austral… Other … "STL(… 1982 Sep    5.08  5.06 -0.0493  0.0686    5.13
## # … with abbreviated variable names ¹​Industry, ²​Turnover_transformed,
## #   ³​season_year, ⁴​remainder, ⁵​season_adjust

We can fit an ETS model on the season_adjust column in the data. We’ll use the methodology outlined in section 8.6 - Estimation and model selection, to select a model by minimizing the AICc.

myseries_season_adjust <- myseries_stl_components %>%
  select(season_adjust)

myseries_season_adjust_test <- myseries_stl_components_test %>%
  select(season_adjust)

myseries_season_adjust_train <- myseries_stl_components_train %>%
  select(season_adjust)

myseries_season_adjust_fit <- myseries_season_adjust_train %>%
  model(ETS(season_adjust))

myseries_forecasts <- myseries_season_adjust_fit |>
  forecast(h = nrow(myseries_test))

report(myseries_season_adjust_fit)
## Series: season_adjust 
## Model: ETS(A,A,N) 
##   Smoothing parameters:
##     alpha = 0.5755745 
##     beta  = 0.001337449 
## 
##   Initial states:
##      l[0]        b[0]
##  5.104853 0.009949651
## 
##   sigma^2:  0.0079
## 
##      AIC     AICc      BIC 
## 352.4649 352.6419 371.6826
accuracy(myseries_forecasts, myseries_season_adjust_test)
## # A tibble: 1 × 10
##   .model             .type    ME  RMSE   MAE   MPE  MAPE  MASE RMSSE  ACF1
##   <chr>              <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 ETS(season_adjust) Test  0.268 0.304 0.269  2.88  2.89   NaN   NaN 0.826
myseries_forecasts %>%
  autoplot(myseries_season_adjust)

The output above shows us that the RMSE is 0.304, which is much, much lower than any of the RMSEs for the previous parts. The final model is an ETS(A, A, N) model(Holt’s linear method), which is to be expected since we have taken the seasonality out of the original data. This model performs the best overall.