library(ggplot2)
library(tidyverse)
library(fpp3)Nuevos modelos suavizacion exponencial
Librerías necesarias
Ejemplo de los modelos anteriores (MEAN, NAIVE, SEASONAL NAIVE Y DRIFT)
Ej: producción de cerveza
recent_production <- aus_production %>% filter(year(Quarter) >= 1992)
beer_train <- recent_production %>% filter(year(Quarter) <= 2007)
beer_fit <- beer_train %>%
model(
Mean = MEAN(Beer),
`Naïve` = NAIVE(Beer),
`Seasonal naïve` = SNAIVE(Beer),
Drift = RW(Beer ~ drift())
)
beer_fc <- beer_fit %>%
forecast(h = 10)
beer_fc %>%
autoplot(filter(aus_production, year(Quarter) >= 1992), level = NULL) +
xlab("Year") + ylab("Megalitres") +
ggtitle("Forecasts for quarterly beer production") +
guides(colour=guide_legend(title="Forecast"))
Métricas de error en el entrenamiento:
beer_accu_train <- accuracy(beer_fit) |>
arrange(MAE)
beer_accu_train# 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 Seasonal naïve Training -2.13e+ 0 16.8 14.3 -0.554 3.31 1 1 -0.288
2 Mean Training 0 43.6 35.2 -0.937 7.89 2.46 2.60 -0.109
3 Naïve Training 4.76e- 1 65.3 54.7 -0.916 12.2 3.83 3.89 -0.241
4 Drift Training -7.22e-15 65.3 54.8 -1.03 12.2 3.83 3.89 -0.241Nos dice que el modelo que mejor se ajusta en el entrenamiento, de acuerdo al MAE es Seasonal naïve.
Métricas de error del pronóstico
beer_accu_fc <- beer_fc |>
accuracy(recent_production) |>
arrange(MAE)
beer_accu_fc# 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 Seasonal naïve Test 5.2 14.3 13.4 1.15 3.17 0.937 0.853 0.132
2 Mean Test -13.8 38.4 34.8 -3.97 8.28 2.44 2.29 -0.0691
3 Naïve Test -51.4 62.7 57.4 -13.0 14.2 4.01 3.74 -0.0691
4 Drift Test -54.0 64.9 58.9 -13.6 14.6 4.12 3.87 -0.0741Esto nos indica que el modelo que tiene un menor MAE en el pronóstico es Seasonal naïve.
Gráficas de los datos de entrenamiento en niveles y en logaritmos:
beer_train |>
autoplot(Beer)
beer_train |>
autoplot(log(Beer))
Hacemos un: Modelo de descomposición
#Con logaritmo
dcmp <- beer_train |>
model(
STL(log(Beer), robust = TRUE)
)
dcmp |>
components() |>
autoplot()
dcmp |>
components() |>
ggplot(aes(x = Quarter, y = season_adjust)) +
geom_line()
#Sin logaritmo
beer_train |>
model(
STL(Beer, robust = TRUE)
) |>
components() |>
autoplot()
Ajuste de modelos
beer_fit <- beer_train |>
model(
snaive = SNAIVE(Beer),
ets_ANA = ETS(Beer ~ error("A") + trend("N") + season("A")),
ets_AAdA = ETS(Beer ~ error("A") + trend("Ad") + season("A")),
ets_MAdM = ETS(Beer ~ error("M") + trend("Ad") + season("M")),
ets_ANA_l = ETS(log(Beer) ~ error("A") + trend("N") + season("A")),
ets_AAdA_l = ETS(log(Beer) ~ error("A") + trend("Ad") + season("A")),
ets_MAdM_l = ETS(log(Beer) ~ error("M") + trend("Ad") + season("M")),
stl_ets_A = decomposition_model(
STL(log(Beer), robust = TRUE),
ETS(season_year ~ error("A") + trend("N") + season("A")),
ETS(season_adjust ~ error("A") + trend("Ad") + season("N"))
),
stl_ets_M = decomposition_model(
STL(log(Beer), robust = TRUE),
ETS(season_year ~ error("M") + trend("N") + season("M")),
ETS(season_adjust ~ error("A") + trend("Ad") + season("N"))
)
)
beer_fit# A mable: 1 x 9
snaive ets_ANA ets_AAdA ets_MAdM ets_ANA_l ets_AAdA_l
<model> <model> <model> <model> <model> <model>
1 <SNAIVE> <ETS(A,N,A)> <ETS(A,Ad,A)> <ETS(M,Ad,M)> <ETS(A,N,A)> <ETS(A,Ad,A)>
# ℹ 3 more variables: ets_MAdM_l <model>, stl_ets_A <model>, stl_ets_M <model>Vamos a utilizar el MAPE como métrica de error:
accuracy(beer_fit) |>
arrange(MAPE)# A tibble: 9 × 10
.model .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
<chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 stl_ets_A Training 0.293 11.5 8.75 0.0389 2.01 0.612 0.684 -0.123
2 ets_MAdM_l Training -0.163 12.3 9.35 -0.0756 2.15 0.654 0.731 -0.265
3 ets_MAdM Training -0.276 12.4 9.45 -0.0936 2.18 0.661 0.736 -0.262
4 ets_AAdA_l Training -0.375 12.3 9.44 -0.128 2.18 0.660 0.731 -0.239
5 ets_AAdA Training -0.533 12.4 9.62 -0.161 2.22 0.673 0.741 -0.243
6 ets_ANA_l Training -2.71 12.8 9.59 -0.652 2.22 0.671 0.760 -0.275
7 ets_ANA Training -2.63 12.8 9.73 -0.636 2.25 0.681 0.763 -0.260
8 snaive Training -2.13 16.8 14.3 -0.554 3.31 1 1 -0.288
9 stl_ets_M Training 62.8 130. 75.8 12.7 15.9 5.30 7.72 -0.0633beer_fit |>
augment() |>
features(.innov, ljung_box, lag = 8)# A tibble: 9 × 3
.model lb_stat lb_pvalue
<chr> <dbl> <dbl>
1 ets_AAdA 6.71 0.568
2 ets_AAdA_l 12.3 0.140
3 ets_ANA 8.88 0.353
4 ets_ANA_l 10.7 0.222
5 ets_MAdM 9.80 0.279
6 ets_MAdM_l 11.0 0.200
7 snaive 32.9 0.0000630
8 stl_ets_A 7.72 0.461
9 stl_ets_M 32.2 0.0000850Hacemos los pronósticos para los siguientes 10 periodos (dos años y medio):
beer_fc <- beer_fit |>
forecast(h = "2 years 6 months")
beer_fc# A fable: 90 x 4 [1Q]
# Key: .model [9]
.model Quarter Beer .mean
<chr> <qtr> <dist> <dbl>
1 snaive 2008 Q1 N(427, 282) 427
2 snaive 2008 Q2 N(383, 282) 383
3 snaive 2008 Q3 N(394, 282) 394
4 snaive 2008 Q4 N(473, 282) 473
5 snaive 2009 Q1 N(427, 563) 427
6 snaive 2009 Q2 N(383, 563) 383
7 snaive 2009 Q3 N(394, 563) 394
8 snaive 2009 Q4 N(473, 563) 473
9 snaive 2010 Q1 N(427, 845) 427
10 snaive 2010 Q2 N(383, 845) 383
# ℹ 80 more rowsbeer_fc |>
autoplot(recent_production, level = NULL, size = 1)
beer_fc |>
autoplot(recent_production |> filter_index("2005 Q1" ~ .), level = NULL, size = 1)
beer_fc |>
autoplot(recent_production |> filter_index("2005 Q1" ~ .), size = 1) +
facet_wrap(~ .model, ncol = 3) +
theme(legend.position = "none")
beer_fc |>
filter(.model != "stl_ets_M") |>
autoplot(recent_production |> filter_index("2005 Q1" ~ .), size = 1) +
facet_wrap(~ .model, ncol = 3) +
theme(legend.position = "none")
Errores del test
beer_fc |>
accuracy(recent_production) |>
arrange(MAPE)# A tibble: 9 × 10
.model .type ME RMSE MAE MPE MAPE MASE RMSSE ACF1
<chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 ets_MAdM Test 2.57 9.68 8.80 0.568 2.12 0.616 0.577 0.220
2 stl_ets_A Test 7.20 11.4 8.99 1.71 2.15 0.629 0.682 0.0805
3 ets_ANA_l Test 0.771 9.88 9.10 0.128 2.20 0.636 0.589 0.200
4 ets_AAdA Test 2.12 10.0 9.29 0.453 2.23 0.650 0.597 0.184
5 ets_MAdM_l Test 2.14 10.1 9.38 0.449 2.25 0.656 0.599 0.185
6 ets_ANA Test 0.857 10.3 9.64 0.143 2.32 0.674 0.616 0.175
7 ets_AAdA_l Test 1.15 10.7 10.0 0.200 2.40 0.702 0.636 0.152
8 snaive Test 5.2 14.3 13.4 1.15 3.17 0.937 0.853 0.132
9 stl_ets_M Test 15.0 43.8 25.2 2.83 5.50 1.76 2.61 -0.154 Pronóstico hacia el futuro, con el mejor modelo:
beer_fut <- recent_production |>
model(
ets_MAdM = ETS(Beer ~ error("M") + trend("Ad") + season("M"))
) |>
forecast(h = "2 years 6 months")
beer_fut# A fable: 10 x 4 [1Q]
# Key: .model [1]
.model Quarter Beer .mean
<chr> <qtr> <dist> <dbl>
1 ets_MAdM 2010 Q3 N(402, 133) 402.
2 ets_MAdM 2010 Q4 N(489, 198) 489.
3 ets_MAdM 2011 Q1 N(418, 145) 418.
4 ets_MAdM 2011 Q2 N(385, 122) 385.
5 ets_MAdM 2011 Q3 N(401, 133) 401.
6 ets_MAdM 2011 Q4 N(489, 198) 489.
7 ets_MAdM 2012 Q1 N(418, 145) 418.
8 ets_MAdM 2012 Q2 N(384, 123) 384.
9 ets_MAdM 2012 Q3 N(401, 133) 401.
10 ets_MAdM 2012 Q4 N(488, 198) 488.beer_fut |>
autoplot(recent_production)