⌨️

[論文読み] 混合整数計画法を用いたブラックボックス最適化

2022/10/16に公開約10,100字

ニューラルネットワークと MILP を用いたブラックボックス最適化 (NN+MILP) の論文読み。

元論文

Papalexopoulos, Theodore P., Christian Tjandraatmadja, Ross Anderson, Juan Pablo Vielma, and David Belanger. “Constrained Discrete Black-Box Optimization Using Mixed-Integer Programming.” In Proceedings of the 39th International Conference on Machine Learning, 17295–322. PMLR, 2022. https://proceedings.mlr.press/v162/papalexopoulos22a.html.

概要

モデルベース最適化 / Model-Based Optimization

ある関数 f を最大化 / 最小化したいような状況を考えます。f について、何らかの情報が得られる場合・・・例えば、「入出力関係は二次計画問題として書ける」とか、「勾配 df / dx が計算できる」とか言う場合には、それぞれの場合に適した計算を利用することができます。

しかし、f がブラックボックスであるような場合、つまり

  • 入出力関係を数式で書くことはできない
  • 勾配も計算できない
  • ...etc.

というような場合、少し異なるアプローチを取る必要があります。このような場合には、次のような モデルベース最適化 を用いることが考えられます。

Model-Based Optimization

  • 入力
    学習回数 N
    初期データセット \mathcal D_n = \{(x_1, y_1), (x_2, y_2), \cdots, (x_n, y_n)\}
    定義域 \mathcal X
  1. for t \in \{n+1, n+2, \cdots, N\}:
  2. \qquad \hat f \leftarrow \mathrm{fit} (\mathcal D_{t-1})
  3. \qquad a \leftarrow \mathrm{get \_ acquisition \_ function}(\hat f)
  4. \qquad x_t \leftarrow \arg \min_{x \in \mathcal X} a(x)
  5. \qquad y_t \leftarrow f(x_t)
  6. \qquad \mathcal D_t \leftarrow \mathcal D_{t-1} \cup \{(x_t, y_t)\}
  7. return \arg \min_{(x, y) \in \mathcal D_N} y
  • 出力
    探索されたデータ \mathcal D_N に含まれる最適解 (x, y)

正体不明の関数 f(x) を直接最適化するのではなく、手元に「代理関数」\hat f(x) を置いておき、それを逐次改良しながら最適化していくのがポイントです。各イテレーションにおいては、獲得関数 a(x) を最適化しています。

NN+MILP

論文で提案された方法 (NN+MILP) では、代理関数 \hat f として ReLU を用いたニューラルネットワークを、獲得関数 a(x) として MILP (混合整数線形計画問題) を使用しています。

Neural Network + Mixed-Integer Linear Programming

  • 入力
    学習回数 N
    初期データセット \mathcal D_n = \{(x_1, y_1), (x_2, y_2), \cdots, (x_n, y_n)\}
    定義域 \mathcal X
  1. for t \in \{n+1, n+2, \cdots, N\}:
  2. \qquad \hat f \leftarrow \mathrm{train \_ neural \_ network} (\mathcal D_{t-1})
  3. \qquad a \leftarrow \mathrm{build \_ MILP \_ formulation}(\hat f)
  4. \qquad x_t \leftarrow \mathrm{MILP \_ solver}(a)
  5. \qquad y_t \leftarrow f(x_t)
  6. \qquad \mathcal D_t \leftarrow \mathcal D_{t-1} \cup \{(x_t, y_t)\}
  7. return \arg \min_{(x, y) \in \mathcal D_N} y
  • 出力
    探索されたデータ \mathcal D_N に含まれる最適解 (x, y)

ただし、定義域 \mathcal X は離散で、しかも x = [x_1, x_2, \cdots, x_d] の各要素はすべて最大値と最小値を持っているとします。

ニューラルネットワーク

単一中間層のニューラルネットワークを考えます。数式で書くと次のようになります。

\begin{aligned} u &= W^{(1)} x + b^{(1)} \\ v &= \mathrm{ReLU} (u) \\ y &= W^{(2)} v + b^{(2)} \end{aligned}

ただし、d_1 を中間層の層数として、

\begin{aligned} x &\in \mathbb R^d \\ W^{(1)} &\in \mathbb R^{d_1 \times d} \\ W^{(2)} &\in \mathbb R^{1 \times d_1} \\ b^{(1)} &\in \mathbb R^{d_1} \\ b^{(1)} &\in \mathbb R \\ \end{aligned}

とします。

\mathrm{train \_ neural \_ network} (\mathcal D_{t-1}) では、\mathcal D_{t-1} を教師データとして、W^{(1)}, W^{(2)}, b^{(1)}, b^{(2)} を学習します。

MILP 定式化

ニューラルネットワークを MILP で扱えるように定式化します。モデルベース最適化の文脈に合わせて、以下のような制約条件を盛り込みます。

  1. 定義域 / Domain
  2. NG制約 / No-Good Constraints
  3. ニューラルネットワーク / Neural Network

目的関数は次のように設定します。

\text{minimize} \quad \hat f(x) = (\text{output of NN})

定義域 についてはかなり問題依存が大きいため、この記事では基本的に No-Good ConstraintsNeural Network について見ていきます。

No-Good Constraints

ノイズなしに対する制約

関数 f はノイズがないことを仮定します:

y = f(x)

つまり、ある点 x に対して f(x) を何度クエリしても、出てくる結果 y は変動しないということを仮定します。

ノイズありの場合には、出力値がランダムに変動することを考慮に入れてモデルを作成する必要があります。典型的には y = f(x) + \varepsilon として確率モデルを作成します。ノイズなしの場合、逆に、出力値が変動しないことを考慮に入れて、同じ点を2度クエリすることがないように 獲得関数を設計する必要が生じます。

  • ノイズありの場合:
    ノイズが生じることを考慮に入れる (y = f(x) + \varepsilon など)
  • ノイズなしの場合:
    同じ点 x に対して2回以上 f(x) を計算しないようにする

ハミング距離

No-Good Constraints は、ハミング距離 という概念を利用して同じ点を2回以上クエリしないようにしています。ここで、2つの文字列 \cal A\cal A のハミング距離とは、\cal A\cal A' の要素のうち異なるものの個数のことを指しています。つまり、例えば次のようなものです。

\begin{darray}{llllllllcr} \cal A &=& [1 & 2 & 3 & 4 & 5] \\ \cal A'&=& [1 & 2 & 3 & 4 & 5] &\rightarrow& 0\\ \cal A'&=& [1 & 2 & 3 & \color{#c03} 6 & 5] &\rightarrow& 1\\ \cal A'&=& [1 & \color{#c03} 4 & 3 & \color{#c03} 6 & 5] &\rightarrow& 2 \end{darray}

同じ点を2度クエリしないようにするには、すべての探索済みの点 \bar x \in \mathcal D_{t-1} と、次にクエリする点 x \in \mathcal X のハミング距離が 1 以上になれば良いはずです。すなわち、x\bar x のハミング距離を D_\mathrm{Ham}(x \| \bar x) と書くとき、

D_\mathrm{Ham}(x \| \bar x) \ge 1 \quad \forall \bar x \in \mathcal D_{t-1}

という制約を追加すればよいのです。

ハミング距離と二値変数

x = [x_1, x_2, \cdots, x_d]\bar x = [\bar x_1, \bar x_2, \cdots, \bar x_d] の要素がすべて \textrm{0-1} の二値変数である場合には、次のような式でハミング距離を計算することができます。

\begin{align*} D_\mathrm{Ham}(x \| \bar x) &= \sum_i \left| x_i - \bar x_i \right| \\ &= \sum_{i: \bar x_i = 0} x_i + \sum_{i: \bar x_i = 1} (1-x_i) \end{align*}

具体的に見てみましょう。例えば、

\begin{darray}{} \bar x &=& [0 & \textcolor{#039} 1 & \textcolor{#039} 1] \\ x &=& [1 & \textcolor{#c03} 0 & \textcolor{#c03} 1] \end{darray}

の場合、ハミング距離は 2 となるはずです。これを実際に上の式で計算してみると、

\begin{align*} D_\mathrm{Ham}(x \| \bar x) &= \sum_{i: \bar x_i = 0} x_i + \sum_{i: \bar x_i = 1} (1-x_i) \\ &= 1 + (\textcolor{#039}1-\textcolor{#c03}0) + (\textcolor{#039}1-\textcolor{#c03}1) \\ &= 2 \end{align*}

となり、たしかに計算結果は 2 となっています。

ハミング距離とワンホットエンコーディング

x\bar x が2値でない場合には、ワンホットエンコーディングを使用することで、同じような計算を行うことができます。

\begin{align*} z &= [z_{ij}] \\ \bar z &= [\bar z_{ij}] \\ z_{ij} &= \mathbb I (x_i=j) = \begin{cases} 1 & x_i = j \\ 0 & \text{otherwise} \end{cases} \\ \bar z_{ij} &= \mathbb I (\bar x_i=j) = \begin{cases} 1 & \bar x_i = j \\ 0 & \text{otherwise} \end{cases} \end{align*}

たとえば x \in \mathcal X については、

\begin{align*} x = \begin{bmatrix} 3 \\ 3 \\ 4\end{bmatrix} \quad \longrightarrow \quad z = \begin{bmatrix} 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 \end{bmatrix} \end{align*}

といった具合。\bar x \in \mathcal D_{t-1} も同様に、

\begin{align*} \bar x = \begin{bmatrix} 3 \\ 4 \\ 4\end{bmatrix} \quad \longrightarrow \quad \bar z = \begin{bmatrix} 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 \end{bmatrix} \end{align*}

とします。すると、ハミング距離は次式で計算できます。

D_\mathrm{Ham}(x \| \bar x) = \sum_{i, j: \bar z_{ij} = 0} z_{ij} + \sum_{i, j: \bar z_{ij} = 1} (1-z_{ij}) \\

ワンホットエンコーディングを MILP で

あとは MILP で変数 x \in \mathcal X のワンホットエンコーディングを表現できれば完成です。いろいろな方法が考えられると思いますが、シンプルには次の方法を利用できます。

\begin{align*} \sum_{j} j z_{ij} &= x_i \\ \sum_{j} z_{ij} &= 1 \end{align*}

既存の値である、定数 \bar x \in \mathcal D_{t-1} については、MILP ソルバ内でワンホットエンコーディングする必要はありません。予め次のようにワンホットエンコーディングして、これを \mathcal D_{t-1} にでも保存しておけば OK です。

\bar z_{ij} = \mathbb I \{\bar x_i = j\} = \begin{cases} 1 & \bar x_i = j \\ 0 & \text{otherwise} \end{cases} \quad \forall \bar x \in \mathcal D_{t-1}

まとめ

以上をまとめると、次のような制約条件を追加すればよいことになります。

\left| \begin{aligned} \sum_{i, j: \bar z_{ij} = 0} z_{ij} + \sum_{i, j: \bar z_{ij} = 1} (1-z_{ij}) &\ge 1 & \forall \bar z \in \mathcal D_{t-1} \\ \sum_{j} j z_{ij} &= x_i & \forall i \in [\mathrm{len}(x)] \\ \sum_{j} z_{ij} &= 1 & \forall i \in [\mathrm{len}(x)] \\ \\ z_{ij} &\in \{0, 1\} \end{aligned} \right.

Neural Network

\begin{aligned} u &= W^{(1)} x + b^{(1)} \\ v &= \mathrm{ReLU} (u) \\ y &= W^{(2)} v + b^{(2)} \end{aligned}

で表される単一中間層ニューラルネットワークを MILP で記述します。詳細は次の記事で。

https://zenn.dev/mory22k/articles/c2773eb131c1f7

全結合層 / Fully Connected Layer

全結合層はそのまま MILP の形で記述できます。

u = W^{(1)} x + b^{(1)}

正規化線形関数 / Rectifier

続いて ReLU 関数を MILP で記述します。

v = \mathrm{ReLU} (u)

まず、要素ごとに分解して記述します。

v_i = \max(0, u_i)

これを MILP で表現するために、十分大きな定数 M と、二値変数 \alpha \in \{0, 1\} を追加して、次のようにします。

\begin{align*} 0 &\le v_i \le M \alpha \\ u_i &\le v_i \le u_i + M(1-\alpha) \end{align*}

なお、この制約によって、二値変数 \alpha は次のような振る舞いをすることになります。

\alpha = \mathbb I \{u_i \ge 0\} = \begin{cases} 1 & u_i \ge 0 \\ 0 & u_i \le 0 \end{cases} >

出力層 / Output Layer

出力層は全結合層と同様の処理を行います。

y = W^{(2)} v + b^{(2)}

まとめ

まとめると、ニューラルネットワークは次のような制約条件によって記述されることになります。

\left| \begin{aligned} u & = W^{(1)} x + b^{(1)} \\ 0 &\le v_i \le M \alpha \\ u_i &\le v_i \le u_i + M(1-\alpha) & \forall i \in [\mathrm{len}(u)] \\ y &= W^{(2)} v + b^{(2)} \\ \\ u_i &\in \mathbb R \\ v_i &\in \mathbb R \\ \alpha &\in \{0, 1\} \end{aligned} \right.

MILP の全体像

MILP の全体像は、次のようになります。

\left| \begin{aligned} & \text{minimize} && y \in \mathbb R \\ \\ & \text{subject to} && \text{\underline{Domain}} \\ &&& \vdots \\ \\ &&& \text{\underline{No-Good Constraints}} \\ &&& \sum_{i, j: \bar z_{ij} = 0} z_{ij} + \sum_{i, j: \bar z_{ij} = 1} (1-z_{ij}) \ge 1 & \forall \bar z \in \mathcal D_{t-1} \\ &&& \sum_{j} j z_{ij} = x_i & \forall i \in [\mathrm{len}(x)] \\ &&& \sum_{j} z_{ij} = 1 & \forall i \in [\mathrm{len}(x)] \\ \\ &&& \text{\underline{Neural Network}} \\ &&& u = W^{(1)} x + b^{(1)} \\ &&& y = W^{(2)} v + b^{(2)} \\ &&& 0 \le v_i \le M \alpha & \forall i \in [\mathrm{len}(u)] \\ &&& u_i \le v_i \le u_i + M(1-\alpha) & \forall i \in [\mathrm{len}(u)] \\ \\ & \text{varibles:} && z_{ij} \in \{0, 1\} \\ &&& u_i \in \mathbb R \\ &&& v_i \in \mathbb R \\ &&& \alpha \in \{0, 1\} \\ \\ & \text{constants:} && \bar z: \text{one-hot encoding of }\bar x \in \mathcal D_{t-1} \\ &&& W^{(1)}, W^{(2)}, b^{(1)}, b^{(2)}: \text{params} \\ &&& M: \text{a sufficiently large number} \end{aligned} \right.

あとは、これをそのまんま MILP ソルバに入力して、最適化してもらいます。

まとめ

  • 離散ブラックボックス最適化を扱うよ
  • モデルベース最適化だよ
  • 関数はノイズなしだよ
  • 代理モデルはニューラルネットワークだよ
  • 獲得関数は代理モデルとおなじものだよ
    (Thompson Sampling にちょっとだけ似てる)
  • 獲得関数の最適化は MILP ソルバを使うよ

Discussion

ログインするとコメントできます