Qiita Teams that are logged in
You are not logged in to any team

Log in to Qiita Team
Community
OrganizationEventAdvent CalendarQiitadon (β)
Service
Qiita JobsQiita ZineQiita Blog
79
Help us understand the problem. What are the problem?

LightGBMのパラメータチューニングまとめ

はじめに

本記事は、下記のハイパーパラメータチューニングに関する記事の、LightGBMにおける実装例を紹介する記事となります。

LightGBMとパラメータチューニング

LightGBMは分類や回帰に用いられる機械学習アルゴリズムで、その性能の高さや使い勝手の良さ(特徴量重要度などが出せる)から、特に回帰においてはXGBoostと並ぶメジャーなアルゴリズムです。

XGBoostと比較したとき、LightGBMは動作が軽量で大規模データの学習に強いと言われています。

一方でLightGBMは多くのハイパーパラメータを持つため、その性能を十分に発揮するためにはパラメータチューニングが重要となります。

チューニング対象のパラメータ

LightGBMの主なパラメータは、こちらの記事で分かりやすく解説されています。
全てのパラメータは、こちらの公式ドキュメントで解説されています

LightGBMは他の回帰アルゴリズム(例:ラッソ回帰(1種類)SVR(3種類))と比べてパラメータの数が多く、また使用するAPI(Scikit-learn API or Training API)によってパラメータの数が変わるなど、複雑なパラメータ構成を持っています。

そのため、全てのパラメータをチューニング対象とすると制御が難しく、一部の重要なパラメータのみを選んでチューニング対象とすることが一般的です。

本記事で使用するパラメータ

本記事では、Optuna公式サンプルコードのチューニング対象パラメータを使用します。

またチューニングに先立ち、学習に使用するAPIおよびboosterのタイプ(boosting type)を選択する必要があります。今回は

選択内容 選択する方式 選択理由
API Scikit-Learn API Scikit-Learnの各種チューニング用メソッドを利用するため
boosting type gbdt 最もよく使われる方式(他方式との差はこちら参照

を採用します

この際、Scikit-learn APIとLearning APIで一部パラメータの名称が異なるのでご注意ください
今回使用するパラメータ一覧と、対応する名称を下記します

Scikit-Learn
APIでの名称
Training
APIでの
名称
型 (スケール) 範囲 デフォルト値 内容
reg_alpha lambda_l1 float (log) 0-∞ 0 L1正則化項の係数(小さいほど過学習寄り)     
reg_lambda lambda_l2 float (log) 0-∞ 0 L2正則化項の係数(小さいほど過学習寄り)
num_leaves num_leaves int (linear) 2-131072 31 1本の木の最大葉枚数(大きいほど過学習寄り)
colsample _bytree feature_fraction float (linear) 0-1 1.0 各決定木においてランダムに抽出される列の割合(大きいほど過学習寄り)
subsample bagging_fraction float (linear) 0-1 1.0 各決定木においてランダムに抽出される標本の割合(大きいほど過学習寄り)
subsample_freq bagging_freq int (linear) 0-∞ 0 ここで指定したイテレーション毎にバギング実施(大きいほど過学習寄りだが、0のときバギング非実施となるので注意)
min_child_samples min_data_in_leaf int (linear) 0-∞ 20 1枚の葉に含まれる最小データ数(小さいほど過学習寄り)

実装

下図のフロー(こちらの記事と同じ)に基づき、LightGBM回帰におけるチューニングを実装します
image.png
コードはこちらのGitHub(lgbm_tuning_tutorials.py)にもアップロードしております。

また、希望があればLightGBM分類の記事も作成しますので、コメント欄に記載いただければと思います。

実装に必要な環境整備

ライブラリ

チューニングの実装に必要なライブラリ等をインストールします。
私の環境では、以下のライブラリを使用しています

・Python本体 (動作確認時は3.9.5を使用)
・Matplotlibで表示可能な環境 (Jupyter等)
・下記ライブラリ (参考として動作確認時のライブラリバージョンも記載)
 lightgbm (3.2.1)
 scikit-learn (0.24.2)
 bayesian-optimization (1.2.0)
 optuna (2.7.0)
 numpy (1.20.3)
 pandas (1.2.4)
 seaborn (0.11.0)
 matplotlib (3.4.2)
こちらの回帰可視化ライブラリ(チューニング前後のモデル可視化にのみ使用するので、なくともチューニング自体はできます)

使用するサンプルデータ

こちらの記事の、大阪都構想選挙データを使用します

ハードウェア

チューニング速度はCPU等の性能に大きく依存します
参考までに本記事で使用したPC(Lenovoのモバイル)のスペックを下記します

デバイス スペック
CPU Intel Core i5-8250U@1.6GHz
メモリ 8GB(4+4) DDR4 2400MHz SODIMM
ストレージ SSD PCIe M.2
OS Windows 10 Pro 64bit

コードの実行にはVSCodeのJupyter拡張を使用しております。

CPUやメモリの使用状況にもチューニングは左右されますが、余計なプロセスを終了してからチューニングを開始することで、繰り返してもほぼ同じ所要時間となり、再現性のとれる結果となりました。

実際の実装例

前述のフローに基づき実装を進めていきます

データの読込

こちらのデータセットをCSVからPandasのDataFrame(Jupyterでの表示用)およびNumpyのndarray(LightGBMモデルへの読込用)に読み込みます

データの読込
import pandas as pd
import time
df_osaka = pd.read_csv(f'./sample_data/osaka_metropolis_english.csv')
OBJECTIVE_VARIALBLE = 'approval_rate'  # 目的変数
USE_EXPLANATORY = ['2_between_30to60', '3_male_ratio', '5_household_member', 'latitude']  # 説明変数
y = df_osaka[OBJECTIVE_VARIALBLE].values  # 目的変数をndarray化
X = df_osaka[USE_EXPLANATORY].values  # 説明変数をndarray化
# データを表示
df_osaka[USE_EXPLANATORY + [OBJECTIVE_VARIALBLE]]
    2_between_30to60    3_male_ratio    5_household_member  latitude    approval_rate
0   0.499141    0.485807    1.83    42.325  0.562728
1   0.425080    0.481006    2.10    42.078  0.526104
2   0.467986    0.473394    2.06    41.535  0.536686
3   0.411267    0.473798    2.20    40.197  0.496048
4   0.390663    0.473722    2.16    43.273  0.471734
5   0.419868    0.474112    2.25    42.190  0.504592
6   0.422815    0.475847    2.52    42.257  0.511019
7   0.342897    0.577474    1.75    38.090  0.492851
8   0.506583    0.468446    1.65    40.870  0.509798
9   0.496780    0.467432    1.96    40.575  0.538094
10  0.374845    0.486567    2.34    39.018  0.465384
11  0.462086    0.505250    1.61    39.560  0.516693
12  0.383394    0.478802    2.36    36.575  0.486810
13  0.392461    0.465524    2.15    36.222  0.450049
14  0.390136    0.471774    2.38    37.270  0.456477
15  0.450164    0.464464    2.06    39.467  0.478386
16  0.368153    0.480087    2.27    39.218  0.480658
17  0.409648    0.460738    2.22    38.316  0.450814
18  0.386937    0.469701    2.23    37.322  0.473037
19  0.402482    0.488469    2.29    40.980  0.468385
20  0.400149    0.486460    2.18    39.835  0.429843
21  0.406196    0.490601    2.35    42.683  0.473805
22  0.433488    0.496974    1.97    43.260  0.550678
23  0.409827    0.491403    1.99    44.470  0.502997

チューニング前のモデルの確認

こちらのライブラリで、チューニング前モデルの特徴量と目的変数の関係を可視化します

チューニング前のモデル
from seaborn_analyzer import regplot
from lightgbm import LGBMRegressor
from sklearn.model_selection import KFold
# 乱数シード
seed = 42
# モデル作成
model = LGBMRegressor(boosting_type='gbdt', objective='regression',
                      random_state=seed, n_estimators=10000)  # チューニング前のモデル
# 学習時fitパラメータ指定
fit_params = {'verbose': 0,  # 学習中のコマンドライン出力
              'early_stopping_rounds': 10,  # 学習時、評価指標がこの回数連続で改善しなくなった時点でストップ
              'eval_metric': 'rmse',  # early_stopping_roundsの評価指標
              'eval_set': [(X, y)]  # early_stopping_roundsの評価指標算出用データ
              }
# クロスバリデーションして決定境界を可視化
cv = KFold(n_splits=3, shuffle=True, random_state=seed)  # KFoldでクロスバリデーション分割指定
regplot.regression_heat_plot(model, USE_EXPLANATORY, OBJECTIVE_VARIALBLE, df_osaka,
                             pair_sigmarange = 0.5, rounddigit_x1=3, rounddigit_x2=3,
                             cv=cv, display_cv_indices=0,
                             fit_params=fit_params)

lgbm_osaka.png

全ての条件で同じ予測値(色の濃さ)となっており、うまく回帰できていないことが分かります。
チューニングでの改善が必要だと一目で分かります。

ちなみに、LGBMRegressorはScikit-Learn APIにおけるLightGBM回帰を実行するクラスで、objectiveが学習時に使用する評価指標、random_stateが使用する乱数シードです。

・early_stopping_roundsについて

LightGBMにはearly_stopping_roundsという便利な機能があります。
XGBoostやLightGBMは学習を繰り返すことで性能を上げていくアルゴリズムですが、学習回数を増やしすぎると性能向上が止まって横ばいとなり、無意味な学習を繰り返して学習時間増加の原因となってしまいます(参考

early_stopping_roundsは、この学習回数を適切なタイミングで打ち止めるための仕組みで、「○回連続で評価指標が向上しなくなったら、学習を打ち止める」というアルゴリズムとなっています。
後述の学習打ち切り問題もあるので、基本的にはearly_stopping_roundsを指定した方が良いです

Scikit-Learn APIのLightGBMでearly_stopping_roundsを利用する場合、fit_params引数にdict形式で'early_stopping_rounds'、'eval_metric'および'eval_set'を指定します。
また、連続条件に至る前に学習が打ち切られないよう、n_estimatorsに大きな値(例:10000)を指定する必要もあります。

fit_paramsの引数名 意味
early_stopping_rounds 上記の「○回連続」を指定する引数
eval_metric 「評価指標」
eval_set 上記「評価指標」算出のためのデータセットを指定する引数
(verbose) 学習中のコマンドライン出力(0から2に大きくなるにつれ経過表示が細かくなる)

eval_setは本来であれば検証用データを入れる事が望ましいですが、cross_val_scoreメソッドの外側で検証用データを分けることができないので、本記事ではCV分割前のデータをそのまま入力します。
(厳密性にこだわる方は、cross_val_scoreやGridSearchCVクラスを使わず、スクラッチでクロスバリデーションを実装してください)
(参考記事でもtest_loglossが収束するタイミングとeval_loglossが収束するタイミングに大きな差はなく、私の経験上からもどちらのデータをeval_setに入れてもスコアに大きな差はないので、実用上は問題がないと考えています)

またeval_metricおよび前述のobjectiveは、後述の「手順1」で指定する評価指標に合わせることが望ましいです。

評価指標の名称 手順1の評価指標 objective eval_metric
RMSE neg_root_mean_squared_errorまたは'neg_mean_squared_error' 'regression' 'rmse'
MAE 'neg_mean_absolute_error' 'regression_l1' 'mae'

ざっと見た限りでは、RMSLEに相当するobjectiveやeval_metricは存在しなさそうです

参考までに、分類での評価指標も列記します
APIドキュメントを見る限りobjectiveはLogLossしか存在しないので、全ての指標を揃えようとすると実質LogLoss一択となりそうです。手順1の評価指標eval_metricはAUCなども選択できます。詳しくはリンクを参照ください)

評価指標の名称 手順1の評価指標 objective eval_metric
LogLoss (2値分類) 'neg_log_loss' 'binary' 'binary_logloss'
LogLoss (多クラス分類) 'neg_log_loss' 'multiclass' 'multi_logloss'

手順1 最大化したい評価指標の定義

今回は
・回帰
・外れ値が多いデータではない
ので、こちらのフローに基づきRMSE(Scikit-Learnでは'neg_mean_squared_error'または'neg_root_mean_squared_error')を採用します。

チューニング前の評価指標をクロスバリデーションで算出します。

チューニング前の評価指標算出
# %% 手順1) チューニング前の評価指標算出
from sklearn.model_selection import cross_val_score
import numpy as np
X = df_osaka[USE_EXPLANATORY].values  
y = df_osaka[OBJECTIVE_VARIALBLE]  # 目的変数をndarray化
scoring = 'neg_mean_squared_error'  # 評価指標をRMSEに指定
# クロスバリデーションで評価指標算出
scores = cross_val_score(model, X, y, cv=cv,
                         scoring=scoring, n_jobs=-1, fit_params=fit_params)
print(f'scores={scores}')
print(f'average_score={np.mean(scores)}')
scores=[-0.00164951 -0.00086281 -0.00165816]
average_score=-0.0013901602481669621

チューニング前のスコア(neg_mean_squared_error = -MSE)は-0.001390であることが分かりました。
前述の可視化の時点で予見はできましたが、XGBoostのときと比べるとかなり悪い値となっています。

チューニングでこの値を改善していきます

手順2 パラメータ種類と範囲の選択

・種類

Optunaサンプルコードに基づき、前述の7パラメータを使用します

・範囲

適当な初期範囲を定めて検証曲線をプロットし、その結果に基づいて範囲を調整していきます。

適当な初期範囲で検証曲線プロット
from sklearn.model_selection import validation_curve
import matplotlib.pyplot as plt
cv_params = {'reg_alpha': [0, 0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1, 3, 10],
             'reg_lambda': [0, 0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1, 3, 10],
             'num_leaves': [2, 4, 8, 16, 32, 64, 96, 128, 160, 192, 224, 256],
             'colsample_bytree': [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
             'subsample': [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
             'subsample_freq': [0, 1, 2, 3, 4, 5, 6, 7],
             'min_child_samples': [0, 2, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100]
             }
param_scales = {'reg_alpha': 'log',
                'reg_lambda': 'log',
                'num_leaves': 'linear',
                'colsample_bytree': 'linear',
                'subsample': 'linear',
                'subsample_freq': 'linear',
                'min_child_samples': 'linear'
                }
# 検証曲線のプロット(パラメータ毎にプロット)
for i, (k, v) in enumerate(cv_params.items()):
    train_scores, valid_scores = validation_curve(estimator=model,
                                                  X=X, y=y,
                                                  param_name=k,
                                                  param_range=v,
                                                  fit_params=fit_params,
                                                  cv=cv, scoring=scoring,
                                                  n_jobs=-1)
    # 学習データに対するスコアの平均±標準偏差を算出
    train_mean = np.mean(train_scores, axis=1)
    train_std  = np.std(train_scores, axis=1)
    train_center = train_mean
    train_high = train_mean + train_std
    train_low = train_mean - train_std
    # テストデータに対するスコアの平均±標準偏差を算出
    valid_mean = np.mean(valid_scores, axis=1)
    valid_std  = np.std(valid_scores, axis=1)
    valid_center = valid_mean
    valid_high = valid_mean + valid_std
    valid_low = valid_mean - valid_std
    # training_scoresをプロット
    plt.plot(v, train_center, color='blue', marker='o', markersize=5, label='training score')
    plt.fill_between(v, train_high, train_low, alpha=0.15, color='blue')
    # validation_scoresをプロット
    plt.plot(v, valid_center, color='green', linestyle='--', marker='o', markersize=5, label='validation score')
    plt.fill_between(v, valid_high, valid_low, alpha=0.15, color='green')
    # スケールをparam_scalesに合わせて変更
    plt.xscale(param_scales[k])
    # 軸ラベルおよび凡例の指定
    plt.xlabel(k)  # パラメータ名を横軸ラベルに
    plt.ylabel(scoring)  # スコア名を縦軸ラベルに
    plt.legend(loc='lower right')  # 凡例
    # グラフを描画
    plt.show()

image.png

ほとんどのパラメータで変化が見られませんが、min_child_samplesのみ変化が見られます。
これはサンプル数23に対し、min_child_samplesのデフォルト値20が大きすぎることが原因です。
同様にnum_leavesのデフォルト値31も大きすぎるので、

min_child_samples=5
num_leaves=2

に変更し、再度検証曲線をプロットします。

検証曲線再プロット(min_child_samplesのデフォルト値を20→5に、num_leavesを31→2に変更)
from sklearn.model_selection import validation_curve
import matplotlib.pyplot as plt
model.set_params(min_child_samples=5, num_leaves=2)
# 検証曲線のプロット(パラメータ毎にプロット)
for i, (k, v) in enumerate(cv_params.items()):
    train_scores, valid_scores = validation_curve(estimator=model,
                                                  X=X, y=y,
                                                  param_name=k,
                                                  param_range=v,
                                                  fit_params=fit_params,
                                                  cv=cv, scoring=scoring,
                                                  n_jobs=-1)
    # 学習データに対するスコアの平均±標準偏差を算出
    train_mean = np.mean(train_scores, axis=1)
    train_std  = np.std(train_scores, axis=1)
    train_center = train_mean
    train_high = train_mean + train_std
    train_low = train_mean - train_std
    # テストデータに対するスコアの平均±標準偏差を算出
    valid_mean = np.mean(valid_scores, axis=1)
    valid_std  = np.std(valid_scores, axis=1)
    valid_center = valid_mean
    valid_high = valid_mean + valid_std
    valid_low = valid_mean - valid_std
    # training_scoresをプロット
    plt.plot(v, train_center, color='blue', marker='o', markersize=5, label='training score')
    plt.fill_between(v, train_high, train_low, alpha=0.15, color='blue')
    # validation_scoresをプロット
    plt.plot(v, valid_center, color='green', linestyle='--', marker='o', markersize=5, label='validation score')
    plt.fill_between(v, valid_high, valid_low, alpha=0.15, color='green')
    # スケールを'log'に(線形なパラメータは'linear'にするので注意)
    plt.xscale(param_scales[k])
    # 軸ラベルおよび凡例の指定
    plt.xlabel(k)  # パラメータ名を横軸ラベルに
    plt.ylabel(scoring)  # スコア名を縦軸ラベルに
    plt.legend(loc='lower right')  # 凡例
    # グラフを描画
    plt.show()

image.png
他のパラメータも変化が見えてきました。

SVMのときほど顕著な傾向は出ていませんが、過学習または未学習になり過ぎない範囲を基準に、以下のようにチューニング範囲を定めます
image.png

手順3&4 パラメータ選択&クロスバリデーション

SVMのときと同様に、下記A、B、C (BayesianOptimization)、C (Optuna)の4種類の実装法をそれぞれ解説します。

名称 概要 メリット デメリット ライブラリ
A グリッドサーチ       パラメータを格子状に総当たり シンプルで結果解釈性が高い 計算時間がかかる Scikit-Learn
B ランダムサーチ パラメータをランダムに決定 平均的にはグリッドサーチより速い 運任せの要素あり Scikit-Learn
C ベイズ最適化 前回結果に基づきパラメータ決定 ランダムサーチより速い ライブラリ操作がやや難 BayesianOptimization,
Optuna

また、SVMと比べLightGBMは処理が重い(XGBoostほどではありませんが)ため、チューニング時間が課題となります。
スコアと共にチューニングにかかった時間も測定し、各アルゴリズムで比較してみます。

処理時間は試行数にも依存するため、後述のearly_stopping_roundなしのとき(=パラメータによる学習の収束速度に左右されない)のチューニング時間が等しくなる、という基準で、各アルゴリズムの試行数を定めました

名称 試行数 early_stopping_roundなし時のチューニング所要時間
グリッドサーチ       1152
48秒
ランダムサーチ 1120 43秒
ベイズ最適化
(Bayesian Optimization)
100
(初期20+ベイズ80)
50秒
ベイズ最適化
(Optuna)
400 49秒

グリッドサーチの場合

GridSearchCVクラスで、グリッドサーチによる最適化を実行します。

下記のようにパラメータの組合せ(3x2x4x3x2x2x4=1152通り)を設定し、グリッドサーチを実行します

グリッドサーチ&クロスバリデーションの実装
from sklearn.model_selection import GridSearchCV
start = time.time()
# 最終的なパラメータ範囲(1152通り)
cv_params = {'reg_alpha': [0.0001, 0.003, 0.1],
             'reg_lambda': [0.0001, 0.1],
             'num_leaves': [2, 3, 4, 6],
             'colsample_bytree': [0.4, 0.7, 1.0],
             'subsample': [0.4, 1.0],
             'subsample_freq': [0, 7],
             'min_child_samples': [0, 2, 5, 10]
             }
# グリッドサーチのインスタンス作成
gridcv = GridSearchCV(model, cv_params, cv=cv,
                      scoring=scoring, n_jobs=-1)
# グリッドサーチ実行(学習実行)
gridcv.fit(X, y, **fit_params)
# 最適パラメータの表示と保持
best_params = gridcv.best_params_
best_score = gridcv.best_score_
print(f'最適パラメータ {best_params}\nスコア {best_score}')
print(f'所要時間{time.time() - start}秒')
結果
最適パラメータ {'colsample_bytree': 0.4, 'min_child_samples': 0, 'num_leaves': 4, 'reg_alpha': 0.0001, 'reg_lambda': 0.0001, 'subsample': 0.4, 'subsample_freq': 0}
スコア -0.00037051319661631437
所要時間63.93486523628235秒

reg_alpha=0.0001
reg_lambda=0.0001
num_leaves=4
colsample_bytree=0.4
subsample=0.4
subsample_freq=0
min_child_samples=0
のときに評価指標neg_mean_squared_errorが最大となり、その時の評価指標は-0.00370‥であることが分かります。

パラメータの値と学習速度について

XGBoostのときほどではありませんが、後述のearly_stopping_roundを設定しないときと比べてやや所要時間が増えています。自作ライブラリ(後日公開&記事化します)でパラメータと所要時間の関係を測定したところ、下図のように一部のパラメータの値に応じて学習時間が増える(収束速度が遅くなる)ことが分かりました。
image.png
image.png

early_stopping_roundを指定しなければこの速度低下が防げるように見えますが、この場合は収束していないのに学習を打ち切っているため、学習性能ををフルに発揮できていないという意味で、やはり避けるべき状況です。

とはいえXGBoostのときと比べるとパラメータの値による学習時間の変化は小さいので、LightGBMに関してはパラメータの値による速度低下にはそれほど神経質にならなくともよいかと思います。

ランダムサーチの場合

RondomizedSearchCVクラスで、ランダムサーチによる最適化を実行します

ランダムサーチ&クロスバリデーションの実装
from sklearn.model_selection import RandomizedSearchCV
start = time.time()
# パラメータの種類と密度をグリッドサーチのときより増やす
cv_params = {'reg_alpha': [0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1],
             'reg_lambda': [0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1],
             'num_leaves': [2, 3, 4, 5, 6],
             'colsample_bytree': [0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
             'subsample': [0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
             'subsample_freq': [0, 1, 2, 3, 4, 5, 6, 7],
             'min_child_samples': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
             }
# ランダムサーチのインスタンス作成
randcv = RandomizedSearchCV(model, cv_params, cv=cv,
                            scoring=scoring, random_state=seed,
                            n_iter=1120, n_jobs=-1)
# ランダムサーチ実行(学習実行)
randcv.fit(X, y, **fit_params)
# 最適パラメータの表示と保持
best_params = randcv.best_params_
best_score = randcv.best_score_
print(f'最適パラメータ {best_params}\nスコア {best_score}')
print(f'所要時間{time.time() - start}秒')
結果
最適パラメータ {'subsample_freq': 6, 'subsample': 0.6, 'reg_lambda': 0.0001, 'reg_alpha': 0.0003, 'num_leaves': 2, 'min_child_samples': 1, 'colsample_bytree': 0.9}
スコア -0.0003321068306584331
所要時間61.21291470527649秒

reg_alpha=0.0003
reg_lambda=0.0001
num_leaves=2
colsample_bytree=0.9
subsample=0.6
subsample_freq=6
min_child_samples=1
のときに評価指標neg_mean_squared_errorが最大となり、その時の評価指標は-0.00332‥であることが分かります。

ベイズ最適化(BayesianOptimization)の場合

BayesianOptimizationによるベイズ最適化の実装法を解説します。
3.1の教訓に基づき、今回は最初から対数スケールに対応した実装法を採用します

ベイズ最適化&クロスバリデーションの実装(Logスケールに変数変換)
from bayes_opt import BayesianOptimization
start = time.time()
# パラメータ範囲(Tupleで範囲選択)
bayes_params = {'reg_alpha': (0.0001, 0.1),
                'reg_lambda': (0.0001, 0.1),
                'num_leaves': (2, 6),
                'colsample_bytree': (0.4, 1.0),
                'subsample': (0.4, 1.0),
                'subsample_freq': (0, 7),
                'min_child_samples': (0, 10)
                }
# 対数スケールパラメータを対数化
param_scales = {'reg_alpha': 'log',
                'reg_lambda': 'log',
                'num_leaves': 'linear',
                'colsample_bytree': 'linear',
                'subsample': 'linear',
                'subsample_freq': 'linear',
                'min_child_samples': 'linear'
                }
bayes_params_log = {k: (np.log10(v[0]), np.log10(v[1])) if param_scales[k] == 'log' else v for k, v in bayes_params.items()}
# 整数型パラメータを指定
int_params = ['num_leaves', 'subsample_freq', 'min_child_samples']

# ベイズ最適化時の評価指標算出メソッド(引数が多いので**kwargsで一括読込)
def bayes_evaluate(**kwargs):
    params = kwargs
    # 対数スケールパラメータは10のべき乗をとる
    params = {k: np.power(10, v) if param_scales[k] == 'log' else v for k, v in params.items()}
    # 整数型パラメータを整数化
    params = {k: round(v) if k in int_params else v for k, v in params.items()}
    # モデルにパラメータ適用
    model.set_params(**params)
    # cross_val_scoreでクロスバリデーション
    scores = cross_val_score(model, X, y, cv=cv,
                             scoring=scoring, fit_params=fit_params, n_jobs=-1)
    val = scores.mean()
    return val

# ベイズ最適化を実行
bo = BayesianOptimization(bayes_evaluate, bayes_params_log, random_state=seed)
bo.maximize(init_points=20, n_iter=80, acq='ei')
# 最適パラメータとスコアを取得
best_params = bo.max['params']
best_score = bo.max['target']
# 対数スケールパラメータは10のべき乗をとる
best_params = {k: np.power(10, v) if param_scales[k] == 'log' else v for k, v in best_params.items()}
# 整数型パラメータを整数化
best_params = {k: round(v) if k in int_params else v for k, v in best_params.items()}
# 最適パラメータを表示
print(f'最適パラメータ {best_params}\nスコア {best_score}')
print(f'所要時間{time.time() - start}秒')

※bo.maximizeの引数init_points, n_iter, acqは、こちらを参考に適宜変更してください

結果
最適パラメータ {'colsample_bytree': 0.7595706008386077, 'min_child_samples': 2, 'num_leaves': 2, 'reg_alpha': 0.0011012364917553066, 'reg_lambda': 0.0028477019767619194, 'subsample': 0.4422343314262108, 'subsample_freq': 3}
スコア -0.000299637259324912
所要時間28.630770444869995秒

reg_alpha=0.00110
reg_lambda=0.0285
num_leaves=2
colsample_bytree=0.760
subsample=0.442
subsample_freq=3
min_child_samples=2
のときに評価指標neg_mean_squared_errorが最大となり、その時の評価指標は-0.00299‥であることが分かります。

評価指標、速度共にグリッドサーチ・ランダムサーチより良好な結果です

ベイズ最適化(Optuna)の場合

Optunaによるベイズ最適化の実装法を解説します。

Optunaベイズ最適化&クロスバリデーションの実装
import optuna
start = time.time()
# ベイズ最適化時の評価指標算出メソッド
def bayes_objective(trial):
    params = {
        'reg_alpha': trial.suggest_float('reg_alpha', 0.0001, 0.1, log=True),
        'reg_lambda': trial.suggest_float('reg_lambda', 0.0001, 0.1, log=True),
        'num_leaves': trial.suggest_int('num_leaves', 2, 6),
        'colsample_bytree': trial.suggest_float('colsample_bytree', 0.4, 1.0),
        'subsample': trial.suggest_float('subsample', 0.4, 1.0),
        'subsample_freq': trial.suggest_int('subsample_freq', 0, 7),
        'min_child_samples': trial.suggest_int('min_child_samples', 0, 10)
    }
    # モデルにパラメータ適用
    model.set_params(**params)
    # cross_val_scoreでクロスバリデーション
    scores = cross_val_score(model, X, y, cv=cv,
                             scoring=scoring, fit_params=fit_params, n_jobs=-1)
    val = scores.mean()
    return val

# ベイズ最適化を実行
study = optuna.create_study(direction='maximize',
                            sampler=optuna.samplers.TPESampler(seed=seed))
study.optimize(bayes_objective, n_trials=400)

# 最適パラメータの表示と保持
best_params = study.best_trial.params
best_score = study.best_trial.value
print(f'最適パラメータ {best_params}\nスコア {best_score}')
print(f'所要時間{time.time() - start}秒')
結果
最適パラメータ {'reg_alpha': 0.003766446430110217, 'reg_lambda': 0.022623094412872694, 'num_leaves': 3, 'colsample_bytree': 0.886490082071905, 'subsample': 0.50263070478602, 'subsample_freq': 4, 'min_child_samples': 1}
スコア -0.0002684894224768026
所要時間45.109469175338745秒

reg_alpha=0.00377
reg_lambda=0.0226
num_leaves=3
colsample_bytree=0.886
subsample=0.503
subsample_freq=4
min_child_samples=1
のときに評価指標neg_mean_squared_errorが最大となり、その時の評価指標は-0.00268‥であることが分かります。

以上、4種類のアルゴリズムでのチューニング結果をまとめると、下の表のようになります。

名称 試行数 スコア(大きいほどGood) 所要時間
グリッドサーチ       1124 -0.000371 64秒
ランダムサーチ 1120 -0.000332 61秒
ベイズ最適化
(Bayesian Optimization)
100 -0.000299 29秒
ベイズ最適化
(Optuna)
400 -0.000268 45秒
※ベイズ最適化
(Optuna)
250 -0.000271 28秒

※Bayesian Optimizationの所要時間を下回るよう試行数を減らして再実行

Optunaが速度、評価指標共にトップであり、効率良い探索ができていることが分かります。

※参考:early_stopping_round不使用時のチューニング結果

上記のパターンではグリッドサーチが多くの時間が掛かっていますが、early_stopping_roundを指定しない場合(常にデフォルトの学習回数100回を使用)では、他のアルゴリズムとほぼ同時間で処理が完了できます(前述の収束前に学習が終了する問題はありますが‥)

参考までに、early_stopping_roundを指定しない場合のスコアと所要時間を下記します

グリッドサーチ
最適パラメータ {'colsample_bytree': 0.4, 'min_child_samples': 0, 'num_leaves': 4, 'reg_alpha': 0.0001, 'reg_lambda': 0.0001, 'subsample': 0.4, 'subsample_freq': 0}
スコア -0.0004023054268113945
所要時間48.24397921562195秒
ランダムサーチ
最適パラメータ {'subsample_freq': 3, 'subsample': 0.5, 'reg_lambda': 0.001, 'reg_alpha': 0.001, 'num_leaves': 4, 'min_child_samples': 3, 'colsample_bytree': 0.7}
スコア -0.00036616522117512536
所要時間43.027225971221924秒
ベイズ最適化(BayesianOptimization)
最適パラメータ {'colsample_bytree': 0.9430083037075744, 'min_child_samples': 3, 'num_leaves': 2, 'reg_alpha': 0.0020869096632570856, 'reg_lambda': 0.0020382643364640087, 'subsample': 0.5853280069552931, 'subsample_freq': 2}
スコア -0.0003700799026914537
所要時間50.49155378341675秒
ベイズ最適化(Optuna)
最適パラメータ {'reg_alpha': 0.0023317132133850566, 'reg_lambda': 0.0006717623390065321, 'num_leaves': 5, 'colsample_bytree': 0.926079384200066, 'subsample': 0.4948413802207167, 'subsample_freq': 4, 'min_child_samples': 0}
スコア -0.00023752018179486973
所要時間48.97835183143616秒

いずれにせよ、Optunaが最も効率よくチューニングできていることが分かります。Optunaの性能恐るべしですね!

手順5 学習・検証曲線の確認

学習、検証曲線で過学習の有無を確認します。
今回はOptunaで調整したパラメータ(colsample_bytree=0.911‥, learning_rate=0.150‥)を例に描画してみます。

学習曲線

学習曲線をプロットし、
・目標性能を達成しているか
・過学習していないか
を確認します

学習曲線のプロット
from sklearn.model_selection import learning_curve
import matplotlib.pyplot as plt
# 最適パラメータを学習器にセット
model.set_params(**best_params)

# 学習曲線の取得
train_sizes, train_scores, valid_scores = learning_curve(estimator=model,
                                                         X=X, y=y,
                                                         train_sizes=np.linspace(0.1, 1.0, 10),
                                                         fit_params=fit_params,
                                                         cv=cv, scoring=scoring, n_jobs=-1)
# 学習データ指標の平均±標準偏差を計算
train_mean = np.mean(train_scores, axis=1)
train_std  = np.std(train_scores, axis=1)
train_center = train_mean
train_high = train_mean + train_std
train_low = train_mean - train_std
# 検証データ指標の平均±標準偏差を計算
valid_mean = np.mean(valid_scores, axis=1)
valid_std  = np.std(valid_scores, axis=1)
valid_center = valid_mean
valid_high = valid_mean + valid_std
valid_low = valid_mean - valid_std
# training_scoresをプロット
plt.plot(train_sizes, train_center, color='blue', marker='o', markersize=5, label='training score')
plt.fill_between(train_sizes, train_high, train_low, alpha=0.15, color='blue')
# validation_scoresをプロット
plt.plot(train_sizes, valid_center, color='green', linestyle='--', marker='o', markersize=5, label='validation score')
plt.fill_between(train_sizes, valid_high, valid_low, alpha=0.15, color='green')
# 最高スコアの表示
best_score = valid_center[len(valid_center) - 1]
plt.text(np.amax(train_sizes), valid_low[len(valid_low) - 1], f'best_score={best_score}',
                color='black', verticalalignment='top', horizontalalignment='right')
# 軸ラベルおよび凡例の指定
plt.xlabel('training examples')  # 学習サンプル数を横軸ラベルに
plt.ylabel(scoring)  # スコア名を縦軸ラベルに
plt.legend(loc='lower right')  # 凡例

image.png

こちらのの確認方法に基づいて確認すると、下記のようにOKと判断できます
・目的の性能を達成しているか → 目標性能は定めていないが、チューニング前より大幅向上(-0.001390 → -0.000268)したのでOKとする
・過学習していないか → データ数が少ない部分では逆転も見られるが、最終的に検証データ指標と学習データ指標が収束しているのでOK

検証曲線

検証曲線をプロットし、
・性能の最大値を捉えられているか
・過学習していないか
を確認します

検証曲線のプロット
from sklearn.model_selection import validation_curve
# 検証曲線描画対象パラメータ
valid_curve_params = {'reg_alpha': [0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1],
                      'reg_lambda': [0.0001, 0.0003, 0.001, 0.003, 0.01, 0.03, 0.1],
                      'num_leaves': [2, 3, 4, 5, 6],
                      'colsample_bytree': [0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
                      'subsample': [0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
                      'subsample_freq': [0, 1, 2, 3, 4, 5, 6, 7],
                      'min_child_samples': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
                      }
param_scales = {'reg_alpha': 'log',
                'reg_lambda': 'log',
                'num_leaves': 'linear',
                'colsample_bytree': 'linear',
                'subsample': 'linear',
                'subsample_freq': 'linear',
                'min_child_samples': 'linear'
                }
# 最適パラメータを上記描画対象に追加
for k, v in valid_curve_params.items():
    if best_params[k] not in v:
        v.append(best_params[k])
        v.sort()
for i, (k, v) in enumerate(valid_curve_params.items()):
    # モデルに最適パラメータを適用
    model.set_params(**best_params)
    # 検証曲線を描画
    train_scores, valid_scores = validation_curve(estimator=model,
                                                  X=X, y=y,
                                                  param_name=k,
                                                  param_range=v,
                                                  fit_params=fit_params,
                                                  cv=cv, scoring=scoring,
                                                  n_jobs=-1)
    # 学習データに対するスコアの平均±標準偏差を算出
    train_mean = np.mean(train_scores, axis=1)
    train_std  = np.std(train_scores, axis=1)
    train_center = train_mean
    train_high = train_mean + train_std
    train_low = train_mean - train_std
    # テストデータに対するスコアの平均±標準偏差を算出
    valid_mean = np.mean(valid_scores, axis=1)
    valid_std  = np.std(valid_scores, axis=1)
    valid_center = valid_mean
    valid_high = valid_mean + valid_std
    valid_low = valid_mean - valid_std
    # training_scoresをプロット
    plt.plot(v, train_center, color='blue', marker='o', markersize=5, label='training score')
    plt.fill_between(v, train_high, train_low, alpha=0.15, color='blue')
    # validation_scoresをプロット
    plt.plot(v, valid_center, color='green', linestyle='--', marker='o', markersize=5, label='validation score')
    plt.fill_between(v, valid_high, valid_low, alpha=0.15, color='green')
    # 最適パラメータを縦線表示
    plt.axvline(x=best_params[k], color='gray')
    # スケールをparam_scalesに合わせて変更
    plt.xscale(param_scales[k])
    # 軸ラベルおよび凡例の指定
    plt.xlabel(k)  # パラメータ名を横軸ラベルに
    plt.ylabel(scoring)  # スコア名を縦軸ラベルに
    plt.legend(loc='lower right')  # 凡例
    # グラフを描画
    plt.show()

image.png
image.png
image.png
image.png
image.png
image.png
image.png

こちらの確認方法に基づいて確認すると、下記のようにOKと判断できます
・性能の最大値を捉えられているか → どのパラメータもだいたい捉えられていそうなのでOKとする
・過学習していないか → データ数が少ないために多少の乖離はあるが、最適パラメータにおいて学習データ指標と検証データ指標が近づいているのでOK

※パラメータreg_lambdaに関してはチューニング後のパラメータ値よりも性能が高くなる部分がありますが、これはチューニング時のパラメータの組合せとは無関係に、検証曲線のパラメータの組合せを選んでいることに起因します。
ベイズ最適化においてはこのような事は往々に起こりますが、まだ改善の余地があると言えるので、気になる方はn_trialsを増やして再実行するのが良いかと思います

チューニング後のモデルを可視化

最後に、チューニング後のモデルを可視化します(チューニング前と比較)

チューニング後のモデルを可視化
regplot.regression_heat_plot(model, USE_EXPLANATORY, OBJECTIVE_VARIALBLE, df_osaka,
                             pair_sigmarange = 0.5, rounddigit_x1=3, rounddigit_x2=3,
                             cv=cv, display_cv_indices=0,
                             fit_params=fit_params, estimator_params=best_params)

lgbm_optuna_osaka.png

チューニング前は予測値が一定でうまく回帰できていませんでしたが、

チューニング後は男性比率(3_male_ratio)が高いほど、30-60代(2_between_30to60)が多いほど、北に(latitude)行くほど賛成率(approval_rate)が高い=色が濃いことが分かり、傾向が捉えられていることが分かります。
チューニング成功と言って良いのではと思います。

まとめ

・LightGBMのパラメータチューニングを、7パラメータを選択して実施
・4種類のアルゴリズムのうち、Optunaが速度、評価指標両方で最も優れていた
・一部のパラメータは学習時間に影響を与えるため、範囲を広くとるとグリッドサーチの速度低下に繋がる

おまけ:OptunaのLightGBM専用チューニング

Optuna1.5.0より、LightGBMをクロスバリデーションでチューニングするLightGBMTunerCVというクラスが追加されたようです。
パラメータを細かく指定する必要がないので簡単にチューニングでき、かつ「重要なハイパーパラメータから順に最適なハイパーパラメータをチューニングする」ことで、効率的にチューニングができるそうです。
詳しくは公式APIリファレンス公式の紹介記事を参照ください

使用法は以下の記事が分かりやすいです

公式のサンプルコードも使用法の参考になるかと思います

Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
79
Help us understand the problem. What are the problem?