Proyecto final e IA responsable

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_openml

# La primera vez se descarga (unos 4 MB); después sale de la caché
adult = fetch_openml("adult", version=2, as_frame=True).frame

# Nombres en español; se dejan fuera fnlwgt (un peso de muestreo
# del censo) y education (repite anios_estudio en palabras)
datos = adult.rename(columns={
    "age": "edad", "workclass": "empleador",
    "education-num": "anios_estudio", "marital-status": "estado_civil",
    "occupation": "ocupacion", "relationship": "parentesco",
    "capital-gain": "ganancia_capital", "capital-loss": "perdida_capital",
    "hours-per-week": "horas", "native-country": "pais", "race": "raza",
}).drop(columns=["fnlwgt", "education"])
sexos = {"Female": "mujer", "Male": "hombre"}
datos["sexo"] = datos["sex"].astype(str).map(sexos)
datos["ingreso_alto"] = (datos["class"] == ">50K").astype(int)
datos = datos.drop(columns=["sex", "class"])

print(datos.shape)
print(f"Ganan más de 50 000 dólares: {datos['ingreso_alto'].mean():.1%}")
# tasa_base: la fracción del grupo que realmente es positiva
print(datos.groupby("sexo")["ingreso_alto"]
      .agg(personas="size", tasa_base="mean").round(3))

from sklearn.model_selection import train_test_split

NUMERICAS = ["edad", "anios_estudio", "ganancia_capital",
             "perdida_capital", "horas"]
CATEGORICAS = ["empleador", "estado_civil", "ocupacion", "parentesco",
               "pais"]
X = datos[NUMERICAS + CATEGORICAS]
y = datos["ingreso_alto"]
sexo = datos["sexo"]

# Se estratifica por ingreso y sexo a la vez: las tres partes tienen
# las mismas proporciones de cada combinación
estratos = y.astype(str) + sexo
X_ent, X_resto, y_ent, y_resto, s_ent, s_resto = train_test_split(
    X, y, sexo, test_size=0.4, stratify=estratos, random_state=0)
X_val, X_pru, y_val, y_pru, s_val, s_pru = train_test_split(
    X_resto, y_resto, s_resto, test_size=0.5,
    stratify=estratos[X_resto.index], random_state=0)
print(f"Entrenamiento: {len(X_ent)}   validación: {len(X_val)}"
      f"   prueba: {len(X_pru)}")

# ¿Quitar la columna sexo basta para que el modelo no la vea?
print(pd.crosstab(datos["parentesco"], datos["sexo"]))

from sklearn.compose import ColumnTransformer
from sklearn.dummy import DummyClassifier
from sklearn.ensemble import HistGradientBoostingClassifier
from sklearn.impute import SimpleImputer
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score, roc_auc_score
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import OneHotEncoder, StandardScaler

def crear_pipeline(modelo, numericas=NUMERICAS, categoricas=CATEGORICAS):
    """Prepara las columnas de la tabla y termina en el modelo dado."""
    prep_categoricas = Pipeline([
        ("imputar", SimpleImputer(strategy="most_frequent")),
        ("onehot", OneHotEncoder(handle_unknown="ignore",
                                 sparse_output=False))])
    preparar = ColumnTransformer([
        ("num", StandardScaler(), numericas),
        ("cat", prep_categoricas, categoricas)])
    return Pipeline([("preparar", preparar), ("modelo", modelo)])

candidatos = {
    "línea base": DummyClassifier(strategy="most_frequent"),
    "logística": LogisticRegression(max_iter=1000),
    "boosting": HistGradientBoostingClassifier(random_state=0),
}
exactitudes = {}
for nombre, estimador in candidatos.items():
    pipe = crear_pipeline(estimador).fit(X_ent, y_ent)
    prob = pipe.predict_proba(X_val)[:, 1]
    exactitudes[nombre] = accuracy_score(y_val, prob >= 0.5)
    print(f"{nombre:10s}  exactitud = {exactitudes[nombre]:.3f}"
          f"   AUC = {roc_auc_score(y_val, prob):.3f}")
modelo = pipe  # el último: el boosting

nombres = list(exactitudes)
colores = ["C1" if n == "boosting" else "C0" for n in nombres]
plt.barh(nombres, list(exactitudes.values()), color=colores)
for i, v in enumerate(exactitudes.values()):
    plt.text(v - 0.01, i, f"{v:.3f}", ha="right", va="center",
             color="white")
plt.xlim(0, 1)
plt.xlabel("Exactitud en validación")
plt.title("Primero, ganarle a la línea base")
plt.show()

prob_pru = modelo.predict_proba(X_pru)[:, 1]
pred_pru = (prob_pru >= 0.5).astype(int)
print(f"Exactitud en prueba: {accuracy_score(y_pru, pred_pru):.3f}")

def metricas(tabla):
    """Métricas de un grupo. tabla tiene las columnas y y pred."""
    real, pred = tabla["y"], tabla["pred"]
    return pd.Series({
        "n": len(tabla),
        "tasa_base": real.mean(),
        "aprobacion": pred.mean(),
        "tvp": pred[real == 1].mean(),
        "tfp": pred[real == 0].mean(),
        "precision": real[pred == 1].mean(),
        "exactitud": (real == pred).mean(),
    })

auditoria = pd.DataFrame({"y": y_pru.to_numpy(), "pred": pred_pru,
                          "sexo": s_pru.to_numpy()})
por_sexo = auditoria.groupby("sexo")[["y", "pred"]].apply(metricas)
por_sexo["n"] = por_sexo["n"].astype(int)
print(por_sexo.round(3).to_string())

rng = np.random.default_rng(0)

def intervalos(tabla, repeticiones=1000):
    """Intervalo del 95 % de cada métrica, por bootstrap."""
    muestras = []
    for _ in range(repeticiones):
        filas = rng.integers(0, len(tabla), size=len(tabla))
        muestras.append(metricas(tabla.iloc[filas]))
    muestras = pd.DataFrame(muestras)
    return muestras.quantile([0.025, 0.975])

MOSTRAR = ["aprobacion", "tvp", "exactitud"]
grupos = {g: t for g, t in auditoria.groupby("sexo")}
for g, tabla in grupos.items():
    iv = intervalos(tabla)
    texto = "   ".join(f"{m}: {iv.loc[0.025, m]:.3f}–"
                       f"{iv.loc[0.975, m]:.3f}" for m in MOSTRAR)
    print(f"{g:7s} {texto}")
    grupos[g] = iv

ancho = 0.38
posiciones = np.arange(len(MOSTRAR))
estilos = [("hombre", "hombres", "C2", ""),
           ("mujer", "mujeres", "C3", "//")]
for k, (g, rotulo, color, rayas) in enumerate(estilos):
    valores = por_sexo.loc[g, MOSTRAR].to_numpy(dtype=float)
    iv = grupos[g][MOSTRAR].to_numpy()
    errores = [valores - iv[0], iv[1] - valores]
    plt.bar(posiciones + (k - 0.5) * ancho, valores, ancho, yerr=errores,
            capsize=4, color=color, edgecolor="C2", hatch=rayas,
            ecolor="C0", label=rotulo)
plt.xticks(posiciones, ["Aprobación", "TVP", "Exactitud"])
plt.ylim(0, 1)
plt.ylabel("Tasa en el conjunto de prueba")
plt.title("Métricas por grupo")
plt.legend(loc="upper left")
plt.show()

prueba = X_pru.assign(sexo=s_pru, real=y_pru, prob=prob_pru.round(2))
errores_fn = prueba[(prueba["real"] == 1) & (prueba["prob"] < 0.5)]
errores_fp = prueba[(prueba["real"] == 0) & (prueba["prob"] >= 0.5)]
print("Falsos negativos por grupo:",
      errores_fn["sexo"].value_counts().to_dict())
print("Falsos positivos por grupo:",
      errores_fp["sexo"].value_counts().to_dict())

columnas = ["sexo", "edad", "anios_estudio", "estado_civil", "horas",
            "prob"]
# Tres falsos negativos al azar de cada grupo
muestra = errores_fn.groupby("sexo").sample(3, random_state=0)
print(muestra[columnas].to_string())

from sklearn.inspection import permutation_importance

imp = permutation_importance(modelo, X_pru, y_pru, n_repeats=5,
                             random_state=0)
orden = np.argsort(imp.importances_mean)
for j in orden[::-1][:5]:
    print(f"{X.columns[j]:17s} {imp.importances_mean[j]:+.3f}")

colores = ["C1" if j == orden[-1] else "C0" for j in orden]
plt.barh(X.columns[orden], imp.importances_mean[orden], color=colores)
plt.xlabel("Exactitud que se pierde")
plt.title("Importancia por permutación")
plt.show()

prob_val = modelo.predict_proba(X_val)[:, 1]
UMBRALES = np.round(np.arange(0.05, 0.96, 0.01), 2)

def tasas(real, prob, umbral):
    """Tasas de aprobación, de verdaderos positivos (TVP) y de falsos
    positivos (TFP)."""
    pred = prob >= umbral
    return pred.mean(), pred[real == 1].mean(), pred[real == 0].mean()

def umbral_parecido(real, prob, objetivo, cual):
    """El umbral cuya tasa (0: aprobación, 1: TVP) queda más cerca
    del objetivo."""
    distancias = [abs(tasas(real, prob, u)[cual] - objetivo)
                  for u in UMBRALES]
    return UMBRALES[int(np.argmin(distancias))]

# En validación: los hombres quedan con 0,5 y se busca el de las mujeres
h, m = (s_val == "hombre").to_numpy(), (s_val == "mujer").to_numpy()
yv = y_val.to_numpy()
meta = tasas(yv[h], prob_val[h], 0.5)
politicas = {
    "un solo umbral": 0.5,
    "igualar aprobación": umbral_parecido(yv[m], prob_val[m], meta[0], 0),
    "igualar oportunidades": umbral_parecido(yv[m], prob_val[m],
                                             meta[1], 1),
}

# En prueba: cada política se mide una sola vez
hp, mp = (s_pru == "hombre").to_numpy(), (s_pru == "mujer").to_numpy()
yp = y_pru.to_numpy()
print("política               umbral M  aprob H/M  TVP H/M    TFP H/M"
      "    exactitud")
for nombre, u_m in politicas.items():
    umbral = np.where(mp, u_m, 0.5)
    pred = prob_pru >= umbral
    ah, th, fh = tasas(yp[hp], prob_pru[hp], 0.5)
    am, tm, fm = tasas(yp[mp], prob_pru[mp], u_m)
    print(f"{nombre:22s} {u_m:.2f}      {ah:.2f}/{am:.2f}  "
          f"{th:.2f}/{tm:.2f}  {fh:.2f}/{fm:.2f}  "
          f"{(pred == yp).mean():.3f}")

# Todas las parejas de umbrales (hombres, mujeres), en validación
rejilla = np.round(np.arange(0.05, 0.96, 0.01), 2)
puntos = []
for u_h in rejilla:
    for u_m in rejilla:
        pred = prob_val >= np.where(m, u_m, u_h)
        dif_aprob = abs(pred[h].mean() - pred[m].mean())
        dif_tvp = abs(pred[h & (yv == 1)].mean()
                      - pred[m & (yv == 1)].mean())
        puntos.append((dif_aprob, dif_tvp, (pred == yv).mean()))
puntos = np.array(puntos)

def frontera(diferencias, exactitud):
    """Para cada diferencia, la mejor exactitud que no la supera."""
    orden = np.argsort(diferencias)
    return diferencias[orden], np.maximum.accumulate(exactitud[orden])

for col, nombre in [(0, "aprobación"), (1, "TVP")]:
    for tope in [0.10, 0.05, 0.02]:
        ok = puntos[:, col] <= tope
        print(f"Diferencia de {nombre} <= {tope:.2f}: mejor exactitud "
              f"{puntos[ok, 2].max():.3f}")

x_a, y_a = frontera(puntos[:, 0], puntos[:, 2])
x_t, y_t = frontera(puntos[:, 1], puntos[:, 2])
plt.step(x_a, y_a, where="post", color="C0", label="aprobación")
plt.step(x_t, y_t, where="post", color="C2", linestyle="--",
         label="oportunidades (TVP)")
ok = puntos[:, 0] <= 0.02
plt.scatter([0.02], [puntos[ok, 2].max()], s=80, color="C1", zorder=3)
plt.xlim(0, 0.3)
plt.ylim(0.84, 0.88)
plt.xlabel("Diferencia máxima entre grupos")
plt.ylabel("Mejor exactitud (validación)")
plt.title("Cuánto cuesta la equidad")
plt.legend(loc="lower right")
plt.show()