「ベイズ最適化」を読む　～第6章多目的ベイズ最適化 ②パレート最適解～

はじめに

データを使って仮説の生成と検証を行なうための方法であるベイズ最適化を学ぶために，今村秀明・松井孝太著「ベイズ最適化　ー適応的実験計画の基礎と実践ー」を読むことにした。

本記事は，「第6章多目的ベイズ最適化」のうち，パレート最適解に関する読書メモである。

本書の紹介ページ

www.kindaikagaku.co.jp

6.1 多目的最適化とは

関連記事(「ベイズ最適化」を読む～第6章多目的ベイズ最適化 ①多目的最適化とは～ - jiku log)において，弱優越とパレート最適解について説明があったが，直感的に理解するのが大変だったので，具体例で確認してみた。

弱優越の図示

あらためて，弱優越の定義を確認する。

定義(解の優越)
2点 $\mathbfit{x}, \mathbfit{x}' \in \mathcal{X}$ に対して，すべての $m=1, \cdots, M$ で $f^{(m)}(\mathbfit{x}) \leq f^{(m)}(\mathbfit{x}')$ が成り立つとき， $\mathbfit{f}_{\mathbfit{x}}$ は $\mathbfit{f}_{\mathbfit{x}'}$ は弱優越する(weakly dominate)といい，記号 $f^{(m)}(\mathbfit{x}) \preceq f^{(m)}(\mathbfit{x}')$ で表す。

弱優越の具体例を考える。 $M=2$ のとき，2つの点 $x, x'$ が与えられると， $(f^{(1)}(x), f^{(2)}(x))^T$ と $(f^{(1)}(x'), f^{(2)}(x'))^T$ の2つのベクトル値が得られる。

今回，目的関数を最小化したいので，「小さい方が優越している」とみなせる。「 $\mathbfit{f}_{x}$ が， $\mathbfit{f}_{x'}$ を弱優越している」という状態は，定義から

$\begin{align} \begin{cases} f^{(1)}(x) \leq f^{(1)}(x') \\ f^{(2)}(x) \leq f^{(2)}(x') \\ \end{cases} \end{align}$

という状態のことである。この状態を， $f^{(1)}-f^{(2)}$ 平面に図示すると以下のようになる。さらに「 $\mathbfit{f}_{x}$ が， $\mathbfit{f}_{x'}$ を弱優越している領域」も以下のように図示できる。

またこのことから，「 $\mathbfit{f}_{x^*}$ が， $\mathbfit{f}_{x}$ に弱優越されていない状態」というのは，以下のようになる。

パレート最適解の定義は，以下の通りであった。

定義(パレート最適解)
実行可能解 $\mathbfit{x}^* \in \mathcal{X}$ における関数値ベクトル $\mathbfit{f}_{\mathbfit{x}^*}$ が任意の $\mathbfit{x} \in \mathcal{X}$ に対して $\mathbfit{f}_{\mathbfit{x}}$ に弱優越されないとき， $\mathbfit{x}^*$ をパレート最適解 (Pareto optimal solution)と呼ぶ。

したがって， $x^*$ に対して，どんな $x \in \mathcal{X}$ を持ってきても図2のような状態になるとき， $x^*$ はパレート最適解になる。
逆に $x^*$ に対して，どれか1つ $x \in \mathcal{X}$ を持ってくると図1のような状態になるとき， $x^*$ はパレート最適解にはならない。

数値例

具体的な目的関数を持ってきて確認してみる。

$\begin{align} &f^{(1)}(x) = x^2 \\ &f^{(2)}(x) = (x-2)^2 \\ &-5 \leq x \leq 5 \end{align}$

とする。

$x$ を変化させて $(f^{(1)}(x), f^{(2)}(x))$ を計算し， $f^{(1)}-f^{(2)}$ 平面に図示すると以下のようになる。

上図の描画用コードは以下の通り。

クリックで展開

import matplotlib.pyplot as plt
import numpy as np

# 目的関数の定義 (例: 2つの目的関数)
def objective1(x):
  """目的関数1: xの2乗"""
  return x**2

def objective2(x):
  """目的関数2: (x-2)^2"""
  return (x-2)**2

# 解候補の生成 (例: -5から5まで0.1刻み)
x_values = np.arange(-5, 5, 0.1)

# 各目的関数の値を計算
obj1_values = [objective1(x) for x in x_values]
obj2_values = [objective2(x) for x in x_values]

# 弱優越解の判定と可視化
plt.figure(figsize=(8, 6))  # 図のサイズを設定
plt.xlabel("f^(1)")  # x軸ラベル
plt.ylabel("f^(2)")  # y軸ラベル
plt.title("Pareto Front Visualization")  # タイトル

# 全ての解をプロット
plt.scatter(obj1_values, obj2_values, label="All Solutions", alpha=0.5, color='gray')

plt.legend()  # 凡例を表示
plt.grid(True) # グリッドを表示
plt.show()  # プロットを表示

パレート最適解の定義にしたがってパレートフロントを計算すると以下のようになる。なおパレート最適解は， $0 \leq x \leq 2$ である。
この図から，パレート最適解は，解の集合のうち「左下」の部分になることが分かる。

上図の描画用コードは以下の通り。

クリックで展開

import matplotlib.pyplot as plt
import numpy as np

# 目的関数の定義 (例: 2つの目的関数)
def objective1(x):
  """目的関数1: xの2乗"""
  return x**2

def objective2(x):
  """目的関数2: (x-2)^2"""
  return (x-2)**2

# 解候補の生成 (例: -5から5まで0.1刻み)
x_values = np.arange(-5, 5, 0.1)

# 各目的関数の値を計算
obj1_values = [objective1(x) for x in x_values]
obj2_values = [objective2(x) for x in x_values]

# 全ての解をプロット
plt.figure(figsize=(8, 6))  # 図のサイズを設定
plt.xlabel("f^(1)")  # x軸ラベル
plt.ylabel("f^(2)")  # y軸ラベル
plt.title("Pareto Front Visualization")  # タイトル

plt.scatter(obj1_values, obj2_values, label="All Solutions", alpha=0.5, color='gray')

# 弱優越解の判定
pareto_front = []
for i in range(len(x_values)):
  is_dominated = False  # 弱優越されていないかどうかのフラグ
  for j in range(len(x_values)):
    if i != j:
      # 他の解と比較し、全てにおいて劣っている場合
      if obj1_values[i] >= obj1_values[j] and obj2_values[i] >= obj2_values[j] and (obj1_values[i] > obj1_values[j] or obj2_values[i] > obj2_values[j]):
        is_dominated = True  # 弱優越されている
        break
  if not is_dominated: # 弱優越されていない解はパレートフロンティアに属する
    pareto_front.append((obj1_values[i], obj2_values[i]))

# パレートフロンティアの点をプロット
pareto_front_x, pareto_front_y = zip(*pareto_front)
plt.scatter(pareto_front_x, pareto_front_y, label="Pareto Front", color='red') # パレートフロンティアを赤色でプロット

# 例としてx=0とx=2の解をプロットし、弱優越の関係を説明
#x = 0の場合
plt.scatter(objective1(0), objective2(0), label="x=0", color='blue',marker='x', s=100)
#x = 2の場合
plt.scatter(objective1(2), objective2(2), label="x=2", color='green', marker='*', s=100)

plt.legend()  # 凡例を表示
plt.grid(True) # グリッドを表示
plt.show()  # プロットを表示