Regresión lineal y logística

El problema

load_diabetes contiene 442 pacientes con diabetes. Para cada paciente hay diez mediciones iniciales y un objetivo numérico: la progresión de la enfermedad un año después. Estas son tres filas con seis de las diez columnas:

age está en años y bmi es el índice de masa corporal. sex solo toma los valores 1 y 2. s1 (colesterol total) va de 97 a 301 en el conjunto completo; s5 (logaritmo de los triglicéridos) va de 3.26 a 6.11.

Necesitas un modelo que prediga el objetivo y una forma de decir qué columnas empujan la predicción hacia arriba o hacia abajo. Después cambias de tarea: con load_breast_cancer hay que marcar tumores como malignos, y dejar pasar un tumor maligno cuesta más que revisar de nuevo uno benigno.

Coeficientes y predicción

Una regresión lineal predice con una suma ponderada: un intercepto más un coeficiente por columna multiplicado por el valor de esa columna. LinearRegression elige el intercepto y los coeficientes que minimizan el error cuadrático en los datos de entrenamiento.

Separamos una partición de validación para medir el error con pacientes que el modelo no vio.

Pregunta. ¿Qué signo esperas para el coeficiente de bmi?

import numpy as np
import pandas as pd
from sklearn.datasets import load_diabetes
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split

diabetes = load_diabetes(as_frame=True, scaled=False)
X, y = diabetes.data, diabetes.target
X_ent, X_val, y_ent, y_val = train_test_split(
    X, y, test_size=0.25, random_state=0
)

lineal = LinearRegression().fit(X_ent, y_ent)

print(X_ent.shape, X_val.shape)
print(round(lineal.intercept_, 2))
print(pd.Series(lineal.coef_, index=X.columns).round(2).to_string())

(331, 10) (111, 10)
-305.97
age    -0.16
sex   -19.89
bmi     6.40
bp      1.04
s1     -0.77
s2      0.41
s3     -0.03
s4      5.02
s5     64.10
s6      0.12

Reconstruimos a mano la predicción del primer paciente de validación con esa suma. Después comparamos el error absoluto medio (MAE) contra un DummyRegressor, que predice siempre la media de entrenamiento.

Pregunta. ¿Coincidirá la suma manual con predict? ¿El MAE del modelo lineal será menor o mayor que el del baseline?

from sklearn.dummy import DummyRegressor
from sklearn.metrics import mean_absolute_error

paciente = X_val.iloc[0]
manual = lineal.intercept_ + (lineal.coef_ * paciente).sum()
modelo = lineal.predict(X_val.iloc[[0]])[0]

base = DummyRegressor().fit(X_ent, y_ent)
mae_base = mean_absolute_error(y_val, base.predict(X_val))
mae_lineal = mean_absolute_error(y_val, lineal.predict(X_val))

print(round(manual, 2), round(modelo, 2), y_val.iloc[0])
print(round(mae_base, 2), round(mae_lineal, 2))

241.85 241.85 321.0
58.31 45.12

Escala de las características

Un coeficiente mide el cambio en la predicción por una unidad de su columna. Una unidad de sex es pasar de 1 a 2, todo su rango. Una unidad de s1 es una fracción pequeña de un rango que va de 97 a 301. Con unidades tan distintas, el tamaño de un coeficiente crudo depende de la unidad elegida.

StandardScaler transforma cada columna para que tenga media 0 y desviación estándar 1 en entrenamiento. Después de escalar, cada coeficiente mide el cambio por una desviación estándar de su columna.

Pregunta. Después de escalar, ¿cambiará el MAE de validación? ¿s5 seguirá teniendo el coeficiente más grande en valor absoluto? ¿Qué pasará con sex y con s1?

from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

escalado = make_pipeline(StandardScaler(), LinearRegression()).fit(X_ent, y_ent)
mae_escalado = mean_absolute_error(y_val, escalado.predict(X_val))

coeficientes = pd.DataFrame(
    {
        "crudo": lineal.coef_,
        "desv_est": X_ent.std(ddof=0),
        "escalado": escalado[-1].coef_,
    },
    index=X.columns,
)

print(coeficientes.round(2).to_string())
print(round(mae_escalado, 2))

     crudo  desv_est  escalado
age  -0.16     13.25     -2.08
sex -19.89      0.50     -9.94
bmi   6.40      4.51     28.86
bp    1.04     13.75     14.34
s1   -0.77     34.40    -26.52
s2    0.41     30.25     12.38
s3   -0.03     12.99     -0.42
s4    5.02      1.31      6.60
s5   64.10      0.54     34.46
s6    0.12     12.05      1.41
45.12

Regularización L2 y L1

s1 es el colesterol total y s2 el colesterol LDL, que forma parte del total. Cuando dos columnas están muy correlacionadas, el modelo puede subir una y bajar la otra casi sin cambiar las predicciones. Sus coeficientes pueden salir grandes y de signo opuesto.

La regularización suma una penalización al error que el modelo minimiza. Ridge (L2) penaliza alpha por la suma de los coeficientes al cuadrado. Lasso (L1) penaliza alpha por la suma de sus valores absolutos. Con alpha más grande, los coeficientes grandes cuestan más. La penalización trata igual a todas las columnas, así que se aplica sobre columnas escaladas.

Pregunta. ¿Cuál de los dos modelos regularizados tendrá coeficientes exactamente en cero? ¿Qué pasará con los coeficientes de s1 y s2?

from sklearn.linear_model import Lasso, Ridge

ridge = make_pipeline(StandardScaler(), Ridge(alpha=100)).fit(X_ent, y_ent)
lasso = make_pipeline(StandardScaler(), Lasso(alpha=5)).fit(X_ent, y_ent)

comparacion = pd.DataFrame(
    {
        "lineal": escalado[-1].coef_,
        "ridge_100": ridge[-1].coef_,
        "lasso_5": lasso[-1].coef_,
    },
    index=X.columns,
)

print(round(X_ent["s1"].corr(X_ent["s2"]), 2))
print((comparacion.round(1) + 0.0).to_string())

0.89
     lineal  ridge_100  lasso_5
age    -2.1       -0.5      0.0
sex    -9.9       -7.0     -1.4
bmi    28.9       23.3     28.0
bp     14.3       12.3      8.4
s1    -26.5       -1.5      0.0
s2     12.4       -4.6      0.0
s3     -0.4       -9.5     -7.3
s4      6.6        5.5      0.0
s5     34.5       19.1     22.4
s6      1.4        3.9      0.0

Un coeficiente en cero saca la columna del modelo. Con alpha creciente, Lasso elimina más columnas. La tabla siguiente compara el mismo modelo Lasso con StandardScaler, ajustado con seis valores de alpha, de 0.1 a 100.

Pregunta. ¿Cuántos coeficientes serán cero con alpha=100? ¿A qué valor conocido se acercará el MAE?

Probabilidad con sigmoide

load_breast_cancer contiene 569 tumores con 30 medidas calculadas a partir de imágenes. En el objetivo original, 0 significa maligno y 1 significa benigno. Aquí lo invertimos para que 1 sea la clase que queremos detectar: maligno.

La regresión logística calcula la misma suma ponderada, z = intercepto + coeficientes · x, y la pasa por la función sigmoide:

sigmoide(z) = 1 / (1 + exp(-z))

Pregunta. Calcula sigmoide(0). Decide si sigmoide(-4) y sigmoide(4) quedan cerca de 0, de 0.5 o de 1.

def sigmoide(z):
    return 1 / (1 + np.exp(-z))


print(np.round(sigmoide(np.array([-4.0, -2.0, 0.0, 2.0, 4.0])), 3))

[0.018 0.119 0.5   0.881 0.982]

Ahora entrenamos LogisticRegression con columnas escaladas y comparamos la suma z de cinco tumores de validación con la probabilidad del modelo.

Pregunta. Si z de un tumor es negativa, ¿su probabilidad de malignidad será mayor o menor que 0.5? ¿Coincidirá sigmoide(z) con predict_proba?

from sklearn.datasets import load_breast_cancer
from sklearn.linear_model import LogisticRegression

tumores = load_breast_cancer(as_frame=True)
X_t = tumores.data
y_t = (tumores.target == 0).astype(int)
Xt_ent, Xt_val, yt_ent, yt_val = train_test_split(
    X_t, y_t, test_size=0.25, stratify=y_t, random_state=0
)

logistica = make_pipeline(StandardScaler(), LogisticRegression())
logistica.fit(Xt_ent, yt_ent)

z = logistica.decision_function(Xt_val.iloc[:5])
proba = logistica.predict_proba(Xt_val.iloc[:5])[:, 1]

print(np.round(z, 2))
print(np.round(sigmoide(z), 3))
print(np.round(proba, 3))
print(yt_val.iloc[:5].tolist())

[ 9.42 -4.11 -3.61 -4.69 -1.11]
[1.    0.016 0.026 0.009 0.248]
[1.    0.016 0.026 0.009 0.248]
[1, 0, 0, 0, 0]

Umbral de decisión

Una probabilidad todavía no es una clase. predict asigna la clase 1 cuando la probabilidad supera 0.5, pero ese umbral se puede cambiar sin volver a entrenar. La precisión y el recall se calculan para la clase maligno.

Pregunta. Compara los umbrales 0.9, 0.5 y 0.1. ¿Con cuál se marcan más tumores como malignos? ¿Qué le pasa al recall al bajar el umbral? ¿Y a la precisión?

from sklearn.metrics import precision_score, recall_score

proba_val = logistica.predict_proba(Xt_val)[:, 1]

print(int(yt_val.sum()), len(yt_val))
print(bool((logistica.predict(Xt_val) == (proba_val > 0.5)).all()))
for umbral in [0.9, 0.5, 0.1]:
    pred = (proba_val >= umbral).astype(int)
    print(
        umbral,
        int(pred.sum()),
        round(precision_score(yt_val, pred), 3),
        round(recall_score(yt_val, pred), 3),
    )

53 143
True
0.9 47 1.0 0.887
0.5 51 1.0 0.962
0.1 64 0.828 1.0

El umbral se elige según el error que más cuesta. En esta tarea, dejar pasar un maligno cuesta más, así que un umbral bajo es defendible.

Práctica

Traza guiada

Un modelo lineal tiene intercepto b y coeficientes w. Para la observación x, calcula a mano la suma z. Después escribe tres respuestas: la predicción si el modelo fuera una regresión lineal, la probabilidad si fuera una regresión logística y la clase con umbral 0.5 y con umbral 0.8.

b = 2.0
w = np.array([0.5, -1.0])
x = np.array([4.0, 3.0])

Comprueba con b + np.dot(w, x) y sigmoide(...). Si z no coincide, compara tu cálculo con w * x, que muestra cada producto por separado.

Umbral por recall mínimo

Modifica la celda de umbrales para recorrer np.arange(0.05, 0.96, 0.05) y quedarte con el umbral más alto que logra un recall de al menos 0.98 en validación. Antes de procesar, predice si ese umbral quedará por encima o por debajo de 0.5. Comprueba con assert que el recall del umbral elegido cumple la condición y anota cuántos benignos se marcan como malignos con él.

Coeficientes recuperados

Genera datos donde conoces los coeficientes verdaderos:

from sklearn.datasets import make_regression

X_s, y_s, coef_real = make_regression(
    n_samples=200,
    n_features=8,
    n_informative=3,
    noise=5.0,
    coef=True,
    random_state=0,
)

Solo tres de las ocho columnas influyen en y_s. Ajusta LinearRegression y Lasso(alpha=1), ambos con StandardScaler. Construye una tabla con coef_real y los dos vectores de coeficientes. Antes de procesar, predice cuál de los dos modelos dará ceros exactos. Comprueba con assert que las columnas con coeficiente cero en Lasso son las mismas que tienen coeficiente cero en coef_real. Para eso, compara modelo[-1].coef_ == 0 con coef_real == 0 usando np.array_equal, donde modelo es tu pipeline con Lasso.

Reto opcional

Añade a X_ent una copia exacta de bmi con X_ent.assign(bmi_copia=X_ent["bmi"]). Ajusta Ridge(alpha=10) y Lasso(alpha=1) con StandardScaler. Antes de procesar, predice cómo reparte cada modelo el coeficiente entre bmi y bmi_copia. Explica el reparto de cada modelo con la forma de su penalización: para dos coeficientes positivos a y b con la misma suma a + b, compara cómo cambian a² + b² y |a| + |b| al mover peso de uno al otro.

Errores que revisar

• Si comparas coeficientes de columnas con unidades distintas, escala primero con StandardScaler. Un coeficiente crudo grande puede venir de una unidad grande, no de una influencia grande.
• Si LogisticRegression muestra ConvergenceWarning, revisa si escalaste las columnas. Sin escalar, las 30 medidas de load_breast_cancer tienen rangos muy distintos y el optimizador no converge en las iteraciones por defecto.
• Si Lasso deja todos los coeficientes en cero, alpha es demasiado grande. El modelo predice la media de entrenamiento y empata con el baseline.
• Si predict_proba[:, 1] parece la probabilidad opuesta, revisa classes_ y la codificación del objetivo. En load_breast_cancer original, 1 significa benigno.
• Si eliges alpha o el umbral con los datos de prueba, el puntaje de prueba deja de estimar el desempeño con datos nuevos. Elige con validación.

Cierre

Con los números de esta lección, explica qué significa el coeficiente 28.86 de bmi en el modelo escalado y por qué el recall pasó de 0.887 a 1.0 sin volver a entrenar la regresión logística.

Lista de comprobación

• ¿Puedes reconstruir una predicción lineal con el intercepto y los coeficientes?
• ¿Puedes explicar por qué escalar cambia los coeficientes de LinearRegression pero no su MAE?
• ¿Sabes qué penalización produce ceros exactos y qué predice el modelo cuando todos los coeficientes son cero?
• ¿Puedes decir qué pasa con la precisión y el recall al bajar el umbral, y por qué?

Documentación oficial

• Modelos lineales en scikit-learn: mínimos cuadrados, Ridge, Lasso y regresión logística con sus funciones objetivo.
• Interpretación de coeficientes de modelos lineales: efecto de la escala y de la correlación entre columnas en los coeficientes.
• Ajuste del umbral de decisión: umbral de predict, probabilidades y cambio de umbral según la métrica.
• Datasets de juguete: descripción de las columnas de load_diabetes y load_breast_cancer.