TPEによるハイパーパラメータ最適化#

本稿では、「ハイパーパラメータ探索の基本的な手法」でのハイパーパラメータおよびその探索に関する基本的な解説、および「ベイズ最適化によるハイパーパラメータ最適化」でのハイパーパラメータ最適化に関する解説を前提に、tree-structured Parzen estimator (TPE) を用いたベイズ最適化によるハイパーパラメータ最適化について解説します。

TPEによるベイズ最適化#

所与の訓練データの下でハイパーパラメータ \(x \in \mathcal{X}\) によって指定されるプロセスで学習を実行し、その結果として得られる適当な指標 \(y\) を最小化したいものとします。このような場合によく用いられるのがsequential model-based optimization (SMBO) と呼ばれる逐次最適化の方法で、これまでに観測した \(( x, y )\) の組の履歴 \(\mathcal{H}\) を元にモデル \(M_t\) を作成し、このモデルから導かれるサロゲート関数 \(S ( x, M )\) を最小化する点として次に評価する点 \(x^*\) を提案します。

smbo

(図は [Bergstra et al., 2011] より引用)

例えばSMBOの一種であるガウス過程によるベイズ最適化では

\[\begin{split} \begin{align} f & \sim \mathrm{GP} (0, K), & & \\ \epsilon_i & \sim \mathrm{N} (0, \sigma^2_\epsilon), & & \text{i.i.d. for } i = 1, 2, \ldots \\ y_i & = f ( x_i ) + \epsilon_i, & & \text{for } i = 1, 2, \ldots \end{align} \end{split}\]

のようなモデリングを行い、この下で、履歴 \(\mathcal{H}\) で条件付けた事後分布 \(p ( y | x )\) を計算します。この \(p ( y | x )\) から期待改善量 (expected improvement; EI)

\[ \mathrm{EI}_{y^*} ( x ) = \int_{- \infty}^\infty \max ( y^* - y, 0 ) p ( y | x ) d y \]

などの値を計算して、次に探索する \(x^*\) としてはこの \(\mathrm{EI}_{y^*} ( x )\) を最大化する点(即ち、サロゲート関数 \(S ( x ) = - \mathrm{EI}_{y^*} ( x )\) を最小化する点)を提案します。ここで \(y^*\) は一般に履歴 \(\mathcal{H}\) の中に含まれる最良の(最小の) \(y\) の値を用います。

この節では、同様の問題設定の下で用いられるまた別の確率モデルとしてtree-structured Parzen estimator (TPE) を導入し、そのモデルの下でのEIなどのサロゲート関数の計算方法を簡単に説明します。加えてハイパーパラメータ最適化におけるTPEの利用について解説し、TPEを実装した代表的なモジュールであるOptunaを利用したハイパーパラメータ探索の実装例を示します。

Tree-structured Parzen Estimator#

Tree-structured Parzen estimator (TPE) ではガウス過程によるベイズ最適化の場合と異なって \(p ( x | y )\)\(p ( y )\) をモデリングし、これらを用いて期待改善量を計算します。具体的には前者の条件付き分布

\[\begin{split} p ( x | y ) = \begin{cases} l ( x ), & \text{if } y < y^* \\ g ( x ), & \text{if } y \ge y^* \end{cases} \end{split}\]

とおき、関数 \(l ( x )\) および \(g ( x )\) は履歴 \(\mathcal{H}\) からのカーネル密度推定(Parzen推定)によって与えます。ここで閾値 \(y^*\) は適当な \(0 < \gamma < 1\) について確率 \(p ( y < y^* ) = \gamma\) となるように与えるものとします。一方で、後者の \(p ( y )\) については明示的なモデリングを行いません。これは実際の逐次最適化のアルゴリズムにおいて \(y^*\) の値を \(\mathcal{H}\) の経験分布から求めているのと、後述するEIの最大化において \(p ( y )\) の具体的な形を考慮する必要がないためです。

このモデリングにおけるEIの最大化について、

\[\begin{split} \begin{align} p ( x ) & = \int_{- \infty}^{\infty} p ( x | y ) p ( y ) d y \\ & = \gamma l ( x ) + ( 1 - \gamma ) g ( x ) \end{align} \end{split}\]

であることに注意すると、EIの式は

\[\begin{split} \begin{align} \mathrm{EI}_{y^*} ( x ) & = \int_{- \infty}^{y^*} ( y^* - y ) p ( y | x ) d y \\ & = \int_{- \infty}^{y^*} ( y^* - y ) \frac{p ( x | y ) p ( y )}{p ( x )} d y \\ & = \frac{l ( x )}{p ( x )} \int_{- \infty}^{y^*} ( y^* - y )p ( y ) d y \\ & \propto \left( \gamma + \frac{g ( x )}{l ( x )} ( 1 - \gamma ) \right)^{- 1} \end{align} \end{split}\]

と変形できます。即ちEIを最大化する \(x\)\(g ( x ) / l ( x )\) の値を最小化する \(x\) ですから、 \(l ( x )\) の値が大きく \(g ( x )\) の値が小さいほど好ましいことになります。いま \(l ( x )\) はカーネル密度推定で求めた確率密度関数だったことを思い出すと、EIの最大値を与える \(x\)

  1. \(l ( x )\) に従う \(x\) を大量にサンプリングする

  2. サンプリングした \(x\) の中から \(g ( x ) / l ( x )\) の値を最小にするものを選ぶ

という方法で近似的に求めることができます。実際のアルゴリズムでもこの方法で次に探索する \(x^*\) を提案しています。

ハイパーパラメータ最適化への応用#

TPEを用いてハイパーパラメータ最適化を行う上でのテクニックや実用上の注意点について簡単に説明します。TPEを用いる場合、ガウス過程によるベイズ最適化を用いる場合と比べてこうした注意点は概ね少なくて済む傾向があります。

入力変数の変換#

ガウス過程によるベイズ最適化の場合と同様、手動での探索やgrid searchの際に指数的に変化させるハイパーパラメータ(多層パーセプトロンの各層におけるユニット数や学習率など)にはしばしば対数変換が行われます。

一方で、各ハイパーパラメータの値(や、それを対数変換した値)に対して標準化正規化の処理が行われることはほとんどありません。これはTPEで \(l ( x )\)\(g ( x )\) を求めるときに行われるカーネル密度推定が適応的なもの、具体的には履歴 \(\mathcal{H}\) 中の各 \(x^{(n)}\) に対応するカーネルの幅が近くの他の \(x^{(n')} \in \mathcal{H}\) との距離に応じて決められるものであることが多いためです。このためTPEは各ハイパーパラメータの線形なスケーリングの影響を比較的受けにくくなっています。

履歴の初期化#

これもガウス過程によるベイズ最適化の場合と同様、実際にはTPEを用いて \(x^*\) を提案するループ処理の開始前に履歴 \(\mathcal{H}\) の初期化が行われます。この際、初期化において \(\mathcal{H}\) を構成する \(x^{( n )}\) はランダムサンプリングで決定することが多く、各パラメータごとに適当な上限・下限を定めた上で独立な(連続または離散の)一様分布がよく利用されます。なお、対数変換を施すハイパーパラメータについては対数変換後の空間における一様分布(即ち、変換前の空間における対数一様分布)からのサンプリングが行われます。

TPE自身のハイパーパラメータの選択#

TPEによるハイパーパラメータ最適化のアルゴリズム自身も、 \(l (x)\) を構成するデータの比率を決める \(\gamma\) や近似的な最適化に際して \(l ( x )\) からサンプリングする個数などのハイパーパラメータを持ちます。これらの値は事前に適切に決定するしかありませんが、一例として \(\gamma\) の値は [Bergstra et al. 2011] の例では \(\gamma = 0.15\) が、 [Bergstra et al. 2013] や後述するOptunaのデフォルト値としては \(\gamma = 0.25\) が採用されています。

Optunaでの実装#

PythonでTPEによるベイズ最適化を実装しているライブラリとしてはHyperoptOptunaが有名で、よく利用されています。本節ではこのうちOptunaを利用した実装例について簡単に解説します。

OptunaはPreferred Networks (PFN) 社が中心となって開発しているオープンソースのハイパーパラメータ最適化のライブラリで、TPEのほかにも複数の最適化アルゴリズムが搭載され、しかも容易に分散並列化ができるように実装されています。本節では

  • 引数にtrialを取る関数オブジェクトを定義する方法、および

  • scikit-learnのestimator向けのラッパーである OptunaSearchCV を利用する方法

の2種類の実装をこの順に見ていきます。

import numpy as np
import optuna
from sklearn.base import clone
from sklearn.datasets import load_diabetes
from sklearn.decomposition import PCA
from sklearn.dummy import DummyRegressor
from sklearn.linear_model import Ridge, Lasso
from sklearn.model_selection import cross_validate
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
import warnings


warnings.filterwarnings("ignore")
optuna.logging.disable_default_handler()  # notebook上で使うのでログ表示を停止
feature, target = load_diabetes(return_X_y=True)  # 糖尿病データセットを利用
feature_train, target_train = feature[:300], target[:300]  # 訓練データ
feature_test, target_test = feature[300:], target[300:]  # 試験データ

基本的な実装#

Optunaで最適化を行う場合の基本的な方法は optuna.trial.Trial クラスのインスタンスを引数に取る関数を実装することです。関数の中ではこのインスタンスのメソッドである

  • Trial.suggest_categorical() : カテゴリ値を取るハイパーパラメータを提案する

  • Trial.suggest_int() : 整数値を取るハイパーパラメータを提案する

  • Trial.suggest_float() : 実数値を取るハイパーパラメータを提案する

などを用いて最適化の対象となるハイパーパラメータを定義していき(いわゆるdefine-by-run方式)、返り値としてはそのハイパーパラメータを入力したときのスコア(例えば回帰問題ならMSEのCV平均など)が返るようにします。

def objective(trial):
    x, y = feature_train, target_train

    pca__n_components = trial.suggest_int(  # 整数値のハイパーパラメータを提案
        "pca__n_components",  # ハイパーパラメータの名前
        1,
        9,  # ハイパーパラメータの取りうる範囲(最小値、最大値の順)
    )
    transformer = PCA(n_components=pca__n_components)
    x_pca = transformer.fit_transform(x)

    regressor_name = trial.suggest_categorical(  # カテゴリ値のハイパーパラメータを提案
        "regressor",  # ハイパーパラメータの名前
        ["Ridge", "Lasso"],  # ハイパーパラメータの取りうる値のリスト
    )
    if regressor_name == "Ridge":
        ridge__alpha = trial.suggest_float(  # 実数値のハイパーパラメータを提案
            "ridge__alpha",  # ハイパーパラメータの名前
            1e-4,
            1.0,  # ハイパーパラメータの取りうる範囲(最小値、最大値の順)
            log=True,  # 対数変換したスケールでカーネル密度推定
        )
        regressor = Ridge(alpha=ridge__alpha)
    else:  # regressor_name == 'Lasso'
        lasso__alpha = trial.suggest_float(  # 実数値のハイパーパラメータを提案
            "lasso__alpha",  # ハイパーパラメータの名前
            1e-4,
            1.0,  # ハイパーパラメータの取りうる範囲(最小値、最大値の順)
            log=True,  # 対数変換したスケールでカーネル密度推定
        )
        regressor = Lasso(alpha=lasso__alpha)

    scores = cross_validate(
        estimator=regressor,
        X=x_pca,
        y=y,
        cv=3,  # 3-fold CVで評価
        scoring="neg_mean_squared_error",
    )
    return scores["test_score"].mean()  # スコアのCV平均のみを返す

このように作成した関数を create_study() 関数に渡して optuna.study.Study クラスのインスタンスを作成し、そのインスタンスの Study.optimize() メソッドを実行することで最適化が実行されます。最適化のアルゴリズムは create_study() 関数の sampler 引数で指定でき、特に指定しない場合はTPEによるベイズ最適化(ただし初期化のランダムサンプリングを10回)が行われます。探索した中での最良の結果は同じインスタンスの Study.best_trial 属性に格納されています。

study = optuna.create_study(
    sampler=optuna.samplers.TPESampler(seed=42),  # ランダムシードを固定する場合
    direction="maximize",  # 目的関数は負のMSEのCV平均なので、最大化
)
study.optimize(objective, n_trials=100)
print(study.best_trial)
FrozenTrial(number=84, values=[-3067.5055000831853], datetime_start=datetime.datetime(2022, 12, 15, 17, 55, 30, 480643), datetime_complete=datetime.datetime(2022, 12, 15, 17, 55, 30, 489856), params={'pca__n_components': 6, 'regressor': 'Lasso', 'lasso__alpha': 0.019872362794542697}, distributions={'pca__n_components': IntUniformDistribution(high=9, low=1, step=1), 'regressor': CategoricalDistribution(choices=('Ridge', 'Lasso')), 'lasso__alpha': LogUniformDistribution(high=1.0, low=0.0001)}, user_attrs={}, system_attrs={}, intermediate_values={}, trial_id=84, state=TrialState.COMPLETE, value=None)

scikit-learn用のラッパーを利用した実装#

また、scikit-learnのハイパーパラメータを最適化する場合には専用のラッパーである OptunaSearchCV クラスを利用することも可能です。このクラスはscikit-learnに実装されている GridSearchCVRandomizedSearchCV などの基本的なハイパーパラメータ探索のクラスと概ね同様に利用できます。ただし OptunaSearchCV の場合はこれらのクラスと違って引数 parameter_distributions に辞書のリストを取ることができないので、上の例のように複数の回帰モデルにわたってハイパーパラメータの最適化をしたい場合はやや煩雑な実装が必要になります。

pipeline = Pipeline(  # ハイパーパラメータを最適化するestimatorを定義
    steps=[
        ("pca", PCA()),
        ("ridge", Ridge()),
    ]
)
param_distributions = (
    {  # ハイパーパラメータの従う分布を指定(初期化時などに利用、分布の台を探索)
        "pca__n_components": optuna.distributions.IntUniformDistribution(1, 9),
        "ridge__alpha": optuna.distributions.LogUniformDistribution(1e-4, 1.0),
    }
)
optuna_search = optuna.integration.OptunaSearchCV(
    estimator=pipeline,
    param_distributions=param_distributions,
    cv=3,
    random_state=42,  # ランダムシードの値を固定する場合
    scoring="neg_mean_squared_error",  # 負のMSE(のCV平均)を最大化するように探索
)

optuna_search.fit(feature_train, target_train)
print(optuna_search.study_.best_trial)
FrozenTrial(number=2, values=[-3077.3934746929845], datetime_start=datetime.datetime(2022, 12, 15, 17, 55, 30, 678085), datetime_complete=datetime.datetime(2022, 12, 15, 17, 55, 30, 686089), params={'pca__n_components': 5, 'ridge__alpha': 0.0015597195848142719}, distributions={'pca__n_components': IntUniformDistribution(high=9, low=1, step=1), 'ridge__alpha': LogUniformDistribution(high=1.0, low=0.0001)}, user_attrs={'mean_fit_time': 0.0011289119720458984, 'std_fit_time': 6.443920970792569e-05, 'mean_score_time': 0.00039776166280110675, 'std_score_time': 3.60642198171362e-05, 'split0_test_score': -2871.2731967341915, 'split1_test_score': -2904.6304792130377, 'split2_test_score': -3456.2767481317233, 'mean_test_score': -3077.3934746929845, 'std_test_score': 268.2568150638913}, system_attrs={}, intermediate_values={}, trial_id=2, state=TrialState.COMPLETE, value=None)

小括#

本稿ではTPEを利用したハイパーパラメータ最適化の手法について簡単に解説し、あわせてOptunaを利用した実装について紹介しました。ガウス過程によるベイズ最適化ではパラメータ空間が大きく、複雑な構造を持つようになるとカーネル関数の構成で頭を悩ませることになりましたが、TPEによる場合はそのような探索空間でもである程度機械的に処理できます。このため、現在ではハイパーパラメータの数が多い場合(文献にもよりますが、5個から10個程度より多い場合)にはTPEを利用してハイパーパラメータの最適化を行うことが一般的になっています。

なお、機械学習全般のハイパーパラメータ探索としてはTPEやガウス過程によるベイズ最適化、あるいはCMA-ESなどを利用した進化計算による最適化がよく用いられますが、ニューラルネットワークの層数・ユニット数といったハイパーパラメータの場合はNeural Architecture Search (NAS) などの全く異なる方法で探索されることもあります。いずれにせよ、探索すべきハイパーパラメータの性質や探索に利用できる計算資源も考慮しながら適切な手法を選択することが重要です。

参考文献#