【Python】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と「ベイズ最適化を活用したスクリプト」を使って、複数の目標を同時に達成する「多目的最適化」の考え方と、その中心にある「ベイズ最適化」の賢さをご紹介しました。

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

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

↓例(AIコーディングにより作成したベイズ最適化アプリ)
この内容が理解できるようになりたいですね。

import streamlit as st
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from skopt import gp_minimize
from skopt.space import Real
from skopt.plots import plot_convergence, plot_objective, plot_evaluations
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.metrics import r2_score
import warnings

warnings.filterwarnings("ignore", category=FutureWarning)

# Set a default font for Matplotlib to prevent garbled text
plt.rcParams['font.family'] = 'Verdana'  # Adjust as necessary for your environment

# --- Session State Initialization ---
if 'optimization_result' not in st.session_state:
    st.session_state.optimization_result = None
if 'last_settings' not in st.session_state:
    st.session_state.last_settings = {}

# --- 1. Surrogate Model Training ---
def train_surrogate_model(df: pd.DataFrame, input_columns: list[str], output_columns: list[str], output_weights: dict):
    """
    Trains a surrogate model using a weighted average of the output columns.
    """
    X = df[input_columns]
    # Calculate weighted average of outputs
    y = np.average(df[output_columns], axis=1, weights=[output_weights[col] for col in output_columns])

    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)
    
    model = RandomForestRegressor(n_estimators=100, random_state=42)
    model.fit(X_train, y_train)

    score = r2_score(y_test, model.predict(X_test))
    
    final_model = RandomForestRegressor(n_estimators=100, random_state=42)
    final_model.fit(X, y)

    return final_model, score

# --- 2. Objective Function Creation ---
def create_objective_function(model, is_maximize=False, **kwargs):
    """
    Creates the objective function for Bayesian optimization.
    """
    def objective(X: list[float]) -> float:
        X_array = np.array(X).reshape(1, -1)
        pred_value = model.predict(X_array)[0]
        return -pred_value if is_maximize else pred_value

    return objective

# --- 3. Bayesian Optimization Execution ---
def run_bayesian_optimization(dimensions, objective_func, n_calls, base_estimator, acq_func):
    """
    Wrapper to run the Bayesian optimization.
    """
    return gp_minimize(
        func=objective_func,
        dimensions=dimensions,
        acq_func=acq_func,
        base_estimator=base_estimator,
        n_calls=n_calls,
        n_initial_points=10,
        random_state=0
    )

# --- 4. Streamlit UI Helper Functions ---
def setup_sidebar():
    """
    Configures the Streamlit sidebar.
    """
    st.sidebar.header("設定")
    settings = {}
    
    uploaded_file = st.sidebar.file_uploader("CSVファイルをアップロード", type=["csv"])
    settings['uploaded_file'] = uploaded_file
    
    data, input_columns, output_columns, dimensions = None, [], [], []
    output_weights = {}

    if uploaded_file:
        data = pd.read_csv(uploaded_file)
        st.sidebar.subheader("CSVデータプレビュー")
        st.sidebar.write(data)
        all_columns = data.columns.tolist()

        st.sidebar.subheader("入力・出力列の選択")
        default_inputs = [col for col in all_columns if 'input' in col.lower() or 'x' in col.lower()] or all_columns[:1]
        default_outputs = [col for col in all_columns if 'output' in col.lower() or 'y' in col.lower() or 'target' in col.lower()] or all_columns[1:2]

        input_columns = st.sidebar.multiselect("入力変数を選択", all_columns, default=default_inputs)
        output_columns = st.sidebar.multiselect("出力変数を選択 (最大5つ)", all_columns, default=default_outputs, max_selections=5)

        if len(output_columns) > 1:
            st.sidebar.subheader("出力変数の重み設定")
            total_weight = 0
            for col in output_columns:
                weight = st.sidebar.slider(f"`{col}` の重み", 0.0, 1.0, 1.0 / len(output_columns), 0.01)
                output_weights[col] = weight
                total_weight += weight
            if not np.isclose(total_weight, 1.0):
                st.sidebar.warning(f"重みの合計が1.0になるように調整してください。(現在: {total_weight:.2f})")
        elif len(output_columns) == 1:
             output_weights[output_columns[0]] = 1.0

        if input_columns:
            st.sidebar.subheader("入力変数の範囲を調整")
            for col in input_columns:
                st.sidebar.write(f"**{col}**")
                cols = st.sidebar.columns(2)
                min_val_csv = float(data[col].min())
                max_val_csv = float(data[col].max())
                
                min_val = cols[0].number_input(f"最小値", value=min_val_csv, step=0.001, format="%.3f", key=f"min_{col}")
                max_val = cols[1].number_input(f"最大値", value=max_val_csv, step=0.001, format="%.3f", key=f"max_{col}")
                
                dimensions.append(Real(min_val, max_val, name=col))
    else:
        st.sidebar.info("最適化を開始するには、CSVファイルをアップロードしてください。")
            
    settings.update({
        'data': data, 'input_columns': input_columns, 'output_columns': output_columns, 
        'output_weights': output_weights, 'dimensions': dimensions
    })

    st.sidebar.subheader("最適化設定")
    settings['is_maximize'] = st.sidebar.radio("最適化の方向", ("最大化", "最小化")) == "最大化"
    settings['n_calls'] = st.sidebar.slider("評価回数 (n_calls)", 10, 200, 50)

    st.sidebar.subheader("高度な設定")
    settings['base_estimator'] = st.sidebar.selectbox("サロゲートモデル", ("GP", "RF", "GBRT", "ET"))
    settings['acq_func'] = st.sidebar.selectbox("獲得関数", ("LCB", "EI", "PI", "gp_hedge"))
    
    return settings

def validate_settings(settings):
    if not settings['uploaded_file']:
        st.error("エラー: 最適化にはCSVファイルのアップロードが必要です。")
        return False
    if not settings['input_columns'] or not settings['output_columns']:
        st.error("エラー: 入力変数と出力変数を1つ以上選択してください。")
        return False
    if not settings['dimensions']:
        st.error("エラー: 入力変数が設定されていません。")
        return False
    if any(dim.low >= dim.high for dim in settings['dimensions']):
        st.error("エラー: 最小値は最大値より小さい必要があります。")
        return False
    if len(settings['output_columns']) > 1 and not np.isclose(sum(settings['output_weights'].values()), 1.0):
        st.error("エラー: 出力変数の重みの合計が1.0になるようにしてください。")
        return False
    return True

def display_results(result, settings):
    st.subheader("最適化結果")
    is_maximize = settings['is_maximize']
    best_val = -result.fun if is_maximize else result.fun
    st.metric(f"最適値 ({'最大化' if is_maximize else '最小化'})", f"{best_val:.4f}")
    
    param_names = [d.name for d in settings['dimensions']]
    st.dataframe(pd.DataFrame([result.x], columns=param_names))

    with st.expander("詳細な最適化プロセスの可視化"):
        # Plot 1: Convergence
        st.subheader("目的関数値の収束履歴")
        fig, ax = plt.subplots(figsize=(10, 4))
        plot_convergence(result, ax=ax)
        ax.set_title("Convergence Plot")
        st.pyplot(fig)

        # Plot 2: Objective Landscape
        if len(param_names) <= 2:
            st.subheader("目的関数と探索点の関係")
            fig, ax = plt.subplots(figsize=(10, 6))
            plot_objective(result, ax=ax, dimensions=param_names, n_points=50)
            st.pyplot(fig)

        # Plot 3: Evaluations
        st.subheader("入力変数の探索履歴")
        fig, ax = plt.subplots(figsize=(10, 6))
        plot_evaluations(result, ax=ax, dimensions=param_names)
        st.pyplot(fig)

        # Download Link
        st.subheader("最適化履歴のダウンロード")
        history_df = pd.DataFrame(result.x_iters, columns=param_names)
        history_df['objective_value'] = -np.array(result.func_vals) if is_maximize else result.func_vals
        csv = history_df.to_csv(index=False).encode('utf-8')
        st.download_button("履歴をCSVでダウンロード", csv, 'optimization_history.csv', 'text/csv')

# --- Main Application Logic ---
def main():
    st.set_page_config(layout="wide")
    st.title("ベイズ最適化 Professional")
    
    settings = setup_sidebar()

    if st.sidebar.button("最適化開始", type="primary"):
        if validate_settings(settings):
            with st.spinner("最適化を実行中..."):
                try:
                    model, score = train_surrogate_model(settings['data'], settings['input_columns'], settings['output_columns'], settings['output_weights'])
                    st.success(f"サロゲートモデル訓練完了 (R²スコア: {score:.3f})")
                    objective_func = create_objective_function(model=model, is_maximize=settings['is_maximize'])
                    
                    if settings['base_estimator'] == 'Dummy':
                        base_estimator = DummyRegressor(strategy='mean')
                    else:
                        base_estimator = settings['base_estimator']

                    result = run_bayesian_optimization(
                        settings['dimensions'], objective_func, settings['n_calls'],
                        base_estimator,
                        settings['acq_func']
                    )
                    st.session_state.optimization_result = result
                    st.session_state.last_settings = settings

                except Exception as e:
                    st.error(f"エラーが発生しました: {e}")
    
    if st.session_state.optimization_result:
        display_results(st.session_state.optimization_result, st.session_state.last_settings)

if __name__ == "__main__":
    main()

コメント