DATA624 Spring 2021 HW5 - Exponential Smoothing

7.8.1

Consider the pigs series – the number of pigs slaughtered in Victoria each month.

a.

Use the ses() function in R to find the optimal values of \(\alpha\) and \(\ell_0\) and generate forecasts for the next four months.

First, we begin by inspecting the pigs time series by looking at its information, its head and tail. We see that the time series is captured on a monthly scale. An autoplot of the time series shows no consistent long term trend.

help(pigs)

## starting httpd help server ... done

Number of pigs slaughtered Description Monthly total number of pigs slaughtered in Victoria, Australia (Jan 1980 to Aug 1995).

Usage pigs Format Time series data

Source Makridakis, Wheelwright and Hyndman (1998) Forecasting: methods and applications, John Wiley & Sons: New York. Chapter 7.

head(pigs)

##         Jan    Feb    Mar    Apr    May    Jun
## 1980  76378  71947  33873  96428 105084  95741

tail(pigs)

##         Mar    Apr    May    Jun    Jul    Aug
## 1995 106723  84307 114896 106749  87892 100506

autoplot(pigs) +
  ylab("Monthly No. of Pigs slaughtered in Victoria, AUS") + xlab("Jan 1980 to  Aug 1995")

Next, R’s simple exponential series, “ses”, function is applied to the pigs time series to forecast the next four months. The results are also plotted below.

fc <- ses(pigs, h=4)

autoplot(fc) +
  autolayer(fitted(fc), series="Fitted") +
   ylab("Monthly No. of Pigs slaughtered in Victoria, AUS") + xlab("Dates") + ggtitle("Forecast No. of Pigs slaughtered in Victoria, AUS")

Looking its summary, we see that the optimal level of alpha is, alpha = 0.2971 and the level is, l = 77260.0561.

summary(fc)

## 
## Forecast method: Simple exponential smoothing
## 
## Model Information:
## Simple exponential smoothing 
## 
## Call:
##  ses(y = pigs, h = 4) 
## 
##   Smoothing parameters:
##     alpha = 0.2971 
## 
##   Initial states:
##     l = 77260.0561 
## 
##   sigma:  10308.58
## 
##      AIC     AICc      BIC 
## 4462.955 4463.086 4472.665 
## 
## Error measures:
##                    ME    RMSE      MAE       MPE     MAPE      MASE       ACF1
## Training set 385.8721 10253.6 7961.383 -0.922652 9.274016 0.7966249 0.01282239
## 
## Forecasts:
##          Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## Sep 1995       98816.41 85605.43 112027.4 78611.97 119020.8
## Oct 1995       98816.41 85034.52 112598.3 77738.83 119894.0
## Nov 1995       98816.41 84486.34 113146.5 76900.46 120732.4
## Dec 1995       98816.41 83958.37 113674.4 76092.99 121539.8

	ME	RMSE	MAE	MPE	MAPE	MASE	ACF1
Training set	385.87	10253.6	7961.38	-0.92	9.27	0.8	0.01

b.

Compute a 95% prediction interval for the first forecast using where is the standard deviation of the residuals. Compare your interval with the interval produced by R.

The custom function below, “calc_95_CI”, calculates the 95% Confidence Interval.

calc_95_CI <- function(fc){
  std_dev <- sd(residuals(fc))
  first_fc <- fc$mean[1]
  lower <- first_fc - 1.96*std_dev
  upper <- first_fc + 1.96*std_dev
  Confi_95 <- tibble(lower,upper)
  colnames(Confi_95) <- c('lower', 'upper')
  return (Confi_95)
  
}

The manually calculated 95% Confidence Interval is:

lower	upper
78679.97	118952.8

The R generated confidence interval is:

lower	upper
78611.97	119020.8

The R generated 95% confidence interval seems to be wider than the manually calculated 95% confidence interval. The difference is greater than a simple rounding error.

7.8.5

Data set books contains the daily sales of paperback and hardcover books at the same store. The task is to forecast the next four days’ sales for paperback and hardcover books.

a.

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

From the help function, head and tail, we see that the time series is a daily time series split between Paperback and Hardcover. The autoplot of the series shows a long term upward trend for both Paperback and Hardcover book sales along with cyclical patterns.

help(books)

Sales of paperback and hardcover books Description Daily sales of paperback and hardcover books at the same store.

Usage books Format Bivariate time series containing the following columns:

Paperback Number of paperback sales each day

Hardcover Number of hardcover sales each day

Source Makridakis, Wheelwright and Hyndman (1998) Forecasting: methods and applications, John Wiley & Sons: New York. Exercise 4.5.

head(books)

## Time Series:
## Start = 1 
## End = 6 
## Frequency = 1 
##   Paperback Hardcover
## 1       199       139
## 2       172       128
## 3       111       172
## 4       209       139
## 5       161       191
## 6       119       168

tail(books)

## Time Series:
## Start = 25 
## End = 30 
## Frequency = 1 
##    Paperback Hardcover
## 25       190       214
## 26       182       200
## 27       222       201
## 28       217       283
## 29       188       220
## 30       247       259

autoplot(books) +
  ylab("Daily sales of paperback and hardcover books at the same store") + xlab("Daily Sales")

b.

Use the ses() function to forecast each series, and plot the forecasts.

Below, R’s “ses” function is applied to the hardcover and paperback book sales time series separately. As per instructions, the time forecast is for the next 4 days.

hardcover <- books[,2]
paperback <- books[,1]
hardcover_fc <- ses(hardcover, h=4)
paperback_fc <- ses(paperback, h=4)

summary(hardcover_fc)

## 
## Forecast method: Simple exponential smoothing
## 
## Model Information:
## Simple exponential smoothing 
## 
## Call:
##  ses(y = hardcover, h = 4) 
## 
##   Smoothing parameters:
##     alpha = 0.3283 
## 
##   Initial states:
##     l = 149.2861 
## 
##   sigma:  33.0517
## 
##      AIC     AICc      BIC 
## 315.8506 316.7737 320.0542 
## 
## Error measures:
##                    ME     RMSE      MAE      MPE     MAPE      MASE       ACF1
## Training set 9.166735 31.93101 26.77319 2.636189 13.39487 0.7987887 -0.1417763
## 
## Forecasts:
##    Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 31       239.5601 197.2026 281.9176 174.7799 304.3403
## 32       239.5601 194.9788 284.1414 171.3788 307.7414
## 33       239.5601 192.8607 286.2595 168.1396 310.9806
## 34       239.5601 190.8347 288.2855 165.0410 314.0792

summary(paperback_fc)

## 
## Forecast method: Simple exponential smoothing
## 
## Model Information:
## Simple exponential smoothing 
## 
## Call:
##  ses(y = paperback, h = 4) 
## 
##   Smoothing parameters:
##     alpha = 0.1685 
## 
##   Initial states:
##     l = 170.8271 
## 
##   sigma:  34.8183
## 
##      AIC     AICc      BIC 
## 318.9747 319.8978 323.1783 
## 
## Error measures:
##                    ME     RMSE     MAE       MPE     MAPE      MASE       ACF1
## Training set 7.175981 33.63769 27.8431 0.4736071 15.57784 0.7021303 -0.2117522
## 
## Forecasts:
##    Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 31       207.1097 162.4882 251.7311 138.8670 275.3523
## 32       207.1097 161.8589 252.3604 137.9046 276.3147
## 33       207.1097 161.2382 252.9811 136.9554 277.2639
## 34       207.1097 160.6259 253.5935 136.0188 278.2005

autoplot(hardcover_fc) +
  autolayer(fitted(hardcover_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Hardcover Sales Projection over 4 days")

autoplot(paperback_fc) +
  autolayer(fitted(paperback_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Paperback Sales Projection over 4 days")

c.

Compute the RMSE values for the training data in each case.

My custom function below, calc_RMSE, calculates the RMSE of the forecast for the hardcover book sales in the variable HC_RMSE.

calc_RMSE <- function(fc){
  MSE <- mean(fc$residuals^2, na.rm=TRUE)
  return (sqrt(MSE))
  
}

The manually calculated RMSE for hardcover book sales training data using simple exponential smoothing is: 31.93101

HC_RMSE <-  calc_RMSE(hardcover_fc)
HC_RMSE

## [1] 31.93101

I confirmed that this is correct by taking the R generated RMSE. The same process is repeated for the paperback forecast.

The R RMSE for the hardcover training data is: 31.93101

accuracy(hardcover_fc)[1,2]

## [1] 31.93101

The manually calculated RMSE for paperback book sales training data using simple exponential smoothing is: 33.63769

Paperback_RMSE <- calc_RMSE(paperback_fc)
Paperback_RMSE 

## [1] 33.63769

The R RMSE for the paperback training data is: 33.63769

accuracy(paperback_fc)[1,2]

## [1] 33.63769

7.8.6

a.

Apply Holt’s linear method to the paperback and hardback series and compute four-day forecasts in each case.

The R function “holt” applies Holt’s linear forecast for the next four days.

hardcover_holt_fc <- holt(hardcover, h=4)
paperback_holt_fc <- holt(paperback, h=4)

The output from the summary from the hardcover forecast shows a very small alpha and beta as the two smoothing parameters.

summary(hardcover_holt_fc)

## 
## Forecast method: Holt's method
## 
## Model Information:
## Holt's method 
## 
## Call:
##  holt(y = hardcover, h = 4) 
## 
##   Smoothing parameters:
##     alpha = 1e-04 
##     beta  = 1e-04 
## 
##   Initial states:
##     l = 147.7935 
##     b = 3.303 
## 
##   sigma:  29.2106
## 
##      AIC     AICc      BIC 
## 310.2148 312.7148 317.2208 
## 
## Error measures:
##                      ME     RMSE      MAE       MPE    MAPE      MASE
## Training set -0.1357882 27.19358 23.15557 -2.114792 12.1626 0.6908555
##                     ACF1
## Training set -0.03245186
## 
## Forecasts:
##    Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 31       250.1739 212.7390 287.6087 192.9222 307.4256
## 32       253.4765 216.0416 290.9113 196.2248 310.7282
## 33       256.7791 219.3442 294.2140 199.5274 314.0308
## 34       260.0817 222.6468 297.5166 202.8300 317.3334

The output from the summary from the paperback forecast also shows a very small alpha and beta as the two smoothing parameters.

summary(paperback_holt_fc)

## 
## Forecast method: Holt's method
## 
## Model Information:
## Holt's method 
## 
## Call:
##  holt(y = paperback, h = 4) 
## 
##   Smoothing parameters:
##     alpha = 1e-04 
##     beta  = 1e-04 
## 
##   Initial states:
##     l = 170.699 
##     b = 1.2621 
## 
##   sigma:  33.4464
## 
##      AIC     AICc      BIC 
## 318.3396 320.8396 325.3456 
## 
## Error measures:
##                     ME     RMSE      MAE       MPE     MAPE      MASE
## Training set -3.717178 31.13692 26.18083 -5.508526 15.58354 0.6602122
##                    ACF1
## Training set -0.1750792
## 
## Forecasts:
##    Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 31       209.4668 166.6035 252.3301 143.9130 275.0205
## 32       210.7177 167.8544 253.5811 145.1640 276.2715
## 33       211.9687 169.1054 254.8320 146.4149 277.5225
## 34       213.2197 170.3564 256.0830 147.6659 278.7735

We can see clearly from the plots that the Holt method alone, with no other parameters set, captures the trend of the time series, but not much else. The forecast should not be relied upon since it does not capture seasonality in the series and that the upward trend continues on into infinity.

autoplot(hardcover_fc) +
  autolayer(fitted(hardcover_holt_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Holt Linear Hardcover Sales Projection over 4 days")

autoplot(paperback_fc) +
  autolayer(fitted(paperback_holt_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Holt Linear Paperback Sales Projection over 4 days")

b.

Using the custom function, calc_RMSE, the RMSE of the simple exponential smoothing forecast are compared to the Holt linear forecast. We see that by adding the linear trend to the forecast, even though the Beta parameter is quite small, did improve the RSME of the forecast.

hardcover_holt_RMSE <- calc_RMSE(hardcover_holt_fc)
paperback_holt_RMSE <- calc_RMSE(paperback_holt_fc)

ses_hardcover_rmse	Holt_linear_hardcover_rmse	ses_paperback_rmse	Holt_linear_paperback_rmse
31.93101	27.19358	33.63769	31.13692

c.

The plots below compare the two sets of forecasts, the SES and Holt Linear. Neither forecasts fit well to the series, but each does capture the trend. The Holt linear forecast smooths out the fluctuations and shows only the trend. As a result the 4 day forecasts for the Holt linear trend is higher than the simple exponential smoothing projection.

ses_hardcover <- autoplot(hardcover_fc) +
  autolayer(fitted(hardcover_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Hardcover Sales Projection over 4 days")

holt_linear_hardcover <- autoplot(hardcover_fc) +
  autolayer(fitted(hardcover_holt_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Holt Linear Hardcover Sales Projection over 4 days")

ses_paperback <- autoplot(paperback_fc) +
  autolayer(fitted(paperback_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Paperback Sales Projection over 4 days")

holt_linear_paperback <- autoplot(paperback_fc) +
  autolayer(fitted(paperback_holt_fc), series="Fitted") +
  ylab("Daily") + xlab("Sales") + ggtitle("Holt Linear Paperback Sales Projection over 4 days")

grid.arrange(ses_hardcover, holt_linear_hardcover, ses_paperback, holt_linear_paperback, ncol=2)

d.

Calculate a 95% prediction interval for the first forecast for each series, using the RMSE values and assuming normal errors. Compare your intervals with those produced using ses and holt.

Using the generate 95% CI custom function and comparing it to the R generated confidence intervals, we see that the R intervals are wider for each series.

Within each series, the Holt confidence intervals are tighter.

	lower 95 CI	upper 95 CI
ses hardcover	178.5848	300.5354
Holt hardcover	195.9640	304.3838
ses paperback	141.5964	272.6230
Holt paperback	147.8390	271.0945

	lower 95 CI	upper 95 CI
ses hardcover	174.7799	304.3403
Holt hardcover	192.9222	307.4256
ses paperback	138.8670	275.3523
Holt paperback	143.9130	275.0205

7.8.7

For this exercise use data set eggs, the price of a dozen eggs in the United States from 1900-1993. Experiment with the various options in the holt() 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 argument is doing to the forecasts.

[Hint: use h=100 when calling holt() so you can clearly see the differences between the various options when plotting the forecasts.]

Which model gives the best RMSE?

help(eggs)

Price of eggs Description Price of dozen eggs in US, 1900-1993, in constant dollars.

Usage eggs Format Time series data

Source Makridakis, Wheelwright and Hyndman (1998) Forecasting: methods and applications, John Wiley & Sons: New York. Chapter 9.

Examples plot(eggs)

head(eggs)

## Time Series:
## Start = 1900 
## End = 1905 
## Frequency = 1 
## [1] 276.79 315.42 314.87 321.25 314.54 317.92

tail(eggs)

## Time Series:
## Start = 1988 
## End = 1993 
## Frequency = 1 
## [1] 64.55 80.36 79.79 74.79 64.86 62.27

Below is the plot of the time series. We see a long term downward trend. We also see outlier spikes between ~1910 to 1920, and between ~1930 to 1940.

autoplot(eggs) +
  ylab("Price") + xlab("1990 to 1993") + ggtitle("Price of dozen eggs in US, 1900-1993, in constant dollars")

This plot show the Holt linear trend with no parameters. We see that the forecasted price goes to zero and below. Although Holt captures the trend, it’s not a reliable forecast since the trend turns negative.

eggs_linear <- holt(eggs, h=100)
autoplot(eggs_linear) + 
  autolayer(fitted(eggs_linear)) +
  ylab("Price") + xlab("1990 to 1993") + ggtitle("Holt Linear Trend with no Parameters")

Setting the damped parameter to TRUE improves on the forecast where it levels off before zero.

eggs_damped_linear <- holt(eggs, damped=TRUE, phi = 0.9, h=100)
autoplot(eggs_damped_linear) + 
  autolayer(eggs_damped_linear$fitted) +
  ylab("Price") + xlab("1990 to 1993") + ggtitle("Holt Damped")

The Holt linear trend with the lambda set to BoXCox lambda captures the trend correctly, but the trend does not go to zero or below.

eggs_box_cox <- holt(eggs, lambda = BoxCox.lambda(eggs), h=100)
autoplot(eggs_box_cox) + 
  autolayer(fitted(eggs_box_cox)) +
  ylab("Price") + xlab("1990 to 1993") + ggtitle("Holt Box Cox")

Setting the damped parameter to TRUE along with setting the lambda parameter to Box Cox for the Holt linear trend shows a trailing off of the trend well above zero.

eggs_box_cox_damped <- holt(eggs, damped=TRUE, lambda = BoxCox.lambda(eggs),  h=100)
autoplot(eggs_box_cox_damped) + 
  autolayer(fitted(eggs_box_cox_damped)) +
  ylab("Price") + xlab("1990 to 1993") + ggtitle("Holt Box Cox Damped")

The manually calculated RMSE values show that the Holt trend with damped parameter has the lowest RMSE. The R generated RMSE shows that the Holt linear trend with the BoX Cox parameter alone has the lowest RMSE.

Manually calculated RMSE comparisons:

	eggs_linear	eggs_damped_linear	eggs_box_cox	eggs_box_cox_damped
RMSE	31.93101	27.19358	33.63769	31.13692

R generated RMSE comparisons:

	eggs_linear	eggs_damped_linear	eggs_box_cox	eggs_box_cox_damped
RMSE	26.58219	26.62317	26.39376	26.53321

7.8.8

retaildata <- read_excel("retail.xlsx", skip=1)
liquor_sales <- ts(retaildata[,"A3349414R"],
  frequency=12, start=c(1982,4))

a.

Why is multiplicative seasonality necessary for this series?

We see from the plots below that multiplicative seasonality is necessary because variations are changing proportional to the level of the series.

liquor_sales  %>% seas(x11="") -> fit
autoplot(fit) +
  ggtitle("X11 decomposition of Australian Liquor Sales")

b.

Apply Holt-Winters’ multiplicative method to the data. Experiment with making the trend damped.

The Holt-Winters model add seasonality to the Holt linear model. Doing this captures the overall trend perfectly for this time series. Damping the trend does not add anything.

liquor_hw <- hw(liquor_sales,seasonal="multiplicative")
autoplot(liquor_sales) +
  autolayer(liquor_hw, series="HW multiplicative forecasts",
    PI=FALSE) +
  xlab("Year") +
  ylab("Sales (millions)") +
  ggtitle("Holt Winters Multiplicative Forecast in Australia") +
  guides(colour=guide_legend(title="Forecast"))

liquor_hw_damped <- hw(liquor_sales,seasonal="multiplicative")
autoplot(liquor_sales) +
  autolayer(liquor_hw_damped, series="HW multiplicative damped forecasts",
    PI=FALSE) +
  xlab("Year") +
  ylab("Sales (millions)") +
  ggtitle("Holt Winters Multiplicative Damped Forecast in Australia") +
  guides(colour=guide_legend(title="Forecast"))

c.

Compare the RMSE of the one-step forecasts from the two methods. Which do you prefer?

The manually generated RMSE values are far lower than the R generated values. In each instance the Holt Winters with the damped option have the lower RMSE although not by much. However, it is the better method.

liquor_hw_one_step <- hw(liquor_sales,seasonal="multiplicative", h=1)
liquor_hw_damped_one_step <- hw(liquor_sales,seasonal="multiplicative",damped = TRUE,  h=1)

The manually generated RMSE comparisons:

	liquor_hw_one_step	liquor_hw_damped_one_step
RMSE	0.0625853	0.0613796

The R generated RMSE comparisons:

	liquor_hw_one_step	liquor_hw_damped_one_step
RMSE	4.237051	4.147425

d.

Check that the residuals from the best method look like white noise.

The residuals from the best method, Holt Winters with the damped parameter set to TRUE, do not show that it is White noise based on the p value of the Ljung-Box test which is less than 0.05, as well as the ACF plot which shows correlations between the lags.

checkresiduals(liquor_hw_damped_one_step)

## 
##  Ljung-Box test
## 
## data:  Residuals from Damped Holt-Winters' multiplicative method
## Q* = 84.18, df = 7, p-value = 1.887e-15
## 
## Model df: 17.   Total lags used: 24

e.

Now find the test set RMSE, while training the model to the end of 2010. Can you beat the naive approach?

liquor_sales_train <- window(liquor_sales, end=c(2010,12))
liquor_sales_test <- window(liquor_sales, start=2011)
autoplot(liquor_sales) +
  autolayer(liquor_sales_train, series="Training") +
  autolayer(liquor_sales_test, series="Test")

The naive approach shows that the RMSE is for the test set is 22.956217

fc_liquor_sales_naive <- snaive(liquor_sales_train)
autoplot(fc_liquor_sales_naive)

accuracy(fc_liquor_sales_naive,liquor_sales_test)

##                     ME      RMSE       MAE      MPE      MAPE     MASE
## Training set  4.455255  8.699864  5.818619  6.15400  9.948117 1.000000
## Test set     19.170833 22.956217 19.520833 11.59039 11.813322 3.354891
##                   ACF1 Theil's U
## Training set 0.7261600        NA
## Test set     0.5801161 0.7479721

The Holt Winters approach shows that the RMSE is for the test set is 12.152124

liquor_sales_train_HW_Multi <- hw(liquor_sales_train, h = 36, seasonal = "multiplicative")

autoplot(liquor_sales_train_HW_Multi )

	ME	RMSE	MAE	MPE	MAPE	MASE	ACF1	Theil’s U
Training set	0.3379808	3.814454	2.454955	0.5733706	4.504932	0.4219137	0.0076021	NA
Test set	-4.8791198	12.152124	9.850873	-3.4254451	6.222885	1.6929918	0.6721037	0.4188652

The Holt Winters Damped approach shows that the RMSE is for the test set is 10.386188

liquor_sales_train_HW_Multi_damped <- hw(liquor_sales_train, h = 36, seasonal = "multiplicative", damped = TRUE)

autoplot(liquor_sales_train_HW_Multi_damped )

I was able to greatly improve on the RMSE of the naive model (22.956217) with the Holt Winters Damped model (10.386188), cutting down the error rate by more than half from the naive method.

	ME	RMSE	MAE	MPE	MAPE	MASE	ACF1	Theil’s U
Training set	0.4223714	3.717405	2.456063	0.3612467	4.543397	0.4221041	-0.0015053	NA
Test set	6.8376905	10.386188	7.953292	3.8069337	4.597764	1.3668694	0.4206074	0.3118086

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

I found that the RMSE for the forecast from the ETS, seasonally adjusted, Box Cox transformed data is (9.873769) which is lower than the RMSE from the Holt Winters Damped model (10.386188)

Below are the steps that I took to arrive at the answer:

Step One: Transform the liquor_sales_train using Box-Cox transformation. The best lambda is -0.07455321. The next two plots show the original training data and the Box Cox transformed data,

head(liquor_sales_train)

##       Apr  May  Jun  Jul  Aug  Sep
## 1982 17.3 18.1 18.1 18.9 19.0 18.4

autoplot(liquor_sales_train)

lambda <- BoxCox.lambda(liquor_sales_train)
lambda

## [1] -0.07455321

train_box_cox_transformed <- BoxCox(liquor_sales_train, lambda)
train_box_cox_transformed %>%  autoplot()

Step 2; Use STL decomposition on the Box Cox transformed data

train_box_cox_transformed2 <- ts(
  train_box_cox_transformed, 
  start = c(1982,1), 
  end = c(2010,12), 
  frequency = 12)

  train_box_cox_tran_stl <- stl(train_box_cox_transformed2,t.window=12, s.window='periodic', robust = TRUE)
  train_box_cox_tran_stl %>%  autoplot()

Step 3: ETS on the seasonally adjusted data

First, as per Hyndman’s blog:

As with other methods of decomposition, it is easy enough to remove the seasonal component to get the seasonally adjusted data.

adjust_df <- df - components[,‘season’] autoplot(df, series=“Data”) + autolayer(adjust_df, series=“Seasonally adjusted”)

https://robjhyndman.com/hyndsight/tslm-decomposition/

train_box_cox_trans_stl_seasonal_adj <- liquor_sales_train -  train_box_cox_tran_stl$time.series[, 'seasonal']
train_box_cox_trans_stl_seasonal_adj_ets <- ets(train_box_cox_trans_stl_seasonal_adj)
summary(train_box_cox_trans_stl_seasonal_adj_ets)

## ETS(M,A,M) 
## 
## Call:
##  ets(y = train_box_cox_trans_stl_seasonal_adj) 
## 
##   Smoothing parameters:
##     alpha = 0.5171 
##     beta  = 0.0025 
##     gamma = 0.2477 
## 
##   Initial states:
##     l = 20.2743 
##     b = 0.1729 
##     s = 1.0222 0.9479 1.0174 1.3945 1.0661 0.9953
##            0.9255 0.9497 0.9247 0.8984 0.9257 0.9327
## 
##   sigma:  0.0611
## 
##      AIC     AICc      BIC 
## 2740.354 2742.226 2805.694 
## 
## Training set error measures:
##                     ME     RMSE      MAE          MPE     MAPE     MASE
## Training set 0.2937121 3.825107 2.453263 -0.009510912 4.412601 0.421623
##                    ACF1
## Training set 0.06158662

	ME	RMSE	MAE	MPE	MAPE	MASE	ACF1	Theil’s U
Training set	0.2937121	3.825107	2.453263	-0.0095109	4.412601	0.421623	0.0615866	NA
Test set	2.1311634	9.873769	8.360848	0.8148709	4.979854	1.436913	0.5532785	0.3151752

train_box_cox_trans_stl_seasonal_adj_ets %>% forecast(h=36) %>% autoplot()