# ══════════════════════════════════════════════════
# T14 · Caso práctico Financiero — IBEX 35 (Decision Desk)
# Abre primero mSeriesTemporales.Rproj en RStudio
# ══════════════════════════════════════════════════

library(tidyverse)
library(rugarch)

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

kupiec_test <- function(viol, alpha) {
  x <- sum(viol); Tt <- length(viol); pihat <- x / Tt
  LR <- -2 * ((Tt - x) * log(1 - alpha) + x * log(alpha) -
              (Tt - x) * log(1 - pihat) - x * log(pihat))
  list(x = x, T = Tt, LR = LR, p = 1 - pchisq(LR, df = 1))
}

christoffersen_test <- function(viol) {
  n00 <- n01 <- n10 <- n11 <- 0
  for (t in 2:length(viol)) {
    if (viol[t - 1] == 0 && viol[t] == 0) n00 <- n00 + 1
    if (viol[t - 1] == 0 && viol[t] == 1) n01 <- n01 + 1
    if (viol[t - 1] == 1 && viol[t] == 0) n10 <- n10 + 1
    if (viol[t - 1] == 1 && viol[t] == 1) n11 <- n11 + 1
  }
  pi01 <- n01 / (n00 + n01); pi11 <- n11 / (n10 + n11); pi <- (n01 + n11) / (n00 + n01 + n10 + n11)
  LL0 <- (n00 + n10) * log(1 - pi) + (n01 + n11) * log(pi)
  LL1 <- n00 * log(1 - pi01) + n01 * log(pi01) + n10 * log(1 - pi11) + n11 * log(pi11)
  LR_ind <- -2 * (LL0 - LL1)
  list(LR = LR_ind, p = 1 - pchisq(LR_ind, df = 1))
}

# ── 1. CARGA Y VaR AL 99%: HISTORICO Y NORMAL ─────
source("scripts/deskR.R")
ret_ibex <- desk_returns("_IBEX", type = "log")
r <- ret_ibex$return
alpha <- 0.01
mu_r <- mean(r); sd_r <- sd(r)

VaR_hist <- as.numeric(-quantile(r, alpha))
VaR_norm <- -(mu_r + qnorm(alpha) * sd_r)
cat(sprintf("VaR historico (99%%): %.2f%% | VaR normal (99%%): %.2f%%\n", VaR_hist * 100, VaR_norm * 100))
pausa()

# ── 2. VaR DINAMICO CON EGARCH(1,1)-t (CAP. 9) ────
spec <- ugarchspec(variance.model = list(model = "eGARCH", garchOrder = c(1, 1)),
                    mean.model = list(armaOrder = c(0, 0), include.mean = TRUE),
                    distribution.model = "std")
fit <- ugarchfit(spec, r, solver = "hybrid")
nu_hat <- coef(fit)["shape"]
sigma_t <- as.numeric(sigma(fit)); mu_t <- as.numeric(fitted(fit))
q_alpha <- qdist("std", p = alpha, shape = nu_hat)
VaR_garch <- -(mu_t + sigma_t * q_alpha)
cat(sprintf("VaR GARCH ultimo dia: %.2f%% | VaR GARCH medio muestra: %.2f%%\n",
            tail(VaR_garch, 1) * 100, mean(VaR_garch) * 100))
pausa()

# ── 3. MONTE CARLO: VALIDACION DEL VaR GARCH ──────
set.seed(2026)
Nsim <- 1e5
mu_T <- tail(mu_t, 1); sigma_T <- tail(sigma_t, 1)
r_sim <- mu_T + sigma_T * rdist("std", n = Nsim, shape = nu_hat)
VaR_mc <- -as.numeric(quantile(r_sim, alpha))
ES_mc  <- -mean(r_sim[r_sim <= -VaR_mc])
cat(sprintf("VaR Monte Carlo: %.2f%% (analitico: %.2f%%) | ES Monte Carlo: %.2f%%\n",
            VaR_mc * 100, -(mu_T + sigma_T * q_alpha) * 100, ES_mc * 100))
pausa()

# ── 4. BACKTESTING: KUPIEC Y CHRISTOFFERSEN ───────
viol_hist  <- as.numeric(r < -VaR_hist)
viol_norm  <- as.numeric(r < -VaR_norm)
viol_garch <- as.numeric(r < -VaR_garch)

for (nm in c("hist", "norm", "garch")) {
  viol <- get(paste0("viol_", nm))
  ku <- kupiec_test(viol, alpha); ch <- christoffersen_test(viol)
  cat(sprintf("%-6s violaciones=%.2f%% | Kupiec p=%.4f | Christoffersen p=%.4f\n",
              nm, mean(viol) * 100, ku$p, ch$p))
}
cat("\nSolo el VaR GARCH (dinamico) supera con holgura el contraste de independencia de Christoffersen.\n")

# === FIN Caso Práctico Financiero — Tema 14 (IBEX 35, VaR y Expected Shortfall) ===