🍣
Gpyoptによるベイズ最適化でフィッティング
GPyOptによるベイズ最適化+実験計画 (1) ~2変数~
参考にしたリンク
- 【GPyOpt】ベイズ最適化で次の実験条件を導出する!
https://qiita.com/oki_kosuke/items/e68d5f68c972c3af32cc
ライブラリ
- python 3.8.13
- matplotlib 3.1.3
- gpy 1.10.0
- gpyopt 1.2.6
- numpy 1.22.3
- pandas 1.2.5
(注) GPyOptのplot_aquisition()でエラーが発生。Matplotlib 3.1.3以前では回避できるとのことから各ライブラリをダウングレードして使用。
やりたいこと
-
の放物線プロットを作成(模擬実験データ)ax^2+b -
と仮定してある(a,b)の組に対し計算上のyとデータのyの差の二乗を計算。これを目的変数とする。ax^2+b - ベイズ最適から(a,b)の組に対し、最も目的変数の小さい組み合わせを算出する。すなわち(a, b)のフィティングを行う
コード
# 必要なライブラリをインポート
import GPy
import GPyOpt
import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
%matplotlib inline
# 元となる関数を定義
def func_test(a,b,x):
y = a*x**2+b
return y
# a=2.0, b = 1.3として、xとyの関係を求める。
a= 2.0
b= 1.3
# 下記x_expとy_expを実験データとして扱う。
x_exp = np.linspace(-1.0,1.0,100)
y_exp = func_test(a,b,x_exp)
plt.scatter(x_exp, y_exp)
# 誤差を算出する関数を作る
def func_def(a,b,x_exp,y_exp):
y_calc = func_test(a,b,x_exp) #x_expがndarrayなのでここもndarrayになる。
delta_sq = (y_calc - y_exp)**2
sum_delta = np.sum(delta_sq)
return sum_delta
# 初点を(a,b)=(0,0)とする。2D arrayであることに注意。
ab_pred = np.asarray([[0.0, 0.0]])
# (a,b)=(0,0)のときのδ^2を算出する。2D arrayのデータを取り出して計算後、再度2D arrayとして保存
delta_calc = np.asarray([[func_def(ab_pred[0,0],ab_pred[0,1],x_exp,y_exp)]])
まずは、初期データから次の予測ができることを確認。
# 保有データから次のデータの予測を1回 行う
bounds = [{'name': 'a', 'type': 'continuous', 'domain': (-5.0, 5.0)},
{'name': 'b', 'type': 'continuous', 'domain': (-5.0, 5.0)}]
params = {'acquisition_type':'EI',#獲得関数としてEIを指定
'f':None, #最適化する関数の設定(実験結果は分からないので設定しない).
'domain':bounds, #パラメータの探索範囲の指定
'model_type':'GP', #標準的なガウシアンプロセスを指定.
'X':ab_pred, #既知データの説明変数(インプットするx)
'Y':delta_calc, #既知データの目的変数(インプットするy)
'de_duplication':True, #重複したデータをサンプルしないように設定.
"normalize_Y":True, #defaltはTrue
"exact_feval":False #defaltはFalse
}
myBopt = GPyOpt.methods.BayesianOptimization(**params)
myBopt.plot_acquisition()
x_suggest = myBopt.suggest_next_locations(ignored_X=ab_pred)
y_predict = myBopt.model.model.predict(x_suggest) #y_predictは(予測平均,予測分散)がタプルで返ってくる
y_mean=y_predict[0]
y_variance=y_predict[1]
print(f'つぎは{x_suggest}、yの予測平均は{y_mean}、分散は{y_variance}')
つぎは[[-3.02910399 1.09674807]]、yの予測平均は[[0.]]、分散は[[1.e-09]]
と出ればOK。
# 更新用の2d array作成。初期データに推奨a,bのペア、および新規目的変数を作成
ab_new_array = np.vstack((ab_pred,x_suggest)) # np.stack(())が二重になる点注意。
delta_new = np.asarray([[func_def(x_suggest[0,0],x_suggest[0,1],x_exp,y_exp)]])
delta_new_array = np.vstack((delta_calc,delta_new)) # np.stack(())が二重になる点注意。
# 更新データから次のデータの予測を行う
# 新しいパラメータ
num = 10 #追加回数
for n in range(num):
params_2 = {'acquisition_type':'EI',#獲得関数としてEIを指定
'f':None, #最適化する関数の設定(実験結果は分からないので設定しない).
'domain':bounds, #パラメータの探索範囲の指定
'model_type':'GP', #標準的なガウシアンプロセスを指定.
'X':ab_new_array, #既知データの説明変数(インプットするx)ここが変更
'Y':delta_new_array, #既知データの目的変数(インプットするy)ここが変更
'de_duplication':True, #重複したデータをサンプルしないように設定.
"normalize_Y":True, #defaltはTrue
"exact_feval":False, #defaltはFalse
'maximize':False
}
myBopt_2 = GPyOpt.methods.BayesianOptimization(**params_2)
ab_suggest_2 = myBopt_2.suggest_next_locations(ignored_X=ab_new_array)
delta_predict_2 = myBopt_2.model.model.predict(ab_suggest_2)
delta_mean_2 =delta_predict_2[0]
delta_variance_2=delta_predict_2[1]
ab_new_array = np.vstack((ab_new_array, ab_suggest_2)) # np.stack(())が二重になる点注意。
delta_new = np.asarray([[func_def(ab_suggest_2[0,0],ab_suggest_2[0,1],x_exp,y_exp)]])
delta_new_array = np.vstack((delta_new_array,delta_new)) # np.stack(())が二重になる点注意。
print(f'つぎは{ab_suggest_2}、deltaの予測平均は{delta_mean_2}、分散は{delta_variance_2},実際は{delta_new}')
myBopt_2.plot_acquisition()
Discussion