3.5 - Notas de Luis - Regresión logística con scikit-learn¶
En el cuaderno Regresión logística desde cero predijimos si un paciente de Fisioterapia faltará a su cita escribiendo nosotros la sigmoide, el coste y el descenso del gradiente. Aquí resolvemos el mismo problema con scikit-learn y añadimos las dos herramientas de este capítulo sin las que un proyecto real no puede funcionar:
- La regularización, que impide que el modelo memorice los datos de entrenamiento (sobreajuste).
- La validación cruzada de k grupos, que permite elegir la intensidad de la regularización sin tocar el conjunto de prueba.
Pasos:
- Importar las bibliotecas y crear los datos.
- Entrenar una regresión logística con scikit-learn.
- Evaluarla: matriz de confusión, informe de clasificación y curva ROC.
- El sobreajuste: qué ocurre con pocos pacientes y muchas variables.
- Regularización L2 y L1.
- Elegir $\lambda$ con validación cruzada de k grupos.
- Elastic Net y variables correlacionadas.
- Clases desbalanceadas.
Paso 1: Bibliotecas y datos¶
Usamos las mismas citas sintéticas de Fisioterapia que en el cuaderno desde cero (mismo generador y misma semilla, para poder comparar los números):
| Columna | Significado |
|---|---|
edad |
Edad del paciente (18–85 años) |
dias_desde_ultima_cita |
Días desde la cita anterior del paciente |
distancia_km |
Distancia de casa al hospital, en km (Minermont atiende una zona rural) |
hora |
Hora de la cita (de 8 a 19 h) |
primera_visita |
1 si es una primera visita, 0 si es una revisión |
ausencias_previas |
Número de citas a las que el paciente faltó anteriormente |
no_asiste |
Objetivo: 1 si el paciente no acudió, 0 si acudió |
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split, StratifiedKFold, cross_validate
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import make_pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import (accuracy_score, confusion_matrix, ConfusionMatrixDisplay,
classification_report, roc_auc_score, roc_curve, RocCurveDisplay)
np.set_printoptions(precision=3, suppress=True)
def generate_appointments(n, seed):
"""Citas sintéticas de Fisioterapia en Minermont: una fila por cita."""
rng = np.random.default_rng(seed)
age = rng.integers(18, 86, n) # años
days = rng.integers(1, 181, n) # días desde la cita anterior
distance = np.round(rng.gamma(2.0, 6.0, n).clip(0.5, 60), 1) # km hasta el hospital (zona rural)
hour = rng.integers(8, 20, n) # hora de la cita, de 8:00 a 19:00
first = rng.binomial(1, 0.35, n) # 1 = primera visita, 0 = revisión
previous = rng.poisson(0.6, n).clip(0, 6) # ausencias en el historial del paciente
late = np.where(hour >= 17, hour - 16, 0) # a última hora de la tarde se falta más
# La regla oculta detrás de los datos. Solo la usamos para generar las citas:
# el modelo nunca la ve y tiene que descubrir el patrón a partir de los ejemplos.
z = (-1.55 - 0.069 * (age - 50) + 0.026 * (days - 90) + 0.11 * (distance - 12)
+ 1.5 * (previous - 0.6) + 0.56 * (late - 0.75) + 0.63 * (first - 0.35))
no_show = rng.binomial(1, 1 / (1 + np.exp(-z))) # 1 = no acudió
return pd.DataFrame({
'edad': age, 'dias_desde_ultima_cita': days, 'distancia_km': distance, 'hora': hour,
'primera_visita': first, 'ausencias_previas': previous, 'no_asiste': no_show,
})
df = generate_appointments(1200, seed=42)
feature_names = ['edad', 'dias_desde_ultima_cita', 'distancia_km', 'hora', 'primera_visita', 'ausencias_previas']
X = df[feature_names]
y = df['no_asiste']
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=0.25, random_state=42, stratify=y)
print(f"Entrenamiento: {len(y_train)} citas, prueba: {len(y_test)}")
print(f"Tasa de ausencias: {y.mean():.1%}")
Entrenamiento: 900 citas, prueba: 300 Tasa de ausencias: 26.2%
Paso 2: Entrenamiento con scikit-learn¶
Un pipeline encadena el escalado y el modelo en un único objeto. Cuando llamamos a fit, el escalador aprende la media y la desviación típica solo con los datos de entrenamiento, y cuando llamamos a predict, las aplica a los datos nuevos. Así se evita, por construcción, el error de dejar que el conjunto de prueba se filtre en el entrenamiento.
LogisticRegression no usa el descenso del gradiente simple, sino un optimizador más rápido (L-BFGS), y aplica por defecto una regularización L2 ligera (C=1). Con 900 citas y seis variables, sus pesos salen muy parecidos a los que obtuvimos desde cero, y algo más pequeños por esa regularización ligera.
model = make_pipeline(StandardScaler(), LogisticRegression())
model.fit(X_train, y_train)
logreg = model[-1]
print("Pesos:", pd.Series(logreg.coef_[0], index=feature_names).round(3).to_dict())
print(f"Sesgo: {logreg.intercept_[0]:.3f}")
Pesos: {'edad': -1.235, 'dias_desde_ultima_cita': 1.136, 'distancia_km': 0.992, 'hora': 0.354, 'primera_visita': 0.28, 'ausencias_previas': 1.0}
Sesgo: -1.776
Paso 3: Evaluación¶
predict aplica el umbral de 0,5; predict_proba devuelve las probabilidades, que necesitamos para probar otros umbrales. El informe de clasificación muestra, para cada clase, la precision (lo que llamamos VPP para la clase de ausencias), el recall (la sensibilidad) y el F1-score, que combina ambas.
Atención a la terminología en castellano: accuracy se traduce como exactitud (la proporción total de aciertos), mientras que precision es la precisión en el sentido del VPP. Son métricas distintas.
y_pred = model.predict(X_test)
p_test = model.predict_proba(X_test)[:, 1]
print(f"Exactitud: {accuracy_score(y_test, y_pred):.3f} (predecir siempre 'acude': {1 - y_test.mean():.3f})\n")
print(classification_report(y_test, y_pred, target_names=['acude', 'no acude'], digits=3))
ConfusionMatrixDisplay.from_predictions(y_test, y_pred, display_labels=['acude', 'no acude'], cmap='Blues')
plt.xlabel('Predicción')
plt.ylabel('Valor real')
plt.title('Matriz de confusión (conjunto de prueba, umbral 0,5)')
plt.show()
Exactitud: 0.827 (predecir siempre 'acude': 0.737)
precision recall f1-score support
acude 0.869 0.900 0.884 221
no acude 0.690 0.620 0.653 79
accuracy 0.827 300
macro avg 0.780 0.760 0.769 300
weighted avg 0.822 0.827 0.824 300
La curva ROC muestra todos los umbrales a la vez: para cada uno representa la sensibilidad (tasa de verdaderos positivos) frente a 1 − especificidad (tasa de falsos positivos). Un modelo inútil sigue la diagonal; cuanto más se acerca la curva a la esquina superior izquierda, mejor. El área bajo la curva (AUC) la resume en un solo número: 0,5 es el azar y 1 es perfecto. También puede leerse como la probabilidad de que el modelo dé una puntuación más alta a una ausencia elegida al azar que a una asistencia elegida al azar.
print(f"AUC: {roc_auc_score(y_test, p_test):.3f}")
fpr, tpr, thresholds = roc_curve(y_test, p_test)
RocCurveDisplay(fpr=fpr, tpr=tpr, roc_auc=roc_auc_score(y_test, p_test)).plot()
for t, marker in [(0.5, 'o'), (0.3, 's')]:
i = np.argmin(np.abs(thresholds - t))
plt.plot(fpr[i], tpr[i], marker, markersize=10, label=f'umbral {t}')
plt.plot([0, 1], [0, 1], 'k--', label='azar')
plt.xlabel('Tasa de falsos positivos (1 - especificidad)')
plt.ylabel('Tasa de verdaderos positivos (sensibilidad)')
plt.title('Curva ROC')
plt.legend()
plt.show()
AUC: 0.880
Pasar del umbral de 0,5 al de 0,3 es subir por la curva: más ausencias detectadas, a cambio de más falsas alarmas.
Paso 4: El sobreajuste¶
En el libro, el primer modelo del hospital era excelente con los datos de entrenamiento y mucho peor con pacientes nuevos: había aprendido las particularidades de los ejemplos de entrenamiento, no el patrón general. Es fácil reproducirlo. Tomamos una muestra pequeña de 300 citas y le añadimos 50 columnas extra que no tienen nada que ver con la asistencia (piensa en todos los demás campos que contiene una historia clínica). La complejidad del modelo crece con el número de variables y, con pocos ejemplos, puede memorizarlas.
Para ver el sobreajuste en estado puro entrenamos con una regularización casi nula (C=10000).
small = generate_appointments(300, seed=3)
X_small = small[feature_names].copy()
rng = np.random.default_rng(7)
for j in range(50):
X_small[f'extra_{j}'] = rng.normal(size=len(small)) # variables irrelevantes
y_small = small['no_asiste']
Xs_train, Xs_test, ys_train, ys_test = train_test_split(
X_small, y_small, test_size=0.3, random_state=42, stratify=y_small)
overfit = make_pipeline(StandardScaler(), LogisticRegression(C=10000, max_iter=5000))
overfit.fit(Xs_train, ys_train)
print(f"Variables: {X_small.shape[1]}, citas de entrenamiento: {len(ys_train)}")
print(f"Exactitud en entrenamiento: {overfit.score(Xs_train, ys_train):.3f}")
print(f"Exactitud en prueba: {overfit.score(Xs_test, ys_test):.3f}")
print(f"Predecir siempre 'acude': {1 - ys_test.mean():.3f}")
Variables: 56, citas de entrenamiento: 210 Exactitud en entrenamiento: 1.000 Exactitud en prueba: 0.633 Predecir siempre 'acude': 0.722
Perfecto con el conjunto de entrenamiento y, con citas nuevas, peor que suponer simplemente que todo el mundo acude. Eso es el sobreajuste.
Paso 5: Regularización L2 y L1¶
La regularización añade al coste una penalización que crece con el tamaño de los pesos, de modo que el modelo solo puede usar pesos grandes cuando de verdad reducen el error:
$$J_{L2}(w, b) = \text{log loss} + \lambda \sum_j w_j^2 \qquad\qquad J_{L1}(w, b) = \text{log loss} + \lambda \sum_j |w_j|$$
- L2 (Ridge) encoge todos los pesos hacia cero, sin anularlos.
- L1 (Lasso) puede dejar pesos exactamente en cero, así que además selecciona variables.
- Cuanto mayor es $\lambda$, más fuerte es la regularización. El sesgo $b$ no se penaliza.
Cuidado con la convención de scikit-learn. LogisticRegression no recibe $\lambda$, sino su inversa, C: una C grande significa una regularización débil. Para las fórmulas de arriba, con $m$ ejemplos de entrenamiento, la equivalencia es $C = \frac{1}{2 \lambda m}$ para L2 y $C = \frac{1}{\lambda m}$ para L1.
m = len(ys_train)
def l2_model(lam):
return make_pipeline(StandardScaler(), LogisticRegression(C=1 / (2 * lam * m), max_iter=5000))
def l1_model(lam):
return make_pipeline(StandardScaler(),
LogisticRegression(C=1 / (lam * m), penalty='l1', solver='liblinear', max_iter=5000, random_state=0))
for lam in [0.001, 0.1]:
mdl = l2_model(lam).fit(Xs_train, ys_train)
w = mdl[-1].coef_[0]
print(f"L2, lambda = {lam:<5}: tamaño de los pesos (suma de |w|) = {np.abs(w).sum():6.2f}, "
f"entrenamiento = {mdl.score(Xs_train, ys_train):.3f}, prueba = {mdl.score(Xs_test, ys_test):.3f}")
L2, lambda = 0.001: tamaño de los pesos (suma de |w|) = 26.14, entrenamiento = 0.943, prueba = 0.733 L2, lambda = 0.1 : tamaño de los pesos (suma de |w|) = 3.81, entrenamiento = 0.857, prueba = 0.744
Un $\lambda$ mayor encoge los pesos y cierra la distancia entre entrenamiento y prueba. Pero ¿qué $\lambda$ usamos? No podemos elegirlo mirando el conjunto de prueba: entonces dejaría de ser "datos no vistos" y nuestra estimación final sería optimista.
Paso 6: Elegir $\lambda$ con validación cruzada de k grupos (las notas de Luis)¶
La validación cruzada de k grupos responde a la pregunta usando solo los datos de entrenamiento:
- Dividir el conjunto de entrenamiento en $k$ grupos (folds) de tamaño parecido. Con grupos estratificados, cada grupo conserva la proporción de ausencias.
- Para cada valor candidato de $\lambda$, entrenar $k$ veces: cada vez con $k-1$ grupos, validando con el restante.
- Promediar las $k$ puntuaciones de validación. Esa media estima cómo se comportará el modelo con datos nuevos.
- Elegir el $\lambda$ con la mejor media y volver a entrenar con todo el conjunto de entrenamiento.
- Evaluar una sola vez en el conjunto de prueba, que no se ha usado para nada más.
Como recomienda el capítulo, usamos $k = 5$ y una rejilla logarítmica de valores: $\lambda \in \{0.001, 0.01, 0.1, 1, 10, 100, 1000\}$.
cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
lambdas = [0.001, 0.01, 0.1, 1, 10, 100, 1000]
rows = []
for lam in lambdas:
scores = cross_validate(l2_model(lam), Xs_train, ys_train, cv=cv, return_train_score=True)
rows.append({'lambda': lam, 'entrenamiento': scores['train_score'].mean(), 'validación': scores['test_score'].mean()})
results = pd.DataFrame(rows)
print(results.round(3).to_string(index=False))
best_lambda = results.loc[results['validación'].idxmax(), 'lambda']
print(f"\nPredecir siempre 'acude' en el conjunto de entrenamiento: {1 - ys_train.mean():.3f}")
print(f"Mejor lambda (mayor exactitud media de validación): {best_lambda:g}")
plt.figure(figsize=(8, 4))
plt.semilogx(results['lambda'], results['entrenamiento'], 'o-', label='entrenamiento (media de los 5 grupos)')
plt.semilogx(results['lambda'], results['validación'], 's-', label='validación (media de los 5 grupos)')
plt.axvline(best_lambda, color='gray', linestyle='--')
plt.xlabel('lambda (intensidad de la regularización)')
plt.ylabel('Exactitud')
plt.title('Validación cruzada de 5 grupos, regularización L2')
plt.legend()
plt.grid(True)
plt.show()
lambda entrenamiento validación 0.001 0.977 0.752 0.010 0.940 0.790 0.100 0.845 0.752 1.000 0.714 0.714 10.000 0.714 0.714 100.000 0.714 0.714 1000.000 0.714 0.714 Predecir siempre 'acude' en el conjunto de entrenamiento: 0.714 Mejor lambda (mayor exactitud media de validación): 0.01
La gráfica muestra los dos peligros. A la izquierda ($\lambda$ pequeño), la exactitud de entrenamiento es muy alta, pero la de validación es pobre: sobreajuste. A la derecha ($\lambda$ grande), los pesos quedan aplastados hacia cero y el modelo ya no aprende nada: subajuste. A partir de $\lambda = 1$, las dos curvas se quedan clavadas justo en la exactitud de predecir que todos los pacientes acuden. El mejor valor está en medio.
El capítulo aconseja además refinar la búsqueda alrededor del mejor valor. Probamos algunos valores cercanos y después evaluamos el modelo elegido, una sola vez, en el conjunto de prueba:
fine = [best_lambda * f for f in [0.3, 0.5, 1, 2, 3]]
fine_scores = [cross_validate(l2_model(lam), Xs_train, ys_train, cv=cv)['test_score'].mean() for lam in fine]
best_lambda = fine[int(np.argmax(fine_scores))]
print("Búsqueda refinada:", ", ".join(f"lambda {l:g} -> {s:.3f}" for l, s in zip(fine, fine_scores)))
print(f"Lambda elegido: {best_lambda:g}\n")
regularized = l2_model(best_lambda).fit(Xs_train, ys_train)
print(f"{'':<24}{'entrenamiento':>15}{'prueba':>8}")
for name, mdl in [('Sin regularización', overfit), (f'L2, lambda = {best_lambda:g}', regularized)]:
print(f"{name:<24}{mdl.score(Xs_train, ys_train):>15.3f}{mdl.score(Xs_test, ys_test):>8.3f}")
Búsqueda refinada: lambda 0.003 -> 0.767, lambda 0.005 -> 0.767, lambda 0.01 -> 0.790, lambda 0.02 -> 0.786, lambda 0.03 -> 0.786
Lambda elegido: 0.01
entrenamiento prueba
Sin regularización 1.000 0.633
L2, lambda = 0.01 0.929 0.811
El modelo regularizado cede un poco de exactitud en el entrenamiento y gana mucha con citas nuevas.
L1 como selector de variables. Con la misma validación cruzada elegimos $\lambda$ para L1 y después miramos cuáles de los 56 pesos sobreviven:
l1_grid = [0.001, 0.003, 0.01, 0.03, 0.1, 0.3]
l1_scores = [cross_validate(l1_model(lam), Xs_train, ys_train, cv=cv)['test_score'].mean() for lam in l1_grid]
best_l1 = l1_grid[int(np.argmax(l1_scores))]
lasso = l1_model(best_l1).fit(Xs_train, ys_train)
weights = pd.Series(lasso[-1].coef_[0], index=X_small.columns)
kept = weights[weights != 0]
print(f"L1, lambda = {best_l1}: exactitud en prueba = {lasso.score(Xs_test, ys_test):.3f}")
print(f"Variables con peso distinto de cero: {len(kept)} de {len(weights)}")
print(kept.round(3).sort_values(key=np.abs, ascending=False).to_string())
L1, lambda = 0.03: exactitud en prueba = 0.867 Variables con peso distinto de cero: 10 de 56 dias_desde_ultima_cita 0.934 edad -0.759 ausencias_previas 0.749 distancia_km 0.500 primera_visita 0.224 extra_5 0.197 extra_3 -0.057 extra_32 -0.057 extra_30 0.046 extra_43 -0.028
L1 ha eliminado la mayoría de las variables. Los pesos más grandes corresponden a variables reales (días desde la última cita, edad, ausencias previas, distancia) y las pocas columnas de ruido que sobreviven tienen pesos pequeños: con 50 columnas aleatorias y solo 210 citas, algunas se correlacionan con la asistencia por pura casualidad. También ha desaparecido la hora, porque su efecto solo aparece en las últimas citas de la tarde y 210 citas son pocas para detectarlo. La regularización ayuda mucho, pero no sustituye al buen criterio sobre qué variables tienen sentido.
Paso 7: Elastic Net y variables correlacionadas¶
Una historia clínica suele guardar la misma información más de una vez. Supongamos que, además de la distancia de casa al hospital, el registro tiene la distancia por carretera y el tiempo de trayecto en coche. Las tres columnas dicen casi lo mismo. El libro recomienda Elastic Net para esta situación, cuando hay variables correlacionadas. Añade al coste las dos penalizaciones:
$$J_{\text{Elastic Net}}(w, b) = \text{log loss} + \lambda_1 \sum_j |w_j| + \lambda_2 \sum_j w_j^2$$
La parte L1 puede seguir eliminando variables, y la parte L2 reparte el peso entre las variables que llevan la misma información. En scikit-learn es penalty='elasticnet' con el optimizador saga, y las dos lambdas se convierten en C y l1_ratio: $C = \frac{1}{m(\lambda_1 + 2\lambda_2)}$ y $\text{l1\_ratio} = \frac{\lambda_1}{\lambda_1 + 2\lambda_2}$.
Volvemos a las 900 citas de entrenamiento del paso 1 y les añadimos las dos columnas nuevas y 20 irrelevantes:
rng = np.random.default_rng(11)
def add_columns(X_part):
X_part = X_part.copy()
X_part['km_carretera'] = np.round(1.25 * X_part['distancia_km'] + rng.normal(0, 1.5, len(X_part)), 1) # distancia por carretera
X_part['minutos_trayecto'] = np.round(1.5 * X_part['distancia_km'] + rng.normal(0, 2.5, len(X_part))) # minutos en coche
for j in range(20):
X_part[f'extra_{j}'] = rng.normal(size=len(X_part)) # variables irrelevantes
return X_part
Xg_train, Xg_test = add_columns(X_train), add_columns(X_test)
distance_group = ['distancia_km', 'km_carretera', 'minutos_trayecto']
print(f"Variables: {Xg_train.shape[1]}, citas de entrenamiento: {len(y_train)}")
print("Correlaciones entre las tres columnas de distancia:")
print(Xg_train[distance_group].corr().round(3).to_string())
Variables: 28, citas de entrenamiento: 900
Correlaciones entre las tres columnas de distancia:
distancia_km km_carretera minutos_trayecto
distancia_km 1.000 0.991 0.983
km_carretera 0.991 1.000 0.975
minutos_trayecto 0.983 0.975 1.000
Elegimos las lambdas de L1, L2 y Elastic Net con la misma validación cruzada de 5 grupos del paso 6. Elastic Net tiene dos, así que su rejilla tiene 12 combinaciones en lugar de 5 valores: es el ajuste adicional del que advierte el libro. Después entrenamos cada modelo con las 900 citas y miramos los pesos de las tres columnas de distancia, cuántas columnas irrelevantes conservan un peso y la exactitud en el conjunto de prueba.
m = len(y_train) # ahora son 900 citas: l1_model, l2_model y la nueva función lo usan para convertir lambda en C
def elastic_net_model(lam1, lam2):
return make_pipeline(StandardScaler(),
LogisticRegression(C=1 / (m * (lam1 + 2 * lam2)), l1_ratio=lam1 / (lam1 + 2 * lam2),
penalty='elasticnet', solver='saga', max_iter=20000, tol=1e-5, random_state=0))
def mean_validation(model):
return cross_validate(model, Xg_train, y_train, cv=cv)['test_score'].mean()
grid = [0.001, 0.003, 0.01, 0.03, 0.1]
lam_l1 = max(grid, key=lambda lam: mean_validation(l1_model(lam)))
lam_l2 = max(grid, key=lambda lam: mean_validation(l2_model(lam)))
en_grid = [(lam1, lam2) for lam1 in [0.003, 0.01, 0.03] for lam2 in [0.003, 0.01, 0.03, 0.1]]
lam1, lam2 = max(en_grid, key=lambda lams: mean_validation(elastic_net_model(*lams)))
l1_name, en_name = f'L1, lambda = {lam_l1:g}', f'Elastic Net, lambdas = {lam1:g}, {lam2:g}'
models = {'Sin regularización': make_pipeline(StandardScaler(), LogisticRegression(C=10000, max_iter=5000)),
l1_name: l1_model(lam_l1),
f'L2, lambda = {lam_l2:g}': l2_model(lam_l2),
en_name: elastic_net_model(lam1, lam2)}
rows = []
for name, mdl in models.items():
weights = pd.Series(mdl.fit(Xg_train, y_train)[-1].coef_[0], index=Xg_train.columns)
rows.append({'modelo': name, **weights[distance_group].round(3).to_dict(),
'irrelevantes': int((weights.filter(like='extra_') != 0).sum()),
'exactitud prueba': round(mdl.score(Xg_test, y_test), 3)})
print(pd.DataFrame(rows).to_string(index=False))
modelo distancia_km km_carretera minutos_trayecto irrelevantes exactitud prueba
Sin regularización 0.997 0.520 -0.415 20 0.830
L1, lambda = 0.01 0.463 0.334 0.000 9 0.840
L2, lambda = 0.001 0.527 0.553 -0.025 20 0.823
Elastic Net, lambdas = 0.01, 0.003 0.338 0.335 0.103 9 0.840
Las tres columnas de distancia tienen correlaciones de entre 0,975 y 0,991. Sin regularización, el modelo reparte la distancia entre ellas de una forma que no tiene sentido: 0,997 para la distancia, 0,520 para la carretera y −0,415 para el tiempo de trayecto, un peso negativo que, por sí solo, significaría que los trayectos largos hacen que los pacientes acudan más. Además, da un peso a las 20 columnas irrelevantes. L2, con la lambda elegida por validación cruzada, 0,001, es muy ligera y apenas cambia ese panorama: el tiempo de trayecto sigue teniendo un pequeño peso negativo, −0,025, y las 20 columnas irrelevantes conservan el suyo. L1 elimina 11 de las 20 columnas irrelevantes y también el tiempo de trayecto, y se queda con 0,463 para la distancia y 0,334 para la carretera. Elastic Net elimina las mismas columnas irrelevantes, quedan 9, y conserva las tres columnas de distancia con pesos positivos: 0,338, 0,335 y 0,103. En el conjunto de prueba, los cuatro modelos están entre 0,823 y 0,840, una diferencia de 5 citas de 300: con 900 citas de entrenamiento, la diferencia está en los pesos, no en la exactitud.
¿Hasta qué punto es estable la elección de L1? Volvemos a entrenar L1 y Elastic Net con 30 muestras bootstrap, cada una formada por 900 citas extraídas al azar, con reemplazo, del conjunto de entrenamiento, y contamos cuántas veces conserva un peso cada columna de distancia:
boot_rng = np.random.default_rng(0)
samples = [boot_rng.integers(0, m, m) for _ in range(30)]
for name in [l1_name, en_name]:
kept = pd.Series(0, index=distance_group)
for idx in samples:
mdl = models[name].fit(Xg_train.iloc[idx], y_train.iloc[idx])
kept += (pd.Series(mdl[-1].coef_[0], index=Xg_train.columns)[distance_group] != 0).astype(int)
print(f"{name:<32}", ", ".join(f"{col} {n} de 30" for col, n in kept.items()))
L1, lambda = 0.01 distancia_km 23 de 30, km_carretera 22 de 30, minutos_trayecto 8 de 30
Elastic Net, lambdas = 0.01, 0.003 distancia_km 30 de 30, km_carretera 30 de 30, minutos_trayecto 26 de 30
Qué columna de distancia conserva L1 depende de la muestra: la distancia en 23 de las 30 muestras, la carretera en 22 y el tiempo de trayecto solo en 8. Dos modelos entrenados con casi las mismas citas cuentan historias distintas sobre qué variables importan. Elastic Net conserva la distancia y la carretera en las 30 muestras y el tiempo de trayecto en 26: su parte L2 hace que las columnas correlacionadas se repartan el peso, así que tienden a quedarse juntas en el modelo. Eso es lo que dice el libro: Elastic Net puede seleccionar variables, como L1, y reducir proporcionalmente los pesos de las demás, como L2, y por eso es especialmente útil cuando hay variables correlacionadas.
Paso 8: Clases desbalanceadas¶
Solo alrededor de una cuarta parte de las citas son ausencias, así que el modelo ve muchos más ejemplos de asistencia y tiende a predecirla. Además de mover el umbral, scikit-learn puede ponderar las clases para que cada error en la clase minoritaria cuente más (class_weight='balanced'). Volvemos al conjunto de datos completo:
def sensitivity_specificity(y_true, y_hat):
tn, fp, fn, tp = confusion_matrix(y_true, y_hat).ravel()
return tp / (tp + fn), tn / (tn + fp)
balanced = make_pipeline(StandardScaler(), LogisticRegression(class_weight='balanced'))
balanced.fit(X_train, y_train)
for name, y_hat in [('Modelo por defecto, umbral 0,5', y_pred),
('Modelo por defecto, umbral 0,3', (p_test >= 0.3).astype(int)),
("class_weight='balanced', umbral 0,5", balanced.predict(X_test))]:
sens, spec = sensitivity_specificity(y_test, y_hat)
print(f"{name:<38} exactitud {accuracy_score(y_test, y_hat):.3f} | sensibilidad {sens:.3f} | especificidad {spec:.3f}")
Modelo por defecto, umbral 0,5 exactitud 0.827 | sensibilidad 0.620 | especificidad 0.900 Modelo por defecto, umbral 0,3 exactitud 0.787 | sensibilidad 0.772 | especificidad 0.792 class_weight='balanced', umbral 0,5 exactitud 0.777 | sensibilidad 0.797 | especificidad 0.769
Ponderar las clases tiene un efecto parecido a bajar el umbral: se detectan más ausencias a cambio de algo de exactitud. Ambas son formas legítimas de decirle al modelo que no detectar una ausencia nos importa más que hacer una llamada de recordatorio.
Conclusión¶
- Con scikit-learn, un pipeline de
StandardScaleryLogisticRegressionreproduce en pocas líneas el modelo que escribimos desde cero. - Un clasificador se juzga con la matriz de confusión, el informe de clasificación y la curva ROC, no solo con la exactitud.
- Con pocos ejemplos y muchas variables el modelo se sobreajusta. La regularización penaliza los pesos grandes: L2 los encoge y L1 puede eliminar variables por completo.
- Elastic Net combina las dos penalizaciones: elimina variables irrelevantes como L1 y reparte el peso entre variables correlacionadas como L2, a cambio de una segunda $\lambda$ que ajustar.
- La intensidad de la regularización se elige con validación cruzada de k grupos sobre los datos de entrenamiento. El conjunto de prueba se usa una sola vez, al final.
- Con clases desbalanceadas, el umbral y los pesos de las clases permiten situar el equilibrio entre sensibilidad y especificidad donde lo necesita el hospital.
Los siguientes cuadernos abordan el mismo problema con el otro modelo de este capítulo, las máquinas de vectores de soporte.