【Python実践】GBDT（Gradient Boosting Decision Tree）で回帰問題と分類問題を解く（コピペで使えるサンプルコード付き）

pip install matplotlib
pip install pandas
pip install scikit-learn

from sklearn.datasets import make_classification
import pandas as pd

# ダミーデータ作成
X, y = make_classification(
    n_samples=1000,
    n_features=7,
    n_informative=5,
    n_redundant=2,
    n_classes=4,
    n_clusters_per_class=1,
    flip_y=0.05,
    random_state=42
)

# データフレーム化
df = pd.DataFrame(X, columns=[f"特徴量_{i+1}" for i in range(7)])
df['正解データ'] = y

# ファイルに保存（任意）
df.to_csv("classification_data.csv", index=False)

import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.ensemble import GradientBoostingClassifier
from sklearn.metrics import accuracy_score, classification_report, confusion_matrix

def train_gbdt_classifier(df, target, columns, params={}, test_size=0.2, random_state=42):
    """
    指定されたデータフレームとパラメータに基づいてGDBTを学習させ、レポートを出力する関数。

    引数:
    df : pandas.DataFrame
        データを含むデータフレーム。
    target : str
        予測対象となるターゲット列名。
    columns : list of str
        使用する特徴量の列名のリスト。Targetが含まれていても無視する。
    params : dict, optional (default={})
        GDBTのハイパーパラメータ。
    test_size : float, optional (default=0.2)
        テストデータの割合。
    random_state : int, optional (default=42)
        乱数シード。

    戻り値:
    model : GradientBoostingClassifier
        訓練されたGDBTモデルのオブジェクト。
    cm : numpy.ndarray
        混同行列を表す2次元配列。
    """
    # 特徴量とターゲットに分ける（target 列を除外）
    columns = [col for col in columns if col != target]
    X = df[columns]
    y = df[target]

    # 訓練データとテストデータに分割する
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=random_state)

    # GBDT分類器のインスタンスを作成する
    model = GradientBoostingClassifier(random_state=random_state, **params)

    # モデルを訓練する
    model.fit(X_train, y_train)

    # テストデータで予測する
    y_pred = model.predict(X_test)

    # モデルの精度を評価する
    accuracy = accuracy_score(y_test, y_pred)
    report = classification_report(y_test, y_pred, target_names=[str(i) for i in range(len(df[target].unique()))])
    cm = confusion_matrix(y_test, y_pred)

    print(f'正解率: {accuracy}')
    print('混同行列:')
    print(cm)
    print('分類レポート:')
    print(report)

    return model, cm

if __name__ == "__main__":
    # データを読み込む
    df = pd.read_csv('classification_data.csv')
    # 関数の呼び出し（ハイパーパラメータはデフォルト値を使用）
    model, cm = train_gbdt_classifier(df, target='正解データ', columns=df.columns)
    # クラス名を設定して混同行列をプロット（プロット部分が必要なら別途追加）
    class_names = df['正解データ'].unique().astype(str)

# 関数の呼び出し（ハイパーパラメータを指定）
train_gbdt_classifier(df, target='正解データ', columns=df.columns,params={'n_estimators': 100})

#---------------------------------------------------------------
# この部分に、前述した train_gbdt_classifier関数 を張り付けて下さい
#---------------------------------------------------------------

import matplotlib.pyplot as plt
import seaborn as sns
from matplotlib import rcParams
rcParams['font.family'] = 'Meiryo'

def plot_confusion_matrix(cm, class_names=[]):
    """
    混同行列をプロットする関数。

    引数:
    cm : numpy.ndarray
        混同行列を表す2次元配列。
    class_names : list of str, optional
        クラス名のリスト。省略すると連番が表示される
    """
    plt.figure(figsize=(10, 7))
    sns.heatmap(cm, annot=True, fmt='d', cmap='Blues', xticklabels=class_names, yticklabels=class_names)
    plt.xlabel('予測ラベル')
    plt.ylabel('実際のラベル')
    plt.title('混同行列')
    plt.show()

if __name__ == "__main__":
    # データを読み込む
    df = pd.read_csv('classification_data.csv')

    # 関数の呼び出し（ハイパーパラメータはデフォルト値を使用）
    model, cm = train_gbdt_classifier(df, target='正解データ', columns=df.columns)

    # クラス名を設定して混同行列をプロット
    class_names = df['正解データ'].unique().astype(str)
    plot_confusion_matrix(cm, class_names=class_names)

import numpy as np
import pandas as pd

# 特徴量の数
n_features = 7

# 各特徴量の標準偏差
std_devs = [1, 1, 1, 1, 1, 1, 1]

# 相関行列の設定（例: 1番目と2番目の特徴量が強い相関を持つ）
correlation_matrix = np.array([
    [1, 0.8, 0.2, 0.1, 0.1, 0.1, 0.1],
    [0.8, 1, 0.3, 0.2, 0.2, 0.2, 0.2],
    [0.2, 0.3, 1, 0.4, 0.3, 0.3, 0.3],
    [0.1, 0.2, 0.4, 1, 0.5, 0.5, 0.5],
    [0.1, 0.2, 0.3, 0.5, 1, 0.4, 0.4],
    [0.1, 0.2, 0.3, 0.5, 0.4, 1, 0.3],
    [0.1, 0.2, 0.3, 0.5, 0.4, 0.3, 1]
])

# 共分散行列の計算
covariance_matrix = np.outer(std_devs, std_devs) * correlation_matrix

# ダミーデータの生成
np.random.seed(42)
n_samples = 1000
data = np.random.multivariate_normal(np.zeros(n_features), covariance_matrix, size=n_samples)

# DataFrameに変換
df = pd.DataFrame(data, columns=[f'特徴量_{i+1}' for i in range(n_features)])

# 結果の確認
print(df.head())

# CSVファイルに保存
df.to_csv('regression_data.csv', index=False)

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.metrics import mean_squared_error, mean_absolute_error, r2_score, mean_absolute_percentage_error

def train_gbdt_regressor(df, target, columns, params={}, test_size=0.2, random_state=42):
    """
    指定されたデータフレームとパラメータに基づいてGBDTモデルを学習させ、レポートを出力する関数。

    引数:
    df : pandas.DataFrame
        データを含むデータフレーム。
    target : str
        予測対象となるターゲット列名。
    columns : list of str
        使用する特徴量の列名のリスト。ターゲットが含まれていても無視する。
    params : dict, optional (default={})
        GBDTのハイパーパラメータ。
    test_size : float, optional (default=0.2)
        テストデータの割合。
    random_state : int, optional (default=42)
        乱数シード。

    戻り値:
    model : GradientBoostingRegressor
        訓練されたGBDTモデルのオブジェクト。
    y_test : list
        テストに使った正解データ。
    y_pred: list
        予測した結果。
    """
    # 特徴量とターゲットに分ける（target 列を除外）
    columns = [col for col in columns if col != target]
    X = df[columns]
    y = df[target]

    # 訓練データとテストデータに分割する
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=random_state)

    # GBDTモデルのインスタンスを作成する
    model = GradientBoostingRegressor(random_state=random_state, **params)

    # モデルを訓練する
    model.fit(X_train, y_train)

    # テストデータで予測する
    y_pred = model.predict(X_test)

    # モデルの精度を評価する
    mse = mean_squared_error(y_test, y_pred)
    rmse = np.sqrt(mse)  # RMSE
    mae = mean_absolute_error(y_test, y_pred)
    r2 = r2_score(y_test, y_pred)
    mape = mean_absolute_percentage_error(y_test, y_pred)

    # 統計情報を表示
    print(f'{target}の統計情報:')
    print(f'最小値(Min): {y.min()}')
    print(f'最大値(Max): {y.max()}')
    print(f'平均値(Mean): {y.mean()}')
    print(f'標準偏差(Std): {y.std()}')
    print(f'分散(Var): {y.var()}')

    # 精度指標を出力する
    print(f'平均絶対誤差 (Mean Absolute Error): {mae}')
    print(f'平均二乗誤差 (Mean Squared Error): {mse}')
    print(f'ルート平均二乗誤差 (Root Mean Squared Error): {rmse}')
    print(f'平均絶対パーセンテージ誤差 (MAPE): {mape}')
    print(f'決定係数 (R^2 Score): {r2}')

    # 特徴量の重要度を表示する
    feature_importances = model.feature_importances_
    importance_df = pd.DataFrame({'特徴量 (Feature)': columns, '重要度 (Importance)': feature_importances})
    importance_df = importance_df.sort_values(by='重要度 (Importance)', ascending=False)
    print("\n特徴量の重要度:")
    print(importance_df)

    return model, y_test, y_pred

if __name__ == "__main__":
    # データの読み込み
    df = pd.read_csv('regression_data.csv')
    
    # 特徴量_1 の予測モデルを作成する
    model, y_test, y_pred = train_gbdt_regressor(df, '特徴量_1', df.columns)

#------------------------------------------------------------------------
# この部分に、前述した train_gbdt_regressor関数 を張り付けて下さい
#------------------------------------------------------------------------

from scipy.stats import norm
from sklearn.metrics import r2_score
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import rcParams
rcParams['font.family'] = 'Meiryo'


def plot_figures(y_test, y_pred):
    """
    グラフを描画し、指定されたファイル名で保存する関数

    引数:
    y_test : pandas.Series
        実際の値
    y_pred : numpy.ndarray
        予測値
    """
    # 残差の計算
    residuals = y_test - y_pred

    fig, axs = plt.subplots(1, 4, figsize=(18, 5))

    # 実際の値と予測値の折れ線グラフ
    axs[0].plot(y_test.values, label='実際の値 (Actual)', linestyle='-', marker=None)
    axs[0].plot(y_pred, label='予測値 (Predicted)', color='red', linestyle='-', marker=None)
    axs[0].set_xlabel('サンプルインデックス (Sample Index)')
    axs[0].set_ylabel('値 (Value)')
    axs[0].set_title('実際の値 VS 予測値\n(Comparison of Actual and Predicted Values)')
    axs[0].legend()

    # 実際の値と予測値の散布図
    axs[1].scatter(y_test, y_pred, alpha=0.3, label='予測値 (Predicted)')
    axs[1].plot([y_test.min(), y_test.max()], [y_test.min(), y_test.max()], 'r--', label='対角線 (y=x)')

    # 回帰線の計算
    slope, intercept = np.polyfit(y_test, y_pred, 1)
    r2 = r2_score(y_test, y_pred)
    axs[1].plot(y_test, slope * y_test + intercept, color='blue', linestyle='-', linewidth=2, label=f'回帰線 (y={slope:.2f}x + {intercept:.2f}, $R^2$={r2:.2f})')
    
    axs[1].set_xlabel('実際の値 (Actual)')
    axs[1].set_ylabel('予測値 (Predicted)')
    axs[1].set_title('実際の値 vs 予測値\n(Actual vs Predicted)')
    axs[1].legend()

    # 残差プロット
    axs[2].scatter(y_pred, residuals, alpha=0.3, label='残差 (Residuals)')
    axs[2].hlines(y=0, xmin=y_pred.min(), xmax=y_pred.max(), linestyles='dashed', label='y=0')
    
    # 残差の平均と標準偏差を計算
    residuals_mean = np.mean(residuals)
    residuals_std = np.std(residuals)
    
    # y_pred をソートし、それに対応する範囲を再計算
    sorted_indices = np.argsort(y_pred)
    sorted_y_pred = y_pred[sorted_indices]
    sorted_upper = residuals_mean + residuals_std
    sorted_lower = residuals_mean - residuals_std

    axs[2].axhline(y=residuals_mean, color='blue', linestyle='--', label=f'平均 (Mean): {residuals_mean:.2f}')
    # 標準偏差の範囲を水色で塗りつぶし
    axs[2].fill_between(
        sorted_y_pred,  # ソートされた予測値
        sorted_lower,   # 下限 (固定値)
        sorted_upper,   # 上限 (固定値)
        color='lightblue',
        alpha=0.4,
        label=f'標準偏差 (Std): {residuals_std:.2f}'
    )

    axs[2].set_xlabel('予測値 (Predicted)')
    axs[2].set_ylabel('残差 (Residuals)')
    axs[2].set_title('残差プロット\n(Residual Plot)')
    axs[2].legend()

    # 残差のヒストグラム
    axs[3].hist(residuals, bins=30, edgecolor='k', density=True, alpha=0.6, color='g')

    # 残差のヒストグラムに正規分布曲線を追加
    xmin, xmax = axs[3].get_xlim()
    x = np.linspace(xmin, xmax, 100)
    p = norm.pdf(x, residuals_mean, residuals_std)
    axs[3].plot(x, p, 'k', linewidth=2)
    
    axs[3].set_xlabel('残差 (Residuals)')
    axs[3].set_ylabel('頻度 (Frequency)')
    axs[3].set_title('残差のヒストグラム\n(Histogram of Residuals)')

    plt.tight_layout()
    plt.show()


if __name__ == "__main__":
    # データを読み込む
    df = pd.read_csv('regression_data.csv')

    # 関数の呼び出し（ハイパーパラメータはデフォルト値を使用）
    model, y_test, y_pred = train_gbdt_regressor(df, target='特徴量_1', columns=df.columns)

    # グラフを描画
    plot_figures(y_test, y_pred)

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split, GridSearchCV
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.metrics import mean_squared_error, mean_absolute_error, r2_score, mean_absolute_percentage_error

def grid_search_gbdt(df, target, columns, param_grid, test_size=0.2, random_state=42):
    """
    グリッドサーチを使用してGBDTモデルのハイパーパラメータをチューニングし、最適なモデルを見つける関数。

    引数:
    df : pandas.DataFrame
        データを含むデータフレーム。
    target : str
        予測対象となるターゲット列名。
    columns : list of str
        使用する特徴量の列名のリスト。ターゲットが含まれていても無視する。
    param_grid : dict
        グリッドサーチで試すハイパーパラメータの辞書。
    test_size : float, optional (default=0.2)
        テストデータの割合。
    random_state : int, optional (default=42)
        乱数シード。

    戻り値:
    best_model : GradientBoostingRegressor
        最適なハイパーパラメータで訓練されたGBDTモデルのオブジェクト。
    best_params : dict
        最適なハイパーパラメータの辞書。
    best_score : float
        グリッドサーチで得られた最適なスコア。
    """
    # 特徴量とターゲットに分ける（target 列を除外）
    columns = [col for col in columns if col != target]
    X = df[columns]
    y = df[target]

    # 訓練データとテストデータに分割する
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=test_size, random_state=random_state)

    # GBDTモデルのインスタンスを作成する
    model = GradientBoostingRegressor(random_state=random_state)

    # グリッドサーチを行う
    grid_search = GridSearchCV(estimator=model, param_grid=param_grid, cv=5, scoring='neg_mean_squared_error', n_jobs=-1)
    grid_search.fit(X_train, y_train)

    # 最適なモデルとパラメータを取得する
    best_model = grid_search.best_estimator_
    best_params = grid_search.best_params_
    best_score = grid_search.best_score_

    # テストデータで予測する
    y_pred = best_model.predict(X_test)

    # モデルの精度を評価する
    mse = mean_squared_error(y_test, y_pred)
    rmse = np.sqrt(mse)  # RMSE
    mae = mean_absolute_error(y_test, y_pred)
    r2 = r2_score(y_test, y_pred)
    mape = mean_absolute_percentage_error(y_test, y_pred)

    # 精度指標を出力する
    print(f'最適なハイパーパラメータ: {best_params}')
    print(f'グリッドサーチの最適なスコア: {best_score}')
    print(f'平均絶対誤差 (Mean Absolute Error): {mae}')
    print(f'平均二乗誤差 (Mean Squared Error): {mse}')
    print(f'ルート平均二乗誤差 (Root Mean Squared Error): {rmse}')
    print(f'平均絶対パーセンテージ誤差 (MAPE): {mape}')
    print(f'決定係数 (R^2 Score): {r2}')

    return best_model, best_params, best_score, y_test, y_pred

if __name__ == "__main__":
    # 回帰問題用テストデータの読み込み
    df = pd.read_csv('regression_data.csv')
    param_grid = {
        'n_estimators': [100, 200],
        'learning_rate': [0.01, 0.1],
        'max_depth': [3, 5, 10],
        'min_samples_split': [2, 5],
        'min_samples_leaf': [1, 2],
        'subsample': [0.8, 1.0]
    }
    best_model, best_params, best_score, y_test, y_pred = grid_search_gbdt(df, target='特徴量_1', columns=df.columns, param_grid=param_grid)