| Institución | Tecnológico de Monterrey, Campus Monterrey |
| Programa | Licenciatura en Inteligencia de Negocios |
| Materia | Generación de escenarios futuros con analítica (CD-3001B). Módulo 1 |
| Alumno | Patricio Hermosillo Sepúlveda (A01567036) |
| Profesor | Raúl Alejandro Cantú |
Los cuatro ejercicios de esta actividad se resuelven con la misma caja de herramientas del Tema 2, y sirven para una idea que los recorre a todos: el método de pronóstico correcto depende de qué tiene adentro la serie. Si la serie no tiene estructura, el método más simple gana. Si tiene tendencia y estación, cualquier método que las ignore va a fallar siempre en la misma dirección.
Ejercicio 1: ventas semanales. Las doce semanas
oscilan alrededor de 19.25 sin tendencia ni patrón repetido. La prueba
de Ljung-Box no encuentra autocorrelación (p = 0.081) y
auto.arima() selecciona un ARIMA(0,0,0) con media, que es
la forma técnica de decir que la serie es ruido alrededor de un nivel
fijo. El mejor pronóstico posible es esa media. El modelo ingenuo es el
peor de los cinco justamente por lo contrario: persigue el ruido.
Ejercicio 2: leche saborizada de Hershey’s. Aquí sí hay estructura. La serie crece 12.3 por ciento anual y repite un patrón de doce meses. El ARIMA estacional con deriva es el único de los cinco modelos que puede seguir esa tendencia; los otros cuatro se quedan cortos mes tras mes, con un sesgo positivo que nunca se corrige. Ahora bien, al reestimar el modelo en cada periodo (y no solo ajustarlo una vez sobre toda la muestra) su error real resulta ser aproximadamente el doble del que sugiere la tabla de ajuste. Reporto las dos cifras, porque la segunda es la que un gerente debería usar para planear.
Ejercicio 3: Vintage Restaurant. La serie tiene una estación muy marcada: enero vale 44 por ciento más que un mes promedio y septiembre 37 por ciento menos, lo cual coincide con la temporada de turismo de invierno y con el cierre de temporada en la costa de Florida. Al quitarle la estación aparece una tendencia limpia de crecimiento. Los dos métodos que pide la actividad, descomposición y regresión con variables ficticias, producen pronósticos parecidos para el cuarto año. Con la venta real de enero (295 mil dólares) el error es de 4 mil por descomposición y de 8.2 mil por regresión: menos del 3 por ciento en los dos casos, así que el procedimiento no está fallando.
Ejercicio 4: Carlson Department Store. Estimo que la tienda habría vendido 12.36 millones de dólares entre septiembre y diciembre. Además, las ventas reales del condado en esos cuatro meses superaron en 28 por ciento a lo que se esperaba sin huracán, y trasladando ese mismo empuje a Carlson resultan 3.51 millones adicionales. El reclamo defendible es de 12.36 millones por ventas perdidas y 3.51 millones más por el exceso de actividad, siempre que la póliza reconozca el segundo concepto.
Una advertencia sobre las cifras. En los cuatro ejercicios sometí el método elegido a datos que no había visto: origen móvil en los ejercicios 1 y 2, y reserva del último año en los ejercicios 3 y 4. Las conclusiones que siguen son las que sobrevivieron esa prueba, no las que se veían mejor en la tabla de ajuste.
La actividad pide cuatro cosas. En los ejercicios 1 y 2 hay que aplicar los cinco modelos vistos en clase a una serie, resumir los resultados en una tabla, elegir el mejor y pronosticar seis periodos. En los ejercicios 3 y 4 hay que entregar un informe gerencial: no basta el número, hay que decirle a alguien qué hacer con él.
Los junto en un solo documento porque comparten método y porque se responden entre sí. Los ejercicios 1 y 2 son el mismo problema con y sin estructura, y sirven para ver qué pasa cuando se aplica un método que no corresponde. Los ejercicios 3 y 4 usan descomposición clásica sobre datos mensuales, y el 4 agrega el paso que le interesa a un negocio: comparar contra un grupo de control.
Los cuatro primeros son modelos de nivel: ninguno tiene un término de tendencia ni de estación, así que su pronóstico a varios periodos tiende a una línea plana. Eso no es un defecto mientras la serie tampoco tenga tendencia. El ejercicio 2 es donde se vuelve un problema, y ahí lo mido.
Una aclaración sobre el modelo 4. Las instrucciones dicen “suavizador exponencial (k = 3) con alfa = 0.2”. El parámetro k pertenece a los promedios móviles, no al suavizamiento exponencial, que usa toda la historia con pesos que decaen. Aplico entonces el suavizamiento con alfa = 0.2 e inicializo con F(2) = y(1), que es la convención del capítulo 18.
Tres decisiones que valen para los cuatro ejercicios.
Misma ventana para todos. El modelo ingenuo puede pronosticar desde el periodo 2, los promedios móviles desde el 4 y el ARIMA ajusta desde el 1. Si a cada uno se le calcula el error sobre los periodos que le tocan, se están comparando promedios de muestras distintas y el resultado no significa nada. Todas las tablas de este documento comparan sobre la ventana común, que empieza en el periodo 4.
Tres métricas y el sesgo. El MAE (error absoluto medio) dice de qué tamaño es el error típico; el RMSE castiga más los errores grandes; el MAPE lo expresa en por ciento y permite comparar entre series de escala distinta. Agrego el sesgo, que es el promedio del error con signo: un modelo puede tener un MAE aceptable y aun así equivocarse siempre hacia el mismo lado, y eso es lo que le cuesta dinero a una empresa que planea inventario.
Origen móvil. Las métricas de ajuste tienen un problema: el ARIMA estima sus parámetros con los mismos datos sobre los que se le mide, así que parte con ventaja. Para saber si esa ventaja es real, en los ejercicios 1 y 2 repito la comparación reestimando todo en cada periodo con la información disponible hasta entonces. Cuando las dos tablas coinciden, la conclusión es sólida. Cuando no, lo digo.
semana <- 1:12
valor <- c(17, 21, 19, 23, 18, 16, 20, 18, 22, 20, 15, 22)
ts1 <- ts(valor, start = c(2025, 1), frequency = 52)
n1 <- length(valor)Doce observaciones semanales. Antes de ajustar cualquier cosa conviene ver qué hay.
Figura 1. Ventas semanales de las 12 semanas observadas. La línea punteada es la media de la serie. Las observaciones cruzan la media doce veces sin quedarse de un lado, que es la firma de una serie sin tendencia.
La serie sube y baja sin irse a ningún lado. No hay tendencia visible y con doce semanas tampoco se puede hablar de estación, porque un ciclo anual necesita al menos 52 datos. Vale la pena confirmarlo antes de modelar.
c(media = mean(valor), desv = sd(valor),
pendiente = coef(lm(valor ~ semana))[2],
p_pendiente = summary(lm(valor ~ semana))$coefficients[2, 4])## media desv pendiente.semana p_pendiente
## 19.25000000 2.52712557 0.02447552 0.91420311
##
## Box-Ljung test
##
## data: valor
## X-squared = 9.8035, df = 5, p-value = 0.081
La pendiente no es distinta de cero y la prueba de Ljung-Box no rechaza la hipótesis de que las observaciones sean independientes. Es decir: nada de lo que pasó la semana pasada ayuda a saber qué va a pasar la próxima. Con esto ya se puede anticipar cuál de los cinco modelos va a ganar, pero los corro todos porque la actividad los pide y porque el contraste es justo lo interesante.
f1_naive <- f_naive(valor)
f1_pm <- f_pm(valor, k = 3)
f1_pmp <- f_pmp(valor)
f1_ses <- head(f_ses(valor, alpha = 0.2), n1)
arima1 <- auto.arima(ts1)
arima1## Series: ts1
## ARIMA(0,0,0) with non-zero mean
##
## Coefficients:
## mean
## 19.2500
## s.e. 0.6985
##
## sigma^2 = 6.386: log likelihood = -27.63
## AIC=59.26 AICc=60.59 BIC=60.23
auto.arima() devuelve un ARIMA(0,0,0) con media distinta
de cero. Conviene detenerse en lo que significa: cero términos
autorregresivos, cero diferencias y cero términos de media móvil. El
modelo dice que la serie es una constante más ruido, y su pronóstico
para cualquier periodo futuro es la media, 19.25. No falló al no
encontrar estructura; encontró que no la hay, que es un resultado.
idx1 <- 4:n1 # ventana comun: desde la semana 4 los cinco pueden pronosticar
res1 <- rbind(
`Ingenuo` = metricas(valor, f1_naive, idx1),
`PM(3)` = metricas(valor, f1_pm, idx1),
`PMP(3)` = metricas(valor, f1_pmp, idx1),
`SES 0.2` = metricas(valor, f1_ses, idx1),
`ARIMA` = metricas(valor, f1_arima, idx1))
round(res1, 3)## MAE RMSE MAPE Sesgo
## Ingenuo 3.889 4.203 20.235 0.333
## PM(3) 2.667 3.197 14.357 0.000
## PMP(3) 2.981 3.390 15.992 0.056
## SES 0.2 2.596 3.007 13.563 0.636
## ARIMA 2.306 2.626 12.371 0.083
El orden es el mismo en las tres métricas, lo cual da confianza: ARIMA, luego el suavizamiento exponencial, luego los promedios móviles y al final el ingenuo. El modelo ingenuo se equivoca en promedio 1.7 veces más que el ARIMA.
La razón es directa. Cuando una serie es ruido alrededor de un nivel, cualquier valor alto va seguido con alta probabilidad de uno más bajo. El modelo ingenuo toma el dato alto y lo proyecta, así que se equivoca por partida doble. Promediar los últimos tres cancela parte de ese ruido, y el suavizamiento exponencial lo cancela un poco mejor porque usa toda la historia. El ARIMA lleva esa lógica al extremo: promedia las doce.
La tabla anterior favorece al ARIMA por construcción, porque su media se estimó con las mismas doce semanas sobre las que se le está midiendo. La prueba honesta es reestimarlo en cada semana usando solo lo anterior.
## Ingenuo PM(3) PMP(3) SES ARIMA
## MAE 3.889 2.667 2.981 2.596 2.447
## RMSE 4.203 3.197 3.390 3.007 2.736
## MAPE 20.235 14.357 15.992 13.563 13.206
## Sesgo 0.333 0.000 0.056 0.636 -0.030
El orden se mantiene y el MAE del ARIMA apenas sube de 2.31 a 2.45. La conclusión no dependía del ajuste.
h <- 6
p1_arima <- forecast(arima1, h = h, level = 95)
pron1 <- data.frame(
Semana = 13:18,
Ingenuo = rep(valor[n1], h),
PM = proyecta(valor, f_pm, h),
PMP = proyecta(valor, f_pmp, h),
SES = rep(f_ses(valor, alpha = 0.2)[n1 + 1], h),
ARIMA = as.numeric(p1_arima$mean))
round(pron1, 2)| Semana | Ingenuo | PM | PMP | SES | ARIMA |
|---|---|---|---|---|---|
| 13 | 22 | 19.00 | 19.33 | 19.18 | 19.25 |
| 14 | 22 | 18.67 | 19.50 | 19.18 | 19.25 |
| 15 | 22 | 19.89 | 19.86 | 19.18 | 19.25 |
| 16 | 22 | 19.19 | 19.65 | 19.18 | 19.25 |
| 17 | 22 | 19.25 | 19.70 | 19.18 | 19.25 |
| 18 | 22 | 19.44 | 19.71 | 19.18 | 19.25 |
El modelo elegido es el ARIMA(0,0,0) con media, y su pronóstico para las seis semanas es constante en 19.25. Gana en las tres métricas, en la tabla de ajuste y en la validación con origen móvil, y sus residuos pasan la prueba de Ljung-Box.
Vale la pena decir qué significa ese resultado, porque es fácil leerlo como que el modelo no sirvió. Al contrario: la serie no contiene información que permita anticipar la semana siguiente, y el mejor pronóstico posible en ese caso es el promedio histórico. Cualquier modelo que produzca una línea con forma estaría inventando estructura que los datos no respaldan.
Lo que sí es accionable es el intervalo. Con 95 por ciento de confianza las ventas de cualquiera de las próximas seis semanas van a caer entre 14.3 y 24.2. Para planear inventario, ese rango es más útil que el punto: dice que hay que estar listo para vender 24 aunque se espere vender 19.
lech <- read.csv("datos/lechitas.csv") # 36 meses, enero 2017 a diciembre 2019
y2 <- lech$ventas
n2 <- length(y2)
ts2 <- ts(y2, start = c(2017, 1), frequency = 12)
str(lech)## 'data.frame': 36 obs. of 3 variables:
## $ anio : int 2017 2017 2017 2017 2017 2017 2017 2017 2017 2017 ...
## $ mes : int 1 2 3 4 5 6 7 8 9 10 ...
## $ ventas: num 25521 23740 26254 25868 27073 ...
Ventas mensuales de leche saborizada de Hershey’s México, en miles de
dólares. El archivo original es
Ventas_Históricas_Lechitas.xlsx, que trae los tres años en
columnas separadas y las fechas como número de serie de Excel; en el
anexo se documenta la conversión a un archivo plano de 36 renglones.
Figura 4. Ventas mensuales de leche saborizada, enero 2017 a diciembre 2019. La línea recta es la tendencia lineal ajustada. La serie sube de forma sostenida y repite un patrón dentro de cada año.
tend <- lm(y2 ~ seq_len(n2))
c(pendiente_mensual = coef(tend)[2],
p_valor = summary(tend)$coefficients[2, 4],
R2 = summary(tend)$r.squared)## pendiente_mensual.seq_len(n2) p_valor
## 298.3734043759332053014077246189 0.0000000000000000000002578471
## R2
## 0.9396981589668018663985549210
## 2017 2018 2019
## 26838.3 30567.7 33837.8
La pendiente es de 298.4 miles de dólares por mes y explica el 94 por ciento de la variación. Los promedios anuales confirman lo mismo desde otro ángulo: el negocio crece alrededor de 12.3 por ciento al año. Esta es la diferencia de fondo con el ejercicio 1, y es la que va a decidir el resultado.
f2_naive <- f_naive(y2)
f2_pm <- f_pm(y2, k = 3)
f2_pmp <- f_pmp(y2)
f2_ses <- head(f_ses(y2, alpha = 0.2), n2)
arima2 <- auto.arima(ts2)
arima2## Series: ts2
## ARIMA(1,0,0)(1,1,0)[12] with drift
##
## Coefficients:
## ar1 sar1 drift
## 0.6383 -0.5517 288.8979
## s.e. 0.1551 0.2047 14.5026
##
## sigma^2 = 202701: log likelihood = -181.5
## AIC=371 AICc=373.11 BIC=375.72
Aquí auto.arima() sí encuentra estructura. La parte
estacional, con una diferencia de orden 12, recoge el patrón que se
repite cada año; el término autorregresivo recoge la inercia de un mes
al siguiente; y la deriva recoge el crecimiento. Es el único de los
cinco modelos con un término que apunta hacia arriba.
idx2 <- 4:n2
res2 <- rbind(
`Ingenuo` = metricas(y2, f2_naive, idx2),
`PM(3)` = metricas(y2, f2_pm, idx2),
`PMP(3)` = metricas(y2, f2_pmp, idx2),
`SES 0.2` = metricas(y2, f2_ses, idx2),
`ARIMA` = metricas(y2, f2_arima, idx2))
round(res2, 2)## MAE RMSE MAPE Sesgo
## Ingenuo 827.45 1012.67 2.68 260.38
## PM(3) 870.47 1031.26 2.81 591.01
## PMP(3) 814.19 958.98 2.62 481.55
## SES 0.2 1400.32 1573.21 4.50 1352.43
## ARIMA 245.59 359.08 0.76 25.28
La columna que hay que leer con atención no es el MAE sino el sesgo. Los cuatro modelos de nivel tienen sesgo positivo, y grande: el suavizamiento exponencial se queda corto en promedio 1352 miles de dólares cada mes, el promedio móvil 591. No es mala suerte, es aritmética. Un modelo que pronostica con el promedio de los tres meses anteriores siempre va a quedar por debajo de una serie que sube, porque promedia un pasado que ya fue superado. Y el error no se corrige nunca: se repite mes tras mes en la misma dirección.
El caso del suavizamiento exponencial es el más claro. Con alfa = 0.2 el modelo solo corrige el 20 por ciento de su error en cada periodo, así que se queda todavía más atrás que el promedio móvil. Un alfa más alto reduciría el retraso, pero la actividad fija el valor y ese es justamente el punto: el parámetro se eligió sin ver los datos.
El ARIMA sale con un MAE de 246 contra 814 del mejor de los otros. Es una diferencia demasiado grande como para aceptarla sin revisarla, porque los parámetros del ARIMA se estimaron con los 36 meses completos.
# desde el mes 25: se le exigen dos años completos de historia al ARIMA, que es
# el minimo para que pueda identificar un patron estacional de 12 meses
om2 <- origen_movil(y2, frecuencia = 12, desde = 25)
round(om2, 2)## Ingenuo PM(3) PMP(3) SES ARIMA
## MAE 912.67 886.91 878.43 1303.55 526.19
## RMSE 1061.02 1048.07 996.52 1491.26 670.96
## MAPE 2.70 2.60 2.58 3.81 1.55
## Sesgo 238.45 470.19 390.47 1265.36 -255.01
El ARIMA sigue ganando con claridad, pero su MAE real es de 526, aproximadamente el doble del 246 que sugería la tabla de ajuste. Esa segunda cifra es la que hay que usar para planear. La conclusión sobre cuál modelo elegir no cambia; la expectativa sobre qué tan bien va a funcionar, sí.
Un detalle metodológico que conviene dejar escrito. Si el origen móvil se arranca en el mes 4, como en el ejercicio 1, el ARIMA queda empatado con el promedio móvil ponderado. No porque el modelo sea malo, sino porque en los primeros periodos se le está pidiendo que identifique un ciclo de doce meses con tres o cuatro datos, cosa imposible. Por eso la validación arranca en el mes 25.
##
## Box-Ljung test
##
## data: residuals(arima2)
## X-squared = 19.372, df = 12, p-value = 0.07994
Figura 6. Autocorrelación de los residuos del ARIMA estacional. Las barras que salen de la banda indican estructura que el modelo no alcanzó a capturar.
La prueba de Ljung-Box rechaza al 5 por ciento: queda algo de estructura en los residuos. Con 36 observaciones y un modelo que gasta una diferencia estacional de orden 12 no hay mucho margen para hacerlo mejor, así que lo reporto como está y anoto la consecuencia práctica: los intervalos de pronóstico que salen abajo son probablemente más angostos de lo que deberían.
p2_arima <- forecast(arima2, h = 6, level = 95)
pron2 <- data.frame(
Mes = MESES[1:6],
Ingenuo = rep(y2[n2], 6),
PM = proyecta(y2, f_pm, 6),
PMP = proyecta(y2, f_pmp, 6),
SES = rep(f_ses(y2, alpha = 0.2)[n2 + 1], 6),
ARIMA = as.numeric(p2_arima$mean))
round(pron2[, -1], 1)| Ingenuo | PM | PMP | SES | ARIMA |
|---|---|---|---|---|
| 34846.2 | 35259.7 | 35045.2 | 34308.3 | 35498.9 |
| 34846.2 | 34968.6 | 34938.0 | 34308.3 | 34202.2 |
| 34846.2 | 35024.8 | 34958.4 | 34308.3 | 36703.0 |
| 34846.2 | 35084.4 | 34966.1 | 34308.3 | 36271.9 |
| 34846.2 | 35025.9 | 34958.9 | 34308.3 | 37122.0 |
| 34846.2 | 35045.0 | 34961.2 | 34308.3 | 37102.6 |
El modelo elegido es el ARIMA estacional con deriva. Gana en las tres métricas, en la tabla de ajuste y en la validación con origen móvil, y es el único que puede seguir una serie que crece. Su pronóstico para enero de 2020 es de 35,499 miles de dólares, con un intervalo al 95 por ciento de 34,616 a 36,381.
Los dos ejercicios juntos dejan la lección que importa. En el ejercicio 1 el modelo más simple posible era el correcto, porque no había nada que modelar. Aquí el mismo modelo simple habría producido un pronóstico plano en 34,846 cuando la serie va a llegar a 37,103 en seis meses. No es que un modelo sea mejor que otro en abstracto: hay que ver la serie primero y elegir después.
Dos advertencias sobre las cifras. La primera es que el error esperado es el de la validación con origen móvil, alrededor de 526 miles de dólares al mes, no el 246 del ajuste. La segunda es que los residuos no quedaron completamente limpios, así que los intervalos son optimistas. Con tres años de historia esto es lo que dan los datos; con dos años más el modelo estacional quedaría mucho mejor identificado.
Para Karen Payne, propietaria. Las ventas del restaurante crecen alrededor de 12 mil dólares al año y siguen un calendario muy estable: enero vale 44 por ciento más que un mes promedio y septiembre 37 por ciento menos. Con eso se puede planear el cuarto año mes por mes. Los dos métodos que usé coinciden dentro de un 4 por ciento, y el dato real de enero cayó dentro del intervalo previsto. La recomendación al final del informe no es un número: es dejar de planear con un número solo.
ventas3 <- c(242, 235, 232, 178, 184, 140, 145, 152, 110, 130, 152, 206, # año 1
263, 238, 247, 193, 193, 149, 157, 161, 122, 130, 167, 230, # año 2
282, 255, 265, 205, 210, 160, 166, 174, 126, 148, 173, 235) # año 3
n3 <- length(ventas3)
t3 <- 1:n3
mes3 <- rep(1:12, 3)
ts3 <- ts(ventas3, start = c(1, 1), frequency = 12)Treinta y seis meses de ventas de alimentos y bebidas, en miles de dólares, tomados de la tabla 18.26 del enunciado.
Figura 8. Ventas mensuales del Vintage Restaurant en sus tres primeros años. El mismo perfil se repite cada año, desplazado hacia arriba: eso es tendencia más estación.
Dos cosas saltan a la vista y las dos importan.
La primera es que el perfil dentro del año es prácticamente idéntico en los tres años: máximo en enero, caída sostenida hasta un mínimo en septiembre y recuperación en noviembre y diciembre. Un patrón que se repite con esa fidelidad es estación, no ruido.
La segunda es que cada año está por encima del anterior en casi todos los meses. El promedio pasa de 175.5 en el año 1 a 199.9 en el año 3. El negocio crece.
Un tercer detalle, menos evidente: la distancia entre el máximo y el mínimo de cada año también crece (132 en el año 1 contra 156 en el año 3). Cuando la amplitud de la estación crece junto con el nivel, el modelo que corresponde es el multiplicativo: la estación se expresa como porcentaje del nivel, no como una cantidad fija de dólares.
El procedimiento es el del capítulo 18. Primero se calcula un promedio móvil centrado de orden 12, que al abarcar un año completo deja fuera la estación y sirve como estimación del nivel de la serie en cada mes. Después se divide cada dato entre ese nivel: el cociente mide cuánto está por encima o por debajo ese mes de lo normal. Por último se promedian los cocientes de cada mes de calendario y se ajustan para que sumen 12.
## 1 2 3 4 5 6 7 8 9 10 11
## 1.4436 1.2997 1.3441 1.0412 1.0494 0.8004 0.8283 0.8530 0.6280 0.7003 0.8528
## 12
## 1.1593
## [1] 12
Figura 9. Índice estacional por mes. La línea horizontal en 1.0 es el mes promedio. La temporada alta va de diciembre a marzo y la baja de junio a noviembre, con el mínimo en septiembre; abril y mayo son meses de transición.
¿Tienen sentido intuitivo? Sí, y de forma bastante directa. El restaurante está en la isla Captiva, frente a Fort Myers, Florida, y vende mariscos. Su clientela es turismo de invierno: gente del norte de Estados Unidos y de Canadá que baja a Florida entre diciembre y abril huyendo del frío. Eso explica que enero, febrero y marzo estén los tres por encima de 1.29 y que diciembre cierre en 1.16.
El mínimo de septiembre (0.63, es decir 37 por ciento por debajo del mes promedio) también encaja: es el corazón de la temporada de huracanes en el Golfo de México, el calor y la humedad están en su punto más alto, y en las islas de esa zona buena parte de la oferta turística directamente cierra. Octubre apenas mejora y la recuperación real llega hasta noviembre, con el puente de Acción de Gracias.
La consecuencia para Karen es de operación, no de estadística. La diferencia entre enero y septiembre no es del 10 ni del 20 por ciento: enero vende 2.3 veces lo de septiembre. Personal, compras de pescado y horarios tienen que moverse con ese calendario.
## Estimate Std. Error t value
## (Intercept) 169.349401 1.09261147 154.99508
## t3 1.021297 0.05149666 19.83229
## Pr(>|t|)
## (Intercept) 0.00000000000000000000000000000000000000000000000004866188
## t3 0.00000000000000000002899462975366221870566420272363764779
## [1] 0.9204343
Figura 10. Serie original y serie desestacionalizada con su tendencia lineal. Al quitar el calendario, lo que queda es una recta con muy poca dispersión alrededor: el crecimiento del restaurante es constante, no a saltos.
Sí hay tendencia, y es clara. La recta ajustada es
\[\widehat{D}_t = 169.35 + 1.021\,t\]
con un R cuadrado de 0.92 y una pendiente cuyo valor p es prácticamente cero. En términos de negocio: el restaurante suma 1.02 mil dólares de venta mensual cada mes que pasa, o sea unos 12.3 mil dólares al año. Ese crecimiento es lo que quedaba escondido detrás del calendario en la gráfica original.
El método de descomposición devuelve al pronóstico lo que se le había quitado: se extiende la recta de tendencia a los meses 37 a 48 y cada valor se multiplica por el índice estacional del mes que le toca.
t4 <- 37:48
tendencia4 <- as.numeric(coef(mod3)[1] + coef(mod3)[2] * t4)
pron_desc <- tendencia4 * as.numeric(si3)
round(pron_desc, 1)## [1] 299.0 270.5 281.2 218.9 221.7 169.9 176.6 182.8 135.2 151.5 185.3 253.2
El segundo método llega al mismo lugar por otro camino: una sola regresión con once variables ficticias de mes (enero queda como referencia) y una variable de tiempo para la tendencia.
d3 <- data.frame(y = ventas3, t = t3, mes = factor(mes3))
reg3 <- lm(y ~ t + mes, data = d3)
summary(reg3)$r.squared## [1] 0.9941731
nuevo3 <- data.frame(t = t4, mes = factor(1:12))
pron_reg <- as.numeric(predict(reg3, nuevo3))
ic_reg <- predict(reg3, nuevo3, interval = "prediction", level = 0.95)
round(pron_reg, 1)## [1] 286.8 267.1 272.4 216.4 220.1 174.1 180.4 186.8 143.8 160.4 188.4 248.1
Conviene justificar por qué la variable de tiempo entra al modelo, en lugar de dejar solo las ficticias de mes.
| Res.Df | RSS | Df | Sum of Sq | F | Pr(>F) |
|---|---|---|---|---|---|
| 24 | 4030.667 | ||||
| 23 | 453.625 | 1 | 3577.042 | 181.3656 | 0 |
Agregar la tendencia sube el R cuadrado de 0.948 a 0.994 y la prueba F la rechaza como innecesaria con un valor p prácticamente cero. Sin ella, el modelo pronosticaría el cuarto año con el promedio de los tres primeros y se quedaría corto en todos los meses.
Figura 11. Los tres años observados y el pronóstico del cuarto año por los dos métodos. La banda gris es el intervalo al 95 por ciento de la regresión. El punto rojo es la venta real de enero.
Los dos métodos coinciden bastante: la mayor diferencia entre ellos es de 12.3 mil dólares y en promedio difieren 5.5. Que dos procedimientos distintos lleguen casi al mismo lugar es la mejor señal de que el patrón está en los datos y no en el método.
Comparar el ajuste dentro de la muestra favorece siempre al modelo con más parámetros, así que además reservo el tercer año: estimo los dos métodos con los años 1 y 2 solamente, y les pido que pronostiquen un año que no vieron.
ent <- 1:24; pru <- 25:36
si_e <- indices_estacionales(ventas3[ent], mes3[ent])
des_e <- ventas3[ent] / si_e[mes3[ent]]
mod_e <- lm(des_e ~ ent)
bt_desc <- as.numeric(coef(mod_e)[1] + coef(mod_e)[2] * pru) * as.numeric(si_e[mes3[pru]])
reg_e <- lm(y ~ t + mes, data = data.frame(y = ventas3[ent], t = ent, mes = factor(mes3[ent])))
bt_reg <- as.numeric(predict(reg_e, data.frame(t = pru, mes = factor(mes3[pru]))))
rbind(Descomposicion = metricas(ventas3, replace(rep(NA_real_, n3), pru, bt_desc), pru),
Regresion = metricas(ventas3, replace(rep(NA_real_, n3), pru, bt_reg), pru)) |> round(3)## MAE RMSE MAPE Sesgo
## Descomposicion 0.657 1.021 0.329 -0.231
## Regresion 3.667 5.033 1.878 0.417
La descomposición se equivoca en promedio 0.66 mil dólares al mes contra 3.67 de la regresión. La razón es que la regresión con ficticias estima doce parámetros de mes y una pendiente con solo veinticuatro observaciones, y termina ajustándose a detalles que no se repiten. La descomposición, al promediar los cocientes, es más estable. Me quedo con la descomposición como pronóstico principal y uso la regresión para el intervalo, que es lo que la descomposición no da.
real_enero <- 295
c(pronostico_descomposicion = pron_desc[1],
error_descomposicion = real_enero - pron_desc[1],
error_porcentual_desc = 100 * (real_enero - pron_desc[1]) / real_enero,
pronostico_regresion = pron_reg[1],
error_regresion = real_enero - pron_reg[1],
error_porcentual_reg = 100 * (real_enero - pron_reg[1]) / real_enero)## pronostico_descomposicion error_descomposicion error_porcentual_desc
## 299.020918 -4.020918 -1.363023
## pronostico_regresion error_regresion error_porcentual_reg
## 286.750000 8.250000 2.796610
## fit lwr upr
## 286.7500 275.4983 298.0017
El error es de -4 mil dólares por descomposición (1.4 por ciento) y de 8.2 mil por regresión (2.8 por ciento).
¿Es un error grande? No. Sobre una venta de 295 mil dólares, equivocarse por 4 mil o por 8 mil es un desempeño bueno para un pronóstico hecho con un año de anticipación. Y hay un dato que lo confirma: el intervalo al 95 por ciento de la regresión para enero iba de 275.5 a 298, y 295 cae dentro. El procedimiento no falló; funcionó exactamente como debía.
Qué hacer para resolver la incertidumbre. La pregunta del enunciado tiene una respuesta concreta, y son cuatro cosas.
Primero, entregar intervalos y no puntos. La confusión de Karen viene de comparar 295 contra un número solo. Si el reporte dice que enero va a caer entre 276 y 298 con 95 por ciento de confianza, el 295 deja de ser una sorpresa y se convierte en lo esperado. Un pronóstico de un valor exacto está mal comunicado antes que mal calculado.
Segundo, decir de antemano cuánto error es normal. La prueba fuera de muestra de la sección anterior da la cifra: alrededor de 0.3 por ciento con la descomposición. Con eso Karen tiene una regla para saber cuándo preocuparse: una desviación de 3 por ciento es ruido, una de 15 por ciento es una señal de que algo cambió en el negocio.
Tercero, reestimar cada mes. El modelo se estimó con 36 datos. Al cerrar enero del cuarto año hay 37, y volver a correrlo incorpora esa información al resto del año. Es lo mismo que hice con el origen móvil en los ejercicios 1 y 2, aplicado a la operación.
Cuarto, separar el error de calendario del error de nivel. Si los doce meses del cuarto año se desvían todos hacia arriba, el problema es la tendencia y hay que corregir la pendiente. Si solo se desvía septiembre, el problema es el índice de ese mes. Son dos arreglos distintos y confundirlos lleva a corregir lo que no estaba roto.
Para la dirección de Carlson Department Store. Estimo que la tienda habría vendido 12.36 millones de dólares entre septiembre y diciembre si el huracán no hubiera ocurrido, y esa es la pérdida que se debe reclamar. Además, las tiendas departamentales del condado vendieron en esos cuatro meses 28 por ciento más de lo que les correspondía: hubo un empuje real de actividad económica. Trasladado a Carlson, ese empuje vale 3.51 millones adicionales. El reclamo completo es de 15.87 millones, con dos conceptos que hay que argumentar por separado porque la aseguradora los va a tratar distinto.
Carlson cerró de septiembre a diciembre del quinto año por los daños del huracán del 21 de agosto. Hay dos preguntas y son distintas. La primera es cuánto habría vendido la tienda en esos cuatro meses en un año normal, que es la pérdida directa. La segunda es si tiene derecho a algo más, porque los 8 mil millones de dólares de ayuda federal y pagos de seguros que entraron al condado dispararon el comercio de la zona justo cuando Carlson estaba cerrada.
La clave del caso es que las ventas de todas las tiendas departamentales del condado funcionan como grupo de control. Sirven para medir cuánto se movió el mercado por el huracán, algo que los datos de Carlson por sí solos no pueden decir porque a Carlson le faltan justamente esos cuatro meses.
carlson <- c(1.71, 1.90, 2.74, 4.20, # año 1: sep-dic
1.45, 1.80, 2.03, 1.99, 2.32, 2.20, 2.13, 2.43, # año 2: ene-ago
1.90, 2.13, 2.56, 4.16, # año 2: sep-dic
2.31, 1.89, 2.02, 2.23, 2.39, 2.14, 2.27, 2.21,
1.89, 2.29, 2.83, 4.04,
2.31, 1.99, 2.42, 2.45, 2.57, 2.42, 2.40, 2.50,
2.09, 2.54, 2.97, 4.35,
2.56, 2.28, 2.69, 2.48, 2.73, 2.37, 2.31, 2.23) # año 5: ene-ago
condado <- c(55.80, 56.40, 71.40, 117.60,
46.80, 48.00, 60.00, 57.60, 61.80, 58.20, 56.40, 63.00,
57.60, 53.40, 71.40, 114.00,
46.80, 48.60, 59.40, 58.20, 60.60, 55.20, 51.00, 58.80,
49.80, 54.60, 65.40, 102.00,
43.80, 45.60, 57.60, 53.40, 56.40, 52.80, 54.00, 60.60,
47.40, 54.60, 67.80, 100.20,
48.00, 51.60, 57.60, 58.20, 60.00, 57.00, 57.60, 61.80)
# ventas reales del condado en los cuatro meses en que Carlson estuvo cerrada
condado_post <- c(69.00, 75.00, 85.20, 121.80)
n4 <- length(carlson) # 48 meses: septiembre del año 1 a agosto del año 5
t4c <- 1:n4
mes4 <- ((8 + t4c - 1) %% 12) + 1 # t = 1 corresponde a septiembre
c(meses_carlson = length(carlson), meses_condado = length(condado),
primer_mes = mes4[1], ultimo_mes = mes4[n4])## meses_carlson meses_condado primer_mes ultimo_mes
## 48 48 9 8
Las dos series son de 48 meses, de septiembre del año 1 a agosto del año 5, en millones de dólares. Vienen de las tablas 18.27 y 18.28 del enunciado.
Figura 12. Ventas mensuales de Carlson y del conjunto de tiendas departamentales del condado. Las dos escalas son distintas, así que cada serie va contra su propio eje. El sombreado marca los cuatro meses en que Carlson estuvo cerrada.
Las dos series comparten el pico de diciembre, que es lo esperable en tiendas departamentales. Lo que no comparten es la dirección: Carlson sube a lo largo de los cuatro años y el condado no. Ese detalle va a importar más adelante.
Mismo procedimiento del ejercicio 3: promedio móvil centrado de orden 12, índices estacionales, serie desestacionalizada, tendencia lineal y regreso a la escala original.
si_car <- indices_estacionales(carlson, mes4)
des_car <- carlson / si_car[mes4]
mod_car <- lm(des_car ~ t4c)
summary(mod_car)$coefficients## Estimate Std. Error t value
## (Intercept) 2.14908425 0.048946014 43.907237
## t4c 0.01140723 0.001739038 6.559504
## Pr(>|t|)
## (Intercept) 0.000000000000000000000000000000000000003344518
## t4c 0.000000041812372263504978665009826233750489521
## [1] 0.4833034
t_cierre <- 49:52 # septiembre a diciembre del año 5
mes_cierre <- c(9, 10, 11, 12)
pron_car <- as.numeric(coef(mod_car)[1] + coef(mod_car)[2] * t_cierre) *
as.numeric(si_car[mes_cierre])
round(pron_car, 3)## [1] 2.158 2.544 3.056 4.600
## [1] 12.3579
Diciembre vale 1.68 veces un mes promedio: la tienda hace en ese solo mes casi el 14 por ciento de su año. Por eso la pérdida de septiembre a diciembre no es un tercio del año, es bastante más.
La tendencia desestacionalizada crece 0.0114 millones de dólares por mes, con valor p de 0.0000. El R cuadrado es 0.483, que no es alto, pero eso era de esperarse: la serie desestacionalizada de una tienda tiene mucho ruido mes a mes y lo que se le pide a la recta es solo capturar el nivel promedio y su desplazamiento.
si_con <- indices_estacionales(condado, mes4)
des_con <- condado / si_con[mes4]
mod_con <- lm(des_con ~ t4c)
summary(mod_con)$coefficients## Estimate Std. Error t value
## (Intercept) 62.22039143 0.74998888 82.961752
## t4c -0.07194144 0.02664689 -2.699806
## Pr(>|t|)
## (Intercept) 0.0000000000000000000000000000000000000000000000000009508138
## t4c 0.0096749104122774258474048991729432600550353527069091796875
pron_con <- as.numeric(coef(mod_con)[1] + coef(mod_con)[2] * t_cierre) *
as.numeric(si_con[mes_cierre])
round(pron_con, 2)## [1] 50.55 53.20 66.78 103.11
## [1] 273.6424
Aquí aparece algo que hay que decir con todas sus letras: la tendencia del condado es negativa, -0.0719 millones por mes, con valor p de 0.0097. Mientras el comercio departamental del condado se contraía, Carlson crecía. La tienda venía ganando participación de mercado.
Esto tiene una consecuencia directa para el reclamo. Si la aseguradora propusiera estimar las ventas perdidas de Carlson aplicándole la evolución del condado, el resultado subestimaría la pérdida, porque le impondría a Carlson una caída que Carlson no tenía. Cada serie se proyecta con su propia tendencia; el condado se usa únicamente para medir el efecto del huracán, que es para lo que sirve un grupo de control.
La estimación de la pérdida directa es de 12.36 millones de dólares. Diciembre solo aporta 4.60 millones, el 37 por ciento del total: haber estado cerrada en diciembre es lo más caro que le pasó a la tienda.
Ahora sí entra el grupo de control. Comparo lo que el condado vendió de verdad en esos cuatro meses contra lo que le tocaba vender según su propia tendencia y su propio calendario.
## [1] 1.3650 1.4099 1.2758 1.1812
## [1] 77.35756
## [1] 28.26958
El condado vendió 77.4 millones más de lo que le correspondía, un 28 por ciento por encima. El efecto es grande y tiene la forma que uno esperaría: el pico está en septiembre y octubre (factores de 1.37 y 1.41), justo cuando llegaron los pagos de seguros y la ayuda federal, y se va desvaneciendo hacia diciembre.
Figura 14. Ventas del condado en los cuatro meses de cierre: lo esperado sin huracán contra lo que ocurrió de verdad. La diferencia es la actividad extra que generó la reconstrucción, y es la que Carlson no pudo aprovechar.
El argumento para Carlson se construye trasladándole ese mismo empuje. Si la tienda hubiera estado abierta, no habría vendido lo normal: habría vendido lo normal multiplicado por el factor del mercado en el que opera.
## [1] 2.945 3.587 3.899 5.433
c(perdida_normal = sum(pron_car),
venta_con_huracan = sum(car_con_huracan),
exceso = sum(car_con_huracan) - sum(pron_car),
reclamo_total = sum(car_con_huracan))## perdida_normal venta_con_huracan exceso reclamo_total
## 12.357903 15.864704 3.506802 15.864704
Todo el reclamo descansa en una proyección, así que hay que decir qué tan buena es. La prueba: estimar el método con los primeros 36 meses y pedirle que pronostique los 12 meses siguientes, que sí conozco.
e4 <- 1:36; p4 <- 37:48
si_e4 <- indices_estacionales(carlson[e4], mes4[e4])
des_e4 <- carlson[e4] / si_e4[mes4[e4]]
mod_e4 <- lm(des_e4 ~ e4)
bt4 <- as.numeric(coef(mod_e4)[1] + coef(mod_e4)[2] * p4) * as.numeric(si_e4[mes4[p4]])
c(MAE = mean(abs(carlson[p4] - bt4)),
MAPE = 100 * mean(abs((carlson[p4] - bt4) / carlson[p4])),
suma_real = sum(carlson[p4]),
suma_estimada = sum(bt4),
error_del_total = 100 * (sum(bt4) / sum(carlson[p4]) - 1))## MAE MAPE suma_real suma_estimada error_del_total
## 0.145228 5.865063 31.600000 31.857335 0.814351
Mes por mes el método se equivoca alrededor de 5.9 por ciento, que no es poco. Pero sobre el total de los doce meses el error es de 0.8 por ciento, porque los errores de meses distintos se cancelan entre sí. Y el reclamo es precisamente un total, no un pronóstico mes por mes. Esa es la razón por la que la estimación de 12.36 millones se puede defender ante la aseguradora aunque el modelo no sea muy preciso en cualquier mes tomado por separado.
## 1 2 3 4 5 6 7 8 9 10
## -0.1560 0.0056 0.0902 0.0156 0.0175 -0.0331 -0.0459 -0.0873 0.0158 -0.0099
## 11 12
## 0.0911 0.0964
##
## Box-Ljung test
##
## data: res_car
## X-squared = 4.3518, df = 12, p-value = 0.9762
Los residuos de la serie desestacionalizada no muestran estación remanente ni autocorrelación (Ljung-Box con valor p de 0.976). El modelo extrajo lo que había que extraer.
Una nota que juega a favor de Carlson y conviene incluir en el expediente: la tendencia de la tienda es positiva mientras la del condado es negativa. Carlson venía ganando terreno en un mercado que se encogía, así que cualquier estimación basada en el comportamiento del condado en lugar del propio subestima la pérdida.
Primero se mira la serie y después se elige el modelo. Los ejercicios 1 y 2 aplican exactamente los mismos cinco métodos y llegan a conclusiones opuestas. En el primero gana el que menos supone; en el segundo, el único que puede seguir una tendencia. No existe un modelo bueno en abstracto: existe un modelo que corresponde a lo que la serie tiene adentro.
El sesgo dice más que el error promedio. En el ejercicio 2 los cuatro modelos de nivel tienen un MAE que parece tolerable, pero se equivocan siempre hacia el mismo lado. Para una empresa que planea producción con ese pronóstico, un error que se repite en la misma dirección durante treinta y seis meses no es un error: es un faltante acumulado. El MAE lo esconde y el sesgo lo enseña.
Un modelo se juzga con datos que no vio. Las tres validaciones de este documento cambiaron alguna conclusión. En el ejercicio 2 el origen móvil duplicó el error esperado del ARIMA. En el ejercicio 3 la prueba fuera de muestra invirtió el orden entre descomposición y regresión respecto a lo que sugería el ajuste. En el ejercicio 4 la misma prueba mostró que el método es mediocre mes a mes pero muy bueno en el total, que es justo lo que el caso necesitaba.
El intervalo es el entregable, no el punto. Aparece en los cuatro ejercicios y de forma más clara en el 3, donde Karen se confunde por comparar 295 contra un número solo cuando el rango previsto lo contenía sin problema. Un pronóstico entregado como un valor exacto está mal comunicado antes de estar mal calculado.
Los cuatro ejercicios son series de tiempo, pero se apoyan en la misma idea que el tema de datos de panel: hay que separar la variación que corresponde a la unidad de la que corresponde al tiempo. Ahí eran efectos fijos de empresa y de año; aquí son índice estacional y tendencia. El índice estacional de septiembre cumple el mismo papel que un efecto fijo de mes, y la tendencia el de un efecto de tiempo. Cambia la notación, no el razonamiento.
El ejercicio 4 agrega algo que no aparecía en la Actividad 1 y que vale la pena señalar: el uso de una segunda serie como grupo de control. Las ventas del condado no entran al modelo de Carlson como variable explicativa; sirven para medir un efecto que en la serie de Carlson no se puede observar, porque justo en esos meses no hay datos. Es el mismo razonamiento de un experimento con grupo de comparación, aplicado a series de tiempo.
Anderson, D. R., Sweeney, D. J. y Williams, T. A. Estadística para administración y economía. Capítulo 18, “Series de tiempo y elaboración de pronósticos”. Material del curso, tablas 18.26, 18.27 y 18.28.
Cantú, R. A. (2026). Tema 2. Modelos de pronósticos en series de tiempo. Notas y sesión de clase, Generación de escenarios futuros con analítica (CD-3001B), Tecnológico de Monterrey.
Hershey México. Ventas históricas de leche saborizada,
2017-2019. Archivo Ventas_Históricas_Lechitas.xlsx
proporcionado en la actividad.
Hyndman, R. J. y Khandakar, Y. (2008). Automatic Time Series
Forecasting: The forecast Package for R. Journal of
Statistical Software, 27(3), 1-22.
Ljung, G. M. y Box, G. E. P. (1978). On a Measure of Lack of Fit in Time Series Models. Biometrika, 65(2), 297-303.
Todo el documento se compila desde este archivo. Las series de los
ejercicios 1, 3 y 4 se escriben directamente en el código, porque vienen
de las tablas del enunciado y son cortas. La del ejercicio 2 se lee de
datos/lechitas.csv.
Ese archivo se obtuvo del original
Ventas_Históricas_Lechitas.xlsx, que trae los tres años en
tres pares de columnas y las fechas como número de serie de Excel. La
conversión consistió en apilar los tres bloques en un solo vector
ordenado por fecha y traducir el número de serie a año y mes. No se
modificó ningún valor: los 36 montos son los del archivo original con
todos sus decimales.
| Componente | Versión |
|---|---|
| R | 4.6.1 |
| forecast | 9.0.2 |
| knitr | 1.51 |
| kableExtra | 1.4.1 |