T11_mST_CP_Macro.R

Ver crudo
# ══════════════════════════════════════════════════
# T11 · Caso práctico Macro — IPC de España (INE)
# Abre primero mSeriesTemporales.Rproj en RStudio
# Datos: data/ipc_mensual_espana.RData (INE, serie IPC290751)
# ══════════════════════════════════════════════════

library(tidyverse)
library(forecast)

pausa <- function(msg = "\n  [Pulsa ENTER para continuar...]") {
  if (interactive()) { cat(msg); invisible(readline()) }
}

recm_f <- function(e) sqrt(mean(e^2, na.rm = TRUE))
eam_f  <- function(e) mean(abs(e), na.rm = TRUE)

# ── 1. CARGA Y ESCALA DEL MASE ────────────────────
load("data/ipc_mensual_espana.RData")
ipc <- ipc_mensual_espana
ipc$fecha <- as.Date(ipc$fecha)
ipc_ts <- ts(ipc$ipc, start = c(2002, 1), frequency = 12)

mae_naive <- mean(abs(diff(ipc_ts, lag = 12)))
mase_f <- function(e) eam_f(e) / mae_naive
cat(sprintf("MAE naive estacional (escala MASE) = %.4f\n", mae_naive))
pausa()

# ── 2. VALIDACION CRUZADA TEMPORAL: SARIMA, ETS, HOLT ─
f_sarima <- function(x, h) forecast(Arima(x, order = c(0, 1, 1),
                                           seasonal = list(order = c(0, 1, 1), period = 12),
                                           method = "CSS"), h = h)
f_ets  <- function(x, h) forecast(ets(x, model = "MAA"), h = h)
f_holt <- function(x, h) forecast::holt(x, h = h)

e_sarima <- tsCV(ipc_ts, f_sarima, h = 1, window = 180)
e_ets    <- tsCV(ipc_ts, f_ets,    h = 1, window = 180)
e_holt   <- tsCV(ipc_ts, f_holt,   h = 1, window = 180)

resultado <- tibble(
  modelo = c("SARIMA", "ETS", "Holt"),
  RECM = c(recm_f(e_sarima), recm_f(e_ets), recm_f(e_holt)),
  EAM  = c(eam_f(e_sarima), eam_f(e_ets), eam_f(e_holt)),
  MASE = c(mase_f(e_sarima), mase_f(e_ets), mase_f(e_holt))
) |> arrange(RECM)
print(resultado)
pausa()

# ── 3. CONTRASTE DE DIEBOLD-MARIANO POR PARES ─────
comun_se <- !is.na(e_sarima) & !is.na(e_ets)
comun_sh <- !is.na(e_sarima) & !is.na(e_holt)
comun_he <- !is.na(e_holt)   & !is.na(e_ets)

dm_se <- dm.test(e_sarima[comun_se], e_ets[comun_se],  h = 1, power = 2)
dm_sh <- dm.test(e_sarima[comun_sh], e_holt[comun_sh], h = 1, power = 2)
dm_he <- dm.test(e_holt[comun_he],   e_ets[comun_he],  h = 1, power = 2)

cat(sprintf("DM SARIMA vs ETS:  stat=%.3f p=%.4f\n", dm_se$statistic, dm_se$p.value))
cat(sprintf("DM SARIMA vs Holt: stat=%.3f p=%.4f\n", dm_sh$statistic, dm_sh$p.value))
cat(sprintf("DM Holt vs ETS:    stat=%.3f p=%.4f\n", dm_he$statistic, dm_he$p.value))
pausa()

# ── 4. COMBINACION SARIMA-ETS: PESO OPTIMO VS. PESOS IGUALES ─
actual    <- as.numeric(ipc_ts)
fc_sarima <- actual - as.numeric(e_sarima)
fc_ets    <- actual - as.numeric(e_ets)

idx  <- which(comun_se)
mitad <- floor(length(idx) / 2)
idx1 <- idx[1:mitad]; idx2 <- idx[(mitad + 1):length(idx)]

s1 <- sd(e_sarima[idx1]); s2 <- sd(e_ets[idx1]); rho <- cor(e_sarima[idx1], e_ets[idx1])
w_star <- (s2^2 - rho * s1 * s2) / (s1^2 + s2^2 - 2 * rho * s1 * s2)
cat(sprintf("Peso optimo (SARIMA) = %.3f  (sigma1=%.3f sigma2=%.3f rho=%.3f)\n", w_star, s1, s2, rho))

e_opt <- actual[idx2] - (w_star * fc_sarima[idx2] + (1 - w_star) * fc_ets[idx2])
e_eq  <- actual[idx2] - (0.5    * fc_sarima[idx2] + 0.5           * fc_ets[idx2])
cat(sprintf("RMSE fuera de muestra: SARIMA=%.4f ETS=%.4f Combo optimo=%.4f Combo igual=%.4f\n",
            recm_f(e_sarima[idx2]), recm_f(e_ets[idx2]), recm_f(e_opt), recm_f(e_eq)))
cat("Nota: la combinacion con pesos iguales suele ser casi tan buena como la 'optima':\n")
cat("      el error de estimacion de sigma1, sigma2 y rho erosiona la ventaja teorica.\n")

# === FIN Caso Práctico Macro — Tema 11 (IPC, evaluacion y combinacion) ===