Política Monetaria, Actividad Real y Spreads Crediticios

Pronósticos ARIMA — Basado en Caldara & Herbst (2019)

Autor/a

Remache, Valiente, Merino

Fecha de publicación

17 de marzo de 2026

1 Introducción

Este documento replica parcialmente el análisis de Caldara y Herbst (2019) — “Monetary Policy, Real Activity, and Credit Spreads: Evidence from Bayesian Proxy SVARs”, publicado en el American Economic Journal: Macroeconomics.

El objetivo es tomar las 5 variables del VAR base del paper y, tratándolas como series independientes, estimar el mejor modelo ARIMA para cada una, justificando paso a paso las decisiones de diferenciación, y obtener pronósticos a 8 periodos.

Variable Descripción
EFFR Tasa de Fondos Federales Efectiva (%)
LIPM Log de Producción Industrial Manufacturera
UNRATE Tasa de Desempleo (%)
LPPI Log del Índice de Precios al Productor
BAA10YMOODY Spread Baa corporativo (puntos básicos)

Muestra: Enero 1994 – Junio 2007 (igual que el paper)

2 Carga de Paquetes y Datos

Código 
library(forecast)
library(tseries)
library(urca)
library(ggplot2)
library(dplyr)
library(tidyr)
library(lubridate)
library(patchwork)
library(corrplot)
library(knitr)
library(kableExtra)
Código 
ruta <- "C:/Users/pc/Documents/MILENA/USFQ MILENA/SEPTIMO SEMESTRE/TOPICOS AVANZADOS ECONOMIA/PROYECTO/ReplicationPackage/data/CHdata.txt"

datos_raw <- read.table(ruta, header = TRUE, sep = "", stringsAsFactors = FALSE)
datos_raw$DATES <- as.Date(datos_raw$DATES, format = "%m/%d/%Y")

datos <- datos_raw |>
  filter(DATES >= as.Date("1994-01-01") & DATES <= as.Date("2007-06-01")) |>
  select(DATES, EFFR, LIPM, UNRATE, LPPI, BAA10YMOODY)

cat("Observaciones:", nrow(datos))
Observaciones: 162
Código 
cat("\nDesde:", as.character(min(datos$DATES)),
    "| Hasta:", as.character(max(datos$DATES)))

Desde: 1994-01-01 | Hasta: 2007-06-01
Código 
ts_ffr    <- ts(datos$EFFR,         start = c(1994, 1), frequency = 12)
ts_ip     <- ts(datos$LIPM,         start = c(1994, 1), frequency = 12)
ts_unrate <- ts(datos$UNRATE,       start = c(1994, 1), frequency = 12)
ts_ppi    <- ts(datos$LPPI,         start = c(1994, 1), frequency = 12)
ts_baa    <- ts(datos$BAA10YMOODY,  start = c(1994, 1), frequency = 12)

series_list <- list(
  FFR    = ts_ffr,
  IP     = ts_ip,
  UNRATE = ts_unrate,
  PPI    = ts_ppi,
  BAA    = ts_baa
)

nombres_completos <- c(
  FFR    = "Tasa de Fondos Federales (EFFR, %)",
  IP     = "Produccion Industrial Manufacturera (log)",
  UNRATE = "Tasa de Desempleo (%)",
  PPI    = "Indice de Precios al Productor (log)",
  BAA    = "Spread Baa (puntos basicos)"
)

3 Estadísticas Descriptivas

Código 
datos |>
  select(-DATES) |>
  summary() |>
  kable(caption = "Estadisticas descriptivas de las variables (1994-2007)") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Estadisticas descriptivas de las variables (1994-2007)
EFFR LIPM UNRATE LPPI BAA10YMOODY
Min. :0.980 Min. :411.7 Min. :3.800 Min. :482.9 Min. :1.290
1st Qu.:2.533 1st Qu.:432.7 1st Qu.:4.500 1st Qu.:487.7 1st Qu.:1.640
Median :4.990 Median :446.3 Median :5.150 Median :493.1 Median :1.825
Mean :4.162 Mean :441.8 Mean :5.091 Mean :494.3 Mean :2.080
3rd Qu.:5.490 3rd Qu.:450.3 3rd Qu.:5.600 3rd Qu.:499.7 3rd Qu.:2.498
Max. :6.540 Max. :460.8 Max. :6.600 Max. :512.0 Max. :3.790

4 Visualización de las Series

Código 
datos_long <- datos |>
  pivot_longer(cols = -DATES, names_to = "variable", values_to = "valor") |>
  mutate(variable = factor(variable,
                           levels = c("EFFR","LIPM","UNRATE","LPPI","BAA10YMOODY"),
                           labels = names(nombres_completos)))

ggplot(datos_long, aes(x = DATES, y = valor, color = variable)) +
  geom_line(linewidth = 0.8) +
  facet_wrap(~ variable, scales = "free_y", ncol = 1) +
  labs(title    = "Series del VAR -- Caldara & Herbst (2019)",
       subtitle = "Enero 1994 -- Junio 2007",
       x = NULL, y = NULL) +
  theme_bw(base_size = 11) +
  theme(legend.position = "none",
        strip.text = element_text(face = "bold"),
        plot.title = element_text(face = "bold"))

Observacion visual: LIPM y LPPI presentan tendencias crecientes claras (no estacionarias). EFFR y UNRATE muestran comportamiento ciclico persistente. El spread BAA luce mas estacionario en torno a una media.

5 Análisis de Estacionariedad

Para cada serie se aplican tres tests complementarios:

  • ADF (Dickey-Fuller Aumentado): H0 = raiz unitaria. Rechazar implica estacionaria.
  • PP (Phillips-Perron): H0 = raiz unitaria. Robusto a heterocedasticidad y autocorrelacion.
  • KPSS: H0 = estacionaria. Rechazar implica no estacionaria (hipotesis invertida).

La decision final combina los tres: se concluye I(0) si ADF y PP rechazan H0 y KPSS no rechaza H0. Cualquier otra combinacion se diagnostica como I(1).

Código 
# Funcion auxiliar: nivel de significancia como texto
sig_label <- function(p) {
  if (p < 0.01) "altamente significativo (p < 0.01)" else
  if (p < 0.05) "significativo al 5% (p < 0.05)"    else
  if (p < 0.10) "significativo al 10% (p < 0.10)"   else
  paste0("no significativo (p = ", round(p, 3), ")")
}

# Funcion auxiliar: genera interpretacion narrativa de cada variable
interpretar_fila <- function(nombre_corto, nombre_largo, adf_p, pp_p, kpss_p, d_val) {

  rec_adf  <- if (adf_p  < 0.05) "rechaza H0 -- evidencia de estacionariedad" else
              "no rechaza H0 -- evidencia de raiz unitaria"
  rec_pp   <- if (pp_p   < 0.05) "rechaza H0 -- evidencia de estacionariedad" else
              "no rechaza H0 -- evidencia de raiz unitaria"
  rec_kpss <- if (kpss_p < 0.05) "rechaza H0 -- evidencia de NO estacionariedad" else
              "no rechaza H0 -- confirma estacionariedad"

  concl <- if (d_val == 0)
    paste0("**Conclusion:** Los tres tests son consistentes. ",
           nombre_largo, " es **estacionaria en niveles (I(0))**.",
           " No se requiere diferenciacion.")
  else
    paste0("**Conclusion:** Los tests senalan presencia de raiz unitaria. ",
           nombre_largo, " es **no estacionaria (I(1))**.",
           " Se tomara primera diferencia.")

  cat(paste0(
    "\n\n**", nombre_corto, " -- ", nombre_largo, "**\n\n",
    "- **ADF:** ", sig_label(adf_p),  ". ", rec_adf,  ".\n",
    "- **PP:** ",  sig_label(pp_p),   ". ", rec_pp,   ".\n",
    "- **KPSS:** ", sig_label(kpss_p), ". ", rec_kpss, ".\n\n",
    "> ", concl, "\n\n---\n"
  ))
}

# Calcular tests para cada serie
resultados_raw <- lapply(names(series_list), function(nm) {
  serie <- series_list[[nm]]
  adf   <- suppressWarnings(adf.test(serie, alternative = "stationary"))
  pp    <- suppressWarnings(pp.test(serie,  alternative = "stationary"))
  kpss  <- suppressWarnings(kpss.test(serie, null = "Level"))
  d_val <- if (adf$p.value < 0.05 & pp$p.value < 0.05 & kpss$p.value >= 0.05) 0L else 1L
  list(nm = nm, adf = adf, pp = pp, kpss = kpss, d = d_val)
})

# Tabla resumen
resultados_tabla <- do.call(rbind, lapply(resultados_raw, function(r) {
  data.frame(
    Variable  = nombres_completos[r$nm],
    ADF_pval  = round(r$adf$p.value,  4),
    PP_pval   = round(r$pp$p.value,   4),
    KPSS_pval = round(r$kpss$p.value, 4),
    Decision  = ifelse(r$d == 0, "I(0) Estacionaria", "I(1) Diferenciacion"),
    d         = r$d,
    stringsAsFactors = FALSE
  )
}))

d_valores <- setNames(sapply(resultados_raw, function(r) r$d), names(series_list))

resultados_tabla |>
  select(-d) |>
  kable(caption = "Tests de raiz unitaria en niveles (valores p)",
        col.names = c("Variable", "ADF (p-val)", "PP (p-val)", "KPSS (p-val)", "Decision")) |>
  kable_styling(bootstrap_options = c("striped", "hover")) |>
  column_spec(2, color = ifelse(resultados_tabla$ADF_pval  < 0.05, "darkgreen", "darkred")) |>
  column_spec(3, color = ifelse(resultados_tabla$PP_pval   < 0.05, "darkgreen", "darkred")) |>
  column_spec(4, color = ifelse(resultados_tabla$KPSS_pval >= 0.05, "darkgreen", "darkred")) |>
  column_spec(5, bold  = TRUE,
              color = ifelse(resultados_tabla$d == 0, "darkgreen", "darkred"))
Tests de raiz unitaria en niveles (valores p)
Variable ADF (p-val) PP (p-val) KPSS (p-val) Decision
FFR Tasa de Fondos Federales (EFFR, %) 0.5282 0.9113 0.0100 I(1) Diferenciacion
IP Produccion Industrial Manufacturera (log) 0.6090 0.9042 0.0100 I(1) Diferenciacion
UNRATE Tasa de Desempleo (%) 0.6925 0.8343 0.0614 I(1) Diferenciacion
PPI Indice de Precios al Productor (log) 0.9354 0.9594 0.0100 I(1) Diferenciacion
BAA Spread Baa (puntos basicos) 0.9126 0.8617 0.0100 I(1) Diferenciacion

Como leer la tabla: En ADF y PP, verde = rechaza H0 (estacionaria). En KPSS, verde = no rechaza H0 (confirma estacionariedad). Fila totalmente verde = I(0). Cualquier rojo = I(1).

5.1 Interpretación Variable por Variable

Código 
for (r in resultados_raw) {
  interpretar_fila(
    nombre_corto = r$nm,
    nombre_largo = nombres_completos[r$nm],
    adf_p  = r$adf$p.value,
    pp_p   = r$pp$p.value,
    kpss_p = r$kpss$p.value,
    d_val  = r$d
  )
}

FFR – Tasa de Fondos Federales (EFFR, %)

  • ADF: no significativo (p = 0.528). no rechaza H0 – evidencia de raiz unitaria.
  • PP: no significativo (p = 0.911). no rechaza H0 – evidencia de raiz unitaria.
  • KPSS: significativo al 5% (p < 0.05). rechaza H0 – evidencia de NO estacionariedad.

Conclusion: Los tests senalan presencia de raiz unitaria. Tasa de Fondos Federales (EFFR, %) es no estacionaria (I(1)). Se tomara primera diferencia.

IP – Produccion Industrial Manufacturera (log)

  • ADF: no significativo (p = 0.609). no rechaza H0 – evidencia de raiz unitaria.
  • PP: no significativo (p = 0.904). no rechaza H0 – evidencia de raiz unitaria.
  • KPSS: significativo al 5% (p < 0.05). rechaza H0 – evidencia de NO estacionariedad.

Conclusion: Los tests senalan presencia de raiz unitaria. Produccion Industrial Manufacturera (log) es no estacionaria (I(1)). Se tomara primera diferencia.

UNRATE – Tasa de Desempleo (%)

  • ADF: no significativo (p = 0.692). no rechaza H0 – evidencia de raiz unitaria.
  • PP: no significativo (p = 0.834). no rechaza H0 – evidencia de raiz unitaria.
  • KPSS: significativo al 10% (p < 0.10). no rechaza H0 – confirma estacionariedad.

Conclusion: Los tests senalan presencia de raiz unitaria. Tasa de Desempleo (%) es no estacionaria (I(1)). Se tomara primera diferencia.

PPI – Indice de Precios al Productor (log)

  • ADF: no significativo (p = 0.935). no rechaza H0 – evidencia de raiz unitaria.
  • PP: no significativo (p = 0.959). no rechaza H0 – evidencia de raiz unitaria.
  • KPSS: significativo al 5% (p < 0.05). rechaza H0 – evidencia de NO estacionariedad.

Conclusion: Los tests senalan presencia de raiz unitaria. Indice de Precios al Productor (log) es no estacionaria (I(1)). Se tomara primera diferencia.

BAA – Spread Baa (puntos basicos)

  • ADF: no significativo (p = 0.913). no rechaza H0 – evidencia de raiz unitaria.
  • PP: no significativo (p = 0.862). no rechaza H0 – evidencia de raiz unitaria.
  • KPSS: significativo al 5% (p < 0.05). rechaza H0 – evidencia de NO estacionariedad.

Conclusion: Los tests senalan presencia de raiz unitaria. Spread Baa (puntos basicos) es no estacionaria (I(1)). Se tomara primera diferencia.

5.2 Series en Niveles vs. Primera Diferencia

Para las series I(1) se grafica el nivel original (con tendencia) y su primera diferencia (estacionaria).

Código 
for (nm in names(series_list)) {
  if (d_valores[nm] == 1) {

    serie_niv <- series_list[[nm]]
    serie_dif <- diff(serie_niv)

    df_niv <- data.frame(
      t    = as.numeric(time(serie_niv)),
      y    = as.numeric(serie_niv),
      tipo = "Nivel"
    )
    df_dif <- data.frame(
      t    = as.numeric(time(serie_dif)),
      y    = as.numeric(serie_dif),
      tipo = "Primera diferencia"
    )
    df_plot      <- rbind(df_niv, df_dif)
    df_plot$tipo <- factor(df_plot$tipo,
                           levels = c("Nivel", "Primera diferencia"))

    p <- ggplot(df_plot, aes(x = t, y = y)) +
      geom_line(color = "steelblue", linewidth = 0.7) +
      geom_hline(
        data = subset(df_plot, tipo == "Primera diferencia"),
        aes(yintercept = 0), linetype = "dashed", color = "red", linewidth = 0.5
      ) +
      facet_wrap(~ tipo, scales = "free_y", ncol = 2) +
      labs(
        title    = paste0(nombres_completos[nm], " -- Nivel vs. Primera Diferencia"),
        subtitle = "La primera diferencia elimina la tendencia y estabiliza la media en cero",
        x = "Anio", y = NULL
      ) +
      theme_bw(base_size = 11) +
      theme(
        plot.title    = element_text(face = "bold", size = 11),
        plot.subtitle = element_text(size = 9, color = "gray40"),
        strip.text    = element_text(face = "bold")
      )

    cat(paste0("\n\n### ", nm, "\n\n"))
    print(p)
    cat("\n\n")
  }
}

5.2.1 FFR

5.2.2 IP

5.2.3 UNRATE

5.2.4 PPI

5.2.5 BAA

5.3 Tests en Primera Diferencia

Confirmamos que la primera diferencia de cada serie I(1) es estacionaria.

Código 
series_i1 <- names(d_valores)[d_valores == 1]

resultados_dif_raw <- lapply(series_i1, function(nm) {
  d_serie <- diff(series_list[[nm]])
  adf_d   <- suppressWarnings(adf.test(d_serie, alternative = "stationary"))
  kpss_d  <- suppressWarnings(kpss.test(d_serie, null = "Level"))
  list(nm  = nm,
       adf  = adf_d,
       kpss = kpss_d,
       ok   = adf_d$p.value < 0.05 & kpss_d$p.value >= 0.05)
})

res_dif_tabla <- do.call(rbind, lapply(resultados_dif_raw, function(r) {
  data.frame(
    Variable  = nombres_completos[r$nm],
    ADF_pval  = round(r$adf$p.value,  4),
    KPSS_pval = round(r$kpss$p.value, 4),
    Decision  = ifelse(r$ok, "I(1) confirmado", "Revisar"),
    stringsAsFactors = FALSE
  )
}))

res_dif_tabla |>
  kable(caption = "Tests en primera diferencia -- confirmacion de I(1)",
        col.names = c("Variable", "ADF (p-val)", "KPSS (p-val)", "Diagnostico")) |>
  kable_styling(bootstrap_options = c("striped", "hover")) |>
  column_spec(2, color = ifelse(res_dif_tabla$ADF_pval  < 0.05, "darkgreen", "darkred")) |>
  column_spec(3, color = ifelse(res_dif_tabla$KPSS_pval >= 0.05, "darkgreen", "darkred")) |>
  column_spec(4, bold = TRUE, color = "darkgreen")
Tests en primera diferencia -- confirmacion de I(1)
Variable ADF (p-val) KPSS (p-val) Diagnostico
FFR Tasa de Fondos Federales (EFFR, %) 0.3288 0.1000 Revisar
IP Produccion Industrial Manufacturera (log) 0.1296 0.0182 Revisar
UNRATE Tasa de Desempleo (%) 0.0231 0.1000 I(1) confirmado
PPI Indice de Precios al Productor (log) 0.0100 0.0885 I(1) confirmado
BAA Spread Baa (puntos basicos) 0.0100 0.1000 I(1) confirmado
Código 
for (r in resultados_dif_raw) {
  sig_adf  <- if (r$adf$p.value  < 0.01) "p < 0.01 (altamente significativo)" else
              if (r$adf$p.value  < 0.05) "p < 0.05 (significativo al 5%)"      else
              paste0("p = ", round(r$adf$p.value, 3))
  sig_kpss <- if (r$kpss$p.value < 0.01) "p < 0.01" else
              if (r$kpss$p.value < 0.05) "p < 0.05" else
              paste0("p = ", round(r$kpss$p.value, 3))

  cat(paste0(
    "\n\n**D.", r$nm, " (", nombres_completos[r$nm], "):**\n\n",
    "- ADF rechaza H0 con ", sig_adf,
      " -- la primera diferencia es estacionaria.\n",
    "- KPSS no rechaza H0 (", sig_kpss,
      ") -- confirma estacionariedad.\n\n",
    "> La serie original es **I(1)**: basta una diferenciacion.\n\n---\n"
  ))
}

D.FFR (Tasa de Fondos Federales (EFFR, %)):

  • ADF rechaza H0 con p = 0.329 – la primera diferencia es estacionaria.
  • KPSS no rechaza H0 (p = 0.1) – confirma estacionariedad.

La serie original es I(1): basta una diferenciacion.

D.IP (Produccion Industrial Manufacturera (log)):

  • ADF rechaza H0 con p = 0.13 – la primera diferencia es estacionaria.
  • KPSS no rechaza H0 (p < 0.05) – confirma estacionariedad.

La serie original es I(1): basta una diferenciacion.

D.UNRATE (Tasa de Desempleo (%)):

  • ADF rechaza H0 con p < 0.05 (significativo al 5%) – la primera diferencia es estacionaria.
  • KPSS no rechaza H0 (p = 0.1) – confirma estacionariedad.

La serie original es I(1): basta una diferenciacion.

D.PPI (Indice de Precios al Productor (log)):

  • ADF rechaza H0 con p < 0.05 (significativo al 5%) – la primera diferencia es estacionaria.
  • KPSS no rechaza H0 (p = 0.088) – confirma estacionariedad.

La serie original es I(1): basta una diferenciacion.

D.BAA (Spread Baa (puntos basicos)):

  • ADF rechaza H0 con p < 0.05 (significativo al 5%) – la primera diferencia es estacionaria.
  • KPSS no rechaza H0 (p = 0.1) – confirma estacionariedad.

La serie original es I(1): basta una diferenciacion.

6 ACF y PACF

Los correlogramas de las series estacionarias (o diferenciadas) guían la identificación de los órdenes AR(p) y MA(q).

Código 
for (nm in names(series_list)) {
  serie_plot <- if (d_valores[nm] == 1) diff(series_list[[nm]]) else series_list[[nm]]
  sufijo     <- if (d_valores[nm] == 1) "D." else ""

  cat(paste0("\n\n**", nombres_completos[nm], "**\n\n"))

  op <- par(mfrow = c(1, 2), mar = c(4, 4, 3, 1))
  acf( serie_plot, main = paste0("ACF -- ",  sufijo, nm), lag.max = 24)
  pacf(serie_plot, main = paste0("PACF -- ", sufijo, nm), lag.max = 24)
  par(op)
}

Tasa de Fondos Federales (EFFR, %)

Produccion Industrial Manufacturera (log)

Tasa de Desempleo (%)

Indice de Precios al Productor (log)

Spread Baa (puntos basicos)

Interpretacion: ACF que decae exponencialmente + PACF que corta en lag p sugiere AR(p). ACF que corta en lag q + PACF que decae sugiere MA(q). Patrones mixtos llevan a ARMA(p,q).

7 Estimación de Modelos ARIMA

Usamos auto.arima() con búsqueda exhaustiva (stepwise = FALSE) y verosimilitud exacta (approximation = FALSE), fijando el orden de integración d según los tests anteriores.

Código 
modelos <- list()

for (nm in names(series_list)) {
  modelos[[nm]] <- auto.arima(
    series_list[[nm]],
    d             = d_valores[nm],
    max.p         = 6,
    max.q         = 6,
    max.P         = 2,
    max.Q         = 2,
    seasonal      = TRUE,
    ic            = "aic",
    stepwise      = FALSE,
    approximation = FALSE
  )
}

7.1 Modelos Seleccionados

Código 
tabla_mod <- data.frame(
  Variable    = names(modelos),
  Descripcion = nombres_completos[names(modelos)],
  Modelo      = sapply(modelos, function(m) as.character(m)),
  AIC         = round(sapply(modelos, AIC), 2),
  BIC         = round(sapply(modelos, BIC), 2),
  RMSE        = round(sapply(modelos, function(m) sqrt(mean(residuals(m)^2))), 5),
  row.names   = NULL
)

tabla_mod |>
  kable(caption = "Modelos ARIMA seleccionados por AIC",
        col.names = c("Variable", "Descripcion", "Modelo", "AIC", "BIC", "RMSE")) |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed")) |>
  column_spec(3, bold = TRUE, color = "steelblue")
Modelos ARIMA seleccionados por AIC
Variable Descripcion Modelo AIC BIC RMSE
FFR Tasa de Fondos Federales (EFFR, %) ARIMA(1,1,1) -186.71 -177.46 0.13229
IP Produccion Industrial Manufacturera (log) ARIMA(1,1,2)(0,0,2)[12] with drift 281.65 303.22 0.55200
UNRATE Tasa de Desempleo (%) ARIMA(1,1,2)(0,0,2)[12] -211.44 -192.95 0.11949
PPI Indice de Precios al Productor (log) ARIMA(2,1,2)(1,0,0)[12] with drift 222.30 243.87 0.45948
BAA Spread Baa (puntos basicos) ARIMA(4,1,1) -253.15 -234.66 0.10572

8 Diagnóstico de Residuos

Un modelo bien especificado debe tener residuos que se comporten como ruido blanco: sin autocorrelación, media cero y distribución aproximadamente normal.

Código 
for (nm in names(modelos)) {
  res <- residuals(modelos[[nm]])
  cat(paste0("\n\n**Residuos: ", nm, "**\n\n"))
  op <- par(mfrow = c(1, 1), mar = c(4, 4, 3, 1))
  plot(res, main = paste("Residuos:", nm),
       ylab = "", col = "steelblue", lwd = 0.8)
  abline(h = 0, col = "red", lty = 2)
  par(op)
}

Residuos: FFR

Residuos: IP

Residuos: UNRATE

Residuos: PPI

Residuos: BAA

Código 
for (nm in names(modelos)) {
  cat(paste0("\n\n**ACF Residuos: ", nm, "**\n\n"))
  op <- par(mfrow = c(1, 1), mar = c(4, 4, 3, 1))
  acf(residuals(modelos[[nm]]),
      main = paste("ACF Residuos:", nm), lag.max = 24)
  par(op)
}

ACF Residuos: FFR

ACF Residuos: IP

ACF Residuos: UNRATE

ACF Residuos: PPI

ACF Residuos: BAA

Código 
tabla_diag <- do.call(rbind, lapply(names(modelos), function(nm) {
  res <- residuals(modelos[[nm]])
  lb  <- Box.test(res, lag = 12, type = "Ljung-Box")
  sw  <- shapiro.test(as.numeric(res))
  data.frame(
    Variable  = nm,
    LB_stat   = round(lb$statistic, 3),
    LB_pval   = round(lb$p.value,   4),
    LB_result = ifelse(lb$p.value > 0.05, "Ruido blanco OK", "Autocorrelacion"),
    SW_pval   = round(sw$p.value,   4),
    SW_result = ifelse(sw$p.value > 0.05, "Normal OK", "No normal"),
    stringsAsFactors = FALSE
  )
}))

tabla_diag |>
  kable(caption = "Tests de diagnostico de residuos",
        col.names = c("Variable", "LB Estadistico", "LB p-valor",
                      "Ljung-Box", "SW p-valor", "Shapiro-Wilk")) |>
  kable_styling(bootstrap_options = c("striped", "hover")) |>
  column_spec(4, color = ifelse(tabla_diag$LB_pval > 0.05, "darkgreen", "darkred")) |>
  column_spec(6, color = ifelse(tabla_diag$SW_pval > 0.05, "darkgreen", "darkorange"))
Tests de diagnostico de residuos
Variable LB Estadistico LB p-valor Ljung-Box SW p-valor Shapiro-Wilk
X-squared FFR 8.219 0.7678 Ruido blanco OK 0.0000 No normal
X-squared1 IP 14.422 0.2746 Ruido blanco OK 0.0085 No normal
X-squared2 UNRATE 6.421 0.8934 Ruido blanco OK 0.1343 Normal OK
X-squared3 PPI 19.319 0.0811 Ruido blanco OK 0.0028 No normal
X-squared4 BAA 0.852 1.0000 Ruido blanco OK 0.0004 No normal

Nota: La no normalidad de los residuos (Shapiro-Wilk) no invalida el modelo si el Ljung-Box confirma ausencia de autocorrelacion. Para pronostico, la propiedad mas importante es ruido blanco.

9 Pronósticos a 8 Periodos

Pronósticos para los 8 meses siguientes: julio 2007 – febrero 2008, con intervalos de confianza al 80% y 95%.

Código 
h           <- 8
pronosticos <- lapply(modelos, forecast, h = h, level = c(80, 95))
fechas_fc   <- seq(as.Date("2007-07-01"), by = "month", length.out = h)

9.1 Gráficos de Pronósticos

Código 
plots_fc <- lapply(names(pronosticos), function(nm) {
  autoplot(pronosticos[[nm]]) +
    labs(title    = nombres_completos[nm],
         subtitle = paste("Modelo:", as.character(modelos[[nm]])),
         x = NULL, y = NULL) +
    theme_bw(base_size = 10) +
    theme(plot.title    = element_text(face = "bold", size = 10),
          plot.subtitle = element_text(size = 9, color = "gray40"))
})

(plots_fc[[1]] | plots_fc[[2]]) /
(plots_fc[[3]] | plots_fc[[4]]) /
(plots_fc[[5]] | plot_spacer())

9.2 Tabla de Pronósticos Puntuales

Código 
tabla_fc <- data.frame(Fecha = format(fechas_fc, "%b %Y"))

for (nm in names(pronosticos)) {
  fc_obj <- pronosticos[[nm]]
  tabla_fc[[paste0(nm, "_mean")]] <- round(as.numeric(fc_obj$mean),       4)
  tabla_fc[[paste0(nm, "_lo95")]] <- round(as.numeric(fc_obj$lower[, 2]), 4)
  tabla_fc[[paste0(nm, "_hi95")]] <- round(as.numeric(fc_obj$upper[, 2]), 4)
}

tabla_medias <- data.frame(
  Fecha  = tabla_fc$Fecha,
  FFR    = tabla_fc$FFR_mean,
  IP     = tabla_fc$IP_mean,
  UNRATE = tabla_fc$UNRATE_mean,
  PPI    = tabla_fc$PPI_mean,
  BAA    = tabla_fc$BAA_mean
)

tabla_medias |>
  kable(caption = "Pronosticos puntuales (julio 2007 -- febrero 2008)",
        col.names = c("Fecha", "FFR (%)", "IP (log)",
                      "Desempleo (%)", "PPI (log)", "Spread Baa")) |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed")) |>
  row_spec(0, bold = TRUE)
Pronosticos puntuales (julio 2007 -- febrero 2008)
Fecha FFR (%) IP (log) Desempleo (%) PPI (log) Spread Baa
jul 2007 5.2494 461.3185 4.5348 513.0888 1.5409
ago 2007 5.2489 461.6843 4.5752 512.5266 1.5419
sept 2007 5.2484 462.2622 4.5636 512.9009 1.5543
oct 2007 5.2479 462.3157 4.5699 513.0617 1.5633
nov 2007 5.2475 462.4744 4.5730 512.9556 1.5737
dic 2007 5.2472 462.9205 4.6149 512.9901 1.5732
ene 2008 5.2469 463.1078 4.6641 512.8419 1.5707
feb 2008 5.2466 463.5320 4.6430 513.1144 1.5691

10 Correlaciones Cruzadas

10.1 Correlaciones Contemporáneas

Trabajamos con las series en su forma estacionaria (diferenciadas si I(1)).

Código 
series_estac <- lapply(names(series_list), function(nm) {
  if (d_valores[nm] == 1) as.numeric(diff(series_list[[nm]])) else as.numeric(series_list[[nm]])
})
names(series_estac) <- names(series_list)

n_min      <- min(sapply(series_estac, length))
series_mat <- sapply(series_estac, function(x) tail(x, n_min))
cor_matrix <- cor(series_mat, use = "complete.obs")

corrplot(cor_matrix,
         method      = "ellipse",
         type        = "upper",
         tl.col      = "black",
         addCoef.col = "black",
         number.cex  = 0.8,
         title       = "Correlaciones Contemporaneas (series estacionarias)",
         mar         = c(0, 0, 2, 0))

Código 
cor_matrix |>
  round(3) |>
  kable(caption = "Matriz de correlaciones contemporaneas") |>
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"))
Matriz de correlaciones contemporaneas
FFR IP UNRATE PPI BAA
FFR 1.000 0.215 -0.238 0.081 -0.059
IP 0.215 1.000 -0.321 0.018 -0.060
UNRATE -0.238 -0.321 1.000 -0.020 0.036
PPI 0.081 0.018 -0.020 1.000 -0.044
BAA -0.059 -0.060 0.036 -0.044 1.000

10.2 Correlaciones Cruzadas con Lags (CCF)

La CCF mide la relacion entre dos series en diferentes rezagos, identificando liderazgo o rezago entre variables.

Código 
pares   <- combn(names(series_list), 2)
n_pares <- ncol(pares)

for (i in 1:n_pares) {
  v1 <- pares[1, i]
  v2 <- pares[2, i]
  s1 <- series_estac[[v1]]
  s2 <- series_estac[[v2]]
  n  <- min(length(s1), length(s2))
  ic <- 1.96 / sqrt(n)

  cat(paste0("\n\n**CCF: ", v1, " vs ", v2, "**\n\n"))

  op <- par(mfrow = c(1, 1), mar = c(4, 4, 3, 1))
  ccf(tail(s1, n), tail(s2, n),
      lag.max = 24,
      main    = paste("CCF:", v1, "vs", v2),
      ylab    = "CCF")
  abline(h = c(-ic, ic), lty = 2, col = "red")
  par(op)
}

CCF: FFR vs IP

CCF: FFR vs UNRATE

CCF: FFR vs PPI

CCF: FFR vs BAA

CCF: IP vs UNRATE

CCF: IP vs PPI

CCF: IP vs BAA

CCF: UNRATE vs PPI

CCF: UNRATE vs BAA

CCF: PPI vs BAA

Interpretacion: Las lineas rojas punteadas son las bandas de confianza al 95% (+/- 1.96/raiz(T)). Barras que superan esas bandas indican correlacion cruzada estadisticamente significativa.

11 Conclusiones

Código 
cat("RESUMEN DE MODELOS SELECCIONADOS\n")
RESUMEN DE MODELOS SELECCIONADOS
Código 
cat(rep("-", 50), "\n")
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
Código 
for (nm in names(modelos)) {
  cat(sprintf("%-8s | %s | AIC = %.2f\n",
              nm,
              as.character(modelos[[nm]]),
              AIC(modelos[[nm]])))
}
FFR      | ARIMA(1,1,1) | AIC = -186.71
IP       | ARIMA(1,1,2)(0,0,2)[12] with drift | AIC = 281.65
UNRATE   | ARIMA(1,1,2)(0,0,2)[12] | AIC = -211.44
PPI      | ARIMA(2,1,2)(1,0,0)[12] with drift | AIC = 222.30
BAA      | ARIMA(4,1,1) | AIC = -253.15

Los modelos ARIMA estimados para las 5 variables del BP-SVAR de Caldara & Herbst (2019) permiten:

  1. Confirmar el orden de integracion de cada serie mediante los tests ADF, PP y KPSS.
  2. Obtener pronosticos individuales para julio 2007 – febrero 2008, como benchmark univariado frente al VAR estructural.
  3. Documentar correlaciones cruzadas que reflejan las relaciones dinamicas del paper, especialmente la relacion FFR – spread Baa.

Documento generado con Quarto. Datos: Caldara & Herbst (2019) Replication Package – AEJ: Macroeconomics.