Regresión lineal y descenso de gradiente

1. El modelo: la forma de la regla que usamos para predecir.
2. Los parámetros: los números del modelo que se ajustan con los datos.
3. La función de pérdida: un solo número que dice qué tan mal predice el modelo con esos parámetros.
4. El entrenamiento: buscar los parámetros que hacen mínima la pérdida.
5. Las horas son la característica (la entrada, x) y la nota es el objetivo (lo que queremos predecir, y). Como cada ejemplo trae su respuesta, es aprendizaje supervisado.

Las predicciones son 3, 4, 5 y 6. Los residuos (real menos predicho) son 0, 1, −1 y 2.

• MSE = (0 + 1 + 1 + 4) / 4 = 1,5
• RMSE = √1,5 ≈ 1,22
• MAE = (0 + 1 + 1 + 2) / 4 = 1

x = np.array([1, 2, 3, 4])
y = np.array([3, 5, 4, 8])
p = x + 2
print("Residuos:", y - p)
print(f"MSE = {mse(y, p)}   RMSE = {np.sqrt(mse(y, p)):.2f}"
      f"   MAE = {mae(y, p)}")

Residuos: [ 0  1 -1  2]
MSE = 1.5   RMSE = 1.22   MAE = 1.0

a) La derivada de una suma es la suma de las derivadas, y la constante 1/n sale afuera. Para cada término usamos la regla de la cadena: la derivada de u² es 2u por la derivada de u. Aquí u = w·xᵢ + b − yᵢ, y su derivada respecto de b es 1 (w, xᵢ e yᵢ no dependen de b). Entonces

∂MSE/∂b = (1/n) · Σ 2 · (w·xᵢ + b − yᵢ) · 1 = (2/n) · Σ (ŷᵢ − yᵢ),

es decir, el doble del error promedio (predicho menos real).

b) Los errores ŷ − y son 0, −1, 1 y −2.

• ∂MSE/∂b = 2 · (0 − 1 + 1 − 2) / 4 = −1
• ∂MSE/∂w = 2 · (0·1 − 1·2 + 1·3 − 2·4) / 4 = 2 · (−7/4) = −3,5

c) w = 1 − 0,1 · (−3,5) = 1,35 y b = 2 − 0,1 · (−1) = 2,1. El MSE baja de 1,5 a 1,279. Las dos derivadas eran negativas, así que el paso aumentó w y b: la recta se inclinó más y subió, porque se quedaba corta con el último punto.

dw, db = gradiente(x, y, 1, 2)
print("dw =", dw, "  db =", db)
w, b = 1 - 0.1 * dw, 2 - 0.1 * db
print(f"Después del paso: w = {w:.2f}   b = {b:.2f}")
print(f"MSE antes = {mse(y, predecir(x, 1, 2)):.3f}"
      f"   después = {mse(y, predecir(x, w, b)):.3f}")

dw = -3.5   db = -1.0
Después del paso: w = 1.35   b = 2.10
MSE antes = 1.500   después = 1.279

pasos_tabla = [0, 1, 2, 5, 10, 50, 100, 300]
curvas = {}
for nombre, tasa in [("A", 0.03), ("B", 0.039), ("C", 0.005)]:
    _, _, h = descenso_gradiente(horas, nota, tasa=tasa, pasos=300)
    curvas[nombre] = h[pasos_tabla, 2]

tabla = pd.DataFrame({"Paso": pasos_tabla, **curvas})
print(tabla.to_string(index=False, float_format=lambda v: f"{v:.4g}"))

• A, adecuada: baja rápido y en el paso 300 ya está en 22,57, casi el mínimo (22,47, el de la teoría). Se puede parar.
• B, demasiado grande: el error sube desde el primer paso y cada vez más rápido. El descenso diverge: hay que bajar la tasa. Fíjate en que al principio sube despacio: con una tasa apenas por encima del límite, la divergencia tarda en notarse.
• C, lenta: baja todo el tiempo, pero en 300 pasos todavía está en 74,43. Funciona, pero conviene subir la tasa o dar muchos más pasos.

Las tasas eran A = 0,03, B = 0,039 y C = 0,005.

pred = modelo.predict(horas.reshape(-1, 1))
antes = mse(nota, pred), mae(nota, pred)
despues = mse(nota_mala, pred), mae(nota_mala, pred)
print(f"MSE: {antes[0]:.1f} -> {despues[0]:.1f}"
      f"   ({despues[0] / antes[0]:.1f} veces)")
print(f"MAE: {antes[1]:.2f} -> {despues[1]:.2f}"
      f"   ({despues[1] / antes[1]:.1f} veces)")
print(f"Residuo del dato malo: {nota_mala[5] - pred[5]:.1f}")

modelo_malo = LinearRegression().fit(horas.reshape(-1, 1), nota_mala)
print(f"Antes:   w = {modelo.coef_[0]:.2f}   b = {modelo.intercept_:.2f}")
print(f"Después: w = {modelo_malo.coef_[0]:.2f}"
      f"   b = {modelo_malo.intercept_:.2f}")

MSE: 22.5 -> 261.6   (11.6 veces)
MAE: 3.86 -> 6.48   (1.7 veces)
Residuo del dato malo: -84.9
Antes:   w = 6.89   b = 31.18
Después: w = 5.02   b = 36.51

a) Con la misma recta, el MSE crece casi 12 veces y el MAE apenas 1,7 veces. Un solo residuo de unos 85 puntos, elevado al cuadrado, pesa más que los otros 29 juntos.

b) Un solo dato malo entre treinta le baja la pendiente en casi dos puntos por hora y le sube el intercepto en más de cinco. Como LinearRegression minimiza el MSE, hereda su sensibilidad a los atípicos. Por eso siempre conviene mirar los datos antes de entrenar (tema 4) y, si hay atípicos que no se pueden corregir, considerar el MAE.

for tasa in [0.036, 0.037, 0.038, 0.0385, 0.039]:
    _, _, h = descenso_gradiente(horas, nota, tasa=tasa, pasos=2000)
    print(f"tasa {tasa}:  MSE final = {h[-1, 2]:.4g}")

tasa 0.036:  MSE final = 22.47
tasa 0.037:  MSE final = 22.47
tasa 0.038:  MSE final = 22.47
tasa 0.0385:  MSE final = 1.146e+36
tasa 0.039:  MSE final = 1.576e+80

Con 0,038 todavía converge al mínimo (22,47) y con 0,0385 ya diverge: el límite está entre 0,038 y 0,0385. Ese límite depende de los datos: cuanto más grandes son los valores de x, más empinadas son las paredes del valle y más chico tiene que ser el paso. Es otra razón para estandarizar.

w_300, b_300 = descenso_gradiente_varias(X_est, y_arr, tasa=0.1, pasos=300)
error_300 = np.mean((X_est @ w_300 + b_300 - y_arr) ** 2)
error_min = np.mean((modelo_diabetes.predict(X_est) - y_arr) ** 2)
print(f"MSE con 300 pasos: {error_300:.1f}   mínimo: {error_min:.1f}"
      f"   ({100 * (error_300 / error_min - 1):.1f} % más)")

pesos = pd.DataFrame({"300 pasos": w_300,
                      "scikit-learn": modelo_diabetes.coef_},
                     index=X.columns)
print(pesos.loc[["bmi", "s1", "s2"]].round(2))

MSE con 300 pasos: 2867.7   mínimo: 2859.7   (0.3 % más)
     300 pasos  scikit-learn
bmi      24.98         24.73
s1      -15.98        -37.68
s2        5.44         22.68

El MSE ya está a menos de medio punto porcentual del mínimo y el peso de bmi está casi bien, pero los pesos de s1 y s2 siguen lejos de los de scikit-learn.

La explicación es que s1 (colesterol total) y s2 (LDL) están muy correlacionadas: subir una y bajar la otra casi no cambia las predicciones. En el valle eso es un fondo largo y casi plano, y el descenso avanza por ahí muy despacio. Moraleja: que la pérdida ya casi no baje no garantiza que cada peso haya llegado a su valor final, y los pesos de variables muy relacionadas entre sí hay que interpretarlos con cuidado.

print("Predicciones para 0 y 12 horas:",
      modelo.predict([[0], [12]]).round(1))
print(f"Horas en los datos: de {horas.min():.2f} a {horas.max():.2f}")

Predicciones para 0 y 12 horas: [ 31.2 113.9]
Horas en los datos: de 0.35 a 7.80

La de 12 horas no tiene sentido: el examen es sobre 100. El modelo solo vio estudiantes de 0,35 a 7,80 horas; fuera de ese rango no sabe nada, y una recta sigue subiendo para siempre. Predecir fuera del rango de los datos se llama extrapolar. Las 0 horas también están fuera del rango, pero apenas: 0,35 horas por debajo del que menos estudió. Es una extrapolación pequeña y la predicción es razonable. Las 12 horas, en cambio, están más de 4 horas por encima del que más estudió, y ahí la recta ya inventa.

• 0,95: el modelo explica el 95 % de la variación del objetivo; sus errores son chicos comparados con cuánto varían los datos.
• 0: predecir siempre el promedio daría el mismo error. El modelo no aporta nada.
• −0,3: es peor que predecir el promedio. Con una regresión lineal bien entrenada y medida en sus propios datos de entrenamiento eso no puede pasar, así que algo falló: por ejemplo, un descenso que divergió o que paró demasiado pronto.

Un R² alto en los datos de entrenamiento todavía no dice cómo le irá al modelo con datos nuevos. Eso lo ves en el tema 7.

Este código hace los pasos 1 a 5 sin gráficos; tus números deberían coincidir.

from sklearn.metrics import r2_score

ingreso = X["MedInc"].to_numpy()
print("Zonas con el valor máximo:", (y == y.max()).sum(),
      "  valor máximo:", y.max())

# Paso 2: una variable, desde cero, con MedInc estandarizado
m, s = ingreso.mean(), ingreso.std()
w_e, b_e, _ = descenso_gradiente((ingreso - m) / s, y, tasa=0.1, pasos=500)
print(f"Descenso, estandarizado:  w = {w_e:.3f}   b = {b_e:.3f}")
print(f"Descenso, original:       w = {w_e / s:.3f}"
      f"   b = {b_e - w_e * m / s:.3f}")

# Paso 3: una variable, con scikit-learn
una = LinearRegression().fit(ingreso.reshape(-1, 1), y)
print(f"scikit-learn:             w = {una.coef_[0]:.3f}"
      f"   b = {una.intercept_:.3f}")

# Paso 4: las ocho variables, estandarizadas
X8 = ((X - X.mean()) / X.std(ddof=0)).to_numpy()
todas = LinearRegression().fit(X8, y)
print(pd.Series(todas.coef_, index=X.columns).round(2).to_string())

# Paso 5: comparar
for nombre, p in [("1 variable ", una.predict(ingreso.reshape(-1, 1))),
                  ("8 variables", todas.predict(X8))]:
    print(f"{nombre}:  RMSE = {np.sqrt(mse(y, p)):.3f}"
          f"   MAE = {mae(y, p):.3f}   R² = {r2_score(y, p):.3f}")

Zonas con el valor máximo: 965   valor máximo: 5.00001
Descenso, estandarizado:  w = 0.794   b = 2.069
Descenso, original:       w = 0.418   b = 0.451
scikit-learn:             w = 0.418   b = 0.451
MedInc        0.83
HouseAge      0.12
AveRooms     -0.27
AveBedrms     0.31
Population   -0.00
AveOccup     -0.04
Latitude     -0.90
Longitude    -0.87
1 variable :  RMSE = 0.837   MAE = 0.626   R² = 0.473
8 variables:  RMSE = 0.724   MAE = 0.531   R² = 0.606

• Muchas zonas tienen exactamente el valor máximo: los valores del censo se recortaron arriba («500 000 dólares o más»). Eso deja una línea recta de residuos en el gráfico, porque el modelo predice valores por encima de 5 que nunca aparecen en los datos.
• Con solo MedInc, cada 10 000 dólares de ingreso mediano suman unos 41 800 dólares al valor de la vivienda. El descenso con MedInc estandarizado llega a lo mismo que scikit-learn, una vez que se vuelve a las unidades originales.
• Con las 8 variables, el R² sube de 0,47 a 0,61. Los pesos más grandes del modelo estandarizado son los de Latitude, Longitude y MedInc: la ubicación importa casi tanto como el ingreso.
• Todo esto está medido con los mismos datos del entrenamiento. En el tema 7 vas a repetir el mini-proyecto midiendo con datos que el modelo no vio.