Implementation shortfall, VWAP, TWAP y Almgren-Chriss.

La cartera de papel y la real

En 1988, André Perold propuso medir la ejecución comparando dos carteras: la de papel, que compra o vende todo al precio del momento de la decisión y sin costes, y la real, que paga lo que paga y quizá no completa la orden. La diferencia es el implementation shortfall (IS). Sus partes cuentan una historia que el precio medio de ejecución no cuenta:

flowchart TB
  D[Decisión] -- retraso --> L[Llegada]
  L -- ejecución --> E[Precio medio]
  L -- oportunidad --> C[Cierre: no hecho]

Para una orden de $Q$ acciones con lado $s$ (+1 compra, −1 venta), precio de decisión $P_d$, precio de llegada al mercado $P_0$, ejecución de $q$ acciones a un precio medio $\overline{P}$ y cierre $P_c$, en fracción del valor de decisión:

Los tres sumandos son, por orden, el retraso, la ejecución y la oportunidad. El retraso es lo que se mueve el precio entre la decisión y el envío (responsabilidad de quien decide o del proceso); la ejecución, lo que cuesta ejecutar; la oportunidad, lo que se deja de ganar con lo que no se ejecutó. La lección 1.1.5 introdujo la idea; aquí se mide.

Los datos son sintéticos: cuatro acciones ficticias (ARVA, muy líquida, a DUNE, poco líquida y volátil) en 60 sesiones de la NYSE de 2025, con 727 órdenes padre de una mesa (TWAP, VWAP, POV e IS), sus 92.409 hijas y las barras de 1 minuto del mercado, que ya llevan el impacto de esas órdenes. El impacto tiene parámetros conocidos (ejecucion4_acciones.csv): temporal $\eta\,\sigma_d\,\rho^{0{,}5}$ con $\eta = 0{,}40$ y $\rho$ la participación de la orden en el volumen del minuto, y permanente lineal $\gamma\,\sigma_d\,q/\text{ADV}$ con $\gamma = 1{,}0$.

import numpy as np
import pandas as pd

acc = pd.read_csv("ejecucion4_acciones.csv", index_col="ticker")
ord_ = pd.read_csv("ejecucion4_ordenes.csv",
                   parse_dates=["inicio", "fin_real"])
ord_["lado_s"] = np.where(ord_["lado"] == "compra", 1, -1)
ord_["sigma"] = ord_["ticker"].map(acc["sigma_diaria"])
ord_["adv"] = ord_["ticker"].map(acc["adv_acciones"])
ord_["minutos"] = (ord_["fin_real"] - ord_["inicio"]
                   ).dt.total_seconds() / 60
print(acc[["adv_acciones", "sigma_diaria"]])
ord_.groupby("algoritmo").agg(
    ordenes=("id_orden", "size"),
    pct_adv_mediana=("pct_adv", "median"),
    minutos_mediana=("minutos", "median"),
    completas=("cantidad_ejecutada",
               lambda q: (q == ord_.loc[q.index, "cantidad"]).mean()),
).round(3)

        adv_acciones  sigma_diaria
ticker
ARVA        11000000         0.014
BELO         6500000         0.018
CIRA         1750000         0.023
DUNE          640000         0.030
           ordenes  pct_adv_mediana  minutos_mediana  completas
algoritmo
IS             154            0.010             75.5      1.000
POV            198            0.012             39.5      0.808
TWAP           171            0.012            106.0      1.000
VWAP           204            0.016             92.5      1.000

Las órdenes van del 0,1 % al 20 % del volumen medio diario (ADV), con medianas del 1 % al 1,6 % según el algoritmo, y duran de media hora a una sesión entera (medianas de 40 a 106 minutos). Las de TWAP, VWAP e IS se completan siempre; las de POV toman una participación fija del volumen y el 19 % no termina antes del cierre. Como precio de decisión se usa el cierre de la sesión anterior (las órdenes del primer día no lo tienen y se dejan fuera):

def pb(x, ref):
    return 1e4 * x / ref


o = ord_.dropna(subset=["cierre_anterior"]).copy()
s, Q, q = o["lado_s"], o["cantidad"], o["cantidad_ejecutada"]
dec, lleg = o["cierre_anterior"], o["precio_llegada"]
o["retraso"] = pb(s * (lleg - dec), dec)
o["ejecucion"] = pb(s * (o["precio_medio"] - lleg), dec) * q / Q
o["oportunidad"] = pb(s * (o["cierre"] - lleg), dec) * (Q - q) / Q
o["is_total"] = o[["retraso", "ejecucion", "oportunidad"]].sum(axis=1)
print(len(o), "órdenes con cierre anterior")
o.groupby("algoritmo")[["retraso", "ejecucion", "oportunidad",
                        "is_total"]].mean().round(1)

710 órdenes con cierre anterior
           retraso  ejecucion  oportunidad  is_total
algoritmo
IS             1.8       40.6          0.0      42.5
POV           17.2       33.1          2.4      52.7
TWAP          -8.0       47.2          0.0      39.2
VWAP          -2.5       44.8          0.0      42.3

El IS medio es de 39 a 53 pb según el algoritmo, y casi todo es ejecución (33 a 47 pb). El retraso es ruido del mercado entre el cierre y el envío: +17,2 pb en POV, −8,0 en TWAP, sin relación con el algoritmo. La oportunidad solo aparece en POV (2,4 pb), el único que deja órdenes sin completar. Es un recordatorio útil: un algoritmo que completa menos parece más barato si solo se mira lo ejecutado.

Cuatro precios de referencia

La parte de ejecución se puede medir frente a varias referencias, y cada una cuenta una cosa distinta:

barras = pd.read_parquet("ejecucion4_barras_1min.parquet")
por_ticker = {t: g.set_index("ts") for t, g in barras.groupby("ticker")}


def twap_intervalo(r):
    b = por_ticker[r.ticker].loc[r.inicio:r.fin_real
                                 - pd.Timedelta(seconds=1)]
    return b["vwap"].mean()


ord_["twap_intervalo"] = ord_.apply(twap_intervalo, axis=1)
ref = {"llegada": "precio_llegada", "vwap": "vwap_intervalo",
       "twap": "twap_intervalo", "cierre": "cierre"}
for nombre, col in ref.items():
    ord_[f"frente_{nombre}"] = pb(ord_["lado_s"] * (
        ord_["precio_medio"] - ord_[col]), ord_[col])
cols = [f"frente_{n}" for n in ref]
print(ord_.groupby("algoritmo")[cols].mean().round(1))
ord_[cols].std().round(1)

           frente_llegada  frente_vwap  frente_twap  frente_cierre
algoritmo
IS                   40.2         29.1         29.5            9.5
POV                  37.5         29.8         30.1           39.4
TWAP                 46.6         26.9         26.0           22.5
VWAP                 44.2         27.9         26.7           24.6
frente_llegada     72.6
frente_vwap        17.0
frente_twap        18.1
frente_cierre     140.3
dtype: float64

Frente a la llegada, las órdenes cuestan de 37,5 a 46,6 pb; frente al VWAP o al TWAP del intervalo, de 26,0 a 30,1: unos 10 a 20 pb menos. Las desviaciones típicas explican por qué los gestores prefieren el VWAP: 17,0 pb frente a 72,6 de la llegada, porque el VWAP se mueve con el mercado durante la orden y elimina casi todo su ruido. El cierre es la referencia más ruidosa (140,3 pb): no dice nada de la calidad de una ejecución que terminó horas antes.

Lo que el VWAP elimina no es solo ruido:

def volumen_intervalo(r):
    b = por_ticker[r.ticker].loc[r.inicio:r.fin_real
                                 - pd.Timedelta(seconds=1)]
    return b["volumen"].sum()


ord_["participacion_media"] = (ord_["cantidad_ejecutada"]
                               / ord_.apply(volumen_intervalo, axis=1))
tramos = pd.cut(ord_["participacion_media"],
                [0, 0.02, 0.05, 0.10, 0.20, 1.0])
ord_.groupby(tramos, observed=True)[
    ["frente_llegada", "frente_vwap"]].agg(["mean", "size"]).round(1)

                    frente_llegada      frente_vwap
                              mean size        mean size
participacion_media
(0.0, 0.02]                   25.0  160        20.8  160
(0.02, 0.05]                  38.6  126        23.2  126
(0.05, 0.1]                   29.2  113        28.1  113
(0.1, 0.2]                    38.6  203        32.0  203
(0.2, 1.0]                    84.9  125        38.1  125

Con menos de un 2 % del volumen, la llegada y el VWAP casi coinciden (25,0 y 20,8 pb). Con más de un 20 %, la orden cuesta 84,9 pb frente a la llegada y 38,1 frente al VWAP: la mitad del coste ha desaparecido de la medida. Las propias ejecuciones son una parte grande del volumen del intervalo y su impacto sube el VWAP; la referencia se mueve con la orden. Un operador evaluado frente al VWAP no tiene incentivo para ejecutar con menos impacto; uno evaluado frente a la llegada, sí.

Qué hay dentro del coste

El generador guarda, para cada orden, la descomposición verdadera de su shortfall frente a la llegada (ejecucion4_ordenes_verdad.csv): media horquilla, impacto temporal, impacto permanente que pagan sus propias hijas y movimiento del medio sin impacto (el mercado). En la realidad esa descomposición no se conoce; aquí sirve para ver de qué está hecho el número:

verdad = pd.read_csv("ejecucion4_ordenes_verdad.csv")
ord_ = ord_.merge(verdad, on="id_orden")
partes = ["coste_horquilla_pb", "coste_temporal_pb",
          "coste_permanente_pb", "mercado_pb"]
tabla = ord_.groupby("ticker")[partes + ["frente_llegada"]].mean()
tabla.columns = ["horquilla", "temporal", "permanente", "mercado",
                 "IS"]
print(tabla.round(1))
ord_[partes + ["frente_llegada"]].std().round(1)

        horquilla  temporal  permanente  mercado    IS
ticker
ARVA          0.3      18.0         1.9     -3.8  16.3
BELO          0.9      24.1         2.5      3.2  30.7
CIRA          4.8      32.0         3.5      2.8  43.0
DUNE         13.6      47.0         5.3     10.3  76.1
coste_horquilla_pb      5.5
coste_temporal_pb      18.7
coste_permanente_pb     4.7
mercado_pb             66.5
frente_llegada         72.6
dtype: float64

El coste crece al pasar de ARVA a DUNE por las tres partes: la horquilla (0,3 a 13,6 pb), el temporal (18,0 a 47,0) y el permanente (1,9 a 5,3). El temporal es el grueso del coste en las cuatro acciones. Pero las desviaciones típicas cambian la lectura: 18,7 pb la del temporal, 4,7 la del permanente y 66,5 la del mercado. Orden a orden, lo que decide el resultado es el movimiento del mercado durante la ejecución, que el operador no controla. Para medir el impacto hacen falta muchas órdenes.

La lección 2.4.7 usó la ley de la raíz cuadrada como supuesto, con su coeficiente declarado: el impacto de una orden de $Q$ acciones ejecutada a lo largo de un día es aproximadamente $Y\,\sigma_d \sqrt{Q/\text{ADV}}$, con $Y$ del orden de 1 en los estudios empíricos (Tóth y otros, 2011, la explican y la recogen). Con las ejecuciones se puede estimar. Cada hija paga, frente al medio vigente justo antes, la media horquilla más el impacto temporal del minuto, $\eta\,\sigma_d\,\rho^{\beta}$; restando la media horquilla, una regresión sin constante sobre $\sigma_d\sqrt{\rho}$ da $\eta$, y los costes medios por tramo de participación en escala doble logarítmica dan el exponente $\beta$:

hijas = pd.read_parquet("ejecucion4_hijas.parquet")
hijas = hijas.merge(ord_[["id_orden", "lado_s", "sigma"]],
                    on="id_orden")
hijas["minuto"] = hijas["ts"].dt.floor("min")
propias = hijas.groupby(["ticker", "minuto"])["cantidad"].sum()
vol = barras.set_index(["ticker", "ts"])["volumen"]
vol.index = vol.index.set_names(["ticker", "minuto"])
hijas = hijas.join(propias.rename("propias"), on=["ticker", "minuto"])
hijas = hijas.join(vol, on=["ticker", "minuto"])
hijas["rho"] = hijas["propias"] / hijas["volumen"]
hijas["coste"] = (hijas["lado_s"] * (hijas["precio"]
                                     - hijas["medio_antes"])
                  - hijas["horquilla"] / 2) / hijas["medio_antes"]

x = hijas["sigma"] * np.sqrt(hijas["rho"])
eta = (x * hijas["coste"]).sum() / (x**2).sum()
print(f"{len(hijas)} hijas; participación mediana "
      f"{hijas['rho'].median():.3f}; eta = {eta:.3f}")

hijas["tramo"] = pd.qcut(np.log(hijas["rho"]), 20)
cubos = hijas.groupby("tramo", observed=True).agg(
    rho=("rho", "mean"), coste_sigma=("coste", "sum"),
    n=("coste", "size"), sigma=("sigma", "sum"))
cubos["coste_sigma"] = cubos["coste_sigma"] / cubos["sigma"]
pendiente, corte = np.polyfit(np.log(cubos["rho"]),
                              np.log(cubos["coste_sigma"]), 1)
print(f"Exponente {pendiente:.3f}; coeficiente {np.exp(corte):.3f}")

92409 hijas; participación mediana 0.098; eta = 0.400
Exponente 0.500; coeficiente 0.400

Con 92.409 hijas la recuperación es exacta en la práctica: $\eta = 0{,}400$, exponente 0,500 y coeficiente del ajuste logarítmico 0,400, los parámetros del generador. La figura de la cabecera enseña los veinte tramos sobre la recta de pendiente 1/2. La alternativa lineal (pendiente 1), la que supone la fórmula cerrada de Almgren-Chriss, solo coincide con la raíz en un punto: exagera el coste de participar mucho y lo infravalora con participaciones bajas.

En real, esta estimación pide los datos de cada hija con el medio vigente justo antes de enviarla, que un bróker no siempre da. Con datos de órdenes padre:

hs = hijas.assign(h=1e4 * hijas["horquilla"] / 2
                  / hijas["medio_antes"] * hijas["cantidad"])
hs = hs.groupby("id_orden")[["h", "cantidad"]].sum()
ord_ = ord_.join((hs["h"] / hs["cantidad"]).rename("horq_obs"),
                 on="id_orden")
ord_["impacto"] = ord_["frente_llegada"] - ord_["horq_obs"]


def sin_constante(y, x):
    b = (x * y).sum() / (x**2).sum()
    e = np.sqrt(((y - b * x)**2).sum() / (len(y) - 1)
                / (x**2).sum())
    return b, e


filas = {}
for alg, d in ord_.groupby("algoritmo"):
    x = 1e4 * d["sigma"] * np.sqrt(d["participacion_media"])
    b, e = sin_constante(d["impacto"], x)
    b2, _ = sin_constante(d["impacto"] - d["mercado_pb"], x)
    filas[alg] = {"Y": b, "error": e, "Y_sin_mercado": b2}
pd.DataFrame(filas).T.round(3)

          Y  error  Y_sin_mercado
IS    0.570  0.046          0.538
POV   0.442  0.055          0.446
TWAP  0.584  0.087          0.497
VWAP  0.573  0.052          0.494

Con las órdenes padre, frente a la llegada y restando la media horquilla observable, el coeficiente sale entre 0,44 y 0,58 según el algoritmo, con errores de 0,05 a 0,09. Está sesgado hacia arriba por dos razones: la orden padre paga también la mitad de su propio impacto permanente, que crece con el tamaño, y el ruido del mercado. Restando el movimiento verdadero del mercado (algo que en la realidad no se puede hacer) el sesgo se reduce (0,45 a 0,54) pero no desaparece, porque el permanente sigue dentro. POV sale más cerca de 0,4 porque sus órdenes son cortas (40 minutos de mediana) y acumulan menos permanente.

El impacto permanente es lo que queda en el precio después de la orden. Con lo observable, se estima restando al coste de cada orden la horquilla y el temporal (este con el $\eta$ de las hijas) y comparando el resto con $\sigma_d\,X/2$, la mitad del impacto permanente total de una orden de $X = q/\text{ADV}$:

t = hijas.assign(t=1e4 * eta * hijas["sigma"] * np.sqrt(hijas["rho"])
                 * hijas["cantidad"])
t = t.groupby("id_orden")["t"].sum() / hs["cantidad"]
ord_["temporal_est"] = ord_["id_orden"].map(t)
X = ord_["cantidad_ejecutada"] / ord_["adv"]
z = 1e4 * ord_["sigma"] * X / 2
resto = ord_["impacto"] - ord_["temporal_est"]
g, e = sin_constante(resto, z)
g2, e2 = sin_constante(resto - ord_["mercado_pb"], z)
print(f"gamma con lo observable: {g:.2f} ± {e:.2f}")
print(f"gamma restando el mercado (verdad): {g2:.3f} ± {e2:.3f}")
n_necesarias = len(ord_) * (e / 0.10)**2
print(f"Órdenes para un error de ±0,10: {n_necesarias:,.0f}")

gamma con lo observable: 1.76 ± 0.42
gamma restando el mercado (verdad): 0.983 ± 0.008
Órdenes para un error de ±0,10: 12,881

Con las 727 órdenes, $\gamma = 1{,}76 \pm 0{,}42$, cuando la verdad es 1,0: a menos de dos errores típicos, pero inútil para decidir. Restando el movimiento verdadero del mercado sale $0{,}983 \pm 0{,}008$: el permanente está ahí, tapado por un ruido decenas de veces mayor. Para un error de ±0,10 harían falta unas 12.900 órdenes como estas. Es la razón por la que las estimaciones del impacto permanente se hacen con cientos de miles de órdenes de muchas mesas, y por la que en un modelo de costes propio el permanente es un supuesto que hay que declarar y variar, no una cifra medida.

Urgencia: el problema de Almgren y Chriss

Ejecutar deprisa sube el impacto temporal (más participación); ejecutar despacio deja la orden expuesta más tiempo al movimiento del mercado. Almgren y Chriss (2000) lo plantearon como una cartera: elegir el calendario de ejecución $x_k$ (lo que queda por ejecutar al final de cada tramo) que minimiza

con $\lambda$ la aversión al riesgo. Con impacto temporal lineal en la velocidad de ejecución y liquidez constante, la solución es cerrada: $x(t) = X \sinh(\kappa (T - t)) / \sinh(\kappa T)$, con $\kappa^2 = \lambda \sigma^2 / \tilde\eta$. Con $\lambda = 0$ sale TWAP; cuanto mayor $\lambda$, más se adelanta la ejecución. Las órdenes IS de la mesa siguen esa trayectoria con urgencias $\kappa \in \{0{,}5;\ 2;\ 5\}$. Con la descomposición verdadera se ve el intercambio:

d = ord_[ord_["algoritmo"] == "IS"].copy()
d["mercado_en_sigmas"] = d["mercado_pb"] / (1e4 * d["sigma"])
d.groupby("urgencia_kappa").agg(
    ordenes=("id_orden", "size"), minutos=("minutos", "median"),
    temporal=("coste_temporal_pb", "mean"),
    riesgo=("mercado_en_sigmas", "std")).round(3)

                ordenes  minutos  temporal  riesgo
urgencia_kappa
0.5                  40     95.5    25.028   0.237
2.0                  53     98.0    31.611   0.233
5.0                  61     60.0    38.690   0.136

Al pasar de $\kappa = 0{,}5$ a 5, el impacto temporal medio sube de 25,0 a 38,7 pb y el riesgo, medido como la desviación típica del movimiento del mercado durante la orden en desviaciones diarias, baja de 0,237 a 0,136. Las órdenes urgentes duran menos (60 minutos de mediana frente a 96) y pagan más por participar más. Entre 0,5 y 2 apenas cambia el riesgo: las dos urgencias más bajas tienen duraciones parecidas.

Los datos no cumplen el supuesto de la fórmula cerrada: el impacto temporal es $\eta\,\sigma_d\,\rho^{0{,}5}$, en raíz de la participación, y el volumen sigue una curva en U. El coste de un tramo $k$ con $n_k$ acciones y un volumen de mercado $V_k$ es $\eta\,\sigma_d\, n_k^{1{,}5} / V_k^{0{,}5}$: convexo en $n_k$, así que el problema sigue siendo convexo y cvxpy lo resuelve con cp.power(u, 1.5) (un cono de potencia que CLARABEL admite). El caso: vender un 10 % del ADV de DUNE en una sesión, en 13 tramos de 30 minutos, con la curva de volumen media de los datos, $\eta$ y $\gamma$ del generador y la varianza de Almgren y Chriss con la volatilidad repartida por igual en la sesión (los datos la tienen en U: es una simplificación). El permanente no depende del calendario con impacto lineal permanente: suma siempre $\gamma\,\sigma_d\,X/2$.

import cvxpy as cp

T_ = "DUNE"
b = por_ticker[T_]
minuto = (b.index.hour * 60 + b.index.minute) - (13 * 60 + 30)
curva = b.groupby([b.index.date, minuto // 30])["volumen"].sum()
V = curva.groupby(level=1).mean().to_numpy()   # 13 tramos de 30 min
N, tau = len(V), 1 / len(V)                    # tau en días
ADV, sig = acc.at[T_, "adv_acciones"], acc.at[T_, "sigma_diaria"]
Xq = 0.10 * ADV                                # 10 % del ADV
ETA, GAMMA = 0.40, 1.0                         # los del generador


def coste(u):
    """Coste esperado y desviación (pb) del calendario u (fracción
    de la orden en cada tramo). Temporal en raíz de la participación,
    permanente lineal, varianza de lo que queda por ejecutar."""
    temporal = 1e4 * ETA * sig * np.sum(u**1.5 * np.sqrt(Xq / V))
    permanente = 1e4 * GAMMA * sig * Xq / ADV / 2
    resto = 1 - np.cumsum(u)
    desv = 1e4 * sig * np.sqrt(tau * np.sum(resto**2))
    return temporal + permanente, desv


def optimo(lam):
    u = cp.Variable(N, nonneg=True)
    temporal = 1e4 * ETA * sig * cp.sum(
        cp.multiply(np.sqrt(Xq / V), cp.power(u, 1.5)))
    varianza = (1e4 * sig)**2 * tau * cp.sum_squares(1 - cp.cumsum(u))
    cp.Problem(cp.Minimize(temporal + lam * varianza),
               [cp.sum(u) == 1]).solve(solver="CLARABEL")
    return np.clip(u.value, 0, None)


print("Volumen medio por tramo (miles):", np.round(V / 1e3, 1))
cal = {"TWAP": np.full(N, 1 / N), "VWAP": V / V.sum(),
       "optimo_lam0": optimo(0.0)}
pd.DataFrame({k: coste(u) for k, u in cal.items()},
             index=["coste_pb", "desv_pb"]).T.round(2)

Volumen medio por tramo (miles): [98.9 68.2 51.2 45.5 44.3 43.1 44.  43.4 44.9 44.5 47.4 52.5 92. ]
             coste_pb  desv_pb
TWAP            51.89   163.18
VWAP            50.78   153.48
optimo_lam0     50.78   153.48

VWAP cuesta 50,8 pb y TWAP 51,9: con impacto en raíz, el óptimo sin aversión al riesgo es VWAP, no TWAP. El coste total es proporcional a $\sum_k n_k^{1,5} / V_k^{0,5}$, y minimizarlo con $\sum n_k$ fijo da $n_k \propto V_k$: participación constante. La solución numérica con $\lambda = 0$ lo confirma (los mismos 50,78 pb). La desviación típica del coste, 153 pb, es tres veces el coste esperado: una orden del 10 % del ADV de una acción con un 3 % de volatilidad diaria es, sobre todo, una apuesta involuntaria sobre el precio de la sesión.

Para usar la fórmula cerrada hay que linealizar el impacto: elegir un $\tilde\eta$ lineal que dé el mismo coste por acción que la raíz con el calendario TWAP, y aplicar la fórmula del seno hiperbólico. Se comparan los dos calendarios para varias $\lambda$, evaluados los dos con el impacto en raíz de los datos:

# eta lineal que da el mismo coste por acción que la raíz con TWAP
n_tw = Xq / N
eta_lin = ETA * sig * tau / np.sqrt(n_tw * V.mean())


def cerrado(lam):
    """Calendario de Almgren-Chriss (2000) con impacto lineal:
    x(t) = sinh(k (T - t)) / sinh(k T), T = 1 día."""
    k = np.sqrt(lam * 1e4 * sig**2 / (eta_lin * Xq))
    t = np.arange(N + 1) / N
    x = np.sinh(k * (1 - t)) / np.sinh(k) if k > 1e-9 else 1 - t
    return -np.diff(x)


filas = {}
for lam in (1e-4, 1e-3, 3e-3, 1e-2, 3e-2):
    for nombre, u in (("raiz_cvxpy", optimo(lam)),
                      ("lineal_cerrado", cerrado(lam))):
        c, sd = coste(u)
        filas[(lam, nombre)] = {"coste_pb": c, "desv_pb": sd,
                                "objetivo": c + lam * sd**2,
                                "primer_tramo": u[0]}
frontera = pd.DataFrame(filas).T
frontera.round(3)

                       coste_pb  desv_pb  objetivo  primer_tramo
0.0001 raiz_cvxpy        50.911  144.054    52.986         0.160
       lineal_cerrado    51.809  160.198    54.375         0.083
0.0010 raiz_cvxpy        55.999   96.594    65.329         0.306
       lineal_cerrado    52.101  139.021    71.428         0.126
0.0030 raiz_cvxpy        65.556   61.825    77.024         0.465
       lineal_cerrado    54.772  112.042    92.432         0.192
0.0100 raiz_cvxpy        80.118   32.253    90.520         0.669
       lineal_cerrado    62.446   77.098   121.888         0.320
0.0300 raiz_cvxpy        93.072   15.553   100.329         0.826
       lineal_cerrado    72.325   49.674   146.349         0.487

Con la misma $\lambda$, el calendario de la fórmula lineal es mucho menos urgente que el óptimo con raíz: con $\lambda = 0{,}01$ ejecuta un 32 % en la primera media hora frente al 67 %, y su objetivo (coste más $\lambda$ por la varianza) es de 121,9 frente a 90,5. La razón está en la forma del impacto: con el lineal, el coste total de un tramo crece con $n^2$; con la raíz, con $n^{1,5}$. La fórmula lineal exagera lo que cuesta concentrar la ejecución y la retrasa.

Como familia, en cambio, los calendarios de seno hiperbólico no están lejos de la frontera: quedan entre 1 y 2 pb por encima para el mismo riesgo. El error de la linealización no es tanto la forma del calendario como la urgencia que asigna a cada $\lambda$: quien calibra su aversión al riesgo con la fórmula lineal y ejecuta con impacto en raíz obtiene un calendario más lento y más arriesgado del que quería.

Errores típicos

• Medir solo lo ejecutado. Sin el coste de oportunidad, un algoritmo que no completa parece barato.
• Evaluar frente al VWAP del intervalo las órdenes grandes: la referencia incluye su propio impacto.
• Confundir ruido con coste. Con 73 pb de desviación por orden, comparar dos algoritmos con unas decenas de órdenes no dice nada; hace falta la media de muchas y su error típico.
• Estimar el impacto con las órdenes padre sin separar el permanente: el coeficiente sale inflado.
• Tratar el permanente como medido. Con una sola mesa es un supuesto.
• Usar la fórmula cerrada de Almgren-Chriss con impacto en raíz sin decirlo: la urgencia sale equivocada.
• Olvidar la curva de volumen. Con liquidez en U, el calendario neutral al riesgo es VWAP, no TWAP.

Ejercicios

1. Estima el coeficiente $Y$ de la ley de la raíz cuadrada en su forma clásica, $Y\,\sigma_d\sqrt{Q/\text{ADV}}$, solo con las órdenes VWAP que duran la sesión entera (390 minutos o casi). ¿Cuántas hay y qué error tiene la estimación?

dia = ord_[(ord_["algoritmo"] == "VWAP") & (ord_["minutos"] >= 380)]
x = 1e4 * dia["sigma"] * np.sqrt(dia["cantidad_ejecutada"]
                                 / dia["adv"])
Y, err = sin_constante(dia["impacto"], x)
Y2, _ = sin_constante(dia["impacto"] - dia["mercado_pb"], x)
print(f"{len(dia)} órdenes de un día: Y = {Y:.2f} ± {err:.2f};"
      f" sin el mercado, {Y2:.2f}")

18 órdenes de un día: Y = 0.94 ± 0.81; sin el mercado, 0.59

2. Repite el óptimo de Almgren-Chriss con $\lambda = 0{,}01$ para un 10 % del ADV de ARVA, la acción más líquida y menos volátil. ¿Cuánto ejecuta en la primera media hora y cuánto cuesta frente a DUNE?

def caso(ticker, pct=0.10, lam=1e-2):
    global V, ADV, sig, Xq
    b = por_ticker[ticker]
    m = (b.index.hour * 60 + b.index.minute) - (13 * 60 + 30)
    curva = b.groupby([b.index.date, m // 30])["volumen"].sum()
    V = curva.groupby(level=1).mean().to_numpy()
    ADV = acc.at[ticker, "adv_acciones"]
    sig = acc.at[ticker, "sigma_diaria"]
    Xq = pct * ADV
    u = optimo(lam)
    c, sd = coste(u)
    return {"primer_tramo": u[0], "coste_pb": c, "desv_pb": sd,
            "vwap_pb": coste(V / V.sum())[0]}


pd.DataFrame({t: caso(t) for t in ("ARVA", "DUNE")}).T.round(3)

      primer_tramo  coste_pb  desv_pb  vwap_pb
ARVA         0.547    33.311   23.138   23.912
DUNE         0.669    80.118   32.253   50.776

3. Un gestor decide comprar al cierre de ayer y su operador recibe la orden a mediodía. El retraso medio de la mesa es de varios puntos básicos, positivo o negativo según el día. ¿Quién debería responder de ese coste y cómo lo medirías para separarlo de la ejecución?