Ejercicio de Pronóstico Modelos SARIMA y Holt-Winters

PIB trimestral precios corrientes SLV 2009-2025

Con las series adjuntas en los enlaces, se solicita:

1. Generar pronósticos para 2026 completo (diciembre) usando un modelo SARIMA:

library(readxl)
PIB_trimestral <- read_excel("PIB a precios corrientes  trimestral.xls", 
    col_types = c("skip", "numeric"), skip = 5)
names(PIB_trimestral)<-c("PIB_trimestral")

# convertir en objeto de Serie Temporal (ts)
PIB_trimestral.ts<-ts(data = PIB_trimestral$PIB_trimestral,start = c(2009,1),frequency = 4) |> print()
##         Qtr1    Qtr2    Qtr3    Qtr4
## 2009      NA 4227.33 4436.39 4312.85
## 2010 4625.05 4375.00 4579.71 4498.45
## 2011 4994.77 4727.87 5163.09 5032.40
## 2012 5360.43 5176.52 5401.29 5225.97
## 2013 5582.38 5208.33 5550.62 5428.02
## 2014 5803.99 5434.08 5701.62 5589.75
## 2015 5868.02 5617.57 5870.22 5828.24
## 2016 6122.20 5695.10 6181.80 5962.28
## 2017 6352.24 5934.25 6271.06 6150.90
## 2018 6622.97 6158.89 6529.02 6480.43
## 2019 6852.52 6402.12 6714.13 6691.76
## 2020 7073.12 6463.49 5377.69 6119.98
## 2021 6960.03 6748.64 7208.52 7239.80
## 2022 7846.18 7556.38 7957.48 7930.34
## 2023 8425.91 7913.51 8422.35 8372.77
## 2024 8856.80 8301.40 8728.83 8539.41
## 2025 9310.10 8657.02 9134.58 9089.93
## 2026 9826.57

a) Usando la metodología paso a paso mostrada en Clase. (manual)

library(forecast)
library(kableExtra)
d<-ndiffs(PIB_trimestral.ts)
D<-nsdiffs(PIB_trimestral.ts)
ordenes_integracion<-c(d,D)
names(ordenes_integracion)<-c("d","D")
ordenes_integracion %>% kable(caption = "Ordenes de Integración")  |>  kable_material()
Ordenes de Integración
x
d 1
D 1
library(forecast)
PIB_trimestral.arima.manual<-Arima(y = PIB_trimestral.ts,order = c(0,1,0), seasonal = c(2,1,0)) |> print()
## Series: PIB_trimestral.ts 
## ARIMA(0,1,0)(2,1,0)[4] 
## 
## Coefficients:
##          sar1     sar2
##       -0.5723  -0.2602
## s.e.   0.1181   0.1152
## 
## sigma^2 = 76982:  log likelihood = -443.54
## AIC=893.07   AICc=893.48   BIC=899.5

Pronostico manual

library(forecast)

#Tabla
pronostico.arima.manual<-forecast(object = PIB_trimestral.arima.manual,h = 4,level = c(0.95)) |> print()
##         Point Forecast    Lo 95     Hi 95
## 2026 Q2       9240.580 8696.777  9784.383
## 2026 Q3       9710.637 8941.583 10479.691
## 2026 Q4       9619.529 8677.635 10561.424
## 2027 Q1      10301.057 9213.451 11388.664
#Gráfica
PIB_trimestral.arima.manual |> forecast(h = 4,level = c(.95)) |> autoplot()

b) Usando la generación automática disponible en forecast.

library(forecast)
library(TSstudio)
PIB_trimestral.arima.automatico<-auto.arima(y = PIB_trimestral.ts)
summary(PIB_trimestral.arima.automatico)
## Series: PIB_trimestral.ts 
## ARIMA(1,0,0)(0,1,1)[4] with drift 
## 
## Coefficients:
##          ar1     sma1    drift
##       0.8176  -0.7784  72.0990
## s.e.  0.0870   0.1507  10.9996
## 
## sigma^2 = 60944:  log likelihood = -443.83
## AIC=895.66   AICc=896.34   BIC=904.29
## 
## Training set error measures:
##                    ME     RMSE      MAE         MPE     MAPE     MASE
## Training set -2.09545 233.8164 124.3587 -0.09474471 1.923808 0.346164
##                     ACF1
## Training set -0.07426543
library(forecast)

#Tabla
pronostico.arima.automatico<-forecast(object = PIB_trimestral.arima.automatico,h = 4,level = c(.95)) |> print()
##         Point Forecast    Lo 95     Hi 95
## 2026 Q2       9302.107 8818.233  9785.981
## 2026 Q3       9540.216 8915.204 10165.228
## 2026 Q4       9440.967 8737.221 10144.714
## 2027 Q1       9977.914 9226.123 10729.705
#Gráfica
PIB_trimestral.arima.automatico |> forecast(h = 4,level = c(.95)) |> autoplot()

2. Generar pronósticos para 2026 completo (diciembre) usando un modelo de Holt-Winters:

Usando Stats y forecast

library(forecast)

Yt <- PIB_trimestral.ts #Serie original

Usando Forecast (Aproximación por Espacios de los Estados ETS )

library(forecast)

Pronosticos_HW2<-hw(y = Yt,
                   h = 12,
                   level = c(0.95),
                   seasonal = "multiplicative",
                   initial = "optimal")
Pronosticos_HW2
##         Point Forecast    Lo 95     Hi 95
## 2026 Q2       9160.425 8517.738  9803.111
## 2026 Q3       9508.059 8632.189 10383.929
## 2026 Q4       9390.505 8361.202 10419.808
## 2027 Q1      10027.395 8778.623 11276.168
## 2027 Q2       9364.693 8074.831 10654.555
## 2027 Q3       9718.904 8264.307 11173.501
## 2027 Q4       9597.595 8056.082 11139.107
## 2028 Q1      10247.318 8497.376 11997.260
## 2028 Q2       9568.961 7843.929 11293.992
## 2028 Q3       9929.748 8050.851 11808.646
## 2028 Q4       9804.684 7866.433 11742.936
## 2029 Q1      10467.241 8313.754 12620.727

Gráfico de la serie original y del pronóstico

Pronosticos_HW2 %>% autoplot()

3. Realice una validación cruzada para los modelos estimados. Interprete el MAPE

library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)

Yt<-Yt %>% as_tsibble() %>% rename(IVAE=value)

data.cross.validation<-Yt %>% 
  as_tsibble() %>% 
  stretch_tsibble(.init = 40,.step = 1)

TSCV1<-data.cross.validation %>% 
model(ARIMA(IVAE ~ pdq(0, 1, 1) + PDQ(1, 0, 0))) %>% 
forecast(h=1) %>% accuracy(Yt)
print(TSCV1)
## # 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(IVAE ~ pdq(0, 1, … Test   135.  503.  325.  1.55  4.58 0.899  1.06 0.141

Indice de precios al consumidor general SLV 2009-2026 (mayo)

library(readxl)

IPCg_mensual <- read_excel("IPC_SLV_2009_2026.xls", 
    col_types = c("skip", "numeric"), skip = 5)
## Warning: Expecting numeric in B7 / R7C2: got 'IPC general'
names(IPCg_mensual)<-c("IPCg_mensual")

# convertir en objeto de Serie Temporal (ts)
IPCg_mensual.ts<-ts(data = IPCg_mensual$IPCg_mensual,start = c(2009,1),frequency = 12) |> print()
##         Jan    Feb    Mar    Apr    May    Jun    Jul    Aug    Sep    Oct
## 2009     NA  99.65  99.53 100.00 100.08 100.36 100.40  99.99  99.79  99.73
## 2010 100.00 100.44 100.57 100.89 100.71 100.50 100.96 100.98 100.81 101.11
## 2011 102.13 102.77 102.96 103.63 106.71 107.23 107.29 107.57 107.69 107.40
## 2012 107.29 107.65 108.00 108.16 108.83 108.53 107.94 107.61 107.80 108.24
## 2013 108.13 108.59 109.05 109.54 108.85 108.69 108.91 108.78 108.87 109.06
## 2014 108.98 109.51 109.72 109.99 109.47 109.72 110.15 110.77 111.04 110.91
## 2015 109.50 108.69 108.56 109.10 109.11 109.33 109.24 109.16 108.82 108.41
## 2016 110.61 110.67 110.37 110.32 110.05 110.13 110.24 110.12 109.85 109.51
## 2017 109.58 110.39 110.69 110.92 111.00 111.19 111.26 111.24 111.10 111.22
## 2018 111.81 111.96 112.05 111.93 111.97 112.11 112.26 112.42 112.71 112.76
## 2019 112.30 112.24 112.44 112.69 112.87 113.01 112.85 112.56 112.16 111.99
## 2020 112.29 112.15 112.00 112.09 111.69 111.94 112.59 112.49 111.82 111.56
## 2021 112.20 112.49 113.19 114.08 114.81 114.84 115.51 116.36 116.63 117.10
## 2022 119.06 119.79 120.73 121.71 122.32 123.43 124.47 124.99 125.56 125.87
## 2023 127.77 128.21 128.97 129.08 128.97 128.88 129.18 129.17 129.44 129.67
## 2024 129.34 129.76 130.00 130.08 130.45 130.71 131.09 131.47 130.96 130.42
## 2025 129.71 130.16 130.08 130.26 130.30 130.44 130.86 131.29 130.81 130.89
## 2026 130.89 131.00 131.60 132.18 133.12 133.74                            
##         Nov    Dec
## 2009  99.03 100.41
## 2010 101.78 102.24
## 2011 107.32 107.48
## 2012 108.33 108.18
## 2013 108.91 109.00
## 2014 110.97 110.40
## 2015 110.76 110.69
## 2016 109.79 109.78
## 2017 111.36 111.62
## 2018 113.02 112.82
## 2019 112.04 112.17
## 2020 111.81 111.98
## 2021 117.95 118.92
## 2022 126.76 127.63
## 2023 130.14 130.32
## 2024 130.04 129.93
## 2025 131.26 131.41
## 2026

1. Generar pronósticos para 2026 completo (diciembre) usando un modelo SARIMA:

a) Usando la metodología paso a paso mostrada en Clase. (manual)

Calculo de las diferencias

library(forecast)
library(kableExtra)

d<-ndiffs(IPCg_mensual.ts)
D<-nsdiffs(IPCg_mensual.ts)

ordenes_integracion<-c(d,D)
names(ordenes_integracion)<-c("d","D")
ordenes_integracion %>% kable(caption = "Ordenes de Integración")  |>  kable_material()
Ordenes de Integración
x
d 1
D 0

Encontramos p,q y P,Q

library(TSstudio)

IPCg_mensual.ts_diff<-IPCg_mensual.ts |>  
  diff(lag = 12,diffences=D) |> 
  diff(diffences=d)
library(forecast)
IPCg_mensual.arima.manual<-Arima(y = IPCg_mensual.ts,order = c(1,1,0), seasonal = c(1,0,0)) |> print()
## Series: IPCg_mensual.ts 
## ARIMA(1,1,0)(1,0,0)[12] 
## 
## Coefficients:
##          ar1    sar1
##       0.3387  0.2506
## s.e.  0.0653  0.0705
## 
## sigma^2 = 0.2048:  log likelihood = -129.65
## AIC=265.3   AICc=265.41   BIC=275.31

Pronostico manual

library(forecast)

#Tabla
pronostico.arima.manual<-forecast(object = IPCg_mensual.arima.manual,h = 12,level = c(0.95)) |> print()
##          Point Forecast    Lo 95    Hi 95
## Jul 2026       134.0434 133.1564 134.9303
## Aug 2026       134.2182 132.7361 135.7003
## Sep 2026       134.1206 132.1563 136.0850
## Oct 2026       134.1484 131.7797 136.5170
## Nov 2026       134.2437 131.5246 136.9628
## Dec 2026       134.2822 131.2512 137.3132
## Jan 2027       134.1522 130.8380 137.4664
## Feb 2027       134.1798 130.6046 137.7551
## Mar 2027       134.3302 130.5117 138.1487
## Apr 2027       134.4756 130.4284 138.5228
## May 2027       134.7112 130.4476 138.9748
## Jun 2027       134.8666 130.3970 139.3362
#Gráfica
IPCg_mensual.arima.manual |> forecast(h = 12,level = c(.95)) |> autoplot()

b) Usando la generación automática disponible en forecast.

library(forecast)
library(TSstudio)
IPCg_mensual.arima.automatico<-auto.arima(y = IPCg_mensual.ts)
summary(IPCg_mensual.arima.automatico)
## Series: IPCg_mensual.ts 
## ARIMA(0,1,1)(1,0,0)[12] with drift 
## 
## Coefficients:
##          ma1    sar1   drift
##       0.2663  0.1957  0.1639
## s.e.  0.0650  0.0719  0.0477
## 
## sigma^2 = 0.1996:  log likelihood = -126.26
## AIC=260.52   AICc=260.72   BIC=273.87
## 
## Training set error measures:
##                      ME      RMSE       MAE          MPE     MAPE     MASE
## Training set 0.00190199 0.4424161 0.2993996 0.0005398963 0.266636 0.138356
##                    ACF1
## Training set 0.01670461
library(forecast)

#Tabla
pronostico.arima.automatico<-forecast(object = IPCg_mensual.arima.automatico,h = 12,level = c(.95)) |> print()
##          Point Forecast    Lo 95    Hi 95
## Jul 2026       134.1894 133.3138 135.0649
## Aug 2026       134.4054 132.9927 135.8181
## Sep 2026       134.4433 132.6475 136.2390
## Oct 2026       134.5908 132.4803 136.7012
## Nov 2026       134.7950 132.4111 137.1790
## Dec 2026       134.9562 132.3271 137.5853
## Jan 2027       134.9863 132.1330 137.8396
## Feb 2027       135.1397 132.0786 138.2008
## Mar 2027       135.3890 132.1333 138.6447
## Apr 2027       135.6343 132.1950 139.0736
## May 2027       135.9502 132.3366 139.5637
## Jun 2027       136.2034 132.4236 139.9832
#Gráfica
IPCg_mensual.arima.automatico |> forecast(h = 12,level = c(.95)) |> autoplot()

2. Generar pronósticos para 2026 completo (diciembre) usando un modelo de Holt-Winters:

Usando Stats y forecast

library(forecast)

Yt <- IPCg_mensual.ts #Serie original

Usando Forecast (Aproximación por Espacios de los Estados ETS )

library(forecast)
#Generar el pronóstico:
Pronosticos_HW2<-hw(y = Yt,
                   h = 7,
                   level = c(0.95),
                   seasonal = "multiplicative",
                   initial = "optimal")
Pronosticos_HW2
##          Point Forecast    Lo 95    Hi 95
## Jul 2026       133.9023 132.3632 135.4414
## Aug 2026       134.2595 132.4936 136.0254
## Sep 2026       134.2632 132.1999 136.3265
## Oct 2026       134.4506 132.0221 136.8791
## Nov 2026       135.0704 132.2120 137.9288
## Dec 2026       135.4902 132.1564 138.8239
## Jan 2027       135.4626 131.6230 139.3022

Gráfico de la serie original y del pronóstico.

Pronosticos_HW2 %>% autoplot()

3. Realice una validación cruzada para los modelos estimados. Interprete el MAPE

library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)

Yt<-Yt %>% as_tsibble() %>% rename(IVAE=value)

data.cross.validation<-Yt %>% 
  as_tsibble() %>% 
  stretch_tsibble(.init = 40,.step = 1)

TSCV1<-data.cross.validation %>% 
model(ARIMA(IVAE ~ pdq(0, 1, 1) + PDQ(1, 0, 0))) %>% 
forecast(h=1) %>% accuracy(Yt)
print(TSCV1)
## # 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(IVAE ~ pdq(0… Test  -0.0379 0.411 0.275 -0.0354 0.248 0.127 0.127 0.0906

Interpretación del MAPE: 0.2610455 = en promedio, los pronósticos del modelo tienen un error absoluto del 0.26% en relación con los valores reales de la serie temporal.

Importación de Hidrocarburos SLV 2017-2026 (mayo)

Carga de datos y conversión a serie temporal

library(readxl)

ImpHidro_m<- read_excel("Importación de hidrocarburos.xls", 
    col_types = c("skip", "numeric"), skip = 5)
## Warning: Expecting numeric in B7 / R7C2: got 'Hidrocarburos valor importado'
names(ImpHidro_m)<-c("ImpHidro_m")

# convertir en objeto de Serie Temporal (ts)
ImpHidro_m.ts<-ts(data = ImpHidro_m$ImpHidro_m,start = c(2017,1),frequency = 12) |> print()
##         Jan    Feb    Mar    Apr    May    Jun    Jul    Aug    Sep    Oct
## 2017     NA 101.71  93.97 149.53  96.28 112.07 128.51 128.47 167.53 139.48
## 2018 140.00 152.20 169.09 136.08 120.74 112.82  97.10 131.95 143.68 162.51
## 2019 127.23 106.52 129.90 127.78 114.39 125.77 100.54 111.09  65.87  29.79
## 2020  57.25  71.73  74.83  91.58  81.36 124.39 108.04 145.58 113.48 138.09
## 2021 146.58 146.68 174.96 172.97 200.09 150.78 180.91 265.99 275.49 188.09
## 2022 200.63 224.95 220.88 247.63 147.73 223.36 202.00 203.09 185.03 207.53
## 2023 229.80 240.41 205.84 186.69 190.69 176.43 189.96 215.51 227.91 227.68
## 2024 211.79 147.99 198.96 162.72 179.72 209.13 158.97 192.69 200.53 189.13
## 2025 175.30 177.24 206.81 180.84 207.55 179.03 157.20 230.00 250.91       
##         Nov    Dec
## 2017 143.59 116.46
## 2018 107.81 132.04
## 2019  53.36  75.73
## 2020 164.20 170.94
## 2021 264.20 305.72
## 2022 225.29 188.08
## 2023 203.01 177.83
## 2024 177.54 262.91
## 2025

1. Generar pronósticos para 2026 completo (diciembre) usando un modelo SARIMA:

a) Usando la metodología paso a paso mostrada en Clase. (manual)

library(forecast)
library(kableExtra)

d<-ndiffs(ImpHidro_m.ts)
D<-nsdiffs(ImpHidro_m.ts)

ordenes_integracion<-c(d,D)
names(ordenes_integracion)<-c("d","D")
ordenes_integracion %>% kable(caption = "Ordenes de Integración")  |>  kable_material()
Ordenes de Integración
x
d 1
D 0
library(TSstudio)
ImpHidro_m.ts_diff<-ImpHidro_m.ts |>  
  diff(lag = 12,diffences=D) |> 
  diff(diffences=d)
library(forecast)
ImpHidro_m.arima.manual<-Arima(y = ImpHidro_m.ts,order = c(1,0,1), seasonal = c(0,0,1)) |> print()
## Series: ImpHidro_m.ts 
## ARIMA(1,0,1)(0,0,1)[12] with non-zero mean 
## 
## Coefficients:
##          ar1      ma1     sma1      mean
##       0.9668  -0.5690  -0.0555  163.8082
## s.e.  0.0264   0.0819   0.1085   29.8276
## 
## sigma^2 = 971.2:  log likelihood = -503.95
## AIC=1017.91   AICc=1018.52   BIC=1031.13

pronostico manual

library(forecast)

#Tabla
pronostico.arima.manual<-forecast(object = ImpHidro_m.arima.manual,h = 12,level = c(0.95)) |> print()
##          Point Forecast    Lo 95    Hi 95
## Oct 2025       216.9221 155.8425 278.0016
## Nov 2025       215.8233 150.0875 281.5592
## Dec 2025       209.4132 139.6050 279.2215
## Jan 2026       212.4701 139.0590 285.8811
## Feb 2026       210.9171 134.2912 287.5429
## Mar 2026       207.5388 128.0253 287.0524
## Apr 2026       207.5623 125.4412 289.6835
## May 2026       204.5452 120.0593 289.0312
## Jun 2026       204.6067 117.9684 291.2449
## Jul 2026       204.5856 115.9827 293.1885
## Aug 2026       199.0962 108.6954 289.4970
## Sep 2026       196.6141 104.5644 288.6638
#Gráfica
ImpHidro_m.arima.manual |> forecast(h = 12,level = c(.95)) |> autoplot()

b) Usando la generación automática disponible en forecast.

library(forecast)
library(TSstudio)
ImpHidro_m.arima.automatico<-auto.arima(y = ImpHidro_m.ts)
summary(ImpHidro_m.arima.automatico)
## Series: ImpHidro_m.ts 
## ARIMA(0,1,1) 
## 
## Coefficients:
##           ma1
##       -0.5891
## s.e.   0.0747
## 
## sigma^2 = 960.8:  log likelihood = -499.55
## AIC=1003.1   AICc=1003.22   BIC=1008.37
## 
## Training set error measures:
##                    ME    RMSE      MAE      MPE     MAPE      MASE        ACF1
## Training set 2.738386 30.6967 24.01432 2.191632 16.04681 0.5207591 -0.03305953
ImpHidro_m.arima.automatico |> forecast(h = 12,level = c(.95)) |> autoplot()

2. Generar pronósticos para 2026 completo (Diciembre) usando un modelo de Holt-Winters:

Usando Stats y forecast

library(forecast)

Yt <- ImpHidro_m.ts #Serie original

Usando Forecast (Aproximación por Espacios de los Estados ETS )

library(forecast)

#Generar el pronóstico:
Pronosticos_HW2<-hw(y = Yt,
                   h = 7,
                   level = c(0.95),
                   seasonal = "multiplicative",
                   initial = "optimal")
Pronosticos_HW2
##          Point Forecast    Lo 95    Hi 95
## Oct 2025       209.5064 129.2995 289.7134
## Nov 2025       221.8555 129.8986 313.8124
## Dec 2025       236.7065 131.6619 341.7511
## Jan 2026       211.4210 111.8120 311.0299
## Feb 2026       205.9395 103.6071 308.2719
## Mar 2026       221.3433 105.9548 336.7318
## Apr 2026       218.3097  99.4309 337.1885

Gráfico de la serie original y del pronóstico.

Pronosticos_HW2 %>% autoplot()

3. Realice una validación cruzada para los modelos estimados. Interprete el MAPE

library(forecast)
library(dplyr)
library(tsibble)
library(fable)
library(fabletools)

Yt<-Yt %>% as_tsibble() %>% rename(IVAE=value)

data.cross.validation<-Yt %>% 
  as_tsibble() %>% 
  stretch_tsibble(.init = 40,.step = 1)

TSCV1<-data.cross.validation %>% 
model(ARIMA(IVAE ~ pdq(0, 0, 0) + PDQ(1, 0, 0))) %>% 
forecast(h=1) %>% accuracy(Yt)
## Warning: The future dataset is incomplete, incomplete out-of-sample data will be treated as missing. 
## 1 observation is missing at 2025 Oct
print(TSCV1)
## # 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(IVAE ~ pdq(0, 0, … Test   38.5  57.9  42.5  16.8  20.2 0.920 0.994 0.445

Interpretación del MAPE: 21.11083 = en promedio, los pronósticos del modelo tienen un error absoluto del 21.11% en relación con los valores reales de la serie temporal.