Assignment UT3-10: PCA y Feature Selection¶

🔁 CONTINUACIÓN DE TAREAS ANTERIORES: Este assignment es la tercera parte de una serie usando Ames Housing: - UT3-08: Feature Engineering - Creaste features derivadas, ratios, transformaciones - UT3-09: Encoding Avanzado - Aplicaste diferentes técnicas de encoding (aunque usaste Adult Income) - UT3-10 (ESTA): PCA y Feature Selection - Reducción dimensional y selección de features

🎯 Objetivos¶

Implementar PCA y analizar componentes principales
Aplicar Feature Selection con múltiples métodos
Comparar PCA vs Feature Selection
Evaluar trade-offs en contexto de negocio

⏱️ Tiempo Estimado: 90-110 minutos¶

💡 OPTIMIZACIÓN: Esta tarea asume que ya conoces el dataset Ames Housing de tareas anteriores. El preprocesamiento es mínimo (una función de ~30 líneas) para enfocarnos directamente en PCA y Feature Selection, que son los temas principales.

⚠️ NOTA IMPORTANTE: Los wrapper methods (Forward, Backward, RFE) son lentos (~2-3 min cada uno). Ten paciencia, es normal.

📊 Dataset: Ames Housing (Resumen Rápido)¶

Fuente: Kaggle / Dean De Cock (Iowa State University) URL: https://www.kaggle.com/c/house-prices-advanced-regression-techniques

📁 Si ya tienes el dataset: Perfecto, úsalo directamente.
📥 Si no lo tienes: Descarga train.csv del link de arriba o usa Kaggle API: kaggle competitions download -c house-prices-advanced-regression-techniques

Contexto de negocio (resumen): Eres data scientist en una empresa de bienes raíces (real estate) que necesita predecir precios de casas con precisión. La empresa tiene 80+ características de cada propiedad (desde calidad de cocina hasta año de construcción), y necesitas:

Identificar qué características realmente importan para el precio de venta
Reducir la complejidad del modelo para que sea más rápido y mantenible
Explicar a agentes inmobiliarios qué factores considerar al tasar una propiedad
Evitar overfitting eliminando features redundantes o irrelevantes

Descripción del Dataset:

~2900 casas vendidas en Ames, Iowa (2006-2010)
Target: SalePrice - Precio de venta en dólares (regresión continua)
~80 features divididas en:
Dimensiones: LotArea, GrLivArea, TotalBsmtSF, GarageArea, 1stFlrSF, 2ndFlrSF
Calidad: OverallQual, OverallCond, KitchenQual, ExterQual, BsmtQual
Temporales: YearBuilt, YearRemodAdd, GarageYrBlt, YrSold
Categóricas: Neighborhood, HouseStyle, RoofStyle, Exterior, Foundation (~40 categóricas)
Numéricas discretas: BedroomAbvGr, FullBath, Fireplaces, GarageCars
Booleanas: CentralAir, PavedDrive, Street

Parte 1: Setup Rápido - Dataset Ames Housing (10 min)¶

⏱️ Tiempo esperado: ~10 minutos

💡 NOTA: Este dataset ya lo usaste en tareas anteriores (Feature Engineering, Encoding). Si ya tienes el dataset descargado y procesado, esta parte te tomará solo 5 minutos.

Paso 1.1: Setup Rápido - Cargar y Preprocesar¶

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import LabelEncoder
from sklearn.impute import SimpleImputer
import warnings
warnings.filterwarnings('ignore')

# Configuración
np.random.seed(42)
plt.style.use('default')
sns.set_palette("husl")

print("=" * 80)
print("ASSIGNMENT UT3-10: PCA y Feature Selection - Ames Housing Dataset")
print("=" * 80)

# ========== FUNCIÓN DE CARGA Y PREPROCESAMIENTO RÁPIDO ==========
def quick_load_and_preprocess_ames(filepath='train.csv'):
    """
    Carga y preprocesa Ames Housing en un solo paso
    (Ya hiciste esto en tareas anteriores, aquí es versión simplificada)
    """
    print("\n🏠 CARGANDO Y PREPROCESANDO AMES HOUSING...")

    # Cargar dataset
    df = pd.read_csv(filepath)
    print(f"✅ Dataset cargado: {df.shape[0]:,} casas, {df.shape[1]} columnas")

    # Eliminar 'Id' (no predictivo)
    df = df.drop('Id', axis=1, errors='ignore')

    # Identificar tipos de variables
    numerical_cols = df.select_dtypes(include=['int64', 'float64']).columns.tolist()
    categorical_cols = df.select_dtypes(include=['object']).columns.tolist()

    # Remover target de features
    if 'SalePrice' in numerical_cols:
        numerical_cols.remove('SalePrice')

    print(f"✅ Features numéricas: {len(numerical_cols)}")
    print(f"✅ Features categóricas: {len(categorical_cols)}")

    # Imputar valores faltantes
    num_imputer = SimpleImputer(strategy='median')
    df[numerical_cols] = num_imputer.fit_transform(df[numerical_cols])

    cat_imputer = SimpleImputer(strategy='most_frequent')
    df[categorical_cols] = cat_imputer.fit_transform(df[categorical_cols])

    print(f"✅ Missing values imputados")

    # Label encoding para categóricas
    le = LabelEncoder()
    for col in categorical_cols:
        df[col] = le.fit_transform(df[col].astype(str))

    print(f"✅ Categóricas encoded")

    # Separar X y y
    X = df.drop('SalePrice', axis=1)
    y = df['SalePrice']

    print(f"\n✅ DATASET LISTO:")
    print(f"   X shape: {X.shape} ({X.shape[1]} features)")
    print(f"   y shape: {y.shape}")
    print(f"   Precio promedio: ${y.mean():,.0f}")
    print(f"   Precio mediana: ${y.median():,.0f}")

    return X, y, X.columns.tolist()

# ========== EJECUTAR CARGA RÁPIDA ==========
# TODO: Completa con el path correcto
X, y, feature_names = quick_load_and_preprocess_ames(______)

print(f"\n📊 RESUMEN DEL DATASET:")
print(f"   Total features: {X.shape[1]}")
print(f"   Total casas: {X.shape[0]:,}")
print(f"   Ejemplos de features: {feature_names[:10]}")

💡 Pista: Pasa el nombre del archivo CSV como string. Debe estar en el mismo directorio que tu notebook/script.
📖 Documentación: pandas.read_csv()

📚 RECORDATORIO DE TAREAS ANTERIORES:

Si necesitas refrescar conceptos de preprocesamiento: - Tarea 08 (Feature Engineering): Creación de features derivadas, ratios, transformaciones - Tarea 09 (Encoding Avanzado): Label, One-Hot, Target Encoding - Esta tarea (10): Nos enfocamos en reducción dimensional (PCA) y feature selection

🤔 Pregunta Reflexiva 1:

¿Con 80+ features, esperarías que todas sean igualmente importantes para predecir precio?
¿Qué problemas puede causar tener tantas features? (Piensa en: overfitting, velocidad, interpretabilidad)
¿Conoces la diferencia entre PCA (transformar features) y Feature Selection (seleccionar features)?

Parte 2: PCA - Análisis de Componentes Principales (35 min)¶

⏱️ Tiempo esperado: ~35 minutos (PCA con 80 features toma más tiempo)

💡 Contexto de Negocio: Con 80+ features, el modelo es complejo, lento de entrenar, y difícil de interpretar. PCA te permitirá reducir a 10-15 componentes mientras mantienes 80-90% de la información. Esto es crítico para: - ⚡ Velocidad: Modelos más rápidos de entrenar y predecir - 🧠 Prevenir overfitting: Menos features = menos riesgo de sobreajuste - 💰 Costos computacionales: Menos features = menos recursos necesarios

Paso 2.1: Estandarización (Crítico para PCA)¶

from sklearn.preprocessing import StandardScaler

# ========== ESTANDARIZACIÓN ==========
print("=== ESTANDARIZACIÓN DE FEATURES ===")
print("⚠️ PCA es sensible a escala. SIEMPRE estandarizar antes de PCA.")

# TODO: Estandarizar datos
scaler = StandardScaler()
X_scaled = ______

# Verificar estandarización: mean ≈ 0, std ≈ 1
print(f"\n✅ Mean después de scaling: {X_scaled.mean():.6f} (esperado: ~0)")
print(f"✅ Std después de scaling: {X_scaled.std():.6f} (esperado: ~1)")

# Verificar shape
print(f"✅ X_scaled shape: {______}")

# Comparar antes vs después
print(f"\n=== COMPARACIÓN ANTES vs DESPUÉS ===")
print(f"Antes - Mean GrLivArea: {X['GrLivArea'].mean():.0f}, Std: {X['GrLivArea'].std():.0f}")
print(f"Después - Mean GrLivArea: {X_scaled[:, X.columns.get_loc('GrLivArea')].mean():.6f}, Std: {X_scaled[:, X.columns.get_loc('GrLivArea')].std():.6f}")

💡 Pistas:

Primer blank: Usa el método del scaler que transforma los datos (combina fit y transform).
Segundo blank: Accede al atributo que contiene las dimensiones del array resultante.

📖 Documentación: StandardScaler

🤔 Pregunta Reflexiva 2:

¿Por qué PCA requiere estandarización? ¿Qué pasaría si no estandarizas?
¿Qué features crees que dominarían PC1 si NO estandarizamos (hint: piensa en escala)?
¿Hay algún caso donde NO deberías estandarizar antes de PCA?

Paso 2.2: Aplicar PCA Completo (80 Componentes)¶

from sklearn.decomposition import PCA
import time

# ========== APLICAR PCA SIN RESTRICCIONES ==========
print("\n=== APLICANDO PCA ===")
print("⏱️ Esto puede tomar 10-20 segundos con 80 features...")

start_time = time.time()

# TODO: Aplicar PCA sin restricción de componentes
pca = PCA()  # Sin n_components = todos los componentes posibles
X_pca = ______

elapsed_time = time.time() - start_time
print(f"✅ PCA completado en {elapsed_time:.2f} segundos")

# ========== ANALIZAR VARIANZA EXPLICADA ==========
explained_variance = pca.explained_variance_ratio_
cumulative_variance = np.cumsum(explained_variance)

print(f"\n=== ANÁLISIS DE COMPONENTES PRINCIPALES ===")
print(f"Total de componentes generados: {______}")
print(f"\nVarianza explicada por componentes principales:")
print(f"  PC1: {explained_variance[0]:.3%} (¡la más importante!)")
print(f"  PC2: {explained_variance[1]:.3%}")
print(f"  PC3: {explained_variance[2]:.3%}")
print(f"  PC4: {explained_variance[3]:.3%}")
print(f"  PC5: {explained_variance[4]:.3%}")

# TODO: Top 10 componentes
print("\n=== TOP 10 COMPONENTES ===")
for i in range(min(10, len(explained_variance))):
    print(f"PC{i+1:2d}: Individual {explained_variance[i]:6.3%} | Acumulada {cumulative_variance[i]:6.3%}")

💡 Pistas:

Primer blank: Usa el método del objeto pca que ajusta y transforma en un solo paso.
Segundo blank: Accede al atributo del objeto pca que guarda el número de componentes resultantes.

📖 Documentación: PCA en scikit-learn

🤔 Pregunta Reflexiva 3:

¿PC1 captura más varianza que en datasets típicos? ¿Por qué crees que pasa esto?
¿Cuántos componentes capturan ~50% de la varianza? ¿Te sorprende?
Si tuvieras que explicar PC1 a un agente inmobiliario, ¿cómo lo harías?

Paso 2.3: Scree Plot y Decisión de Dimensionalidad¶

# ========== CREAR SCREE PLOT ==========
print("\n=== SCREE PLOT: VISUALIZACIÓN DE VARIANZA ===")

# TODO: Crear scree plot con 80 componentes
plt.figure(figsize=(16, 6))

# Subplot 1: Varianza individual (primeros 30 componentes para claridad)
plt.subplot(1, 2, 1)
n_to_show = min(30, len(explained_variance))
plt.bar(range(1, n_to_show + 1), explained_variance[:n_to_show], alpha=0.7, color='steelblue')
plt.xlabel('Componente Principal', fontsize=12)
plt.ylabel('Varianza Explicada (Individual)', fontsize=12)
plt.title(f'Scree Plot - Primeros {n_to_show} Componentes', fontsize=14)
plt.grid(axis='y', alpha=0.3)

# Subplot 2: Varianza acumulada (TODOS los componentes)
plt.subplot(1, 2, 2)
plt.plot(range(1, len(cumulative_variance) + 1), cumulative_variance, 'o-', 
         color='steelblue', markersize=4, linewidth=2)

# Líneas de referencia
plt.axhline(y=______, color='r', linestyle='--', label='80% varianza', linewidth=2)
plt.axhline(y=______, color='g', linestyle='--', label='90% varianza', linewidth=2)
plt.axhline(y=______, color='orange', linestyle='--', label='95% varianza', linewidth=2)

plt.xlabel('Número de Componentes', fontsize=12)
plt.ylabel('Varianza Acumulada', fontsize=12)
plt.title('Varianza Acumulada por Número de Componentes', fontsize=14)
plt.legend(fontsize=11)
plt.grid(alpha=0.3)
plt.ylim([0, 1.05])

plt.tight_layout()
plt.show()

# ========== DECISIÓN DE DIMENSIONALIDAD ==========
print("\n=== DECISIÓN: ¿CUÁNTOS COMPONENTES NECESITAMOS? ===")

# TODO: Calcular componentes necesarios para diferentes umbrales
n_components_80 = ______
n_components_90 = ______
n_components_95 = ______

print(f"📊 Para 80% de varianza: {n_components_80} componentes")
print(f"📊 Para 90% de varianza: {n_components_90} componentes")
print(f"📊 Para 95% de varianza: {n_components_95} componentes")

# Análisis de reducción dimensional
original_features = X.shape[1]
reduction_80 = (1 - n_components_80 / original_features) * 100
reduction_90 = (1 - n_components_90 / original_features) * 100
reduction_95 = (1 - n_components_95 / original_features) * 100

print(f"\n🎯 IMPACTO DE REDUCCIÓN DIMENSIONAL:")
print(f"   Original: {original_features} features")
print(f"   80% varianza: {original_features} → {n_components_80} ({reduction_80:.1f}% reducción)")
print(f"   90% varianza: {original_features} → {n_components_90} ({reduction_90:.1f}% reducción)")
print(f"   95% varianza: {original_features} → {n_components_95} ({reduction_95:.1f}% reducción)")

print(f"\n💡 RECOMENDACIÓN PRÁCTICA:")
print(f"   Para este assignment, usaremos {n_components_80} componentes (80% varianza)")
print(f"   Esto balancea reducción dimensional con retención de información")

💡 Pistas:

Para las líneas de referencia: Usa valores decimales (0.80 = 80%, 0.90 = 90%, 0.95 = 95%).
Para calcular componentes necesarios: Usa np.argmax() para encontrar el índice donde la varianza acumulada alcanza o supera el umbral, y suma 1 (los índices empiezan en 0).

📖 Documentación: matplotlib.pyplot.axhline | numpy.argmax

🤔 Pregunta Reflexiva 4:

¿Con 80% de varianza reduces a cuántas dimensiones? ¿Es una reducción significativa?
¿Preferirías 80%, 90% o 95% de varianza para un modelo de producción? ¿Por qué?
¿El trade-off entre simplicidad y información retenida vale la pena en este caso?

Paso 2.4: Interpretación de Loadings (¿Qué representa cada PC?)¶

# ========== OBTENER LOADINGS ==========
print("\n=== INTERPRETACIÓN DE COMPONENTES PRINCIPALES ===")
print("Los loadings te dicen qué features originales contribuyen a cada componente")

# TODO: Obtener loadings de PC1 y PC2
loadings = pca.components_.T * np.sqrt(pca.explained_variance_)

# Crear DataFrame de loadings para PC1 y PC2
loadings_df = pd.DataFrame(
    loadings[:, :2],
    columns=['PC1', 'PC2'],
    index=X.columns
)

# ========== FEATURES MÁS IMPORTANTES PARA PC1 ==========
print("\n=== PC1: COMPONENTE PRINCIPAL #1 ===")
print(f"Explica {explained_variance[0]:.1%} de la varianza total")
print(f"\nTop 10 features más importantes para PC1:")
pc1_top = loadings_df['PC1'].abs().nlargest(10)
for i, (feature, loading_abs) in enumerate(pc1_top.items(), 1):
    loading_val = loadings_df.loc[feature, 'PC1']
    direction = "↑ positivo" if loading_val > 0 else "↓ negativo"
    print(f"  {i:2d}. {feature:20s}: {loading_val:+7.3f} ({direction})")

# ========== FEATURES MÁS IMPORTANTES PARA PC2 ==========
print("\n=== PC2: COMPONENTE PRINCIPAL #2 ===")
print(f"Explica {explained_variance[1]:.1%} de la varianza total")
print(f"\nTop 10 features más importantes para PC2:")
pc2_top = loadings_df['PC2'].abs().nlargest(10)
for i, (feature, loading_abs) in enumerate(pc2_top.items(), 1):
    loading_val = loadings_df.loc[feature, 'PC2']
    direction = "↑ positivo" if loading_val > 0 else "↓ negativo"
    print(f"  {i:2d}. {feature:20s}: {loading_val:+7.3f} ({direction})")

# ========== VISUALIZAR LOADINGS (solo top features para claridad) ==========
print("\n=== VISUALIZACIÓN DE LOADINGS ===")

# Seleccionar top features para visualizar (top 15 de PC1 o PC2)
top_features_pc1 = set(loadings_df['PC1'].abs().nlargest(15).index)
top_features_pc2 = set(loadings_df['PC2'].abs().nlargest(15).index)
top_features = list(top_features_pc1.union(top_features_pc2))

loadings_df_viz = loadings_df.loc[top_features]

plt.figure(figsize=(14, 10))
plt.scatter(loadings_df_viz['PC1'], loadings_df_viz['PC2'], alpha=0.7, s=150, c='steelblue', edgecolors='black')

# Anotar features
for feature in loadings_df_viz.index:
    plt.annotate(feature, 
                (loadings_df_viz.loc[feature, 'PC1'], loadings_df_viz.loc[feature, 'PC2']), 
                fontsize=10, alpha=0.9, ha='center')

plt.xlabel(f'PC1 Loadings ({explained_variance[0]:.1%} varianza explicada)', fontsize=13)
plt.ylabel(f'PC2 Loadings ({explained_variance[1]:.1%} varianza explicada)', fontsize=13)
plt.title('Loadings Plot - Top Features en PC1 y PC2\n(Features más influyentes en componentes principales)', fontsize=14)
plt.axhline(y=0, color='k', linestyle='-', alpha=0.3, linewidth=1)
plt.axvline(x=0, color='k', linestyle='-', alpha=0.3, linewidth=1)
plt.grid(alpha=0.3)
plt.tight_layout()
plt.show()

# ========== INTERPRETACIÓN DE NEGOCIO ==========
print("\n=== 💼 INTERPRETACIÓN DE NEGOCIO ===")
print("\nPC1 representa principalmente (top 3):")
for i, (feature, loading_abs) in enumerate(loadings_df['PC1'].abs().nlargest(3).items(), 1):
    loading_val = loadings_df.loc[feature, 'PC1']
    direction = "positivamente" if loading_val > 0 else "negativamente"
    print(f"  {i}. {feature}: Contribuye {direction}")

print("\n💡 Interpretación: PC1 probablemente captura el 'tamaño/calidad general' de la casa")

print("\nPC2 representa principalmente (top 3):")
for i, (feature, loading_abs) in enumerate(loadings_df['PC2'].abs().nlargest(3).items(), 1):
    loading_val = loadings_df.loc[feature, 'PC2']
    direction = "positivamente" if loading_val > 0 else "negativamente"
    print(f"  {i}. {feature}: Contribuye {direction}")

print("\n💡 Interpretación: PC2 probablemente captura otra dimensión (edad, ubicación, etc.)")

🤔 Pregunta Reflexiva 5:

¿Las features que dominan PC1 tienen sentido para predecir precio de una casa?
¿PC1 y PC2 son interpretables o son "cajas negras"? ¿Por qué importa esto?
¿Cómo le explicarías a un agente inmobiliario qué significa PC1 en términos de negocio?

Paso 2.5: Feature Selection Basada en PCA Loadings (CRÍTICO)¶

💡 CONCEPTO CLAVE: En lugar de usar PC1, PC2... como features, vamos a identificar las features ORIGINALES que más contribuyen a los PCs importantes, y usar ESAS features originales.

# ========== FEATURE SELECTION BASADA EN PCA LOADINGS ==========
print("\n=== FEATURE SELECTION BASADA EN PCA LOADINGS ===")
print("💡 En lugar de usar PC1, PC2... usaremos las features ORIGINALES")
print("   que tienen mayor loading (peso) en los componentes principales")

# Decidir cuántos componentes considerar
n_top_components = ______

# Obtener loadings absolutos de todos los componentes importantes
print(f"\n🔍 Analizando loadings de los primeros {n_top_components} componentes...")

# Para cada componente, obtener las features con mayor loading absoluto
all_loadings = pca.components_[:n_top_components, :]  # Primeros n componentes

# Crear DataFrame con loadings de todos los componentes
loadings_all = pd.DataFrame(
    all_loadings.T,
    columns=[f'PC{i+1}' for i in range(n_top_components)],
    index=X.columns
)

# ========== ESTRATEGIA: SUMAR LOADINGS ABSOLUTOS ==========
# Para cada feature, sumar su importancia (loading absoluto) en todos los componentes
print("\n📊 ESTRATEGIA: Ranking de features por suma de loadings absolutos")

# TODO: Calcular importancia total de cada feature
feature_importance_from_pca = ______

# Ordenar por importancia
feature_importance_from_pca = feature_importance_from_pca.sort_values(ascending=False)

print(f"\n🔝 TOP 20 FEATURES POR IMPORTANCIA EN PCA:")
for i, (feature, importance) in enumerate(feature_importance_from_pca.head(20).items(), 1):
    print(f"  {i:2d}. {feature:25s}: {importance:.4f}")

# ========== SELECCIONAR TOP-K FEATURES ==========
k = n_components_80  # Mismo número que usamos con PCA reducido

print(f"\n✅ Seleccionando top {k} features originales basadas en loadings de PCA")

# TODO: Seleccionar features
selected_features_pca = ______

print(f"\n📋 Features seleccionadas ({k}):")
for i, feat in enumerate(selected_features_pca, 1):
    print(f"  {i:2d}. {feat}")

# ========== PREPARAR DATASET CON FEATURES SELECCIONADAS ==========
X_pca_selected = X_scaled[:, X.columns.isin(selected_features_pca)]

print(f"\n✅ Dataset con features seleccionadas por PCA:")
print(f"   Shape: {X_pca_selected.shape}")
print(f"   Reducción: {X.shape[1]} → {X_pca_selected.shape[1]} features")

# ========== VISUALIZAR COMPARACIÓN ==========
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(16, 6))

# Top 20 features por importancia PCA
ax1.barh(range(20), feature_importance_from_pca.head(20).values, color='steelblue', alpha=0.7)
ax1.set_yticks(range(20))
ax1.set_yticklabels(feature_importance_from_pca.head(20).index, fontsize=9)
ax1.set_xlabel('Importancia Total (Suma de Loadings Absolutos)', fontsize=11)
ax1.set_title('Top 20 Features por Importancia en PCA\n(Features originales con mayor peso en componentes)', fontsize=12)
ax1.invert_yaxis()
ax1.grid(axis='x', alpha=0.3)

# Distribución de importancias
ax2.hist(feature_importance_from_pca, bins=30, alpha=0.7, color='lightgreen', edgecolor='black')
ax2.set_xlabel('Importancia', fontsize=11)
ax2.set_ylabel('Frecuencia', fontsize=11)
ax2.set_title('Distribución de Importancia de Features', fontsize=12)
ax2.axvline(feature_importance_from_pca.iloc[k-1], color='red', linestyle='--', 
            label=f'Umbral (top {k})', linewidth=2)
ax2.legend(fontsize=11)
ax2.grid(alpha=0.3)

plt.tight_layout()
plt.show()

print("\n💡 INTERPRETACIÓN:")
print("   Estas features originales son las que 'explican' los componentes principales")
print("   Ventaja: Mantienen interpretabilidad (puedes decir 'GrLivArea importa')")
print("   Diferencia con PCA: Usas features originales, no combinaciones lineales")

💡 Pistas:

Primer blank: Usa la variable que calculaste anteriormente para 80% de varianza.
Segundo blank: Calcula la suma de valores absolutos del DataFrame de loadings a lo largo del eje de componentes (.abs().sum(axis=1)).
Tercer blank: Usa el método para obtener los k elementos más grandes y extrae sus índices como lista (.nlargest(k).index.tolist()).

📖 Documentación: pandas.DataFrame.abs | pandas.Series.nlargest

🤔 Pregunta Reflexiva 5.5 (Nueva):

¿Las features seleccionadas por loadings de PCA coinciden con tu intuición?
¿Por qué este enfoque es diferente a usar directamente PC1, PC2, etc.?
¿Cuál método preferirías para un negocio: PCA completo o features basadas en loadings? ¿Por qué?

Paso 2.6: Evaluación de Performance con REGRESIÓN (RMSE y R²)¶

⏱️ Tiempo esperado: ~3-5 minutos (cross-validation puede ser lento con 80 features)

from sklearn.ensemble import RandomForestRegressor  # ⚠️ REGRESSOR, no Classifier
from sklearn.model_selection import cross_val_score
from sklearn.metrics import mean_squared_error, r2_score, mean_absolute_error
import warnings
warnings.filterwarnings('ignore')  # Suprimir warnings de convergencia

# ========== MODELO BASELINE: TODAS LAS FEATURES ORIGINALES ==========
print("\n=== EVALUACIÓN DE PERFORMANCE: PCA vs ORIGINAL ===")
print("⏱️ Esto puede tomar 1-2 minutos (cross-validation con 80 features)...\n")

# TODO: Evaluar modelo con datos originales
print("🔄 Evaluando modelo con features originales...")
rf_original = RandomForestRegressor(
    random_state=42, 
    n_estimators=______,
    max_depth=______,
    n_jobs=-1             # Usar todos los cores
)

# Usar neg_mean_squared_error y neg_mean_absolute_error para CV
scores_mse_original = -cross_val_score(rf_original, X_scaled, y, cv=5, 
                                        scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_original = cross_val_score(rf_original, X_scaled, y, cv=5, 
                                     scoring='r2', n_jobs=-1)

rmse_original = np.sqrt(scores_mse_original)

print(f"\n✅ BASELINE - Features Originales ({X.shape[1]} features):")
print(f"   RMSE: ${rmse_original.mean():,.0f} ± ${rmse_original.std():,.0f}")
print(f"   R²:   {scores_r2_original.mean():.4f} ± {scores_r2_original.std():.4f}")
print(f"   Scores RMSE: {[f'${x:,.0f}' for x in rmse_original]}")

# ========== MODELO CON PCA (80% VARIANZA) ==========
print(f"\n🔄 Evaluando modelo con PCA ({n_components_80} componentes)...")

# TODO: Aplicar PCA reducido
pca_reduced = PCA(n_components=______)
X_pca_reduced = pca_reduced.fit_transform(X_scaled)

print(f"✅ PCA transformado: {X.shape[1]} → {X_pca_reduced.shape[1]} features")

# Evaluar con PCA
rf_pca = RandomForestRegressor(
    random_state=42, 
    n_estimators=100, 
    max_depth=15,
    n_jobs=-1
)

scores_mse_pca = -cross_val_score(rf_pca, X_pca_reduced, y, cv=5, 
                                   scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_pca = cross_val_score(rf_pca, X_pca_reduced, y, cv=5, 
                                scoring='r2', n_jobs=-1)

rmse_pca = np.sqrt(scores_mse_pca)

print(f"\n✅ PCA - Componentes Reducidos ({n_components_80} componentes):")
print(f"   RMSE: ${rmse_pca.mean():,.0f} ± ${rmse_pca.std():,.0f}")
print(f"   R²:   {scores_r2_pca.mean():.4f} ± {scores_r2_pca.std():.4f}")
print(f"   Scores RMSE: {[f'${x:,.0f}' for x in rmse_pca]}")

# ========== MODELO CON FEATURES SELECCIONADAS POR PCA LOADINGS ==========
print(f"\n🔄 Evaluando modelo con features originales seleccionadas por PCA loadings...")

rf_pca_selected = RandomForestRegressor(
    random_state=42, 
    n_estimators=100, 
    max_depth=15,
    n_jobs=-1
)

scores_mse_pca_selected = -cross_val_score(rf_pca_selected, X_pca_selected, y, cv=5, 
                                             scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_pca_selected = cross_val_score(rf_pca_selected, X_pca_selected, y, cv=5, 
                                          scoring='r2', n_jobs=-1)

rmse_pca_selected = np.sqrt(scores_mse_pca_selected)

print(f"\n✅ PCA Loadings - Features Originales Seleccionadas ({len(selected_features_pca)} features):")
print(f"   RMSE: ${rmse_pca_selected.mean():,.0f} ± ${rmse_pca_selected.std():,.0f}")
print(f"   R²:   {scores_r2_pca_selected.mean():.4f} ± {scores_r2_pca_selected.std():.4f}")
print(f"   Scores RMSE: {[f'${x:,.0f}' for x in rmse_pca_selected]}")

# ========== ANÁLISIS COMPARATIVO ==========
print(f"\n" + "="*80)
print(f"{'COMPARACIÓN: ORIGINAL vs PCA vs PCA LOADINGS':^80}")
print(f"="*80)

reduction_pct = (1 - n_components_80 / X.shape[1]) * 100
rmse_diff_pca = rmse_pca.mean() - rmse_original.mean()
rmse_diff_pca_selected = rmse_pca_selected.mean() - rmse_original.mean()
r2_diff_pca = scores_r2_pca.mean() - scores_r2_original.mean()
r2_diff_pca_selected = scores_r2_pca_selected.mean() - scores_r2_original.mean()

print(f"\n📊 REDUCCIÓN DIMENSIONAL:")
print(f"   Original: {X.shape[1]} features")
print(f"   PCA: {X.shape[1]} → {n_components_80} componentes ({reduction_pct:.1f}% reducción)")
print(f"   PCA Loadings: {X.shape[1]} → {len(selected_features_pca)} features originales ({reduction_pct:.1f}% reducción)")
print(f"   Varianza retenida (PCA): {pca_reduced.explained_variance_ratio_.sum():.1%}")

print(f"\n📊 PERFORMANCE COMPARATIVO:")
print(f"\n   {'Método':<25s} {'RMSE':>15s} {'R²':>10s} {'Features':>10s}")
print(f"   {'-'*25} {'-'*15} {'-'*10} {'-'*10}")
print(f"   {'Original':<25s} ${rmse_original.mean():>14,.0f} {scores_r2_original.mean():>10.4f} {X.shape[1]:>10d}")
print(f"   {'PCA Componentes':<25s} ${rmse_pca.mean():>14,.0f} {scores_r2_pca.mean():>10.4f} {n_components_80:>10d}")
print(f"   {'PCA Loadings (Originales)':<25s} ${rmse_pca_selected.mean():>14,.0f} {scores_r2_pca_selected.mean():>10.4f} {len(selected_features_pca):>10d}")

print(f"\n📊 DIFERENCIAS VS ORIGINAL:")
print(f"   PCA Componentes:  RMSE {rmse_diff_pca:+,.0f} ({(rmse_diff_pca/rmse_original.mean())*100:+.1f}%) | R² {r2_diff_pca:+.4f}")
print(f"   PCA Loadings:     RMSE {rmse_diff_pca_selected:+,.0f} ({(rmse_diff_pca_selected/rmse_original.mean())*100:+.1f}%) | R² {r2_diff_pca_selected:+.4f}")

# Interpretación
print(f"\n💡 INTERPRETACIÓN:")
print(f"\n   🔵 PCA Componentes (PC1, PC2...):")
if rmse_pca.mean() < rmse_original.mean() * 1.05:
    print(f"      ✅ Mantiene performance similar con {reduction_pct:.0f}% reducción")
    print(f"      ⚠️ Pero: Componentes son combinaciones lineales (menos interpretables)")
else:
    print(f"      ⚠️ Pierde precisión significativa ({(rmse_diff_pca/rmse_original.mean())*100:.1f}%)")

print(f"\n   🟢 PCA Loadings (Features originales):")
if rmse_pca_selected.mean() < rmse_original.mean() * 1.05:
    print(f"      ✅ Mantiene performance similar con {reduction_pct:.0f}% reducción")
    print(f"      ✅ Plus: Usa features originales (interpretables)")
else:
    print(f"      ⚠️ Pierde precisión ({(rmse_diff_pca_selected/rmse_original.mean())*100:.1f}%)")

print(f"\n   💼 PARA NEGOCIO:")
print(f"      - PCA Componentes: Mejor para modelos 'black box' donde solo importa precisión")
print(f"      - PCA Loadings: Mejor para negocio (puedes decir 'GrLivArea es importante')")

💡 Pistas:

n_estimators: Número de árboles en el Random Forest. Valores típicos: 50-200. Más árboles = más lento pero más estable.
max_depth: Profundidad máxima de cada árbol. Valores típicos: 10-20. Mayor profundidad = más complejo.
Para PCA reducido: Usa la variable que calculaste para 80% de varianza.

📖 Documentación: RandomForestRegressor | PCA | cross_val_score

🤔 Pregunta Reflexiva 6:

¿El RMSE aumentó o disminuyó con PCA? ¿Es aceptable la diferencia?
¿R² bajó significativamente? ¿Qué te dice esto sobre la información perdida?
¿Vale la pena perder {reduction_pct:.0f}% de features si RMSE sube ${rmse_diff:,.0f}?
¿Para qué casos de negocio preferirías PCA sobre mantener todas las features?
¿Cómo afecta PCA al tiempo de entrenamiento e inferencia en producción?

Parte 3: Feature Selection - Filter Methods (25 min)¶

⏱️ Tiempo esperado: ~25 minutos

💡 Contexto de Negocio: A diferencia de PCA (que crea nuevas features abstractas), Feature Selection mantiene las features originales. Esto es CRÍTICO para: - 🗣️ Explicabilidad: Puedes decir "GrLivArea y OverallQual son las más importantes" - 📊 Monitoreo: Detectar feature drift es más fácil - 🤝 Confianza: Stakeholders entienden qué features usa el modelo

Paso 3.1: Filter Method - F-test (ANOVA F-value para Regresión)¶

from sklearn.feature_selection import SelectKBest, f_regression  # ⚠️ f_regression, no f_classif

# ========== F-TEST PARA REGRESIÓN ==========
print("\n=== FILTER METHOD: F-TEST (ANOVA) ===")
print("F-test mide la relación lineal entre cada feature y el target (SalePrice)")

# TODO: Seleccionar top-k features por F-test
k = n_components_80  # Mismo número que PCA para comparación justa

print(f"\nSeleccionando top {k} features con F-test...")

selector_f = SelectKBest(_______, k=k)
X_filter_f = selector_f.fit_transform(X_scaled, y)

# Identificar features seleccionadas
selected_features_f = X.columns[selector_f.get_support()]
print(f"\n✅ Features seleccionadas por F-test ({k}):")
for i, feat in enumerate(selected_features_f, 1):
    print(f"  {i:2d}. {feat}")

# ========== SCORES DE F-TEST ==========
scores_f = pd.Series(selector_f.scores_, index=X.columns).sort_values(ascending=False)
print(f"\n=== TOP 15 F-SCORES (Mayor correlación con SalePrice) ===")
for i, (feat, score) in enumerate(scores_f.head(15).items(), 1):
    print(f"  {i:2d}. {feat:20s}: {score:,.0f}")

# TODO: Visualizar scores (top 30 para claridad)
plt.figure(figsize=(14, 10))
scores_f.head(30).sort_values(ascending=True).plot(kind='barh', color='steelblue')
plt.xlabel('F-Score (ANOVA)', fontsize=12)
plt.title('Top 30 Features por F-test\n(Mayor F-score = Mayor relación lineal con SalePrice)', fontsize=14)
plt.grid(axis='x', alpha=0.3)
plt.tight_layout()
plt.show()

💡 Pista: Usa la función f_regression que calcula los F-scores de ANOVA para problemas de regresión.

📖 Documentación: SelectKBest | f_regression

🤔 Pregunta Reflexiva 7:

¿Las features con mayor F-score tienen sentido intuitivamente para predecir precio?
¿Esperabas que esas features fueran las más importantes? ¿Por qué sí/no?
¿F-test captura relaciones no-lineales? ¿Cuál es su limitación principal?

Paso 3.2: Filter Method - Mutual Information (Captura relaciones no-lineales)¶

from sklearn.feature_selection import mutual_info_regression  # ⚠️ mutual_info_regression, no classif

# ========== MUTUAL INFORMATION PARA REGRESIÓN ==========
print("\n=== FILTER METHOD: MUTUAL INFORMATION ===")
print("MI captura relaciones LINEALES Y NO-LINEALES (más flexible que F-test)")
print("⏱️ Esto puede tomar 30-60 segundos...")

# TODO: Seleccionar top-k features por Mutual Information
selector_mi = SelectKBest(_______, k=k)
X_filter_mi = selector_mi.fit_transform(X_scaled, y)

selected_features_mi = X.columns[selector_mi.get_support()]
print(f"\n✅ Features seleccionadas por Mutual Information ({k}):")
for i, feat in enumerate(selected_features_mi, 1):
    print(f"  {i:2d}. {feat}")

# Scores
scores_mi = pd.Series(selector_mi.scores_, index=X.columns).sort_values(ascending=False)
print(f"\n=== TOP 15 MI SCORES ===")
for i, (feat, score) in enumerate(scores_mi.head(15).items(), 1):
    print(f"  {i:2d}. {feat:20s}: {score:.4f}")

# ========== COMPARACIÓN: F-TEST vs MUTUAL INFORMATION ==========
common_features = set(selected_features_f) & set(selected_features_mi)
print(f"\n" + "="*70)
print(f"{'COMPARACIÓN: F-TEST vs MUTUAL INFORMATION':^70}")
print(f"="*70)
print(f"\n📊 Features en común: {len(common_features)} de {k} ({len(common_features)/k*100:.1f}% coincidencia)")

print(f"\n✅ Features comunes (ambos métodos las eligieron):")
for i, feat in enumerate(sorted(common_features), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n🔵 Features SOLO en F-test:")
only_f = set(selected_features_f) - set(selected_features_mi)
for i, feat in enumerate(sorted(only_f), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n🟢 Features SOLO en Mutual Information:")
only_mi = set(selected_features_mi) - set(selected_features_f)
for i, feat in enumerate(sorted(only_mi), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n💡 INTERPRETACIÓN:")
if len(common_features) / k > 0.7:
    print(f"   Alta coincidencia ({len(common_features)/k*100:.0f}%) → Ambos métodos están de acuerdo")
else:
    print(f"   Baja coincidencia ({len(common_features)/k*100:.0f}%) → MI captura relaciones no-lineales diferentes")

💡 Pista: Usa la función mutual_info_regression que estima la información mutua entre cada feature y el target para problemas de regresión.

📖 Documentación: mutual_info_regression

🤔 Pregunta Reflexiva 8:

¿Por qué MI y F-test eligieron features diferentes? ¿Qué significa esto?
¿Cuáles features "solo en MI" crees que tienen relación no-lineal con el precio?
¿Preferirías usar F-test o MI para seleccionar features? ¿Por qué?

Paso 3.3: Evaluar Performance de Filter Methods¶

# ========== EVALUACIÓN: F-TEST ==========
print("\n=== EVALUANDO FILTER METHODS ===")
print("⏱️ Esto puede tomar 1-2 minutos...\n")

print("🔄 Evaluando F-test...")
rf_filter_f = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_filter_f = -cross_val_score(rf_filter_f, X_filter_f, y, cv=5, 
                                        scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_filter_f = cross_val_score(rf_filter_f, X_filter_f, y, cv=5, 
                                     scoring='r2', n_jobs=-1)
rmse_filter_f = np.sqrt(scores_mse_filter_f)

print(f"✅ F-test ({k} features):")
print(f"   RMSE: ${rmse_filter_f.mean():,.0f} ± ${rmse_filter_f.std():,.0f}")
print(f"   R²:   {scores_r2_filter_f.mean():.4f} ± {scores_r2_filter_f.std():.4f}")

# ========== EVALUACIÓN: MUTUAL INFORMATION ==========
print(f"\n🔄 Evaluando Mutual Information...")
rf_filter_mi = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_filter_mi = -cross_val_score(rf_filter_mi, X_filter_mi, y, cv=5, 
                                         scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_filter_mi = cross_val_score(rf_filter_mi, X_filter_mi, y, cv=5, 
                                      scoring='r2', n_jobs=-1)
rmse_filter_mi = np.sqrt(scores_mse_filter_mi)

print(f"✅ Mutual Information ({k} features):")
print(f"   RMSE: ${rmse_filter_mi.mean():,.0f} ± ${rmse_filter_mi.std():,.0f}")
print(f"   R²:   {scores_r2_filter_mi.mean():.4f} ± {scores_r2_filter_mi.std():.4f}")

# ========== COMPARACIÓN ACTUALIZADA ==========
print(f"\n" + "="*70)
print(f"{'COMPARACIÓN: BASELINE vs PCA vs FILTER METHODS':^70}")
print(f"="*70)

comparison_data = {
    'Método': ['Original', 'PCA', 'F-test', 'MI'],
    'N_Features': [X.shape[1], n_components_80, k, k],
    'RMSE': [rmse_original.mean(), rmse_pca.mean(), rmse_filter_f.mean(), rmse_filter_mi.mean()],
    'R²': [scores_r2_original.mean(), scores_r2_pca.mean(), scores_r2_filter_f.mean(), scores_r2_filter_mi.mean()]
}
comparison_df = pd.DataFrame(comparison_data)
comparison_df['Reducción%'] = (1 - comparison_df['N_Features'] / X.shape[1]) * 100

print(f"\n{comparison_df.to_string(index=False)}")

print(f"\n💡 OBSERVACIONES:")
print(f"   - PCA: Reduce a componentes abstractos (no interpretables)")
print(f"   - F-test & MI: Mantienen features originales (interpretables)")
print(f"   - ¿Cuál tiene mejor RMSE? ¿Vale la pena sacrificar interpretabilidad?")

🤔 Pregunta Reflexiva 9:

¿Cuál método (PCA, F-test, MI) tuvo mejor RMSE? ¿Por qué crees que fue así?
¿La diferencia en RMSE es suficiente para justificar perder interpretabilidad (PCA)?
Si tuvieras que presentar resultados al equipo de bienes raíces, ¿qué método elegirías?

Parte 4: Feature Selection - Wrapper Methods (20 min)¶

⏱️ Tiempo esperado: ~20 minutos (reducido usando two-stage selection)

💡 Contexto: Wrapper Methods evalúan subconjuntos de features entrenando el modelo. Son más lentos pero más precisos que Filter Methods.

🚀 ESTRATEGIA OPTIMIZADA: En lugar de aplicar Forward/Backward sobre las ~79 features originales (muy lento), las aplicaremos sobre las features pre-seleccionadas por PCA Loadings. Esto es una estrategia "two-stage" común en la industria: 1. Stage 1: PCA Loadings reduce 79 → ~36-45 features (rápido) 2. Stage 2: Wrapper methods refinan esa selección (más rápido)

Paso 4.1: Forward Selection (Sequential Feature Selection)¶

from sklearn.feature_selection import SequentialFeatureSelector

# ========== TWO-STAGE SELECTION: PCA + FORWARD ==========
print("\n=== WRAPPER METHOD: FORWARD SELECTION ===")
print("💡 ESTRATEGIA: Forward Selection sobre features pre-seleccionadas por PCA")
print(f"   Stage 1 (ya hecho): PCA Loadings → {len(selected_features_pca)} features")
print(f"   Stage 2 (ahora): Forward Selection → refinar a menos features")
print("⏱️ Esto tomará ~30-60 segundos (mucho más rápido que sobre 79 features)...\n")

# Decidir cuántas features seleccionar con wrapper
k_wrapper = max(15, k // 2)  # Aproximadamente la mitad de las features PCA, o mínimo 15
print(f"🎯 Target: Seleccionar {k_wrapper} features con Forward Selection")

# TODO: Forward Selection sobre features pre-seleccionadas
estimator_forward = RandomForestRegressor(
    random_state=42, 
    n_estimators=______,
    max_depth=______,
    n_jobs=-1
)

selector_forward = SequentialFeatureSelector(
    estimator=estimator_forward, 
    n_features_to_select=k_wrapper,
    direction=______,
    cv=3,              # 3 folds para rapidez
    n_jobs=-1
)

print(f"🔄 Ejecutando Forward Selection sobre {len(selected_features_pca)} features pre-seleccionadas...")
import time
start_time = time.time()
X_forward = selector_forward.fit_transform(X_pca_selected, y)
elapsed_time = time.time() - start_time

print(f"✅ Forward Selection completado en {elapsed_time:.1f} segundos")

# Features seleccionadas (mapear índices de vuelta a nombres)
selected_indices_forward = selector_forward.get_support()
selected_features_forward = [selected_features_pca[i] for i, sel in enumerate(selected_indices_forward) if sel]

print(f"\n✅ Features seleccionadas por Forward Selection ({len(selected_features_forward)}):")
for i, feat in enumerate(selected_features_forward, 1):
    print(f"  {i:2d}. {feat}")

Paso 4.2: Backward Elimination (Sequential Feature Selection)¶

# ========== TWO-STAGE SELECTION: PCA + BACKWARD ==========
print("\n=== WRAPPER METHOD: BACKWARD ELIMINATION ===")
print("💡 ESTRATEGIA: Backward Elimination sobre features pre-seleccionadas por PCA")
print(f"   Stage 1 (ya hecho): PCA Loadings → {len(selected_features_pca)} features")
print(f"   Stage 2 (ahora): Backward Elimination → refinar a {k_wrapper} features")
print("⏱️ Esto tomará ~30-60 segundos...\n")

# TODO: Backward Selection sobre features pre-seleccionadas
estimator_backward = RandomForestRegressor(
    random_state=42, 
    n_estimators=50,
    max_depth=10,
    n_jobs=-1
)

selector_backward = SequentialFeatureSelector(
    estimator=estimator_backward, 
    n_features_to_select=k_wrapper,
    direction=______,
    cv=3,
    n_jobs=-1
)

print(f"🔄 Ejecutando Backward Elimination sobre {len(selected_features_pca)} features...")
start_time = time.time()
X_backward = selector_backward.fit_transform(X_pca_selected, y)
elapsed_time_backward = time.time() - start_time

print(f"✅ Backward Elimination completado en {elapsed_time_backward:.1f} segundos")

# Features seleccionadas (mapear índices de vuelta a nombres)
selected_indices_backward = selector_backward.get_support()
selected_features_backward = [selected_features_pca[i] for i, sel in enumerate(selected_indices_backward) if sel]

print(f"\n✅ Features seleccionadas por Backward Elimination ({len(selected_features_backward)}):")
for i, feat in enumerate(selected_features_backward, 1):
    print(f"  {i:2d}. {feat}")

# ========== COMPARAR FORWARD VS BACKWARD ==========
print(f"\n" + "="*70)
print(f"{'COMPARACIÓN: FORWARD vs BACKWARD':^70}")
print(f"="*70)

common_features_fb = set(selected_features_forward) & set(selected_features_backward)
print(f"\n📊 Features en común: {len(common_features_fb)} de {k_wrapper} ({len(common_features_fb)/k_wrapper*100:.1f}% coincidencia)")

print(f"\n✅ Features comunes (ambos métodos las eligieron):")
for i, feat in enumerate(sorted(common_features_fb), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n🔵 Features SOLO en Forward:")
only_forward = set(selected_features_forward) - set(selected_features_backward)
for i, feat in enumerate(sorted(only_forward), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n🟢 Features SOLO en Backward:")
only_backward = set(selected_features_backward) - set(selected_features_forward)
for i, feat in enumerate(sorted(only_backward), 1):
    print(f"  {i:2d}. {feat}")

print(f"\n⏱️ TIEMPO DE EJECUCIÓN:")
print(f"   Forward:  {elapsed_time:.1f}s")
print(f"   Backward: {elapsed_time_backward:.1f}s")

print(f"\n💡 OBSERVACIÓN:")
if len(common_features_fb) / k_wrapper > 0.7:
    print(f"   Alta coincidencia ({len(common_features_fb)/k_wrapper*100:.0f}%) → Ambos métodos convergen")
else:
    print(f"   Baja coincidencia ({len(common_features_fb)/k_wrapper*100:.0f}%) → Orden de selección importa")

💡 Pistas:

Para Forward: n_estimators=50 y max_depth=10 (valores más pequeños para rapidez), direction='forward'.
Para Backward: Igual que Forward pero direction='backward'.
Two-Stage Strategy: Aplicamos wrapper methods sobre las features pre-seleccionadas por PCA (~36-45 features) en lugar de las 79 originales. Esto reduce el tiempo de 2-3 minutos a 30-60 segundos.
k_wrapper se calcula automáticamente como la mitad de k (o mínimo 15) para refinar aún más la selección.

📖 Documentación: SequentialFeatureSelector

🤔 Pregunta Reflexiva (Forward vs Backward):

¿Forward o Backward fue más rápido? ¿Por qué?
¿Por qué los dos métodos eligieron features diferentes si ambos usan el mismo modelo?
¿Cuándo preferirías usar uno sobre el otro?

Paso 4.3: Recursive Feature Elimination (RFE)¶

from sklearn.feature_selection import RFE

# ========== TWO-STAGE SELECTION: PCA + RFE ==========
print("\n=== WRAPPER METHOD: RFE (Recursive Feature Elimination) ===")
print("💡 ESTRATEGIA: RFE sobre features pre-seleccionadas por PCA")
print(f"   Stage 1 (ya hecho): PCA Loadings → {len(selected_features_pca)} features")
print(f"   Stage 2 (ahora): RFE → refinar a {k_wrapper} features")
print("⏱️ Esto tomará ~45-90 segundos...\n")

# TODO: RFE con Random Forest Regressor sobre features pre-seleccionadas
estimator = RandomForestRegressor(
    random_state=42, 
    n_estimators=______,
    max_depth=______,
    n_jobs=-1
)
selector_rfe = RFE(estimator=estimator, n_features_to_select=k_wrapper, step=______)

print(f"🔄 Ejecutando RFE sobre {len(selected_features_pca)} features...")
import time
start_time = time.time()
X_rfe = selector_rfe.fit_transform(X_pca_selected, y)
elapsed_time = time.time() - start_time

print(f"✅ RFE completado en {elapsed_time:.1f} segundos")

# Features seleccionadas (mapear índices de vuelta a nombres)
selected_indices_rfe = selector_rfe.get_support()
selected_features_rfe = [selected_features_pca[i] for i, sel in enumerate(selected_indices_rfe) if sel]

print(f"\n✅ Features seleccionadas por RFE ({len(selected_features_rfe)}):")
for i, feat in enumerate(selected_features_rfe, 1):
    print(f"  {i:2d}. {feat}")

# Ranking de features (solo sobre las pre-seleccionadas por PCA)
ranking = pd.Series(selector_rfe.ranking_, index=selected_features_pca).sort_values()
print(f"\nRanking de features (1 = seleccionada, solo mostrando top 20):")
print(ranking.head(20))

# Visualizar ranking (top 30 para claridad)
plt.figure(figsize=(12, 8))
ranking.head(30).sort_values(ascending=False).plot(kind='barh')
plt.xlabel('Ranking (1 = mejor, números mayores = eliminadas antes)')
plt.title(f'RFE Feature Ranking - Top 30 de {len(selected_features_pca)} features pre-seleccionadas')
plt.tight_layout()
plt.show()

💡 Pistas:

n_estimators y max_depth: Igual que Forward/Backward (50 y 10 respectivamente para rapidez).
step: Número de features a eliminar en cada iteración. step=2 elimina 2 features por iteración (más rápido que step=1).

📖 Documentación: RFE

Paso 4.4: Comparación de Features Seleccionadas (Todos los Wrapper Methods)¶

# Comparar features seleccionadas por diferentes métodos
print("\n=== COMPARACIÓN DE FEATURES SELECCIONADAS ===")

# Crear conjunto de features por método
features_dict = {
    'F-test': set(selected_features_f),
    'Mutual Info': set(selected_features_mi),
    'Forward': set(selected_features_forward),
    'Backward': set(selected_features_backward),
    'RFE': set(selected_features_rfe)
}

# Features en al menos 2 métodos
all_features = set()
for features in features_dict.values():
    all_features.update(features)

feature_counts = {}
for feature in all_features:
    count = sum(1 for features in features_dict.values() if feature in features)
    feature_counts[feature] = count

# Features consistentes (en todos los métodos)
consistent_features = [f for f, count in feature_counts.items() if count == 3]
print(f"\nFeatures consistentes (en todos los métodos): {len(consistent_features)}")
print(consistent_features)

# Features en al menos 2 métodos
robust_features = [f for f, count in feature_counts.items() if count >= 2]
print(f"\nFeatures robustas (≥2 métodos): {len(robust_features)}")
print(robust_features)

print("\n💡 OBSERVACIÓN:")
print(f"   Forward, Backward y RFE son todos wrapper methods, pero usan estrategias diferentes")

🤔 Pregunta Reflexiva 10:

¿Cuánto tiempo tomaron los wrapper methods (Forward, Backward, RFE) vs Filter Methods?
¿Los 3 wrapper methods seleccionaron features similares? ¿Qué significa esto?
¿Cuándo usarías Forward vs Backward vs RFE en un proyecto real?

Paso 4.5: Evaluar Performance de Wrapper Methods¶

# ========== EVALUACIÓN: FORWARD SELECTION ==========
print("\n=== EVALUANDO WRAPPER METHODS ===")
print("⏱️ Cross-validation con features de cada método...\n")

print("🔄 Evaluando Forward Selection...")
rf_forward = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_forward = -cross_val_score(rf_forward, X_forward, y, cv=5, 
                                       scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_forward = cross_val_score(rf_forward, X_forward, y, cv=5, 
                                    scoring='r2', n_jobs=-1)
rmse_forward = np.sqrt(scores_mse_forward)

print(f"✅ Forward Selection ({len(selected_features_forward)} features):")
print(f"   RMSE: ${rmse_forward.mean():,.0f} ± ${rmse_forward.std():,.0f}")
print(f"   R²:   {scores_r2_forward.mean():.4f} ± {scores_r2_forward.std():.4f}")

# ========== EVALUACIÓN: BACKWARD ELIMINATION ==========
print(f"\n🔄 Evaluando Backward Elimination...")
rf_backward = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_backward = -cross_val_score(rf_backward, X_backward, y, cv=5, 
                                         scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_backward = cross_val_score(rf_backward, X_backward, y, cv=5, 
                                      scoring='r2', n_jobs=-1)
rmse_backward = np.sqrt(scores_mse_backward)

print(f"✅ Backward Elimination ({len(selected_features_backward)} features):")
print(f"   RMSE: ${rmse_backward.mean():,.0f} ± ${rmse_backward.std():,.0f}")
print(f"   R²:   {scores_r2_backward.mean():.4f} ± {scores_r2_backward.std():.4f}")

# ========== EVALUACIÓN: RFE ==========
print(f"\n🔄 Evaluando RFE...")
rf_rfe = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_rfe = -cross_val_score(rf_rfe, X_rfe, y, cv=5, 
                                   scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_rfe = cross_val_score(rf_rfe, X_rfe, y, cv=5, 
                                scoring='r2', n_jobs=-1)
rmse_rfe = np.sqrt(scores_mse_rfe)

print(f"✅ RFE ({len(selected_features_rfe)} features):")
print(f"   RMSE: ${rmse_rfe.mean():,.0f} ± ${rmse_rfe.std():,.0f}")
print(f"   R²:   {scores_r2_rfe.mean():.4f} ± {scores_r2_rfe.std():.4f}")

# ========== COMPARACIÓN ACTUALIZADA ==========
print(f"\n" + "="*80)
print(f"{'COMPARACIÓN: TODOS LOS MÉTODOS HASTA AHORA':^80}")
print(f"="*80)

comparison_updated = {
    'Método': ['Original', 'PCA Componentes', 'PCA Loadings', 'F-test', 'MI', 'Forward', 'Backward', 'RFE'],
    'N_Features': [X.shape[1], n_components_80, k, k, k, k, k, k],
    'RMSE': [rmse_original.mean(), rmse_pca.mean(), rmse_pca_selected.mean(), rmse_filter_f.mean(), 
             rmse_filter_mi.mean(), rmse_forward.mean(), rmse_backward.mean(), rmse_rfe.mean()],
    'R²': [scores_r2_original.mean(), scores_r2_pca.mean(), scores_r2_pca_selected.mean(), scores_r2_filter_f.mean(), 
           scores_r2_filter_mi.mean(), scores_r2_forward.mean(), scores_r2_backward.mean(), scores_r2_rfe.mean()]
}
comparison_updated_df = pd.DataFrame(comparison_updated)
comparison_updated_df['Reducción%'] = (1 - comparison_updated_df['N_Features'] / X.shape[1]) * 100
comparison_updated_df = comparison_updated_df.sort_values('RMSE')

print(f"\n{comparison_updated_df.to_string(index=False)}")

print(f"\n💡 OBSERVACIÓN:")
best_method = comparison_updated_df.iloc[0]['Método']
best_rmse = comparison_updated_df.iloc[0]['RMSE']
print(f"   🏆 Mejor RMSE: {best_method} (${best_rmse:,.0f})")

Parte 5: Feature Selection - Embedded Methods (15 min)¶

Paso 5.1: Random Forest Feature Importance¶

# TODO: Entrenar Random Forest y obtener importances
print("\n=== EMBEDDED METHODS: Random Forest ===")
from sklearn.ensemble import RandomForestClassifier

rf_embedded = RandomForestClassifier(random_state=42, n_estimators=200, max_depth=10)
rf_embedded.fit(X_scaled, y)

# Feature importances
importances = pd.Series(rf_embedded.feature_importances_, index=X.columns).sort_values(ascending=False)
print("Top 10 features por importancia:")
print(importances.head(10))

# Visualizar importances
plt.figure(figsize=(12, 8))
importances.sort_values(ascending=True).plot(kind='barh')
plt.xlabel('Feature Importance')
plt.title('Random Forest Feature Importances')
plt.tight_layout()
plt.show()

# TODO: Seleccionar top-k features
top_k_features = importances.nlargest(k).index
X_rf_importance = X_scaled[:, X.columns.isin(top_k_features)]

print(f"\nFeatures seleccionadas por RF Importance ({k}):")
print(list(top_k_features))

# Evaluar
rmse_rf_importance = np.sqrt(-cross_val_score(rf_embedded, X_rf_importance, y, cv=5, scoring='neg_mean_squared_error'))
print(f"\nRMSE RF Importance ({k} features): ${rmse_rf_importance.mean():,.0f} ± ${rmse_rf_importance.std():,.0f}")
scores_rf_importance = cross_val_score(rf_embedded, X_rf_importance, y, cv=5, scoring='accuracy')
print(f"\nRF Importance ({k} features):")
print(f"  Mean: {scores_rf_importance.mean():.3f}")
print(f"  Std: {scores_rf_importance.std():.3f}")

🤔 Pregunta Reflexiva 11:

¿Qué features crees que son "realmente importantes" basándote en los 3 métodos?
¿Por qué RF Importance podría diferir de F-test/MI/RFE?
¿Confiarías más en un método o en el consenso de varios?

Paso 5.2: Lasso (L1 Regularization para Regresión)¶

from sklearn.linear_model import LassoCV

# ========== LASSO PARA FEATURE SELECTION ==========
print("\n=== EMBEDDED METHOD: Lasso (L1 Regularization) ===")
print("Lasso penaliza coeficientes, forzando a 0 features no importantes")
print("⏱️ Esto puede tomar 30-60 segundos...\n")

# TODO: Lasso para regresión
lasso = LassoCV(cv=5, random_state=42, max_iter=______)
lasso.fit(X_scaled, y)

print(f"✅ Lasso alpha seleccionado: {lasso.alpha_:.4f}")

# Features seleccionadas (coef != 0)
lasso_nonzero = X.columns[lasso.coef_ != 0]
print(f"\n📊 Features con coeficiente no-cero: {len(lasso_nonzero)} de {X.shape[1]}")

# Seleccionar top-k por magnitud de coeficiente
coef_abs = pd.Series(np.abs(lasso.coef_), index=X.columns).sort_values(ascending=False)
lasso_features = coef_abs.nlargest(k).index

print(f"\n✅ Top {k} features por magnitud de coeficiente Lasso:")
for i, (feat, coef) in enumerate(coef_abs.nlargest(k).items(), 1):
    print(f"  {i:2d}. {feat:20s}: |{coef:.6f}|")

# TODO: Visualizar coeficientes (top 30)
plt.figure(figsize=(14, 10))
coef_abs.head(30).sort_values(ascending=True).plot(kind='barh', color='purple')
plt.xlabel('|Coeficiente Lasso|', fontsize=12)
plt.title('Top 30 Features por Magnitud de Coeficiente Lasso\n(Mayor magnitud = Mayor importancia)', fontsize=14)
plt.grid(axis='x', alpha=0.3)
plt.tight_layout()
plt.show()

# Preparar features para evaluación
X_lasso = X_scaled[:, X.columns.isin(lasso_features)]

# ========== EVALUAR CON RANDOM FOREST ==========
print(f"\n🔄 Evaluando Lasso selection...")
rf_lasso = RandomForestRegressor(random_state=42, n_estimators=100, max_depth=15, n_jobs=-1)

scores_mse_lasso = -cross_val_score(rf_lasso, X_lasso, y, cv=5, 
                                     scoring='neg_mean_squared_error', n_jobs=-1)
scores_r2_lasso = cross_val_score(rf_lasso, X_lasso, y, cv=5, 
                                  scoring='r2', n_jobs=-1)
rmse_lasso = np.sqrt(scores_mse_lasso)

print(f"✅ Lasso selection ({k} features):")
print(f"   RMSE: ${rmse_lasso.mean():,.0f} ± ${rmse_lasso.std():,.0f}")
print(f"   R²:   {scores_r2_lasso.mean():.4f} ± {scores_r2_lasso.std():.4f}")

💡 Pista: max_iter es el número máximo de iteraciones para convergencia. Usa 10000 para asegurar que el algoritmo converja completamente (valores muy bajos pueden causar warnings de convergencia).

📖 Documentación: LassoCV

🤔 Pregunta Reflexiva (Lasso):

¿Lasso forzó muchas features a 0? ¿Qué te dice esto sobre la redundancia de features?
¿Las features que Lasso considera importantes coinciden con los otros métodos?
¿Lasso es más útil para interpretabilidad o para performance?

🤔 Pregunta Reflexiva 12 (Final de Parte 6):

¿Qué método dio mejor RMSE? ¿Es el que recomendarías para producción? ¿Por qué sí o no?
¿Por qué Feature Selection es preferible a PCA en el contexto de bienes raíces?
Si tuvieras que presentar estos resultados al CEO de la empresa inmobiliaria (no técnico), ¿qué 3 puntos clave destacarías?
¿Cómo comunicarías el trade-off entre reducir features (velocidad) y mantener precisión (RMSE)?
¿El método ganador es interpretable? ¿Por qué esto es crítico para esta industria?

Parte 7: Investigación Libre (opcional, +15 min)¶

Explora uno de estos temas avanzados:

Opción A: Forward/Backward Selection¶

Implementa Sequential Feature Selection y compara con RFE:

from sklearn.feature_selection import SequentialFeatureSelector

# Forward selection
sfs_forward = SequentialFeatureSelector(
    RandomForestClassifier(random_state=42),
    n_features_to_select=k,
    direction='forward',
    cv=3
)
# Implementa y evalúa

Opción B: SHAP Values para Feature Importance¶

Usa SHAP para una interpretación más robusta:

import shap

# Train model
model = RandomForestClassifier(random_state=42)
model.fit(X_scaled, y)

# SHAP values
explainer = shap.TreeExplainer(model)
shap_values = explainer.shap_values(X_scaled)

# Visualize
shap.summary_plot(shap_values, X_scaled, feature_names=X.columns)

Opción C: PCA Incremental¶

Para datasets que no caben en memoria:

from sklearn.decomposition import IncrementalPCA

# PCA incremental (mini-batches)
ipca = IncrementalPCA(n_components=k, batch_size=1000)
# Implementa con batches

Opción D: Feature Engineering + Selection¶

Crea interaction features y luego aplica selection:

from sklearn.preprocessing import PolynomialFeatures

# Crear interacciones
poly = PolynomialFeatures(degree=2, include_bias=False, interaction_only=True)
X_interactions = poly.fit_transform(X_scaled)

# Aplicar feature selection en features expandidas
# Evaluar si mejora

📝 Reflexión Final Integradora (OBLIGATORIO)¶

⏱️ Tiempo estimado: 15-20 minutos

Responde las siguientes preguntas en una celda Markdown al final de tu notebook. Estas preguntas te ayudarán a consolidar tu comprensión y demostrar pensamiento crítico.

A. Sobre PCA¶

Interpretabilidad: ¿Puedes explicar en términos simples qué representa PC1 en el contexto de precios de casas? ¿Es esto útil para un agente inmobiliario?
Varianza Explicada: Si PC1 captura 40% de varianza, ¿qué significa esto exactamente? ¿Qué información se "pierde" en el 60% restante?
Cuándo usar PCA: Menciona 3 escenarios reales donde PCA sería MÁS útil que Feature Selection (no en este dataset, sino en general).
Limitaciones: ¿Cuál es la mayor desventaja de PCA para este problema de bienes raíces?

B. Sobre Feature Selection¶

Consistencia: Si F-test, MI, RFE y Lasso eligieron features diferentes, ¿cómo decides cuál conjunto usar? ¿Qué estrategia aplicarías?
Features Redundantes: Si GarageArea y GarageCars están altamente correlacionadas, ¿cuál deberías eliminar? ¿Cómo decides?
Métodos Filter vs Wrapper: ¿Por qué RFE es más lento que F-test? ¿Cuándo justifica el tiempo extra?
Lasso Shrinkage: Si Lasso forzó 40 features a coeficiente 0, ¿qué te dice esto sobre la redundancia en el dataset?

🏠 Trabajos Domiciliarios (Opcionales pero Recomendados)¶

⏱️ Tiempo estimado: 2-4 horas por trabajo

Estos trabajos te permitirán profundizar en aspectos específicos de PCA y Feature Selection. Puedes elegir uno, varios, o todos.

Trabajo Domiciliario 1: Comparación de Umbrales de Varianza en PCA¶

Objetivo: Entender cómo el umbral de varianza afecta performance y complejidad.

Tarea: 1. Implementa PCA con diferentes umbrales de varianza explicada: 70%, 80%, 90%, 95%, 99% 2. Para cada umbral: - Registra número de componentes - Calcula RMSE y R² con cross-validation - Mide tiempo de entrenamiento e inferencia 3. Grafica: - Gráfico 1: Varianza explicada vs Número de componentes - Gráfico 2: RMSE vs Número de componentes - Gráfico 3: Tiempo de entrenamiento vs Número de componentes 4. Responde: - ¿Dónde está el "punto óptimo" (elbow) en el trade-off varianza-performance? - ¿A partir de cuántos componentes el RMSE deja de mejorar significativamente? - ¿Cuánto tiempo ahorras en inferencia con 70% vs 99% varianza? - ¿Qué umbral recomendarías para una app móvil que necesita respuestas instantáneas?

Entregable: Notebook con código, gráficos y análisis (max 3 páginas).

Trabajo Domiciliario 2: Feature Selection - Comparación de Métodos Filter¶

Objetivo: Comparar diferentes métricas de importancia en Filter Methods.

Tarea: 1. Implementa los siguientes Filter Methods y selecciona top-k=15 features: - F-test (ANOVA F-value): f_regression - Mutual Information: mutual_info_regression - Correlación de Pearson: df.corr()['SalePrice'].abs() - Chi-squared: chi2 (requiere discretización de features) 2. Para cada método: - Lista las 15 features seleccionadas - Calcula RMSE y R² con Random Forest 3. Crea un diagrama de Venn mostrando overlap entre los 4 métodos 4. Identifica: - Features "robustas" (seleccionadas por los 4 métodos) - Features "controversiales" (seleccionadas solo por 1-2 métodos) 5. Responde: - ¿Por qué Correlación de Pearson eligió features diferentes a F-test si ambos miden relaciones lineales? - ¿Cuál método fue mejor en RMSE? ¿Coincide con tu intuición? - ¿Cómo decidirías qué features usar si los métodos no están de acuerdo? - ¿Preferirías usar un "ensemble" de métodos? ¿Cómo lo implementarías?

Entregable: Notebook con código, Venn diagram, tabla comparativa y análisis (max 3 páginas).

Trabajo Domiciliario 3: Wrapper Methods - Sequential Feature Selection¶

Objetivo: Comparar RFE con Sequential Feature Selection (SFS) y entender trade-offs de tiempo vs precisión.

Tarea: 1. Implementa los siguientes Wrapper Methods y selecciona top-k=15 features: - RFE (Recursive Feature Elimination): elimina features de peor a mejor - Forward Selection: agrega features de a una (empieza con 0) - Backward Selection: elimina features de a una (empieza con 80) - Bi-Directional Selection: combina forward y backward

Usa SequentialFeatureSelector de sklearn:

from sklearn.feature_selection import SequentialFeatureSelector

# Forward Selection
sfs_forward = SequentialFeatureSelector(
    RandomForestRegressor(random_state=42, n_estimators=50),
    n_features_to_select=15,
    direction='forward',
    cv=3,
    n_jobs=-1
)

Para cada método:
Registra features seleccionadas
Mide tiempo de ejecución (crítico para comparar)
Calcula RMSE y R² con Random Forest (100 estimators, cv=5)
Compara:
¿Qué método es más rápido? ¿Por qué?
¿Qué método tiene mejor RMSE?
¿El overlap de features es alto o bajo?
Experimenta con Early Stopping:
Implementa Forward Selection con early stopping (para si RMSE no mejora en 3 iteraciones)
Compara tiempo y features seleccionadas vs Forward Selection completo
Responde:
¿Vale la pena el tiempo extra de RFE/SFS vs Filter Methods (F-test)?
¿En qué escenarios usarías Wrapper Methods en lugar de Filter Methods?
¿Cómo podrías acelerar Wrapper Methods sin sacrificar mucha precisión?

Entregable: Notebook con código, tabla de tiempos, comparación de features y análisis (max 4 páginas).

📚 Recursos Adicionales¶

🎉 ¡Assignment completado! Sube tu notebook al portfolio y prepárate para discutir tus hallazgos en clase.