Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters
1. IVAE SLV 2009 - Diciembre 2022
1.1. Pronóstico para 2022 completo usando un modelo SARIMA
Importación de datos
convertir en objeto de Serie Temporal (ts)
## Qtr1 Qtr2 Qtr3 Qtr4
## 2009 4227.33 4436.39 4312.85 4625.05
## 2010 4375.00 4579.71 4498.45 4994.77
## 2011 4727.87 5163.09 5032.40 5360.43
## 2012 5176.52 5401.29 5225.97 5582.38
## 2013 5208.33 5550.62 5428.02 5803.99
## 2014 5434.08 5701.62 5589.75 5868.02
## 2015 5617.57 5870.22 5828.24 6122.21
## 2016 5695.11 6181.81 5962.28 6352.22
## 2017 5934.20 6271.01 6150.91 6623.07
## 2018 6159.12 6529.22 6480.41 6852.11
## 2019 6401.11 6713.27 6691.85 7074.90
## 2020 6462.08 5377.99 6118.58 6962.54
## 2021 6736.71 7203.08 7236.92 7866.431. Cálculo automático del modelo ARIMA
library(forecast)
library(TSstudio)
PIB_trim.arima.automatico<-auto.arima(y = PIB_trim.ts)
summary(PIB_trim.arima.automatico)## Series: PIB_trim.ts
## ARIMA(1,0,0)(2,1,0)[4] with drift
##
## Coefficients:
## ar1 sar1 sar2 drift
## 0.5874 -0.8878 -0.3696 53.5071
## s.e. 0.1216 0.1435 0.2155 9.6498
##
## sigma^2 = 62490: log likelihood = -332.83
## AIC=675.66 AICc=677.09 BIC=685.02
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 1.201683 229.9488 105.3604 -0.05167754 1.814617 0.3266285
## ACF1
## Training set 0.04734491.1. Pronóstico usando el modelo automático
library(forecast)
#Tabla
pronostico.arima.automatico<-forecast(object = PIB_trim.arima.automatico,h = 4,level = c(.95)) |> print()## Point Forecast Lo 95 Hi 95
## 2022 Q1 7190.394 6700.441 7680.347
## 2022 Q2 6698.502 6130.286 7266.718
## 2022 Q3 7020.770 6427.947 7613.593
## 2022 Q4 7636.601 7035.523 8237.680
2. Cálculo “manual” del modelo
Calcular las Diferencias:
library(kableExtra)
d<-ndiffs(PIB_trim.ts)
D<-nsdiffs(PIB_trim.ts)
ordenes_integracion<-c(d,D)
names(ordenes_integracion)<-c("d","D")
ordenes_integracion %>% kable(caption = "Ordenes de Integración") |> kable_material()| x | |
|---|---|
| d | 1 |
| D | 1 |
library(TSstudio)
PIB_trim.ts_diff<-PIB_trim.ts |>
diff(lag = 4,diffences=D) |>
diff(diffences=d)
PIB_trim.ts_diff |>
ts_cor(lag.max = 3*4)Lectura del Correlograma:
Parte regular p=0, q=0 Parte estacional P=2,Q=0
library(forecast)
PIB_trim.arima.manual<-Arima(y = PIB_trim.ts,order = c(0,1,0), seasonal = c(2,1,0)) |> print()## Series: PIB_trim.ts
## ARIMA(0,1,0)(2,1,0)[4]
##
## Coefficients:
## sar1 sar2
## -0.8752 -0.3604
## s.e. 0.1439 0.2166
##
## sigma^2 = 75961: log likelihood = -331.39
## AIC=668.77 AICc=669.33 BIC=674.322.1. Pronóstico usando el modelo manual
library(forecast)
#Tabla
pronostico.arima.manual<-forecast(object = PIB_trim.arima.manual,h = 4,level = c(0.95)) |> print()## Point Forecast Lo 95 Hi 95
## 2022 Q1 7360.215 6820.028 7900.403
## 2022 Q2 6972.775 6208.835 7736.716
## 2022 Q3 7350.562 6414.930 8286.194
## 2022 Q4 8001.659 6921.284 9082.034
2. Generar pronósticos para 2022 completo (diciembre) usando un modelo de Holt-Winters:
Carga serie
library(readxl)
library(forecast)
library(magrittr)
SERIEIVAE <- read_excel("C:/Users/USER/Desktop/Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters/tableIMAE.xlsx",col_types = c("skip", "numeric"), skip = 5)
SERIEIVAE.ts <- ts(data = SERIEIVAE, start = c(2009, 1), frequency = 12)
SERIEIVAE.ts %>% autoplot(main = "IVAE, El Salvador 2009-2022[Diciembre]", xlab = "Años/Meses", ylab = "Indice")
Estimación de modelo HW
## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
##
## Call:
## HoltWinters(x = Yt, seasonal = "multiplicative", optim.start = c(0.9, 0.9, 0.9))
##
## Smoothing parameters:
## alpha: 0.8470194
## beta : 0
## gamma: 1
##
## Coefficients:
## [,1]
## a 119.1738548
## b 0.1586990
## s1 0.9743058
## s2 0.9631373
## s3 1.0017599
## s4 0.9523792
## s5 1.0217635
## s6 1.0392864
## s7 0.9958466
## s8 1.0060704
## s9 0.9836083
## s10 0.9573668
## s11 1.0145028
## s12 1.0873190## Point Forecast Lo 95 Hi 95
## Jan 2023 116.2664 110.91216 121.6206
## Feb 2023 115.0865 108.10318 122.0698
## Mar 2023 119.8605 111.33639 128.3847
## Apr 2023 114.1033 104.77479 123.4317
## May 2023 122.5783 111.64638 133.5101
## Jun 2023 124.8454 112.84839 136.8423
## Jul 2023 119.7852 107.40058 132.1697
## Aug 2023 121.1746 107.87257 134.4766
## Sep 2023 118.6253 104.83919 132.4113
## Oct 2023 115.6124 101.43343 129.7914
## Nov 2023 122.6732 107.01024 138.3362
## Dec 2023 131.6507 91.42329 171.8781Aproximación de espacio de los estados
PronosticosHW2<-hw(y = Yt,
h = 12,
level = c(0.95),
seasonal = "multiplicative",
initial = "optimal")
PronosticosHW2## Point Forecast Lo 95 Hi 95
## Jan 2023 115.7968 110.15569 121.4378
## Feb 2023 113.2495 105.80792 120.6912
## Mar 2023 120.0137 110.44819 129.5793
## Apr 2023 114.7658 104.19767 125.3340
## May 2023 122.4225 109.75882 135.0861
## Jun 2023 123.4600 109.37441 137.5455
## Jul 2023 117.8473 103.20859 132.4861
## Aug 2023 121.0112 104.80173 137.2206
## Sep 2023 119.8008 102.62464 136.9770
## Oct 2023 117.2151 99.33405 135.0962
## Nov 2023 124.0884 104.04481 144.1320
## Dec 2023 132.4478 109.88642 155.0091
Identificando orden de integración
library(kableExtra)
d<-ndiffs(Yt)
D<-nsdiffs(Yt)
ordenes_integracion<-c(d,D)
names(ordenes_integracion)<-c("d","D")
ordenes_integracion %>% kable(caption = "Ordenes de Integración") %>% kable_material()| x | |
|---|---|
| d | 1 |
| D | 1 |
Grafico de la serie estacionaria
Valores para (p,q) y (P,Q)
En caso observando el PACF los Componentes autoregresivos el valor de p=0 mientras que P=1 mientras que los Componentes de media movil los valores de Q=1 y q=0
Estimando SARIMA
\(SARIMA(0,1,0)(1,1,1)[12]\)
library(ggthemes)
modelo_estimado <- Yt %>%
Arima(order = c(0, 1, 0),
seasonal = c(1, 1, 1))
summary(modelo_estimado)## Series: .
## ARIMA(0,1,0)(1,1,1)[12]
##
## Coefficients:
## sar1 sma1
## -0.0799 -0.7630
## s.e. 0.1119 0.0955
##
## sigma^2 = 6.507: log likelihood = -370.05
## AIC=746.1 AICc=746.26 BIC=755.23
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.04316329 2.434352 1.620592 0.02056879 1.590072 0.4138902
## ACF1
## Training set 0.02134069
Para el valor de MAPE es de 1.59% lo que quiere decir que de cada 100 valores pronosticados el error cometido sera de 1.59% y ante los residuos estamos en presencia de ruido blanco.
3. PIB trimestral precios corrientes SLV 2009-2021
3.1. Pronóstico para 2022 completo usando un modelo SARIMA
a) Usando la metodología paso a paso
Cargar serie
library(readxl)
library(forecast)
PIB_SLV<- read_excel("C:/Users/USER/Desktop/Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters/PIB a precios corrientes trimestral.xlsx", col_types = c("skip","numeric"),
skip = 5)
serie.PIB.ts<- ts(data = PIB_SLV, start = c(2009,1),frequency = 4)
Yt_1 <- serie.PIB.ts
serie.PIB.ts %>%
autoplot(main = "PIB trimestral precios corrientes SLV 2009-2021",
xlab = "años/ trimestres",
ylab = "Indice")
1.Gráfico de la Serie (Metodología Box-Jenkins)
2. Identificacion
2.2.Valores para (p,q) y (P,Q)
En caso observando el PACF los Componentes autoregresivos el valor de p=0 mientra que P=1 ya que el retardo 8 no es tan significativo mientras que los Componentes de media movil los valores de Q=1 y q=0
2.3 Estimacion del modelo
2.3.1 Usando forecast
library(forecast)
library(ggthemes)
modelo_estimado_PIB<-Yt_1 %>%
Arima(order = c(0, 1, 0),
seasonal = c(1, 1, 1))
summary(modelo_estimado_PIB)## Series: .
## ARIMA(0,1,0)(1,1,1)[4]
##
## Coefficients:
## sar1 sma1
## -0.4041 -0.4656
## s.e. 0.3684 0.3588
##
## sigma^2 = 77017: log likelihood = -331.62
## AIC=669.24 AICc=669.8 BIC=674.79
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 8.149971 258.1647 120.7587 0.01976588 2.03297 0.3743648 -0.1371202#grafico de las raices (estabilidad/convertibilidad)
modelo_estimado_PIB %>% autoplot(type="both")+theme_solarized()
Para el valor de MAPE es de 2.03% lo que quiere decir que de cada 100 valores pronosticados el error cometido sera de 2.03% y ante los residuos estamos en presencia de ruido blanco.
library(forecast)
#Estimar el modelo
modeloHW_PIB<-HoltWinters(x = Yt_1,
seasonal = "multiplicative",
optim.start = c(0.9,0.9,0.9))
modeloHW_PIB## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
##
## Call:
## HoltWinters(x = Yt_1, seasonal = "multiplicative", optim.start = c(0.9, 0.9, 0.9))
##
## Smoothing parameters:
## alpha: 0.8241646
## beta : 0
## gamma: 0
##
## Coefficients:
## [,1]
## a 7531.8419181
## b 48.1826250
## s1 0.9754034
## s2 1.0055201
## s3 0.9784087
## s4 1.0406678#Generar el pronostico:
pronosticosMW_PIB<-forecast(object = modeloHW_PIB,h=4, level=c(0.95))
pronosticosMW_PIB## Point Forecast Lo 95 Hi 95
## 2022 Q1 7393.582 6979.770 7807.393
## 2022 Q2 7670.315 7075.995 8264.636
## 2022 Q3 7510.646 6799.544 8221.749
## 2022 Q4 8038.714 7282.309 8795.1193. Verificacion de sobre ajuste/sub ajuste
library(tsibble)
library(feasts)
library(fable)
library(fabletools)
library(tidyr)
library(dplyr)
a<-Yt_1 %>% as_tsibble() %>%
model(arima_original=ARIMA(value ~ pdq(0, 1, 0) + PDQ(1, 1, 1)),
arima_010_011 = ARIMA(value ~ pdq(0, 1, 0) + PDQ(0, 1, 1)),
arima_010_110 = ARIMA(value ~ pdq(0, 1, 0) + PDQ(1, 1, 0)),
arima_automatico=ARIMA(value,ic="bic",stepwise = FALSE)
)
print(a)## # A mable: 1 x 4
## arima_original arima_010_011 arima_010_110
## <model> <model> <model>
## 1 <ARIMA(0,1,0)(1,1,1)[4]> <ARIMA(0,1,0)(0,1,1)[4]> <ARIMA(0,1,0)(1,1,0)[4]>
## # ℹ 1 more variable: arima_automatico <model>a %>% pivot_longer(everything(), names_to = "Model name",
values_to = "Orders") %>% glance() %>%
arrange(AICc) ->tabla
tabla## # A tibble: 4 × 9
## `Model name` .model sigma2 log_lik AIC AICc BIC ar_roots ma_roots
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <list> <list>
## 1 arima_010_011 Orders 75717. -332. 668. 669. 672. <cpl [0]> <cpl [4]>
## 2 arima_original Orders 77017. -332. 669. 670. 675. <cpl [4]> <cpl [4]>
## 3 arima_010_110 Orders 78488. -333. 670. 670. 673. <cpl [4]> <cpl [0]>
## 4 arima_automatico Orders 61662. -334. 675. 676. 683. <cpl [1]> <cpl [4]>Validacion Cruzada (Cross Validated)
library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)
#convertimos
Yt<-Yt_1 %>%
as_tsibble() %>%
rename(PIB=value)
#generador de la dta (roja)
data.cross.validation.PIB<-Yt %>%
as_tsibble() %>%
stretch_tsibble(.init = 20,.step = 1)
#generador de pronostico (verde) y simulacion
TSCV<-data.cross.validation.PIB %>%
model(ARIMA(PIB ~ pdq(0, 1, 0) + PDQ(1, 1, 1))) %>%
forecast(h=1) %>% accuracy(Yt)
print(TSCV) ## # 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 ARIMA(PIB ~ pdq(0, 1… Test 35.4 381. 198. 0.296 3.18 0.614 0.825 -0.00416El MAPE anterior era de 2.03% y ahora con la data desconocida el MAPE es de 3.18%, podemos concluir entonces que sigue teniendo un buen poder predictivo de nuestro modelo.
4. IPC general SLV 2009-2022 (mayo)
4.1. Pronóstico para 2022 completo usando un modelo SARIMA
a) Usando la metodología paso a paso
Cargar serie
library(readxl)
library(forecast)
IPC_SLV<- read_excel("C:/Users/USER/Desktop/Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters/IPC general.xlsx", col_types = c("skip","numeric"),
skip = 6)
serie.IPC.ts<- ts(data = IPC_SLV, start = c(2009,1),frequency = 12)
Yt_2 <- serie.IPC.ts
serie.IPC.ts %>%
autoplot(main = "IPC general SLV 2009-2022 (mayo)",
xlab = "años/ trimestres",
ylab = "Indice")
1. Gráfico de la Serie (Metodología Box-Jenkins)
2. Identificacion
2.2.Valores para (p,q) y (P,Q)
En caso observando el PACF los Componentes autoregresivos el valor de p=4 mientras que P=2 ya que el retardo 36 no es tan significativo mientras que los Componentes de media movil los valores de Q=1 y q=2
2.3 Estimacion del modelo
2.3.1 Usando forecast
library(forecast)
library(ggthemes)
modelo_estimado_IPC<-Yt_2 %>%
Arima(order = c(1, 1, 1),
seasonal = c(1, 0, 2))
summary(modelo_estimado_IPC)## Series: .
## ARIMA(1,1,1)(1,0,2)[12]
##
## Coefficients:
## ar1 ma1 sar1 sma1 sma2
## 0.9406 -0.7934 0.0571 0.0615 -0.0683
## s.e. 0.0859 0.1447 0.6621 0.6595 0.1106
##
## sigma^2 = 0.2251: log likelihood = -104.87
## AIC=221.74 AICc=222.29 BIC=240.15
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.05582688 0.4654425 0.3091819 0.04952956 0.2834401 0.1668998
## ACF1
## Training set 0.08746749#grafico de las raices (estabilidad/convertibilidad)
modelo_estimado_IPC %>% autoplot(type="both")+theme_solarized()
Para el valor de MAPE es de 0.28% lo que quiere decir que de cada 100 valores pronosticados el error cometido sera de 0.28% y ante los residuos estamos en presencia de ruido blanco.
library(forecast)
#Estimar el modelo
modeloHW_IPC<-HoltWinters(x = Yt_2,
seasonal = "multiplicative",
optim.start = c(0.9,0.9,0.9))
modeloHW_IPC## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
##
## Call:
## HoltWinters(x = Yt_2, seasonal = "multiplicative", optim.start = c(0.9, 0.9, 0.9))
##
## Smoothing parameters:
## alpha: 0.8419073
## beta : 0.0743462
## gamma: 1
##
## Coefficients:
## [,1]
## a 122.7484331
## b 0.4728397
## s1 1.0043101
## s2 1.0017791
## s3 0.9980880
## s4 0.9970592
## s5 0.9984875
## s6 0.9981783
## s7 0.9959355
## s8 0.9977953
## s9 1.0010455
## s10 1.0045211
## s11 1.0052312
## s12 1.0055526#Generar el pronostico:
pronosticosMW_IPC<-forecast(object = modeloHW_IPC,h=12, level=c(0.95))
pronosticosMW_IPC## Point Forecast Lo 95 Hi 95
## May 2022 123.7524 122.6054 124.8994
## Jun 2022 123.9142 122.3694 125.4590
## Jul 2022 123.9295 122.0324 125.8267
## Aug 2022 124.2733 122.0409 126.5056
## Sep 2022 124.9234 122.3615 127.4853
## Oct 2022 125.3567 122.4738 128.2396
## Nov 2022 125.5459 122.3495 128.7424
## Dec 2022 126.2522 122.7300 129.7744
## Jan 2023 127.1368 123.2817 130.9918
## Feb 2023 128.0532 123.8606 132.2457
## Mar 2023 128.6190 124.0961 133.1419
## Apr 2023 129.1356 123.8244 134.44683. Verificacion de sobre ajuste/sub ajuste
library(tsibble)
library(feasts)
library(fable)
library(fabletools)
library(tidyr)
library(dplyr)
#generando modelo automatico
a_2<-Yt_2 %>% as_tsibble() %>%
model(arima_original=ARIMA(value ~ pdq(1, 1, 1) + PDQ(1, 0, 2)),
arima_010_011 = ARIMA(value ~ pdq(1, 1, 1) + PDQ(0, 0, 2)),
arima_010_110 = ARIMA(value ~ pdq(1, 1, 1) + PDQ(1, 0, 1)),
arima_automatico=ARIMA(value,ic="bic",stepwise = FALSE)
)
print(a_2)## # A mable: 1 x 4
## arima_original arima_010_011
## <model> <model>
## 1 <ARIMA(1,1,1)(1,0,2)[12] w/ drift> <ARIMA(1,1,1)(0,0,2)[12] w/ drift>
## # ℹ 2 more variables: arima_010_110 <model>, arima_automatico <model>#ordenar
a_2 %>% pivot_longer(everything(), names_to = "Model name",
values_to = "Orders") %>% glance() %>%
arrange(AICc) ->tabla_2
tabla_2## # A tibble: 4 × 9
## `Model name` .model sigma2 log_lik AIC AICc BIC ar_roots ma_roots
## <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <list> <list>
## 1 arima_010_011 Orders 0.219 -102. 217. 218. 235. <cpl [1]> <cpl [25]>
## 2 arima_automatico Orders 0.218 -102. 217. 218. 235. <cpl [26]> <cpl [0]>
## 3 arima_010_110 Orders 0.219 -103. 217. 218. 236. <cpl [13]> <cpl [13]>
## 4 arima_original Orders 0.220 -102. 219. 220. 240. <cpl [13]> <cpl [25]>Validacion Cruzada (Cross Validated)
library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)
#convertimos
Yt<-Yt_2 %>% as_tsibble() %>% rename(IPC=value)
#generador de la dta (roja)
data.cross.validation.IPC<-Yt %>%
as_tsibble() %>%
stretch_tsibble(.init = 60,.step = 1)
#generador de pronostico (verde) y simulacion
TSCV_2<-data.cross.validation.IPC %>%
model(ARIMA(IPC ~ pdq(1, 1, 1) + PDQ(1, 0, 2))) %>%
forecast(h=1) %>% accuracy(Yt)
print(TSCV_2) ## # 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 ARIMA(IPC ~ pdq(1, 1… Test 0.0695 0.437 0.295 0.0579 0.261 0.159 0.160 0.0775El MAPE anterior era de 0.2834% y ahora con la data desconocida el MAPE es de 0.260%, podemos concluir entonces que sigue teniendo un buen poder predictivo de nuestro modelo.
5. Importación de Hidrocarburos SLV 2017-2022 (mayo)
5.1. Pronóstico para 2022 completo usando un modelo SARIMA
a) Usando la metodología paso a paso
Cargar serie
library(readxl)
library(forecast)
M_Hidrocarburos<- read_excel("C:/Users/USER/Desktop/Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters/hidrocarburos.xlsx", col_types = c("skip","numeric"),
skip = 7)
serie.M_HC.ts<- ts(data = M_Hidrocarburos, start = c(2017,1),frequency = 12)
Yt_3 <- serie.M_HC.ts
serie.M_HC.ts %>%
autoplot(main = "Importación de Hidrocarburos SLV 2017-2022 (mayo)",
xlab = "años/ trimestres",
ylab = "Indice")
1. Gráfico de la Serie (Metodología Box-Jenkins)
2. Identificacion
2.2.Valores para (p,q) y (P,Q)
En caso observando el PACF los Componentes autoregresivos el valor de p=2 ya que el retardo 4 no es significativo mientras que P=0 ya que el retardo 12 no es tan significativo mientras que los Componentes de media movil los valores de Q=1 y q=1 ya que el retardo 5 no es tan significativo.
2.3 Estimacion del modelo
2.3.1 Usando forecast
library(forecast)
library(ggthemes)
modelo_estimado_M_HC<-Yt_3 %>%
Arima(order = c(1, 0, 1),
seasonal = c(1, 0, 0))
summary(modelo_estimado_M_HC)## Series: .
## ARIMA(1,0,1)(1,0,0)[12] with non-zero mean
##
## Coefficients:
## ar1 ma1 sar1 mean
## 0.7675 -0.6480 -0.1675 215250.480
## s.e. 0.2240 0.2506 0.1519 5963.638
##
## sigma^2 = 1.374e+09: log likelihood = -750.32
## AIC=1510.64 AICc=1511.7 BIC=1521.36
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -369.5025 35877.48 28571.92 -3.554277 14.51947 0.6523367
## ACF1
## Training set -0.02814368#grafico de las raices (estabilidad/convertibilidad)
modelo_estimado_M_HC %>% autoplot(type="both")+theme_solarized()
Para el valor de MAPE es de 14.51% lo que quiere decir que de cada 100 valores pronosticados el error cometido sera de 14.51% y ante los residuos estamos en presencia de ruido blanco.
library(forecast)
#Estimar el modelo
modeloHW_M_HC<-HoltWinters(x = Yt_3,
seasonal = "multiplicative",
optim.start = c(0.9,0.9,0.9))
modeloHW_M_HC## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
##
## Call:
## HoltWinters(x = Yt_3, seasonal = "multiplicative", optim.start = c(0.9, 0.9, 0.9))
##
## Smoothing parameters:
## alpha: 0.1661945
## beta : 0
## gamma: 0.4642478
##
## Coefficients:
## [,1]
## a 2.573334e+05
## b 1.683628e+03
## s1 1.018236e+00
## s2 9.556789e-01
## s3 8.400515e-01
## s4 8.933522e-01
## s5 9.081175e-01
## s6 1.022234e+00
## s7 1.003639e+00
## s8 9.295311e-01
## s9 9.323749e-01
## s10 1.099022e+00
## s11 1.001407e+00
## s12 7.715030e-01#Generar el pronostico:
pronosticosMW_M_HC<-forecast(object = modeloHW_M_HC,h=12, level=c(0.95))
pronosticosMW_M_HC## Point Forecast Lo 95 Hi 95
## Apr 2022 263740.5 214576.5 312904.5
## May 2022 249146.2 198065.7 300226.6
## Jun 2022 220416.4 168128.3 272704.4
## Jul 2022 235905.7 180931.4 290880.0
## Aug 2022 241333.7 184169.9 298497.4
## Sep 2022 273381.2 211970.8 334791.7
## Oct 2022 270098.1 207387.3 332809.0
## Nov 2022 251719.3 189236.0 314202.6
## Dec 2022 254059.2 189773.1 318345.3
## Jan 2023 301318.4 229877.6 372759.2
## Feb 2023 276241.4 206751.9 345730.9
## Mar 2023 214120.6 174612.5 253628.73. Verificacion de sobre ajuste/sub ajuste
library(tsibble)
library(feasts)
library(fable)
library(fabletools)
library(tidyr)
library(dplyr)
#generando modelo automatico
a_3<-Yt_3 %>% as_tsibble() %>%
model(arima_original=ARIMA(value ~ pdq(1, 0, 1) + PDQ(1, 0, 0)),
arima_010_011 = ARIMA(value ~ pdq(1, 1, 1) + PDQ(1, 0, 0)),
arima_010_110 = ARIMA(value ~ pdq(1, 0, 1) + PDQ(1, 1, 0)),
arima_automatico=ARIMA(value,ic="bic",stepwise = FALSE) #modelo automatico
)
print(a_3)## # A mable: 1 x 4
## arima_original arima_010_011
## <model> <model>
## 1 <ARIMA(1,0,1)(1,0,0)[12] w/ mean> <ARIMA(1,1,1)(1,0,0)[12]>
## # ℹ 2 more variables: arima_010_110 <model>, arima_automatico <model>Validacion Cruzada (Cross Validated)
library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)
#convertimos
Yt<-Yt_3 %>%
as_tsibble() %>%
rename(M_HC= value)
#generador de la dta (roja)
data.cross.validation.M_HC<-Yt %>%
as_tsibble() %>%
stretch_tsibble(.init = 60,.step = 1)
#generador de pronostico (verde) y simulacion
TSCV_3<-data.cross.validation.M_HC %>%
model(ARIMA(M_HC ~ pdq(1, 0, 1) + PDQ(1, 0, 1))) %>%
forecast(h=1) %>% accuracy(Yt)
print(TSCV_3) ## # 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 ARIMA(M_HC ~ pdq(1… Test 17780. 65601. 64297. 0.727 27.9 1.47 1.18 -0.0682El MAPE anterior era de 14.51% y ahora con la data desconocida el MAPE es de 27.87%, podemos concluir entonces que no tiene un buen poder predictivo de nuestro modelo.