Unidad 3.3 — Modelos de Regresion en Machine Learning

Unidad 3.3 — Modelos de Regresion en Machine Learning#

Esta unidad cubre de forma teorica y practica todos los modelos de regresion, desde los fundamentos de correlacion hasta tecnicas avanzadas como Gradient Boosting y SVR. Cada seccion incluye la explicacion conceptual seguida de codigo ejecutable en Python.

0. Instalacion de dependencias y preparacion del entorno#

Antes de comenzar, instalamos y cargamos todas las librerias necesarias.

# Instalacion (descomentar si es necesario)
# !pip install numpy pandas matplotlib seaborn scikit-learn xgboost lightgbm catboost

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.model_selection import train_test_split, cross_val_score, KFold
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
from sklearn.metrics import (
    mean_squared_error, mean_absolute_error, r2_score,
    mean_absolute_percentage_error
)

import warnings
warnings.filterwarnings('ignore')

# Configuracion global de graficos
plt.rcParams['figure.figsize'] = (10, 5)
plt.rcParams['figure.dpi'] = 100
sns.set_style('whitegrid')
print("Entorno listo.")

0.1 Dataset de trabajo#

Usaremos un dataset sintetico de precios de viviendas con multiples variables. Esto permite ilustrar todos los modelos sin depender de archivos externos.

np.random.seed(42)
n = 500

# Variables independientes
metros = np.random.uniform(40, 300, n)
habitaciones = np.random.randint(1, 7, n)
antiguedad = np.random.uniform(0, 50, n)
dist_centro = np.random.uniform(0.5, 30, n)
piso = np.random.randint(1, 15, n)
tiene_garaje = np.random.choice([0, 1], n, p=[0.4, 0.6])

# Variable dependiente con relacion no perfectamente lineal
precio = (
    20000
    + 1500 * metros
    + 8000 * habitaciones
    - 400 * antiguedad
    - 2500 * dist_centro
    + 300 * piso
    + 15000 * tiene_garaje
    + 0.8 * metros**1.3          # componente no lineal
    + np.random.normal(0, 15000, n)  # ruido
)

df = pd.DataFrame({
    'metros': metros,
    'habitaciones': habitaciones,
    'antiguedad': antiguedad,
    'dist_centro': dist_centro,
    'piso': piso,
    'tiene_garaje': tiene_garaje,
    'precio': precio
})

print(f"Shape: {df.shape}")
df.head(10)

df.describe().round(2)

1. Correlacion: el paso previo obligatorio#

La correlacion mide la fuerza y direccion de la relacion entre dos variables. Su valor va de -1 a +1. Antes de construir cualquier modelo de regresion, este analisis es imprescindible por cinco razones:

  1. Validacion de relacion: si no hay correlacion significativa entre X e Y, un modelo lineal sera inutil.

  2. Deteccion de multicolinealidad: si dos predictores tienen correlacion > 0.80 entre si, incluir ambos genera inestabilidad en los coeficientes.

  3. Seleccion de variables: ayuda a priorizar que variables incluir.

  4. Direccion de la relacion: indica si es directa (+) o inversa (-), lo cual debe coincidir con la logica del dominio.

  5. Identificacion de no linealidad: si el scatter plot muestra relacion clara pero Pearson da valor bajo, la relacion no es lineal.

Recordar siempre: correlacion no implica causalidad.

1.1 Coeficiente de Pearson#

Mide la relacion lineal entre dos variables continuas. Requiere normalidad aproximada y ausencia de outliers extremos.

\[r = \frac{\sum (x_i - \bar{x})(y_i - \bar{y})}{\sqrt{\sum (x_i - \bar{x})^2 \cdot \sum (y_i - \bar{y})^2}}\]
# Matriz de correlacion de Pearson
corr_pearson = df.corr(method='pearson')

fig, ax = plt.subplots(figsize=(9, 7))
mask = np.triu(np.ones_like(corr_pearson, dtype=bool))
sns.heatmap(
    corr_pearson, mask=mask, annot=True, fmt='.3f',
    cmap='RdBu_r', center=0, vmin=-1, vmax=1,
    square=True, linewidths=0.5, ax=ax
)
ax.set_title('Matriz de Correlacion de Pearson', fontsize=14, fontweight='bold')
plt.tight_layout()
plt.show()

print("\nCorrelacion de cada variable con el precio:")
print(corr_pearson['precio'].drop('precio').sort_values(ascending=False).round(4))

1.2 Coeficiente de Spearman#

Trabaja con los rangos de los datos. Mide relaciones monotonas, no necesariamente lineales. Es mas robusto frente a outliers.

\[\rho = 1 - \frac{6 \sum d_i^2}{n(n^2 - 1)}\]
# Correlacion de Spearman
corr_spearman = df.corr(method='spearman')

fig, ax = plt.subplots(figsize=(9, 7))
mask = np.triu(np.ones_like(corr_spearman, dtype=bool))
sns.heatmap(
    corr_spearman, mask=mask, annot=True, fmt='.3f',
    cmap='RdBu_r', center=0, vmin=-1, vmax=1,
    square=True, linewidths=0.5, ax=ax
)
ax.set_title('Matriz de Correlacion de Spearman', fontsize=14, fontweight='bold')
plt.tight_layout()
plt.show()

1.3 Coeficiente Tau de Kendall#

Mas robusto que Spearman con muestras pequenas. Se basa en pares concordantes y discordantes.

\[\tau = \frac{\text{pares concordantes} - \text{pares discordantes}}{n(n-1)/2}\]
# Correlacion de Kendall
corr_kendall = df.corr(method='kendall')

fig, ax = plt.subplots(figsize=(9, 7))
mask = np.triu(np.ones_like(corr_kendall, dtype=bool))
sns.heatmap(
    corr_kendall, mask=mask, annot=True, fmt='.3f',
    cmap='RdBu_r', center=0, vmin=-1, vmax=1,
    square=True, linewidths=0.5, ax=ax
)
ax.set_title('Matriz de Correlacion de Kendall', fontsize=14, fontweight='bold')
plt.tight_layout()
plt.show()

1.4 Comparacion visual de los tres coeficientes#

# Comparar los tres coeficientes respecto al precio
comparacion = pd.DataFrame({
    'Pearson': corr_pearson['precio'].drop('precio'),
    'Spearman': corr_spearman['precio'].drop('precio'),
    'Kendall': corr_kendall['precio'].drop('precio')
})

comparacion.plot(kind='bar', figsize=(10, 5), width=0.75)
plt.title('Comparacion de Coeficientes de Correlacion con Precio', fontsize=13, fontweight='bold')
plt.ylabel('Coeficiente')
plt.xlabel('Variable')
plt.xticks(rotation=0)
plt.legend(loc='lower right')
plt.axhline(y=0, color='black', linewidth=0.5)
plt.tight_layout()
plt.show()

1.5 Scatter plots: visualizar la relacion antes de modelar#

fig, axes = plt.subplots(2, 3, figsize=(15, 9))
variables = ['metros', 'habitaciones', 'antiguedad', 'dist_centro', 'piso', 'tiene_garaje']

for ax, var in zip(axes.flat, variables):
    ax.scatter(df[var], df['precio'], alpha=0.3, s=10, color='steelblue')
    ax.set_xlabel(var)
    ax.set_ylabel('precio')
    r = df[var].corr(df['precio'])
    ax.set_title(f'{var} (r = {r:.3f})')

plt.suptitle('Relacion de cada variable con el precio', fontsize=14, fontweight='bold')
plt.tight_layout()
plt.show()

1.6 Guia de interpretacion#

Rango

Fuerza

Implicacion

0.00 – 0.10

Despreciable

No vale la pena modelar

0.10 – 0.30

Debil

Bajo poder predictivo individual

0.30 – 0.50

Moderada

Util combinada con otras variables

0.50 – 0.70

Fuerte

Buena candidata para regresion

0.70 – 1.00

Muy fuerte

Excelente. Si es entre predictores, cuidado con multicolinealidad

2. Preparacion de datos para los modelos#

Separamos variables, dividimos en entrenamiento/prueba y estandarizamos.

X = df.drop('precio', axis=1)
y = df['precio']

X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.2, random_state=42
)

# Estandarizar (necesario para Ridge, Lasso, SVR, KNN)
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

print(f"Entrenamiento: {X_train.shape[0]} muestras")
print(f"Prueba:        {X_test.shape[0]} muestras")

Funciones auxiliares para evaluar modelos#

Definimos funciones reutilizables para evaluar y comparar modelos de forma consistente.

def evaluar_modelo(nombre, y_real, y_pred):
    """Calcula y retorna metricas de regresion."""
    mse  = mean_squared_error(y_real, y_pred)
    rmse = np.sqrt(mse)
    mae  = mean_absolute_error(y_real, y_pred)
    r2   = r2_score(y_real, y_pred)
    mape = mean_absolute_percentage_error(y_real, y_pred) * 100
    return {
        'Modelo': nombre,
        'R2': round(r2, 4),
        'RMSE': round(rmse, 2),
        'MAE': round(mae, 2),
        'MAPE (%)': round(mape, 2)
    }

def graficar_prediccion(y_real, y_pred, titulo):
    """Grafica valores reales vs predichos y residuos."""
    fig, axes = plt.subplots(1, 2, figsize=(13, 5))

    # Real vs Predicho
    axes[0].scatter(y_real, y_pred, alpha=0.4, s=15, color='steelblue')
    lim_min = min(y_real.min(), y_pred.min())
    lim_max = max(y_real.max(), y_pred.max())
    axes[0].plot([lim_min, lim_max], [lim_min, lim_max], 'r--', linewidth=1.5)
    axes[0].set_xlabel('Valor Real')
    axes[0].set_ylabel('Valor Predicho')
    axes[0].set_title(f'{titulo}: Real vs Predicho')

    # Residuos
    residuos = y_real - y_pred
    axes[1].scatter(y_pred, residuos, alpha=0.4, s=15, color='coral')
    axes[1].axhline(y=0, color='black', linewidth=1)
    axes[1].set_xlabel('Valor Predicho')
    axes[1].set_ylabel('Residuo')
    axes[1].set_title(f'{titulo}: Residuos')

    plt.tight_layout()
    plt.show()

# Almacen de resultados
resultados = []
print("Funciones auxiliares definidas.")

3. Tipos de regresion#

3.1 Regresion lineal simple#

El modelo mas basico. Modela la relacion entre una variable independiente y la variable dependiente mediante una linea recta.

\[Y = \beta_0 + \beta_1 X + \varepsilon\]

Se ajusta por Minimos Cuadrados Ordinarios (OLS), minimizando \(\sum (y_i - \hat{y}_i)^2\).

Usamos metros como unico predictor para ilustrar el concepto.

from sklearn.linear_model import LinearRegression

# Regresion simple: solo metros
X_simple_train = X_train[['metros']]
X_simple_test = X_test[['metros']]

modelo_simple = LinearRegression()
modelo_simple.fit(X_simple_train, y_train)
y_pred_simple = modelo_simple.predict(X_simple_test)

print(f"Intercepto (beta_0): {modelo_simple.intercept_:,.2f}")
print(f"Coeficiente metros (beta_1): {modelo_simple.coef_[0]:,.2f}")
print(f"\nInterpretacion: cada metro cuadrado adicional agrega ~{modelo_simple.coef_[0]:,.0f} al precio.")

res = evaluar_modelo('Lineal Simple', y_test, y_pred_simple)
resultados.append(res)
print(f"\nR2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_simple, 'Regresion Lineal Simple')

# Visualizar la linea de regresion
fig, ax = plt.subplots(figsize=(10, 5))
ax.scatter(df['metros'], df['precio'], alpha=0.3, s=10, color='steelblue', label='Datos')

x_line = np.linspace(df['metros'].min(), df['metros'].max(), 100).reshape(-1, 1)
y_line = modelo_simple.predict(x_line)
ax.plot(x_line, y_line, color='red', linewidth=2, label='Regresion Lineal')

ax.set_xlabel('Metros cuadrados')
ax.set_ylabel('Precio')
ax.set_title('Regresion Lineal Simple: Metros vs Precio', fontweight='bold')
ax.legend()
plt.tight_layout()
plt.show()

3.2 Regresion lineal multiple#

Usa dos o mas variables independientes. Cada coeficiente \(\beta_i\) indica el cambio en Y por unidad de \(X_i\) manteniendo constantes las demas.

\[Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_n X_n + \varepsilon\]

Es crucial verificar la multicolinealidad mediante el VIF (Factor de Inflacion de Varianza). Un VIF > 10 es problematico.

# Regresion multiple con todas las variables
modelo_multiple = LinearRegression()
modelo_multiple.fit(X_train, y_train)
y_pred_multiple = modelo_multiple.predict(X_test)

print("Coeficientes del modelo:")
for var, coef in zip(X.columns, modelo_multiple.coef_):
    print(f"  {var:15s}: {coef:>12,.2f}")
print(f"  {'intercepto':15s}: {modelo_multiple.intercept_:>12,.2f}")

res = evaluar_modelo('Lineal Multiple', y_test, y_pred_multiple)
resultados.append(res)
print(f"\nR2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_multiple, 'Regresion Lineal Multiple')

Verificacion de multicolinealidad con VIF

from statsmodels.stats.outliers_influence import variance_inflation_factor

vif_data = pd.DataFrame()
vif_data['Variable'] = X.columns
vif_data['VIF'] = [
    variance_inflation_factor(X.values, i) for i in range(X.shape[1])
]
vif_data = vif_data.sort_values('VIF', ascending=False)

print("Factor de Inflacion de Varianza (VIF)")
print("VIF > 5: preocupante  |  VIF > 10: grave")
print("-" * 35)
print(vif_data.to_string(index=False))

3.3 Regresion polinomial#

Agrega terminos de potencia (\(X^2\), \(X^3\), …) para capturar curvaturas. Tecnicamente sigue siendo lineal en los parametros.

\[Y = \beta_0 + \beta_1 X + \beta_2 X^2 + \cdots + \beta_n X^n + \varepsilon\]

A mayor grado, mayor riesgo de overfitting. Comparamos grados 2 y 3.

from sklearn.pipeline import Pipeline

resultados_poly = {}

for grado in [2, 3]:
    pipe = Pipeline([
        ('poly', PolynomialFeatures(degree=grado, include_bias=False)),
        ('scaler', StandardScaler()),
        ('reg', LinearRegression())
    ])
    pipe.fit(X_train, y_train)
    y_pred = pipe.predict(X_test)

    res = evaluar_modelo(f'Polinomial grado {grado}', y_test, y_pred)
    resultados.append(res)
    resultados_poly[grado] = res
    print(f"Grado {grado}  ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Grafico del mejor
pipe_best = Pipeline([
    ('poly', PolynomialFeatures(degree=2, include_bias=False)),
    ('scaler', StandardScaler()),
    ('reg', LinearRegression())
])
pipe_best.fit(X_train, y_train)
y_pred_poly = pipe_best.predict(X_test)
graficar_prediccion(y_test.values, y_pred_poly, 'Regresion Polinomial (grado 2)')

3.4 Regresion Ridge (Regularizacion L2)#

Agrega una penalizacion proporcional al cuadrado de los coeficientes. Los encoge hacia cero pero nunca los hace exactamente cero.

\[\min \sum (y_i - \hat{y}_i)^2 + \alpha \sum \beta_j^2\]

Es especialmente util cuando hay multicolinealidad o mas variables que observaciones.

from sklearn.linear_model import Ridge, RidgeCV

# Buscar mejor alpha con validacion cruzada
alphas = np.logspace(-2, 4, 100)
ridge_cv = RidgeCV(alphas=alphas, cv=5, scoring='r2')
ridge_cv.fit(X_train_scaled, y_train)

print(f"Mejor alpha: {ridge_cv.alpha_:.4f}")

y_pred_ridge = ridge_cv.predict(X_test_scaled)
res = evaluar_modelo('Ridge (L2)', y_test, y_pred_ridge)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Comparar coeficientes Ridge vs OLS
coef_comp = pd.DataFrame({
    'Variable': X.columns,
    'OLS': modelo_multiple.coef_,
    'Ridge': ridge_cv.coef_
})
print("\nComparacion de coeficientes (datos estandarizados para Ridge):")
print(coef_comp.to_string(index=False))

graficar_prediccion(y_test.values, y_pred_ridge, 'Ridge (L2)')

3.5 Regresion Lasso (Regularizacion L1)#

Penaliza con el valor absoluto de los coeficientes. A diferencia de Ridge, Lasso puede hacer coeficientes exactamente cero, eliminando variables automaticamente.

\[\min \sum (y_i - \hat{y}_i)^2 + \alpha \sum |\beta_j|\]
from sklearn.linear_model import Lasso, LassoCV

lasso_cv = LassoCV(alphas=np.logspace(-2, 4, 100), cv=5, random_state=42, max_iter=10000)
lasso_cv.fit(X_train_scaled, y_train)

print(f"Mejor alpha: {lasso_cv.alpha_:.4f}")

y_pred_lasso = lasso_cv.predict(X_test_scaled)
res = evaluar_modelo('Lasso (L1)', y_test, y_pred_lasso)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Variables seleccionadas (coef != 0)
coef_lasso = pd.DataFrame({
    'Variable': X.columns,
    'Coeficiente': lasso_cv.coef_,
    'Seleccionada': ['Si' if abs(c) > 1e-6 else 'ELIMINADA' for c in lasso_cv.coef_]
})
print("\nSeleccion automatica de variables por Lasso:")
print(coef_lasso.to_string(index=False))

graficar_prediccion(y_test.values, y_pred_lasso, 'Lasso (L1)')

3.6 Regresion Elastic Net#

Combina Ridge (L2) y Lasso (L1). Util cuando hay grupos de variables correlacionadas.

\[\min \sum (y_i - \hat{y}_i)^2 + \alpha \left[ \rho \sum |\beta_j| + (1 - \rho) \sum \beta_j^2 \right]\]

Tiene dos hiperparametros: \(\alpha\) (fuerza total) y \(\rho\) (proporcion L1 vs L2).

from sklearn.linear_model import ElasticNetCV

enet_cv = ElasticNetCV(
    l1_ratio=[0.1, 0.3, 0.5, 0.7, 0.9, 0.95],
    alphas=np.logspace(-2, 4, 50),
    cv=5, random_state=42, max_iter=10000
)
enet_cv.fit(X_train_scaled, y_train)

print(f"Mejor alpha:    {enet_cv.alpha_:.4f}")
print(f"Mejor l1_ratio: {enet_cv.l1_ratio_:.2f}")

y_pred_enet = enet_cv.predict(X_test_scaled)
res = evaluar_modelo('Elastic Net', y_test, y_pred_enet)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_enet, 'Elastic Net')

3.7 Regresion logistica#

A pesar de su nombre, es un modelo de clasificacion. Predice la probabilidad de pertenecer a una clase usando la funcion sigmoide.

\[p = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X_1 + \cdots + \beta_n X_n)}}\]

Creamos una variable binaria: el precio es caro (1) si esta por encima de la mediana, barato (0) si esta por debajo.

from sklearn.linear_model import LogisticRegression
from sklearn.metrics import (
    accuracy_score, precision_score, recall_score,
    f1_score, confusion_matrix, classification_report,
    roc_curve, roc_auc_score
)

# Crear variable binaria
mediana_precio = y.median()
y_bin = (y > mediana_precio).astype(int)
y_bin_train = y_bin.loc[y_train.index]
y_bin_test = y_bin.loc[y_test.index]

log_reg = LogisticRegression(max_iter=1000, random_state=42)
log_reg.fit(X_train_scaled, y_bin_train)

y_pred_log = log_reg.predict(X_test_scaled)
y_prob_log = log_reg.predict_proba(X_test_scaled)[:, 1]

print("=== Regresion Logistica ===")
print(f"Accuracy:  {accuracy_score(y_bin_test, y_pred_log):.4f}")
print(f"Precision: {precision_score(y_bin_test, y_pred_log):.4f}")
print(f"Recall:    {recall_score(y_bin_test, y_pred_log):.4f}")
print(f"F1-Score:  {f1_score(y_bin_test, y_pred_log):.4f}")
print(f"AUC-ROC:   {roc_auc_score(y_bin_test, y_prob_log):.4f}")

print("\nReporte de clasificacion:")
print(classification_report(y_bin_test, y_pred_log, target_names=['Barato', 'Caro']))

# Curva ROC y Matriz de Confusion
fig, axes = plt.subplots(1, 2, figsize=(13, 5))

# Curva ROC
fpr, tpr, _ = roc_curve(y_bin_test, y_prob_log)
auc = roc_auc_score(y_bin_test, y_prob_log)
axes[0].plot(fpr, tpr, color='steelblue', linewidth=2, label=f'AUC = {auc:.4f}')
axes[0].plot([0, 1], [0, 1], 'r--', linewidth=1)
axes[0].set_xlabel('Tasa de Falsos Positivos')
axes[0].set_ylabel('Tasa de Verdaderos Positivos')
axes[0].set_title('Curva ROC', fontweight='bold')
axes[0].legend()

# Matriz de confusion
cm = confusion_matrix(y_bin_test, y_pred_log)
sns.heatmap(cm, annot=True, fmt='d', cmap='Blues', ax=axes[1],
            xticklabels=['Barato', 'Caro'], yticklabels=['Barato', 'Caro'])
axes[1].set_xlabel('Predicho')
axes[1].set_ylabel('Real')
axes[1].set_title('Matriz de Confusion', fontweight='bold')

plt.tight_layout()
plt.show()

3.8 Regresion con arboles de decision#

Divide el espacio de datos en regiones y asigna a cada una el promedio de los valores observados. No requiere supuestos sobre la distribucion ni escalar variables.

Hiperparametros clave: max_depth, min_samples_split, min_samples_leaf.

from sklearn.tree import DecisionTreeRegressor, plot_tree

arbol = DecisionTreeRegressor(max_depth=5, min_samples_leaf=10, random_state=42)
arbol.fit(X_train, y_train)
y_pred_arbol = arbol.predict(X_test)

res = evaluar_modelo('Arbol de Decision', y_test, y_pred_arbol)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_arbol, 'Arbol de Decision')

# Visualizar las primeras ramas del arbol
fig, ax = plt.subplots(figsize=(20, 8))
plot_tree(arbol, feature_names=X.columns, filled=True, rounded=True,
          max_depth=3, fontsize=8, ax=ax)
ax.set_title('Arbol de Decision (primeras 3 capas)', fontweight='bold', fontsize=14)
plt.tight_layout()
plt.show()

# Importancia de variables
importancia = pd.DataFrame({
    'Variable': X.columns,
    'Importancia': arbol.feature_importances_
}).sort_values('Importancia', ascending=True)

fig, ax = plt.subplots(figsize=(8, 4))
ax.barh(importancia['Variable'], importancia['Importancia'], color='steelblue')
ax.set_title('Importancia de Variables - Arbol de Decision', fontweight='bold')
ax.set_xlabel('Importancia')
plt.tight_layout()
plt.show()

3.9 Random Forest para regresion#

Ensemble de multiples arboles entrenados con muestras bootstrap y subconjuntos aleatorios de variables. La prediccion es el promedio de todos los arboles. Reduce la varianza sin aumentar mucho el sesgo.

from sklearn.ensemble import RandomForestRegressor

rf = RandomForestRegressor(
    n_estimators=200, max_depth=10,
    min_samples_leaf=5, random_state=42, n_jobs=-1
)
rf.fit(X_train, y_train)
y_pred_rf = rf.predict(X_test)

res = evaluar_modelo('Random Forest', y_test, y_pred_rf)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Importancia de variables
importancia_rf = pd.DataFrame({
    'Variable': X.columns,
    'Importancia': rf.feature_importances_
}).sort_values('Importancia', ascending=True)

fig, ax = plt.subplots(figsize=(8, 4))
ax.barh(importancia_rf['Variable'], importancia_rf['Importancia'], color='steelblue')
ax.set_title('Importancia de Variables - Random Forest', fontweight='bold')
ax.set_xlabel('Importancia')
plt.tight_layout()
plt.show()

graficar_prediccion(y_test.values, y_pred_rf, 'Random Forest')

3.10 Gradient Boosting para regresion#

Construye arboles secuencialmente. Cada arbol corrige los errores del anterior. Es consistentemente uno de los mejores algoritmos para datos tabulares.

Comparamos tres implementaciones: sklearn GradientBoosting, XGBoost y LightGBM.

from sklearn.ensemble import GradientBoostingRegressor

# sklearn Gradient Boosting
gb_sklearn = GradientBoostingRegressor(
    n_estimators=200, learning_rate=0.1, max_depth=5,
    subsample=0.8, random_state=42
)
gb_sklearn.fit(X_train, y_train)
y_pred_gb = gb_sklearn.predict(X_test)

res = evaluar_modelo('Gradient Boosting (sklearn)', y_test, y_pred_gb)
resultados.append(res)
print(f"sklearn GB  ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_gb, 'Gradient Boosting (sklearn)')

# XGBoost
try:
    from xgboost import XGBRegressor

    xgb = XGBRegressor(
        n_estimators=200, learning_rate=0.1, max_depth=5,
        subsample=0.8, colsample_bytree=0.8,
        random_state=42, verbosity=0
    )
    xgb.fit(X_train, y_train)
    y_pred_xgb = xgb.predict(X_test)

    res = evaluar_modelo('XGBoost', y_test, y_pred_xgb)
    resultados.append(res)
    print(f"XGBoost  ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")
except ImportError:
    print("XGBoost no instalado. Ejecutar: pip install xgboost")

# LightGBM
try:
    from lightgbm import LGBMRegressor

    lgbm = LGBMRegressor(
        n_estimators=200, learning_rate=0.1, max_depth=5,
        subsample=0.8, colsample_bytree=0.8,
        random_state=42, verbosity=-1
    )
    lgbm.fit(X_train, y_train)
    y_pred_lgbm = lgbm.predict(X_test)

    res = evaluar_modelo('LightGBM', y_test, y_pred_lgbm)
    resultados.append(res)
    print(f"LightGBM ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")
except ImportError:
    print("LightGBM no instalado. Ejecutar: pip install lightgbm")

# CatBoost
try:
    from catboost import CatBoostRegressor

    cat = CatBoostRegressor(
        iterations=200, learning_rate=0.1, depth=5,
        random_state=42, verbose=0
    )
    cat.fit(X_train, y_train)
    y_pred_cat = cat.predict(X_test)

    res = evaluar_modelo('CatBoost', y_test, y_pred_cat)
    resultados.append(res)
    print(f"CatBoost ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")
except ImportError:
    print("CatBoost no instalado. Ejecutar: pip install catboost")

3.11 Support Vector Regression (SVR)#

Busca un «tubo» (epsilon-tube) que contenga la mayor cantidad de datos. Solo los puntos fuera del tubo (vectores de soporte) definen el modelo. Requiere escalar las variables.

Kernels: lineal, RBF (el mas comun) y polinomial.

from sklearn.svm import SVR

# SVR con kernel RBF
svr = SVR(kernel='rbf', C=100, epsilon=0.1, gamma='scale')
svr.fit(X_train_scaled, y_train)
y_pred_svr = svr.predict(X_test_scaled)

res = evaluar_modelo('SVR (RBF)', y_test, y_pred_svr)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")
print(f"Vectores de soporte: {svr.n_support_} de {len(y_train)} muestras")

graficar_prediccion(y_test.values, y_pred_svr, 'SVR (RBF)')

3.12 K-Nearest Neighbors para regresion (KNN)#

Predice como el promedio de los K vecinos mas cercanos. No aprende una funcion: memoriza los datos. Requiere estandarizar las variables.

Probamos diferentes valores de K para encontrar el optimo.

from sklearn.neighbors import KNeighborsRegressor

# Buscar mejor K
k_range = range(1, 31)
scores_knn = []
for k in k_range:
    knn_temp = KNeighborsRegressor(n_neighbors=k, weights='distance')
    cv_score = cross_val_score(knn_temp, X_train_scaled, y_train, cv=5, scoring='r2')
    scores_knn.append(cv_score.mean())

mejor_k = list(k_range)[np.argmax(scores_knn)]
print(f"Mejor K: {mejor_k} (R2 CV: {max(scores_knn):.4f})")

# Grafico de K vs R2
fig, ax = plt.subplots(figsize=(10, 4))
ax.plot(k_range, scores_knn, marker='o', markersize=4, color='steelblue')
ax.axvline(x=mejor_k, color='red', linestyle='--', label=f'Mejor K = {mejor_k}')
ax.set_xlabel('K (numero de vecinos)')
ax.set_ylabel('R2 (Validacion Cruzada)')
ax.set_title('Seleccion de K para KNN', fontweight='bold')
ax.legend()
plt.tight_layout()
plt.show()

# Modelo final
knn = KNeighborsRegressor(n_neighbors=mejor_k, weights='distance')
knn.fit(X_train_scaled, y_train)
y_pred_knn = knn.predict(X_test_scaled)

res = evaluar_modelo('KNN', y_test, y_pred_knn)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

graficar_prediccion(y_test.values, y_pred_knn, f'KNN (K={mejor_k})')

3.13 Regresion bayesiana#

Estima una distribucion de probabilidad para cada coeficiente, no un valor unico. Incorpora conocimiento previo (prior) y lo actualiza con los datos (posterior). Scikit-learn ofrece BayesianRidge como implementacion practica.

from sklearn.linear_model import BayesianRidge

bayes = BayesianRidge(compute_score=True)
bayes.fit(X_train_scaled, y_train)
y_pred_bayes, y_std_bayes = bayes.predict(X_test_scaled, return_std=True)

res = evaluar_modelo('Bayesiana (BayesianRidge)', y_test, y_pred_bayes)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Visualizar predicciones con intervalos de confianza
idx_sort = np.argsort(y_pred_bayes)
fig, ax = plt.subplots(figsize=(12, 5))
ax.scatter(range(len(y_test)), y_test.values[idx_sort], s=10, alpha=0.5, label='Real', color='steelblue')
ax.plot(range(len(y_test)), y_pred_bayes[idx_sort], color='red', linewidth=1, label='Prediccion')
ax.fill_between(
    range(len(y_test)),
    y_pred_bayes[idx_sort] - 2 * y_std_bayes[idx_sort],
    y_pred_bayes[idx_sort] + 2 * y_std_bayes[idx_sort],
    alpha=0.2, color='red', label='Intervalo 95%'
)
ax.set_xlabel('Muestra (ordenada por prediccion)')
ax.set_ylabel('Precio')
ax.set_title('Regresion Bayesiana: Predicciones con Incertidumbre', fontweight='bold')
ax.legend()
plt.tight_layout()
plt.show()

3.14 Regresion cuantilica#

Predice cuantiles especificos en lugar de la media. Util para estimar peores escenarios, intervalos de prediccion o cuando la distribucion de Y no es simetrica.

Comparamos la mediana (cuantil 0.5), el percentil 10 y el percentil 90.

from sklearn.linear_model import QuantileRegressor

cuantiles = [0.10, 0.50, 0.90]
predicciones_q = {}

for q in cuantiles:
    qr = QuantileRegressor(quantile=q, alpha=0.1, solver='highs')
    qr.fit(X_train_scaled, y_train)
    predicciones_q[q] = qr.predict(X_test_scaled)
    if q == 0.5:
        res = evaluar_modelo('Cuantilica (mediana)', y_test, predicciones_q[q])
        resultados.append(res)
        print(f"Mediana  ->  R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Visualizar los tres cuantiles
idx_sort = np.argsort(predicciones_q[0.5])
fig, ax = plt.subplots(figsize=(12, 5))
ax.scatter(range(len(y_test)), y_test.values[idx_sort], s=10, alpha=0.4, label='Real', color='steelblue')
ax.plot(range(len(y_test)), predicciones_q[0.50][idx_sort], color='red', linewidth=1.5, label='Mediana (q=0.5)')
ax.fill_between(
    range(len(y_test)),
    predicciones_q[0.10][idx_sort],
    predicciones_q[0.90][idx_sort],
    alpha=0.2, color='orange', label='Intervalo 10%-90%'
)
ax.set_xlabel('Muestra (ordenada por prediccion mediana)')
ax.set_ylabel('Precio')
ax.set_title('Regresion Cuantilica: Mediana e Intervalo 10%-90%', fontweight='bold')
ax.legend()
plt.tight_layout()
plt.show()

3.15 Regresion isotonica#

Ajusta una funcion escalonada monotona a los datos. No asume forma parametrica, solo que la relacion tiene una direccion consistente. Solo funciona con una variable predictora.

from sklearn.isotonic import IsotonicRegression

# Usamos metros como unica variable
x_iso_train = X_train['metros'].values
x_iso_test = X_test['metros'].values

iso = IsotonicRegression(increasing=True)
iso.fit(x_iso_train, y_train.values)

# Para prediccion, clipeamos al rango de entrenamiento
x_iso_test_clip = np.clip(x_iso_test, x_iso_train.min(), x_iso_train.max())
y_pred_iso = iso.predict(x_iso_test_clip)

res = evaluar_modelo('Isotonica', y_test, y_pred_iso)
resultados.append(res)
print(f"R2: {res['R2']}  |  RMSE: {res['RMSE']}  |  MAE: {res['MAE']}")

# Visualizacion
fig, ax = plt.subplots(figsize=(10, 5))
ax.scatter(df['metros'], df['precio'], alpha=0.2, s=10, color='steelblue', label='Datos')

x_plot = np.linspace(x_iso_train.min(), x_iso_train.max(), 500)
y_plot = iso.predict(x_plot)
ax.plot(x_plot, y_plot, color='red', linewidth=2, label='Isotonica')

ax.set_xlabel('Metros cuadrados')
ax.set_ylabel('Precio')
ax.set_title('Regresion Isotonica: Metros vs Precio', fontweight='bold')
ax.legend()
plt.tight_layout()
plt.show()

4. Metricas de evaluacion#

Cada metrica captura un aspecto diferente del rendimiento:

Metrica

Formula

Interpretacion

\(R^2\)

\(1 - \frac{\sum(y_i - \hat{y}_i)^2}{\sum(y_i - \bar{y})^2}\)

Proporcion de varianza explicada (0 a 1)

\(R^2_{adj}\)

Penaliza por numero de variables

Para comparar modelos con distinto numero de variables

MSE

\(\frac{1}{n}\sum(y_i - \hat{y}_i)^2\)

Error cuadratico medio. Penaliza errores grandes

RMSE

\(\sqrt{MSE}\)

Misma escala que Y. La mas usada

MAE

\(\frac{1}{n}\sum|y_i - \hat{y}_i|\)

Menos sensible a outliers que RMSE

MAPE

\(\frac{100}{n}\sum\frac{|y_i - \hat{y}_i|}{|y_i|}\)

Error como porcentaje. Facil de comunicar

 
# Ejemplo practico de como cada metrica reacciona diferente
np.random.seed(0)
y_ejemplo = np.array([100, 200, 300, 400, 500])
y_pred_bueno = y_ejemplo + np.array([5, -8, 12, -3, 7])
y_pred_outlier = y_ejemplo + np.array([5, -8, 12, -3, 200])  # un error enorme

metricas = ['MSE', 'RMSE', 'MAE']
for nombre_pred, y_p in [('Sin outlier', y_pred_bueno), ('Con outlier', y_pred_outlier)]:
    mse = mean_squared_error(y_ejemplo, y_p)
    rmse = np.sqrt(mse)
    mae = mean_absolute_error(y_ejemplo, y_p)
    print(f"{nombre_pred:12s} -> MSE: {mse:>10.1f}  |  RMSE: {rmse:>8.1f}  |  MAE: {mae:>6.1f}")

print("\nObservar como MSE y RMSE se disparan con el outlier, mientras MAE es mas estable.")

5. Trade-off sesgo-varianza#

El error total de un modelo se descompone en:

\[\text{Error total} = \text{Sesgo}^2 + \text{Varianza} + \text{Error irreducible}\]

  • Sesgo alto (underfitting): el modelo es demasiado simple para capturar la realidad.

  • Varianza alta (overfitting): el modelo se ajusta al ruido y falla con datos nuevos.

El objetivo es encontrar el punto donde la suma de ambos es minima.

# Demostrar sesgo-varianza con arboles de distinta profundidad
profundidades = range(1, 25)
train_scores = []
test_scores = []

for d in profundidades:
    tree = DecisionTreeRegressor(max_depth=d, random_state=42)
    tree.fit(X_train, y_train)
    train_scores.append(r2_score(y_train, tree.predict(X_train)))
    test_scores.append(r2_score(y_test, tree.predict(X_test)))

fig, ax = plt.subplots(figsize=(10, 5))
ax.plot(profundidades, train_scores, 'o-', label='Entrenamiento', color='steelblue', markersize=4)
ax.plot(profundidades, test_scores, 'o-', label='Prueba', color='coral', markersize=4)
ax.axvline(x=list(profundidades)[np.argmax(test_scores)], color='green',
           linestyle='--', label=f'Mejor profundidad = {list(profundidades)[np.argmax(test_scores)]}')
ax.set_xlabel('Profundidad del arbol')
ax.set_ylabel('R2')
ax.set_title('Trade-off Sesgo-Varianza: Arbol de Decision', fontweight='bold')
ax.legend()
ax.set_ylim(0, 1.05)
plt.tight_layout()
plt.show()

print("A la izquierda del optimo: underfitting (alto sesgo)")
print("A la derecha del optimo: overfitting (alta varianza)")

6. Comparacion final de todos los modelos#

# Tabla comparativa
df_resultados = pd.DataFrame(resultados)
df_resultados = df_resultados.sort_values('R2', ascending=False).reset_index(drop=True)
df_resultados.index += 1
print("=== RANKING DE MODELOS POR R2 ===\n")
print(df_resultados.to_string())

# Grafico comparativo
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# R2
df_plot = df_resultados.sort_values('R2', ascending=True)
colors = ['#2ecc71' if r > 0.9 else '#3498db' if r > 0.8 else '#e74c3c' for r in df_plot['R2']]
axes[0].barh(df_plot['Modelo'], df_plot['R2'], color=colors)
axes[0].set_xlabel('R2')
axes[0].set_title('R2 por Modelo (mayor es mejor)', fontweight='bold')
axes[0].set_xlim(0, 1)

# RMSE
df_plot2 = df_resultados.sort_values('RMSE', ascending=False)
colors2 = ['#2ecc71' if r < df_resultados['RMSE'].median() else '#e74c3c' for r in df_plot2['RMSE']]
axes[1].barh(df_plot2['Modelo'], df_plot2['RMSE'], color=colors2)
axes[1].set_xlabel('RMSE')
axes[1].set_title('RMSE por Modelo (menor es mejor)', fontweight='bold')

plt.tight_layout()
plt.show()

7. Validacion cruzada y buenas practicas#

7.1 Validacion cruzada K-Fold#

Se dividen los datos en K partes. Se entrena con K-1 y se valida con la restante. Se rota K veces y se promedian los resultados. Valores tipicos de K: 5 o 10.

from sklearn.ensemble import RandomForestRegressor, GradientBoostingRegressor

modelos_cv = {
    'Lineal': LinearRegression(),
    'Ridge': Ridge(alpha=1.0),
    'Lasso': Lasso(alpha=1.0, max_iter=10000),
    'Arbol': DecisionTreeRegressor(max_depth=5, random_state=42),
    'Random Forest': RandomForestRegressor(n_estimators=100, max_depth=8, random_state=42),
    'Gradient Boosting': GradientBoostingRegressor(n_estimators=100, max_depth=5, random_state=42),
}

cv_resultados = []
kf = KFold(n_splits=5, shuffle=True, random_state=42)

for nombre, modelo in modelos_cv.items():
    # Usar datos escalados para los que lo necesitan
    if nombre in ['Ridge', 'Lasso']:
        scores = cross_val_score(modelo, X_train_scaled, y_train, cv=kf, scoring='r2')
    else:
        scores = cross_val_score(modelo, X_train, y_train, cv=kf, scoring='r2')

    cv_resultados.append({
        'Modelo': nombre,
        'R2 Media': round(scores.mean(), 4),
        'R2 Std': round(scores.std(), 4),
        'R2 Min': round(scores.min(), 4),
        'R2 Max': round(scores.max(), 4)
    })

df_cv = pd.DataFrame(cv_resultados).sort_values('R2 Media', ascending=False)
df_cv.index = range(1, len(df_cv) + 1)
print("=== VALIDACION CRUZADA 5-FOLD ===\n")
print(df_cv.to_string())

# Boxplot de los scores de CV
fig, ax = plt.subplots(figsize=(10, 5))

cv_data = []
for nombre, modelo in modelos_cv.items():
    if nombre in ['Ridge', 'Lasso']:
        scores = cross_val_score(modelo, X_train_scaled, y_train, cv=kf, scoring='r2')
    else:
        scores = cross_val_score(modelo, X_train, y_train, cv=kf, scoring='r2')
    cv_data.append(scores)

bp = ax.boxplot(cv_data, labels=list(modelos_cv.keys()), patch_artist=True)
for patch in bp['boxes']:
    patch.set_facecolor('steelblue')
    patch.set_alpha(0.6)
ax.set_ylabel('R2')
ax.set_title('Distribucion de R2 por Modelo (5-Fold CV)', fontweight='bold')
plt.xticks(rotation=15)
plt.tight_layout()
plt.show()

7.2 Preprocesamiento: escalamiento de variables#

Obligatorio para SVR, KNN, Ridge, Lasso. No necesario para arboles.

from sklearn.preprocessing import StandardScaler, MinMaxScaler, RobustScaler

escaladores = {
    'StandardScaler': StandardScaler(),
    'MinMaxScaler': MinMaxScaler(),
    'RobustScaler': RobustScaler()
}

muestra = X_train[['metros', 'antiguedad', 'dist_centro']].head(5)
print("Datos originales:")
print(muestra.round(2).to_string())

for nombre, esc in escaladores.items():
    transformado = pd.DataFrame(
        esc.fit_transform(X_train[['metros', 'antiguedad', 'dist_centro']]),
        columns=['metros', 'antiguedad', 'dist_centro']
    ).head(5)
    print(f"\n{nombre}:")
    print(transformado.round(4).to_string())

7.3 Diagnostico de overfitting vs underfitting#

  • Overfitting: error de entrenamiento muy bajo, error de prueba alto. Soluciones: mas datos, regularizacion, reducir complejidad.

  • Underfitting: error alto en ambos. Soluciones: modelo mas complejo, mas features, reducir regularizacion.

# Curvas de aprendizaje para diagnosticar
from sklearn.model_selection import learning_curve

fig, axes = plt.subplots(1, 3, figsize=(16, 4.5))
modelos_lc = [
    ('Lineal (alto sesgo)', LinearRegression()),
    ('Random Forest', RandomForestRegressor(n_estimators=50, max_depth=8, random_state=42)),
    ('Arbol sin poda (alta varianza)', DecisionTreeRegressor(random_state=42))
]

for ax, (nombre, modelo) in zip(axes, modelos_lc):
    train_sizes, train_scores, test_scores = learning_curve(
        modelo, X_train, y_train, cv=5,
        train_sizes=np.linspace(0.1, 1.0, 10),
        scoring='r2', n_jobs=-1
    )
    ax.plot(train_sizes, train_scores.mean(axis=1), 'o-', label='Entrenamiento', color='steelblue', markersize=4)
    ax.plot(train_sizes, test_scores.mean(axis=1), 'o-', label='Validacion', color='coral', markersize=4)
    ax.fill_between(train_sizes,
                     train_scores.mean(axis=1) - train_scores.std(axis=1),
                     train_scores.mean(axis=1) + train_scores.std(axis=1),
                     alpha=0.1, color='steelblue')
    ax.fill_between(train_sizes,
                     test_scores.mean(axis=1) - test_scores.std(axis=1),
                     test_scores.mean(axis=1) + test_scores.std(axis=1),
                     alpha=0.1, color='coral')
    ax.set_xlabel('Muestras de entrenamiento')
    ax.set_ylabel('R2')
    ax.set_title(nombre, fontweight='bold')
    ax.legend(fontsize=8)
    ax.set_ylim(0, 1.05)

plt.suptitle('Curvas de Aprendizaje: Diagnostico de Sesgo-Varianza', fontweight='bold', fontsize=13)
plt.tight_layout()
plt.show()

print("Lineal: las curvas convergen pero a un nivel bajo (alto sesgo / underfitting)")
print("Random Forest: buen balance entre entrenamiento y validacion")
print("Arbol sin poda: gran brecha entre entrenamiento y validacion (alta varianza / overfitting)")

8. Guia de seleccion de modelo#

Escenario

Modelo recomendado

Razon

Primer modelo (baseline)

Regresion lineal

Rapido, interpretable

Relacion no lineal

Polinomial o arboles

Capturan curvaturas

Muchas variables con multicolinealidad

Ridge o Elastic Net

Estabilizan coeficientes

Muchas variables, pocas relevantes

Lasso o Elastic Net

Seleccion automatica

Maximo rendimiento predictivo

Gradient Boosting

Generalmente el mejor en datos tabulares

Robustez sin mucho tuning

Random Forest

Buen rendimiento por defecto

Dataset pequeno

SVR o Bayesiana

Funcionan bien con pocas muestras

Cuantificar incertidumbre

Bayesiana o cuantilica

Dan distribuciones, no solo puntos

Clasificacion binaria

Logistica

Probabilidades como salida

Flujo de trabajo recomendado#

  1. Analisis exploratorio y correlacion.

  2. Empezar con regresion lineal como baseline.

  3. Probar regularizacion si hay multicolinealidad o muchas variables.

  4. Probar modelos de ensemble si la relacion parece no lineal.

  5. Comparar modelos con validacion cruzada.

  6. Ajustar hiperparametros del mejor candidato.

  7. Evaluar en conjunto de prueba final (una sola vez).

Fin de la Unidad 3.3 — Modelos de Regresion en Machine Learning

Contenido