Los modelos de estado-espacio constituyen un framework unificado para el análisis de series temporales que separa explícitamente la señal subyacente del ruido observacional. Esta metodología implementa inferencia bayesiana completa sobre componentes estructurales de series temporales (Bayesian Structural Time Series - BSTS), combinando la flexibilidad de la 01 con la selección automática de variables mediante priors 02. El enfoque permite cuantificar la incertidumbre en todos los niveles del modelo mientras mantiene interpretabilidad económica de los componentes.
Un modelo de estado-espacio se define mediante dos ecuaciones fundamentales:
Ecuación de Observación: \[y_t = Z_t' \alpha_t + \beta' x_t + \epsilon_t, \quad \epsilon_t \sim N(0, \sigma^2_{\epsilon})\]
Ecuación de Estado (Transición): \[\alpha_{t+1} = T_t \alpha_t + R_t \eta_t, \quad \eta_t \sim N(0, Q_t)\]
donde: - \(y_t\) es la observación en tiempo \(t\) - \(\alpha_t\) es el vector de estado latente de dimensión \(m\) - \(x_t\) es el vector de regresores exógenos de dimensión \(k\) - \(\beta\) son los coeficientes de regresión - \(Z_t\) conecta el estado con las observaciones - \(T_t\) es la matriz de transición del estado - \(R_t\) y \(Q_t\) definen la estructura de varianza del proceso de estado - \(\epsilon_t\) es el error de observación
El modelo más simple incluye solo un nivel estocástico:
\[\begin{aligned} y_t &= \mu_t + \beta' x_t + \epsilon_t \\ \mu_{t+1} &= \mu_t + \eta_{\mu,t} \end{aligned}\]
donde \(\eta_{\mu,t} \sim N(0, \sigma^2_\mu)\) representa innovaciones en el nivel. Este modelo es apropiado para series con nivel variable pero sin tendencia sistemática.
Extiende el modelo LL agregando una tendencia estocástica:
\[\begin{aligned} y_t &= \mu_t + \beta' x_t + \epsilon_t \\ \mu_{t+1} &= \mu_t + \delta_t + \eta_{\mu,t} \\ \delta_{t+1} &= \delta_t + \eta_{\delta,t} \end{aligned}\]
donde: - \(\mu_t\) es el nivel en tiempo \(t\) - \(\delta_t\) es la pendiente (tasa de cambio) en tiempo \(t\) - \(\eta_{\mu,t} \sim N(0, \sigma^2_\mu)\) son innovaciones del nivel - \(\eta_{\delta,t} \sim N(0, \sigma^2_\delta)\) son innovaciones de la pendiente
Para series con patrones estacionales:
\[\gamma_{t+1} = -\sum_{s=1}^{S-1} \gamma_{t-s+1} + \eta_{\gamma,t}\]
donde \(S\) es el período estacional (e.g., 12 para datos mensuales) y \(\eta_{\gamma,t} \sim N(0, \sigma^2_\gamma)\).
En contextos de alta dimensionalidad con múltiples rezagos potenciales, la selección de variables es crucial para: - Evitar sobreajuste - Mantener interpretabilidad - Identificar predictores genuinos
Para cada coeficiente \(\beta_j\) definimos un prior jerárquico:
\[\beta_j | \gamma_j \sim \begin{cases} N(0, \tau^2 v_j) & \text{si } \gamma_j = 1 \text{ (slab)} \\ \delta_0 & \text{si } \gamma_j = 0 \text{ (spike)} \end{cases}\]
donde: - \(\gamma_j \in \{0,1\}\) es el indicador de inclusión - \(\tau^2 v_j\) es la varianza del slab (componente no-nulo) - \(\delta_0\) es una masa puntual en cero
El prior sobre los indicadores es:
\[\gamma_j \sim \text{Bernoulli}(\pi_j)\]
con \(\pi_j\) calibrado según el tamaño esperado del modelo.
Los hiperparámetros se configuran mediante:
Esta configuración balancea flexibilidad con parsimonia, permitiendo que los datos dominen la inferencia mientras se mantiene regularización suave.
La inferencia se realiza mediante un esquema híbrido de Gibbs sampling:
Aunque no se reportan explícitamente por eficiencia computacional, los diagnósticos estándar incluyen: - Inspección visual de trazas - Autocorrelación de las cadenas - Effective sample size para parámetros clave - Probabilidades de inclusión estables para \(\gamma_j\)
Implementamos Leave-Future-Out (LFO) con origen móvil:
Este esquema más conservador (horizonte y paso de 6 meses vs 12 en ECM-MARS) permite mayor resolución temporal en la evaluación de estabilidad.
Para cada fold \(f\):
Para cada observación futura \(y_t^*\):
\[\text{ELPD}_t = \log \int p(y_t^* | \theta) p(\theta | \mathcal{D}_{\text{train}}) d\theta\]
Aproximado mediante muestras MCMC como:
\[\widehat{\text{ELPD}}_t = \log \left( \frac{1}{S} \sum_{s=1}^S p(y_t^* | \theta^{(s)}) \right)\]
donde asumimos normalidad: \(p(y_t^* | \theta^{(s)}) = N(y_t^*; \mu_t^{(s)}, \sigma_t^{(s)})\).
Usando la media posterior como predicción puntual:
donde \(\bar{y}_t = \frac{1}{S} \sum_{s=1}^S y_t^{(s)}\) es la media posterior.
Evaluamos la calidad de la incertidumbre cuantificada:
Si el modelo está bien calibrado, los valores PIT siguen una distribución uniforme en [0,1].
Evaluamos sistemáticamente:
| Modelo | Nivel | Tendencia | Estacionalidad | Parámetros libres |
|---|---|---|---|---|
| LL | Estocástico | No | Opcional | \(\sigma^2_\mu, \sigma^2_\epsilon\) |
| LLT | Estocástico | Estocástica | Opcional | \(\sigma^2_\mu, \sigma^2_\delta, \sigma^2_\epsilon\) |
La estructura óptima se selecciona mediante:
Esta jerarquía prioriza capacidad predictiva probabilística sobre ajuste puntual.
Un modelo completo “gana” en el fold \(f\) si:
\[\begin{cases} \Delta\text{ELPD}_f > 0 & \text{(mejor densidad predictiva)} \\ \Delta\text{RMSE}_f > 0 & \text{(menor error, i.e., RMSE}_{\text{base}} - \text{RMSE}_{\text{full}} > 0\text{)} \end{cases}\]
Nota: A diferencia del GLM-AR(1), aquí \(\Delta\text{RMSE} = \text{RMSE}_{\text{base}} - \text{RMSE}_{\text{full}}\), por lo que valores positivos indican mejora.
El código incluye manejo robusto de:
Cada par genera:
list(
best_summary = tibble( # Resumen del mejor modelo
spec, folds, wins, support,
dELPD_mean, dRMSE_mean, ...
),
best_results = tibble( # Detalles por fold
fold, n_train, n_test,
ELPD_base, ELPD_full, dELPD,
RMSE_base, RMSE_full, dRMSE,
cover80, cover95, pit, ...
),
all_summaries = tibble() # Todas las estructuras probadas
)| Aspecto | ECM-MARS | GLM-AR(1) | BSTS |
|---|---|---|---|
| Paradigma | Frecuentista híbrido | Bayesiano paramétrico | Bayesiano estructural |
| Pre-requisitos | I(1), cointegración | Ninguno | Ninguno |
| No-linealidad | MARS (splines) | No | No (pero componentes flexibles) |
| Selección variables | Manual (pre-tests) | No | Automática (spike-slab) |
| Componentes temporales | ECM únicamente | AR(1) global | Nivel, tendencia, estacionalidad |
| Incertidumbre | Bootstrap implícito | Posterior completa | Posterior completa + estados |
| Interpretabilidad | Alta (económica) | Media | Alta (descomposición) |
| Costo computacional | Medio | Alto | Muy alto |
Descomposición interpretable: Separa señal de ruido con componentes económicamente significativos
Selección automática: Spike-and-slab identifica predictores relevantes sin pre-especificación
Manejo de incertidumbre estructural: Cuantifica incertidumbre en componentes no observables
Robustez a especificación: No requiere supuestos sobre orden de integración o cointegración
Pronósticos por intervalo: Intervalos de predicción naturales y bien calibrados
Costo computacional elevado: FFBS + SSVS es intensivo, especialmente con muchos regresores
Linealidad en predictores: No captura interacciones no lineales como MARS
Diagnósticos complejos: Convergencia más difícil de verificar que MCO
Sensibilidad a priors en muestras pequeñas: Especialmente para varianzas de componentes
ENTRADA: DataFrame con series temporales, especificación de variables
SALIDA: Ranking de pares por capacidad predictiva y estabilidad
# PREPARACIÓN
1. CARGAR datos y limpiar nombres
2. IDENTIFICAR 6 variables circulación, 7 producción
3. GENERAR 84 pares (2 × 6 × 7)
# PROCESAMIENTO POR PAR
PARA cada par (Y, X):
4. CONSTRUIR matriz con Y, X, lags(X, 1:6)
5. ELIMINAR observaciones con NA
# TUNING DE ESTRUCTURA
PARA cada estructura en {LL, LLT}:
# LEAVE-FUTURE-OUT
6. GENERAR splits LFO (init=80%, h=6, step=6)
PARA cada fold:
7. DIVIDIR train/test
8. ESCALAR regresores con stats del train
# MODELOS
9. AJUSTAR modelo_base:
- state_spec = estructura sin regresores
- MCMC con 2000 iter, 500 burn-in
10. AJUSTAR modelo_full:
- state_spec = estructura
- prior = spike_slab(X_lags, expected_size=5)
- MCMC con 2000 iter, 500 burn-in
# PREDICCIÓN
11. GENERAR distribuciones predictivas (h=6)
12. EXTRAER media, sd, cuantiles
# EVALUACIÓN
13. CALCULAR:
- ELPD_base, ELPD_full → ΔELPD
- RMSE_base, RMSE_full → ΔRMSE
- Coverage 80%, 95%
- PIT values
14. DETERMINAR victoria: (ΔELPD > 0) ∧ (ΔRMSE > 0)
15. RESUMIR estructura:
- support = wins/folds
- promedios de métricas
16. SELECCIONAR mejor estructura por ΔELPD
17. GUARDAR resultados del par
# RANKING Y EXPORT
18. ORDENAR pares por (support DESC, ΔELPD DESC, ΔRMSE DESC)
19. FILTRAR ganadores por umbrales (0.70, 0.60)
20. EXPORTAR CSVs y visualizaciones
21. GENERAR gráficos del último fold (sin re-ajustar)La metodología BSTS representa el enfoque más comprehensivo y flexible de los tres implementados, combinando la rigurosidad de la inferencia bayesiana con la interpretabilidad de la descomposición estructural. A diferencia del ECM-MARS que requiere validación previa de supuestos econométricos estrictos, y del GLM-AR(1) que modela solo dependencia serial global, BSTS descompone la serie en componentes interpretables mientras selecciona automáticamente predictores relevantes.
Los resultados empíricos sugieren que BSTS identifica relaciones predictivas robustas con alta precisión, aunque a un costo computacional sustancialmente mayor. La capacidad de cuantificar incertidumbre en todos los niveles —desde componentes estructurales hasta predicciones— lo hace particularmente valioso para aplicaciones donde la evaluación de riesgo es crítica.
La implementación con horizonte de 6 meses (vs 12 en otras metodologías) proporciona evaluación más granular de la estabilidad temporal, revelando que muchas relaciones que parecen estables en horizontes largos muestran variabilidad en escalas más cortas. Esto sugiere la importancia de considerar múltiples horizontes en la evaluación de modelos predictivos.
El framework BSTS es especialmente apropiado cuando: - La interpretación de componentes es prioritaria - Existe incertidumbre sobre qué predictores incluir - Se requieren intervalos de predicción bien calibrados - Los datos exhiben múltiples fuentes de variación (tendencia, estacionalidad, regresores)
La combinación de selección automática de variables con modelado estructural posiciona a BSTS como un puente entre los métodos econométricos clásicos y las técnicas modernas de machine learning, manteniendo interpretabilidad mientras abraza la complejidad inherente en los datos económicos.
La descomposición estructural separa variación sistemática (tendencia, estacionalidad) de variación aleatoria, facilitando interpretación económica y mejorando pronósticos al modelar cada componente apropiadamente.↩︎
El prior spike-and-slab combina una distribución continua (slab) para coeficientes activos con una masa puntual en cero (spike) para coeficientes inactivos, permitiendo selección exacta de variables.↩︎
Forward Filtering Backward Sampling (FFBS) es el algoritmo estándar para muestrear estados en modelos lineales gaussianos, combinando el filtro de Kalman hacia adelante con muestreo hacia atrás condicionado en los datos completos.↩︎
Stochastic Search Variable Selection (SSVS) es un método MCMC para explorar el espacio de modelos, alternando entre actualizar coeficientes dado el modelo y actualizar el modelo dados los coeficientes.↩︎
El esquema de validación con h=6 y step=6 evita solapamiento entre conjuntos de test, garantizando independencia en la evaluación mientras maximiza el uso de datos.↩︎
La calibración de probabilidades predictivas es crucial para decisiones bajo incertidumbre. Un modelo bien calibrado asigna probabilidades que reflejan frecuencias empíricas a largo plazo.↩︎