DATA624 - HW 5
Farhana Akther
2025-03-09
library(fpp3)
library(tidyverse)Exercise 8.1
Consider the the number of pigs slaughtered in Victoria, available in the ‘aus_livestock’ dataset. let’s extract the data for Victoria and plot.
pigs_vic <- aus_livestock |>
filter(State == "Victoria", Animal == "Pigs")autoplot(pigs_vic, Count)+
labs(y = "Count"
, x = "Date"
, title = "Pigs Slaughtered, Victoria") 
a.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.
Answer:
fit_pigs <- pigs_vic |>
model(ETS(Count ~ error("A") + trend("N") + season("N")))
report(fit_pigs) ## 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.07Let’s generate Forecasts for the Next 4 Months
fc <- fit_pigs |>
forecast(h = "4 months")Plot Forecast:
pigs_new <- pigs_vic |>
filter(yearmonth(Month) >= yearmonth("2015 Jan")) # last ~5 years
autoplot(pigs_new, Count) +
autolayer(fc, color = "blue") + # Add forecast line
labs(title = "Forecast for Pig Slaughter in Victoria",
y = "Number of pigs slaughtered") +
theme_minimal()
Analysis: We used the ETS(A,N,N) model (Simple Exponential Smoothing) to estimate the pig slaughter data. The optimal values found were \(\alpha = 0.3221\) and \(\ell_0 = 100646.6\). The four-month forecast shows a stable trend, with the widening 80% and 95% confidence intervals indicating increasing uncertainty over time.
b. 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.
Answer:
int <- hilo(fc)
print(int$`95%`[1])## <hilo[1]>
## [1] [76854.79, 113518.3]95resids <- augment(fit_pigs)
s <- sd(resids$.resid)
print(paste("Standard Deviation: ",s))## [1] "Standard Deviation: 9344.66579258967"upper <- fc$.mean[1]+1.96*s
lower <- fc$.mean[1]-1.96*s
print(paste("Upper limit: ", upper))## [1] "Upper limit: 113502.102384467"print(paste("Lower limit: ", lower))## [1] "Lower limit: 76871.0124775157"We computed the 95% prediction interval for the first forecast using \(\hat{y} \pm 1.96s\), where \(\hat{y}\) is the first forecasted value and \(s\) is the standard deviation of residuals.
- Manual Calculation: [76871.01, 113502.1]
- R’s Built-in Interval: [76854.79, 113518.3]
As we can see the the two intervals are nearly identical.This confirms that our manual calculation matches R’s prediction interval.
Exercise 8.5
Data set global_economy contains the annual Exports from many countries. Select one country to analyse.
glimpse(global_economy)## Rows: 15,150
## Columns: 9
## Key: Country [263]
## $ Country <fct> "Afghanistan", "Afghanistan", "Afghanistan", "Afghanistan",…
## $ Code <fct> AFG, AFG, AFG, AFG, AFG, AFG, AFG, AFG, AFG, AFG, AFG, AFG,…
## $ Year <dbl> 1960, 1961, 1962, 1963, 1964, 1965, 1966, 1967, 1968, 1969,…
## $ GDP <dbl> 537777811, 548888896, 546666678, 751111191, 800000044, 1006…
## $ Growth <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
## $ CPI <dbl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…
## $ Imports <dbl> 7.024793, 8.097166, 9.349593, 16.863910, 18.055555, 21.4128…
## $ Exports <dbl> 4.132233, 4.453443, 4.878051, 9.171601, 8.888893, 11.258279…
## $ Population <dbl> 8996351, 9166764, 9345868, 9533954, 9731361, 9938414, 10152…#unique(global_economy$Country)First, let’s select a country to analyze using filter and plot the Exports data over time. let’s choose Australia.
global_economy |>
filter(Country == "Australia") |>
select(Year, GDP, Exports) |>
head()## # A tsibble: 6 x 3 [1Y]
## Year GDP Exports
## <dbl> <dbl> <dbl>
## 1 1960 18573188487. 13.0
## 2 1961 19648336880. 12.4
## 3 1962 19888005376. 13.9
## 4 1963 21501847911. 13.0
## 5 1964 23758539590. 14.9
## 6 1965 25931235301. 13.2a. Plot the Exports series and discuss the main features of the data.
exports_aus <- global_economy |>
filter(Country == "Australia") |>
select(Year, Exports)
autoplot(exports_aus, Exports) +
labs(title = "Annual Exports of Australia",
y = "Exports (% of GDP)", # Updated label
x = "Year") +
theme_minimal()
The plot shows that Australia’s exports have increased over time, with fluctuations becoming more prominent after 2000. A sharp spike and drop around the early 2000s may indicate an economic event. No clear seasonal pattern is visible.
b.Use an ETS(A,N,N) model to forecast the series, and plot the forecasts.
fit_aus <- exports_aus |>
model(ETS(Exports ~ error("A") + trend("N") + season("N")))
report(fit_aus) # estimated alpha and other model parameters## Series: Exports
## Model: ETS(A,N,N)
## Smoothing parameters:
## alpha = 0.5659948
##
## Initial states:
## l[0]
## 12.98943
##
## sigma^2: 1.3621
##
## AIC AICc BIC
## 257.3943 257.8387 263.5756fc_aus <- fit_aus |> forecast(h = "10 years")
autoplot(exports_aus, Exports) +
autolayer(fc_aus, level = c(80, 95)) +
geom_line(aes(y = .fitted), col="#D55E00",
data = augment(fit_aus)) +
labs(title = "ETS(A,N,N) Forecast for Australia's Exports",
y = "Exports (% of GDP)",
x = "Year") +
theme_minimal()
We applied the ETS(A,N,N) model (Simple Exponential Smoothing) to forecast Australia’s Exports (% of GDP) for the next 10 years. The optimal smoothing parameter found was \(\alpha = 0.566\), meaning recent data has a moderate influence on the forecast.The initial level was estimated as \(\ell_0 = 12.99\%\). of GDP. The forecast maintains a stable trend, with 80% and 95% confidence intervals widening over time, indicating increasing uncertainty in future predictions. The orange line represents the fitted values, showing how well the model tracks past data.
c. Compute the RMSE values for the training data.
rmse_ANN <- accuracy(fit_aus) |> select(RMSE)
print(rmse_ANN)## # A tibble: 1 × 1
## RMSE
## <dbl>
## 1 1.15The Root Mean Squared Error (RMSE) for the ETS(A,N,N) model is 1.147. This indicates that the average deviation of the model’s predictions from actual export values in the training data.
d. 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 Both Models Together
mix_aus <- exports_aus |>
model(
ANN = ETS(Exports ~ error("A") + trend("N") + season("N")), # ETS(A,N,N)
AAN = ETS(Exports ~ error("A") + trend("A") + season("N")) # ETS(A,A,N)
)
# Compare RMSE for both models
accuracy(mix_aus)## # A tibble: 2 × 10
## .model .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 ANN Training 0.232 1.15 0.914 1.09 5.41 0.928 0.928 0.0125
## 2 AAN Training -0.000000746 1.12 0.893 -0.387 5.39 0.907 0.904 0.109We compared the ETS(A,N,N) (Simple Exponential Smoothing) and ETS(A,A,N) (Additive Trend) models. The ETS(A,A,N) model had a slightly lower RMSE (1.117 vs. 1.147), suggesting it provides a slightly better fit. However, the difference is small, and the trend parameter is very close to zero, meaning the trend effect is weak. While ETS(A,A,N) is more flexible, ETS(A,N,N) is a simpler model and still performs well if no strong trend is expected
e. Compare the forecasts from both methods. Which do you think is best?
Let’s generate Forecasts for ETS(A,A,N) and plot both forecasts
mix_aus |>
forecast(h = "10 years") |>
autoplot(exports_aus, level = NULL) +
labs(y = "Exports (% of GDP)",
title = "Comparison of ETS(A,N,N) vs. ETS(A,A,N) Forecasts",
x = "Year") +
theme_minimal()
We compared the forecasts from ETS(A,N,N) (Simple Exponential Smoothing) and ETS(A,A,N) (Additive Trend). The ETS(A,N,N) model produces a stable forecast, while the ETS(A,A,N) model captures an upward trend. Since ETS(A,A,N) has a slightly lower RMSE (1.117 vs. 1.147), it is the better model for forecasting, assuming exports continue to rise. However, ETS(A,N,N) is simpler and still performs well if no strong trend is expected.
f. 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.
R’s Built-in Prediction Interval
mix_aus |>
forecast(h=1) |>
mutate(interval = hilo(Exports, 95))## # A fable: 2 x 5 [1Y]
## # Key: .model [2]
## .model Year Exports
## <chr> <dbl> <dist>
## 1 ANN 2018 N(21, 1.4)
## 2 AAN 2018 N(21, 1.3)
## # ℹ 2 more variables: .mean <dbl>, interval <hilo>Manual Calculation
resids <- augment(fit_aus)
s <- sd(resids$.resid) mix_aus |>
forecast(h=1)|>
mutate(low= .mean -1.96*s,
high= .mean +1.96*s)## # A fable: 2 x 6 [1Y]
## # Key: .model [2]
## .model Year Exports
## <chr> <dbl> <dist>
## 1 ANN 2018 N(21, 1.4)
## 2 AAN 2018 N(21, 1.3)
## # ℹ 3 more variables: .mean <dbl>, low <dbl>, high <dbl>We computed the 95% prediction interval for the first forecast using two methods:
- R’s built-in hilo() function, which gave intervals [18.32, 22.89] (ANN) and [18.57, 23.10] (AAN).
- Manual calculation using \(\hat{y}\pm 1.96s\), which produced similar results [18.39, 22.83] (ANN) and [18.62, 23.06] (AAN).
Since both methods gave nearly identical results, our manual calculation is correct and matches R’s built-in prediction intervals.
Overall, both models performed well, but ETS(A,A,N) was slightly better as it captured a small upward trend, while ETS(A,N,N) provided a simpler, stable forecast. If exports are expected to grow, ETS(A,A,N) is a better choice, but if no strong trend is anticipated, ETS(A,N,N) is sufficient. Our manually calculated 95% prediction intervals closely matched R’s, confirming the accuracy of our calculations.
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.]
Answer:
First, let’s extract and Plot the Chinese GDP Data. We’ll first filter the global_economy dataset to isolate China’s GDP and plot the time series to understand its trend and behavior.
head(global_economy)## # A tsibble: 6 x 9 [1Y]
## # Key: Country [1]
## Country Code Year GDP Growth CPI Imports Exports Population
## <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Afghanistan AFG 1960 537777811. NA NA 7.02 4.13 8996351
## 2 Afghanistan AFG 1961 548888896. NA NA 8.10 4.45 9166764
## 3 Afghanistan AFG 1962 546666678. NA NA 9.35 4.88 9345868
## 4 Afghanistan AFG 1963 751111191. NA NA 16.9 9.17 9533954
## 5 Afghanistan AFG 1964 800000044. NA NA 18.1 8.89 9731361
## 6 Afghanistan AFG 1965 1006666638. NA NA 21.4 11.3 9938414china_gdp <- global_economy |>
filter(Country == "China") |>
select(Year, GDP)
autoplot(china_gdp, GDP) +
labs(title = "China's GDP Over Time",
y = "GDP (in USD)",
x = "Year") +
theme_minimal()
Since our data (China’s GDP) shows a strong upward trend and potentially non-linear growth, performing STL decomposition can help us Understand the key components and Decide if a damped trend or Box-Cox transformation is necessary by observing the trend’s behavior.
# STL decomposition
dcmp <- china_gdp |>
model(stl = STL(GDP ~ trend()))
components(dcmp) |>
autoplot() +
labs(title = "STL Decomposition of China's GDP")
From the plot, we can observe that the trend shows a strong exponential growth starting around the early 2000s. The residuals show some fluctuations, especially in the later years, but no clear recurring pattern.This confirms the absence of seasonality. The presence of larger residuals post-2000 could indicate increased variability.
Since the trend appears exponential we can: - Fit an ETS model with a damped trend — this will control the sharp upward projection.
Try a Box-Cox transformation — this can stabilize the variance and improve forecast reliability.
Use a larger forecast horizon (e.g., 20+ years) to clearly visualize the differences.
Estimate the Box-Cox Transformation Parameter
lambda <- china_gdp %>%
features(GDP, features = guerrero) %>%
pull(lambda_guerrero)
lambda## [1] -0.03446284Fit ETS Models
fit <- china_gdp |>
model(
"Holt's Method" = ETS(GDP ~ error('A') + trend('A') + season('N')),
"Damped Holt's Method" = ETS(GDP ~ error('A') + trend('Ad') + season('N')),
"Box-Cox" = ETS(box_cox(GDP, lambda) ~ error('A') + trend('Ad') + season('N')),
"Log-Transformed Model" = ETS(log(GDP) ~ error('A') + trend('Ad') + season('N'))
)fit|>
forecast(h = 20) %>%
autoplot(china_gdp, level = NULL) +
labs(title = "China GDP",
y = "GDP") +
guides(colour = guide_legend(title = "Forecasts"))
fit |> accuracy()## # A tibble: 4 × 10
## .model .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Holt's Method Trai… 2.36e10 1.90e11 9.59e10 1.41 7.62 0.442 0.453 0.00905
## 2 Damped Holt's… Trai… 2.95e10 1.90e11 9.49e10 1.62 7.62 0.438 0.454 -0.00187
## 3 Box-Cox Trai… 6.35e 9 1.96e11 1.02e11 1.77 7.26 0.472 0.468 0.135
## 4 Log-Transform… Trai… -3.05e10 2.43e11 1.12e11 0.794 7.24 0.518 0.579 0.476We explored four ETS models to forecast China’s GDP: Holt’s Method, Damped Holt’s Method, a Box-Cox transformation model, and a Log-Transformed model. Although the ‘guerrero()’ method suggested a log transformation \(\lambda ≈ 0)\), the resulting model performed poorly, producing the highest RMSE and heavily correlated residuals. The Box-Cox model also showed aggressive growth. Among the four models, the Damped Holt’s Method provided the most stable and realistic forecast. It shows the lowest MAE and well-controlled residuals which, making it the best choice for this data.
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?
glimpse(aus_production)## Rows: 218
## Columns: 7
## $ Quarter <qtr> 1956 Q1, 1956 Q2, 1956 Q3, 1956 Q4, 1957 Q1, 1957 Q2, 1957…
## $ Beer <dbl> 284, 213, 227, 308, 262, 228, 236, 320, 272, 233, 237, 313…
## $ Tobacco <dbl> 5225, 5178, 5297, 5681, 5577, 5651, 5317, 6152, 5758, 5641…
## $ Bricks <dbl> 189, 204, 208, 197, 187, 214, 227, 222, 199, 229, 249, 234…
## $ Cement <dbl> 465, 532, 561, 570, 529, 604, 603, 582, 554, 620, 646, 637…
## $ Electricity <dbl> 3923, 4436, 4806, 4418, 4339, 4811, 5259, 4735, 4608, 5196…
## $ Gas <dbl> 5, 6, 7, 6, 5, 7, 7, 6, 5, 7, 8, 6, 5, 7, 8, 6, 6, 8, 8, 7…# Filter Gas data from 1980 Q1 onward
gas_data <- aus_production |>
filter(Quarter >= yearquarter("1980 Q1")) |>
select(Quarter, Gas)
# Plot the filtered data
autoplot(gas_data, Gas) +
labs(title = "Australian Gas Production (from 1980 Q1)",
y = "Gas Production (in PJ)",
x = "Year") +
theme_minimal()
# STL decomposition on filtered data
dcmp2 <- gas_data |>
model(stl = STL(Gas ~ trend(window = 13)))
components(dcmp2) |>
autoplot() +
labs(title = "STL Decomposition of Australian Gas Production",
y = "Gas Production (in PJ)")
The clear increase in seasonal variation supports the need for a multiplicative seasonality model in the ETS framework and the trend’s gradual flattening suggests that adding a damped trend may improve long-term forecast stability.
fit_gas_data <- gas_data |>
model(
"Additive" = ETS(Gas ~ error('A') + trend('A') + season('A')),
"Multiplicative" = ETS(Gas ~ error('M') + trend('A') + season('M')),
"Damped Holt's Method" = ETS(Gas ~ error('A') + trend('Ad') + season('M'))
)
fit_gas_data <- fit_gas_data |>
forecast(h = 10)fit_gas_data |>
autoplot(gas_data, level = NULL) +
labs(title = "Forecasts for Australian Gas Production (from 1980 Q1)",
y = "Gas Production (in PJ)",
x = "Year") +
guides(colour = guide_legend(title = "Forecasts")) +
facet_grid(.model ~ .)
The multiplicative model was necessary because the seasonal variation increased over time, which additive models could not capture effectively. The Damped Holt’s Method provided a realistic alternative if the growth trend is expected to slow down.
Exercise 8.8
Recall your retail time series data (from Exercise 7 in Section 2.10).
set.seed(12345678) #from Exercise 7 in Section 2.10
myseries <- aus_retail |>
filter(`Series ID` == sample(aus_retail$`Series ID`, 1))a. Why is multiplicative seasonality necessary for this series?
Multiplicative seasonality is necessary because the seasonal fluctuations increase in size as turnover grows. This pattern suggests that seasonal effects are proportional to the overall trend, this is why a multiplicative model is more appropriate.
#STL decomposition
dcmp3 <- myseries|>
model(stl = STL(Turnover))
components(dcmp3) |> autoplot()
b.Apply Holt-Winters’ multiplicative method to the data. Experiment with making the trend damped
fit_retail <- myseries |>
model(
"HW Mult" = ETS(Turnover ~ error("M") + trend("A") + season("M")),
"Damped HW" = ETS(Turnover ~ error("M") + trend("Ad") + season("M")),
"HW Damped 95" = ETS(Turnover ~ error("A") + trend("Ad", phi = 0.95) + season("M"))
)
# Generate forecasts
fc_retail <- fit_retail |>
forecast(h = 15)# Plot the forecasts
fc_retail |>
autoplot(myseries, level = NULL) +
labs(title = "Retail Turnover Forecasts Using Holt-Winters’ Methods (Faceted)",
y = "Turnover (in AUD)",
x = "Year") +
guides(colour = guide_legend(title = "Models")) +
facet_grid(.model ~ .) +
theme_minimal()
c. Compare the RMSE of the one-step forecasts from the two methods. Which do you prefer?
accuracy(fit_retail) |>
select(.model, RMSE)## # A tibble: 3 × 2
## .model RMSE
## <chr> <dbl>
## 1 HW Mult 0.613
## 2 Damped HW 0.616
## 3 HW Damped 95 0.603Based on the RMSE values from the table, the HW Damped 95 model has the lowest RMSE (0.6027). This makes it the most accurate among the three models for one-step forecasts. Thus, HW Damped 95 is the preferred model in this case.
d. Check that the residuals from the best method look like white noise.
myseries |>
model(best_model = ETS(Turnover ~ error("A") + trend("Ad", phi = 0.95) + season("M"))) |>
gg_tsresiduals() +
labs(title = "Residuals, Best Model (HW Damped 95)")
The residuals from the HW Damped 95 model resemble white noise. All three elements — random fluctuations around zero, minimal autocorrelation within the significance limits, and a roughly normal distribution — align with the characteristics of white noise. This confirms that the model effectively captures the data’s structure.
e. 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?
Let’s filter the data so that the training set includes data up to the end of 2010 and then fit the models using the training data.
Split Data
myseries_train <- myseries |>
filter(year(Month) <= 2010)
myseries_test <- myseries |>
filter(year(Month) >= 2011)Fit the Models on the Training Data
fit_train <- myseries_train |>
model(
"HW Mult" = ETS(Turnover ~ error("M") + trend("A") + season("M")),
"Damped HW" = ETS(Turnover ~ error("M") + trend("Ad") + season("M")),
"HW Damped 95" = ETS(Turnover ~ error("A") + trend("Ad", phi = 0.95) + season("M")),
"SNAIVE" = SNAIVE(Turnover)
)Forecast for the Test Period
fc_test <- fit_train |>
forecast(new_data = myseries_test)Calculate RMSE on the test set for each model.
test_accuracy <- fc_test |>
accuracy(myseries_test) |>
select(.model, RMSE)Combine the RMSE values for comparison.
test_accuracy <- fc_test |>
accuracy(myseries_test) |>
select(.model, RMSE)
test_accuracy## # A tibble: 4 × 2
## .model RMSE
## <chr> <dbl>
## 1 Damped HW 1.15
## 2 HW Damped 95 1.38
## 3 HW Mult 1.78
## 4 SNAIVE 1.55Yes, the Damped Holt-Winters model successfully outperformed the Seasonal Naïve approach as it achieves a lower RMSE 1.151 compared to 1.552. This improved accuracy shows that the Damped HW model better captures the data’s trend and seasonal patterns which, makes it the preferred choice for forecasting in this case.
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?
Answer:
# Step 1: Box-Cox Transformation
lambda2 <- myseries_train |>
features(Turnover, features = guerrero) |>
pull(lambda_guerrero)
myseries_train_bx <- myseries_train |>
mutate(Turnover = box_cox(Turnover, lambda2))# Step 2: STL Decomposition
dcmp4 <- myseries_train_bx |>
model(stl = STL(Turnover))
components(dcmp4) |> autoplot() 
As we can observe above, the Box-Cox transformation made the seasonality a lot more constant as well as the variability in the remainder.
# Step 3: Fit ETS Model on Box-Cox Transformed Data
fit_train2 <- myseries_train_bx |>
model(
"Damped Holt's Method" = ETS(Turnover ~ error('A') + trend('Ad') + season('M'))
)
# Step 4: Forecast
fc_train2 <- fit_train2 |>
forecast(h = 15)
# Step 5: Visualize the Forecast
fc_train2 |>
autoplot(myseries_train_bx, level = NULL) +
labs(title = "Household Goods Turnover (STL + ETS Forecast)") +
guides(colour = guide_legend(title = "Forecasts")) +
theme_minimal()
accuracy(fit_train) |>
select(.model, RMSE)## # A tibble: 4 × 2
## .model RMSE
## <chr> <dbl>
## 1 HW Mult 0.518
## 2 Damped HW 0.519
## 3 HW Damped 95 0.518
## 4 SNAIVE 1.21accuracy(fit_train2)|>
select(.model, RMSE)## # A tibble: 1 × 2
## .model RMSE
## <chr> <dbl>
## 1 Damped Holt's Method 0.109The STL decomposition combined with an ETS model on the Box-Cox transformed series resulted in a significant improvement in forecast accuracy. The RMSE for the STL + ETS model was 0.1091, which is much lower than the RMSE values from the previous models, including the Damped HW, HW Damped 95, HW Mult, and the Seasonal Naïve model. This improvement indicates that the STL decomposition effectively stabilized the seasonal pattern and reduced variability in the data, enabling the ETS model to better capture the underlying trend and seasonal structure. As a result, the STL + ETS model outperformed all other models and is the preferred choice for forecasting in this scenario.