Pythonで「賢く」複数の目標を叶える!ベイズ最適化の魔法

Python

はじめに:欲張りな願いを叶える「賢い」方法

材料開発やプロセス開発の現場では、「強度を上げつつ、コストは下げたい」「性能を最大にしつつ、環境負荷は最小にしたい」といった、複数の目標を同時に達成したい場面がよくありますよね。しかし、これらの目標は互いに反発し合うことが多く、闇雲に試行錯誤しても、なかなか良い答えは見つかりません。

そんな「欲張りな願い」を効率的に叶えるための強力なツールが「ベイズ最適化」です。本記事では、この賢い最適化手法が、どのように複数の目標を同時に追いかける「多目的最適化」に役立つのかを、分かりやすく解説していきます。

多目的最適化って、どんなこと?

まず、多目的最適化とは、その名の通り、複数の「目的」(ゴール)を同時に達成しようとする問題のことです。

例えば、あなたが新しい素材を開発するとします。 目的1:素材の「硬さ」を最大にしたい 目的2:素材の「軽さ」を最大にしたい

この2つの目的は、一般的に「トレードオフ」の関係にあります。つまり、硬くしようとすると重くなりがちで、軽くしようとすると硬さが犠牲になる、といった具合です。

このように、複数の目的が互いに相反する場合、すべての目的を同時に最高にする「たった一つの完璧な答え」は、ほとんど存在しません。その代わりに、私たちは「パレート最適解」と呼ばれる特別な解のグループを探します。これは、「これ以上、どれか一つの目的を良くしようとすると、必ず他の目的のどれかが悪くなってしまう」という、バランスの取れた解の集まりのことです。このグループの中から、あなたの状況に一番合った「最高のバランス」を見つけるのが、多目的最適化のゴールなんです。

少ない試行で最適な答えを見つける「ベイズ最適化」の賢さ

では、どうやってこの「パレート最適解」を効率的に見つけるのでしょうか?ここで登場するのが「ベイズ最適化」です。

ベイズ最適化は、特に「試行回数をできるだけ少なくしたい」という場合に絶大な威力を発揮する最適化手法です。例えば、新しい材料の実験には時間もコストもかかりますよね。そんな時、手当たり次第に試すのではなく、ベイズ最適化は「次にどこを試せば、一番良い結果が得られそうか?」を賢く予測してくれます。

その賢さの秘密は、これまでの実験データから「目的関数がどんな形をしているか」を統計的に学習し、まだ試していない場所で「どれくらい良い結果が出そうか」と「どれくらい新しい情報が得られそうか」のバランスを考えて、次の実験点を決める点にあります。まるで、経験豊富な研究者が次にどこを試すべきか直感的に判断するようなイメージです。

ベイズ最適化が多目的最適化を「効率化」する仕組み

このベイズ最適化の「賢さ」は、複数の目標を同時に追いかける多目的最適化において、特に力を発揮します。

ベイズ最適化は、それぞれの目的関数について、これまでのデータから「次に試すべき最適な条件」を提案してくれます。これにより、手探りで実験条件を探すよりも、はるかに少ない試行回数で、複数の目標をバランス良く達成できる「パレート最適解」の候補を効率的に見つけ出すことができるのです。

本記事でご紹介する「ベイズ最適化を活用したスクリプト」は、まさにこのベイズ最適化の考え方を使って、複数の目的を同時に最適化するためのツールです。このスクリプトを使うことで、複雑な多目的最適化の問題も、より少ない労力で解決に導くことが可能になります。

さあ、実践!Pythonでベイズ最適化を動かしてみよう

それでは、実際にPythonと「ベイズ最適化を活用したスクリプト」を使って、最適化を進める手順を見ていきましょう。ここでは、まずベイズ最適化の基本的な動きを理解するために、シンプルな単一目的の最適化の例を示します。

準備:必要なライブラリのインストール

ベイズ最適化には、scikit-optimize (skopt) というライブラリが便利です。まだインストールしていない場合は、以下のコマンドでインストールしておきましょう。

pip install scikit-optimize matplotlib

コード例:シンプルなベイズ最適化

以下のPythonコードは、ある関数(目的関数)の最小値を効率的に見つけ出すベイズ最適化の基本的な流れを示しています。

import numpy as np
from skopt import gp_minimize
from skopt.space import Real
import matplotlib.pyplot as plt
from skopt.plots import plot_convergence

# 目的関数の定義(例:最小化したい関数)
# この関数が、例えば「材料の硬さ」や「製品のコスト」など、
# 最適化したい対象の性能を表すと考えられます。
# 実際には、実験結果やシミュレーション結果を返す関数になります。
def objective_function(x):
    # xはリスト形式で渡されるので、最初の要素を使います
    # この関数は、xが0.5付近で最小値を取るような形をしています
    return (x[0] - 0.5)**2 + np.sin(x[0] * 5) * 0.1

# 探索空間の定義
# 今回は、入力変数 'x' が0から1の範囲で探索するように設定します。
# 複数の入力変数がある場合は、ここにRealを複数追加します。
dimensions = [Real(0.0, 1.0, name='x')]

# ベイズ最適化の実行
# gp_minimize: ガウス過程を用いたベイズ最適化を実行する関数
# n_calls: 目的関数を評価する(実験を行う)回数。少ない試行回数で最適解を探します。
# random_state: 乱数シード。結果の再現性を確保するために設定します。
res = gp_minimize(objective_function,  # 最適化する関数
                  dimensions,          # 探索空間
                  n_calls=20,          # 試行回数
                  random_state=0)      # 乱数シード

# 最適化結果の表示
print(f"最適化されたxの値: {res.x[0]:.4f}")
print(f"目的関数の最小値: {res.fun:.4f}")

# 収束プロットの表示
# ベイズ最適化がどのように最適な値に近づいていったかを示します。
# 試行回数を重ねるごとに、目的関数の最小値が改善されていく様子が分かります。
plot_convergence(res)
plt.title("ベイズ最適化の収束過程")
plt.xlabel("試行回数")
plt.ylabel("目的関数の最小値")
plt.show()

# 目的関数の形状と探索点のプロット
# 目的関数がどのような形をしていて、ベイズ最適化がどの点を探索したかを示します。
# 赤い点が探索した点、緑の星が最終的に見つけ出した最適解です。
x_plot = np.linspace(0, 1, 100)
y_plot = [objective_function([val]) for val in x_plot]

plt.figure(figsize=(8, 6))
plt.plot(x_plot, y_plot, label="目的関数")
plt.scatter(res.x_iters, res.func_vals, color='red', label="探索点")
plt.scatter(res.x[0], res.fun, color='green', marker='*', s=200, label="最適解")
plt.title("目的関数とベイズ最適化の探索点")
plt.xlabel("x")
plt.ylabel("目的関数の値")
plt.legend()
plt.grid(True)
plt.show()

多目的最適化への応用

上記のコードは単一の目的関数を最適化する例ですが、これを複数の目的関数に応用することで、多目的最適化が可能になります。

基本的な考え方は以下の通りです。

  1. 複数の目的関数を定義: 最適化したいそれぞれの目標(硬さ、軽さ、コストなど)に対応する目的関数を用意します。
  2. ベイズ最適化の適用: 各目的関数に対して、上記のコードのようにベイズ最適化を適用し、それぞれの目的関数にとっての最適な条件を探索します。
  3. パレート最適解の導出: 各目的関数で得られた探索結果を統合し、複数の目的をバランス良く満たす「パレート最適解」の集合を導出します。この際、各目的関数の重要度に応じて重み付けを行ったり、制約条件を設けたりすることもあります。

このように、ベイズ最適化の「賢く探索する」能力を複数の目的に適用することで、複雑な多目的最適化の問題も、効率的に解決へと導くことができるのです。

まとめ:あなたの研究・開発を「賢く」加速させる

本記事では、Pythonと「ベイズ最適化を活用したスクリプト」を使って、複数の目標を同時に達成する「多目的最適化」の考え方と、その中心にある「ベイズ最適化」の賢さをご紹介しました。

相反する目標を同時に追いかけるのは、一見すると複雑で難しい課題に思えるかもしれません。しかし、ベイズ最適化のような強力なツールを活用することで、効率的に、そして賢く、あなたの研究や開発における「最高のバランス」を見つけ出すことが可能です。

ぜひ、このベイズ最適化の技術を、皆さんの日々の研究や開発に取り入れて、新たな発見やブレークスルーを「賢く」生み出してくださいね!

↓5変数で1つの値を最適化したいときの例(AIコーディングにより作成)
この内容が理解できるようになりたいですね。

import pandas as pd
import numpy as np
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel as C, WhiteKernel
from sklearn.preprocessing import StandardScaler
from skopt import gp_minimize
from skopt.space import Real
from skopt.utils import use_named_args
import matplotlib.pyplot as plt
import warnings
import tkinter as tk
from tkinter import filedialog
import argparse
import os
from datetime import datetime

def main():
    # --- 1. コマンドライン引数の設定 ---
    parser = argparse.ArgumentParser(description='CSVデータに基づいて単一目的最適化(特性Aの最大化)を実行します。')
    parser.add_argument(
        '--n_calls', type=int, default=50, 
        help='ベイズ最適化の試行回数。デフォルトは50(少データ対応)。')
    parser.add_argument(
        '--n_random_starts', type=int, default=10, 
        help='初期ランダム探索の回数。デフォルトは10(少データ対応)。')
    args = parser.parse_args()

    # 警告を非表示にする設定
    warnings.filterwarnings("ignore")

    # --- 2. CSVファイル選択ダイアログの表示 ---
    root = tk.Tk()
    root.withdraw()

    print("最適化に使用するCSVファイルを選択してください...")
    filepath = filedialog.askopenfilename(
        title="最適化の元データとなるCSVファイルを選択",
        filetypes=[("CSVファイル", "*.csv"), ("すべてのファイル", "*.*")]
    )

    if not filepath:
        print("ファイルが選択されませんでした。プログラムを終了します。")
        exit()

    # --- 3. データの読み込みと前処理 ---
    try:
        print(f"'{filepath}' を読み込んでいます...")
        df = pd.read_csv(filepath)
    except Exception as e:
        print(f"ファイルの読み込み中にエラーが発生しました: {e}")
        exit()

    # 必須カラムの定義
    required_columns = ['A', 'B', 'C', 'D', 'E', 'F']
    if not all(col in df.columns for col in required_columns):
        print(f"エラー: CSVファイルには、次のカラムが必要です: {required_columns}")
        exit()

    # データの確認
    print(f"読み込んだデータ点数: {len(df)}")
    if len(df) < 3:
        print("エラー: 最低3つのデータ点が必要です。")
        exit()

    # 変数(X)と特性(y_A)にデータを分割
    X = df[['A', 'B', 'C', 'D', 'E']].values
    y_A = df['F'].values

    # データの統計情報を表示
    print("\n--- データの統計情報 ---")
    print(f"特性F (目的変数) の範囲: {y_A.min():.4f} ~ {y_A.max():.4f}")
    print("各変数の範囲:")
    for i, col in enumerate(['A', 'B', 'C', 'D', 'E']):
        print(f"  {col}: {X[:, i].min():.4f} ~ {X[:, i].max():.4f}")

    # --- 4. データの正規化 ---
    # 少ないデータでも安定して学習できるよう、データを正規化します
    scaler_X = StandardScaler()
    scaler_y = StandardScaler()
    
    X_scaled = scaler_X.fit_transform(X)
    y_A_scaled = scaler_y.fit_transform(y_A.reshape(-1, 1)).flatten()

    # --- 5. 代理モデルの構築 (少データ対応) ---
    # 少ないデータに適したカーネルとパラメータを使用
    kernel = (
        C(1.0, (1e-3, 1e3)) * 
        RBF(length_scale=1.0, length_scale_bounds=(1e-2, 1e2)) + 
        WhiteKernel(noise_level=1e-3, noise_level_bounds=(1e-5, 1e-1))
    )
    
    model_A = GaussianProcessRegressor(
        kernel=kernel, 
        n_restarts_optimizer=20,  # 少データでより多くの再試行
        random_state=42,
        alpha=1e-6,  # 数値安定性のため
        normalize_y=False  # 既に正規化済み
    )
    
    model_A.fit(X_scaled, y_A_scaled)
    print("CSVから読み込んだデータで、特性Aの代理モデルを構築しました。")

    # --- 6. 最適化のための問題設定 ---
    # 実際のデータの範囲に基づいて探索範囲を設定
    margin = 0.1  # 既存データ範囲の10%のマージンを追加
    dimensions = []
    for i, col in enumerate(['A', 'B', 'C', 'D', 'E']):
        min_val = X[:, i].min()
        max_val = X[:, i].max()
        range_val = max_val - min_val
        if range_val == 0:  # 全て同じ値の場合
            range_val = 1.0
        
        low = min_val - margin * range_val
        high = max_val + margin * range_val
        dimensions.append(Real(low=low, high=high, name=col))

    print("\n--- 最適化の探索範囲 ---")
    for dim in dimensions:
        print(f"  {dim.name}: {dim.low:.4f} ~ {dim.high:.4f}")

    # --- 7. 目的関数の定義 (正規化対応) ---
    @use_named_args(dimensions=dimensions)
    def objective(**params):
        x_pred = np.array(list(params.values())).reshape(1, -1)
        
        # 正規化
        x_pred_scaled = scaler_X.transform(x_pred)
        
        # 予測
        pred_A_scaled = model_A.predict(x_pred_scaled)[0]
        
        # 逆正規化
        pred_A = scaler_y.inverse_transform([[pred_A_scaled]])[0, 0]
        
        return -pred_A  # 最大化のため負の値を返す

    # --- 8. ベイズ最適化の実行 ---
    print("\n代理モデルに基づいたベイズ最適化を開始します...")
    
    # 少データの場合は試行回数を調整
    n_calls = min(args.n_calls, len(df) * 10)  # データ数の10倍を上限
    n_random_starts = min(args.n_random_starts, max(5, len(df)))  # 最低5回、データ数まで
    
    print(f"試行回数: {n_calls}, 初期ランダム探索: {n_random_starts}")
    
    result = gp_minimize(
        func=objective,
        dimensions=dimensions,
        n_calls=n_calls,
        random_state=42,
        n_random_starts=n_random_starts,
        acq_func='EI',  # Expected Improvement
        acq_optimizer='auto'
    )
    print("最適化が完了しました。")

    # --- 9. 結果の表示 ---
    best_parameters = result.x
    final_x_pred = np.array(best_parameters).reshape(1, -1)
    final_x_pred_scaled = scaler_X.transform(final_x_pred)
    final_A_pred_scaled = model_A.predict(final_x_pred_scaled)[0]
    final_A_pred = scaler_y.inverse_transform([[final_A_pred_scaled]])[0, 0]

    print("\n--- 最適化結果 ---")
    print("最適な変数の組み合わせ (予測):")
    for dim, value in zip(dimensions, best_parameters):
        print(f"  {dim.name}: {value:.4f}")
    print(f"\nそのときの特性F (最大化) の予測値: {final_A_pred:.4f}")
    
    # 予測の不確実性も表示
    pred_mean, pred_std = model_A.predict(final_x_pred_scaled, return_std=True)
    pred_std_original = pred_std[0] * scaler_y.scale_[0]  # 元のスケールに戻す
    print(f"予測の標準偏差: ±{pred_std_original:.4f}")
    print(f"95%信頼区間: [{final_A_pred - 1.96*pred_std_original:.4f}, {final_A_pred + 1.96*pred_std_original:.4f}]")

    # --- 10. 結果のCSV保存 ---
    output_dir = "optimization_results"
    os.makedirs(output_dir, exist_ok=True)
    timestamp = datetime.now().strftime("%Y%m%d_%H%M%S")
    output_csv_path = os.path.join(output_dir, f"optimal_parameters_{timestamp}.csv")

    result_df = pd.DataFrame({
        dim.name: [value] for dim, value in zip(dimensions, best_parameters)
    })
    result_df['Predicted_Trait_F'] = final_A_pred
    result_df['Prediction_Std'] = pred_std_original
    result_df['Confidence_Interval_Lower'] = final_A_pred - 1.96*pred_std_original
    result_df['Confidence_Interval_Upper'] = final_A_pred + 1.96*pred_std_original
    result_df.to_csv(output_csv_path, index=False)
    print(f"\n最適化結果を '{output_csv_path}' に保存しました。")

    # --- 11. 探索過程の可視化と画像保存 ---
    plt.figure(figsize=(12, 8))
    
    # サブプロット1: 目的関数の推移
    plt.subplot(2, 1, 1)
    plt.plot(result.func_vals, marker='o', linestyle='-', color='blue', alpha=0.7)
    plt.title('Objective Function Value (Negative of Trait F) over Iterations')
    plt.xlabel('Iteration')
    plt.ylabel('Objective Value (Lower is Better)')
    plt.grid(True, alpha=0.3)
    
    # サブプロット2: 累積最良値
    plt.subplot(2, 1, 2)
    best_so_far = np.minimum.accumulate(result.func_vals)
    plt.plot(best_so_far, marker='s', linestyle='-', color='red', alpha=0.7)
    plt.title('Best Value Found So Far')
    plt.xlabel('Iteration')
    plt.ylabel('Best Objective Value')
    plt.grid(True, alpha=0.3)
    
    plt.tight_layout()
    
    output_plot_path = os.path.join(output_dir, f"optimization_plot_{timestamp}.png")
    plt.savefig(output_plot_path, dpi=300, bbox_inches='tight')
    print(f"探索過程のグラフを '{output_plot_path}' に保存しました。")
    plt.show()

    print("\nグラフを生成しました。縦軸の値が小さいほど、特性Fが大きいことを示します。")
    print(f"元データの特性F最大値: {y_A.max():.4f}")
    print(f"予測された最適値: {final_A_pred:.4f}")

    # モデルの性能評価(簡易版)
    print("\n--- モデル性能の簡易評価 ---")
    try:
        # 学習データに対する予測
        X_scaled_all = scaler_X.transform(X)
        y_pred_scaled = model_A.predict(X_scaled_all)
        y_pred = scaler_y.inverse_transform(y_pred_scaled.reshape(-1, 1)).flatten()
        
        # R²スコア
        ss_res = np.sum((y_A - y_pred) ** 2)
        ss_tot = np.sum((y_A - np.mean(y_A)) ** 2)
        r2_score = 1 - (ss_res / ss_tot)
        
        print(f"学習データに対するR²スコア: {r2_score:.4f}")
        print(f"平均絶対誤差: {np.mean(np.abs(y_A - y_pred)):.4f}")
        
    except Exception as e:
        print(f"モデル性能評価中にエラーが発生しました: {e}")

if __name__ == '__main__':
    main()

コメント