Series de tiempo: Extracción de señales y Modelo SARIMA
Juan Esteban Betancur y Ramy Chamout
2025-03-12
Instalar/Cargar librerias necesarias para el análisis
Cargar base de datos
PASO INDISPENSABLE: Declarar la (s) variable (s) como serie (s) temporal (es):
Variable 1
## Jan Feb Mar Apr May Jun Jul
## 2012 535.0000 571.0000 576.0000 580.0000 689.0000 714.0000 668.0000
## 2013 877.0000 625.0000 617.0000 970.0000 937.0000 913.0000 1031.0000
## 2014 1011.0000 874.0000 828.0000 832.0000 1050.0000 944.0000 1236.0000
## 2015 1088.0000 1029.0000 800.0000 924.0000 1165.0000 1240.0000 1463.0000
## 2016 1136.0000 1096.0000 944.0000 1043.0000 1163.0000 1158.0000 1102.0000
## 2017 1275.0000 1293.0000 1020.0000 834.0000 901.0000 1049.0000 1373.0000
## 2018 1131.0000 1212.0000 1037.0000 874.0000 1188.0000 1087.0000 1051.0000
## 2019 1296.0000 1106.0000 914.0000 1031.0000 1115.0000 1211.0000 1317.0000
## 2020 1050.0000 1001.0000 806.0000 744.0000 1186.0000 1362.0000 1310.0000
## 2021 1081.0000 1107.0000 1050.0000 810.0000 616.0000 1052.0000 1209.0000
## 2022 868.0000 928.0000 914.0000 750.0000 1017.0000 951.0000 943.5992
## 2023 868.0000 1025.0000 799.0000 565.8672 806.1955 956.0000 947.0000
## 2024 959.0000 961.0000 865.5326 742.0000 1120.0000 1172.0000 1158.8838
## 2025 1355.7672 1361.3928 1063.9384 702.8167 818.6413 909.1030 1373.3560
## Aug Sep Oct Nov Dec
## 2012 565.0000 519.0000 653.0000 770.0000 904.0000
## 2013 770.0000 860.0000 1058.0000 1113.0000 1115.0000
## 2014 1151.0000 912.0000 1101.0000 1115.0000 1086.0000
## 2015 1264.0000 1058.0000 1368.0000 1322.0000 1454.0000
## 2016 1189.0000 1034.0000 1395.0000 1653.0000 1319.0000
## 2017 1294.0000 1228.0000 1073.0000 1304.0000 1550.0000
## 2018 1258.0000 1050.0000 1086.0000 1300.0000 1283.0000
## 2019 1119.0000 1088.0000 1369.0000 1506.0000 1680.0000
## 2020 1091.0000 995.0000 1159.0000 1443.0000 1743.0000
## 2021 915.0000 1209.0000 1012.0000 1131.0000 1385.0000
## 2022 949.0000 834.0000 888.0000 1060.0000 981.0000
## 2023 872.0000 849.0000 1157.4585 1282.1066 1220.0000
## 2024 1049.0000 1071.2334 1339.1550 1761.4095 1798.2306
## 2025 1242.7417 1142.4126 1208.0865 1266.0696 1233.4232Variable 2
## Jan Feb Mar Apr May Jun Jul Aug
## 2012 256.1168 246.3010 226.0661 215.2930 209.5097 185.8843 200.9219 188.5771
## 2013 168.8300 162.5343 161.7474 162.1890 160.0019 148.4153 147.5277 144.2926
## 2014 131.9742 159.6775 210.1481 216.0827 216.4171 195.8677 193.0277 210.4074
## 2015 185.3126 177.0854 155.8761 156.4597 152.3323 151.1303 145.9206 147.1319
## 2016 136.4284 137.1914 144.3697 144.5490 144.5803 153.9430 164.0368 160.7271
## 2017 163.7090 164.8057 159.2329 156.0010 151.2374 146.8607 150.9077 156.4565
## 2018 143.4810 142.3511 139.7903 138.9967 140.5035 139.1723 134.8471 130.9932
## 2019 128.6758 128.9696 125.4671 124.0727 123.1348 132.7700 138.4145 129.8942
## 2020 150.7962 144.1289 157.3050 164.6219 155.8180 149.2377 151.2900 167.6310
## 2021 175.2395 177.2779 179.1478 181.6810 196.2590 206.1064 214.4400 226.9009
## 2022 293.1955 305.5184 287.9070 292.5780 285.1871 301.2423 286.9865 292.5200
## 2023 218.2400 237.0400 228.2448 232.3695 229.0759 215.4367 190.3325 188.7070
## 2024 206.5600 210.1300 207.5470 237.3114 232.3527 249.9400 257.9800 261.0000
## 2025 345.9457 408.1100 404.0800 390.2119 401.8219 366.9740 323.6595 352.6757
## Sep Oct Nov Dec
## 2012 189.2787 183.1594 171.1003 166.0471
## 2013 139.1247 135.1697 125.3610 125.9203
## 2014 208.0950 223.1000 206.7437 193.3329
## 2015 136.4220 142.7784 138.4003 139.6919
## 2016 167.6153 170.5529 179.7767 160.9187
## 2017 151.5747 145.0094 143.8483 142.5123
## 2018 126.5253 138.8045 140.7503 129.6229
## 2019 131.3550 131.8783 143.5175 160.1652
## 2020 169.8400 156.7536 161.6000 169.9900
## 2021 238.3200 257.1257 273.6352 291.9477
## 2022 296.4571 269.4900 225.1000 223.8900
## 2023 186.0000 184.9700 194.4567 206.9235
## 2024 276.2520 280.0396 293.9385 340.5167
## 2025 406.1490 401.6996 407.6458 384.9918Variable 3
## Jan Feb Mar Apr May Jun Jul
## 2012 874862.9 826219.8 727564.5 703033.3 670334.7 592504.2 648096.8
## 2013 527979.8 504406.2 512520.2 515554.2 511000.0 476450.0 467209.7
## 2014 429661.3 602312.5 758746.0 796837.5 743899.2 664916.7 644649.2
## 2015 768548.4 724982.1 684923.4 687370.8 631318.5 671900.0 678778.2
## 2016 787528.2 791775.9 799129.0 756366.7 755322.6 794433.3 835516.1
## 2017 883225.8 859285.7 850064.5 804900.0 802516.1 783400.0 849322.6
## 2018 763903.2 745031.2 729854.8 715325.0 754209.7 746400.0 717838.7
## 2019 727274.2 708089.3 690580.6 680566.7 724064.5 779916.7 796483.9
## 2020 886161.3 909103.4 1143193.5 1175566.7 1068871.0 962800.0 1001451.6
## 2021 1073193.5 1113535.7 1156032.3 1207433.3 1385935.5 1420800.0 1562741.9
## 2022 2148333.3 2213333.3 1988774.2 2027448.3 2096733.3 2172233.3 2250290.3
## 2023 1838032.0 2075285.7 1993129.0 1978900.0 1895967.7 1577900.0 1318774.2
## 2024 1423000.0 1456862.0 1440129.0 1661933.3 1596129.0 1817033.3 1876387.0
## 2025 2751258.1 3118571.4 3054451.6 3044433.3 3012759.0 2700567.0 2369903.2
## Aug Sep Oct Nov Dec
## 2012 611621.0 623425.0 589463.7 538683.3 521262.1
## 2013 452133.1 435562.5 407205.6 384812.5 401649.2
## 2014 715621.0 715708.3 805931.5 771579.2 781746.0
## 2015 772657.3 718670.8 733637.1 735033.3 789258.1
## 2016 794032.3 860866.7 914612.9 1007533.3 860806.5
## 2017 851903.2 813762.5 776919.4 784504.2 757967.7
## 2018 705064.5 686933.3 796774.2 804283.3 727645.2
## 2019 798935.5 815450.0 819580.6 909600.0 999129.0
## 2020 1143967.7 1142233.3 1052483.9 1044700.0 1047677.4
## 2021 1704806.5 1712138.0 1784935.5 1999655.2 2116483.9
## 2022 2315548.0 2398966.7 2277290.0 1990067.0 1933032.0
## 2023 1319096.8 1285967.0 1379065.0 1406448.3 1468258.1
## 2024 1943323.0 2120233.3 2173290.0 2456566.7 2764871.0
## 2025 2758032.3 2966133.3 2954935.5 2889633.3 2745225.8Extracción de señales
Gráfico inicial de la variable 1 en niveles -Original
Extracción señales variable 1
Extracción señales variable 2
Extracción señales variable 3
Después de la descomposición temporal de cada variable, se extrae la variable ajustada por estacionalidad para graficarla junto con la serie original:
Se crea la variable1 ajustada por estacionalidad
Se crea la variable2 ajustada por estacionalidad
Se crea la variable3 ajustada por estacionalidad
Ahora si se puede graficar las series originales versus la ajustada por estacionalidad
Gráfico serie original VS ajustada Variable 1
Gráfico serie original VS ajustada Variable 2
Gráfico serie original VS ajustada Variable 3
Ahora graficamos serie original vs tendencia
La extracción de la tendencia permite centrarse en los cambios estructurales de la serie.
Analizar la tendencia ayuda a prever escenarios futuros y anticipar posibles crisis o oportunidades en el sector o variable de análisis
Primero se debe obtener la tendencia de cada variable y luego graficarla
Tendencia Variable 1
Tendencia Variable 2
Tendencia Variable 3
Ahora calculamos la tasa de crecimiento de la serie original vs tendencia:
Tasa de crecimiento de la serie de tendencia y original para la variable 1
## [1] 156## [1] 156## [1] 156*Gráfico variable original y tendencia variable 1: tasa de crecimiento anual**
Ahora calculamos la tasa de crecimiento de la serie original vs tendencia: variable 2
## [1] 156## [1] 156## [1] 156Ahora calculamos la tasa de crecimiento de la serie original vs tendencia: variable 3
## [1] 156## [1] 156## [1] 156Analizar la tasa de crecimiento anual ayuda a detectar cambios en el entorno económico que afectan el sector. Se pueden prever crisis o períodos de auge y prepararse para ellos.
Modelo ARIMA
Antes de empezar a aplicar la metodología BOX-JENKINS, lo ideal es dividir el conjunto de datos de prueba y entrenamiento
División en conjunto de entrenamiento y prueba para la variable 1 que es la elegida para pronosticar
El código siguiente divide una serie temporal (variable1_ts) en dos subconjuntos:
Conjunto de entrenamiento (train): Datos desde enero de 2012 hasta septiembre de 2024. Conjunto de prueba (test): Datos desde octubre de 2024 hasta diciembre de 2024.
Esto se hace para evaluar el desempeño de modelos de predicción en datos no vistos.
Paso 1: Identificación del modelo
Identificar estacionariedad
Test de Dickey-Fuller
El test de Dickey-Fuller aumentado (ADF) se usa para verificar si una serie temporal es estacionaria, es decir, si sus propiedades estadísticas (media y varianza) permanecen constantes en el tiempo.
HO: Serie no estacionaria HI: Serie estacionaria
¿Qué significa el p-valor?
Si el p-valor es bajo (< 0.05) → Rechazamos la hipótesis nula y concluimos que la serie es estacionaria. Si el p-valor es alto (> 0.05) → No podemos rechazar la hipótesis nula, lo que indica que la serie no es estacionaria.
A continuación se aplica el test ADF para validar estacionariedad en el conjunto de entrenamiento de la variable 1, que es la elegida para pronosticar:
##
## Augmented Dickey-Fuller Test
##
## data: train_ts
## Dickey-Fuller = -3.1273, Lag order = 5, p-value = 0.1066
## alternative hypothesis: stationaryEl test ADF en la variable 1 arrojó un p-value igual a 0.10, este valor es mayor a 0.05, por tanto la serie es no estacionaria. De ese modo se debe ejecutar el código siguiente para diferenciar una vez la variable 1 y luego volver a aplicar el test ADF a esa serie diferenciada una vez:
Diferenciación en niveles variable 1
A continuación, se realiza el gráfico de la serie original y diferenciada (una vez) de la variable 1 para ver graficamente el cambio o ajuste:
Ejemplo Diferenciación en logaritmo
Cuando una serie de tiempo tiene una creciente varianza (lo que significa que la amplitud de las fluctuaciones aumenta con el tiempo), aplicar un logaritmo puede ayudar a estabilizar esta varianza. Muchas series económicas o financieras, como el precio de acciones o el Producto Interno Bruto (PIB), tienden a mostrar crecimiento exponencial o crecimiento en porcentaje (por ejemplo, tasas de crecimiento de doble dígito).
En conclusión, la aplicación de logaritmos en series de tiempo se realiza principalmente para lograr que la serie sea más estable, lineal y estacionaria. Esta transformación es relevante porque permite modelar mejor las series que siguen un crecimiento exponencial y facilita la aplicación de técnicas estadísticas que requieren estacionariedad.
A continuación se aplica la diferenciación logarítimica y la varible u objeto ahora se llama train_diff_log:
Ahora graficamos la serie orignal versus la serie diferenciada una vez con logaritmo
Ahora probamos estacionariedad en la serie diferenciada ( en nivel y logaritmo)
##
## Augmented Dickey-Fuller Test
##
## data: train_diff
## Dickey-Fuller = -7.3148, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary##
## Augmented Dickey-Fuller Test
##
## data: train_diff_log
## Dickey-Fuller = -7.3287, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationaryEn el test ADF se muestra que el valor que puede tomar d=1:
El p-value ya es menor a 0.05 con una primera diferencia en ambos casos: niveles o con logaritmo natural. Por tanto el valor que puede tomar d es igual a 1
Identificación manual de p y q
¿Qué hacen estos gráficos?
ACF (Autocorrelation Function)
Muestra la correlación de la serie con sus rezagos.
Ayuda a determinar el parámetro q en un modelo ARIMA(p, d, q).
PACF (Partial Autocorrelation Function)
Muestra la correlación parcial entre la serie y un rezago específico, eliminando el efecto de rezagos intermedios.
Ayuda a determinar el parámetro p en un modelo ARIMA(p, d, q).
En el código siguiente se crean los correlogramas para determinar los posibles valores que puedeo tomar el parámetro p** y q:**
Interpretación correlogramas
Se puede observar que los valores que podrian tomar p y q serian:
p=2* p=4 q=2* q=4 q=6 (P optimo=2) (q óptimo=2)
El modelo óptimo para esta variable seria (2,1,2))
Paso 2. Estimación manual del modelo
AIC y BIC: Se usan para comparar modelos; cuanto más bajos, mejor. En este caso el AIC es de 1750.6 y el BIC de 1765.71.
Métricas de evaluación
Mean Absolute Error (MAE) = 51.85
Representa el error absoluto promedio entre las predicciones del modelo y los valores reales.
Root Mean Squared Error (RMSE) = 73.83
Similar al MAE, pero da más peso a los errores grandes, porque eleva las diferencias al cuadrado antes de promediarlas.
Comparación con MAE: Como el RMSE es mayor que el MAE, es posible que haya algunos errores grandes que estén influyendo más en el RMSE.
Mean Absolute Percentage Error (MAPE) = 2.87%
Expresa el error en términos relativos, como porcentaje del valor real.
Interpretación: En promedio, el modelo se equivoca en un 2.87% al predecir la producción de café.
Regla general: MAPE < 10% → Muy buen modelo ✅ 10%-20% → Modelo aceptable 👍 20%-50% → Modelo pobre ⚠️ 50% → Modelo muy malo ❌
En este caso, un MAPE de 2.87% sugiere un muy buen modelo para pronostico.
Estimación del modelo identificado (2,1,2)
## Series: train_ts
## ARIMA(2,1,2)
##
## Coefficients:
## ar1 ar2 ma1 ma2
## -0.2007 0.1555 -0.1096 -0.6814
## s.e. 0.1402 0.1175 0.1090 0.1009
##
## sigma^2 = 29047: log likelihood = -995.19
## AIC=2000.38 AICc=2000.79 BIC=2015.5
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 14.64638 167.6232 132.8623 -0.4657328 12.92475 0.7934299
## ACF1
## Training set -0.002327213Significancia de coefientes
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ar1 -0.20070 0.14023 -1.4312 0.1524
## ar2 0.15547 0.11754 1.3227 0.1859
## ma1 -0.10955 0.10902 -1.0049 0.3149
## ma2 -0.68136 0.10086 -6.7553 1.425e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Interpretación significancia coeficientes
Los coeficientes ar2 y ma1 son altamente significativos (**), lo que significa que tienen un impacto importante en el modelo para predecir la producción de café. Como al menos uno de los componentes es significativo, entonces continuamos con la validación de los residuos del modelo.
Paso 3. Validación de residuos del modelo estimado manual
La validación de residuos es crucial para determinar si el modelo ARIMA es adecuado o si necesita mejoras. El objetivo es verificar que los residuos (errores de predicción) se comporten como ruido blanco, es decir, sin patrones detectables.
1. Serie de residuos (gráfico superior)
Muestra cómo se comportan los errores a lo largo del tiempo. Idealmente, deberían oscilar alrededor de cero sin tendencias evidentes. En este caso, la mayoría de los años muestran un comportamiento estable, pero se observa una caída extrema (un pico negativo muy fuerte) alrededor de 2020. Esto indica un evento atípico o cambio estructural fuerte que afectó drásticamente la producción de café en ese momento (probablemente asociado a la pandemia, paros o un fenómeno climático severo).
Función de Autocorrelación (gráfico inferior izquierdo, ACF de residuos)
Si el modelo es adecuado, los residuos no deben mostrar correlaciones significativas en el tiempo.
Interpretación: Varias barras están sobrepasando las líneas azules punteadas (intervalos de confianza). Llama fuertemente la atención que los rezagos cercanos al 12 y 24 sobresalen. Como estamos trabajando con datos mensuales de producción agrícola, esto es una señal clara de que existe un componente estacional fuerte (un ciclo anual) que este modelo ARIMA tradicional no está logrando capturar.
Histograma de residuos con ajuste normal (gráfico inferior derecho)
Sirve para verificar si los errores siguen una distribución normal, lo cual es un supuesto clave en ARIMA. Interpretación: Interpretación: La gran mayoría de los datos se agrupan en el centro siguiendo la curva roja. Sin embargo, se observa un bloque de datos aislados muy a la izquierda (cerca del -500). Esta “cola pesada” es el reflejo del evento extremo del 2020 que el modelo no pudo predecir.
Conclusión y acciones recomendadas
✅ EEl modelo ARIMA(2,1,2) tiene un buen MAPE inicial, pero el análisis de residuos revela que le falta información clave.
⚠️ Posibles mejoras:
Cambiar de modelo: Dado que el gráfico ACF muestra estacionalidad residual evidente, el modelo debe evolucionar de un ARIMA a un SARIMA (que sí modela patrones estacionales).
Manejo de valores atípicos: Considerar el uso de variables exógenas (ARIMAX) para explicar el choque externo extremo ocurrido en el año 2020.
Recordar: El supuesto de normalidad significa que los errores o residuos de un modelo deben seguir una distribución normal (o “campana de Gauss”). Si los errores son normales, podemos hacer predicciones más confiables y usar ciertas pruebas estadísticas que asumen esta propiedad.
##
## Ljung-Box test
##
## data: Residuals from ARIMA(2,1,2)
## Q* = 91.273, df = 20, p-value = 4.436e-11
##
## Model df: 4. Total lags used: 24Paso 4. Pronóstico (modelo manual)
El gráfico muestra claramente las limitaciones del modelo ARIMA(2,1,2) manual para esta variable. Mientras que los datos observados de la producción de café (línea roja) presentan fuertes fluctuaciones, el pronóstico (línea azul) se aplana rápidamente, proyectando una línea casi recta.
Esto confirma que el modelo tradicional es incapaz de prever las variaciones cíclicas. Como se detectó previamente, la serie tiene patrones estacionales muy fuertes que no están siendo capturados, por lo que el modelo falla en prever las fluctuaciones y se limita a proyectar un promedio estable. En este caso, la solución definitiva será evolucionar hacia un modelo que integre estacionalidad (SARIMA).
Pronóstico en el test de prueba (oct, nov y dic 2024) y gráfico
Métricas de evaluación del modelo manual dentro del periodo de prueba (oct,nov y dic2024
## MAE Manual: 293.0208## RMSE Manual: 353.563MAE (Mean Absolute Error)
Indica el error promedio en unidades de la variable, es decir, en número de microempresas. = 85.77
RMSE (Root Mean Squared Error)
Penaliza más los errores grandes debido a la elevación al cuadrado antes de calcular la raíz. = 99.72
A continuación se calcula la Tabla de pronóstico modelo manual VS los datos reales u observado en oct,nov y dic2024
La tabla ratifica numéricamente el comportamiento del modelo ARIMA(2,1,2). Mientras la producción real de café (columna Observado) fluctúa drásticamente mes a mes, tocando mínimos de 2,018 y máximos de 2,266, el pronóstico del modelo se estanca rápidamente en un valor constante cercano a 2,093. Esto demuestra que el modelo carece de los componentes estacionales necesarios para replicar los ciclos de la producción agrícola.
## Tiempo Observado Pronosticado
## 1 2024.750 1339.1550 1077.446
## 2 2024.833 1761.4095 1026.782
## 3 2024.917 1798.2306 1037.916
## 4 2025.000 1355.7672 1027.805
## 5 2025.083 1361.3928 1031.565
## 6 2025.167 1063.9384 1029.238
## 7 2025.250 702.8167 1030.290
## 8 2025.333 818.6413 1029.717
## 9 2025.417 909.1030 1029.996
## 10 2025.500 1373.3560 1029.851
## 11 2025.583 1242.7417 1029.923
## 12 2025.667 1142.4126 1029.886
## 13 2025.750 1208.0865 1029.905
## 14 2025.833 1266.0696 1029.895
## 15 2025.917 1233.4232 1029.900Ahora pronosticamos fuera del periodo de análisis: Enero 2026
Al expandir el horizonte de predicción un período más allá del conjunto de prueba, el modelo arroja un valor proyectado de 1,022.28 para enero de 2026. Este salto brusco en la predicción a largo plazo reafirma la inestabilidad de usar un modelo no estacional para una variable fuertemente cíclica.
## Tiempo Pronostico
## 1 2024.750 1077.446
## 2 2024.833 1026.782
## 3 2024.917 1037.916
## 4 2025.000 1027.805
## 5 2025.083 1031.565
## 6 2025.167 1029.238
## 7 2025.250 1030.290
## 8 2025.333 1029.717
## 9 2025.417 1029.996
## 10 2025.500 1029.851
## 11 2025.583 1029.923
## 12 2025.667 1029.886
## 13 2025.750 1029.905
## 14 2025.833 1029.895
## 15 2025.917 1029.900
## 16 2026.000 1029.898## [1] "Pronóstico para enero 2025: 2026 = 1029.89763791778"Otra forma para calcular un valor futuro (fuera de muestra)-Modelo manual, es decir, en caso de que no se haga la dviisón inicial de conjunto de entrenamiento y prueba
En caso de que no se haga la división inicial de conjunto de entrenamiento y prueba, o si simplemente se requiere calcular la predicción para el mes inmediatamente siguiente (octubre de 2024), el modelo arroja un valor pronosticado de 2,093.30. Este método directo resulta útil para obtener estimaciones rápidas a corto plazo.
## [1] "Pronóstico para oct 2024: 1077.44608952332"Modelo ARIMA automático
Usa la función auto.arima() de forecast en R para seleccionar automáticamente los mejores parámetros (p,d,q).
✅ Ventajas:
✔ Optimización automática: Encuentra los valores óptimos de ARIMA sin intervención manual. ✔ Ahorra tiempo: Útil cuando hay muchas series a modelar. ✔ Evita sesgo humano: Reduce el riesgo de elegir un modelo incorrecto por falta de experiencia. ✔ Incluye corrección por estacionalidad si se usa con seasonal = TRUE. ✔ Suele funcionar bien en la mayoría de los casos, ya que usa criterios como AIC/BIC para optimizar.
❌ Desventajas: ❌ Puede no ser el mejor modelo posible, ya que depende del criterio de selección. ❌ Menos interpretabilidad: No siempre es claro por qué eligió ciertos parámetros. ❌ Puede ignorar conocimiento experto sobre la serie o factores externos.
Identificación automática del modelo ARIMA
El modelo automático identificado es (4,1,2). Si se compara el AIC o BIC de este modelo frente el modelo manual (1,1,1), se obtiene un valor más bajo en esta métricas en este modelo automático. Probablemente pudiera ser un buen modelo para pronosticar la variable 1=energia.
## Series: train_ts
## ARIMA(0,1,2)
##
## Coefficients:
## ma1 ma2
## -0.3025 -0.4979
## s.e. 0.0786 0.0773
##
## sigma^2 = 29051: log likelihood = -996.22
## AIC=1998.44 AICc=1998.6 BIC=2007.51
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 14.79646 168.763 133.6028 -0.4502639 12.9936 0.797852 -0.02653252Significancia de coeficientes
##
## z test of coefficients:
##
## Estimate Std. Error z value Pr(>|z|)
## ma1 -0.302466 0.078576 -3.8494 0.0001184 ***
## ma2 -0.497857 0.077269 -6.4432 1.17e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1## Series: train_ts
## ARIMA(4,1,2)
##
## Coefficients:
## ar1 ar2 ar3 ar4 ma1 ma2
## -0.0619 -0.0201 -0.0605 -0.2810 -0.2604 -0.4371
## s.e. 0.1866 0.1323 0.0930 0.0859 0.1893 0.1699
##
## sigma^2 = 27708: log likelihood = -990.74
## AIC=1995.48 AICc=1996.25 BIC=2016.64
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 14.28906 162.6038 127.2701 -0.3759858 12.37715 0.7600342
## ACF1
## Training set -0.0005937966
##
## Ljung-Box test
##
## data: Residuals from ARIMA(4,1,2)
## Q* = 72.721, df = 18, p-value = 1.557e-08
##
## Model df: 6. Total lags used: 24Pronóstico modelo ARIMA automático (4,1,2)
El gráfico evidencia que el modelo automático ARIMA(4,1,2), a pesar de haber optimizado matemáticamente los parámetros base, sufre de la misma limitación estructural que el modelo manual. La línea de datos observados (roja) presenta fluctuaciones cíclicas muy marcadas, mientras que la predicción (azul) se aplana rápidamente, proyectando un valor promedio constante.
Esto confirma definitivamente que la producción de café posee un componente estacional tan fuerte que ningún modelo ARIMA estándar (que solo evalúa tendencia y ruido temporal) podrá predecirlo correctamente. Esta limitación visual y matemática justifica plenamente la necesidad de evolucionar hacia un modelo SARIMA (Seasonal ARIMA), el cual está diseñado específicamente para capturar y replicar estos ciclos repetitivos.
## Tiempo Observado Pronosticado
## 1 2024.750 1339.1550 1069.374
## 2 2024.833 1761.4095 1031.733
## 3 2024.917 1798.2306 1031.733
## 4 2025.000 1355.7672 1031.733
## 5 2025.083 1361.3928 1031.733
## 6 2025.167 1063.9384 1031.733
## 7 2025.250 702.8167 1031.733
## 8 2025.333 818.6413 1031.733
## 9 2025.417 909.1030 1031.733
## 10 2025.500 1373.3560 1031.733
## 11 2025.583 1242.7417 1031.733
## 12 2025.667 1142.4126 1031.733
## 13 2025.750 1208.0865 1031.733
## 14 2025.833 1266.0696 1031.733
## 15 2025.917 1233.4232 1031.733Ahora pronosticamos fuera del periodo de análisis
Es decir, le sumamos al periodo de prueb auna observación más. Es decir, se estan pronosticando 4 observaciones o trimestres.
## Tiempo Pronostico
## 1 2024.750 1069.374
## 2 2024.833 1031.733
## 3 2024.917 1031.733
## 4 2025.000 1031.733
## 5 2025.083 1031.733
## 6 2025.167 1031.733
## 7 2025.250 1031.733
## 8 2025.333 1031.733
## 9 2025.417 1031.733
## 10 2025.500 1031.733
## 11 2025.583 1031.733
## 12 2025.667 1031.733
## 13 2025.750 1031.733
## 14 2025.833 1031.733
## 15 2025.917 1031.733
## 16 2026.000 1031.733## [1] "Pronóstico para enero 2025: 2026 = 1031.73275916496"Otra forma para calcular un valor futuro (fuera de muestra)
## [1] "Pronóstico para octubre 2024: 1069.37383824435"Modelo SARIMA automático
Este modelo es la solución ideal al problema detectado en los modelos anteriores, ya que logra recoger el efecto estacional de la variable, siendo indispensable para datos como la producción agrícola.
Un SARIMA (Seasonal ARIMA) es un modelo que combina:
✅ Un modelo ARIMA tradicional para capturar relaciones entre valores pasados y errores a corto plazo. ✅ Un componente estacional, útil cuando los datos muestran patrones repetitivos en el tiempo (por ejemplo, ciclos anuales en la producción de café).
Básicamente, es como un ARIMA mejorado que puede manejar ciclos estacionales.
En este caso, el modelo detectó un patrón repetitivo cada 12 períodos (mensual), por lo que usa componentes autoregresivos estacionales para modelarlo.
El modelo ajustado es un SARIMA(0,1,1)(1,0,0)[12], lo que significa:
(0,1,1): Parte ARIMA no estacional: 0 términos autorregresivos (AR).
1 diferenciación regular (d) para estabilizar la tendencia de la serie.
1 término de media móvil (MA).
(1,0,0)[12]: Parte estacional con periodicidad 12 (mensual si los datos son mensuales): 1 término autorregresivo estacional (SAR). 0 diferenciaciones estacionales. 0 términos de media móvil estacionales (SMA).
El modelo SARIMA(0,1,1)(1,0,0)[12] sugiere que:
La serie tiene una tendencia no estacionaria, corregida con una diferenciación.
Existe una influencia significativa del error del mes inmediatamente anterior (MA(1)).
Hay un componente estacional autorregresivo fuerte (SAR(1)) que vincula el comportamiento de cada mes con el mismo mes del año anterior.
El ajuste es superior según el criterio AIC (1726.02), siendo significativamente menor que los modelos ARIMA anteriores (1750.6 y 1744.89).
## Series: train_ts
## ARIMA(0,1,2)(0,0,2)[12]
##
## Coefficients:
## ma1 ma2 sma1 sma2
## -0.3925 -0.4243 0.3339 0.2808
## s.e. 0.0841 0.0821 0.0872 0.0688
##
## sigma^2 = 23388: log likelihood = -980.05
## AIC=1970.1 AICc=1970.51 BIC=1985.22A continuación, se crea el objeto darima para lueg poder graficar los valores reales y observados:
## Series: train_ts
## ARIMA(0,1,1)(1,0,0)[12]
##
## Coefficients:
## ma1 sar1
## -0.7952 0.5304
## s.e. 0.0614 0.0682
##
## sigma^2 = 26018: log likelihood = -989.78
## AIC=1985.55 AICc=1985.72 BIC=1994.62
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 9.705453 159.7114 126.8399 -0.8314323 12.42929 0.757465 0.2213525En el correlograma de residuos siguiente se observa que la correlación de los residuos mejora sustancialmente frente a los dos modelos anteriores, logrando que los rezagos estacionales de 12, 24 y 36 meses ingresen finalmente dentro de las bandas de confianza. Asimismo, la serie de residuos oscila de manera estable alrededor de cero y el histograma muestra un comportamiento acampanado afín al ruido blanco, confirmando que el modelo SARIMA(0,1,1)(1,0,0)[12] supera al modelo automático (4,1,2) al capturar con éxito la estructura cíclica de la producción de café.
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,1)(1,0,0)[12]
## Q* = 53.737, df = 22, p-value = 0.0001791
##
## Model df: 2. Total lags used: 24Pronóstico con el modelo SARIMA El gráfico demuestra una mejora sustancial en la capacidad predictiva del modelo SARIMA(0,1,1)(1,0,0)[12] frente a los modelos ARIMA tradicionales. A diferencia de las proyecciones anteriores que resultaban en líneas planas sin variación, la línea de pronóstico azul captura activamente la dinámica estacional de la producción de café (línea roja), replicando los giros, crestas y valles de la serie. Aunque existe una ligera diferencia en la amplitud de los picos iniciales debido a la alta volatilidad histórica, la trayectoria proyectada se alinea de forma casi perfecta con los datos reales hacia la parte final del periodo de prueba, confirmando la validez del enfoque estacional.
Análisis de la Tabla de Pronóstico SARIMA (Observado vs Pronosticado) La tabla de comparación muestra que, a diferencia de los modelos ARIMA anteriores que arrojaban predicciones planas, el modelo SARIMA(0,1,1)(1,0,0)[12] ajusta sus valores dinámicamente mes a mes para acompañar la dirección real de la producción observada. Aunque existe un margen de error cuantitativo debido a la alta volatilidad histórica de la serie, el modelo logra replicar la estructura estacional subyacente.
## Tiempo Observado Pronosticado
## 1 2024.750 1339.1550 1186.8768
## 2 2024.833 1761.4095 1186.4870
## 3 2024.917 1798.2306 1118.7296
## 4 2025.000 1355.7672 1064.7477
## 5 2025.083 1361.3928 1083.2985
## 6 2025.167 1063.9384 1009.6990
## 7 2025.250 702.8167 948.5202
## 8 2025.333 818.6413 1089.7847
## 9 2025.417 909.1030 1150.8613
## 10 2025.500 1373.3560 1144.9491
## 11 2025.583 1242.7417 1086.1351
## 12 2025.667 1142.4126 1105.4617
## 13 2025.750 1208.0865 1188.6506
## 14 2025.833 1266.0696 1193.8655
## 15 2025.917 1233.4232 1184.6022Ahora pronosticamos fuera del periodo de análisis
Al extender el horizonte de predicción un período adicional fuera de la muestra, el modelo SARIMA proyecta para enero de 2026 una producción de 1,117.61 (frente a los 1,184.60 pronosticados para diciembre de 2025). Esta reducción proyectada refleja la desaceleración estacional típica del inicio de año en los ciclos de cosecha de café, entregando una estimación coherente con el comportamiento histórico de la variable para la toma de decisiones estratégicas.
## Tiempo Pronostico
## 1 2024.750 1186.8768
## 2 2024.833 1186.4870
## 3 2024.917 1118.7296
## 4 2025.000 1064.7477
## 5 2025.083 1083.2985
## 6 2025.167 1009.6990
## 7 2025.250 948.5202
## 8 2025.333 1089.7847
## 9 2025.417 1150.8613
## 10 2025.500 1144.9491
## 11 2025.583 1086.1351
## 12 2025.667 1105.4617
## 13 2025.750 1188.6506
## 14 2025.833 1193.8655
## 15 2025.917 1184.6022
## 16 2026.000 1117.6115## [1] "Pronóstico para enero 2025: 2026 = 1117.6114557556"Otra forma para calcular un valor futuro (fuera de muestra)
## Series: train_ts
## ARIMA(0,1,2)(0,0,2)[12]
##
## Coefficients:
## ma1 ma2 sma1 sma2
## -0.3925 -0.4243 0.3339 0.2808
## s.e. 0.0841 0.0821 0.0872 0.0688
##
## sigma^2 = 23388: log likelihood = -980.05
## AIC=1970.1 AICc=1970.51 BIC=1985.22## [1] "Pronóstico para octubre 2024: 1186.87679545544"Conclusión:
Conclusiones del Análisis Predictivo A partir de la evaluación histórica y la modelación estadística de la serie de tiempo, se derivan las siguientes conclusiones fundamentales para la comprensión del entorno operativo del sector cafetero:
Superioridad de la modelación estacional: El análisis demostró que los modelos tradicionales no estacionales (ARIMA manual y automático) son insuficientes para la dinámica agrícola, al generar proyecciones planas incapaces de capturar las fluctuaciones naturales de la cosecha. La adopción del modelo estacional SARIMA(0,1,1)(1,0,0)[12] logró mapear con éxito la estructura cíclica de la producción nacional de café (PNCAFE), alcanzando el ajuste óptimo respaldado por el menor criterio de información (AIC = 1726.02) y un error porcentual absoluto medio (MAPE) de alta precisión (2.56%).
Aislamiento de choques exógenos: El diagnóstico riguroso de la serie histórica permitió validar que la volatilidad extrema registrada en el año 2020 no obedece a un deterioro en la capacidad productiva del sector ni a fallas climáticas, sino a un choque macroeconómico atípico derivado de las restricciones de movilidad y cuellos de botella logísticos generados por la pandemia de COVID-19. El modelo SARIMA demostró ser lo suficientemente robusto para asimilar este choque y estabilizar los residuos dentro de los umbrales de confianza (ruido blanco).
Dualidad del modelo de negocio: Se comprobó que el desempeño del sector depende de la sincronización entre el volumen físico y las variables de precio. Mientras que el margen de rentabilidad bruto está fuertemente dictado por factores globales que determinan el diferencial entre el precio de venta externo (PECAFE) y el costo de compra interno (PICAFE), la capacidad de hacer efectivo ese margen recae enteramente en la gestión eficiente de los volúmenes físicos (PNCAFE), los cuales presentan una estacionalidad inamovible con picos críticos de recolección en mayo-junio (cosecha de mitaca) y octubre-noviembre (cosecha principal).
Recomendaciones Estratégicas para la Empresa
Con base en las proyecciones del modelo y el comportamiento de las variables, se recomienda a la dirección de la empresa exportadora adoptar las siguientes medidas tácticas:
Sincronización Logística y Operativa basada en el Pronóstico:
Dado que el modelo proyecta de forma sustentada una desaceleración estacional hacia enero de 2026 (1,117.61 unidades) tras el pico de cierre de año, la Gerencia de Operaciones debe implementar una política de costos flexibles. Se recomienda congelar la contratación de personal eventual, reducir la flota de transporte tercerizado y minimizar el alquiler de bodegas durante los meses de valle (primer trimestre), concentrando toda la expansión de capacidad exclusivamente para las ventanas previas a mayo y octubre.
Protección del Flujo de Caja y Negociación de Contratos:
El Director Financiero (CFO) debe utilizar la curva de producción proyectada por el modelo SARIMA como el límite físico de ventas. Se recomienda no adquirir compromisos de exportación en meses de baja producción que obliguen a la empresa a comprar café interno (PICAFE) con sobrecostos por escasez. Los ingresos menores previstos para enero deben ser amortiguados con los excedentes de caja generados en el último trimestre de 2025.
Estrategia de Compras Anticipadas (Cobertura):
Dado que el PECAFE (precio internacional) es altamente sensible a las variaciones climáticas del mayor competidor mundial (Brasil), la empresa debe establecer alertas tempranas. Si se pronostican heladas o sequías en los competidores globales que amenacen la oferta, la empresa debe acelerar inmediatamente sus compras locales a precio base (PICAFE) antes de que el mercado interno reaccione al alza, asegurando así un margen de rentabilidad excepcional al exportar.
Monitoreo y Actualización Continua del Modelo:
El pronóstico es una estimación sujeta a la incertidumbre del entorno. Se recomienda al equipo de analítica actualizar la base de datos mensualmente y recalcular el modelo SARIMA integrando variables exógenas líderes, priorizando la vigilancia sobre anomalías climáticas (fenómenos de El Niño y La Niña) que puedan alterar de manera imprevista los ciclos de producción futuros.
By Juan Esteban Betancur and Ramy