CA-0415 Series de Tiempo — Semana 8
Al finalizar este capítulo, se espera que el estudiantado pueda:
fable, sin tratar la selección automática como sustituto de la interpretación estadística;En Sección 7.19 establecimos que una comparación predictiva responsable debe fijar el origen, el horizonte, los métodos de referencia y la información disponible antes de observar el período de evaluación. Esta semana aplicaremos ese protocolo a una familia de métodos cuya lógica es distinta de ARIMA.
En un ARIMA describimos la dinámica mediante rezagos de la serie y de sus innovaciones. En suavizamiento exponencial adoptaremos otra perspectiva: resumiremos el pasado mediante componentes que evolucionan en el tiempo, como un nivel local, una pendiente local y un patrón estacional. La guía bibliográfica del curso asigna esta semana a suavizamiento exponencial y modelos ETS, con Hyndman y Athanasopoulos como texto conductor (Hyndman y Athanasopoulos 2021, chap. 8). Petris, Petrone y Campagnoli proporcionan un complemento importante al mostrar cómo los algoritmos clásicos de suavizamiento pueden conectarse con modelos dinámicos (Petris et al. 2009, sec. 3.1).
El hilo conductor será
\[ \boxed{ \text{promedio con memoria decreciente} \longrightarrow \text{nivel} \longrightarrow \text{tendencia} \longrightarrow \text{estacionalidad} \longrightarrow \text{modelo ETS} \longrightarrow \text{distribución predictiva} \longrightarrow \text{evaluación temporal} }. \]
Suponga que una serie no presenta una tendencia sistemática ni un patrón estacional claro. Dos métodos de referencia de Sección 7.4 representan extremos conceptuales.
El método naïve utiliza solamente la observación más reciente:
\[ \widehat y_{T+h\mid T}=y_T, \qquad h\geq 1. \tag{8.1}\]
El método de la media utiliza todas las observaciones con el mismo peso:
\[ \widehat y_{T+h\mid T} = \bar y_T = \frac{1}{T}\sum_{t=1}^{T}y_t. \tag{8.2}\]
Ambas reglas pueden ser inadecuadas cuando la estructura cambia gradualmente. El primer método descarta casi todo el pasado; el segundo conserva observaciones antiguas con el mismo peso que observaciones recientes.
El suavizamiento exponencial propone una solución intermedia:
conservar todo el pasado, pero reducir geométricamente su influencia conforme envejece.
Esta idea es fundamental. El adjetivo exponencial no se refiere a que la serie observada tenga que crecer exponencialmente. Se refiere a la forma de los pesos asignados a observaciones pasadas (Hyndman y Athanasopoulos 2021, sec. 8.1).
El suavizamiento exponencial simple (SES, simple exponential smoothing) está diseñado para series cuyo comportamiento puede resumirse mediante un nivel que evoluciona, sin una tendencia o estacionalidad que deban extrapolarse explícitamente.
Sea \(\ell_t\) el nivel estimado después de observar \(y_t\). La ecuación de actualización es
\[ \ell_t = \alpha y_t+(1-\alpha)\ell_{t-1}, \qquad 0\leq\alpha\leq1. \tag{8.3}\]
El pronóstico para cualquier horizonte futuro es
\[ \widehat y_{t+h\mid t}=\ell_t, \qquad h=1,2,\ldots \tag{8.4}\]
Por tanto, SES genera una función de pronóstico plana. Esto no significa que el nivel histórico haya sido constante: \(\ell_t\) puede cambiar cada vez que llega una nueva observación. Significa que, una vez fijado el origen \(t\), no extrapolamos una pendiente ni un patrón estacional hacia \(t+h\).
El parámetro \(\alpha\) controla la rapidez de adaptación:
La interpretación más importante se obtiene expandiendo recursivamente Ecuación 8.3. Primero,
\[ \ell_t = \alpha y_t+(1-\alpha)\ell_{t-1}. \]
Pero
\[ \ell_{t-1} = \alpha y_{t-1}+(1-\alpha)\ell_{t-2}, \]
por lo que
\[ \ell_t = \alpha y_t + \alpha(1-\alpha)y_{t-1} + (1-\alpha)^2\ell_{t-2}. \]
Continuando hasta el estado inicial \(\ell_0\) obtenemos
\[ \ell_t = \sum_{j=0}^{t-1} \alpha(1-\alpha)^j y_{t-j} + (1-\alpha)^t\ell_0. \tag{8.5}\]
El peso de \(y_{t-j}\) es
\[ w_j=\alpha(1-\alpha)^j. \tag{8.6}\]
Como \(0<1-\alpha<1\) cuando \(0<\alpha<1\), los pesos decaen geométricamente con \(j\). De aquí proviene el nombre del método.
Además,
\[ \sum_{j=0}^{t-1}\alpha(1-\alpha)^j = 1-(1-\alpha)^t. \tag{8.7}\]
Al agregar el peso \((1-\alpha)^t\) asignado al estado inicial \(\ell_0\), los pesos de Ecuación 8.5 suman exactamente uno.
En muestras cortas no debemos omitir el término \((1-\alpha)^t\ell_0\). La representación como promedio infinito es una aproximación útil cuando \(t\) es grande y \((1-\alpha)^t\) es despreciable, pero la recursión real debe comenzar con algún valor inicial.
El siguiente gráfico muestra cuánto peso recibe una observación situada \(j\) períodos en el pasado.
alphas_pesos <- c(0.1, 0.3, 0.6, 0.9)
pesos_ets <- crossing(
alpha = alphas_pesos,
rezago = 0:30
) |>
mutate(
peso = alpha * (1 - alpha)^rezago,
alpha = factor(alpha)
)
ggplot(pesos_ets, aes(x = rezago, y = peso, linetype = alpha)) +
geom_line(linewidth = 0.8) +
labs(
x = "Antigüedad j",
y = expression(alpha * (1 - alpha)^j),
linetype = expression(alpha)
)
Para hacer visible el compromiso entre estabilidad y adaptación, simularemos una serie cuyo nivel cambia abruptamente. No suponemos que SES sea el modelo verdadero; la simulación se utiliza únicamente para estudiar la recursión.
n_nivel <- 90
nivel_verdadero <- c(rep(20, 45), rep(28, 45))
y_nivel <- nivel_verdadero + rnorm(n_nivel, sd = 1.5)
ses_manual <- function(y, alpha, l0 = y[1]) {
n <- length(y)
nivel <- numeric(n)
nivel[1] <- alpha * y[1] + (1 - alpha) * l0
if (n >= 2) {
for (t in 2:n) {
nivel[t] <- alpha * y[t] + (1 - alpha) * nivel[t - 1]
}
}
nivel
}
niveles_alpha <- tibble(
t = seq_len(n_nivel),
observado = y_nivel,
verdadero = nivel_verdadero,
`alpha = 0.10` = ses_manual(y_nivel, 0.10),
`alpha = 0.35` = ses_manual(y_nivel, 0.35),
`alpha = 0.80` = ses_manual(y_nivel, 0.80)
) |>
pivot_longer(
cols = starts_with("alpha"),
names_to = "alpha",
values_to = "nivel"
)niveles_alpha |>
ggplot(aes(x = t)) +
geom_point(aes(y = observado), size = 0.8, alpha = 0.55) +
geom_line(aes(y = verdadero), linewidth = 0.6, linetype = 2) +
geom_line(aes(y = nivel, linetype = alpha), linewidth = 0.8) +
labs(
x = "Tiempo",
y = "Nivel",
linetype = NULL
)
Un valor grande de \(\alpha\) no es universalmente mejor. Permite adaptación rápida, pero también transmite con mayor intensidad el ruido de cada observación al estado actualizado. El valor apropiado debe determinarse en función de los datos y del objetivo predictivo.
La recursión de SES puede reescribirse de una forma que será esencial más adelante. El pronóstico de un paso antes de observar \(y_t\) es
\[ \widehat y_{t\mid t-1}=\ell_{t-1}. \]
Defina el error de pronóstico o innovación
\[ e_t = y_t-\widehat y_{t\mid t-1} = y_t-\ell_{t-1}. \tag{8.8}\]
Entonces
\[ \begin{aligned} \ell_t &= \alpha y_t+(1-\alpha)\ell_{t-1}\\ &= \ell_{t-1} + \alpha(y_t-\ell_{t-1})\\ &= \ell_{t-1}+\alpha e_t. \end{aligned} \tag{8.9}\]
Esta representación permite leer la actualización como
\[ \boxed{ \text{nuevo estado} = \text{estado anterior} + \text{ganancia}\times\text{nueva información} }. \]
La “nueva información” es la parte de \(y_t\) que no estaba anticipada por \(\ell_{t-1}\). Esta estructura aparecerá nuevamente en los modelos de espacio-estado y en el filtro de Kalman, aunque todavía no desarrollaremos formalmente ese filtro.
\(\ell_t\) es un estado actualizado secuencialmente. No es el valor de una curva global estimada una sola vez utilizando todos los tiempos. La distinción será importante cuando estudiemos filtrado y suavizamiento retrospectivo en modelos de espacio-estado.
Para ejecutar SES debemos conocer al menos
\[ \alpha \quad\text{y}\quad \ell_0. \]
Una forma clásica de estimarlos es minimizar la suma de cuadrados de los errores de un paso:
\[ \operatorname{SSE}(\alpha,\ell_0) = \sum_{t=1}^{T} \left(y_t-\widehat y_{t\mid t-1}\right)^2 = \sum_{t=1}^{T}e_t^2. \tag{8.10}\]
La dependencia de \(e_t\) respecto de \(\alpha\) y \(\ell_0\) es recursiva, por lo que el problema se resuelve numéricamente (Hyndman y Athanasopoulos 2021, sec. 8.1).
Dos mensajes son importantes:
Más adelante reemplazaremos la minimización de SSE por la estimación mediante verosimilitud dentro de un modelo probabilístico ETS.
SES produce pronósticos planos. Si una serie muestra una tendencia persistente durante el intervalo de interés, necesitamos un segundo estado que describa la pendiente local.
El método lineal de Holt utiliza dos componentes:
Sus ecuaciones son (Hyndman y Athanasopoulos 2021, sec. 8.2; Petris et al. 2009, sec. 3.1.1)
\[ \widehat y_{t+h\mid t} = \ell_t+h b_t, \tag{8.11}\]
\[ \ell_t = \alpha y_t + (1-\alpha)(\ell_{t-1}+b_{t-1}), \tag{8.12}\]
\[ b_t = \beta^*(\ell_t-\ell_{t-1}) + (1-\beta^*)b_{t-1}, \tag{8.13}\]
con \(0\leq\alpha\leq1\) y \(0\leq\beta^*\leq1\) en la parametrización clásica.
Antes de observar \(y_t\), el nivel esperado a partir del estado anterior sería
\[ \ell_{t-1}+b_{t-1}. \]
La ecuación Ecuación 8.12 combina ese valor extrapolado con la observación nueva. Por ello, \(\ell_t\) no es simplemente una media suavizada de \(y_t\); incorpora la pendiente que se esperaba entre \(t-1\) y \(t\).
La cantidad
\[ \ell_t-\ell_{t-1} \]
es el cambio reciente del nivel estimado. Ecuación 8.13 combina ese cambio con la pendiente previa \(b_{t-1}\).
Por tanto, \(b_t\) es una pendiente que evoluciona. No debe interpretarse como el coeficiente fijo de una regresión lineal global.
Utilizamos la serie WWWusage, incluida en R, para comparar SES, Holt y Holt con tendencia amortiguada. La serie no es estacional y presenta cambios claros en su nivel y pendiente.
www_usage <- as_tsibble(WWWusage)
stopifnot(
nrow(www_usage) == length(WWWusage),
!anyNA(www_usage$value),
all(is.finite(www_usage$value))
)www_usage |>
autoplot(value) +
labs(
x = "Minuto",
y = "Usuarios"
)
Ajustamos tres especificaciones de error aditivo para aislar el efecto del componente de tendencia.
fit_www <- www_usage |>
model(
SES = ETS(value ~ error("A") + trend("N") + season("N")),
Holt = ETS(value ~ error("A") + trend("A") + season("N")),
Amortiguado = ETS(value ~ error("A") + trend("Ad") + season("N"))
)
fc_www <- fit_www |>
forecast(h = 15)fc_www |>
autoplot(www_usage, level = NULL) +
labs(
x = "Minuto",
y = "Usuarios",
colour = "Método"
)
No debemos seleccionar el método mirando solamente cuál extensión gráfica parece más plausible. Siguiendo Sección 7.8, podemos comparar pronósticos de un paso con orígenes sucesivos.
www_cv <- www_usage |>
stretch_tsibble(.init = 30, .step = 2)
acc_www <- www_cv |>
model(
SES = ETS(value ~ error("A") + trend("N") + season("N")),
Holt = ETS(value ~ error("A") + trend("A") + season("N")),
Amortiguado = ETS(value ~ error("A") + trend("Ad") + season("N"))
) |>
forecast(h = 1) |>
accuracy(www_usage) |>
select(.model, RMSE, MAE, MASE) |>
arrange(RMSE)knitr::kable(
acc_www,
digits = 3,
caption = "Evaluación de un paso con origen móvil para tres métodos de suavizamiento no estacionales."
)| .model | RMSE | MAE | MASE |
|---|---|---|---|
| Amortiguado | 3.064 | 2.563 | 0.566 |
| Holt | 3.286 | 2.850 | 0.630 |
| SES | 5.685 | 4.600 | 1.017 |
La tabla debe leerse como evidencia específica para esta serie, este horizonte y este diseño de validación. No establece un orden universal entre los tres métodos.
El pronóstico de Holt en Ecuación 8.11 mantiene la pendiente \(b_t\) indefinidamente. Para un horizonte grande,
\[ \widehat y_{t+h\mid t}=\ell_t+h b_t \]
puede producir extrapolaciones poco plausibles, especialmente cuando \(b_t\) representa solamente una tendencia local reciente.
La tendencia amortiguada introduce un parámetro
\[ 0<\phi<1 \]
y reemplaza \(h b_t\) por una suma geométricamente decreciente (Hyndman y Athanasopoulos 2021, sec. 8.2):
\[ \widehat y_{t+h\mid t} = \ell_t + (\phi+\phi^2+\cdots+\phi^h)b_t. \tag{8.14}\]
Las actualizaciones se escriben como
\[ \ell_t = \alpha y_t + (1-\alpha)(\ell_{t-1}+\phi b_{t-1}), \tag{8.15}\]
\[ b_t = \beta^*(\ell_t-\ell_{t-1}) + (1-\beta^*)\phi b_{t-1}. \tag{8.16}\]
Si \(\phi=1\), recuperamos el método lineal de Holt. Cuando \(\phi<1\), la contribución marginal de la tendencia al pronóstico disminuye con \(h\).
Defina
\[ \phi_h = \sum_{j=1}^{h}\phi^j. \]
Para \(\phi\neq1\),
\[ \phi_h = \frac{\phi(1-\phi^h)}{1-\phi}. \tag{8.17}\]
Como \(0<\phi<1\),
\[ \lim_{h\to\infty}\phi^h=0, \]
y por tanto
\[ \lim_{h\to\infty}\widehat y_{t+h\mid t} = \ell_t + \frac{\phi}{1-\phi}b_t. \tag{8.18}\]
El método amortiguado no obliga a que la serie futura se vuelva constante en realidad. Lo que hace es evitar que la pendiente estimada hoy se extrapole sin reducción para siempre.
h_damp <- 1:40
l0_damp <- 100
b0_damp <- 2
trayectorias_damp <- crossing(
h = h_damp,
phi = c(0.80, 0.90, 0.98, 1.00)
) |>
group_by(phi) |>
mutate(
contribucion = if_else(
phi == 1,
h,
phi * (1 - phi^h) / (1 - phi)
),
pronostico = l0_damp + contribucion * b0_damp,
phi = factor(phi)
) |>
ungroup()
ggplot(trayectorias_damp, aes(x = h, y = pronostico, linetype = phi)) +
geom_line(linewidth = 0.8) +
labs(
x = "Horizonte h",
y = "Pronóstico",
linetype = expression(phi)
)
Sea \(m\) el período estacional. Para datos mensuales con estacionalidad anual, \(m=12\); para datos trimestrales, \(m=4\).
El método de Holt–Winters agrega un estado \(s_t\) que describe el patrón estacional. Existen dos formas principales.
La formulación aditiva es adecuada cuando las oscilaciones estacionales tienen una amplitud aproximadamente constante en las unidades originales de la serie (Hyndman y Athanasopoulos 2021, sec. 8.3).
Las ecuaciones son
\[ \widehat y_{t+h\mid t} = \ell_t+h b_t+s_{t+h-m(k+1)}, \tag{8.19}\]
\[ \ell_t = \alpha(y_t-s_{t-m}) + (1-\alpha)(\ell_{t-1}+b_{t-1}), \tag{8.20}\]
\[ b_t = \beta^*(\ell_t-\ell_{t-1}) + (1-\beta^*)b_{t-1}, \tag{8.21}\]
\[ s_t = \gamma(y_t-\ell_{t-1}-b_{t-1}) + (1-\gamma)s_{t-m}, \tag{8.22}\]
con
\[ k=\left\lfloor\frac{h-1}{m}\right\rfloor. \]
El índice \(k\) asegura que, para horizontes superiores a un ciclo estacional, reutilicemos la posición estacional correspondiente del último ciclo estimado.
La ecuación de nivel utiliza una observación desestacionalizada,
\[ y_t-s_{t-m}, \]
y la ecuación estacional compara la observación nueva con el nivel y tendencia que habrían sido esperados sin estacionalidad.
Cuando la amplitud estacional crece o disminuye aproximadamente en proporción con el nivel, una representación multiplicativa puede ser más natural:
\[ \widehat y_{t+h\mid t} = (\ell_t+h b_t)s_{t+h-m(k+1)}, \tag{8.23}\]
\[ \ell_t = \alpha\frac{y_t}{s_{t-m}} + (1-\alpha)(\ell_{t-1}+b_{t-1}), \tag{8.24}\]
\[ b_t = \beta^*(\ell_t-\ell_{t-1}) + (1-\beta^*)b_{t-1}, \tag{8.25}\]
\[ s_t = \gamma\frac{y_t}{\ell_{t-1}+b_{t-1}} + (1-\gamma)s_{t-m}. \tag{8.26}\]
Aquí \(s_t\) es un factor relativo, no una desviación en las unidades de \(y_t\).
Una formulación multiplicativa no es adecuada cuando los datos pueden ser negativos y puede ser problemática cuando contienen ceros. Además, “la amplitud parece aumentar” no basta para elegirla: una transformación como el logaritmo o Box–Cox puede ofrecer otra representación razonable. La decisión debe considerar interpretación, ajuste y desempeño predictivo.
La distinción se refiere a cómo cambia la amplitud estacional con el nivel.
n_season <- 72
t_season <- seq_len(n_season)
patron <- rep(c(-1.0, -0.3, 0.4, 0.9, 0.5, -0.5), length.out = n_season)
base_season <- 20 + 0.20 * t_season
season_sim <- tibble(
t = t_season,
Aditiva = base_season + 3.2 * patron + rnorm(n_season, sd = 0.5),
Multiplicativa = base_season * (1 + 0.16 * patron) + rnorm(n_season, sd = 0.5)
) |>
pivot_longer(-t, names_to = "Estructura", values_to = "y")
ggplot(season_sim, aes(x = t, y = y)) +
geom_line() +
facet_wrap(~ Estructura, ncol = 1, scales = "free_y") +
labs(
x = "Tiempo",
y = "Valor"
)
Antes de introducir ETS conviene separar dos clasificaciones.
La primera clasifica métodos de pronóstico según tendencia y estacionalidad (Hyndman y Athanasopoulos 2021, sec. 8.4). Utilizaremos:
N: ninguna;A: aditiva;Ad: tendencia aditiva amortiguada;M: multiplicativa, cuando corresponde a estacionalidad.metodos_ets <- tribble(
~Tendencia, ~`Sin estacionalidad (N)`, ~`Estacionalidad aditiva (A)`, ~`Estacionalidad multiplicativa (M)`,
"Ninguna (N)", "SES: (N,N)", "(N,A)", "(N,M)",
"Aditiva (A)", "Holt: (A,N)", "Holt-Winters aditivo: (A,A)", "Holt-Winters multiplicativo: (A,M)",
"Aditiva amortiguada (Ad)", "Holt amortiguado: (Ad,N)", "(Ad,A)", "Holt-Winters amortiguado: (Ad,M)"
)
knitr::kable(
metodos_ets,
caption = "Taxonomía de métodos de suavizamiento exponencial según tendencia y estacionalidad."
)| Tendencia | Sin estacionalidad (N) | Estacionalidad aditiva (A) | Estacionalidad multiplicativa (M) |
|---|---|---|---|
| Ninguna (N) | SES: (N,N) | (N,A) | (N,M) |
| Aditiva (A) | Holt: (A,N) | Holt-Winters aditivo: (A,A) | Holt-Winters multiplicativo: (A,M) |
| Aditiva amortiguada (Ad) | Holt amortiguado: (Ad,N) | (Ad,A) | Holt-Winters amortiguado: (Ad,M) |
Esta clasificación todavía describe algoritmos que generan pronósticos puntuales. No especifica una distribución para \(Y_t\).
El paso conceptual decisivo consiste en preguntar:
¿qué proceso probabilístico podría generar observaciones de manera compatible con estas recursiones?
Los métodos de suavizamiento producen pronósticos puntuales. Un modelo estadístico debe especificar, además, cómo surge la incertidumbre alrededor de esos pronósticos. Esto permite construir verosimilitudes, distribuciones predictivas e intervalos (Hyndman y Athanasopoulos 2021, secs. 8.5–8.7).
La formulación ETS utiliza modelos de espacio-estado de innovaciones. Cada modelo contiene:
La notación es
\[ \operatorname{ETS}(E,T,S), \tag{8.27}\]
con
Por ejemplo, SES es un método de pronóstico. ETS\((A,N,N)\) y ETS\((M,N,N)\) son dos modelos probabilísticos que pueden generar el mismo tipo de pronóstico puntual de SES, pero implican estructuras de incertidumbre diferentes.
Considere
\[ Y_t = \ell_{t-1}+\varepsilon_t, \qquad \varepsilon_t\overset{\text{iid}}{\sim}N(0,\sigma^2), \tag{8.28}\]
\[ \ell_t = \ell_{t-1}+\alpha\varepsilon_t. \tag{8.29}\]
Antes de observar \(Y_t\),
\[ E(Y_t\mid\mathcal F_{t-1})=\ell_{t-1}, \]
de modo que \(\varepsilon_t\) es exactamente la innovación de un paso.
De Ecuación 8.28,
\[ \varepsilon_t=Y_t-\ell_{t-1}. \]
Sustituyendo en Ecuación 8.29,
\[ \begin{aligned} \ell_t &= \ell_{t-1} + \alpha(Y_t-\ell_{t-1})\\ &= \alpha Y_t+(1-\alpha)\ell_{t-1}. \end{aligned} \]
Hemos recuperado exactamente la recursión SES de Ecuación 8.3.
Esta derivación muestra que el suavizamiento puede verse como una actualización dirigida por una innovación aleatoria. La ecuación de observación especifica la distribución de esa innovación y convierte al algoritmo en un modelo probabilístico.
Condicionalmente a los parámetros y al estado \(\ell_t\),
\[ Y_{t+1}\mid\mathcal F_t \sim N(\ell_t,\sigma^2). \tag{8.30}\]
Por eso podemos hablar no sólo de
\[ \widehat y_{t+1\mid t}=\ell_t, \]
sino de una distribución para la observación futura.
Para horizontes mayores, la incertidumbre se propaga porque las innovaciones futuras también modifican los estados futuros. Este principio será central en modelos de espacio-estado.
Una alternativa es
\[ Y_t = \ell_{t-1}(1+\varepsilon_t), \tag{8.31}\]
\[ \ell_t = \ell_{t-1}(1+\alpha\varepsilon_t), \tag{8.32}\]
con \(E(\varepsilon_t)=0\).
En este caso
\[ \varepsilon_t = \frac{Y_t-\ell_{t-1}}{\ell_{t-1}}, \tag{8.33}\]
por lo que la innovación mide un error relativo. Si el nivel crece, una innovación del mismo tamaño relativo genera una desviación absoluta mayor.
ETS\((A,N,N)\) y ETS\((M,N,N)\) pueden producir pronósticos puntuales de tipo SES, pero no producen la misma distribución predictiva. En el modelo multiplicativo, la escala de la incertidumbre depende del nivel.
El mismo principio permite reescribir Holt. Defina el error de un paso como
\[ e_t = y_t-(\ell_{t-1}+b_{t-1}). \tag{8.34}\]
Entonces Ecuación 8.12 implica
\[ \ell_t = \ell_{t-1}+b_{t-1}+\alpha e_t. \tag{8.35}\]
Para la pendiente,
\[ \begin{aligned} b_t &= \beta^*(\ell_t-\ell_{t-1}) +(1-\beta^*)b_{t-1}\\ &= \beta^*(b_{t-1}+\alpha e_t) +(1-\beta^*)b_{t-1}\\ &= b_{t-1}+\alpha\beta^*e_t. \end{aligned} \]
Si definimos
\[ \beta=\alpha\beta^*, \tag{8.36}\]
obtenemos el modelo ETS\((A,A,N)\):
\[ Y_t = \ell_{t-1}+b_{t-1}+\varepsilon_t, \tag{8.37}\]
\[ \ell_t = \ell_{t-1}+b_{t-1}+\alpha\varepsilon_t, \tag{8.38}\]
\[ b_t = b_{t-1}+\beta\varepsilon_t. \tag{8.39}\]
En las ecuaciones clásicas del método de Holt hemos usado \(\beta^*\). En la representación de innovaciones es conveniente usar \(\beta=\alpha\beta^*\). Cuando un programa reporta parámetros ETS, debemos verificar qué parametrización utiliza antes de interpretar directamente el número.
Esta equivalencia es un ejemplo concreto de por qué el lenguaje de estados será útil: nivel y tendencia pueden actualizarse conjuntamente utilizando la misma innovación observada.
Una de las ventajas de la formulación probabilística es que podemos estimar los parámetros mediante máxima verosimilitud (Hyndman y Athanasopoulos 2021, sec. 8.6).
Suponga que, dados los estados iniciales y los parámetros, el modelo produce un pronóstico de un paso
\[ \mu_t = E(Y_t\mid\mathcal F_{t-1}). \]
Con error aditivo gaussiano,
\[ Y_t=\mu_t+\varepsilon_t, \qquad \varepsilon_t\sim N(0,\sigma^2). \]
Por tanto,
\[ f(y_t\mid\mathcal F_{t-1}) = \frac{1}{\sigma} \varphi\left( \frac{y_t-\mu_t}{\sigma} \right), \tag{8.40}\]
con \(\varphi(\cdot)\) la densidad normal estándar. La descomposición temporal de la verosimilitud es
\[ L(\theta;y_{1:T}) = \prod_{t=1}^{T} f(y_t\mid\mathcal F_{t-1};\theta), \tag{8.41}\]
para un vector \(\theta\) que contiene parámetros de suavizamiento, estados iniciales y parámetros de escala.
La log-verosimilitud aditiva contiene el término
\[ -\frac{1}{2\sigma^2} \sum_{t=1}^{T}e_t^2. \]
Por ello, bajo normalidad y para una especificación aditiva, maximizar la verosimilitud está estrechamente relacionado con minimizar SSE (Hyndman y Athanasopoulos 2021, sec. 8.6).
Con error multiplicativo,
\[ Y_t=\mu_t(1+\varepsilon_t), \]
y entonces
\[ \varepsilon_t=\frac{y_t}{\mu_t}-1. \]
Si \(\varepsilon_t\sim N(0,\sigma^2)\) y \(\mu_t\neq0\), el cambio de variable introduce
\[ f(y_t\mid\mathcal F_{t-1}) = \frac{1}{\sigma|\mu_t|} \varphi\left( \frac{y_t/\mu_t-1}{\sigma} \right). \tag{8.42}\]
Aparece el factor \(1/|\mu_t|\). Por ello, la estimación por verosimilitud de un modelo con error multiplicativo no equivale simplemente a minimizar una suma de cuadrados de errores absolutos.
Esta diferencia es una de las razones por las que conviene distinguir el algoritmo de suavizamiento del modelo estadístico que lo sustenta.
Un modelo ETS puede involucrar:
La estimación no consiste, por tanto, en ajustar únicamente dos o tres coeficientes. Los estados iniciales también deben determinarse.
Dentro del marco ETS, modelos diferentes pueden compararse mediante criterios de información como AIC, AICc y BIC (Hyndman y Athanasopoulos 2021, sec. 8.6). La lógica es la misma estudiada en Sección 4.9: recompensar ajuste mediante la verosimilitud, pero penalizar complejidad.
Para AIC,
\[ \operatorname{AIC} = -2\log L+2k, \tag{8.43}\]
mientras que
\[ \operatorname{AICc} = \operatorname{AIC} + \frac{2k(k+1)}{T-k-1}. \tag{8.44}\]
Aquí \(k\) debe contabilizar los parámetros estimados pertinentes, incluidos los estados iniciales que forman parte del ajuste ETS.
AICc compara modelos mediante ajuste penalizado dentro de la muestra utilizada para estimarlos. La semana 7 mostró por qué esto no sustituye una comparación de errores genuinamente fuera de muestra. En este capítulo utilizaremos ambos tipos de evidencia.
fableLa función ETS() permite fijar componentes concretos o dejar que el procedimiento considere combinaciones compatibles con la serie.
Por ejemplo,
# SES con error aditivo
ETS(y ~ error("A") + trend("N") + season("N"))
# Holt con error aditivo
ETS(y ~ error("A") + trend("A") + season("N"))
# Holt amortiguado
ETS(y ~ error("A") + trend("Ad") + season("N"))
# Selección automática
ETS(y)La selección automática es útil, pero debemos interpretar qué estructura fue elegida, verificar las innovaciones y evaluar el pronóstico fuera de muestra.
No toda combinación algebraicamente imaginable es igualmente apropiada. Los modelos con componentes multiplicativos requieren cuidado con ceros, valores negativos y estados cercanos a cero. Además, existen regiones de parámetros que producen dinámicas no admisibles.
Por estas razones, una implementación automática moderna impone restricciones de parámetros y descarta ciertas combinaciones problemáticas. En una aplicación conviene evitar dos errores:
La estabilidad es parte de la definición práctica del modelo.
Utilizaremos el conjunto tourism de tsibbledata, que contiene viajes domésticos trimestrales en Australia (O’Hara-Wild et al. 2022). Siguiendo el ejemplo de FPP3, agregaremos los viajes cuyo propósito es vacaciones (Hyndman y Athanasopoulos 2021, secs. 8.3 y 8.6).
Nuestro objetivo será responder cuatro preguntas:
aus_holidays <- tourism |>
filter(Purpose == "Holiday") |>
summarise(Trips = sum(Trips) / 1e3)
stopifnot(
nrow(aus_holidays) > 40,
!anyNA(aus_holidays$Trips),
all(is.finite(aus_holidays$Trips)),
all(aus_holidays$Trips > 0)
)La división entre \(1000\) expresa los viajes en millones aproximadamente y mejora la escala gráfica sin cambiar la estructura temporal.
aus_holidays |>
autoplot(Trips) +
labs(
x = "Trimestre",
y = "Viajes (millones)"
)
Antes de ajustar un ETS, conviene recuperar las preguntas exploratorias de la semana 1:
No necesitamos contestarlas de forma perfecta mirando un solo gráfico. Servirán para definir candidatos que luego compararemos.
Reservaremos los últimos 12 trimestres como período de prueba. Toda selección de modelos del análisis inicial se realizará con el entrenamiento.
h_test_ets <- 12
n_tourism <- nrow(aus_holidays)
tourism_train <- aus_holidays |>
slice_head(n = n_tourism - h_test_ets)
tourism_test <- aus_holidays |>
slice_tail(n = h_test_ets)
stopifnot(
nrow(tourism_train) + nrow(tourism_test) == nrow(aus_holidays),
max(tourism_train$Quarter) < min(tourism_test$Quarter)
)particion_tourism <- bind_rows(
tourism_train |>
summarise(
conjunto = "Entrenamiento",
inicio = min(Quarter),
fin = max(Quarter),
n = n()
),
tourism_test |>
summarise(
conjunto = "Prueba",
inicio = min(Quarter),
fin = max(Quarter),
n = n()
)
)
knitr::kable(
particion_tourism,
caption = "Partición temporal utilizada en la aplicación de turismo."
)| Quarter | conjunto | inicio | fin | n |
|---|---|---|---|---|
| 1998 Q1 | Entrenamiento | 1998 Q1 | 1998 Q1 | 1 |
| 1998 Q2 | Entrenamiento | 1998 Q2 | 1998 Q2 | 1 |
| 1998 Q3 | Entrenamiento | 1998 Q3 | 1998 Q3 | 1 |
| 1998 Q4 | Entrenamiento | 1998 Q4 | 1998 Q4 | 1 |
| 1999 Q1 | Entrenamiento | 1999 Q1 | 1999 Q1 | 1 |
| 1999 Q2 | Entrenamiento | 1999 Q2 | 1999 Q2 | 1 |
| 1999 Q3 | Entrenamiento | 1999 Q3 | 1999 Q3 | 1 |
| 1999 Q4 | Entrenamiento | 1999 Q4 | 1999 Q4 | 1 |
| 2000 Q1 | Entrenamiento | 2000 Q1 | 2000 Q1 | 1 |
| 2000 Q2 | Entrenamiento | 2000 Q2 | 2000 Q2 | 1 |
| 2000 Q3 | Entrenamiento | 2000 Q3 | 2000 Q3 | 1 |
| 2000 Q4 | Entrenamiento | 2000 Q4 | 2000 Q4 | 1 |
| 2001 Q1 | Entrenamiento | 2001 Q1 | 2001 Q1 | 1 |
| 2001 Q2 | Entrenamiento | 2001 Q2 | 2001 Q2 | 1 |
| 2001 Q3 | Entrenamiento | 2001 Q3 | 2001 Q3 | 1 |
| 2001 Q4 | Entrenamiento | 2001 Q4 | 2001 Q4 | 1 |
| 2002 Q1 | Entrenamiento | 2002 Q1 | 2002 Q1 | 1 |
| 2002 Q2 | Entrenamiento | 2002 Q2 | 2002 Q2 | 1 |
| 2002 Q3 | Entrenamiento | 2002 Q3 | 2002 Q3 | 1 |
| 2002 Q4 | Entrenamiento | 2002 Q4 | 2002 Q4 | 1 |
| 2003 Q1 | Entrenamiento | 2003 Q1 | 2003 Q1 | 1 |
| 2003 Q2 | Entrenamiento | 2003 Q2 | 2003 Q2 | 1 |
| 2003 Q3 | Entrenamiento | 2003 Q3 | 2003 Q3 | 1 |
| 2003 Q4 | Entrenamiento | 2003 Q4 | 2003 Q4 | 1 |
| 2004 Q1 | Entrenamiento | 2004 Q1 | 2004 Q1 | 1 |
| 2004 Q2 | Entrenamiento | 2004 Q2 | 2004 Q2 | 1 |
| 2004 Q3 | Entrenamiento | 2004 Q3 | 2004 Q3 | 1 |
| 2004 Q4 | Entrenamiento | 2004 Q4 | 2004 Q4 | 1 |
| 2005 Q1 | Entrenamiento | 2005 Q1 | 2005 Q1 | 1 |
| 2005 Q2 | Entrenamiento | 2005 Q2 | 2005 Q2 | 1 |
| 2005 Q3 | Entrenamiento | 2005 Q3 | 2005 Q3 | 1 |
| 2005 Q4 | Entrenamiento | 2005 Q4 | 2005 Q4 | 1 |
| 2006 Q1 | Entrenamiento | 2006 Q1 | 2006 Q1 | 1 |
| 2006 Q2 | Entrenamiento | 2006 Q2 | 2006 Q2 | 1 |
| 2006 Q3 | Entrenamiento | 2006 Q3 | 2006 Q3 | 1 |
| 2006 Q4 | Entrenamiento | 2006 Q4 | 2006 Q4 | 1 |
| 2007 Q1 | Entrenamiento | 2007 Q1 | 2007 Q1 | 1 |
| 2007 Q2 | Entrenamiento | 2007 Q2 | 2007 Q2 | 1 |
| 2007 Q3 | Entrenamiento | 2007 Q3 | 2007 Q3 | 1 |
| 2007 Q4 | Entrenamiento | 2007 Q4 | 2007 Q4 | 1 |
| 2008 Q1 | Entrenamiento | 2008 Q1 | 2008 Q1 | 1 |
| 2008 Q2 | Entrenamiento | 2008 Q2 | 2008 Q2 | 1 |
| 2008 Q3 | Entrenamiento | 2008 Q3 | 2008 Q3 | 1 |
| 2008 Q4 | Entrenamiento | 2008 Q4 | 2008 Q4 | 1 |
| 2009 Q1 | Entrenamiento | 2009 Q1 | 2009 Q1 | 1 |
| 2009 Q2 | Entrenamiento | 2009 Q2 | 2009 Q2 | 1 |
| 2009 Q3 | Entrenamiento | 2009 Q3 | 2009 Q3 | 1 |
| 2009 Q4 | Entrenamiento | 2009 Q4 | 2009 Q4 | 1 |
| 2010 Q1 | Entrenamiento | 2010 Q1 | 2010 Q1 | 1 |
| 2010 Q2 | Entrenamiento | 2010 Q2 | 2010 Q2 | 1 |
| 2010 Q3 | Entrenamiento | 2010 Q3 | 2010 Q3 | 1 |
| 2010 Q4 | Entrenamiento | 2010 Q4 | 2010 Q4 | 1 |
| 2011 Q1 | Entrenamiento | 2011 Q1 | 2011 Q1 | 1 |
| 2011 Q2 | Entrenamiento | 2011 Q2 | 2011 Q2 | 1 |
| 2011 Q3 | Entrenamiento | 2011 Q3 | 2011 Q3 | 1 |
| 2011 Q4 | Entrenamiento | 2011 Q4 | 2011 Q4 | 1 |
| 2012 Q1 | Entrenamiento | 2012 Q1 | 2012 Q1 | 1 |
| 2012 Q2 | Entrenamiento | 2012 Q2 | 2012 Q2 | 1 |
| 2012 Q3 | Entrenamiento | 2012 Q3 | 2012 Q3 | 1 |
| 2012 Q4 | Entrenamiento | 2012 Q4 | 2012 Q4 | 1 |
| 2013 Q1 | Entrenamiento | 2013 Q1 | 2013 Q1 | 1 |
| 2013 Q2 | Entrenamiento | 2013 Q2 | 2013 Q2 | 1 |
| 2013 Q3 | Entrenamiento | 2013 Q3 | 2013 Q3 | 1 |
| 2013 Q4 | Entrenamiento | 2013 Q4 | 2013 Q4 | 1 |
| 2014 Q1 | Entrenamiento | 2014 Q1 | 2014 Q1 | 1 |
| 2014 Q2 | Entrenamiento | 2014 Q2 | 2014 Q2 | 1 |
| 2014 Q3 | Entrenamiento | 2014 Q3 | 2014 Q3 | 1 |
| 2014 Q4 | Entrenamiento | 2014 Q4 | 2014 Q4 | 1 |
| 2015 Q1 | Prueba | 2015 Q1 | 2015 Q1 | 1 |
| 2015 Q2 | Prueba | 2015 Q2 | 2015 Q2 | 1 |
| 2015 Q3 | Prueba | 2015 Q3 | 2015 Q3 | 1 |
| 2015 Q4 | Prueba | 2015 Q4 | 2015 Q4 | 1 |
| 2016 Q1 | Prueba | 2016 Q1 | 2016 Q1 | 1 |
| 2016 Q2 | Prueba | 2016 Q2 | 2016 Q2 | 1 |
| 2016 Q3 | Prueba | 2016 Q3 | 2016 Q3 | 1 |
| 2016 Q4 | Prueba | 2016 Q4 | 2016 Q4 | 1 |
| 2017 Q1 | Prueba | 2017 Q1 | 2017 Q1 | 1 |
| 2017 Q2 | Prueba | 2017 Q2 | 2017 Q2 | 1 |
| 2017 Q3 | Prueba | 2017 Q3 | 2017 Q3 | 1 |
| 2017 Q4 | Prueba | 2017 Q4 | 2017 Q4 | 1 |
Ajustaremos cuatro modelos:
ETS().fit_tourism_ets <- tourism_train |>
model(
`ETS(A,A,A)` = ETS(
Trips ~ error("A") + trend("A") + season("A")
),
`ETS(M,A,M)` = ETS(
Trips ~ error("M") + trend("A") + season("M")
),
`ETS(M,Ad,M)` = ETS(
Trips ~ error("M") + trend("Ad") + season("M")
),
`ETS automatico` = ETS(Trips)
)
# Especificación seleccionada automáticamente
fit_tourism_ets |>
select(`ETS automatico`) |>
report()Series: Trips
Model: ETS(M,N,A)
Smoothing parameters:
alpha = 0.252709
gamma = 0.0001000325
Initial states:
l[0] s[0] s[-1] s[-2] s[-3]
9.662824 -0.5621126 -0.6894698 -0.2545154 1.506098
sigma^2: 0.0022
AIC AICc BIC
181.0660 182.9327 196.6026 En lugar de asumir que glance() siempre contiene exactamente las mismas columnas en todas las versiones del paquete, seleccionamos únicamente las columnas disponibles que necesitamos.
ic_tourism_ets <- glance(fit_tourism_ets) |>
select(any_of(c(".model", "log_lik", "AIC", "AICc", "BIC"))) |>
arrange(AICc)
knitr::kable(
ic_tourism_ets,
digits = 3,
caption = "Criterios de información de los candidatos ETS ajustados al período de entrenamiento."
)| .model | log_lik | AIC | AICc | BIC |
|---|---|---|---|---|
| ETS automatico | -83.533 | 181.066 | 182.933 | 196.603 |
| ETS(M,A,M) | -84.120 | 186.240 | 189.343 | 206.216 |
| ETS(M,Ad,M) | -83.360 | 186.720 | 190.580 | 208.915 |
| ETS(A,A,A) | -85.255 | 188.510 | 191.613 | 208.485 |
El menor AICc identifica el candidato preferido según ese criterio y esa muestra de entrenamiento. Todavía no hemos utilizado el período de prueba.
param_tourism_ets <- tidy(fit_tourism_ets) |>
filter(term %in% c("alpha", "beta", "gamma", "phi")) |>
select(.model, term, estimate) |>
arrange(.model, term)
knitr::kable(
param_tourism_ets,
digits = 4,
caption = "Parámetros de suavizamiento estimados para los candidatos ETS."
)| .model | term | estimate |
|---|---|---|
| ETS automatico | alpha | 0.2527 |
| ETS automatico | gamma | 0.0001 |
| ETS(A,A,A) | alpha | 0.2562 |
| ETS(A,A,A) | beta | 0.0254 |
| ETS(A,A,A) | gamma | 0.0001 |
| ETS(M,A,M) | alpha | 0.2466 |
| ETS(M,A,M) | beta | 0.0156 |
| ETS(M,A,M) | gamma | 0.0001 |
| ETS(M,Ad,M) | alpha | 0.2614 |
| ETS(M,Ad,M) | beta | 0.0001 |
| ETS(M,Ad,M) | gamma | 0.0001 |
| ETS(M,Ad,M) | phi | 0.9783 |
La interpretación debe ser cualitativa y coherente con la parametrización ETS:
No debemos interpretar un parámetro aislado sin observar la escala del error, los demás estados y el comportamiento predictivo del modelo.
Para visualizar los tres tipos de estado, utilizaremos el candidato ETS\((M,A_d,M)\), cuya especificación contiene nivel, pendiente y estacionalidad.
comp_tourism <- fit_tourism_ets |>
select(`ETS(M,Ad,M)`) |>
components() |>
as_tibble() |>
select(Quarter, any_of(c("level", "slope", "season"))) |>
pivot_longer(
cols = -Quarter,
names_to = "estado",
values_to = "valor"
)comp_tourism |>
ggplot(aes(x = Quarter, y = valor)) +
geom_line() +
facet_wrap(~ estado, ncol = 1, scales = "free_y") +
labs(
x = "Trimestre",
y = NULL
)
El nivel, la pendiente y la estacionalidad son cantidades internas del modelo. Su significado depende de la especificación aditiva o multiplicativa. En particular, un estado estacional multiplicativo se interpreta como un factor relativo.
Como en Sección 4.11, un modelo de pronóstico debe dejar innovaciones sin estructura temporal evidente. Utilizaremos el ETS automático del período de entrenamiento.
fit_tourism_ets |>
select(`ETS automatico`) |>
gg_tsresiduals()
El diagnóstico tiene dos objetivos distintos:
Una ACF residual con dependencia marcada sugeriría que la estructura de nivel, tendencia y estacionalidad no agota la dinámica temporal de la serie.
Para el holdout compararemos:
La referencia estacional es especialmente importante porque la serie presenta una periodicidad trimestral fuerte. Para el competidor ARIMA no repetiremos todo el proceso de identificación de la semana 6. A partir del período de entrenamiento fijamos \(d=0\) y \(D=1\): preservamos el nivel no estacional y eliminamos la estacionalidad trimestral mediante una diferencia de período cuatro. Condicionalmente a esa decisión, ARIMA() selecciona \(p,q,P,Q\) dentro de rangos moderados mediante AICc. Esta separación entre elección de la diferenciación y selección de órdenes AR/MA es consistente con Sección 6.15.1.
Además, usamos .safely = FALSE. Por defecto, model() reemplaza un ajuste que falla por un modelo nulo para permitir que los demás modelos continúen. En una nota reproducible preferimos que un error del ARIMA detenga el render y muestre su causa en lugar de producir pronósticos faltantes que luego aparezcan como NaN en las métricas.
fit_tourism_holdout <- tourism_train |>
model(
`Naive estacional` = SNAIVE(Trips),
ETS = ETS(Trips),
ARIMA = ARIMA(
Trips ~
pdq(p = 0:3, d = 0, q = 0:3) +
PDQ(P = 0:2, D = 1, Q = 0:2)
),
.safely = FALSE
)
# Hacemos visible la especificación ARIMA seleccionada en entrenamiento.
report(fit_tourism_holdout |> select(ARIMA))Series: Trips
Model: ARIMA(0,0,0)(0,1,2)[4]
Coefficients:
sma1 sma2
-0.6082 -0.2880
s.e. 0.2574 0.1871
sigma^2 estimated as 0.2277: log likelihood=-45.29
AIC=96.59 AICc=96.99 BIC=103.06fc_tourism_holdout <- fit_tourism_holdout |>
forecast(h = h_test_ets)
# Una falla silenciosa del ajuste no debe propagarse a tablas o gráficos.
fc_tourism_holdout_check <- fc_tourism_holdout |>
as_tibble()
stopifnot(
nrow(fc_tourism_holdout_check) == 3 * h_test_ets,
!anyNA(fc_tourism_holdout_check$.mean),
all(is.finite(fc_tourism_holdout_check$.mean))
)La siguiente figura distingue explícitamente el período de entrenamiento y el período de prueba mediante el tipo de línea. Además, las observaciones de prueba se marcan con puntos para evitar que se confundan con las trayectorias pronosticadas.
tourism_periodos <- bind_rows(
tourism_train |> mutate(periodo = "Entrenamiento"),
tourism_test |> mutate(periodo = "Prueba")
)
fc_tourism_holdout_tbl <- fc_tourism_holdout |>
as_tibble() |>
select(.model, Quarter, .mean)
ggplot() +
geom_line(
data = tourism_periodos,
aes(x = Quarter, y = Trips, linetype = periodo),
linewidth = 0.65
) +
geom_point(
data = tourism_test,
aes(x = Quarter, y = Trips),
shape = 21,
fill = "white",
size = 2
) +
geom_line(
data = fc_tourism_holdout_tbl,
aes(x = Quarter, y = .mean, colour = .model),
linewidth = 0.85
) +
scale_linetype_manual(
values = c("Entrenamiento" = "solid", "Prueba" = "longdash")
) +
labs(
x = "Trimestre",
y = "Viajes (millones)",
colour = "Método",
linetype = "Serie observada"
)
acc_tourism_holdout <- fc_tourism_holdout |>
accuracy(aus_holidays) |>
select(.model, RMSE, MAE, MASE) |>
arrange(RMSE)
knitr::kable(
acc_tourism_holdout,
digits = 4,
caption = "Desempeño en los 12 trimestres reservados para prueba."
)| .model | RMSE | MAE | MASE |
|---|---|---|---|
| Naive estacional | 0.9476 | 0.8106 | 2.0347 |
| ETS | 1.0342 | 0.9099 | 2.2840 |
| ARIMA | 1.3796 | 1.2605 | 3.1641 |
Una sola partición es fácil de interpretar, pero puede depender excesivamente de qué 12 observaciones quedaron al final. Por eso repetiremos la comparación en varios orígenes.
Utilizaremos ventanas expansivas. Para que cada origen tenga ocho trimestres futuros observables, eliminamos temporalmente los últimos ocho puntos al construir los orígenes y después pronosticamos \(h=8\) en cada uno.
h_cv_ets <- 8
init_cv_ets <- 48
step_cv_ets <- 2
tourism_cv_base <- aus_holidays |>
slice_head(n = nrow(aus_holidays) - h_cv_ets)
tourism_cv <- tourism_cv_base |>
stretch_tsibble(
.init = init_cv_ets,
.step = step_cv_ets
)
stopifnot(
n_distinct(tourism_cv$.id) >= 5,
all(table(tourism_cv$.id) >= init_cv_ets)
)Se hace el ajuste automático o naive en cada conjunto móvil de entrenamiento:
fit_tourism_cv <- tourism_cv |>
model(
`Naive estacional` = SNAIVE(Trips),
ETS = ETS(Trips),
ARIMA = ARIMA(
Trips ~
pdq(p = 0:3, d = 0, q = 0:3) +
PDQ(P = 0:2, D = 1, Q = 0:2)
),
.safely = FALSE
)
fc_tourism_cv <- fit_tourism_cv |>
forecast(h = h_cv_ets) |>
group_by(.id, .model) |>
mutate(h = row_number()) |>
ungroup() |>
as_fable(
response = "Trips",
distribution = Trips
)
# Verificación explícita antes de calcular métricas.
stopifnot(
!anyNA(as_tibble(fc_tourism_cv)$.mean),
all(is.finite(as_tibble(fc_tourism_cv)$.mean))
)Como en la semana 7, construiremos la escala de MASE usando únicamente la información disponible en cada ventana. Para datos trimestrales, la referencia de escala utiliza rezago cuatro.
escalas_tourism_cv <- tourism_cv |>
as_tibble() |>
group_by(.id) |>
arrange(Quarter, .by_group = TRUE) |>
mutate(dif_snaive = Trips - lag(Trips, 4)) |>
summarise(
escala_MASE = mean(abs(dif_snaive), na.rm = TRUE),
.groups = "drop"
)
errores_tourism_cv <- fc_tourism_cv |>
as_tibble() |>
select(.id, .model, Quarter, h, .mean) |>
left_join(
aus_holidays |>
as_tibble() |>
transmute(Quarter, observado = Trips),
by = "Quarter"
) |>
left_join(escalas_tourism_cv, by = ".id") |>
mutate(
error = observado - .mean,
abs_escalado = abs(error) / escala_MASE
)acc_tourism_cv_h <- errores_tourism_cv |>
group_by(.model, h) |>
summarise(
RMSE = sqrt(mean(error^2)),
MAE = mean(abs(error)),
MASE = mean(abs_escalado),
.groups = "drop"
)
acc_tourism_cv_global <- errores_tourism_cv |>
group_by(.model) |>
summarise(
RMSE = sqrt(mean(error^2)),
MAE = mean(abs(error)),
MASE = mean(abs_escalado),
n_pronosticos = n(),
.groups = "drop"
) |>
arrange(RMSE)knitr::kable(
acc_tourism_cv_global,
digits = 4,
caption = "Desempeño agregado bajo validación temporal expansiva para turismo."
)| .model | RMSE | MAE | MASE | n_pronosticos |
|---|---|---|---|---|
| Naive estacional | 0.6325 | 0.5197 | 1.3296 | 104 |
| ETS | 0.6951 | 0.5726 | 1.4655 | 104 |
| ARIMA | 0.7500 | 0.6052 | 1.5439 | 104 |
acc_tourism_cv_h |>
ggplot(aes(x = h, y = RMSE, linetype = .model)) +
geom_line(linewidth = 0.8) +
geom_point() +
scale_x_continuous(breaks = 1:h_cv_ets) +
labs(
x = "Horizonte h",
y = "RMSE",
linetype = "Método"
)
La figura es más informativa que una sola métrica agregada si la calidad relativa cambia con \(h\). Un método puede responder bien a corto plazo y deteriorarse más rápido que otro.
Evaluaremos intervalos centrales de 80%. Como ya discutimos en Sección 7.10, cobertura y amplitud deben leerse conjuntamente.
intervalos_tourism_cv <- fc_tourism_cv |>
hilo(level = 80) |>
select(.id, .model, Quarter, h, `80%`) |>
fabletools::unpack_hilo(`80%`, names_sep = "_") |>
transmute(
.id,
.model,
Quarter,
h,
inferior = `80%_lower`,
superior = `80%_upper`
) |>
as_tibble() |>
left_join(
aus_holidays |>
as_tibble() |>
transmute(Quarter, observado = Trips),
by = "Quarter"
) |>
mutate(
cubierto = observado >= inferior & observado <= superior,
amplitud = superior - inferior
)resumen_intervalos_tourism <- intervalos_tourism_cv |>
group_by(.model) |>
summarise(
cobertura_80 = mean(cubierto),
amplitud_media = mean(amplitud),
n = n(),
.groups = "drop"
)
knitr::kable(
resumen_intervalos_tourism,
digits = 4,
caption = "Cobertura y amplitud de intervalos predictivos nominales de 80% bajo validación temporal."
)| .model | cobertura_80 | amplitud_media | n |
|---|---|---|---|
| ARIMA | 0.5385 | 1.2949 | 104 |
| ETS | 0.5577 | 1.1324 | 104 |
| Naive estacional | 0.8365 | 1.6502 | 104 |
cobertura_tourism_h <- intervalos_tourism_cv |>
group_by(.model, h) |>
summarise(
cobertura = mean(cubierto),
.groups = "drop"
)
ggplot(cobertura_tourism_h, aes(x = h, y = cobertura, linetype = .model)) +
geom_hline(yintercept = 0.80, linetype = 2) +
geom_line(linewidth = 0.8) +
geom_point() +
scale_x_continuous(breaks = 1:h_cv_ets) +
labs(
x = "Horizonte h",
y = "Cobertura empírica",
linetype = "Método"
)
Un modelo con menor RMSE puede tener intervalos mal calibrados. Otro puede alcanzar alta cobertura porque sus intervalos son demasiado anchos. La conclusión predictiva debe considerar el objetivo de la aplicación y las métricas relevantes.
En Sección 8.19.2 evaluamos los intervalos predictivos mediante dos cantidades complementarias: cobertura empírica y amplitud promedio. La cobertura estudia calibración, mientras que la amplitud estudia agudeza. Sin embargo, analizar ambas cantidades por separado puede dificultar la comparación entre procedimientos.
En la semana 7 introdujimos dos reglas de puntuación que permiten resumir distintos aspectos de la distribución predictiva:
En ambos casos, valores menores indican mejor desempeño predictivo. Por tanto, estas medidas complementan, pero no reemplazan, RMSE, MAE, MASE, cobertura y amplitud (Hyndman y Athanasopoulos 2021, sec. 5.9).
Para un intervalo predictivo central de nivel \(100(1-\alpha)\%\),
\[ I_{t+h\mid t} = [\ell_{\alpha,t+h\mid t},u_{\alpha,t+h\mid t}], \]
la puntuación de Winkler es
\[ W_{\alpha,t+h\mid t} = u_{\alpha,t+h\mid t} - \ell_{\alpha,t+h\mid t} + \begin{cases} \dfrac{2}{\alpha} \left( \ell_{\alpha,t+h\mid t}-y_{t+h} \right), & y_{t+h}<\ell_{\alpha,t+h\mid t}, \\[1.2ex] 0, & \ell_{\alpha,t+h\mid t} \leq y_{t+h} \leq u_{\alpha,t+h\mid t}, \\[1.2ex] \dfrac{2}{\alpha} \left( y_{t+h}-u_{\alpha,t+h\mid t} \right), & y_{t+h}>u_{\alpha,t+h\mid t}. \end{cases} \]
Si la observación futura cae dentro del intervalo, el score coincide con su amplitud. Si queda fuera, se añade una penalización proporcional a la distancia entre la observación y el extremo del intervalo. De esta forma, Winkler evita favorecer automáticamente intervalos muy anchos: aumentar la amplitud puede mejorar cobertura, pero también aumenta directamente el score.
Continuaremos utilizando los intervalos centrales de 80% estudiados en Sección 8.19.2, por lo que en este caso
\[ \alpha=0.20. \]
Calculamos primero la puntuación promedio sobre todos los orígenes y horizontes.
winkler_tourism <- fc_tourism_cv |>
accuracy(
aus_holidays,
measures = list(Winkler = winkler_score),
level = 80,
by = ".model"
) |>
select(.model, Winkler) |>
arrange(Winkler)
winkler_tourism_h <- fc_tourism_cv |>
accuracy(
aus_holidays,
measures = list(Winkler = winkler_score),
level = 80,
by = c("h", ".model")
)knitr::kable(
winkler_tourism,
digits = 4,
caption = "Puntuación de Winkler promedio para intervalos predictivos de 80% bajo validación temporal. Menor es mejor."
)| .model | Winkler |
|---|---|
| Naive estacional | 2.0828 |
| ETS | 2.7336 |
| ARIMA | 2.9359 |
La puntuación agregada resume todos los errores cometidos por los intervalos, pero puede ocultar diferencias importantes entre horizontes. Por ello, examinamos también
\[ \overline W_\alpha(h), \]
es decir, el score promedio para cada horizonte \(h\).
ggplot(
winkler_tourism_h,
aes(x = h, y = Winkler, linetype = .model)
) +
geom_line(linewidth = 0.8) +
geom_point() +
scale_x_continuous(breaks = 1:h_cv_ets) +
labs(
x = "Horizonte h",
y = "Puntuación de Winkler",
linetype = "Método"
)
La interpretación debe hacerse conjuntamente con la cobertura de Figura 8.12. Por ejemplo, un método podría producir intervalos relativamente estrechos y obtener un buen Winkler cuando las observaciones quedan dentro de ellos, pero recibir penalizaciones grandes en horizontes donde la cobertura se deteriora.
En consecuencia, Winkler resume el compromiso entre cobertura y amplitud, pero sigue evaluando solamente un intervalo específico de la distribución predictiva.
Para evaluar la distribución predictiva completa retomamos el Continuous Ranked Probability Score. Si \(F_{t+h\mid t}\) es la función de distribución predictiva y \(y_{t+h}\) es la realización observada,
\[ \operatorname{CRPS} \left( F_{t+h\mid t},y_{t+h} \right) = \int_{-\infty}^{\infty} \left[ F_{t+h\mid t}(z) - \mathbf 1\{y_{t+h}\leq z\} \right]^2 \,dz. \]
El CRPS compara toda la función de distribución predictiva con una distribución degenerada en la observación realizada. Tiene las mismas unidades que la variable pronosticada y, nuevamente, menor es mejor.
Una distribución excesivamente concentrada puede recibir una penalización grande cuando la observación cae lejos de su centro; una distribución demasiado dispersa, por otra parte, pierde agudeza. Por tanto, CRPS evalúa simultáneamente localización e incertidumbre predictiva sin tener que escoger previamente un nivel particular como 80% o 95%.
El cálculo debe realizarse para todas las distribuciones predictivas generadas en cada origen y horizonte. Este es uno de los bloques computacionalmente más costosos del capítulo, por lo que lo almacenamos en caché.
crps_tourism <- fc_tourism_cv |>
accuracy(
aus_holidays,
measures = list(CRPS = CRPS),
by = ".model"
) |>
select(.model, CRPS) |>
arrange(CRPS)
crps_tourism_h <- fc_tourism_cv |>
accuracy(
aus_holidays,
measures = list(CRPS = CRPS),
by = c("h", ".model")
)El objeto almacenado en caché contiene tanto la evaluación agregada como la evaluación por horizonte. De esta manera, cambios posteriores en texto, tablas o figuras no obligan a recalcular el CRPS mientras no cambien los pronósticos producidos en ajustar-tourism-cv.
knitr::kable(
crps_tourism,
digits = 4,
caption = "CRPS promedio bajo validación temporal para ETS, ARIMA y naïve estacional. Menor es mejor."
)| .model | CRPS |
|---|---|
| Naive estacional | 0.3569 |
| ETS | 0.4154 |
| ARIMA | 0.4426 |
Al igual que con RMSE y Winkler, el valor agregado puede ocultar heterogeneidad entre horizontes.
ggplot(
crps_tourism_h,
aes(x = h, y = CRPS, linetype = .model)
) +
geom_line(linewidth = 0.8) +
geom_point() +
scale_x_continuous(breaks = 1:h_cv_ets) +
labs(
x = "Horizonte h",
y = "CRPS",
linetype = "Método"
)
La figura permite investigar si la distribución predictiva de alguno de los métodos se deteriora más rápidamente conforme aumenta \(h\). Una ordenación distinta a la obtenida mediante RMSE no constituye una contradicción: ambas métricas responden a preguntas diferentes. RMSE evalúa el pronóstico puntual; CRPS evalúa la distribución predictiva completa.
Podemos finalmente reunir en una misma tabla algunas de las medidas utilizadas en esta aplicación.
resumen_tourism_integral <- acc_tourism_cv_global |>
select(.model, RMSE, MAE, MASE) |>
left_join(
resumen_intervalos_tourism |>
select(
.model,
cobertura_80,
amplitud_media
),
by = ".model"
) |>
left_join(
winkler_tourism,
by = ".model"
) |>
left_join(
crps_tourism,
by = ".model"
) |>
arrange(CRPS)
knitr::kable(
resumen_tourism_integral,
digits = 4,
caption = "Síntesis del desempeño fuera de muestra para turismo bajo los mismos orígenes y horizontes."
)| .model | RMSE | MAE | MASE | cobertura_80 | amplitud_media | Winkler | CRPS |
|---|---|---|---|---|---|---|---|
| Naive estacional | 0.6325 | 0.5197 | 1.3296 | 0.8365 | 1.6502 | 2.0828 | 0.3569 |
| ETS | 0.6951 | 0.5726 | 1.4655 | 0.5577 | 1.1324 | 2.7336 | 0.4154 |
| ARIMA | 0.7500 | 0.6052 | 1.5439 | 0.5385 | 1.2949 | 2.9359 | 0.4426 |
La tabla debe interpretarse por columnas y no como un concurso en el que un único número determina automáticamente el procedimiento elegido:
En particular, no es necesario que el procedimiento con menor RMSE tenga también el menor Winkler o CRPS. Si ocurre un desacuerdo, la pregunta relevante es qué aspecto del pronóstico explica esa diferencia.
Todas las columnas de Tabla 8.11 corresponden a los mismos orígenes, horizontes y observaciones futuras. Esta condición es fundamental: una diferencia de score sólo puede atribuirse razonablemente al procedimiento de pronóstico cuando el experimento predictivo es comparable.
Además, no deben compararse directamente las magnitudes numéricas de RMSE, Winkler y CRPS entre sí. Aunque todas están expresadas en la escala de la variable, corresponden a funciones de pérdida diferentes. Su utilidad está en comparar los mismos procedimientos dentro de cada métrica.
ETS y ARIMA no representan dos versiones del mismo mecanismo. Ambos pueden producir buenos pronósticos, pero organizan la dependencia temporal de maneras diferentes.
comparacion_ets_arima <- tribble(
~Aspecto, ~ETS, ~ARIMA,
"Representación central", "Estados de nivel, tendencia y estacionalidad", "Rezagos, diferencias e innovaciones",
"No estacionariedad", "Puede representarse mediante estados evolutivos", "Se maneja típicamente mediante diferenciación",
"Estacionalidad", "Estado estacional aditivo o multiplicativo", "Operadores AR y MA estacionales y diferenciación estacional",
"Memoria", "Implícita en recursiones de actualización", "Explícita mediante estructura ARMA",
"Selección automática", "Combinaciones E-T-S", "Órdenes p-d-q y P-D-Q",
"Incertidumbre", "Distribución predictiva del modelo de innovaciones", "Distribución predictiva del ARIMA ajustado",
"Diagnóstico", "Innovaciones", "Innovaciones"
)
knitr::kable(
comparacion_ets_arima,
caption = "Comparación conceptual entre las representaciones ETS y ARIMA."
)| Aspecto | ETS | ARIMA |
|---|---|---|
| Representación central | Estados de nivel, tendencia y estacionalidad | Rezagos, diferencias e innovaciones |
| No estacionariedad | Puede representarse mediante estados evolutivos | Se maneja típicamente mediante diferenciación |
| Estacionalidad | Estado estacional aditivo o multiplicativo | Operadores AR y MA estacionales y diferenciación estacional |
| Memoria | Implícita en recursiones de actualización | Explícita mediante estructura ARMA |
| Selección automática | Combinaciones E-T-S | Órdenes p-d-q y P-D-Q |
| Incertidumbre | Distribución predictiva del modelo de innovaciones | Distribución predictiva del ARIMA ajustado |
| Diagnóstico | Innovaciones | Innovaciones |
Existen relaciones matemáticas importantes entre ciertos modelos de suavizamiento y modelos ARIMA. Por ejemplo, SES está estrechamente vinculado con un ARIMA\((0,1,1)\) bajo una parametrización apropiada. Sin embargo, esa equivalencia parcial no convierte a las familias completas en idénticas.
Para este curso nos interesa una conclusión más operativa:
ETS y ARIMA son candidatos que deben compararse bajo los mismos orígenes, horizontes y datos disponibles, en lugar de escogerse por una regla universal.
La aplicación de Sección 8.19 implementa exactamente esa comparación.
Una ventaja decisiva del modelo ETS respecto del algoritmo aislado es que podemos trabajar con
\[ p(y_{T+h}\mid y_{1:T},\widehat\theta), \]
además del pronóstico puntual.
En modelos puramente aditivos, la propagación de innovaciones suele producir expresiones analíticas manejables para la varianza predictiva. En otros modelos, especialmente cuando aparecen componentes multiplicativos, la distribución futura puede ser más compleja y la simulación de trayectorias constituye una alternativa natural (Hyndman y Athanasopoulos 2021, sec. 8.7).
Es importante distinguir tres objetos:
Con componentes multiplicativos, estos objetos no tienen por qué coincidir exactamente debido a la no linealidad. Por eso conviene utilizar directamente las distribuciones generadas por el modelo cuando la meta es evaluar probabilísticamente el pronóstico.
El ajuste ETS que hemos utilizado es fundamentalmente frecuentista por máxima verosimilitud: los parámetros de suavizamiento y los estados iniciales se estiman y posteriormente se utilizan para construir pronósticos.
La presencia de estados no convierte automáticamente al análisis en bayesiano.
En una formulación bayesiana de un modelo dinámico especificaríamos, por ejemplo,
\[ p(\theta_0), \qquad p(\theta_t\mid\theta_{t-1}), \qquad p(y_t\mid\theta_t), \]
y actualizaríamos distribuciones completas sobre estados y parámetros conforme llega información nueva. Los modelos lineales dinámicos de Prado, Ferreira y West y de Petris, Petrone y Campagnoli desarrollan precisamente esa lógica (Prado et al. 2021, chap. 4; Petris et al. 2009, chaps. 2–4).
La relación conceptual puede resumirse así:
\[ \text{ETS esta semana} \quad\longrightarrow\quad \text{estados + innovación + actualización} \quad\longrightarrow\quad \text{DLM y filtro de Kalman en semanas posteriores}. \]
Petris et al. destacan que algoritmos como el promedio móvil exponencial y Holt pueden verse como procedimientos empíricos que adquieren una justificación probabilística al derivarse de modelos dinámicos (Petris et al. 2009, sec. 3.1.1).
No construiremos un “ETS bayesiano” mediante MCMC sólo para mantener simetría con semanas anteriores. El objetivo es comprender la formulación de innovaciones y la evaluación predictiva. La inferencia bayesiana de estados latentes y parámetros desconocidos se desarrollará con mayor coherencia en el bloque de modelos lineales dinámicos.
Las recursiones necesitan comenzar en algún punto. Para SES aparece \(\ell_0\); con tendencia aparece \(b_0\); con estacionalidad necesitamos un conjunto de estados estacionales iniciales.
Una estrategia ingenua sería fijarlos manualmente y olvidarse de ellos. Sin embargo, en muestras cortas pueden afectar apreciablemente:
Por ello, en el marco ETS los estados iniciales forman parte del problema de estimación (Hyndman y Athanasopoulos 2021, sec. 8.6).
Esto también explica por qué dos implementaciones aparentemente equivalentes pueden producir resultados ligeramente diferentes si utilizan reglas de inicialización distintas.
Un ETS bien especificado debería producir innovaciones sin dependencia temporal importante. Sin embargo, debemos evitar dos interpretaciones excesivas.
Como en ARIMA, un buen diagnóstico indica que no encontramos una estructura residual obvia que el modelo haya omitido. No demuestra que la estructura de estados sea el mecanismo verdadero que generó los datos.
\(\ell_t\), \(b_t\) y \(s_t\) son estados definidos por el modelo. Dos modelos distintos pueden descomponer la misma serie de formas diferentes y aun producir pronósticos similares.
La utilidad de los estados debe evaluarse por:
stopifnot(all(datos$y > 0))No debemos permitir que una función automática elija una formulación multiplicativa sin haber pensado qué significa en la escala del problema.
La siguiente secuencia sería incorrecta:
# NO hacer esto:
# 1. ajustar varios ETS a toda la serie;
# 2. escoger el que se ve mejor en los últimos años;
# 3. declarar esos mismos años como período de prueba.La selección pertenece al procedimiento de entrenamiento.
Cuando evaluamos un procedimiento automático,
datos_cv |>
model(ETS(y)) |>
forecast(h = h)la selección y estimación se repiten dentro de cada ventana. Eso representa mucho mejor lo que habría ocurrido históricamente que seleccionar una única especificación con toda la muestra.
La validación temporal puede requerir decenas o cientos de ajustes. En Quarto con knitr, puede ser útil emplear
#| cache: trueen los chunks cuyo resultado depende únicamente de objetos previos estables. Debemos recordar que, si cambia el código o los datos de los que depende el chunk, la caché debe invalidarse para evitar resultados obsoletos.
Si dos modelos se ajustan en escalas distintas, por ejemplo uno a \(y_t\) y otro a \(\log y_t\), las verosimilitudes no siempre son directamente comparables sin considerar el Jacobiano y la definición exacta utilizada por el software. La evaluación fuera de muestra en una escala común evita muchas comparaciones ambiguas.
La actividad propone reconstruir progresivamente la aplicación de turismo sin depender de una única llamada automática.
y, alpha y l0.tourism_train.ETS(Trips).report() y tidy() para identificar la estructura seleccionada y los parámetros.gg_tsresiduals().step_cv_ets de 2 a 1.Los principales resultados del capítulo pueden resumirse de la siguiente manera.
SES asigna pesos geométricamente decrecientes al pasado:
\[ w_j=\alpha(1-\alpha)^j. \]
La actualización de SES puede escribirse como corrección por innovación:
\[ \ell_t=\ell_{t-1}+\alpha e_t. \]
Holt añade una pendiente local \(b_t\) y produce pronósticos lineales en el horizonte.
El amortiguamiento reemplaza la extrapolación lineal indefinida por
\[ \sum_{j=1}^{h}\phi^j b_t, \qquad 0<\phi<1. \]
Holt–Winters añade un estado estacional. La versión aditiva es natural cuando la amplitud estacional es aproximadamente constante; la multiplicativa representa efectos proporcionales al nivel.
Los métodos de suavizamiento son algoritmos de pronóstico. ETS añade una estructura probabilística mediante ecuaciones de observación y estado.
En ETS\((A,N,N)\),
\[ Y_t=\ell_{t-1}+\varepsilon_t, \qquad \ell_t=\ell_{t-1}+\alpha\varepsilon_t, \]
y la recursión SES aparece como consecuencia directa.
La verosimilitud ETS puede construirse secuencialmente con las densidades predictivas de un paso.
Los parámetros de suavizamiento y los estados iniciales forman parte del problema de estimación.
AICc ayuda a seleccionar dentro de muestra, pero no reemplaza la evaluación temporal.
ETS y ARIMA describen la dinámica desde perspectivas diferentes. Deben compararse bajo el mismo experimento predictivo.
La formulación ETS introduce de manera natural las ideas de estado, innovación y actualización secuencial, que retomaremos cuando estudiemos formalmente espacio-estado y filtro de Kalman.
Pesos SES. Partiendo de Ecuación 8.3, derive Ecuación 8.5 por sustitución sucesiva.
Suma de pesos. Demuestre Ecuación 8.7 y verifique que los pesos de \(y_1,\ldots,y_t\) y \(\ell_0\) suman uno.
Caso límite. Muestre qué ocurre con SES cuando \(\alpha=1\). ¿Qué método de referencia se recupera?
Memoria lenta. Explique qué ocurre cuando \(\alpha\) se aproxima a cero. ¿Significa esto que el modelo no utiliza la observación nueva?
Corrección por innovación. Derive Ecuación 8.9 y explique cada uno de sus términos.
Cambio de nivel. Compare cualitativamente la respuesta esperada de \(\alpha=0.1\) y \(\alpha=0.9\) después de un salto permanente del nivel.
Holt. Partiendo de Ecuación 8.12 y Ecuación 8.13, explique por qué \(b_t\) no es la pendiente de una recta global.
Amortiguamiento. Derive Ecuación 8.17 a partir de la suma de una progresión geométrica.
Límite. Demuestre Ecuación 8.18 para \(0<\phi<1\).
Sin amortiguamiento. Muestre que \(\phi=1\) recupera el pronóstico lineal de Holt.
Estacionalidad aditiva. Explique las unidades de \(s_t\) en Ecuación 8.19.
Estacionalidad multiplicativa. Explique las unidades de \(s_t\) en Ecuación 8.23.
Desestacionalización. Compare \(y_t-s_{t-m}\) con \(y_t/s_{t-m}\). ¿Qué estructura supone cada operación?
Modelo versus método. Explique por qué SES no especifica por sí solo una distribución predictiva completa.
ETS\((A,N,N)\). Derive SES a partir de Ecuación 8.28 y Ecuación 8.29.
ETS\((M,N,N)\). Derive Ecuación 8.33 y explique por qué \(\varepsilon_t\) es adimensional.
Holt y \(\beta\). Demuestre Ecuación 8.36.
Verosimilitud. Derive Ecuación 8.42 mediante un cambio de variable desde una innovación gaussiana.
SSE y verosimilitud. Explique por qué minimizar SSE puede coincidir con máxima verosimilitud bajo error aditivo gaussiano, pero no necesariamente con error multiplicativo.
AICc. ¿Por qué los estados iniciales estimados deben considerarse al medir la complejidad de un modelo ETS?
Serie sin tendencia. Proponga un contexto en el que SES sea más razonable que Holt. Justifique.
Tendencia reciente. Una serie creció rápidamente durante los últimos seis meses después de varios años estables. Explique por qué Holt sin amortiguamiento podría producir pronósticos extremos a largo plazo.
Amplitud proporcional. Una serie de ventas duplica su nivel y aproximadamente duplica también la amplitud de sus oscilaciones estacionales. Compare una especificación aditiva y una multiplicativa.
Ceros. ¿Qué problema conceptual y numérico aparece al utilizar error multiplicativo en una serie que contiene ceros?
Innovaciones correlacionadas. Un ETS deja una ACF residual con picos importantes en rezagos 1 y 2. ¿Qué indica esto? ¿Puede concluir inmediatamente que ARIMA resolverá el problema?
AICc versus prueba. Un ETS tiene menor AICc que ARIMA en entrenamiento, pero mayor RMSE en prueba. Explique por qué no hay contradicción.
Cobertura. ETS produce intervalos 80% con cobertura 92% y ARIMA con cobertura 81%. ¿Puede concluir que ETS cuantifica mejor la incertidumbre sin conocer la amplitud? Explique.
Estados. Dos ETS producen niveles \(\ell_t\) distintos pero pronósticos parecidos. ¿Implica esto que uno de los niveles está mal estimado? Discuta.
Bayes. Explique por qué tener estados latentes no convierte automáticamente a ETS en un modelo ajustado bayesianamente.
Puente a Kalman. Compare verbalmente Ecuación 8.9 con la idea “estado previo + corrección × innovación”. ¿Qué elementos esperaría encontrar en un filtro de Kalman?
SES manual. Escriba una función en R que implemente SES y devuelva niveles, pronósticos de un paso y errores.
Optimización. Use optim() para minimizar Ecuación 8.10 en una serie simulada. Compare el \(\widehat\alpha\) obtenido con ETS().
Sensibilidad al inicio. Repita el ejercicio anterior con tres valores de \(\ell_0\). Compare el resultado para muestras de tamaño 20 y 200.
Uso de internet. Reproduzca Tabla 8.1 utilizando .step = 1. Compare con .step = 2.
Holt largo. Pronostique WWWusage a 50 pasos con Holt y Holt amortiguado. Discuta la diferencia.
Turismo aditivo. Ajuste ETS\((A,A,A)\) a tourism_train y grafique sus estados.
Turismo multiplicativo. Ajuste ETS\((M,A,M)\) y compare su estado estacional con el ejercicio anterior.
Selección automática. Use report() sobre ETS(Trips). Interprete la estructura elegida sin limitarse a repetir sus letras.
Diagnóstico. Use gg_tsresiduals() para el ETS automático. Identifique cualquier rezago que merezca investigación adicional.
Holdout. Cambie el período de prueba de 12 a 8 trimestres. ¿Cambia la comparación ETS–ARIMA–naïve estacional?
MASE. Verifique manualmente el denominador estacional de MASE para una de las ventanas de Sección 8.19.
Origen móvil. Repita Sección 8.19 con .step = 1. Compare RMSE\((h)\).
Ventana móvil. Reemplace stretch_tsibble() por slide_tsibble(.size = 48, .step = 2). Discuta si una ventana móvil tiene sentido económico para esta serie.
Intervalos. Compare cobertura y amplitud de niveles 80% y 95%.
Winkler. Aplique la puntuación de Winkler de la semana 7 a los tres procedimientos de Sección 8.19.
CRPS. Si el tiempo de cómputo lo permite, compare las distribuciones predictivas mediante CRPS.
Transformación. Ajuste modelos ETS a log(Trips) y transforme los pronósticos a la escala original. Discuta cuidadosamente la diferencia entre transformar una media y calcular una media después de la transformación inversa.
Benchmark. Elimine el naïve estacional de la comparación y escriba una crítica metodológica de ese experimento.
Horizonte. Repita la validación para \(h=4\) y \(h=12\). Explique por qué cambiar el horizonte puede cambiar la conclusión.
Informe breve. Redacte una conclusión de máximo 200 palabras sobre la aplicación de turismo que incluya: modelo de referencia, diseño temporal, métrica principal, incertidumbre y una limitación.
La semana 9 cambia de perspectiva y estudia la serie en el dominio de la frecuencia. Volveremos después a las ideas de estado e innovación cuando iniciemos formalmente los modelos de espacio-estado.