ブラックボックス(ベイズ)最適化とハイパラの自動調整
この記事は学ロボアドベントカレンダー2023最終日の記事です。
学ロボアドカレ2023、多数のご参加ありがとうございました!!
さて、2023年も年の瀬が迫って参りました。
皆様いかがお過ごしでしょうか?
私はこの年末が学部生活の冬という現実に戦々恐々としています。
この記事では、ハイパーパラメータの自動調整に用いられるblackbox最適化を概説します。
blackbox最適化の中でも特にベイズ最適化を中心に解説し、optunaを用いた実装についても解説します。
環境としては、python 3.11.6を想定しています。
==============目次==============
===============================
1.ブラックボックス最適化とは
式で陽に表現できない関数をblackbox関数と言い、
blackbox関数に対する最適化問題をblackbox最適化と言います。
以降、blackbox関数をと表す事にします。
特にblackbox最適化では以下の様な状況を考える事が多いです。
・が観測できない
・の評価コストが高い
・の観測にノイズが乗る (が観測される)
勾配やヘッセ行列などが得られないので、一般的な最適化手法を適用する事は出来ません。
blackbox最適化では、この様な関数に対してを最小にするを見つける事を考えます。
(例) blackbox最適化の例としては、何らかの問題のハイパーパラメータの自動調整が挙げられます。
この記事の後半では、ハイパラの自動調整に広く用いられるoptunaで実際にblackbox最適化を実装します。
2. ベイズ最適化
この記事では、特にblackbox最適化の中でも広く用いられているベイズ最適化を扱います。
ベイズ最適化では、blackbox関数に統計モデルを仮定し、を分布からのサンプルとして捉えます。
具体的には、以下の手続きでを最小化する点を求めます。
1. 観測データから統計モデルを更新
2. 統計モデルから獲得関数を求め、獲得関数を最大化する点を次の入力として選択
3. 選択した点を入力して、評価値を観測し、とする。
統計モデルは入力に対するの分布の近似に使用し、
獲得関数は探索と活用のトレードオフを考慮した入力を定めるのに使用します。
強化学習などでもおなじみ探索と活用のトレードオフはここでも登場します。
勿論UCB系のアプローチも使われます。(今回はminimizeなのでLCB)
3. ベイズ最適化の手法
ベイズ最適化には、統計モデルと獲得関数の組み合わせで様々な種類の手法が存在します。
ここでは、ベイズ最適化でよく使われる手法を統計モデル/獲得関数のそれぞれで列挙します。
3-1. 統計モデル
まず統計モデルから概説します。
統計モデルを用いてblackbox関数をモデル化します。
ガウス過程回帰(GPR)
ガウス過程(GP)というのは、入力空間上のランダムな関数を定める確率過程です。
出力の平均値を定める関数 とカーネル関数 から構成され、
入力系列 に対する、関数の出力系列が
となる時、関数はガウス過程に従うと言います。
ただし、、とします。
以降では、が得られている時にGPの平均を表す関数を、分散を表す関数をとして、
関数がガウス過程に従う事をと表す事にします。
更に、観測に乗るノイズの分散をとして、
を定義します。()
ガウス過程回帰(GPR)は、ガウス過程と正規分布の再生性・同時分布の性質を用いて、
ガウス過程のをオンラインで更新する事で関数回帰を行う手法です。
GPRでは、の事前分布をとし、観測ノイズとすると、
が観測された場合の事後分布は以下の様に表されます。
・
・
ただし、としとしています。
上記の様にガウス過程としてモデル化する事で、回帰する対象の分散(不確実さ)も表現出来るのがメリットです。
ただし、観測した点数について更新にかかる計算量がと計算量が大きいです。
TPE(Tree-structured Parzen Estimator)
TPEでは、のみをモデル化してを求めます。
獲得関数にはEI(詳細は後述)を用います。
具体的には、
としてモデル化し、
関数はを用いたParzen推定(カーネル密度推定法の一種)で求めます。
は陽にモデル化しない事に注意して下さい。
ただし、はとなる様に与えた閾値とします。(の上位となるとか)
optunaでは、がデフォルト値として採用されています。
がよりも小さいか否かを考えたいので、この様に分割しています。
このモデルを用いたEIを計算すると以下の様になります。
(TPEは獲得関数までセットなので、ここで解説します)
※ に注意。
以上の結果から、を最大化するにはを小さくする点を考えれば良いと分かります。
特にの定義から、以下の手続きで獲得関数を最大化する点が求まります。
1. を複数サンプリングする
2. サンプリングしたから、最もを小さくする点を次の入力として採用
以上がEIとTPEを用いたベイズ最適化法の手続きとなります。
ガウス過程回帰と比べて計算量が小さく、観測した点数が多い場合にも使えるのがメリットです。
Optunaのdefaultの最適化アルゴリズムです。
こちらの記事が分かりやすいです。
gochikika.ntt.com
3-2. 獲得関数
獲得関数についてもまとめておきます。
統計モデルから求めた獲得関数を最大化する点を次の入力に用います。
EI(Expected Improvement)
EI(期待改善量)はblackbox最適化で用いられる最も代表的な獲得関数です。
単純かつ比較的高速に動作するので、手っ取り早く試すのに向いている手法です。
EIの獲得関数は、閾値に対して]として定義されます。
には様々な決め方が存在しますが、をに含まれる観測値として、
と定める方法が用いられたりします。
LCB(Lower Confidence Bound)
多碗bandit問題や強化学習でもおなじみの考え方です。
LCBは信頼区間の下限を用いた獲得関数であり、
GPRなど明示的に分散も求まる統計モデルを用いて以下の様に定義されます。
ただし、はハイパーパラメータで などとして定めます。
(最大化問題を考える場合はUCBを考えれば良いです。)
LCBは"不確かな時は楽観的に"というヒューリスティックに基づいていて、
標準偏差を期待値に足し合わせた評価値で探索と活用のトレードオフを取ります。
直感的かつ単純な実装で実現できる事や、リグレット等の理論解析のしやすさがメリットとなりますが、
実用上性能がそこまで良くない事も実験的に知られています。
Thompson Sampling
Thompson Samplingは、目的関数の事後分布から関数を一つサンプリングし、
その関数を最小化する点を次の入力として採用する手法です。
こちらもLCB/UCBと合わせて多碗バンディット問題でよく使われます。
※ blackbox最適化は連続腕bandit問題とも解釈出来るので、banditの手法が幾つか登場します。
ガウス過程を用いる場合、Thompson Samplingの獲得関数は以下の様に定義されます。
は確率的に定まる関数であり、直接を最大化するのは難しいので、
今回はを有限個に分割して近似的にを計算します。
まず、として互いに素な集合族を取り、
各集合から代表点を取ります。
この代表点を用いて、として有限個の点でを近似します。
さらに、を用いてに対応するをサンプルし、
によって、関数をとして有限個の点で近似します。
これらの近似によって、以下の単純な計算でを最大化する点を近似的に求められる様になります。
勿論この計算法では予め集合の分割法を定める必要があります。
※ を最大化するは、を最小化する点の事後分布に従う点としても解釈できます。
この観点から捉えると、banditでのThompson Samplingとの関連性が分かりやすいと思います。
Thompson Samplingは多碗bandit問題では有力な手法なのですが、
ベイズ最適化ではそれほど使われていない様です。(高次元で性能が悪い)
※ banditの手法の比較は以下で試せます。(気になる方は試してみて下さい)
github.com
(補足)獲得関数の最小化について
特にGPRを用いた獲得関数の最小化について補足します。
先程定義した、を考えます。
これらの関数のについての勾配は以下の様に求まります。(としています)
これらの勾配はを事前計算しておく事で、計算量をまで削減できます。
EIやLCBの獲得関数の勾配はの勾配を用いて表現されるので、
上記の行列を事前計算しておく事で、獲得関数の最適化計算の計算量を削減できます。
今回の実装では事前計算を行っておりません
4. 実装(GP-LCB)
GPR+LCBを用いた場合のプログラムを実装してみます。
この手法はGP-LCB(GP-UCB)と呼ばれ、連続腕バンディット問題の方策の1つです。
ここからは、を最小化するをblackbox最適化で求める事を考えます。
ただし]とします。
GPRはscikit-learnのGaussianProcessRegressorを使用し、
kernel関数にはmatern kernelにwhite kernel、Constant kernelを組み合わせた物を用います。
matern kernelの定義は以下を参照して下さい。
sklearn.gaussian_process.kernels.Matern — scikit-learn 1.3.2 documentation
また、LCBのにはを使用し、
の観測ノイズはとします。
import numpy as np import matplotlib.pyplot as plt from sklearn.gaussian_process import GaussianProcessRegressor from sklearn.gaussian_process.kernels import Matern, WhiteKernel, ConstantKernel from scipy.optimize import minimize from functools import partial def f(x: np.ndarray) -> np.ndarray: x = x.reshape(-1) return np.exp(x/2.0) * np.sin(2 * np.pi * x) def obs(value: np.ndarray, noise: bool | float = False) -> np.ndarray: return value + np.random.normal(0, noise, size=value.shape) if noise else value def plot_data(xs: np.ndarray, ys: np.ndarray, x_bound: np.ndarray, gpr: GaussianProcessRegressor, n) -> None: x = np.linspace(x_bound[:, 0], x_bound[:, 1], 200).reshape(-1) y, sigma = gpr.predict(x.reshape(-1, 1), return_std=True) plt.plot(x, f(x), label='f(blackbox)') plt.plot(x, y, label='GPR(model)') plt.fill_between(x, y - sigma, y + sigma, alpha=0.2, label='sigma') plt.fill_between(x, y - sigma * np.sqrt(2.0 * np.log(n)), y + sigma * np.sqrt(2.0 * np.log(n)), alpha=0.2, label='UCB/LCB') plt.scatter(xs, ys, marker='x', color='red') plt.xlabel("x") plt.ylabel("f(x)") plt.legend() plt.show() def lcb_policy(x_bound: np.ndarray, gpr: GaussianProcessRegressor, n: int) -> np.ndarray: def lcb_score(x: np.ndarray, gpr: GaussianProcessRegressor, n: int) -> float: mean, sigma = gpr.predict(x.reshape(1, -1), return_std=True) score = (mean - sigma * np.sqrt(2.0 * np.log(n)))[0] return score objective = partial(lcb_score, gpr=gpr, n=n) point_num = 30 x_next_candidates: list = [] for x in np.random.uniform(x_bound[:, 0], x_bound[:, 1], size=(point_num, x_bound.shape[0])): # 多点スタートで最適化 (local minimum対策) ans = minimize(objective, x0=x, bounds=x_bound) x_next_candidates.append([ans.fun, ans.x]) x_next = min(x_next_candidates)[1] return x_next def GP_LCB(x_bound: np.ndarray, trial_num: int): kernel = Matern(nu=1.5) + WhiteKernel() + ConstantKernel(1.0) # カーネルの設定 gpr = GaussianProcessRegressor(kernel=kernel, alpha=0, random_state=0) # モデルの作成 # observed data xs: np.ndarray = np.empty(0) ys: np.ndarray = np.empty(0) # 初期値の設定 x_next: np.ndarray = np.random.uniform(x_bound[:, 0], x_bound[:, 1], size=(x_bound.shape[0],)) for i in range(trial_num): print(f"trial: {i+1}") y_next = obs(f(x_next), noise=0.01) # 観測値の取得 xs = np.append(xs, x_next).reshape(-1, x_next.shape[0]) ys = np.append(ys, y_next) gpr.fit(xs, ys) # モデルの学習 plot_data(xs, ys, x_bound, gpr, i+1) x_next = lcb_policy(x_bound, gpr, i + 1) return min(ys) if __name__ == '__main__': x_bound_array = np.array([[-1.0, 1.0]]) GP_LCB(x_bound_array, 20)
※ GPRのkernelの選び方などGP-LCB自体にもハイパラが多く、上手く最小化出来ない時もあります。
最適化計算等めんどうな部分はscikit-learnのモジュールで誤魔化しています。
本来GPRを用いた場合の獲得関数の最適化計算には逆行列の事前計算テクが使えるのですが、
今回は省略させて頂きました。
以下がプログラムの実行結果となります。
青の帯が分散を表し、オレンジの帯の上端がUCBを下端がLCBを表しています。
赤のxが入力した点を表しています。
反復に従って、を最小化する点が求まっている事が分かります。
5. optunaを用いた最適化
・optuna公式
optuna.org
5-1. install
以下のコマンドでinstall出来ます。
pip install optuna
dashboardで可視化したい人は以下もinstallして下さい。
pip install optuna-dashboard
installで躓く事は無いと思います。
5-2. toy-modelに対する最適化(sample code)
簡単な関数に対して、optunaを用いたblackbox最適化を実装してみます。
先ずは1変数で考えます。
最適化したい変数をtrial objectで記述し、以下のcodeで最適化します。
import optuna def objective(trial): x = trial.suggest_float('x', -10, 10) # xの範囲を[-10, 10]とする return (x - 2) ** 2 study = optuna.create_study() # problem setting study.optimize(objective, n_trials=100) # blackbox最適化の実行 study.best_params # 結果
結果はstudy objectのパラメータに格納されます。
optuna.study.Study — Optuna 3.5.0 documentation
また、optunaでは以下の様にアルゴリズムが選択出来ます。
Efficient Optimization Algorithms — Optuna 3.5.0 documentation
デフォルトではTPEを使用します。
2変数になっても同じ様に実装出来ます。
import optuna def func(x, y): return x**2+y**2+x*y def objective(trial): x = trial.suggest_float('x', -1, 1) # xの範囲指定 y = trial.suggest_float('y', -1, 1) # yの範囲指定 return func(x, y) study = optuna.create_study() study.optimize(objective, n_trials=100) study.best_params
簡単にblackbox最適化を試す事が出来ます。
5-3. 結果の可視化
以下のコードで結果を可視化出来ます。
# Empirical Distribution Function Plot fig = optuna.visualization.plot_edf(study) fig.write_image('edf.png') # Optimization History Plot fig = optuna.visualization.plot_optimization_history(study) fig.write_image('optimization_history.png') # parallel coordinate plot fig = optuna.visualization.plot_parallel_coordinate(study) fig.write_image('parallel_coordinate.png') # slice fig = optuna.visualization.plot_slice(study) fig.write_image('slice.png')
こちらの記事に詳しく説明されています。(コードもこちらの記事から引用しました。)
ibetoge.hatenablog.com
5-4. マルチプロセス化
optunaはRDBを用いたマルチプロセス化にも対応しています。
今回はSQLiteを用いて実装します。
以下のコマンドでデータベースを生成する事が出来ます。(pythonからでも生成出来ます)
optuna create-study --study 'optuna_test' --storage 'sqlite:///optuna_testdb.db'
pythonでも生成出来ます。
DATABASE_URI = 'sqlite:///optuna_testdb.db' study_name = 'optuna_test' optuna.create_study(study_name=study_name, storage=DATABASE_URI)
データベースが生成できたら、次はmulti-processでの実行を実装します。
今回はmultiprocessingを用いて実装します。
multi-processでoptunaを実行するコードは以下の様になります。
簡単かつ直感的に実装出来る様になっています。
import optuna import multiprocessing from multiprocessing import Process def func(x, y): return x ** 2 + y ** 2 + x * y def objective(trial): x = trial.suggest_float('x', -1, 1) # xの範囲指定 y = trial.suggest_float('y', -1, 1) # yの範囲指定 return func(x, y) def optimize(study_name, storage, n_trials): study = optuna.create_study(study_name=study_name, storage=storage, load_if_exists=True) study.optimize(objective, n_trials=n_trials) n_trials = 10 concurrency = multiprocessing.cpu_count() # CPUの数を取得 n_trials_per_cpu = n_trials / concurrency DATABASE_URI = 'sqlite:///optuna_testdb.db' study_name = 'optuna_test' optuna.create_study(study_name=study_name, storage=DATABASE_URI) # dbの生成 # multi-process化 workers = [Process(target=optimize, args=(study_name, DATABASE_URI, n_trials_per_cpu)) for _ in range(concurrency)] for worker in workers: worker.start() for worker in workers: worker.join() study = optuna.load_study(study_name=study_name, storage=DATABASE_URI) # 結果の取得 fig = optuna.visualization.plot_edf(study) fig.write_image('edf.png') # Optimization History Plot fig = optuna.visualization.plot_optimization_history(study) fig.write_image('optimization_history.png') # parallel coordinate plot fig = optuna.visualization.plot_parallel_coordinate(study) fig.write_image('parallel_coordinate.png') # slice fig = optuna.visualization.plot_slice(study) fig.write_image('slice.png')
公式Documentにも記述があります。
Easy Parallelization — Optuna 3.5.0 documentation
参考文献
この記事を書くに当たって参考にした解説/本を列挙しておきます。
・blackbox最適化の幅広いサーベイ
search.ieice.org
・ベイズ最適化
www.kindaikagaku.co.jp
・ガウス過程回帰
www.jstage.jst.go.jp
・bandit問題
www.kspub.co.jp