La prueba \(t\) del capítulo 9 compara medias. Muchas preguntas de un equipo de datos no tienen una media que comparar, porque lo que observas es una categoría: la respuesta fue correcta, incorrecta o venció el tiempo; el request cayó en la réplica 1, 2, 3 o 4; el usuario aceptó la recomendación o no. Con categorías trabajas contando, y la pregunta pasa a ser si los conteos observados se parecen a los que esperarías bajo una hipótesis.
En este capítulo recorres la familia de pruebas que responde esa pregunta. Todas comparten una idea: medir la distancia entre lo observado y lo esperado, y compararla con la que produciría el azar.
El contexto es el servicio de inferencia del capítulo 9: un balanceador que reparte requests, un A/B entre la versión actual del modelo y una cuantizada, y la distribución de clases que predice un clasificador en producción. Los datos son simulados: la mayoría se generan con rng.multinomial o rng.binomial a partir de proporciones que declara cada etapa, y unos pocos ejemplos pequeños son tablas escritas a mano.
scipy.stats) y Matplotlib.Al terminar este capítulo vas a ser capaz de:
Ejecuta las celdas en orden: las etapas comparten variables y el generador aleatorio
rng.
El capítulo 9 presentó \(H_0\), \(H_1\), \(\alpha\), el p-value y los errores tipo I y II; aquí los usas sin volver a definirlos. Toda prueba de este capítulo y del siguiente recorre los mismos cinco pasos, y solo cambian el estadístico y su distribución.
| Paso | Qué haces | En la prueba \(t\) | En este capítulo |
|---|---|---|---|
| 1. Hipótesis | Escribes \(H_0\) y \(H_1\) antes de mirar el resultado | \(H_0:\mu=\mu_0\) | \(H_0\): las proporciones son las indicadas, o las dos variables son independientes |
| 2. Estadístico | Resumes en un número cuánto se alejan los datos de \(H_0\) | \(t\) | \(\chi^2\) |
| 3. Distribución de referencia | Sabes cómo se comporta ese número cuando \(H_0\) es cierta | \(t\) de Student | chi-cuadrado |
| 4. p-value y decisión | Área de cola más allá del valor observado, comparada con \(\alpha\) | dos colas | cola derecha |
| 5. Tamaño del efecto | Dices cuánto, en unidades que se entiendan | diferencia en ms | diferencia de proporciones, V de Cramér |
El quinto paso forma parte de la prueba: un p-value pequeño dice que el azar explica mal los datos, pero no si la diferencia es grande.
Una variable categórica toma valores que son etiquetas, no cantidades: correcta, incorrecta, timeout. Su resumen natural es la tabla de frecuencias, que cuenta cuántas observaciones cayeron en cada categoría, con conteos (frecuencias absolutas) y proporciones (frecuencias relativas). Las pruebas trabajan con los conteos, no con los porcentajes: el 60% de 10 casos y el 60% de 10 000 no pesan lo mismo. Por ejemplo, si de 600 requests 178 fueron a la réplica 1, el conteo es 178 y la proporción \(178/600=0.297\).
La frecuencia observada \(O_i\) es lo que contaste en la categoría \(i\). La frecuencia esperada \(E_i\) es lo que contarías, en promedio, si \(H_0\) fuera cierta: es la referencia contra la que mides la sorpresa, y puede tener decimales aunque un conteo real no los tenga. Si el balanceador reparte por igual entre 4 réplicas, de 600 requests esperas \(600\times 0.25=150\) en cada una.
El estadístico suma, sobre todas las categorías, la distancia entre observado y esperado:
\[ \chi^2=\sum_{i}\frac{(O_i-E_i)^2}{E_i}. \]
Cada parte de la fórmula tiene una razón:
Si \(H_0\) es cierta, \(\chi^2\) tiende a tomar valores moderados, cercanos a sus grados de libertad, y cuanto más se alejan los datos de \(H_0\), más crece. Por eso la prueba mira solo la cola derecha. Una categoría con \(O=178\) y \(E=150\) aporta \((178-150)^2/150=5.23\) al estadístico.
La distribución chi-cuadrado con \(k\) grados de libertad, \(\chi^2_k\), es la distribución de la suma de \(k\) normales estándar independientes elevadas al cuadrado:
\[ Z_1^2+Z_2^2+\dots+Z_k^2\sim\chi^2_k . \]
Solo toma valores positivos, es asimétrica hacia la derecha y su media es \(k\). Los grados de libertad cuentan cuántas discrepancias pueden variar libremente: con 4 categorías y un total fijo de 600, conocidas tres frecuencias la cuarta queda determinada, así que hay \(4-1=3\) grados de libertad.
Aparece aquí porque, con \(n\) grande, cada término \((O_i-E_i)/\sqrt{E_i}\) se comporta aproximadamente como una normal, y el estadístico es una suma de sus cuadrados. Son \(k\) términos, pero no independientes: como los conteos suman \(n\), las discrepancias suman 0 y una queda determinada por las demás. Esa restricción quita un grado de libertad y deja \(\chi^2_{k-1}\). La palabra aproximadamente importa: es una aproximación para muestras grandes, y la etapa 3 muestra cuándo falla.
| Pregunta | Prueba | \(H_0\) | Grados de libertad |
|---|---|---|---|
| ¿Una variable se reparte según las proporciones que esperabas? | Bondad de ajuste | \(p_i=p_{i,0}\) para toda categoría | \(k-1\) |
| ¿Dos variables categóricas están relacionadas? | Independencia | Las dos variables son independientes | \((f-1)(c-1)\) |
| ¿Dos grupos tienen la misma proporción de éxitos? | Dos proporciones | \(p_A=p_B\) | 1 |
En la tabla, \(k\) es el número de categorías, \(f\) el de filas y \(c\) el de columnas. La tercera prueba es el caso \(2\times 2\) de la segunda.
En una tabla de contingencia, si las dos variables fueran independientes, la frecuencia esperada de cada celda sería
\[ E_{ij}=\frac{(\text{total de la fila } i)\times(\text{total de la columna } j)}{n}. \]
Es la regla del producto de probabilidades independientes del capítulo 7, \(P(A\cap B)=P(A)\,P(B)\), escrita con conteos: si el 10% de las respuestas es incorrecto y la versión no influye, esperas un 10% de incorrectas en cada versión. Con 1000 requests en la versión A, 168 respuestas incorrectas en total y \(n=2000\), esperas \(1000\times 168/2000=84\) incorrectas en A.
Con \(n\) grande, discrepancias muy pequeñas terminan dando un p-value pequeño. Estas medidas acompañan a la prueba:
\[ V=\sqrt{\frac{\chi^2}{n\,(m-1)}},\qquad m=\min(f,c). \]
Va de 0 (sin asociación) a 1 (asociación perfecta) y no crece sistemáticamente con \(n\); en muestras pequeñas es imprecisa y tiende a sobrestimar.
\[ r_{ij}=\frac{O_{ij}-E_{ij}}{\sqrt{E_{ij}\,(1-p_{i\cdot})\,(1-p_{\cdot j})}}, \]
donde \(p_{i\cdot}\) y \(p_{\cdot j}\) son las proporciones de la fila y de la columna. Bajo \(H_0\) se comportan aproximadamente como una normal estándar, así que valores más allá de \(\pm 2\) señalan las celdas que más contribuyen.
Dos preguntas ordenan todo el mapa: qué tipo de dato es la respuesta (categoría o número) y cómo se obtuvieron los grupos (independientes o pareados).
Buena parte de lo que produce y consume un sistema de aprendizaje automático es categórico:
Si lo que mides es una categoría, partes de una tabla de conteos y de una pregunta sobre lo esperado; la proporción ocupa el lugar que tenía la media.
«En
[población], con \(n=\)[valor], se cruzó[variable 1]con[variable 2]. La frecuencia esperada mínima fue[valor]. Se probó[H0]: χ²([gl]) =[valor], p =[valor]. A α =[valor]se[rechaza o no se rechaza]H0. El tamaño del efecto es[V o diferencia de proporciones con IC]. La discrepancia se concentra en[celdas]. Limitación:[n, esperados pequeños, varias pruebas o diseño observacional].»
Mira siempre la tabla de esperados antes que el p-value. La etapa 10 llena esta plantilla con el A/B de resultados.
La celda carga NumPy, pandas, Matplotlib y el módulo stats de SciPy, fija un estilo de gráficos y crea el generador rng con la semilla 42, del que salen todas las simulaciones.
%matplotlib inline
%config InlineBackend.figure_format = "retina"
import sys
print("Python:", sys.version.split()[0])
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
SEMILLA = 42 # semilla central del capítulo
rng = np.random.default_rng(SEMILLA)
plt.rcParams["figure.figsize"] = (8, 5)
plt.rcParams["axes.grid"] = True
plt.rcParams["grid.alpha"] = 0.25Python: 3.13.14El entorno queda con Python 3.13.14.
La siguiente celda, plegada en la versión web, define las funciones de gráficos. Solo dibujan: ningún cálculo estadístico ocurre dentro de ellas, y los números se obtienen a la vista en cada etapa.
def graficar_suma_cuadrados(sumas, titulo):
fig, axes = plt.subplots(1, len(sumas), figsize=(11, 3.4))
for ax, (k, q) in zip(axes, sumas.items()):
tope = np.quantile(q, 0.995)
xs = np.linspace(0.02, tope, 300)
ax.hist(q, bins=60, range=(0, tope), density=True, alpha=0.5, label="simulación")
ax.plot(xs, stats.chi2.pdf(xs, k), color="C3", label=fr"$\chi^2_{{{k}}}$")
ax.axvline(k, color="0.3", linestyle=":", label=f"media = {k}")
if k == 1:
ax.set_ylim(0, 1.5)
ax.set_title(f"suma de {k} normales al cuadrado", fontsize=10)
ax.set_xlabel("valor de la suma")
ax.legend(fontsize=8)
axes[0].set_ylabel("densidad")
fig.suptitle(titulo, fontsize=11)
plt.tight_layout()
plt.show()
def graficar_obs_esp(etiquetas, obs, esp, titulo, ylabel="requests"):
x = np.arange(len(etiquetas))
fig, ax = plt.subplots(figsize=(7.6, 4))
ax.bar(x - 0.2, obs, width=0.4, label="observado")
ax.bar(x + 0.2, esp, width=0.4, color="0.7", label="esperado bajo H0")
for xi, o in zip(x, obs):
ax.text(xi - 0.2, o, f"{o:g}", ha="center", va="bottom", fontsize=8)
ax.set_xticks(x)
ax.set_xticklabels(etiquetas)
ax.set_ylabel(ylabel)
ax.set_title(titulo)
ax.legend(fontsize=8)
plt.tight_layout()
plt.show()
def graficar_cola_chi2(estadistico, gl, titulo, alpha=0.05):
critico = stats.chi2.ppf(1 - alpha, gl)
tope = max(estadistico, critico) * 1.25 + 2
xs = np.linspace(0.01, tope, 500)
fx = stats.chi2.pdf(xs, gl)
fig, ax = plt.subplots(figsize=(8, 3.8))
ax.plot(xs, fx, label=fr"$\chi^2_{{{gl}}}$ (si $H_0$ es cierta)")
ax.fill_between(xs, fx, where=xs >= critico, alpha=0.3, color="C3", label=fr"zona de rechazo ($\alpha$ = {alpha:g})")
ax.axvline(critico, color="C3", linestyle=":", label=f"valor crítico = {critico:.2f}")
ax.axvline(estadistico, color="0.15", linestyle="--", label=f"observado = {estadistico:.2f}")
if gl <= 2:
ax.set_ylim(0, min(fx.max(), 0.6))
ax.set_xlabel(r"estadístico $\chi^2$")
ax.set_ylabel("densidad")
ax.set_title(titulo)
ax.legend(fontsize=8)
plt.tight_layout()
plt.show()
def graficar_calibracion(sim_a, sim_b, gl, etiqueta_a, etiqueta_b, titulo, alpha=0.05):
critico = stats.chi2.ppf(1 - alpha, gl)
xs = np.linspace(0.05, 16, 300)
fig, axes = plt.subplots(1, 2, figsize=(11, 3.8), sharey=True)
for ax, sim, etiqueta in zip(axes, (sim_a, sim_b), (etiqueta_a, etiqueta_b)):
ax.hist(sim, bins=48, range=(0, 16), density=True, alpha=0.5, label="simulación bajo H0")
ax.plot(xs, stats.chi2.pdf(xs, gl), color="C3", label=fr"$\chi^2_{{{gl}}}$ teórica")
ax.axvline(critico, color="0.3", linestyle=":", label=f"valor crítico = {critico:.2f}")
ax.set_title(etiqueta, fontsize=10)
ax.set_xlabel(r"estadístico $\chi^2$")
ax.legend(fontsize=8)
axes[0].set_ylabel("densidad")
axes[0].set_ylim(0, 0.45)
fig.suptitle(titulo, fontsize=11)
plt.tight_layout()
plt.show()
def graficar_perfiles(tabla, titulo):
perfiles = tabla.div(tabla.sum(axis=1), axis=0) * 100
x = np.arange(tabla.shape[1])
ancho = 0.8 / tabla.shape[0]
fig, axes = plt.subplots(1, 2, figsize=(11, 3.8))
for i, fila in enumerate(tabla.index):
for ax in axes:
ax.bar(x + i * ancho - 0.4 + ancho / 2, perfiles.loc[fila], width=ancho, label=fila)
for ax in axes:
ax.set_xticks(x)
ax.set_xticklabels(tabla.columns)
ax.legend(fontsize=8)
axes[0].set_ylabel("% de los requests de la versión")
axes[0].set_title("escala completa", fontsize=10)
axes[1].set_ylim(0, perfiles.iloc[:, 1:].values.max() * 1.35)
axes[1].set_title("acercamiento a las categorías poco frecuentes", fontsize=10)
fig.suptitle(titulo, fontsize=11)
plt.tight_layout()
plt.show()
def graficar_residuos(residuos, titulo):
fig, ax = plt.subplots(figsize=(6.6, 2.9))
tope = max(3.0, np.abs(residuos.values).max())
im = ax.imshow(residuos.values, cmap="RdBu_r", vmin=-tope, vmax=tope, aspect="auto")
for i in range(residuos.shape[0]):
for j in range(residuos.shape[1]):
valor = residuos.values[i, j]
ax.text(j, i, f"{valor:+.2f}", ha="center", va="center", fontsize=11,
color="white" if abs(valor) > 0.6 * tope else "black")
ax.set_xticks(range(residuos.shape[1]))
ax.set_xticklabels(residuos.columns)
ax.set_yticks(range(residuos.shape[0]))
ax.set_yticklabels(residuos.index)
ax.grid(False)
fig.colorbar(im, ax=ax, label="residuo ajustado")
ax.set_title(titulo, fontsize=10)
plt.tight_layout()
plt.show()
def graficar_dos_proporciones(p_a, p_b, ic_a, ic_b, dif, ic_dif, titulo):
fig, axes = plt.subplots(1, 2, figsize=(10.5, 3.4))
ax = axes[0]
for y, p, ic, nombre, color in ((1, p_a, ic_a, "A (actual)", "C0"), (0, p_b, ic_b, "B (cuantizado)", "C1")):
ax.hlines(y, ic[0] * 100, ic[1] * 100, color=color, linewidth=4)
ax.plot(p * 100, y, "o", color=color)
ax.text(p * 100, y + 0.15, f"{p * 100:.1f}%", ha="center", fontsize=9)
ax.set_yticks([1, 0])
ax.set_yticklabels(["A (actual)", "B (cuantizado)"])
ax.set_ylim(-0.6, 1.6)
ax.set_xlabel("tasa de aceptación (%) con IC 95%")
ax.set_title("cada versión por separado", fontsize=10)
ax = axes[1]
ax.hlines(0, ic_dif[0] * 100, ic_dif[1] * 100, color="C2", linewidth=4)
ax.plot(dif * 100, 0, "o", color="C2")
ax.axvline(0, color="0.3", linestyle="--", label="sin diferencia")
ax.set_yticks([])
ax.set_xlabel("diferencia B − A (puntos porcentuales) con IC 95%")
ax.set_title("la diferencia", fontsize=10)
ax.legend(fontsize=8)
fig.suptitle(titulo, fontsize=11)
plt.tight_layout()
plt.show()
def graficar_deriva(clases, p_ref, proporciones, pvalores, titulo, alpha=0.05):
fig, axes = plt.subplots(1, 2, figsize=(11.5, 3.8), gridspec_kw={"width_ratios": [1.6, 1]})
ax = axes[0]
x = np.arange(len(clases))
ancho = 0.8 / (len(proporciones) + 1)
ax.bar(x - 0.4 + ancho / 2, p_ref * 100, width=ancho, color="0.6", label="referencia")
for s, prop in enumerate(proporciones, start=1):
ax.bar(x - 0.4 + ancho / 2 + s * ancho, prop * 100, width=ancho, label=f"semana {s}")
ax.set_xticks(x)
ax.set_xticklabels(clases)
ax.set_ylabel("% de las predicciones")
ax.set_title("distribución de clases predichas", fontsize=10)
ax.legend(fontsize=8)
ax = axes[1]
semanas = np.arange(1, len(pvalores) + 1)
ax.plot(semanas, pvalores, "o-")
ax.axhline(alpha, color="C3", linestyle="--", label=fr"$\alpha$ = {alpha:g}")
ax.axhline(alpha / len(pvalores), color="C1", linestyle=":", label=fr"$\alpha$/{len(pvalores)} = {alpha / len(pvalores):g}")
ax.set_yscale("log")
ax.set_xticks(semanas)
ax.set_xlabel("semana")
ax.set_ylabel("p-value (escala logarítmica)")
ax.set_title("prueba de bondad de ajuste por semana", fontsize=10)
ax.legend(fontsize=8)
fig.suptitle(titulo, fontsize=11)
plt.tight_layout()
plt.show()
Las 12 etapas siguen el orden del mapa: primero la distribución de referencia, después una variable (bondad de ajuste), luego dos variables (independencia y proporciones) y tres situaciones que cambian la lectura o la prueba: pocos datos, una muestra enorme y datos pareados. Las tres últimas aplican todo a un reporte, a un cálculo a mano y a la vigilancia de un clasificador en producción.
Antes de usar la distribución conviene verla nacer. Generas 20 000 grupos de \(k\) normales estándar con rng.normal, para \(k=1\), 3 y 8, elevas cada una al cuadrado y sumas las de cada grupo. El histograma de esas sumas debería coincidir con la densidad \(\chi^2_k\).
sumas = {}
for k in (1, 3, 8):
z = rng.normal(size=(20_000, k))
sumas[k] = (z ** 2).sum(axis=1)
print(f"k = {k}: media simulada = {sumas[k].mean():.3f} (teórica {k})",
f" percentil 95 simulado = {np.quantile(sumas[k], 0.95):.2f}",
f"(teórico {stats.chi2.ppf(0.95, k):.2f})")
graficar_suma_cuadrados(sumas, "La chi-cuadrado es una suma de normales estándar al cuadrado")k = 1: media simulada = 1.009 (teórica 1) percentil 95 simulado = 3.87 (teórico 3.84)
k = 3: media simulada = 3.013 (teórica 3) percentil 95 simulado = 7.89 (teórico 7.81)
k = 8: media simulada = 8.009 (teórica 8) percentil 95 simulado = 15.63 (teórico 15.51)
Las tres simulaciones siguen la curva teórica. La media simulada coincide con los grados de libertad (1.009, 3.013 y 8.009 frente a 1, 3 y 8), y el percentil 95 simulado queda cerca del teórico (3.87 frente a 3.84, 7.89 frente a 7.81, 15.63 frente a 15.51). En la figura, la curva de \(k=1\) se dispara cerca de 0 y las de 3 y 8 tienen su pico cada vez más a la derecha.
Dos rasgos se usan en todo el capítulo: la distribución solo toma valores positivos, porque es una suma de cuadrados, y se corre hacia la derecha al aumentar los grados de libertad. Por eso un mismo valor del estadístico, por ejemplo 8, es raro con 3 grados de libertad y típico con 8: un \(\chi^2\) se lee siempre junto con sus grados de libertad.
Un balanceador distribuye los requests entre 4 réplicas del servicio y está configurado para repartir por igual. Registras 600 requests y cuentas cuántos recibió cada réplica. Los conteos se simulan con rng.multinomial a partir de un reparto real de 0.28, 0.26, 0.25 y 0.21, que quien analiza no conoce.
\[ H_0: p_1=p_2=p_3=p_4=0.25,\qquad H_1:\text{al menos una } p_i\neq 0.25,\qquad \alpha=0.05. \]
Bajo \(H_0\) esperas \(600\times 0.25=150\) requests por réplica, con \(4-1=3\) grados de libertad.
replicas = ["réplica 1", "réplica 2", "réplica 3", "réplica 4"]
n_req = 600
p_diseno = np.array([0.28, 0.26, 0.25, 0.21]) # el reparto real, desconocido para quien analiza
obs = rng.multinomial(n_req, p_diseno)
p_h0 = np.full(4, 0.25)
esp = n_req * p_h0
aporte = (obs - esp) ** 2 / esp
chi2_obs = aporte.sum()
gl = len(obs) - 1
p_manual = stats.chi2.sf(chi2_obs, gl)
tabla_ba = pd.DataFrame({"observado": obs, "esperado": esp, "O − E": obs - esp, "(O − E)²/E": aporte.round(3)},
index=replicas)
print(tabla_ba.to_string())
print()
print("χ² a mano =", round(chi2_obs, 3), " gl =", gl)
print("valor crítico (α = 0.05) =", round(stats.chi2.ppf(0.95, gl), 3))
print("p-value a mano =", round(p_manual, 5))
res = stats.chisquare(obs, f_exp=esp)
print("SciPy chisquare: χ² =", round(res.statistic, 3), " p =", round(res.pvalue, 5))
print("¿Rechazar H0 a α = 0.05?", res.pvalue < 0.05)
graficar_obs_esp(replicas, obs, esp, "Requests por réplica: observado frente a esperado (n = 600)")
graficar_cola_chi2(chi2_obs, gl, "Dónde cae el estadístico observado si el reparto fuera uniforme") observado esperado O − E (O − E)²/E
réplica 1 178 150.0 28.0 5.227
réplica 2 155 150.0 5.0 0.167
réplica 3 161 150.0 11.0 0.807
réplica 4 106 150.0 -44.0 12.907
χ² a mano = 19.107 gl = 3
valor crítico (α = 0.05) = 7.815
p-value a mano = 0.00026
SciPy chisquare: χ² = 19.107 p = 0.00026
¿Rechazar H0 a α = 0.05? True
El estadístico vale 19.11 con 3 grados de libertad. El valor crítico al 5% es 7.81 y el p-value es 0.00026: si el reparto fuera uniforme, una discrepancia así o mayor aparecería en menos de 3 de cada 10 000 muestras. Rechazas \(H_0\), y en la segunda figura la línea del valor observado queda muy a la derecha de la zona de rechazo.
La tabla dice además dónde está la discrepancia. La réplica 4 aporta 12.91 y la réplica 1, 5.23: entre las dos suman 18.13 de los 19.11 puntos del estadístico. La réplica 4 recibió 106 requests (17.7%) y la réplica 1, 178 (29.7%), cuando esperabas un 25% en cada una.
La prueba detecta que el reparto observado no es compatible con uno uniforme, pero no dice por qué: la causa puede estar en la configuración, en revisiones de salud que sacan una réplica de rotación o en otra cosa que estos datos no registran.
La distribución \(\chi^2\) es una aproximación para muestras grandes. Lo compruebas simulando 20 000 experimentos en los que \(H_0\) es cierta y contando cuántas veces la prueba rechaza: si la aproximación funciona, debe rechazar cerca del 5% de las veces. Comparas dos escenarios con 4 categorías y 3 grados de libertad:
def chi2_bajo_h0(n, p, n_rep=20_000):
p = np.asarray(p)
conteos = rng.multinomial(n, p, size=n_rep)
esperado = n * p
return ((conteos - esperado) ** 2 / esperado).sum(axis=1)
critico = stats.chi2.ppf(0.95, 3)
sim_grande = chi2_bajo_h0(600, [0.25, 0.25, 0.25, 0.25])
sim_chico = chi2_bajo_h0(20, [0.70, 0.20, 0.08, 0.02])
print("Valor crítico χ²(3) al 5% =", round(critico, 3))
print("Esperados grandes (150 por categoría): P(rechazar | H0) =", round((sim_grande >= critico).mean(), 4))
print("Esperados pequeños (14, 4, 1.6, 0.4): P(rechazar | H0) =", round((sim_chico >= critico).mean(), 4))
print("Valores distintos del estadístico: grandes =", len(np.unique(sim_grande.round(6))),
" pequeños =", len(np.unique(sim_chico.round(6))))
graficar_calibracion(
sim_grande, sim_chico, 3,
"esperados grandes (150 por categoría)", "esperados pequeños (hasta 0.4)",
"Con H0 cierta, ¿el estadístico sigue la curva teórica?",
)Valor crítico χ²(3) al 5% = 7.815
Esperados grandes (150 por categoría): P(rechazar | H0) = 0.0474
Esperados pequeños (14, 4, 1.6, 0.4): P(rechazar | H0) = 0.0663
Valores distintos del estadístico: grandes = 860 pequeños = 198
Con esperados de 150 por categoría, la prueba rechaza en el 4.74% de las simulaciones, muy cerca del 5% nominal, y el histograma sigue la curva. Con esperados de 14, 4, 1.6 y 0.4 rechaza en el 6.63%: la prueba anuncia una tasa de falsas alarmas del 5% y en este escenario comete alrededor de un tercio más.
El panel derecho muestra dos fenómenos a la vez. El estadístico es discreto y avanza a saltos: en 20 000 simulaciones tomó 198 valores distintos, frente a 860 en el otro escenario, y el histograma forma barras altas y separadas. Y la categoría con esperado 0.4 lo dispara: una sola observación allí aporta \((1-0.4)^2/0.4=0.9\), y dos aportan \((2-0.4)^2/0.4=6.4\). Una curva continua describe mal esa cola.
De aquí sale la regla práctica de esperados de al menos 5. Es empírica y no una frontera exacta, porque la aproximación se degrada de forma gradual. Cuando no se cumple tienes cuatro salidas: reunir más datos; juntar categorías con un criterio definido de antemano; calcular el p-value por simulación, generando muchas muestras bajo \(H_0\) como en esta etapa y contando en qué fracción el estadístico iguala o supera al observado; o, en tablas \(2\times 2\), usar una prueba exacta (etapa 7).
En el A/B entre la versión actual (A) y la cuantizada (B), cada request termina en uno de tres resultados: respuesta correcta, respuesta incorrecta o timeout. Asignas 1000 requests al azar a cada versión, y los resultados se simulan con proporciones de diseño de 0.88, 0.07 y 0.05 en A y de 0.86, 0.115 y 0.025 en B.
\[ H_0:\text{el resultado es independiente de la versión},\qquad H_1:\text{la distribución de resultados cambia con la versión},\qquad \alpha=0.05. \]
La tabla tiene 2 filas y 3 columnas, así que hay \((2-1)(3-1)=2\) grados de libertad.
Como el tamaño de cada versión se fijó de antemano (1000 y 1000), el nombre técnico es prueba de homogeneidad: preguntas si las dos versiones tienen la misma distribución de resultados. El cálculo es idéntico al de la prueba de independencia, y por eso suelen presentarse juntas.
resultados = ["correcta", "incorrecta", "timeout"]
n_a = n_b = 1000
conteo_a = rng.multinomial(n_a, [0.88, 0.07, 0.05])
conteo_b = rng.multinomial(n_b, [0.86, 0.115, 0.025])
tabla = pd.DataFrame([conteo_a, conteo_b], index=["A (actual)", "B (cuantizado)"], columns=resultados)
print("Frecuencias observadas")
print(tabla.assign(total=tabla.sum(axis=1)).to_string())
print()
print("Porcentaje por versión")
print((tabla.div(tabla.sum(axis=1), axis=0) * 100).round(1).to_string())
# Esperados a mano: (total fila × total columna) / n
O = tabla.values
n_total = O.sum()
E = np.outer(O.sum(axis=1), O.sum(axis=0)) / n_total
chi2_ind = ((O - E) ** 2 / E).sum()
gl_ind = (O.shape[0] - 1) * (O.shape[1] - 1)
print()
print("Frecuencias esperadas bajo independencia")
print(pd.DataFrame(E, index=tabla.index, columns=resultados).round(2).to_string())
print("Esperado mínimo =", round(E.min(), 2), " (regla: al menos 5)")
print()
print("χ² a mano =", round(chi2_ind, 3), " gl =", gl_ind, " p =", f"{stats.chi2.sf(chi2_ind, gl_ind):.2e}")
chi2_sp, p_sp, gl_sp, E_sp = stats.chi2_contingency(O, correction=False)
print("SciPy chi2_contingency: χ² =", round(chi2_sp, 3), " gl =", gl_sp, " p =", f"{p_sp:.2e}")
print("¿Rechazar H0 a α = 0.05?", p_sp < 0.05)
graficar_perfiles(tabla, "Distribución de resultados por versión (n = 1 000 por versión)")Frecuencias observadas
correcta incorrecta timeout total
A (actual) 874 69 57 1000
B (cuantizado) 880 99 21 1000
Porcentaje por versión
correcta incorrecta timeout
A (actual) 87.4 6.9 5.7
B (cuantizado) 88.0 9.9 2.1
Frecuencias esperadas bajo independencia
correcta incorrecta timeout
A (actual) 877.0 84.0 39.0
B (cuantizado) 877.0 84.0 39.0
Esperado mínimo = 39.0 (regla: al menos 5)
χ² a mano = 21.993 gl = 2 p = 1.68e-05
SciPy chi2_contingency: χ² = 21.993 gl = 2 p = 1.68e-05
¿Rechazar H0 a α = 0.05? True
El esperado mínimo es 39, muy por encima de 5, así que la aproximación es confiable. El estadístico vale 21.99 con 2 grados de libertad y el p-value es \(1.7\times10^{-5}\): rechazas \(H_0\), porque con estos datos la distribución de resultados no es la misma en las dos versiones.
Los porcentajes muestran algo que la prueba global no dice. La proporción de respuestas correctas es parecida (87.4% en A y 88.0% en B); la diferencia está en cómo fallan: B tiene más respuestas incorrectas (9.9% frente a 6.9%) y menos timeouts (2.1% frente a 5.7%). En el panel derecho de la figura, que acerca las categorías poco frecuentes, las barras naranjas de B suben en incorrecta y bajan en timeout. Comparar solo la tasa de respuestas correctas habría escondido ese cambio.
Ver los porcentajes no basta para afirmar qué celdas se apartan de lo esperado más de lo que explica el azar. Eso lo responden los residuos.
La prueba global dice que la tabla no se parece a la esperada, pero no en qué celdas. Dos herramientas completan la lectura: los residuos ajustados, que dicen dónde, y la V de Cramér, que dice cuánto. También calculas el intervalo de la diferencia en la tasa de respuestas correctas, la columna que no se aparta.
p_fila = O.sum(axis=1, keepdims=True) / n_total
p_col = O.sum(axis=0, keepdims=True) / n_total
residuos = pd.DataFrame((O - E) / np.sqrt(E * (1 - p_fila) * (1 - p_col)), index=tabla.index, columns=resultados)
aportes = pd.DataFrame((O - E) ** 2 / E, index=tabla.index, columns=resultados)
print("Residuos ajustados (|r| > 2 señala una celda que se aparta de lo esperado)")
print(residuos.round(2).to_string())
print()
print("Aporte de cada celda al χ²")
print(aportes.round(2).to_string())
print("Aporte por columna:", aportes.sum(axis=0).round(2).to_dict())
m = min(O.shape)
v_cramer = np.sqrt(chi2_ind / (n_total * (m - 1)))
print()
print("V de Cramér =", round(v_cramer, 3))
# La columna que no se aparta: ¿qué tan bien se conoce la diferencia en la tasa de respuestas correctas?
pc_a, pc_b = O[0, 0] / O[0].sum(), O[1, 0] / O[1].sum()
se_pc = np.sqrt(pc_a * (1 - pc_a) / O[0].sum() + pc_b * (1 - pc_b) / O[1].sum())
print("Tasa de correctas, B − A =", round((pc_b - pc_a) * 100, 1), "puntos IC 95%: [",
round((pc_b - pc_a - 1.96 * se_pc) * 100, 1), ",", round((pc_b - pc_a + 1.96 * se_pc) * 100, 1), "]")
graficar_residuos(residuos, "Residuos ajustados: rojo = más de lo esperado, azul = menos")Residuos ajustados (|r| > 2 señala una celda que se aparta de lo esperado)
correcta incorrecta timeout
A (actual) -0.41 -2.42 4.16
B (cuantizado) 0.41 2.42 -4.16
Aporte de cada celda al χ²
correcta incorrecta timeout
A (actual) 0.01 2.68 8.31
B (cuantizado) 0.01 2.68 8.31
Aporte por columna: {'correcta': 0.02, 'incorrecta': 5.36, 'timeout': 16.62}
V de Cramér = 0.105
Tasa de correctas, B − A = 0.6 puntos IC 95%: [ -2.3 , 3.5 ]
Los residuos ajustados de las 6 celdas ordenan las tres columnas:
En la figura, la celda de timeout de A es la más roja y la de B la más azul. En una tabla de dos filas, los residuos de una columna son iguales y de signo opuesto: lo que una versión tiene de más, la otra lo tiene de menos.
La V de Cramér es 0.105. Con las referencias habituales (alrededor de 0.1 pequeño, 0.3 mediano y 0.5 grande) es una asociación pequeña. Esas referencias son orientativas y vienen de otras disciplinas; lo que decide es el costo de cada tipo de fallo. Aquí el p-value es muy pequeño y el efecto también, y las dos cosas son ciertas a la vez porque con 2000 requests se detectan diferencias modestas.
Miras los residuos después de la prueba global y son varios a la vez: sirven para ubicar la discrepancia, no son seis pruebas independientes.
El caso más frecuente en producto tiene dos versiones y un resultado de sí o no. Mides si el usuario aceptó la recomendación del modelo, con 4000 usuarios asignados al azar a cada versión; las aceptaciones se simulan con rng.binomial y tasas de diseño de 11.0% en A y 13.5% en B.
\[ H_0:p_A=p_B,\qquad H_1:p_A\neq p_B,\qquad \alpha=0.05. \]
El estadístico es una \(z\): la diferencia de proporciones dividida por su error estándar. Bajo \(H_0\) las dos versiones comparten una misma proporción, que estimas juntando los datos (proporción combinada, \(\hat p\)):
\[ z=\frac{\hat p_B-\hat p_A}{\sqrt{\hat p\,(1-\hat p)\left(\dfrac1{n_A}+\dfrac1{n_B}\right)}}. \]
Para el intervalo de confianza de la diferencia no supones \(H_0\), así que cada versión usa su propia proporción:
\[ (\hat p_B-\hat p_A)\pm z_{0.975}\sqrt{\frac{\hat p_A(1-\hat p_A)}{n_A}+\frac{\hat p_B(1-\hat p_B)}{n_B}}. \]
n_ua = n_ub = 4000
acepta_a = rng.binomial(n_ua, 0.110)
acepta_b = rng.binomial(n_ub, 0.135)
pa, pb = acepta_a / n_ua, acepta_b / n_ub
dif = pb - pa
# Prueba: error estándar con la proporción combinada (lo que vale bajo H0)
p_comb = (acepta_a + acepta_b) / (n_ua + n_ub)
se_h0 = np.sqrt(p_comb * (1 - p_comb) * (1 / n_ua + 1 / n_ub))
z_obs = dif / se_h0
p_z = 2 * stats.norm.sf(abs(z_obs))
# Intervalo: error estándar sin asumir H0
se_dif = np.sqrt(pa * (1 - pa) / n_ua + pb * (1 - pb) / n_ub)
z975 = stats.norm.ppf(0.975)
ic_dif = (dif - z975 * se_dif, dif + z975 * se_dif)
ic_a = (pa - z975 * np.sqrt(pa * (1 - pa) / n_ua), pa + z975 * np.sqrt(pa * (1 - pa) / n_ua))
ic_b = (pb - z975 * np.sqrt(pb * (1 - pb) / n_ub), pb + z975 * np.sqrt(pb * (1 - pb) / n_ub))
print("A: aceptaron", acepta_a, "de", n_ua, " →", round(pa * 100, 2), "%")
print("B: aceptaron", acepta_b, "de", n_ub, " →", round(pb * 100, 2), "%")
print("Proporción combinada =", round(p_comb, 5))
print("Diferencia B − A =", round(dif * 100, 2), "puntos porcentuales")
print("IC 95% de la diferencia: [", round(ic_dif[0] * 100, 2), ",", round(ic_dif[1] * 100, 2), "] puntos")
print("Cambio relativo =", round(dif / pa * 100, 1), "% sobre la tasa de A")
print("z =", round(z_obs, 3), " p-value =", round(p_z, 4))
tabla_2x2 = np.array([[acepta_a, n_ua - acepta_a], [acepta_b, n_ub - acepta_b]])
chi2_2x2, p_2x2, _, _ = stats.chi2_contingency(tabla_2x2, correction=False)
print()
print("La misma pregunta como tabla 2×2: χ² =", round(chi2_2x2, 3), " p =", round(p_2x2, 4))
print("z² =", round(z_obs ** 2, 3), " ← coincide con χ²")
graficar_dos_proporciones(pa, pb, ic_a, ic_b, dif, ic_dif, "A/B de aceptación de la recomendación (n = 4 000 por versión)")A: aceptaron 472 de 4000 → 11.8 %
B: aceptaron 535 de 4000 → 13.38 %
Proporción combinada = 0.12587
Diferencia B − A = 1.58 puntos porcentuales
IC 95% de la diferencia: [ 0.12 , 3.03 ] puntos
Cambio relativo = 13.3 % sobre la tasa de A
z = 2.123 p-value = 0.0337
La misma pregunta como tabla 2×2: χ² = 4.509 p = 0.0337
z² = 4.509 ← coincide con χ²
La aceptación fue de 11.80% en A (472 de 4000) y de 13.38% en B (535 de 4000). La diferencia es de 1.58 puntos porcentuales, con IC 95% de 0.12 a 3.03 puntos. El estadístico es \(z=2.12\) y el p-value 0.034, así que rechazas la igualdad de proporciones a \(\alpha=0.05\).
Planteada como tabla \(2\times 2\), la misma pregunta da \(\chi^2=4.509\) y el mismo p-value. La igualdad \(z^2=\chi^2\) se cumple siempre que hay 1 grado de libertad, porque una \(\chi^2_1\) es una normal estándar al cuadrado (etapa 1). La prueba \(z\) conserva algo que la \(\chi^2\) pierde: el signo, que dice cuál de las dos versiones es mayor.
Antes del A/B grande se hizo un piloto con 8 usuarios por versión; la tabla está escrita a mano. Es una tabla \(2\times 2\) con conteos pequeños, y la regla de la etapa 3 anticipa el problema: tienes que mirar los esperados.
La prueba exacta de Fisher no usa la aproximación \(\chi^2\). Fija los totales de filas y columnas, enumera todas las tablas posibles con esos totales y suma la probabilidad de las que son tan extremas como la observada o más, calculada con la distribución hipergeométrica. Es un cálculo exacto: vale con cualquier tamaño de muestra, sin exigir esperados grandes.
piloto = pd.DataFrame([[7, 1], [3, 5]], index=["B (cuantizado)", "A (actual)"], columns=["aceptó", "no aceptó"])
print(piloto.to_string())
chi2_p, p_chi, _, E_piloto = stats.chi2_contingency(piloto.values, correction=False)
chi2_y, p_yates, _, _ = stats.chi2_contingency(piloto.values, correction=True)
fisher = stats.fisher_exact(piloto.values)
print()
print("Frecuencias esperadas:")
print(pd.DataFrame(E_piloto, index=piloto.index, columns=piloto.columns).to_string())
print("Esperado mínimo =", E_piloto.min(), " ← por debajo de 5")
print()
print("χ² sin corrección: χ² =", round(chi2_p, 3), " p =", round(p_chi, 4))
print("χ² con corrección de Yates: χ² =", round(chi2_y, 3), " p =", round(p_yates, 4))
print("Prueba exacta de Fisher: p =", round(fisher.pvalue, 4))
print("Prueba exacta de Barnard: p =", round(stats.barnard_exact(piloto.values).pvalue, 4)) aceptó no aceptó
B (cuantizado) 7 1
A (actual) 3 5
Frecuencias esperadas:
aceptó no aceptó
B (cuantizado) 5.0 3.0
A (actual) 5.0 3.0
Esperado mínimo = 3.0 ← por debajo de 5
χ² sin corrección: χ² = 4.267 p = 0.0389
χ² con corrección de Yates: χ² = 2.4 p = 0.1213
Prueba exacta de Fisher: p = 0.1189
Prueba exacta de Barnard: p = 0.0573En el piloto aceptaron 7 de 8 usuarios en B (87.5%) y 3 de 8 en A (37.5%). La diferencia muestral es enorme, pero los esperados son 5 y 3: el supuesto no se cumple, y las pruebas no concuerdan. La \(\chi^2\) sin corrección da \(p=0.039\) y llevaría a rechazar \(H_0\). La prueba exacta de Fisher da \(p=0.119\) y no rechaza. La corrección de Yates, que resta 0.5 al valor absoluto de cada discrepancia antes de elevar al cuadrado, da \(p=0.121\), cerca de la de Fisher.
Con esperados pequeños conviene una prueba exacta, porque no depende de una aproximación. La de Fisher es además conservadora: con tan pocos datos rechaza menos del 5% de las veces aun con \(H_0\) cierta. Otra prueba exacta, la de Barnard, que no fija los totales de columna, da aquí \(p=0.057\). Ninguna de las dos rechaza a \(\alpha=0.05\): con 16 usuarios no puedes descartar que la diferencia se deba al azar. Eso tampoco muestra que no haya diferencia, porque con tan pocos datos la potencia es baja, y el A/B de la etapa 6 sí la detectó con 4000 usuarios por versión.
chi2_contingency aplica la corrección de Yates por defecto en tablas \(2\times 2\). En este capítulo se usa correction=False para que el resultado coincida con la fórmula hecha a mano. Con muestras grandes la corrección casi no cambia nada; con muestras pequeñas, usa Fisher.
El contraste que viste con medias en el capítulo 9 vale igual para conteos. Simulas dos experimentos con el mismo efecto real: la tasa de respuestas incorrectas pasa de 10.0% en A a 10.5% en B, medio punto porcentual.
def experimento_error(n):
a = rng.multinomial(n, [0.900, 0.100])
b = rng.multinomial(n, [0.895, 0.105])
chi2_e, p_e, _, _ = stats.chi2_contingency([a, b], correction=False)
return {
"n por versión": n,
"error A (%)": a[1] / n * 100,
"error B (%)": b[1] / n * 100,
"diferencia (puntos)": (b[1] - a[1]) / n * 100,
"chi2": chi2_e,
"p-value": p_e,
"V de Cramér": np.sqrt(chi2_e / (2 * n)),
}
comparacion = pd.DataFrame([experimento_error(500), experimento_error(500_000)])
print(comparacion.round({"error A (%)": 2, "error B (%)": 2, "diferencia (puntos)": 2, "chi2": 2, "V de Cramér": 4})
.to_string(index=False, formatters={"p-value": "{:.2e}".format})) n por versión error A (%) error B (%) diferencia (puntos) chi2 p-value V de Cramér
500 9.40 11.40 2.00 1.07 3.00e-01 0.0328
500000 9.98 10.57 0.59 93.77 3.54e-22 0.0097En el experimento pequeño la diferencia observada es de 2.00 puntos y el p-value 0.30, así que no rechazas \(H_0\). En el enorme, la diferencia observada es de 0.59 puntos y el p-value \(3.5\times10^{-22}\): el experimento con la diferencia observada más pequeña es el que tiene el p-value más extremo.
La V de Cramér del experimento enorme es 0.0097, una asociación casi nula. El p-value mide qué tan compatible es la diferencia con el azar solo; con un millón de requests el azar explica muy poco, y medio punto se detecta sin dificultad. La V del experimento pequeño (0.033) es mayor solo porque su diferencia observada, 2 puntos, salió alta por azar. Si medio punto de error importa es una pregunta sobre el costo de un error, y la responde el negocio.
El experimento pequeño tampoco contradice al grande. Con 500 requests por versión, medio punto es indetectable, y la diferencia de 2.00 puntos que observaste es compatible con el efecto real de 0.5 y también con ningún efecto.
Para comparar la exactitud de A y B, lo habitual es pasar el mismo conjunto de prueba por los dos modelos. Eso rompe el primer supuesto de la prueba chi-cuadrado: las observaciones ya no son independientes, porque cada caso aparece dos veces. Un caso difícil tiende a fallar en los dos modelos, y uno fácil, a acertar en los dos.
Con datos pareados, la tabla correcta no cruza modelo con acierto. Cruza el resultado de A con el de B, caso por caso:
| B acierta | B falla | |
|---|---|---|
| A acierta | ambos aciertan | \(b\): solo A acierta |
| A falla | \(c\): solo B acierta | ambos fallan |
Los casos en que los dos modelos coinciden no dicen nada sobre cuál es mejor; toda la información está en los discordantes, \(b\) y \(c\). Si los dos modelos fueran igual de buenos, un caso discordante favorecería a A o a B con la misma probabilidad:
\[ H_0: P(\text{solo A acierta})=P(\text{solo B acierta}),\qquad \chi^2_{\text{McNemar}}=\frac{(b-c)^2}{b+c}\ \text{ con 1 grado de libertad}. \]
La versión exacta usa la distribución binomial: bajo \(H_0\), \(b\sim\text{Binomial}(b+c,\,0.5)\). Los conteos de los 400 casos de prueba están escritos a mano.
n_casos = 400
# Conteos diseñados. Cada caso cae en una de cuatro situaciones:
# ambos aciertan, solo A acierta, solo B acierta, ambos fallan
ambos, solo_a, solo_b, ninguno = 318, 28, 12, 42
pareada = pd.DataFrame([[ambos, solo_a], [solo_b, ninguno]],
index=["A acierta", "A falla"], columns=["B acierta", "B falla"])
print("Tabla pareada (cada caso cuenta una sola vez)")
print(pareada.to_string())
acc_a = (ambos + solo_a) / n_casos
acc_b = (ambos + solo_b) / n_casos
print()
print("Exactitud A =", round(acc_a * 100, 2), "% Exactitud B =", round(acc_b * 100, 2), "%")
print("Diferencia A − B =", round((acc_a - acc_b) * 100, 2), "puntos = (b − c)/n =", round((solo_a - solo_b) / n_casos * 100, 2))
se_par = np.sqrt(solo_a + solo_b - (solo_a - solo_b) ** 2 / n_casos) / n_casos # error estándar con datos pareados
print("IC 95% de la diferencia pareada: [", round((acc_a - acc_b - 1.96 * se_par) * 100, 2), ",",
round((acc_a - acc_b + 1.96 * se_par) * 100, 2), "] puntos")
b, c = int(solo_a), int(solo_b)
chi2_mc = (b - c) ** 2 / (b + c)
print()
print("Discordantes: b =", b, " c =", c)
print("McNemar: χ² =", round(chi2_mc, 3), " p =", round(stats.chi2.sf(chi2_mc, 1), 4))
print("McNemar exacta (binomial): p =", round(stats.binomtest(b, b + c, 0.5).pvalue, 4))
# El análisis equivocado: tratar los 800 resultados como si fueran de casos distintos
mal = np.array([[ambos + solo_a, n_casos - ambos - solo_a], [ambos + solo_b, n_casos - ambos - solo_b]])
chi2_mal, p_mal, _, _ = stats.chi2_contingency(mal, correction=False)
print()
print("Chi-cuadrado ignorando el pareo: χ² =", round(chi2_mal, 3), " p =", round(p_mal, 4))Tabla pareada (cada caso cuenta una sola vez)
B acierta B falla
A acierta 318 28
A falla 12 42
Exactitud A = 86.5 % Exactitud B = 82.5 %
Diferencia A − B = 4.0 puntos = (b − c)/n = 4.0
IC 95% de la diferencia pareada: [ 0.93 , 7.07 ] puntos
Discordantes: b = 28 c = 12
McNemar: χ² = 6.4 p = 0.0114
McNemar exacta (binomial): p = 0.0166
Chi-cuadrado ignorando el pareo: χ² = 2.443 p = 0.118La exactitud es 86.5% en A y 82.5% en B: 4 puntos de diferencia sobre los mismos 400 casos.
La diferencia de exactitud es la misma en los dos análisis, porque \((b-c)/n=16/400=4\) puntos; lo que cambia es contra qué variabilidad la comparas. Los resultados de A y B sobre un mismo caso están correlacionados: coinciden en 360 de 400 (318 aciertos comunes y 42 fallos comunes), presumiblemente porque hay casos fáciles y difíciles para ambos. El análisis pareado aprovecha esa correlación; el otro la ignora y sobrestima el error estándar de la diferencia.
En este ejemplo, ignorar el pareo vuelve la prueba menos sensible. No siempre ocurre en esa dirección, pero el p-value calculado siempre deja de corresponder al diseño real. Ten presente también que McNemar compara dos modelos sobre un conjunto de prueba, y no considera la variación que habría al reentrenarlos con otros datos.
El reporte sigue siete pasos, los del capítulo 9 con dos propios de los conteos:
Los aplicas al A/B de resultados de las etapas 4 y 5.
correcta, incorrecta, timeout).En un A/B con 1000 requests por versión, la distribución de resultados difiere entre versiones (\(\chi^2(2)=22.0\), \(p<0.001\); esperado mínimo 39; V de Cramér de 0.11). No se detecta diferencia en la proporción de respuestas correctas (87.4% en A, 88.0% en B). La versión cuantizada tiene menos timeouts (2.1% frente a 5.7%) y más respuestas incorrectas (9.9% frente a 6.9%). La recomendación depende del costo relativo de un timeout y de una respuesta incorrecta, que estos datos no miden.
Con un ejemplo pequeño y números redondos puedes hacer todo a mano y dejar el código como verificación. Quieres saber si la aceptación de la recomendación se relaciona con el tipo de dispositivo, en una muestra escrita a mano de 120 sesiones de usuarios distintos, 60 desde móvil y 60 desde escritorio:
| aceptó | no aceptó | total | |
|---|---|---|---|
| móvil | 18 | 42 | 60 |
| escritorio | 30 | 30 | 60 |
| total | 48 | 72 | 120 |
\[ E_{\text{móvil, aceptó}}=\frac{60\times 48}{120}=24,\qquad E_{\text{móvil, no aceptó}}=\frac{60\times 72}{120}=36, \]
y lo mismo para escritorio, porque las dos filas tienen 60 sesiones. El esperado mínimo es 24, así que el supuesto se cumple.
| celda | \(O\) | \(E\) | \(O-E\) | \((O-E)^2/E\) |
|---|---|---|---|---|
| móvil, aceptó | 18 | 24 | \(-6\) | \(36/24=1.5\) |
| móvil, no aceptó | 42 | 36 | \(+6\) | \(36/36=1.0\) |
| escritorio, aceptó | 30 | 24 | \(+6\) | \(36/24=1.5\) |
| escritorio, no aceptó | 30 | 36 | \(-6\) | \(36/36=1.0\) |
\[ \chi^2=1.5+1.0+1.5+1.0=5.0 . \]
\[ V=\sqrt{\frac{\chi^2}{n\,(m-1)}}=\sqrt{\frac{5.0}{120\times 1}}=0.204 . \]
La celda siguiente verifica cada número.
ejemplo = pd.DataFrame([[18, 42], [30, 30]], index=["móvil", "escritorio"], columns=["aceptó", "no aceptó"])
O_ej = ejemplo.values
n_ej = O_ej.sum()
E_ej = np.outer(O_ej.sum(axis=1), O_ej.sum(axis=0)) / n_ej
print("Esperadas:")
print(pd.DataFrame(E_ej, index=ejemplo.index, columns=ejemplo.columns).to_string())
chi2_ej = ((O_ej - E_ej) ** 2 / E_ej).sum()
print()
print("χ² a mano =", round(chi2_ej, 3), " valor crítico χ²(1) al 5% =", round(stats.chi2.ppf(0.95, 1), 3))
print("p-value =", round(stats.chi2.sf(chi2_ej, 1), 4))
print("V de Cramér =", round(np.sqrt(chi2_ej / n_ej), 3))
sp_sin = stats.chi2_contingency(O_ej, correction=False)
sp_con = stats.chi2_contingency(O_ej, correction=True)
print()
print("SciPy sin corrección: χ² =", round(sp_sin[0], 3), " p =", round(sp_sin[1], 4))
print("SciPy con Yates: χ² =", round(sp_con[0], 3), " p =", round(sp_con[1], 4))
print("Fisher exacta: p =", round(stats.fisher_exact(O_ej).pvalue, 4))
p_mov, p_esc = 18 / 60, 30 / 60
se_ej = np.sqrt(p_mov * (1 - p_mov) / 60 + p_esc * (1 - p_esc) / 60)
print()
print("Diferencia escritorio − móvil =", round((p_esc - p_mov) * 100, 1), "puntos",
" IC 95%: [", round((p_esc - p_mov - 1.96 * se_ej) * 100, 1), ",", round((p_esc - p_mov + 1.96 * se_ej) * 100, 1), "]")Esperadas:
aceptó no aceptó
móvil 24.0 36.0
escritorio 24.0 36.0
χ² a mano = 5.0 valor crítico χ²(1) al 5% = 3.841
p-value = 0.0253
V de Cramér = 0.204
SciPy sin corrección: χ² = 5.0 p = 0.0253
SciPy con Yates: χ² = 4.201 p = 0.0404
Fisher exacta: p = 0.0399
Diferencia escritorio − móvil = 20.0 puntos IC 95%: [ 2.8 , 37.2 ]El código reproduce cada número del cálculo a mano: esperados de 24 y 36, \(\chi^2=5.0\) y \(V=0.204\). El p-value es 0.025: si el dispositivo no se relacionara con la aceptación, una tabla tan alejada de la esperada o más aparecería en alrededor de 25 de cada 1000 muestras.
En 120 sesiones (60 por dispositivo), la aceptación fue de 30% en móvil y 50% en escritorio. La diferencia es de 20 puntos porcentuales (IC 95%: 2.8 a 37.2). \(\chi^2(1)=5.0\), \(p=0.025\), esperado mínimo 24, V de Cramér de 0.20. Los datos son observacionales: muestran asociación entre dispositivo y aceptación, sin establecer la causa.
Con la misma tabla, calcula a mano la prueba \(z\) de dos proporciones de la etapa 6 (con la proporción combinada) y comprueba que \(z^2\) coincide con el \(\chi^2\) de 5.0.
La proporción combinada es \(\hat p=48/120=0.40\), y el error estándar bajo \(H_0\) es \(\sqrt{0.40\times 0.60\times(1/60+1/60)}=\sqrt{0.008}\approx 0.0894\). La diferencia es \(0.50-0.30=0.20\), así que \(z=0.20/0.0894\approx 2.236\) y \(z^2=5.0\), igual que el \(\chi^2\). El p-value bilateral de \(z=2.236\) es 0.025, el mismo de la celda.
Un clasificador se entrena con una distribución de clases, y en producción esa distribución puede cambiar: aparece un producto nuevo, cambia la temporada, cambia el tipo de usuario. Ese cambio se llama deriva (drift) y puede degradar el modelo sin que se note a tiempo, porque en producción casi nunca tienes pronto la etiqueta verdadera.
Lo que sí tienes es lo que el modelo predice, y vigilar la distribución de las clases predichas es una prueba de bondad de ajuste. Las proporciones de referencia son las de las predicciones del modelo sobre un conjunto de validación amplio y reciente; no se usan las etiquetas del entrenamiento, porque un modelo imperfecto no reproduce exactamente su distribución. Los conteos observados son las predicciones de cada semana.
Un clasificador de intención ficticio asigna cada mensaje de soporte a una de cinco clases. En validación, sus predicciones se repartieron en proporciones 0.40, 0.25, 0.15, 0.12 y 0.08. Revisas cuatro semanas de producción, con 1500 mensajes cada una. Las tres primeras se simulan con esa distribución de referencia. En la cuarta, por diseño, queja sube de 8% a 12% y el resto se reacomoda: consulta baja de 40% a 36%, pago de 25% a 24% y devolución sube de 12% a 13%. Simulas directamente los conteos de predicciones con rng.multinomial; no hay un modelo entrenado detrás.
\[ H_0:\text{las predicciones de la semana siguen las proporciones de referencia}. \]
clases = ["consulta", "pago", "envío", "devolución", "queja"]
p_ref = np.array([0.40, 0.25, 0.15, 0.12, 0.08])
p_semana4 = np.array([0.36, 0.24, 0.15, 0.13, 0.12])
n_sem = 1500
filas, proporciones = [], []
for semana, p_real in enumerate([p_ref, p_ref, p_ref, p_semana4], start=1):
conteos = rng.multinomial(n_sem, p_real)
esperado = n_sem * p_ref
prueba = stats.chisquare(conteos, f_exp=esperado)
resid = (conteos - esperado) / np.sqrt(esperado * (1 - p_ref)) # residuo ajustado, como en la sección 1.9
proporciones.append(conteos / n_sem)
filas.append({
"semana": semana,
**dict(zip(clases, conteos)),
"chi2": round(prueba.statistic, 2),
"p-value": prueba.pvalue,
"w": round(np.sqrt(prueba.statistic / n_sem), 3),
"mayor residuo": f"{clases[np.abs(resid).argmax()]} ({resid[np.abs(resid).argmax()]:+.2f})",
})
monitoreo = pd.DataFrame(filas)
print("Esperado por semana:", dict(zip(clases, (n_sem * p_ref).astype(int).tolist())))
print()
print(monitoreo.to_string(index=False, formatters={"p-value": lambda p: f"{p:.4f}" if p >= 0.0001 else f"{p:.1e}"}))
print()
print("Residuos ajustados de la semana 4:", dict(zip(clases, resid.round(2).tolist())))
print("Semanas con p < 0.05: ", monitoreo.loc[monitoreo["p-value"] < 0.05, "semana"].tolist())
print("Semanas con p < 0.05/4 = 0.0125:", monitoreo.loc[monitoreo["p-value"] < 0.05 / 4, "semana"].tolist())
graficar_deriva(clases, p_ref, proporciones, monitoreo["p-value"].values,
"Vigilancia de la distribución de clases predichas (n = 1 500 por semana)")Esperado por semana: {'consulta': 600, 'pago': 375, 'envío': 225, 'devolución': 180, 'queja': 120}
semana consulta pago envío devolución queja chi2 p-value w mayor residuo
1 593 367 219 187 134 2.32 0.6775 0.039 queja (+1.33)
2 613 361 202 210 114 8.46 0.0762 0.075 devolución (+2.38)
3 581 402 233 176 108 4.12 0.3901 0.052 pago (+1.61)
4 546 362 203 217 172 37.60 1.4e-07 0.158 queja (+4.95)
Residuos ajustados de la semana 4: {'consulta': -2.85, 'pago': -0.78, 'envío': -1.59, 'devolución': 2.94, 'queja': 4.95}
Semanas con p < 0.05: [4]
Semanas con p < 0.05/4 = 0.0125: [4]
En las tres primeras semanas, generadas sin deriva, los p-values son 0.68, 0.076 y 0.39: ninguna alarma. En la cuarta, \(\chi^2=37.6\) y \(p=1.4\times10^{-7}\). El mayor residuo corresponde a queja, con 172 predicciones cuando esperabas 120 (11.5% frente a 8%; residuo de +4.95). También se apartan devolución (+2.94) y consulta (−2.85), que cambiaron por diseño, pero la prueba no distingue cuál es el cambio de origen y cuál el reacomodo. La columna w resume la magnitud: de 0.04 a 0.08 en las semanas sin deriva y 0.16 en la cuarta. En el panel derecho de la figura, el p-value de la semana 4 cae por debajo de \(10^{-6}\), muy lejos de las dos líneas de umbral.
La semana 2 merece atención. Su p-value es 0.076 y un residuo llega a +2.38 (devolución), aunque se generó sin deriva. Si miraras solo ese residuo, verías una señal donde no la hay: los residuos son varios, y alguno pasa de 2 por azar.
Varias semanas son varias pruebas. Con \(\alpha=0.05\) por semana y sin deriva alguna, la probabilidad de al menos una falsa alarma en cuatro semanas es \(1-0.95^4=0.185\). Un ajuste simple es exigir \(p<0.05/4=0.0125\), la línea punteada de la figura, y la semana 4 lo pasa con holgura; el capítulo 11 desarrolla este problema.
Este uso de la prueba tiene tres límites:
El ejemplo usa la misma prueba de la etapa 2: las réplicas del balanceador son ahora las clases del modelo, y el reparto uniforme, la distribución de referencia.
chi2_contingency aplica la corrección de Yates salvo que indiques correction=False; di cuál usaste.| Situación | Qué usar | En SciPy |
|---|---|---|
| Una variable categórica frente a proporciones conocidas | Bondad de ajuste, \(k-1\) gl | stats.chisquare(obs, f_exp=esp) |
| Dos variables categóricas, muestras independientes | Independencia, \((f-1)(c-1)\) gl | stats.chi2_contingency(tabla) (en \(2\times 2\) aplica Yates salvo correction=False) |
| Dos grupos, respuesta sí o no | \(z\) de dos proporciones o \(\chi^2\) en tabla \(2\times 2\) | stats.chi2_contingency(tabla, correction=False) |
| Tabla \(2\times 2\) con esperados menores que 5 | Prueba exacta de Fisher | stats.fisher_exact(tabla) |
| Los mismos casos medidos dos veces, respuesta sí o no | McNemar sobre los discordantes | stats.binomtest(b, b + c, 0.5) |
| ¿Qué celdas se apartan? | Residuos ajustados | a mano (etapa 5) |
| ¿Qué tan fuerte es la asociación? | V de Cramér, diferencia de proporciones con IC | a mano (etapas 5 y 6) |
| La respuesta es un número, no una categoría | Pruebas del capítulo 11 | — |
Al reportar, sigue este orden: tabla observada, esperado mínimo, estadístico con sus grados de libertad, p-value y tamaño del efecto.
binomtest(b, b + c, 0.5). La diferencia \(b-c\) es necesariamente \(150-138=12\), pero la significancia depende de cuántos casos discordantes hay.scipy.stats.chisquare, chi2_contingency, fisher_exact y binomtest: https://docs.scipy.org/doc/scipy/reference/stats.html.Una prueba con conteos compara lo observado con lo que esperarías bajo \(H_0\): el estadístico \(\chi^2\) mide esa distancia, la distribución chi-cuadrado dice cuánta distancia produce el azar, los residuos ubican la discrepancia y el tamaño del efecto dice si importa.
Este capítulo cubrió la rama del mapa en la que la respuesta es una categoría. El capítulo 11 recorre la otra rama, en la que la respuesta es un número: qué hacer con mediciones pareadas, con más de dos grupos (ANOVA), con distribuciones muy asimétricas (pruebas no paramétricas y de permutación) y cuando haces muchas pruebas a la vez, un problema que ya asomó en la etapa 12, donde cuatro semanas eran cuatro oportunidades de una falsa alarma.
La distribución chi-cuadrado vuelve en el nivel Avanzado, en los modelos lineales generalizados, donde sirve para comparar modelos y medir qué tan bien se ajustan a los datos. Y las tablas de conteos reaparecen en el capítulo 14 como matrices de confusión.