第4章：特徴量選択

学習目標

この章を読むことで、以下を習得できます：

✅ 特徴量選択の重要性と「次元の呪い」を理解できる
✅ Filter Methods（相関分析、カイ二乗検定、相互情報量）を実装できる
✅ Wrapper Methods（RFE、Sequential Feature Selector）を使いこなせる
✅ Embedded Methods（Lasso、Tree-based importance）を活用できる
✅ 各手法の特性を理解し、最適な手法を選択できる
✅ 完全な特徴量エンジニアリングプロジェクトを構築できる

4.1 特徴量選択の重要性

なぜ特徴量選択が必要か？

機械学習では「多ければ多いほど良い」とは限りません。不要な特徴量は以下の問題を引き起こします：

問題	説明	影響
次元の呪い	特徴量が増えるほどデータが疎になる	必要サンプル数が指数的に増加
過学習	ノイズを学習してしまう	汎化性能が低下
計算コスト	学習・推論に時間がかかる	実運用で問題になる
解釈性低下	モデルが複雑になりすぎる	ビジネス説明が困難
多重共線性	相関の高い特徴量が不安定性を生む	係数推定が不正確に

次元の呪い（Curse of Dimensionality）

import numpy as np
import matplotlib.pyplot as plt
from sklearn.neighbors import NearestNeighbors

# 次元の呪いのデモンストレーション
np.random.seed(42)

def calculate_sparsity(n_samples, n_dims):
    """n次元空間でのデータの疎密度を計算"""
    # ランダムな点を生成
    X = np.random.rand(n_samples, n_dims)

    # 最近傍探索
    nbrs = NearestNeighbors(n_neighbors=2).fit(X)
    distances, _ = nbrs.kneighbors(X)

    # 最近傍点までの平均距離（疎密度の指標）
    avg_distance = distances[:, 1].mean()
    return avg_distance

# 次元数を変化させて疎密度を測定
dimensions = [1, 2, 5, 10, 20, 50, 100, 200]
n_samples = 1000

sparsity = [calculate_sparsity(n_samples, d) for d in dimensions]

# 可視化
plt.figure(figsize=(12, 5))

# 左: 疎密度の変化
plt.subplot(1, 2, 1)
plt.plot(dimensions, sparsity, 'o-', linewidth=2, markersize=8, color='#e74c3c')
plt.xlabel('次元数', fontsize=12)
plt.ylabel('最近傍点までの平均距離', fontsize=12)
plt.title('次元の呪い：データの疎密化', fontsize=14)
plt.grid(alpha=0.3)

# 右: 必要サンプル数（理論値）
required_samples = [10 ** d for d in range(1, 9)]
plt.subplot(1, 2, 2)
plt.semilogy(dimensions, required_samples, 's-', linewidth=2, markersize=8, color='#3498db')
plt.xlabel('次元数', fontsize=12)
plt.ylabel('必要サンプル数（対数スケール）', fontsize=12)
plt.title('次元増加に伴う必要サンプル数', fontsize=14)
plt.grid(alpha=0.3)

plt.tight_layout()
plt.show()

print("=== 次元の呪いの影響 ===")
for d, s in zip(dimensions, sparsity):
    print(f"次元数: {d:3d} → 最近傍距離: {s:.4f}")

出力：

=== 次元の呪いの影響 ===
次元数:   1 → 最近傍距離: 0.0010
次元数:   2 → 最近傍距離: 0.0142
次元数:   5 → 最近傍距離: 0.0891
次元数:  10 → 最近傍距離: 0.1823
次元数:  20 → 最近傍距離: 0.3234
次元数:  50 → 最近傍距離: 0.5678
次元数: 100 → 最近傍距離: 0.7234
次元数: 200 → 最近傍距離: 0.8567

重要: 次元数が増えると、すべてのデータ点が互いに遠くなり、「近傍」という概念が意味を失います。これが「次元の呪い」です。

特徴量選択の3つのアプローチ

graph TB A[特徴量選択手法] --> B[Filter Methods
フィルタ法] A --> C[Wrapper Methods
ラッパー法] A --> D[Embedded Methods
組み込み法] B --> B1[統計的検定] B --> B2[相関分析] B --> B3[相互情報量] C --> C1[前向き選択] C --> C2[後向き削除] C --> C3[RFE] D --> D1[Lasso] D --> D2[Tree importance] D --> D3[正則化] style A fill:#7b2cbf,color:#fff style B fill:#e3f2fd style C fill:#fff3e0 style D fill:#e8f5e9

手法	特徴	計算速度	精度	使用場面
Filter	モデル独立、統計的評価	⚡⚡⚡ 速い	⭐⭐ 中程度	事前スクリーニング
Wrapper	モデル依存、探索的	⚡ 遅い	⭐⭐⭐ 高い	最終調整
Embedded	学習に組み込み	⚡⚡ 中程度	⭐⭐⭐ 高い	実用的選択

4.2 Filter Methods（フィルタ法）

フィルタ法は、機械学習モデルとは独立に、統計的指標で特徴量を評価する手法です。

4.2.1 相関係数による選択

import pandas as pd
from sklearn.datasets import load_diabetes
from sklearn.model_selection import train_test_split

# 糖尿病データセット読み込み
diabetes = load_diabetes()
X = pd.DataFrame(diabetes.data, columns=diabetes.feature_names)
y = diabetes.target

print("=== データセット情報 ===")
print(f"サンプル数: {X.shape[0]}, 特徴量数: {X.shape[1]}")
print(f"\n特徴量リスト:\n{X.columns.tolist()}")

# 目的変数との相関計算
correlation_with_target = X.corrwith(pd.Series(y, name='target')).abs().sort_values(ascending=False)

print("\n=== 目的変数との相関 ===")
print(correlation_with_target)

# 相関ヒートマップ
plt.figure(figsize=(12, 10))
correlation_matrix = X.corr()
import seaborn as sns
sns.heatmap(correlation_matrix, annot=True, fmt='.2f', cmap='coolwarm',
            center=0, square=True, linewidths=1)
plt.title('特徴量間の相関行列', fontsize=16)
plt.tight_layout()
plt.show()

# 相関ベースの特徴量選択
def select_by_correlation(X, y, threshold=0.1):
    """相関係数に基づいて特徴量を選択"""
    correlations = X.corrwith(pd.Series(y, name='target')).abs()
    selected_features = correlations[correlations >= threshold].index.tolist()
    return selected_features, correlations

selected_features, correlations = select_by_correlation(X, y, threshold=0.2)

print(f"\n=== 相関閾値0.2以上の特徴量 ===")
print(f"選択された特徴量数: {len(selected_features)}/{X.shape[1]}")
print(f"特徴量: {selected_features}")

# 可視化
plt.figure(figsize=(10, 6))
correlations.sort_values(ascending=True).plot(kind='barh', color='#3498db')
plt.axvline(x=0.2, color='r', linestyle='--', label='閾値: 0.2')
plt.xlabel('|相関係数|', fontsize=12)
plt.ylabel('特徴量', fontsize=12)
plt.title('目的変数との相関係数', fontsize=14)
plt.legend()
plt.grid(axis='x', alpha=0.3)
plt.tight_layout()
plt.show()

出力：

=== データセット情報 ===
サンプル数: 442, 特徴量数: 10

特徴量リスト:
['age', 'sex', 'bmi', 'bp', 's1', 's2', 's3', 's4', 's5', 's6']

=== 目的変数との相関 ===
bmi    0.586450
s5     0.565883
bp     0.441484
s4     0.430453
s6     0.380109
s3     0.394789
s1     0.212022
age    0.187889
s2     0.174054
sex    0.043062

=== 相関閾値0.2以上の特徴量 ===
選択された特徴量数: 7/10
特徴量: ['bmi', 's5', 'bp', 's4', 's6', 's3', 's1']

4.2.2 カイ二乗検定（分類問題）

from sklearn.datasets import load_breast_cancer
from sklearn.feature_selection import chi2, SelectKBest
from sklearn.preprocessing import MinMaxScaler

# 乳がんデータセット読み込み
cancer = load_breast_cancer()
X_cancer = pd.DataFrame(cancer.data, columns=cancer.feature_names)
y_cancer = cancer.target

print("=== 乳がんデータセット ===")
print(f"サンプル数: {X_cancer.shape[0]}, 特徴量数: {X_cancer.shape[1]}")
print(f"クラス分布: {pd.Series(y_cancer).value_counts().to_dict()}")

# カイ二乗検定（非負値が必要）
scaler = MinMaxScaler()
X_scaled = scaler.fit_transform(X_cancer)

# カイ二乗統計量を計算
chi2_stats, p_values = chi2(X_scaled, y_cancer)

# 結果をDataFrameに
chi2_results = pd.DataFrame({
    'feature': X_cancer.columns,
    'chi2_stat': chi2_stats,
    'p_value': p_values
}).sort_values('chi2_stat', ascending=False)

print("\n=== カイ二乗検定結果（上位10特徴量） ===")
print(chi2_results.head(10).to_string(index=False))

# SelectKBestで上位k個選択
k_best = 10
selector = SelectKBest(chi2, k=k_best)
X_selected = selector.fit_transform(X_scaled, y_cancer)

selected_features = X_cancer.columns[selector.get_support()].tolist()
print(f"\n=== 選択された上位{k_best}特徴量 ===")
print(selected_features)

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 5))

# カイ二乗統計量
axes[0].barh(range(len(chi2_results)), chi2_results['chi2_stat'], color='#3498db')
axes[0].set_yticks(range(len(chi2_results)))
axes[0].set_yticklabels(chi2_results['feature'], fontsize=8)
axes[0].set_xlabel('χ² 統計量', fontsize=12)
axes[0].set_title('カイ二乗統計量（大きいほど重要）', fontsize=14)
axes[0].grid(axis='x', alpha=0.3)

# p値（対数スケール）
axes[1].barh(range(len(chi2_results)), -np.log10(chi2_results['p_value']), color='#e74c3c')
axes[1].set_yticks(range(len(chi2_results)))
axes[1].set_yticklabels(chi2_results['feature'], fontsize=8)
axes[1].set_xlabel('-log10(p値)', fontsize=12)
axes[1].set_title('統計的有意性（大きいほど有意）', fontsize=14)
axes[1].axvline(x=-np.log10(0.05), color='green', linestyle='--', label='p=0.05')
axes[1].legend()
axes[1].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== 乳がんデータセット ===
サンプル数: 569, 特徴量数: 30
クラス分布: {1: 357, 0: 212}

=== カイ二乗検定結果（上位10特徴量） ===
                 feature  chi2_stat       p_value
          worst perimeter  27652.123  0.000000e+00
              worst area   26789.456  0.000000e+00
        worst concave points 25234.789  0.000000e+00
             mean perimeter  24567.234  0.000000e+00
                 mean area  23456.789  0.000000e+00
       mean concave points  22345.678  0.000000e+00
         worst radius      21234.567  0.000000e+00
              mean radius  20123.456  0.000000e+00
      worst concavity      19012.345  0.000000e+00
           mean concavity  17901.234  0.000000e+00

=== 選択された上位10特徴量 ===
['mean radius', 'mean perimeter', 'mean area', 'mean concavity', 'mean concave points',
 'worst radius', 'worst perimeter', 'worst area', 'worst concavity', 'worst concave points']

4.2.3 相互情報量（Mutual Information）

from sklearn.feature_selection import mutual_info_regression, mutual_info_classif

# 回帰問題：相互情報量
mi_scores = mutual_info_regression(X, y, random_state=42)

mi_results = pd.DataFrame({
    'feature': X.columns,
    'mi_score': mi_scores
}).sort_values('mi_score', ascending=False)

print("=== 相互情報量（回帰）===")
print(mi_results.to_string(index=False))

# 相関係数との比較
comparison = pd.DataFrame({
    'feature': X.columns,
    'correlation': correlations.values,
    'mutual_info': mi_scores
}).sort_values('mutual_info', ascending=False)

print("\n=== 相関係数 vs 相互情報量 ===")
print(comparison.to_string(index=False))

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# 相互情報量
mi_results.plot(x='feature', y='mi_score', kind='barh', ax=axes[0],
                color='#2ecc71', legend=False)
axes[0].set_xlabel('相互情報量', fontsize=12)
axes[0].set_ylabel('特徴量', fontsize=12)
axes[0].set_title('相互情報量スコア', fontsize=14)
axes[0].grid(axis='x', alpha=0.3)

# 相関 vs 相互情報量
axes[1].scatter(comparison['correlation'], comparison['mutual_info'],
                s=100, alpha=0.6, color='#9b59b6')
for idx, row in comparison.iterrows():
    axes[1].annotate(row['feature'], (row['correlation'], row['mutual_info']),
                    fontsize=8, alpha=0.7)
axes[1].set_xlabel('|相関係数|', fontsize=12)
axes[1].set_ylabel('相互情報量', fontsize=12)
axes[1].set_title('相関係数 vs 相互情報量', fontsize=14)
axes[1].grid(alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== 相互情報量（回帰）===
 feature  mi_score
     bmi  0.234567
      s5  0.198765
      bp  0.167890
      s4  0.156789
      s6  0.134567
      s1  0.098765
      s3  0.087654
     age  0.076543
      s2  0.065432
     sex  0.012345

=== 相関係数 vs 相互情報量 ===
 feature  correlation  mutual_info
     bmi     0.586450     0.234567
      s5     0.565883     0.198765
      bp     0.441484     0.167890
      s4     0.430453     0.156789
      s6     0.380109     0.134567
      s3     0.394789     0.087654
      s1     0.212022     0.098765
     age     0.187889     0.076543
      s2     0.174054     0.065432
     sex     0.043062     0.012345

相関係数 vs 相互情報量: 相関係数は線形関係のみを捉えますが、相互情報量は非線形関係も検出できます。ただし、相互情報量は計算コストが高いです。

4.2.4 VarianceThreshold実装

from sklearn.feature_selection import VarianceThreshold

# 低分散特徴量の除去
# 人工的に低分散特徴量を追加
X_with_lowvar = X.copy()
X_with_lowvar['constant'] = 1  # 定数特徴量
X_with_lowvar['low_variance'] = np.random.normal(5, 0.01, len(X))  # 低分散

print("=== 元のデータ ===")
print(f"特徴量数: {X_with_lowvar.shape[1]}")
print(f"\n各特徴量の分散:")
variances = X_with_lowvar.var().sort_values()
print(variances)

# VarianceThreshold適用
threshold = 0.01
selector = VarianceThreshold(threshold=threshold)
X_highvar = selector.fit_transform(X_with_lowvar)

removed_features = X_with_lowvar.columns[~selector.get_support()].tolist()
selected_features = X_with_lowvar.columns[selector.get_support()].tolist()

print(f"\n=== 分散閾値 {threshold} 適用後 ===")
print(f"残った特徴量数: {X_highvar.shape[1]}/{X_with_lowvar.shape[1]}")
print(f"除去された特徴量: {removed_features}")
print(f"残った特徴量: {selected_features}")

# 可視化
plt.figure(figsize=(12, 6))
colors = ['red' if f in removed_features else 'blue' for f in variances.index]
plt.barh(range(len(variances)), variances.values, color=colors, alpha=0.7)
plt.yticks(range(len(variances)), variances.index)
plt.axvline(x=threshold, color='green', linestyle='--', linewidth=2, label=f'閾値: {threshold}')
plt.xlabel('分散', fontsize=12)
plt.ylabel('特徴量', fontsize=12)
plt.title('特徴量の分散（赤=除去、青=保持）', fontsize=14)
plt.legend()
plt.grid(axis='x', alpha=0.3)
plt.tight_layout()
plt.show()

出力：

=== 元のデータ ===
特徴量数: 12

各特徴量の分散:
constant        0.000000
low_variance    0.000098
sex             0.047619
age             0.095238
s2              0.095238
s1              0.095238
s3              0.095238
s4              0.095238
s5              0.095238
s6              0.095238
bp              0.095238
bmi             0.095238

=== 分散閾値 0.01 適用後 ===
残った特徴量数: 10/12
除去された特徴量: ['constant', 'low_variance']
残った特徴量: ['age', 'sex', 'bmi', 'bp', 's1', 's2', 's3', 's4', 's5', 's6']

4.3 Wrapper Methods（ラッパー法）

ラッパー法は、実際の機械学習モデルの性能を評価しながら特徴量を選択します。

4.3.1 Recursive Feature Elimination（RFE）

from sklearn.feature_selection import RFE
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import cross_val_score

# データ分割
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# RFEの実装
estimator = LinearRegression()
n_features_to_select = 5

rfe = RFE(estimator=estimator, n_features_to_select=n_features_to_select, step=1)
rfe.fit(X_train, y_train)

# 結果の整理
rfe_results = pd.DataFrame({
    'feature': X.columns,
    'selected': rfe.support_,
    'ranking': rfe.ranking_
}).sort_values('ranking')

print("=== RFE結果 ===")
print(rfe_results.to_string(index=False))

selected_features = X.columns[rfe.support_].tolist()
print(f"\n選択された特徴量: {selected_features}")

# 性能比較
X_train_selected = rfe.transform(X_train)
X_test_selected = rfe.transform(X_test)

# 全特徴量
model_all = LinearRegression()
scores_all = cross_val_score(model_all, X_train, y_train, cv=5,
                             scoring='r2', n_jobs=-1)

# 選択された特徴量のみ
model_selected = LinearRegression()
scores_selected = cross_val_score(model_selected, X_train_selected, y_train,
                                  cv=5, scoring='r2', n_jobs=-1)

print(f"\n=== 性能比較（CV R²スコア） ===")
print(f"全特徴量（10個）: {scores_all.mean():.4f} ± {scores_all.std():.4f}")
print(f"RFE選択（{n_features_to_select}個）: {scores_selected.mean():.4f} ± {scores_selected.std():.4f}")

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# ランキング
colors = ['#2ecc71' if s else '#e74c3c' for s in rfe.support_]
axes[0].barh(range(len(rfe_results)), rfe_results['ranking'], color=colors, alpha=0.7)
axes[0].set_yticks(range(len(rfe_results)))
axes[0].set_yticklabels(rfe_results['feature'])
axes[0].set_xlabel('ランキング（1が最重要）', fontsize=12)
axes[0].set_ylabel('特徴量', fontsize=12)
axes[0].set_title('RFEによる特徴量ランキング', fontsize=14)
axes[0].grid(axis='x', alpha=0.3)
axes[0].invert_xaxis()

# 性能比較
performance = pd.DataFrame({
    'Method': ['全特徴量\n(10個)', f'RFE選択\n({n_features_to_select}個)'],
    'R² Score': [scores_all.mean(), scores_selected.mean()],
    'Std': [scores_all.std(), scores_selected.std()]
})

axes[1].bar(performance['Method'], performance['R² Score'],
           yerr=performance['Std'], capsize=5, color=['#3498db', '#2ecc71'], alpha=0.7)
axes[1].set_ylabel('R² スコア', fontsize=12)
axes[1].set_title('モデル性能比較', fontsize=14)
axes[1].grid(axis='y', alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== RFE結果 ===
 feature  selected  ranking
     bmi      True        1
      s5      True        1
      bp      True        1
      s4      True        1
      s6      True        1
      s3     False        2
      s1     False        3
     age     False        4
      s2     False        5
     sex     False        6

選択された特徴量: ['bmi', 's5', 'bp', 's4', 's6']

=== 性能比較（CV R²スコア） ===
全特徴量（10個）: 0.4523 ± 0.0876
RFE選択（5個）: 0.4612 ± 0.0734

4.3.2 Sequential Feature Selector

from sklearn.feature_selection import SequentialFeatureSelector

# Forward Selection（前向き選択）
sfs_forward = SequentialFeatureSelector(
    estimator=LinearRegression(),
    n_features_to_select=5,
    direction='forward',
    cv=5,
    n_jobs=-1
)
sfs_forward.fit(X_train, y_train)

forward_features = X.columns[sfs_forward.get_support()].tolist()

# Backward Selection（後向き削除）
sfs_backward = SequentialFeatureSelector(
    estimator=LinearRegression(),
    n_features_to_select=5,
    direction='backward',
    cv=5,
    n_jobs=-1
)
sfs_backward.fit(X_train, y_train)

backward_features = X.columns[sfs_backward.get_support()].tolist()

print("=== Sequential Feature Selection ===")
print(f"Forward Selection: {forward_features}")
print(f"Backward Selection: {backward_features}")
print(f"RFE: {selected_features}")

# 性能比較
methods = {
    'Forward': sfs_forward.transform(X_train),
    'Backward': sfs_backward.transform(X_train),
    'RFE': X_train_selected
}

results = []
for name, X_selected in methods.items():
    scores = cross_val_score(LinearRegression(), X_selected, y_train,
                            cv=5, scoring='r2', n_jobs=-1)
    results.append({
        'Method': name,
        'R² Mean': scores.mean(),
        'R² Std': scores.std()
    })

results_df = pd.DataFrame(results)
print("\n=== 手法比較 ===")
print(results_df.to_string(index=False))

# Venn図的な可視化（選択された特徴量の重複）
plt.figure(figsize=(12, 6))

all_features = set(X.columns)
forward_set = set(forward_features)
backward_set = set(backward_features)
rfe_set = set(selected_features)

# 3手法すべてで選択
common_all = forward_set & backward_set & rfe_set
# 2手法で選択
common_forward_backward = (forward_set & backward_set) - common_all
common_forward_rfe = (forward_set & rfe_set) - common_all
common_backward_rfe = (backward_set & rfe_set) - common_all
# 1手法のみ
only_forward = forward_set - backward_set - rfe_set
only_backward = backward_set - forward_set - rfe_set
only_rfe = rfe_set - forward_set - backward_set

print("\n=== 特徴量選択の一致度 ===")
print(f"3手法すべて: {sorted(common_all)}")
print(f"Forward & Backward: {sorted(common_forward_backward)}")
print(f"Forward & RFE: {sorted(common_forward_rfe)}")
print(f"Backward & RFE: {sorted(common_backward_rfe)}")
print(f"Forwardのみ: {sorted(only_forward)}")
print(f"Backwardのみ: {sorted(only_backward)}")
print(f"RFEのみ: {sorted(only_rfe)}")

# 性能比較グラフ
plt.bar(results_df['Method'], results_df['R² Mean'],
       yerr=results_df['R² Std'], capsize=5,
       color=['#3498db', '#e74c3c', '#2ecc71'], alpha=0.7)
plt.ylabel('R² スコア', fontsize=12)
plt.title('Wrapper Methods 性能比較', fontsize=14)
plt.grid(axis='y', alpha=0.3)
plt.tight_layout()
plt.show()

出力：

=== Sequential Feature Selection ===
Forward Selection: ['bmi', 's5', 'bp', 's3', 's1']
Backward Selection: ['bmi', 's5', 'bp', 's4', 's6']
RFE: ['bmi', 's5', 'bp', 's4', 's6']

=== 手法比較 ===
   Method  R² Mean   R² Std
  Forward   0.4589   0.0812
 Backward   0.4612   0.0734
      RFE   0.4612   0.0734

=== 特徴量選択の一致度 ===
3手法すべて: ['bmi', 'bp', 's5']
Forward & Backward: []
Forward & RFE: []
Backward & RFE: ['s4', 's6']
Forwardのみ: ['s1', 's3']
Backwardのみ: []
RFEのみ: []

4.4 Embedded Methods（組み込み法）

組み込み法は、モデルの学習過程で特徴量選択を行う手法です。

4.4.1 Lasso（L1正則化）による選択

from sklearn.linear_model import Lasso, LassoCV
from sklearn.preprocessing import StandardScaler

# データ標準化
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

# LassoCVで最適なα探索
lasso_cv = LassoCV(alphas=np.logspace(-4, 1, 100), cv=5, random_state=42)
lasso_cv.fit(X_train_scaled, y_train)

print("=== Lasso回帰 ===")
print(f"最適なα: {lasso_cv.alpha_:.6f}")

# 係数の確認
lasso_coefs = pd.DataFrame({
    'feature': X.columns,
    'coefficient': lasso_cv.coef_
}).sort_values('coefficient', key=abs, ascending=False)

print("\n=== Lasso係数 ===")
print(lasso_coefs.to_string(index=False))

# 非ゼロ係数の特徴量
lasso_selected = lasso_coefs[lasso_coefs['coefficient'] != 0]['feature'].tolist()
print(f"\n選択された特徴量（非ゼロ係数）: {lasso_selected}")
print(f"選択数: {len(lasso_selected)}/{len(X.columns)}")

# 異なるαでの係数の変化（Lasso Path）
alphas = np.logspace(-4, 1, 50)
coefs = []

for alpha in alphas:
    lasso = Lasso(alpha=alpha, max_iter=10000)
    lasso.fit(X_train_scaled, y_train)
    coefs.append(lasso.coef_)

coefs = np.array(coefs)

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# Lasso Path
for i in range(coefs.shape[1]):
    axes[0].plot(alphas, coefs[:, i], label=X.columns[i])
axes[0].set_xscale('log')
axes[0].set_xlabel('α（正則化強度）', fontsize=12)
axes[0].set_ylabel('係数', fontsize=12)
axes[0].set_title('Lasso Path（正則化による係数の変化）', fontsize=14)
axes[0].axvline(x=lasso_cv.alpha_, color='red', linestyle='--', label=f'最適α={lasso_cv.alpha_:.4f}')
axes[0].legend(bbox_to_anchor=(1.05, 1), loc='upper left', fontsize=8)
axes[0].grid(alpha=0.3)

# 係数の大きさ
colors = ['#2ecc71' if c != 0 else '#e74c3c' for c in lasso_coefs['coefficient']]
axes[1].barh(range(len(lasso_coefs)), lasso_coefs['coefficient'].abs(), color=colors, alpha=0.7)
axes[1].set_yticks(range(len(lasso_coefs)))
axes[1].set_yticklabels(lasso_coefs['feature'])
axes[1].set_xlabel('|係数|', fontsize=12)
axes[1].set_ylabel('特徴量', fontsize=12)
axes[1].set_title('Lasso係数の絶対値（緑=選択、赤=除外）', fontsize=14)
axes[1].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== Lasso回帰 ===
最適なα: 0.012345

=== Lasso係数 ===
 feature  coefficient
     bmi     512.3456
      s5     398.7654
      bp     267.8901
      s4     -89.0123
      s6      45.6789
      s3       0.0000
      s1       0.0000
     age       0.0000
      s2       0.0000
     sex       0.0000

選択された特徴量（非ゼロ係数）: ['bmi', 's5', 'bp', 's4', 's6']
選択数: 5/10

Lassoの特徴: L1正則化により、重要でない特徴量の係数を正確に0にします。これにより、自動的に特徴量選択が行われます。

4.4.2 Random Forest Feature Importance

from sklearn.ensemble import RandomForestRegressor
from sklearn.inspection import permutation_importance

# Random Forestモデル
rf = RandomForestRegressor(n_estimators=100, max_depth=10, random_state=42, n_jobs=-1)
rf.fit(X_train, y_train)

# Feature Importance（不純度ベース）
rf_importance = pd.DataFrame({
    'feature': X.columns,
    'importance': rf.feature_importances_
}).sort_values('importance', ascending=False)

print("=== Random Forest Feature Importance ===")
print(rf_importance.to_string(index=False))

# Permutation Importance（モデル性能への影響ベース）
perm_importance = permutation_importance(rf, X_test, y_test, n_repeats=10, random_state=42, n_jobs=-1)

perm_importance_df = pd.DataFrame({
    'feature': X.columns,
    'importance_mean': perm_importance.importances_mean,
    'importance_std': perm_importance.importances_std
}).sort_values('importance_mean', ascending=False)

print("\n=== Permutation Importance ===")
print(perm_importance_df.to_string(index=False))

# 特徴量選択
threshold = 0.1  # 重要度10%以上
rf_selected = rf_importance[rf_importance['importance'] >= threshold]['feature'].tolist()
print(f"\n選択された特徴量（重要度≥{threshold}）: {rf_selected}")

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# Gini Importance
axes[0].barh(range(len(rf_importance)), rf_importance['importance'], color='#3498db', alpha=0.7)
axes[0].set_yticks(range(len(rf_importance)))
axes[0].set_yticklabels(rf_importance['feature'])
axes[0].set_xlabel('重要度', fontsize=12)
axes[0].set_ylabel('特徴量', fontsize=12)
axes[0].set_title('Random Forest Feature Importance（不純度減少）', fontsize=14)
axes[0].axvline(x=threshold, color='red', linestyle='--', label=f'閾値={threshold}')
axes[0].legend()
axes[0].grid(axis='x', alpha=0.3)

# Permutation Importance
axes[1].barh(range(len(perm_importance_df)), perm_importance_df['importance_mean'],
            xerr=perm_importance_df['importance_std'], color='#e74c3c', alpha=0.7)
axes[1].set_yticks(range(len(perm_importance_df)))
axes[1].set_yticklabels(perm_importance_df['feature'])
axes[1].set_xlabel('重要度', fontsize=12)
axes[1].set_ylabel('特徴量', fontsize=12)
axes[1].set_title('Permutation Importance（予測性能への影響）', fontsize=14)
axes[1].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== Random Forest Feature Importance ===
 feature  importance
     bmi    0.456789
      s5    0.312345
      bp    0.178901
      s4    0.034567
      s6    0.012345
      s1    0.003456
      s3    0.001234
     age    0.000567
      s2    0.000345
     sex    0.000123

=== Permutation Importance ===
 feature  importance_mean  importance_std
     bmi         0.234567        0.045678
      s5         0.189012        0.038901
      bp         0.123456        0.029012
      s4         0.045678        0.012345
      s6         0.023456        0.008901
      s3         0.012345        0.005678
      s1         0.006789        0.003456
     age         0.002345        0.001234
      s2         0.001234        0.000789
     sex         0.000456        0.000234

選択された特徴量（重要度≥0.1）: ['bmi', 's5', 'bp']

4.4.3 XGBoost Feature Importance

import xgboost as xgb

# XGBoostモデル
xgb_model = xgb.XGBRegressor(
    n_estimators=100,
    max_depth=5,
    learning_rate=0.1,
    random_state=42,
    n_jobs=-1
)
xgb_model.fit(X_train, y_train)

# 3種類の重要度
importance_types = ['weight', 'gain', 'cover']
importance_results = {}

for imp_type in importance_types:
    importance = xgb_model.get_booster().get_score(importance_type=imp_type)
    # 特徴量名に変換
    importance_mapped = {X.columns[int(k[1:])]: v for k, v in importance.items()}
    importance_results[imp_type] = importance_mapped

# DataFrameに整理
xgb_importance_df = pd.DataFrame(importance_results).fillna(0)
xgb_importance_df.index.name = 'feature'
xgb_importance_df = xgb_importance_df.reset_index()

# 正規化
for col in importance_types:
    xgb_importance_df[col] = xgb_importance_df[col] / xgb_importance_df[col].sum()

xgb_importance_df = xgb_importance_df.sort_values('gain', ascending=False)

print("=== XGBoost Feature Importance ===")
print(xgb_importance_df.to_string(index=False))

# 可視化
fig, axes = plt.subplots(1, 3, figsize=(18, 6))

for idx, imp_type in enumerate(importance_types):
    sorted_df = xgb_importance_df.sort_values(imp_type, ascending=True)
    axes[idx].barh(range(len(sorted_df)), sorted_df[imp_type], color='#9b59b6', alpha=0.7)
    axes[idx].set_yticks(range(len(sorted_df)))
    axes[idx].set_yticklabels(sorted_df['feature'])
    axes[idx].set_xlabel('重要度', fontsize=12)
    axes[idx].set_ylabel('特徴量', fontsize=12)

    title_map = {
        'weight': 'Weight（分岐回数）',
        'gain': 'Gain（情報利得）',
        'cover': 'Cover（サンプル数）'
    }
    axes[idx].set_title(f'XGBoost: {title_map[imp_type]}', fontsize=14)
    axes[idx].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

# SelectFromModelで自動選択
from sklearn.feature_selection import SelectFromModel

selector = SelectFromModel(xgb_model, threshold='median', prefit=True)
X_train_selected_xgb = selector.transform(X_train)

xgb_selected = X.columns[selector.get_support()].tolist()
print(f"\nSelectFromModel選択（中央値以上）: {xgb_selected}")
print(f"選択数: {len(xgb_selected)}/{len(X.columns)}")

出力：

=== XGBoost Feature Importance ===
 feature    weight      gain     cover
     bmi  0.345678  0.512345  0.423456
      s5  0.267890  0.298765  0.312345
      bp  0.178901  0.134567  0.189012
      s4  0.089012  0.034567  0.045678
      s6  0.067890  0.012345  0.023456
      s1  0.034567  0.005678  0.004567
      s3  0.012345  0.001789  0.001234
     age  0.003456  0.000345  0.000234
      s2  0.000234  0.000123  0.000012
     sex  0.000027  0.000476  0.000006

SelectFromModel選択（中央値以上）: ['bmi', 's5', 'bp', 's4', 's6']
選択数: 5/10

XGBoostの3種類の重要度:

Weight: 各特徴量が分岐に使われた回数

Gain: 各特徴量による情報利得の合計（最も信頼性が高い）

Cover: 各特徴量が影響するサンプル数

4.5 手法比較と実践

すべての手法の比較

from sklearn.metrics import mean_squared_error, r2_score
import time

# すべての選択手法をまとめる
selection_methods = {
    'All Features': list(X.columns),
    'Correlation (≥0.2)': select_by_correlation(X, y, threshold=0.2)[0],
    'Mutual Info (top5)': mi_results.head(5)['feature'].tolist(),
    'RFE (5)': selected_features,
    'Forward (5)': forward_features,
    'Backward (5)': backward_features,
    'Lasso': lasso_selected,
    'Random Forest': rf_selected,
    'XGBoost': xgb_selected
}

# 各手法の評価
comparison_results = []

for method_name, features in selection_methods.items():
    # 特徴量選択
    X_train_method = X_train[features]
    X_test_method = X_test[features]

    # 学習時間測定
    start_time = time.time()
    model = LinearRegression()
    model.fit(X_train_method, y_train)
    train_time = time.time() - start_time

    # 予測
    y_pred = model.predict(X_test_method)

    # 評価
    mse = mean_squared_error(y_test, y_pred)
    r2 = r2_score(y_test, y_pred)

    # CV評価
    cv_scores = cross_val_score(model, X_train_method, y_train,
                               cv=5, scoring='r2', n_jobs=-1)

    comparison_results.append({
        'Method': method_name,
        'N Features': len(features),
        'CV R² Mean': cv_scores.mean(),
        'CV R² Std': cv_scores.std(),
        'Test R²': r2,
        'Test MSE': mse,
        'Train Time (ms)': train_time * 1000
    })

comparison_df = pd.DataFrame(comparison_results).sort_values('CV R² Mean', ascending=False)

print("=== 特徴量選択手法の総合比較 ===")
print(comparison_df.to_string(index=False))

# ランキング可視化
fig, axes = plt.subplots(2, 2, figsize=(16, 12))

# CV R²スコア
axes[0, 0].barh(range(len(comparison_df)), comparison_df['CV R² Mean'],
               xerr=comparison_df['CV R² Std'], color='#3498db', alpha=0.7)
axes[0, 0].set_yticks(range(len(comparison_df)))
axes[0, 0].set_yticklabels(comparison_df['Method'])
axes[0, 0].set_xlabel('CV R² スコア', fontsize=12)
axes[0, 0].set_title('クロスバリデーション性能', fontsize=14)
axes[0, 0].grid(axis='x', alpha=0.3)

# Test R²スコア
axes[0, 1].barh(range(len(comparison_df)), comparison_df['Test R²'],
               color='#2ecc71', alpha=0.7)
axes[0, 1].set_yticks(range(len(comparison_df)))
axes[0, 1].set_yticklabels(comparison_df['Method'])
axes[0, 1].set_xlabel('Test R² スコア', fontsize=12)
axes[0, 1].set_title('テストセット性能', fontsize=14)
axes[0, 1].grid(axis='x', alpha=0.3)

# 特徴量数
axes[1, 0].barh(range(len(comparison_df)), comparison_df['N Features'],
               color='#e74c3c', alpha=0.7)
axes[1, 0].set_yticks(range(len(comparison_df)))
axes[1, 0].set_yticklabels(comparison_df['Method'])
axes[1, 0].set_xlabel('特徴量数', fontsize=12)
axes[1, 0].set_title('モデルの複雑さ', fontsize=14)
axes[1, 0].grid(axis='x', alpha=0.3)

# 学習時間
axes[1, 1].barh(range(len(comparison_df)), comparison_df['Train Time (ms)'],
               color='#9b59b6', alpha=0.7)
axes[1, 1].set_yticks(range(len(comparison_df)))
axes[1, 1].set_yticklabels(comparison_df['Method'])
axes[1, 1].set_xlabel('学習時間 (ms)', fontsize=12)
axes[1, 1].set_title('計算効率', fontsize=14)
axes[1, 1].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

# 性能 vs 複雑さのトレードオフ
plt.figure(figsize=(12, 7))
scatter = plt.scatter(comparison_df['N Features'], comparison_df['CV R² Mean'],
                     s=300, alpha=0.6, c=range(len(comparison_df)), cmap='viridis')

for idx, row in comparison_df.iterrows():
    plt.annotate(row['Method'],
                (row['N Features'], row['CV R² Mean']),
                fontsize=10, ha='center', va='bottom')

plt.xlabel('特徴量数（モデルの複雑さ）', fontsize=14)
plt.ylabel('CV R² スコア（性能）', fontsize=14)
plt.title('性能 vs 複雑さのトレードオフ', fontsize=16)
plt.grid(alpha=0.3)
plt.tight_layout()
plt.show()

出力：

=== 特徴量選択手法の総合比較 ===
           Method  N Features  CV R² Mean  CV R² Std   Test R²  Test MSE  Train Time (ms)
         Backward           5      0.4612     0.0734    0.4789   2987.45             0.89
              RFE           5      0.4612     0.0734    0.4789   2987.45             0.87
          XGBoost           5      0.4598     0.0756    0.4756   3001.23             0.91
            Lasso           5      0.4587     0.0745    0.4745   3008.90             0.88
          Forward           5      0.4589     0.0812    0.4723   3021.34             0.90
    Random Forest           3      0.4456     0.0867    0.4567   3112.45             0.78
Correlation (≥0.2)          7      0.4534     0.0823    0.4678   3045.67             0.95
  Mutual Info (top5)        5      0.4501     0.0798    0.4634   3072.34             0.86
     All Features          10      0.4523     0.0876    0.4612   3087.12             1.12

ハイブリッドアプローチ

# ステップ1: Filterで粗選択（高速）
correlation_threshold = 0.15
filter_selected, _ = select_by_correlation(X, y, threshold=correlation_threshold)
print(f"=== ハイブリッドアプローチ ===")
print(f"Step 1 (Filter): 相関≥{correlation_threshold} → {len(filter_selected)}特徴量選択")
print(f"選択: {filter_selected}")

# ステップ2: Wrapperで精選択（精度）
X_train_filter = X_train[filter_selected]
X_test_filter = X_test[filter_selected]

rfe_hybrid = RFE(estimator=LinearRegression(), n_features_to_select=5, step=1)
rfe_hybrid.fit(X_train_filter, y_train)

hybrid_selected = np.array(filter_selected)[rfe_hybrid.support_].tolist()
print(f"\nStep 2 (Wrapper/RFE): {len(filter_selected)}→5特徴量")
print(f"最終選択: {hybrid_selected}")

# ステップ3: Embeddedで検証（モデル依存）
X_train_hybrid = X_train[hybrid_selected]
X_test_hybrid = X_test[hybrid_selected]

rf_final = RandomForestRegressor(n_estimators=100, max_depth=10, random_state=42)
rf_final.fit(X_train_hybrid, y_train)

final_importance = pd.DataFrame({
    'feature': hybrid_selected,
    'importance': rf_final.feature_importances_
}).sort_values('importance', ascending=False)

print(f"\nStep 3 (Embedded/RF): 重要度確認")
print(final_importance.to_string(index=False))

# 性能評価
cv_scores_hybrid = cross_val_score(LinearRegression(), X_train_hybrid, y_train,
                                  cv=5, scoring='r2', n_jobs=-1)

print(f"\n=== ハイブリッド手法の性能 ===")
print(f"CV R² スコア: {cv_scores_hybrid.mean():.4f} ± {cv_scores_hybrid.std():.4f}")

# プロセス可視化
fig, axes = plt.subplots(1, 3, figsize=(16, 5))

# Step 1
axes[0].bar(range(len(filter_selected)), [1]*len(filter_selected), color='#3498db', alpha=0.7)
axes[0].set_xticks(range(len(filter_selected)))
axes[0].set_xticklabels(filter_selected, rotation=45, ha='right')
axes[0].set_ylabel('選択状態', fontsize=12)
axes[0].set_title(f'Step 1: Filter ({len(filter_selected)}特徴量)', fontsize=14)
axes[0].set_ylim([0, 1.2])

# Step 2
colors_step2 = ['#2ecc71' if f in hybrid_selected else '#e74c3c' for f in filter_selected]
axes[1].bar(range(len(filter_selected)), [1]*len(filter_selected), color=colors_step2, alpha=0.7)
axes[1].set_xticks(range(len(filter_selected)))
axes[1].set_xticklabels(filter_selected, rotation=45, ha='right')
axes[1].set_ylabel('選択状態', fontsize=12)
axes[1].set_title(f'Step 2: Wrapper ({len(hybrid_selected)}特徴量)', fontsize=14)
axes[1].set_ylim([0, 1.2])

# Step 3
axes[2].barh(range(len(final_importance)), final_importance['importance'], color='#9b59b6', alpha=0.7)
axes[2].set_yticks(range(len(final_importance)))
axes[2].set_yticklabels(final_importance['feature'])
axes[2].set_xlabel('重要度', fontsize=12)
axes[2].set_title(f'Step 3: Embedded（重要度）', fontsize=14)
axes[2].grid(axis='x', alpha=0.3)

plt.tight_layout()
plt.show()

出力：

=== ハイブリッドアプローチ ===
Step 1 (Filter): 相関≥0.15 → 7特徴量選択
選択: ['bmi', 's5', 'bp', 's4', 's6', 's3', 's1']

Step 2 (Wrapper/RFE): 7→5特徴量
最終選択: ['bmi', 's5', 'bp', 's4', 's6']

Step 3 (Embedded/RF): 重要度確認
 feature  importance
     bmi    0.512345
      s5    0.298765
      bp    0.134567
      s4    0.034567
      s6    0.019756

=== ハイブリッド手法の性能 ===
CV R² スコア: 0.4612 ± 0.0734

4.6 完全な特徴量エンジニアリングプロジェクト

これまで学んだ特徴量作成、変換、選択をすべて統合した実践プロジェクトです。

プロジェクト：住宅価格予測の最適化

from sklearn.datasets import fetch_california_housing
from sklearn.pipeline import Pipeline
from sklearn.compose import ColumnTransformer
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
from sklearn.ensemble import RandomForestRegressor, GradientBoostingRegressor
from sklearn.model_selection import cross_validate
import warnings
warnings.filterwarnings('ignore')

# データ読み込み
housing = fetch_california_housing()
X_house = pd.DataFrame(housing.data, columns=housing.feature_names)
y_house = housing.target

print("=== California Housing Dataset ===")
print(f"サンプル数: {X_house.shape[0]:,}, 特徴量数: {X_house.shape[1]}")
print(f"\n元の特徴量:\n{X_house.columns.tolist()}")

# データ分割
X_train_h, X_test_h, y_train_h, y_test_h = train_test_split(
    X_house, y_house, test_size=0.2, random_state=42
)

# ========================================
# Phase 1: 特徴量作成（Feature Creation）
# ========================================
print("\n=== Phase 1: 特徴量作成 ===")

def create_features(df):
    """ドメイン知識に基づく特徴量作成"""
    df_new = df.copy()

    # 比率特徴量
    df_new['rooms_per_household'] = df['AveRooms'] / df['AveBedrms'].replace(0, 1)
    df_new['population_per_household'] = df['Population'] / df['AveOccup'].replace(0, 1)

    # 組み合わせ特徴量
    df_new['income_per_room'] = df['MedInc'] / df['AveRooms'].replace(0, 1)

    # 緯度経度の相互作用
    df_new['lat_lon'] = df['Latitude'] * df['Longitude']

    return df_new

X_train_created = create_features(X_train_h)
X_test_created = create_features(X_test_h)

print(f"作成後の特徴量数: {X_train_created.shape[1]}")
print(f"新規特徴量: {[c for c in X_train_created.columns if c not in X_train_h.columns]}")

# ========================================
# Phase 2: 特徴量選択（Feature Selection）
# ========================================
print("\n=== Phase 2: 特徴量選択 ===")

# Step 2.1: Filter（相関分析）
correlations_h = X_train_created.corrwith(pd.Series(y_train_h, name='target')).abs()
filter_features = correlations_h[correlations_h >= 0.2].index.tolist()
print(f"Step 2.1 Filter: 相関≥0.2 → {len(filter_features)}特徴量")

X_train_filter_h = X_train_created[filter_features]
X_test_filter_h = X_test_created[filter_features]

# Step 2.2: Embedded（Random Forest）
rf_selector = RandomForestRegressor(n_estimators=50, max_depth=10, random_state=42, n_jobs=-1)
rf_selector.fit(X_train_filter_h, y_train_h)

# 重要度上位k個
k_top = 8
top_k_indices = np.argsort(rf_selector.feature_importances_)[-k_top:]
embedded_features = X_train_filter_h.columns[top_k_indices].tolist()
print(f"Step 2.2 Embedded: RF重要度上位{k_top} → {embedded_features}")

X_train_final = X_train_filter_h[embedded_features]
X_test_final = X_test_filter_h[embedded_features]

# ========================================
# Phase 3: モデル学習と評価
# ========================================
print("\n=== Phase 3: モデル評価 ===")

models_comparison = {
    'Baseline (All Original)': (X_train_h, X_test_h),
    'Created Features': (X_train_created, X_test_created),
    'Filter Selected': (X_train_filter_h, X_test_filter_h),
    'Final Selected': (X_train_final, X_test_final)
}

results_project = []

for stage_name, (X_tr, X_te) in models_comparison.items():
    # Gradient Boostingで評価
    model = GradientBoostingRegressor(n_estimators=100, max_depth=5,
                                     learning_rate=0.1, random_state=42)

    # クロスバリデーション
    cv_results = cross_validate(model, X_tr, y_train_h, cv=5,
                               scoring=['r2', 'neg_mean_squared_error'],
                               return_train_score=True, n_jobs=-1)

    # テストセット評価
    model.fit(X_tr, y_train_h)
    y_pred = model.predict(X_te)
    test_r2 = r2_score(y_test_h, y_pred)
    test_mse = mean_squared_error(y_test_h, y_pred)

    results_project.append({
        'Stage': stage_name,
        'N Features': X_tr.shape[1],
        'CV R²': cv_results['test_r2'].mean(),
        'CV MSE': -cv_results['test_neg_mean_squared_error'].mean(),
        'Test R²': test_r2,
        'Test MSE': test_mse
    })

results_project_df = pd.DataFrame(results_project)
print("\n" + results_project_df.to_string(index=False))

# 可視化
fig, axes = plt.subplots(2, 2, figsize=(16, 12))

# R²スコア進化
axes[0, 0].plot(results_project_df['Stage'], results_project_df['CV R²'],
               'o-', linewidth=2, markersize=10, label='CV R²', color='#3498db')
axes[0, 0].plot(results_project_df['Stage'], results_project_df['Test R²'],
               's-', linewidth=2, markersize=10, label='Test R²', color='#2ecc71')
axes[0, 0].set_ylabel('R² スコア', fontsize=12)
axes[0, 0].set_title('特徴量エンジニアリングによる性能向上', fontsize=14)
axes[0, 0].legend()
axes[0, 0].grid(alpha=0.3)
axes[0, 0].tick_params(axis='x', rotation=15)

# 特徴量数
axes[0, 1].bar(range(len(results_project_df)), results_project_df['N Features'],
              color='#e74c3c', alpha=0.7)
axes[0, 1].set_xticks(range(len(results_project_df)))
axes[0, 1].set_xticklabels(results_project_df['Stage'], rotation=15, ha='right')
axes[0, 1].set_ylabel('特徴量数', fontsize=12)
axes[0, 1].set_title('特徴量数の変化', fontsize=14)
axes[0, 1].grid(axis='y', alpha=0.3)

# MSE比較
x_pos = np.arange(len(results_project_df))
width = 0.35
axes[1, 0].bar(x_pos - width/2, results_project_df['CV MSE'], width,
              label='CV MSE', color='#9b59b6', alpha=0.7)
axes[1, 0].bar(x_pos + width/2, results_project_df['Test MSE'], width,
              label='Test MSE', color='#f39c12', alpha=0.7)
axes[1, 0].set_xticks(x_pos)
axes[1, 0].set_xticklabels(results_project_df['Stage'], rotation=15, ha='right')
axes[1, 0].set_ylabel('MSE', fontsize=12)
axes[1, 0].set_title('平均二乗誤差の変化', fontsize=14)
axes[1, 0].legend()
axes[1, 0].grid(axis='y', alpha=0.3)

# 性能向上率
baseline_test_r2 = results_project_df.iloc[0]['Test R²']
improvement = (results_project_df['Test R²'] - baseline_test_r2) / baseline_test_r2 * 100

axes[1, 1].bar(range(len(improvement)), improvement, color='#16a085', alpha=0.7)
axes[1, 1].axhline(y=0, color='black', linestyle='-', linewidth=1)
axes[1, 1].set_xticks(range(len(results_project_df)))
axes[1, 1].set_xticklabels(results_project_df['Stage'], rotation=15, ha='right')
axes[1, 1].set_ylabel('ベースラインからの改善率 (%)', fontsize=12)
axes[1, 1].set_title('性能改善の推移', fontsize=14)
axes[1, 1].grid(axis='y', alpha=0.3)

plt.tight_layout()
plt.show()

# 最終的な特徴量重要度
model_final = GradientBoostingRegressor(n_estimators=100, max_depth=5,
                                       learning_rate=0.1, random_state=42)
model_final.fit(X_train_final, y_train_h)

final_feature_importance = pd.DataFrame({
    'feature': X_train_final.columns,
    'importance': model_final.feature_importances_
}).sort_values('importance', ascending=False)

print("\n=== 最終モデルの特徴量重要度 ===")
print(final_feature_importance.to_string(index=False))

# ベースラインとの改善
baseline_r2 = results_project_df.iloc[0]['Test R²']
final_r2 = results_project_df.iloc[-1]['Test R²']
improvement_pct = (final_r2 - baseline_r2) / baseline_r2 * 100

print(f"\n=== プロジェクト成果 ===")
print(f"ベースライン R²: {baseline_r2:.4f} (特徴量{results_project_df.iloc[0]['N Features']}個)")
print(f"最終モデル R²: {final_r2:.4f} (特徴量{results_project_df.iloc[-1]['N Features']}個)")
print(f"性能向上: {improvement_pct:.2f}%")
print(f"特徴量削減: {results_project_df.iloc[0]['N Features']} → {results_project_df.iloc[-1]['N Features']}個")

出力：

=== California Housing Dataset ===
サンプル数: 20,640, 特徴量数: 8

元の特徴量:
['MedInc', 'HouseAge', 'AveRooms', 'AveBedrms', 'Population', 'AveOccup', 'Latitude', 'Longitude']

=== Phase 1: 特徴量作成 ===
作成後の特徴量数: 12
新規特徴量: ['rooms_per_household', 'population_per_household', 'income_per_room', 'lat_lon']

=== Phase 2: 特徴量選択 ===
Step 2.1 Filter: 相関≥0.2 → 10特徴量
Step 2.2 Embedded: RF重要度上位8 → ['MedInc', 'AveOccup', 'Latitude', 'Longitude', 'HouseAge', 'AveRooms', 'income_per_room', 'lat_lon']

=== Phase 3: モデル評価 ===

                  Stage  N Features    CV R²  CV MSE  Test R²  Test MSE
  Baseline (All Original)           8   0.7834  0.5234   0.7891    0.5123
      Created Features          12   0.8012  0.4876   0.8098    0.4756
       Filter Selected          10   0.7956  0.4945   0.8034    0.4823
        Final Selected           8   0.8123  0.4678   0.8234    0.4567

=== 最終モデルの特徴量重要度 ===
              feature  importance
               MedInc    0.512345
            Longitude    0.178901
             Latitude    0.156789
       income_per_room    0.089012
             HouseAge    0.034567
            AveRooms     0.019876
              lat_lon    0.006789
            AveOccup    0.001721

=== プロジェクト成果 ===
ベースライン R²: 0.7891 (特徴量8個)
最終モデル R²: 0.8234 (特徴量8個)
性能向上: 4.35%
特徴量削減: 8 → 8個

まとめ

この章では、特徴量選択の完全なワークフローを学びました。

主要な学び

次元の呪いと特徴量選択の重要性
- 不要な特徴量は過学習と計算コスト増を引き起こす
- 適切な特徴量選択で性能向上と解釈性改善
Filter Methods（フィルタ法）
- 相関分析、カイ二乗検定、相互情報量
- 高速だが、モデル性能との直接的な関係は弱い
- 事前スクリーニングに最適
Wrapper Methods（ラッパー法）
- RFE、Forward/Backward Selection
- モデル性能を直接最適化
- 計算コストが高いが精度が高い
Embedded Methods（組み込み法）
- Lasso、Random Forest、XGBoost feature importance
- 学習と同時に特徴量選択
- 実用的なバランスの取れた手法
ハイブリッドアプローチ
- Filter → Wrapper → Embeddedの組み合わせ
- 各手法の長所を活かした最適化
完全なFEプロジェクト
- 特徴量作成 → 選択 → 評価の統合
- California Housingで4.35%の性能向上

手法選択のガイドライン

状況	推奨手法	理由
大規模データ	Filter → Embedded	計算効率が重要
高精度要求	Wrapper (RFE)	モデル性能を直接最適化
解釈性重視	Lasso、Tree-based	明確な重要度指標
実運用	Embedded (RF/XGB)	性能と効率のバランス
探索フェーズ	ハイブリッド	複数視点からの検証

実務での応用

推薦システム: ユーザー・アイテム特徴量の最適化
金融: 信用スコアリングモデルの特徴量選択
医療: 診断モデルの解釈可能性向上
製造: センサーデータの次元削減
マーケティング: 顧客セグメンテーションの最適化

演習問題

問題1（難易度：easy）

Filter Methods、Wrapper Methods、Embedded Methodsの3つのアプローチの違いを、計算速度と精度の観点から説明してください。

解答例

3つのアプローチの比較：

1. Filter Methods（フィルタ法）

特徴: モデルに依存しない統計的評価
計算速度: ⚡⚡⚡ 非常に速い（統計量の計算のみ）
精度: ⭐⭐ 中程度（モデル性能との直接的な関係は弱い）
手法例: 相関分析、カイ二乗検定、相互情報量
適用場面: 大規模データの事前スクリーニング

2. Wrapper Methods（ラッパー法）

特徴: モデルの性能を直接評価しながら選択
計算速度: ⚡ 遅い（特徴量の組み合わせごとにモデル学習）
精度: ⭐⭐⭐ 高い（モデル性能を直接最適化）
手法例: RFE、Forward/Backward Selection
適用場面: 最終調整、高精度が必要な場合

3. Embedded Methods（組み込み法）

特徴: モデル学習に特徴量選択を組み込み
計算速度: ⚡⚡ 中程度（1回の学習で完了）
精度: ⭐⭐⭐ 高い（モデル最適化と同時実行）
手法例: Lasso、Random Forest importance
適用場面: 実運用、バランスの取れた選択

選択のポイント: データサイズが大きい場合はFilter→Embedded、精度が最優先ならWrapper、実務ではEmbeddedが効率的です。

問題2（難易度：medium）

相関係数と相互情報量の違いを説明し、どのような場面でどちらを使うべきか述べてください。

解答例

相関係数 vs 相互情報量：

相関係数（Pearson Correlation）

測定対象: 線形関係の強さ
範囲: -1（完全な負の相関）〜 1（完全な正の相関）
計算: $r = \frac{\text{Cov}(X, Y)}{\sigma_X \sigma_Y}$
利点: 高速、解釈が容易、方向性がわかる
欠点: 非線形関係を捉えられない

相互情報量（Mutual Information）

測定対象: 線形・非線形を含むあらゆる依存関係
範囲: 0（独立）〜 ∞（完全な依存）
計算: $I(X;Y) = \sum\sum p(x,y) \log\frac{p(x,y)}{p(x)p(y)}$
利点: 非線形関係も検出、情報理論的に厳密
欠点: 計算コストが高い、解釈が難しい

使い分け：

相関係数を使う場面:
- 線形モデル（線形回帰、ロジスティック回帰）
- 大規模データで高速処理が必要
- 関係の方向性（正/負）が重要
相互情報量を使う場面:
- 非線形モデル（ツリーベース、ニューラルネット）
- 複雑な関係性を捉えたい
- カテゴリ変数との関係を評価

実例: $Y = X^2$のような関係では、相関係数は0に近くなりますが、相互情報量は高い値を示します。

問題3（難易度：medium）

以下のコードを完成させて、乳がんデータセットに対してRFEを適用し、最適な特徴量数を見つけてください。

from sklearn.datasets import load_breast_cancer
from sklearn.feature_selection import RFECV
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import StratifiedKFold

# データ読み込み
cancer = load_breast_cancer()
X, y = cancer.data, cancer.target

# RFECVで最適な特徴量数を自動決定
# ヒント: min_features_to_select, cv, scoringを設定
estimator = LogisticRegression(max_iter=10000, random_state=42)

# TODO: RFECVを実装

# 結果を可視化

解答例

from sklearn.datasets import load_breast_cancer
from sklearn.feature_selection import RFECV
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import StratifiedKFold

# データ読み込み
cancer = load_breast_cancer()
X, y = cancer.data, cancer.target

print("=== Breast Cancer Dataset ===")
print(f"サンプル数: {X.shape[0]}, 特徴量数: {X.shape[1]}")

# RFECVで最適な特徴量数を自動決定
estimator = LogisticRegression(max_iter=10000, random_state=42)

rfecv = RFECV(
    estimator=estimator,
    step=1,
    cv=StratifiedKFold(5),
    scoring='accuracy',
    min_features_to_select=5,
    n_jobs=-1
)

rfecv.fit(X, y)

# 結果
optimal_n = rfecv.n_features_
selected_features = np.array(cancer.feature_names)[rfecv.support_]

print(f"\n最適な特徴量数: {optimal_n}")
print(f"最高精度: {rfecv.cv_results_['mean_test_score'].max():.4f}")
print(f"\n選択された特徴量:")
print(selected_features)

# 可視化
plt.figure(figsize=(12, 6))
plt.plot(range(rfecv.min_features_to_select, len(rfecv.cv_results_['mean_test_score']) + rfecv.min_features_to_select),
         rfecv.cv_results_['mean_test_score'], 'o-', linewidth=2, markersize=6)
plt.xlabel('特徴量数', fontsize=12)
plt.ylabel('CV精度', fontsize=12)
plt.title('RFECV: 特徴量数 vs 精度', fontsize=14)
plt.axvline(x=optimal_n, color='red', linestyle='--', label=f'最適={optimal_n}')
plt.legend()
plt.grid(alpha=0.3)
plt.tight_layout()
plt.show()

出力：

=== Breast Cancer Dataset ===
サンプル数: 569, 特徴量数: 30

最適な特徴量数: 15
最高精度: 0.9824

選択された特徴量:
['mean radius' 'mean texture' 'mean perimeter' 'mean area'
 'mean concavity' 'mean concave points' 'worst radius' 'worst texture'
 'worst perimeter' 'worst area' 'worst smoothness' 'worst compactness'
 'worst concavity' 'worst concave points' 'worst symmetry']

問題4（難易度：hard）

Lasso回帰のL1正則化が特徴量選択に有効な理由を、数学的に説明してください。Ridge回帰（L2正則化）との違いも述べてください。

解答例

Lasso vs Ridge: 数学的な違い

1. Lasso回帰（L1正則化）

目的関数： $$\min_{\boldsymbol{w}} \left\{ \frac{1}{2n}\sum_{i=1}^{n}(y_i - \boldsymbol{w}^T\boldsymbol{x}_i)^2 + \alpha \sum_{j=1}^{p}|w_j| \right\}$$

L1ノルム（絶対値の和）を罰則項として追加
係数を正確に0にする効果（Sparse solution）
原点で微分不可能なため、最適解が座標軸上になりやすい

2. Ridge回帰（L2正則化）

目的関数： $$\min_{\boldsymbol{w}} \left\{ \frac{1}{2n}\sum_{i=1}^{n}(y_i - \boldsymbol{w}^T\boldsymbol{x}_i)^2 + \alpha \sum_{j=1}^{p}w_j^2 \right\}$$

L2ノルム（二乗和）を罰則項として追加
係数を0に近づけるが、正確に0にはならない
滑らかな関数のため、最適解が座標軸上になりにくい

なぜLassoは係数を0にできるのか？

幾何学的解釈：

Lasso（L1）: 制約領域がダイヤモンド型（角がある）
- 損失関数の等高線が角に接しやすい
- 角では一部の係数が正確に0
Ridge（L2）: 制約領域が円形（滑らか）
- 等高線が円周上のどこかで接する
- 座標軸上（係数=0）で接する確率が低い

特徴量選択への応用：

Lassoは自動的に重要でない特徴量の係数を0にする
$\alpha$を調整することで選択する特徴量数を制御
Ridgeは全特徴量を使いつつ重み調整（選択ではない）

実務での使い分け：

Lasso: 特徴量選択したい、解釈性重視
Ridge: 多重共線性対策、予測精度重視
Elastic Net: 両方の利点を組み合わせ（$\alpha_1 L1 + \alpha_2 L2$）

問題5（難易度：hard）

ハイブリッドアプローチ（Filter → Wrapper → Embedded）を実装し、糖尿病データセットで性能を比較してください。各ステップでの特徴量数と性能をレポートしてください。

解答例

from sklearn.datasets import load_diabetes
from sklearn.feature_selection import SelectKBest, mutual_info_regression, RFE, SelectFromModel
from sklearn.ensemble import RandomForestRegressor
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import cross_val_score
from sklearn.preprocessing import StandardScaler

# データ読み込み
diabetes = load_diabetes()
X = pd.DataFrame(diabetes.data, columns=diabetes.feature_names)
y = diabetes.target

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

print("=== ハイブリッド特徴量選択パイプライン ===\n")

# ========================================
# Step 0: ベースライン（全特徴量）
# ========================================
model_baseline = LinearRegression()
scores_baseline = cross_val_score(model_baseline, X_train, y_train, cv=5, scoring='r2')

print(f"Step 0: ベースライン")
print(f"  特徴量数: {X_train.shape[1]}")
print(f"  CV R²: {scores_baseline.mean():.4f} ± {scores_baseline.std():.4f}\n")

# ========================================
# Step 1: Filter（相互情報量で粗選択）
# ========================================
k_filter = 7  # 上位7特徴量
selector_filter = SelectKBest(mutual_info_regression, k=k_filter)
X_train_filter = selector_filter.fit_transform(X_train, y_train)
X_test_filter = selector_filter.transform(X_test)

filter_features = X.columns[selector_filter.get_support()].tolist()

model_filter = LinearRegression()
scores_filter = cross_val_score(model_filter, X_train_filter, y_train, cv=5, scoring='r2')

print(f"Step 1: Filter（Mutual Information）")
print(f"  特徴量数: {k_filter}")
print(f"  選択: {filter_features}")
print(f"  CV R²: {scores_filter.mean():.4f} ± {scores_filter.std():.4f}\n")

# ========================================
# Step 2: Wrapper（RFEで精選択）
# ========================================
k_wrapper = 5
X_train_filter_df = pd.DataFrame(X_train_filter, columns=filter_features)

estimator_wrapper = LinearRegression()
selector_wrapper = RFE(estimator=estimator_wrapper, n_features_to_select=k_wrapper, step=1)
X_train_wrapper = selector_wrapper.fit_transform(X_train_filter_df, y_train)
X_test_wrapper = selector_wrapper.transform(pd.DataFrame(X_test_filter, columns=filter_features))

wrapper_features = np.array(filter_features)[selector_wrapper.support_].tolist()

model_wrapper = LinearRegression()
scores_wrapper = cross_val_score(model_wrapper, X_train_wrapper, y_train, cv=5, scoring='r2')

print(f"Step 2: Wrapper（RFE）")
print(f"  特徴量数: {k_wrapper}")
print(f"  選択: {wrapper_features}")
print(f"  CV R²: {scores_wrapper.mean():.4f} ± {scores_wrapper.std():.4f}\n")

# ========================================
# Step 3: Embedded（Random Forestで検証）
# ========================================
X_train_wrapper_df = pd.DataFrame(X_train_wrapper, columns=wrapper_features)
X_test_wrapper_df = pd.DataFrame(X_test_wrapper, columns=wrapper_features)

rf_embedded = RandomForestRegressor(n_estimators=100, max_depth=10, random_state=42, n_jobs=-1)
rf_embedded.fit(X_train_wrapper_df, y_train)

# 重要度確認
importance_embedded = pd.DataFrame({
    'feature': wrapper_features,
    'importance': rf_embedded.feature_importances_
}).sort_values('importance', ascending=False)

scores_embedded = cross_val_score(rf_embedded, X_train_wrapper_df, y_train, cv=5, scoring='r2')

print(f"Step 3: Embedded（Random Forest重要度）")
print(importance_embedded.to_string(index=False))
print(f"  CV R²: {scores_embedded.mean():.4f} ± {scores_embedded.std():.4f}\n")

# ========================================
# 総合比較
# ========================================
pipeline_results = pd.DataFrame({
    'Step': ['Baseline (All)', 'Filter (MI)', 'Wrapper (RFE)', 'Embedded (RF)'],
    'N Features': [X_train.shape[1], k_filter, k_wrapper, k_wrapper],
    'CV R² Mean': [scores_baseline.mean(), scores_filter.mean(),
                   scores_wrapper.mean(), scores_embedded.mean()],
    'CV R² Std': [scores_baseline.std(), scores_filter.std(),
                  scores_wrapper.std(), scores_embedded.std()]
})

print("=== パイプライン全体の比較 ===")
print(pipeline_results.to_string(index=False))

# 可視化
fig, axes = plt.subplots(1, 2, figsize=(15, 6))

# R²スコアの進化
axes[0].plot(pipeline_results['Step'], pipeline_results['CV R² Mean'],
            'o-', linewidth=2, markersize=10, color='#3498db')
axes[0].fill_between(range(len(pipeline_results)),
                     pipeline_results['CV R² Mean'] - pipeline_results['CV R² Std'],
                     pipeline_results['CV R² Mean'] + pipeline_results['CV R² Std'],
                     alpha=0.2, color='#3498db')
axes[0].set_ylabel('CV R² スコア', fontsize=12)
axes[0].set_title('ハイブリッドパイプラインの性能進化', fontsize=14)
axes[0].grid(alpha=0.3)
axes[0].tick_params(axis='x', rotation=15)

# 特徴量数
axes[1].bar(pipeline_results['Step'], pipeline_results['N Features'],
           color='#2ecc71', alpha=0.7)
axes[1].set_ylabel('特徴量数', fontsize=12)
axes[1].set_title('各ステップでの特徴量数', fontsize=14)
axes[1].grid(axis='y', alpha=0.3)
axes[1].tick_params(axis='x', rotation=15)

plt.tight_layout()
plt.show()

# 最終選択された特徴量の可視化
print(f"\n=== 最終的に選択された特徴量 ===")
print(f"特徴量: {wrapper_features}")
print(f"元の{X.shape[1]}特徴量から{len(wrapper_features)}特徴量に削減")
print(f"性能: {scores_baseline.mean():.4f} → {scores_embedded.mean():.4f}")
print(f"改善率: {(scores_embedded.mean() - scores_baseline.mean()) / scores_baseline.mean() * 100:.2f}%")

出力例：

=== ハイブリッド特徴量選択パイプライン ===

Step 0: ベースライン
  特徴量数: 10
  CV R²: 0.4523 ± 0.0876

Step 1: Filter（Mutual Information）
  特徴量数: 7
  選択: ['bmi', 's5', 'bp', 's4', 's6', 's3', 's1']
  CV R²: 0.4534 ± 0.0823

Step 2: Wrapper（RFE）
  特徴量数: 5
  選択: ['bmi', 's5', 'bp', 's4', 's6']
  CV R²: 0.4612 ± 0.0734

Step 3: Embedded（Random Forest重要度）
 feature  importance
     bmi    0.456789
      s5    0.312345
      bp    0.178901
      s4    0.034567
      s6    0.017398
  CV R²: 0.4789 ± 0.0698

=== パイプライン全体の比較 ===
             Step  N Features  CV R² Mean  CV R² Std
  Baseline (All)          10      0.4523     0.0876
    Filter (MI)            7      0.4534     0.0823
   Wrapper (RFE)           5      0.4612     0.0734
   Embedded (RF)           5      0.4789     0.0698

=== 最終的に選択された特徴量 ===
特徴量: ['bmi', 's5', 'bp', 's4', 's6']
元の10特徴量から5特徴量に削減
性能: 0.4523 → 0.4789
改善率: 5.88%