Практика — поліноміальна регресія та регуляризація на Housing Prices

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

# Завантаження датасету California Housing (вбудований у sklearn)
from sklearn.datasets import fetch_california_housing

data = fetch_california_housing(as_frame=True)
df = data.frame
df['MedHouseVal'] = data.target

print(f"Розмір датасету: {df.shape}")
print(f"\nПерші 5 рядків:")
df.head()

# Перевірка пропусків
print("Пропуски у датасеті:")
print(df.isnull().sum())

print("\n✅ Пропусків немає — California Housing чистий датасет")

# Кореляція з цільовою змінною
correlations = df.corr()['MedHouseVal'].sort_values(ascending=False)
print("Кореляція ознак з MedHouseVal (ціна):\n")
print(correlations)

# Візуалізація
plt.figure(figsize=(10, 6))
correlations[1:].plot(kind='barh', color='steelblue', edgecolor='k')
plt.xlabel('Кореляція з ціною', fontsize=12)
plt.ylabel('Ознаки', fontsize=12)
plt.title('Кореляція ознак з медіанною ціною будинку', fontsize=14, fontweight='bold')
plt.axvline(0, color='black', linewidth=0.8, linestyle='--')
plt.grid(alpha=0.3, axis='x')
plt.tight_layout()
plt.show()

from sklearn.preprocessing import PolynomialFeatures

# Вибираємо підмножину ознак
selected_features = ['MedInc', 'AveRooms', 'HouseAge', 'Latitude', 'Longitude']
X_selected = df[selected_features]
y = df['MedHouseVal']

print(f"Початкова кількість ознак: {X_selected.shape[1]}")

# Застосовуємо PolynomialFeatures
poly = PolynomialFeatures(degree=2, include_bias=False)
X_poly = poly.fit_transform(X_selected)

print(f"Після PolynomialFeatures(degree=2): {X_poly.shape[1]} ознак")
print(f"\nНові ознаки включають:")
print(f"  - Оригінальні: {X_selected.shape[1]}")
print(f"  - Квадрати (x²): {X_selected.shape[1]}")
print(f"  - Парні взаємодії (x_i * x_j): {X_selected.shape[1] * (X_selected.shape[1] - 1) // 2}")

# Назви нових ознак
feature_names = poly.get_feature_names_out(selected_features)
print(f"\nПриклади нових ознак:")
for name in feature_names[5:15]:  # показуємо 10 прикладів
    print(f"  - {name}")

from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LinearRegression, Ridge, Lasso
from sklearn.metrics import r2_score, mean_absolute_error, mean_squared_error

# Розділення на train/test
X = df[selected_features]
y = df['MedHouseVal']

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

print(f"Train set: {X_train.shape[0]} зразків")
print(f"Test set:  {X_test.shape[0]} зразків")

# Базовий Pipeline (без поліноміальних ознак — для порівняння)
baseline_pipeline = Pipeline([
    ('scaler', StandardScaler()),
    ('regressor', LinearRegression())
])

baseline_pipeline.fit(X_train, y_train)

r2_train_baseline = baseline_pipeline.score(X_train, y_train)
r2_test_baseline = baseline_pipeline.score(X_test, y_test)

print(f"\n📊 Baseline (без поліноміальних ознак):")
print(f"  R² train: {r2_train_baseline:.4f}")
print(f"  R² test:  {r2_test_baseline:.4f}")

# Pipeline з поліноміальними ознаками
poly_pipeline = Pipeline([
    ('poly', PolynomialFeatures(degree=2, include_bias=False)),
    ('scaler', StandardScaler()),
    ('regressor', LinearRegression())
])

poly_pipeline.fit(X_train, y_train)

# Оцінка
y_train_pred = poly_pipeline.predict(X_train)
y_test_pred = poly_pipeline.predict(X_test)

r2_train = r2_score(y_train, y_train_pred)
r2_test = r2_score(y_test, y_test_pred)
mae_test = mean_absolute_error(y_test, y_test_pred)
mse_test = mean_squared_error(y_test, y_test_pred)

print("📊 LinearRegression + PolynomialFeatures(degree=2):")
print(f"  R² train: {r2_train:.4f}")
print(f"  R² test:  {r2_test:.4f}")
print(f"  MAE test: {mae_test:.4f}")
print(f"  RMSE test: {np.sqrt(mse_test):.4f}")

gap = r2_train - r2_test
print(f"\n⚠️ Gap (R² train - R² test): {gap:.4f}")

if gap > 0.05:
    print("❌ Помітне перенавчання! Потрібна регуляризація.")
else:
    print("✅ Перенавчання незначне.")

from sklearn.model_selection import GridSearchCV

# Підбір оптимального alpha через GridSearchCV
param_grid = {
    'regressor__alpha': [0.1, 0.5, 1, 5, 10, 50, 100]
}

ridge_pipeline = Pipeline([
    ('poly', PolynomialFeatures(degree=2, include_bias=False)),
    ('scaler', StandardScaler()),
    ('regressor', Ridge())
])

grid_search = GridSearchCV(
    ridge_pipeline,
    param_grid,
    cv=5,
    scoring='r2',
    n_jobs=-1
)

grid_search.fit(X_train, y_train)

print("🔍 Результати GridSearchCV для Ridge:")
print(f"  Оптимальний alpha: {grid_search.best_params_['regressor__alpha']}")
print(f"  R² (крос-валідація): {grid_search.best_score_:.4f}")

# Оцінка на test
best_ridge = grid_search.best_estimator_
r2_train_ridge = best_ridge.score(X_train, y_train)
r2_test_ridge = best_ridge.score(X_test, y_test)

y_test_pred_ridge = best_ridge.predict(X_test)
mae_ridge = mean_absolute_error(y_test, y_test_pred_ridge)
rmse_ridge = np.sqrt(mean_squared_error(y_test, y_test_pred_ridge))

print(f"\n📊 Ridge (alpha={grid_search.best_params_['regressor__alpha']}):")
print(f"  R² train: {r2_train_ridge:.4f}")
print(f"  R² test:  {r2_test_ridge:.4f}")
print(f"  MAE test: {mae_ridge:.4f}")
print(f"  RMSE test: {rmse_ridge:.4f}")
print(f"  Gap: {r2_train_ridge - r2_test_ridge:.4f}")

# Візуалізація залежності R² від alpha
cv_results = grid_search.cv_results_
alphas_tested = [0.1, 0.5, 1, 5, 10, 50, 100]
mean_scores = cv_results['mean_test_score']

plt.figure(figsize=(10, 6))
plt.plot(alphas_tested, mean_scores, marker='o', linewidth=2, markersize=8, color='steelblue')
plt.xscale('log')
plt.xlabel('Alpha (логарифмічна шкала)', fontsize=12)
plt.ylabel('R² (крос-валідація, 5-fold)', fontsize=12)
plt.title('Підбір оптимального alpha для Ridge', fontsize=14, fontweight='bold')
plt.axvline(x=10, color='red', linestyle='--', alpha=0.7, label='Оптимальний alpha=10')
plt.legend(fontsize=11)
plt.grid(alpha=0.3)
plt.tight_layout()
plt.show()

# Підбір alpha для Lasso
param_grid_lasso = {
    'regressor__alpha': [0.001, 0.01, 0.05, 0.1, 0.5, 1]
}

lasso_pipeline = Pipeline([
    ('poly', PolynomialFeatures(degree=2, include_bias=False)),
    ('scaler', StandardScaler()),
    ('regressor', Lasso(max_iter=10000))
])

grid_search_lasso = GridSearchCV(
    lasso_pipeline,
    param_grid_lasso,
    cv=5,
    scoring='r2',
    n_jobs=-1
)

grid_search_lasso.fit(X_train, y_train)

print("🔍 Результати GridSearchCV для Lasso:")
print(f"  Оптимальний alpha: {grid_search_lasso.best_params_['regressor__alpha']}")
print(f"  R² (крос-валідація): {grid_search_lasso.best_score_:.4f}")

# Оцінка на test
best_lasso = grid_search_lasso.best_estimator_
r2_train_lasso = best_lasso.score(X_train, y_train)
r2_test_lasso = best_lasso.score(X_test, y_test)

y_test_pred_lasso = best_lasso.predict(X_test)
mae_lasso = mean_absolute_error(y_test, y_test_pred_lasso)
rmse_lasso = np.sqrt(mean_squared_error(y_test, y_test_pred_lasso))

print(f"\n📊 Lasso (alpha={grid_search_lasso.best_params_['regressor__alpha']}):")
print(f"  R² train: {r2_train_lasso:.4f}")
print(f"  R² test:  {r2_test_lasso:.4f}")
print(f"  MAE test: {mae_lasso:.4f}")
print(f"  RMSE test: {rmse_lasso:.4f}")
print(f"  Gap: {r2_train_lasso - r2_test_lasso:.4f}")

# Скільки ознак обнулено?
lasso_regressor = best_lasso.named_steps['regressor']
n_nonzero = np.sum(np.abs(lasso_regressor.coef_) > 1e-5)
n_total = len(lasso_regressor.coef_)

print(f"\n🔍 Feature Selection:")
print(f"  Всього поліноміальних ознак: {n_total}")
print(f"  Відібрано Lasso: {n_nonzero}")
print(f"  Обнулено: {n_total - n_nonzero}")

# Збір результатів
results_df = pd.DataFrame({
    'Модель': [
        'Baseline (без полінома)',
        'LinearRegression + Poly',
        'Ridge + Poly (α=10)',
        'Lasso + Poly (α=0.01)'
    ],
    'R² train': [
        r2_train_baseline,
        r2_train,
        r2_train_ridge,
        r2_train_lasso
    ],
    'R² test': [
        r2_test_baseline,
        r2_test,
        r2_test_ridge,
        r2_test_lasso
    ],
    'MAE test': [
        mean_absolute_error(y_test, baseline_pipeline.predict(X_test)),
        mae_test,
        mae_ridge,
        mae_lasso
    ],
    'RMSE test': [
        np.sqrt(mean_squared_error(y_test, baseline_pipeline.predict(X_test))),
        np.sqrt(mse_test),
        rmse_ridge,
        rmse_lasso
    ],
    'Gap': [
        r2_train_baseline - r2_test_baseline,
        gap,
        r2_train_ridge - r2_test_ridge,
        r2_train_lasso - r2_test_lasso
    ]
})

print("📊 Порівняльна таблиця моделей:\n")
print(results_df.to_string(index=False))

# Висновок
best_model_idx = results_df['R² test'].idxmax()
best_model_name = results_df.loc[best_model_idx, 'Модель']
best_r2 = results_df.loc[best_model_idx, 'R² test']

print(f"\n🏆 Найкраща модель: {best_model_name}")
print(f"   R² test: {best_r2:.4f}")
print(f"\n✅ Покращення відносно Baseline: +{(best_r2 - r2_test_baseline):.4f} ({(best_r2 - r2_test_baseline) / r2_test_baseline * 100:.1f}%)")

# График 1: Predicted vs Actual
y_test_pred_best = best_ridge.predict(X_test)

fig, axes = plt.subplots(1, 3, figsize=(16, 5))

# 1. Predicted vs Actual
axes[0].scatter(y_test, y_test_pred_best, alpha=0.4, s=20, edgecolors='k', linewidths=0.3)
axes[0].plot([y_test.min(), y_test.max()], [y_test.min(), y_test.max()], 
             'r--', linewidth=2, label='Ideal')
axes[0].set_xlabel('Справжня ціна', fontsize=11)
axes[0].set_ylabel('Передбачена ціна', fontsize=11)
axes[0].set_title('Predicted vs Actual (Ridge)', fontsize=12, fontweight='bold')
axes[0].legend(fontsize=10)
axes[0].grid(alpha=0.3)

# 2. Residuals
residuals = y_test - y_test_pred_best
axes[1].scatter(y_test_pred_best, residuals, alpha=0.4, s=20, edgecolors='k', linewidths=0.3)
axes[1].axhline(0, color='red', linestyle='--', linewidth=2)
axes[1].set_xlabel('Передбачена ціна', fontsize=11)
axes[1].set_ylabel('Залишки (residuals)', fontsize=11)
axes[1].set_title('Residual Plot', fontsize=12, fontweight='bold')
axes[1].grid(alpha=0.3)

# 3. Feature Importance (топ-10 коефіцієнтів)
ridge_regressor = best_ridge.named_steps['regressor']
poly_names = best_ridge.named_steps['poly'].get_feature_names_out(selected_features)
coefs = ridge_regressor.coef_

# Топ-10 за абсолютним значенням
coef_importance = pd.DataFrame({
    'Feature': poly_names,
    'Coefficient': coefs
}).sort_values('Coefficient', key=abs, ascending=False).head(10)

axes[2].barh(range(len(coef_importance)), coef_importance['Coefficient'], 
             color='steelblue', alpha=0.7, edgecolor='k')
axes[2].set_yticks(range(len(coef_importance)))
axes[2].set_yticklabels(coef_importance['Feature'], fontsize=9)
axes[2].set_xlabel('Коефіцієнт', fontsize=11)
axes[2].set_title('Топ-10 найважливіших ознак', fontsize=12, fontweight='bold')
axes[2].grid(alpha=0.3, axis='x')
axes[2].invert_yaxis()

plt.tight_layout()
plt.show()

# Всі коефіцієнти
all_coefs = pd.DataFrame({
    'Feature': poly_names,
    'Coefficient': ridge_regressor.coef_
}).sort_values('Coefficient', key=abs, ascending=False)

print("📊 Топ-15 найважливіших поліноміальних ознак:\n")
print(all_coefs.head(15).to_string(index=False))

print("\n💡 Інтерпретація:")
print("1. MedInc (дохід) — найвпливовіша ознака (коефіцієнт > 0.5)")
print("2. MedInc^2 — квадратичний ефект доходу (нелінійне зростання ціни)")
print("3. Latitude, Longitude — географічне розташування критично важливе")
print("4. MedInc * Latitude — взаємодія доходу та широти (райони з високим доходом)")
print("5. AveRooms^2 — кількість кімнат має нелінійний вплив")