Laboratorio: modelos fundacionales frente a GARCH, HAR y LightGBM.

El banco de pruebas

La lección 5.8.4 terminó con una pista: el mismo Chronos-Bolt que no predice el rendimiento explica casi la mitad de la varianza del logaritmo de la varianza realizada. Falta saber si eso es bueno, es decir, si mejora lo que ya hacen los modelos de volatilidad de siempre. Este laboratorio lo mide con cinco familias de modelos:

Los datos son los de vol5 (README de los datos del nivel 5, sección «Volatilidad»), sintéticos. El proceso verdadero es un HAR en el logaritmo de la varianza integrada, con efecto apalancamiento, saltos, cambios de nivel y choques en fechas conocidas. Eso tiene dos consecuencias que hay que tener presentes al leer los resultados: un HAR bien estimado está muy cerca del óptimo, y existe un oráculo, la esperanza exacta de la RV con la información del día anterior, que marca el suelo de lo alcanzable. Un modelo por debajo del oráculo tendría una fuga.

Se usan tres activos con perfiles distintos: el índice (apalancamiento fuerte, γ = 0,18, sin saltos), el bono (sin apenas apalancamiento y con un cambio de nivel de +0,60 en log-varianza el 3 de enero de 2022) y la cripto (12 saltos al año de desviación 7 % y un cambio de nivel de −0,60 en 2019). La prueba cubre 2018-2025.

flowchart TB
  D["vol5: 3 activos"] --> O["Orígenes 2018-2025"]
  O --> E["HAR y GARCH"]
  O --> G["LightGBM"]
  O --> F["Fundacionales"]
  E --> P["QLIKE y ECM"]
  G --> P
  F --> P
  P --> T["DM y MCS"]

El protocolo es de origen móvil sin fugas. Los modelos que se estiman (GARCH, HAR, LightGBM) se reestiman cada 1 de enero con todos los datos anteriores y prevén con esos parámetros cada día del año; sus features en el día $t$ solo usan información hasta $t-1$. Los fundacionales no se estiman: en cada origen reciben las 512 sesiones anteriores. Como cada previsión de TimesFM cuesta del orden de una décima de segundo en CPU, los modelos se comparan en una de cada 20 sesiones (101 orígenes por activo); al final se comprueba con todas las sesiones lo que dicen los modelos baratos.

import numpy as np
import pandas as pd

d = pd.read_parquet("vol5_diario.parquet")
v = pd.read_parquet("vol5_verdad.parquet")
ACT = ["indice", "bono", "cripto"]
# parte nocturna ω, saltos por año y desviación del salto (README)
PAR = {"indice": (0.20, 0, 0.0), "bono": (0.15, 0, 0.0),
       "cripto": (0.0, 12, 0.07)}
fechas = pd.DatetimeIndex(d.loc[d["activo"] == "indice", "fecha"])
n = len(fechas)
corte = fechas.searchsorted(pd.Timestamp("2018-01-02"))
orig = np.arange(corte, n, 20)              # una de cada 20 sesiones
inicio = {y: fechas.searchsorted(pd.Timestamp(f"{y}-01-01"))
          for y in range(2018, 2027)}
r = {a: d.loc[d["activo"] == a, "rendimiento"].to_numpy() for a in ACT}
rv = {a: d.loc[d["activo"] == a, "rv_5min"].to_numpy() for a in ACT}
lrv = {a: np.log(rv[a]) for a in ACT}
print(f"{len(orig)} orígenes por activo de {fechas[orig[0]].date()} "
      f"a {fechas[orig[-1]].date()}; {n - corte} sesiones en 2018-2025")

101 orígenes por activo de 2018-01-02 a 2025-12-16; 2011 sesiones en 2018-2025

Pérdidas, oráculo y HAR

La varianza verdadera no se observa; se compara con la RV de 5 minutos, que la mide con error. Patton (2011) demostró que no todas las pérdidas sirven con un proxy ruidoso: el ECM y la QLIKE sí ordenan los modelos igual que si se conociera la varianza verdadera. La QLIKE, con $h$ la previsión,

vale 0 cuando $h = RV$, castiga más quedarse corto que pasarse y depende solo del cociente $RV/h$, así que un día de varianza extrema no la domina como al ECM. El HAR se estima en logaritmos y su previsión de la media lleva la corrección lognormal $e^{\hat x + \hat\sigma^2/2}$, con $\hat\sigma^2$ la varianza de los residuos del entrenamiento. El oráculo es la fórmula del README: $E[RV_t \mid t-1] = (1-\omega)\,E[IV_t \mid t-1] + \lambda\,\sigma_J^2$, con la parte nocturna $\omega$, la intensidad de saltos $\lambda$ (por sesión) y su varianza $\sigma_J^2$.

flowchart LR
  A["1 de enero"] --> B["Reestimar con < t"]
  B --> C["Prever cada día"]
  C --> D["Año siguiente"]
  D --> A

def qlike(obs, f):
    """QLIKE de Patton: obs/f - log(obs/f) - 1 (0 si f = obs)."""
    return obs / f - np.log(obs / f) - 1

def medias_har(x):
    s = pd.Series(x)
    return pd.DataFrame({"d1": s.shift(1),
                         "s5": s.rolling(5).mean().shift(1),
                         "m22": s.rolling(22).mean().shift(1)})

prev, coef = {a: {} for a in ACT}, {}
for a in ACT:
    om, nj, sj = PAR[a]
    iv_prev = v.loc[v["activo"] == a, "iv_prevista"].to_numpy()
    prev[a]["oraculo"] = iv_prev * (1 - om) + nj / 252 * sj ** 2
    prev[a]["ingenua"] = np.r_[np.nan, rv[a][:-1]]
    X = medias_har(lrv[a]).to_numpy()
    har = np.full(n, np.nan)
    for y in range(2018, 2026):             # reestimado cada enero
        t0, t1 = inicio[y], inicio[y + 1]
        ok = np.isfinite(X[:t0]).all(axis=1)
        A = np.column_stack([np.ones(ok.sum()), X[:t0][ok]])
        b, *_ = np.linalg.lstsq(A, lrv[a][:t0][ok], rcond=None)
        s2 = np.var(lrv[a][:t0][ok] - A @ b)
        har[t0:t1] = np.exp(b[0] + X[t0:t1] @ b[1:] + s2 / 2)
    prev[a]["har"], coef[a] = har, b
pd.DataFrame(coef, index=["c", "d1", "s5", "m22"]).round(3)

     indice   bono  cripto
c    -0.808 -1.215  -0.615
d1    0.373  0.273   0.137
s5    0.307  0.385   0.274
m22   0.238  0.236   0.503

Los coeficientes del último ajuste (enero de 2025) muestran la estructura HAR: el peso se reparte entre la RV de ayer, la media semanal y la mensual, con una suma de 0,92 en el índice, 0,89 en el bono y 0,91 en la cripto, cerca de la persistencia verdadera (0,92). En la cripto pesa más la media mensual (0,503): los saltos meten en la RV diaria un ruido que la media larga filtra.

GARCH con arch

El GARCH(1,1) se estima con los rendimientos diarios, no con la RV. arch_model(...).fit(last_obs=t0) estima con los datos anteriores a t0 y forecast(start=t0 - 1, ...) filtra el resto del año con esos parámetros fijos, de modo que la previsión para el día $t$ solo usa rendimientos hasta $t-1$ (lección 2.1.4). La varianza del GARCH es de cierre a cierre; la RV no ve la parte nocturna, así que se reescala con el cociente $k$ entre la RV media y la varianza de los rendimientos del entrenamiento, como en el README.

from arch import arch_model

for a in ACT:
    gar = np.full(n, np.nan)
    for y in range(2018, 2026):
        t0, t1 = inicio[y], inicio[y + 1]
        am = arch_model(100 * r[a][:t1], mean="Constant", vol="GARCH",
                        p=1, q=1, dist="normal")
        res = am.fit(last_obs=t0, disp="off")
        # la fila de t usa datos hasta t y prevé t + 1
        h = res.forecast(start=t0 - 1, horizon=1, reindex=False)
        h = h.variance.to_numpy()[: t1 - t0, 0] / 1e4
        k = rv[a][:t0].mean() / np.var(res.resid[:t0] / 100)
        gar[t0:t1] = k * h        # a la escala de la RV (sin la noche)
    prev[a]["garch"] = gar
print(res.params.round(4).to_dict(), f"k = {k:.2f}")

{'mu': 0.1113, 'omega': 0.1264, 'alpha[1]': 0.0183, 'beta[1]': 0.9706} k = 1.04

La salida es la del último ajuste, la cripto en enero de 2025: α = 0,018 y β = 0,971, una persistencia de 0,989 típica de un GARCH diario, y $k$ = 1,04 (en la cripto no hay parte nocturna: ω = 0). El GARCH reacciona solo a los rendimientos al cuadrado, una medida de la varianza mucho más ruidosa que la RV; por eso, en datos con RV, el HAR suele ganarle con holgura.

LightGBM con retardos

El modelo de árboles recibe más información que el HAR: cinco retardos del log de la RV, las medias de 5, 22 y 66 sesiones, el rendimiento de ayer y la parte negativa de los rendimientos del último día y de la última semana, que es lo que necesita para captar el efecto apalancamiento (la varianza sube tras las caídas). Configuración pequeña y determinista (lección 5.5.3), reestimada cada enero; la previsión de la media lleva la misma corrección lognormal que el HAR.

import lightgbm as lgb

def features(a):
    s, ret = pd.Series(lrv[a]), pd.Series(r[a])
    F = medias_har(lrv[a])
    for k in range(2, 6):
        F[f"d{k}"] = s.shift(k)
    F["m66"] = s.rolling(66).mean().shift(1)
    F["r1"] = ret.shift(1)
    F["neg1"] = np.minimum(ret, 0).shift(1)
    F["neg5"] = np.minimum(ret, 0).rolling(5).sum().shift(1)
    return F

def lgbm():
    return lgb.LGBMRegressor(
        n_estimators=300, learning_rate=0.03, num_leaves=8,
        min_child_samples=100, subsample=0.8, subsample_freq=1,
        colsample_bytree=0.8, random_state=0, n_jobs=4,
        deterministic=True, force_row_wise=True, verbose=-1)

ult = {}
for a in ACT:
    F = features(a)
    ok = F.notna().all(axis=1).to_numpy()
    p = np.full(n, np.nan)
    for y in range(2018, 2026):
        t0, t1 = inicio[y], inicio[y + 1]
        tr = ok & (np.arange(n) < t0)
        m = lgbm().fit(F[tr], lrv[a][tr])
        s2 = np.var(lrv[a][tr] - m.predict(F[tr]))
        p[t0:t1] = np.exp(m.predict(F[t0:t1]) + s2 / 2)
    prev[a]["lgbm"], ult[a] = p, m
imp = pd.Series(ult["indice"].booster_.feature_importance("gain"),
                index=F.columns)
(imp / imp.sum()).sort_values(ascending=False).head(5).round(3)

s5      0.437
d1      0.334
m22     0.096
r1      0.042
neg1    0.032
dtype: float64

La importancia por ganancia del modelo del índice en 2025 la dominan las mismas variables que el HAR: la media semanal (43,7 %), la RV de ayer (33,4 %) y la media mensual (9,6 %). El rendimiento de ayer (4,2 %) y su parte negativa (3,2 %) pesan poco, pero son información que un HAR-RV no tiene; si LightGBM gana en algún activo, será por ahí.

Los modelos fundacionales sobre el log de la RV

Los dos modelos reciben el logaritmo de la RV, no la RV: en logaritmos la serie es casi gaussiana y su escala no cambia con el nivel (el ejercicio 1 prueba la RV en niveles). Los dos devuelven cuantiles; la mediana del log, exponenciada, es la mediana de la RV, que queda por debajo de su media porque la distribución de la RV tiene cola a la derecha. Para prever la media, que es lo que premian la QLIKE y el ECM, se usa la corrección lognormal con la dispersión que da el propio modelo: si el log es aproximadamente normal, $\sigma \approx (q_{0,9} - q_{0,1}) / (2 \cdot 1{,}2816)$ y la media es $e^{q_{0,5} + \sigma^2/2}$. Las 303 series (101 orígenes × 3 activos) van en un solo lote.

import torch

torch.manual_seed(0)
torch.use_deterministic_algorithms(True)
from chronos import BaseChronosPipeline
import timesfm

CTX = 512
ctx = [lrv[a][t - CTX:t].astype("float32") for a in ACT for t in orig]
Z = 2 * 1.2815515655446004            # distancia entre q10 y q90 en σ

chronos = BaseChronosPipeline.from_pretrained(
    "amazon/chronos-bolt-small", device_map="cpu", dtype=torch.float32)
q, _ = chronos.predict_quantiles(
    inputs=[torch.tensor(c) for c in ctx], prediction_length=1,
    quantile_levels=[0.1, 0.5, 0.9])
q = q[:, 0, :].numpy().astype(float)
ch_med, ch_s = q[:, 1], (q[:, 2] - q[:, 0]) / Z

tfm = timesfm.TimesFM_2p5_200M_torch.from_pretrained(
    "google/timesfm-2.5-200m-pytorch", torch_compile=False)
tfm.compile(timesfm.ForecastConfig(
    max_context=CTX, max_horizon=16, normalize_inputs=True,
    use_continuous_quantile_head=True, per_core_batch_size=len(ctx)))
punto, cuant = tfm.forecast(horizon=1, inputs=ctx)
tf_med = punto[:, 0].astype(float)
tf_s = (cuant[:, 0, 9] - cuant[:, 0, 1]).astype(float) / Z

k = len(orig)
for i, a in enumerate(ACT):
    sl = slice(i * k, (i + 1) * k)
    for nom, med, s in [("chronos", ch_med, ch_s),
                        ("timesfm", tf_med, tf_s)]:
        prev[a][nom + "_mediana"] = np.exp(med[sl])
        prev[a][nom] = np.exp(med[sl] + s[sl] ** 2 / 2)
print("σ log prevista (media):", round(ch_s.mean(), 3),
      round(tf_s.mean(), 3))

σ log prevista (media): 0.409 0.42

La dispersión media prevista del log de la RV es 0,41 en Chronos y 0,42 en TimesFM: la corrección multiplica la mediana por $e^{0{,}41^2/2} \approx 1{,}09$, un 9 % más de varianza prevista.

Resultados: QLIKE y ECM

MOD = ["ingenua", "garch", "har", "lgbm", "chronos_mediana",
       "chronos", "timesfm_mediana", "timesfm", "oraculo"]

def en_orig(a, nom):
    f = prev[a][nom]
    return f if len(f) == len(orig) else f[orig]

perd = {a: {m: qlike(rv[a][orig], en_orig(a, m)) for m in MOD}
        for a in ACT}
tq = pd.DataFrame({a: {m: perd[a][m].mean() for m in MOD} for a in ACT})
ecm = pd.DataFrame({a: {m: np.mean((rv[a][orig] - en_orig(a, m)) ** 2)
                        for m in MOD} for a in ACT})
te = ecm / ecm.loc["oraculo"]
pd.concat({"QLIKE": tq, "ECM / ECM oráculo": te}, axis=1).round(3)

                 QLIKE               ECM / ECM oráculo
                indice   bono cripto            indice   bono cripto
ingenua          0.105  0.082  0.662             1.644  2.160  3.004
garch            0.135  0.119  0.362             1.790  1.843  1.013
har              0.080  0.056  0.383             1.191  0.995  1.029
lgbm             0.074  0.054  0.398             1.155  0.917  1.028
chronos_mediana  0.085  0.060  0.461             1.375  1.071  1.049
chronos          0.082  0.059  0.408             1.275  1.139  1.045
timesfm_mediana  0.089  0.060  0.465             1.324  1.025  1.049
timesfm          0.085  0.059  0.414             1.234  1.030  1.037
oraculo          0.069  0.054  0.350             1.000  1.000  1.000

La figura de la cabecera resume la tabla dividiendo cada QLIKE por la del HAR del mismo activo. Lo que dice:

• El oráculo es el suelo (0,069, 0,054 y 0,350) y nadie lo cruza en QLIKE: no hay fugas.
• El HAR queda cerca del suelo en el índice y el bono, como anuncia el README: es casi el modelo verdadero.
• LightGBM mejora al HAR en el índice (0,074 frente a 0,080) y en el bono (0,054 frente a 0,056), y queda detrás en la cripto (0,398 frente a 0,383).
• Los dos fundacionales, sin entrenar nada, quedan entre la ingenua y el HAR: Chronos 0,082, 0,059 y 0,408; TimesFM 0,085, 0,059 y 0,414. Están un 2-8 % por encima del HAR y muy por debajo de la ingenua (0,105, 0,082 y 0,662). Chronos, con 48 millones de parámetros, queda ligeramente por delante de TimesFM, con 231 millones.
• La corrección de la media importa, sobre todo con saltos: en la cripto, la QLIKE de Chronos baja de 0,461 con la mediana a 0,408 con la media; la de TimesFM, de 0,465 a 0,414.
• El GARCH es el peor en el índice y el bono (0,135 y 0,119) y el mejor no oráculo en la cripto (0,362). La sección «Una muestra corta engaña» vuelve sobre este resultado.

El ECM relativo al oráculo enseña por qué no se usa como pérdida principal: en el bono, LightGBM (0,917) y el HAR (0,995) salen por debajo del oráculo, algo imposible en esperanza. En 101 orígenes el ECM lo deciden unos pocos días de varianza extrema; la QLIKE, que mira el cociente, es mucho más estable.

La figura muestra el bono entre 2021 y 2023. El 3 de enero de 2022 su volatilidad sube de forma permanente (+0,60 en log-varianza, ≈ 35 % más de volatilidad). Ni el HAR, reestimado con todo el pasado, ni TimesFM, que normaliza cada contexto de 512 sesiones, saben que el cambio ha ocurrido; los dos lo absorben a través de la RV reciente y siguen casi el mismo camino. El modelo fundacional no tiene ventaja de adaptación sobre un HAR, que también se alimenta de las medias de 1, 5 y 22 días.

¿Son significativas las diferencias?

Con 101 orígenes por activo, cada diferencia tiene su error. Primero, Diebold-Mariano de cada modelo frente al HAR con la QLIKE, con varianza de Newey-West (4 retardos, aunque con orígenes separados 20 sesiones la autocorrelación de las pérdidas es pequeña) y contraste bilateral. Un $t$ positivo es que el modelo pierde más que el HAR.

from scipy import stats

def dm_hac(l1, l2, rezagos=4):
    """t de Diebold-Mariano con varianza de Newey-West (bilateral)."""
    dd = l1 - l2
    u = dd - dd.mean()
    nn = len(dd)
    var = u @ u / nn
    for j in range(1, rezagos + 1):
        var += 2 * (1 - j / (rezagos + 1)) * (u[j:] @ u[:-j]) / nn
    t = dd.mean() / np.sqrt(var / nn)
    return t, 2 * stats.t.sf(abs(t), nn - 1)

COMP = ["ingenua", "garch", "lgbm", "chronos", "timesfm", "oraculo"]
pd.DataFrame({a: {m: "{:+.2f} ({:.2f})".format(
                     *dm_hac(perd[a][m], perd[a]["har"]))
                  for m in COMP} for a in ACT})

               indice          bono        cripto
ingenua  +2.33 (0.02)  +2.73 (0.01)  +1.93 (0.06)
garch    +3.74 (0.00)  +5.56 (0.00)  -0.19 (0.85)
lgbm     -0.97 (0.33)  -1.04 (0.30)  +0.83 (0.41)
chronos  +0.59 (0.56)  +1.33 (0.19)  +0.67 (0.50)
timesfm  +1.66 (0.10)  +1.14 (0.26)  +0.76 (0.45)
oraculo  -1.94 (0.06)  -1.01 (0.31)  -0.88 (0.38)

Los fundacionales pierden frente al HAR en los tres activos, pero ninguna diferencia es significativa: $t$ entre +0,59 y +1,66 (TimesFM en el índice, p = 0,10). LightGBM mejora al HAR en el índice y el bono sin significación ($t$ = −0,97 y −1,04). Solo el GARCH en el índice y el bono ($t$ = +3,74 y +5,56) y la ingenua en el índice y el bono son claramente peores. La fila del oráculo es la más reveladora: ni siquiera la previsión perfecta bate al HAR con significación en esta muestra ($t$ = −1,94 en el índice, p = 0,06). Con 101 puntos, lo que se puede afirmar es poco.

Diebold-Mariano compara parejas. Para preguntar qué modelos podrían ser el mejor a la vez, sin elegir uno de referencia, está el Model Confidence Set de Hansen, Lunde y Nason (2011): parte de todos los modelos, contrasta si todos tienen la misma pérdida esperada y, si se rechaza, elimina el peor y repite. Los que quedan forman un conjunto que contiene al mejor modelo con probabilidad $1-\alpha$. La versión de la celda usa el estadístico $T_{max}$: para cada modelo, su pérdida media menos la media de los que siguen vivos, dividida por su error típico; la distribución bajo la nula sale de un bootstrap de bloques móviles (bloques de 5, para respetar la dependencia). La librería arch trae la misma prueba (arch.bootstrap.MCS), que sirve para contrastar la implementación propia.

from arch.bootstrap import MCS

def mcs(L, alfa=0.10, B=2000, bloque=5, semilla=0):
    """MCS de Hansen, Lunde y Nason con el estadístico T_max y
    bootstrap de bloques móviles. Devuelve p del MCS por modelo."""
    rng = np.random.default_rng(semilla)
    nn, _ = L.shape
    nb = -(-nn // bloque)
    ini = rng.integers(0, nn - bloque + 1, size=(B, nb))
    idx = (ini[:, :, None] + np.arange(bloque)).reshape(B, -1)[:, :nn]
    vivos, p_mcs, pmax = list(L.columns), {}, 0.0
    while len(vivos) > 1:
        X = L[vivos].to_numpy()
        di = X.mean(0) - X.mean()          # pérdida frente a la media
        Xb = X[idx].mean(1)
        dib = Xb - Xb.mean(1, keepdims=True)
        se = np.sqrt(((dib - di) ** 2).mean(0))
        t = di / se
        Tb = ((dib - di) / se).max(1)       # T_max bajo la nula
        pmax = max(pmax, float(np.mean(Tb >= t.max())))
        if pmax >= alfa:
            break
        peor = vivos[int(np.argmax(t))]
        p_mcs[peor] = pmax
        vivos.remove(peor)
    return {**p_mcs, **{m: 1.0 for m in vivos}}

CAND = ["ingenua", "garch", "har", "lgbm", "chronos", "timesfm"]
for a in ACT:
    L = pd.DataFrame({m: perd[a][m] for m in CAND})
    propio = mcs(L)
    ref = MCS(L, size=0.10, reps=2000, block_size=5, method="max",
              bootstrap="mbb", seed=0)
    ref.compute()
    dentro = sorted(m for m in CAND if propio[m] >= 0.1)
    print(f"{a:7s} propio: {dentro}\n"
          f"        arch:   {sorted(ref.included)}")

indice  propio: ['chronos', 'har', 'lgbm']
        arch:   ['chronos', 'har', 'lgbm']
bono    propio: ['chronos', 'har', 'lgbm', 'timesfm']
        arch:   ['chronos', 'har', 'lgbm', 'timesfm']
cripto  propio: ['chronos', 'garch', 'har', 'ingenua', 'lgbm', 'timesfm']
        arch:   ['chronos', 'garch', 'har', 'ingenua', 'lgbm', 'timesfm']

La implementación propia y la de arch dan los mismos conjuntos al 10 %. En el índice quedan HAR, LightGBM y Chronos; en el bono, esos tres y TimesFM; en la cripto no se puede descartar ninguno, ni siquiera la ingenua, cuya QLIKE es 1,7 veces la del HAR: con saltos de un 7 % de desviación, unos pocos días dominan las pérdidas y 101 orígenes no bastan para separar nada. El MCS no da un ganador; da los modelos que los datos no permiten descartar.

Una muestra corta engaña

Los modelos baratos se pueden evaluar en las 2011 sesiones de 2018-2025 en vez de en 101. La celda lo hace y repite Diebold-Mariano de LightGBM frente al HAR, ahora con 10 retardos de Newey-West porque las sesiones son consecutivas.

BAR = ["ingenua", "garch", "har", "lgbm", "oraculo"]
todo = pd.DataFrame({a: {m: qlike(rv[a][corte:], prev[a][m][corte:])
                         .mean() for m in BAR} for a in ACT})
t_lg = {a: float(round(dm_hac(
            qlike(rv[a][corte:], prev[a]["lgbm"][corte:]),
            qlike(rv[a][corte:], prev[a]["har"][corte:]),
            rezagos=10)[0], 2)) for a in ACT}
print(f"t de LightGBM frente al HAR, {n - corte} sesiones:", t_lg)
todo.round(3)

t de LightGBM frente al HAR, 2011 sesiones: {'indice': -2.41, 'bono': 3.07, 'cripto': 1.33}
         indice   bono  cripto
ingenua   0.126  0.101   0.525
garch     0.149  0.113   0.385
har       0.089  0.067   0.366
lgbm      0.086  0.069   0.377
oraculo   0.080  0.065   0.325

Dos órdenes de la submuestra se invierten. En el bono, LightGBM pasa de 0,96 veces la QLIKE del HAR a 1,03, y con 2011 sesiones la diferencia es significativa en contra ($t$ = +3,07). En la cripto, el GARCH pasa de 0,95 a 1,05. En el índice, en cambio, LightGBM sigue por delante (0,086 frente a 0,089) y ahora con significación ($t$ = −2,41). La diferencia entre los dos activos tiene una explicación plausible: el índice tiene un efecto apalancamiento fuerte (γ = 0,18) que LightGBM capta con la parte negativa de los rendimientos y el HAR-RV no ve; en el bono (γ = 0,04) esa información no existe y la flexibilidad de los árboles solo añade varianza. No es una prueba: las dos features de rendimientos pesan solo un 7 % de la ganancia del modelo. Las QLIKE de la ingenua, el HAR y el oráculo coinciden con las de la tabla del README (0,126, 0,089 y 0,080 en el índice); la del GARCH difiere un poco porque aquí se reestima cada año con arch.

El coste de cada modelo

tam = {
    "HAR": 4, "GARCH(1,1)": 4,
    "LightGBM (hojas)": int(sum(
        t["num_leaves"] for t in
        ult["indice"].booster_.dump_model()["tree_info"])),
    "Chronos-Bolt small": sum(x.numel() for x in
                              chronos.model.parameters()),
    "TimesFM-2.5": sum(x.numel() for x in tfm.model.parameters()),
}
pd.Series(tam).map("{:,}".format)

HAR                             4
GARCH(1,1)                      4
LightGBM (hojas)            2,400
Chronos-Bolt small     47,718,016
TimesFM-2.5           231,289,280
dtype: str

El HAR y el GARCH tienen cuatro parámetros cada uno; LightGBM, 2 400 hojas (300 árboles de 8); Chronos-Bolt small, 47,7 millones de parámetros, y TimesFM-2.5, 231,3 millones. La figura pone el tamaño frente a la QLIKE media relativa al HAR: cincuenta millones de veces más parámetros no compran precisión en una serie cuyo proceso es un HAR.

El tiempo depende de la máquina y de su carga, así que el curso no lo ejecuta; este es el código para medirlo:

import time

t0 = time.perf_counter()
chronos.predict_quantiles(inputs=[torch.tensor(c) for c in ctx],
                          prediction_length=1, quantile_levels=[0.5])
print(f"Chronos, {len(ctx)} previsiones: "
      f"{time.perf_counter() - t0:.1f} s")
t0 = time.perf_counter()
tfm.forecast(horizon=1, inputs=ctx)
print(f"TimesFM, {len(ctx)} previsiones: "
      f"{time.perf_counter() - t0:.1f} s")

En la máquina donde se escribió la lección (16 núcleos compartidos con otros procesos, carga ≈ 12, 4 hilos por proceso), las 24 reestimaciones del HAR (3 activos × 8 años) tardaron 0,02 s; las del GARCH con sus previsiones, 0,9 s; las de LightGBM, 1,9 s; las 303 previsiones de Chronos-Bolt, 3,1 s, y las de TimesFM-2.5, 35 s (62 s con la máquina más cargada). La diferencia de fondo no es el entrenamiento sino la inferencia: el HAR prevé con un producto escalar, y TimesFM necesita una pasada de 231 millones de parámetros por cada previsión, unos 0,1-0,2 s en esta CPU. Prever las 2011 sesiones de los tres activos serían unas 6 000 pasadas, del orden de 10-20 minutos aquí; en GPU (device="cuda") bajaría mucho, pero sigue sin competir con cuatro multiplicaciones.

Qué hacer con esto

Es la misma conclusión que el estudio de Brini (2026) con 50 activos reales, citado en la lección 5.8.4: los modelos fundacionales no dan una mejora uniforme sobre un Log-HAR bien especificado, y el tamaño o la arquitectura importan menos que tener una baseline econométrica seria.

Errores típicos

• Comparar con baselines débiles. Frente a la ingenua, los fundacionales parecen excelentes (del 22 al 38 % menos de QLIKE); frente al HAR, pierden.
• Usar la mediana como previsión de la varianza. Con QLIKE o ECM se prevé la media; la mediana de una distribución con cola a la derecha se queda corta, y más con saltos.
• Fiarse del ECM con pocas observaciones. En el bono, dos modelos «baten al oráculo» en ECM.
• Sacar conclusiones de 100 puntos. Ni el oráculo es significativamente mejor que el HAR; el orden de LightGBM y el GARCH se invierte con la muestra completa.
• Fugas en la evaluación. Features con el día $t$, parámetros estimados con el año de prueba, normalizaciones con toda la serie (lección 5.2.4): el oráculo marca un suelo que ningún modelo honesto cruza.
• Ignorar la inferencia. Un modelo de 231 millones de parámetros cuesta en cada previsión, no solo al entrenar.

Ejercicios

1. Repite la previsión de Chronos-Bolt pasándole la RV en niveles en lugar de su logaritmo (y su mediana como previsión). ¿Qué representación le conviene más?

ctx_n = [rv[a][t - CTX:t].astype("float32") for a in ACT for t in orig]
qn, _ = chronos.predict_quantiles(
    inputs=[torch.tensor(c) for c in ctx_n], prediction_length=1,
    quantile_levels=[0.5])
qn = qn[:, 0, 0].numpy().astype(float)
pd.DataFrame({a: {"niveles": qlike(rv[a][orig],
                                   qn[i * k:(i + 1) * k]).mean(),
                  "log (mediana)": perd[a]["chronos_mediana"].mean(),
                  "log (media)": perd[a]["chronos"].mean()}
              for i, a in enumerate(ACT)}).round(3)

               indice   bono  cripto
niveles         0.096  0.061   0.499
log (mediana)   0.085  0.060   0.461
log (media)     0.082  0.059   0.408

2. La recalibración de Mincer-Zarnowitz regresa la variable sobre la previsión: si la previsión es insesgada, la pendiente es 1 y la constante 0. Calcula, para cada activo, la pendiente de la regresión del log de la RV sobre el log de la previsión del HAR, LightGBM, Chronos y TimesFM en los orígenes. ¿Qué indica una pendiente menor que 1? ¿Se podría usar esa regresión para corregir las previsiones?

mz = {}
for a in ACT:
    yv = np.log(rv[a][orig])
    mz[a] = {m: np.polyfit(np.log(en_orig(a, m)), yv, 1)[0]
             for m in ["har", "lgbm", "chronos", "timesfm"]}
pd.DataFrame(mz).round(2)

         indice  bono  cripto
har        1.07  1.06    1.06
lgbm       1.11  1.11    1.06
chronos    1.01  0.96    0.92
timesfm    1.05  1.01    1.00

3. Sin código: tu equipo usa un HAR para dimensionar posiciones con volatility targeting (lección 4.2.2) en 40 futuros. Alguien propone sustituirlo por TimesFM-2.5 «porque es el estado del arte». Con lo que has medido, ¿qué le responderías y qué prueba pedirías antes de decidir?