Deep Learning界隈でよく目にするReLU関数をQUBO(Quadratic Unconstrained Binary Optimization)で表現する方法についてレビューします。
https://arxiv.org/abs/1811.03829
この手法は上の論文で提案されています。
数学的に間違った表現や計算をしている場合はお知らせ下さい。
(おそらく色々と誤植があると思うので適宜修正中です。)
概要
論文上でReLU関数は次の様に定義されています。
\begin{align}
f(m) &= -\min_{m}\{0,m\}\\
&= \begin{cases}
-m & (m <0) \\
0 & (m \ge 0)
\end{cases}
\end{align}
よく見かけるReLUとは定義が逆な感じがしますが、右肩上がりにしたい場合は$m\to-m$に変換すれば大丈夫です。
しかし、後で説明しますがこの変換を最終的なQUBOの式に適用しても所望の式を得ることができません。なので$m\to -m$に変換した元の式からスタートして一連の操作を行ってください。
QUBO化の流れとしては以下の通りとなります。
- ルジャンドル変換を施す
- 逆ルジャンドル変換を施す
- Wolf変換を施す
- ラグランジュの未定乗数法を用いる
- 補助変数をバイナリ変数で展開する
まずは簡単にルジャンドル変換やWolf変換について説明した後に本題に入ります。
先行研究として紹介されているこちらの論文(https://arxiv.org/abs/1205.1148) でもルジャンドル変換、逆ルジャンドル変換までは行っています。ここまでの操作を行った後のReLU関数は
\begin{align}
F(m) = - \min_{t} \{-mt \}
\end{align}
となります。$t$はルジャンドル変換で導入された新たな変数です。
一緒に最適化したいコスト関数を$C(m)$とすると、$F(m)$に入っている$\min$の前の負符号のせいで、
\begin{align}
\min_{m}\{C(m) + F(m) \} &= \min_{m}\biggl\{ C(m) - \min_{t}\{ -mt\}\biggr\} \\
&\neq \min_{m,t}\{C(m) - (-mt) \}
\end{align}
となり、$m,t$に対して同時に最適化することが出来ません。
提案手法ではWolf変換を施してこれを回避しています。
このWolf変換とやらがポイントのようです。
数学的準備
証明で使われている数学について説明します。
冗長なところもあるので知っていたら飛ばして下さい。
ルジャンドル変換
ルジャンドル変換は新たな変数を導入して非凸関数を凸関数として表す変換です。
下で定義される$f^*(t)$を$f(x)$のルジャンドル変換と呼ぶ。
\begin{align}
f^*(t) = \sup_{x}\{xt - f(x)\} \\
f^*(t) = -\inf_{x}\{f(x) -xt\}
\end{align}
実際に変換するときは$f(x)$が微分可能な区間と不可能な区間に分けて使用します。
これにより$t$の値域が決定されます。
ReLUは原点で微分不可能なので、原点以外と原点でそれぞれ場合分けが必要になってきます。
最後にマージしますが。
(1) $f(x)$が微分可能な区間を含む場合
$f^*(t)$は$f(x) - xt$をxについて最小にすることにより得られる関数なので
\begin{align}
\frac{\partial f^*}{\partial x} = 0 \Leftrightarrow t = \frac{\partial f}{\partial x}
\end{align}
$f(x)$が凸関数であれば$f^{\prime}(x)$は単調増加なので、上式の逆関数$\psi$が存在する。
\begin{align}
x &= \psi(t) \\
f^{*}(t) &= - \biggl(f(\psi(t)) - t\psi(t) \biggr)
\end{align}
が求まる。
$\frac{\partial f^*}{\partial x}$が定数となる場合は$t$の値を1つ選んでも$x$は一意に決まらない。
この場合、逆ルジャンドル変換をした後の$x$の範囲は実数全体になる。
また、
\begin{align}
df = \biggl(\frac{\partial f}{\partial x} \biggr)dx = t dx
\end{align}
より、
\begin{align}
df^* = df - t dx - x dt = -x dt
\end{align}
となり、ルジャンドル変換により独立変数が$x$から$t$に入れ替わったことが確認できる。
(2) $f(x)$が微分不可能な区間を含む場合
劣微分$\partial f(x)= [f^{\prime}_{-}(x), f^{\prime}_{+}(x)]$が動く範囲の全体を$C^{\ast}$とすると、$C^{\ast}$は$f(x)$の接線の傾きとして許される値全ての集まりで、$\mathbb{R}$の中の単一の区間になる。$u\in C^{\ast}$についてルジャンドル変換の最小値が存在する。つまり$C^{\ast}$上の関数$f^{\ast}(t)$が定義される。
これより、
\begin{align}
t_{\pm} = \frac{\partial f}{\partial x}\bigg|_{x \to x_{\pm}}
\end{align}
を用いて$t$の値域が$t\in [t_{-}, t_{+} ]$となる。
不連続な関数をルジャンドル変換すると、不連続点で$f^{\ast}(t)$の値が一意に決まらない。
QUBO化する時にはルジャンドル変換する関数には連続性を仮定した方が良さそう。
- 参考文献
http://www.phys.cs.is.nagoya-u.ac.jp/~tanimura/class/H28-tanimura-mechanics-6.pdf
https://ocw.u-tokyo.ac.jp/lecture_files/gf_06/7/notes/ja/07murota.pdf
田崎さんの統計力学2の付録
Wolf変換
正しくはWolf Duality Theoremと呼ぶらしいです。
原論文はこちらになります。
P.Wolfe, "A DUALITY THEOREM FOR NON-LINEAR PROGRAMMING." Quarterly ofApplied Mathematics 19, no. 3 (1961): 239-44. http://www.jstor.org/stable/43635235
Wolf変換を行う事で不等式制約条件のある最小化問題を最大化問題にマップすることができる。
元問題を次のように定義する。
\begin{align}
\min_x\{ f(x) \}\\
{\rm subject\;to\;} \boldsymbol{g}(x) \le 0
\end{align}
$\boldsymbol{g}(x)$が不等式制約条件の集合です。
Wolf変換による元問題は次の双対問題にマップされる。
\begin{align}
\max_{x,\boldsymbol{z}}&\{ \mathcal{L}(x,\boldsymbol{z}) \}\\
{\rm subject\;to\;} &\nabla_x \mathcal{L}(x,\boldsymbol{z})= 0,\; \boldsymbol{g}\le 0, \; \boldsymbol{z} \ge 0 \\
\mathcal{L}(x,\boldsymbol{z}) &= f(x) + \boldsymbol{z} \cdot \boldsymbol{g}(x)
\end{align}
証明には後のセクションで説明(紹介)しているラグランジュの未定乗数法を利用します。
正確には不等式制約の入ったラグランジュの未定乗数法(KKT条件)です。
証明の気持ちについて説明します。
KKT条件を利用して、ラグランジュの未定乗数を$\boldsymbol{z}$とおくとラグランジアン$\mathcal{L}$は
\begin{align}
\mathcal{L} &= f(x) + \boldsymbol{z} \cdot \boldsymbol{g}(x) = f(x) + \sum_i \boldsymbol{z}_i g_i(x) \\
\boldsymbol{z} &\ge 0
\end{align}
となります。$x$についての微分が0となることより、
\begin{align}
\nabla_x f(x) &= - \sum_i z_i \nabla_x g_i(x)
\end{align}
が得られる。また$\boldsymbol{g}(x) \le 0$を満たす$x$を$x^*$とおく。これらの条件
\begin{align}
\boldsymbol{g}(x^*) &\le 0 \\
\nabla_x f(x) &= - \sum_i z_i \nabla_x g_i(x)
\end{align}
を利用して、$f(x^*) - f(x)$を計算していくと、
\begin{align}
f(x^*) - f(x) &= \nabla_x f(x)(x^* - x) = f(x) + \nabla_{x^*} f(x)(x^*-x) - f(x) + O(x^{*2})\\
& \ge \nabla_{x^*} f(x)(x^*-x) = - \sum_i z_i \nabla_x g_i(x) (x^*-x)\\
&\ge -\sum_i z_i\biggl[g_i(x^*) - g_i(x)\biggr] \ge\sum_i z_i g_i(x)
\end{align}
つまり、
\begin{align}
f(x^*) &\ge f(x) + \sum_i z_i g_i(x) = \mathcal{L} (x,\boldsymbol{z})\\
\boldsymbol{g}\le 0 & , \;\; \boldsymbol{z} \ge 0
\end{align}
$f(x^*)$の最小化にすることは、$\mathcal{L} (x,\boldsymbol{z})$を最大化することと等価である。
証明の雰囲気としてはこんな感じです。
ラグランジュの未定乗数法
こちらの記事がすごく分かりやすいです。丸投げします。
https://www.iwanttobeacat.com/entry/2018/01/28/114614
制約条件が不等式の場合(KKT条件)もこちらで説明されています。
https://www.iwanttobeacat.com/entry/2018/01/30/210446
QUBO化
いよいよ本題に入ります。
ReLU関数を次の様に定義する。
\begin{align}
f(m) &= -\min\{0,m\}\\
&= \begin{cases}
-m & (m <0) \\
0 & (m \ge 0)
\end{cases}
\end{align}
ルジャンドル変換を施す
この$f(m)$にルジャンドル変換を行うと
\begin{align}
f^*(t) &= -\inf_m\{h(m,t)\}\\
h(m,t) &\equiv f(m) - mt
\end{align}
$\min$を行う中身の関数を$h(m,t)$と書き直しています。
$f^*(t)$の中身を計算するために、(a)$m<0$、(b)$m=0$、(c)$m>0$ に場合分けして計算する。
(a) $m<0$
\begin{align}
h(m,t) &= -m(t+1)\\
-\inf_m\{h(m,t)\} &= -\inf_m\{m(t+1)\} =0 \Leftrightarrow t = -1\\
f^*(t) &= 0
\end{align}
この区間では、$f^*(t)=0$となります。また$t=-1$のみが有効となる。
(b) $m=0$
微分が不連続なので片側の極限値をそれぞれ計算します。
\begin{align}
h_-(m) &= -m(t+1)\\
h_+(m) &= -mt\\
t_- = \left.\frac{\partial h_-(m) }{\partial m}\right|_{m=0} = -1 ,\;\;\; t_+ &= \left.\frac{\partial h_+(m=0) }{\partial m} \right|_{m=0}= 0\;\rightarrow t \in [-1,0]\\
f^*(t) &= 0
\end{align}
この区間では、$f^*(t)=0$となります。また$t \in [-1,0]$となる。
(c) $m>0$
最後に、
\begin{align}
h(m,t) &= -mt\\
-\inf_m\{h(m,t)\} &= -\inf_m\{-m t\} = 0 \Leftrightarrow t = 0\\
f^*(t) &= 0
\end{align}
この区間では、$f^*(t)=0$となります。また$t=1$のみが有効となる。
以上より$f(m)$のルジャンドル変換は
\begin{align}
f^*(t) = 0, \; t\in [-1,0]
\end{align}
となります。
逆ルジャンドル変換を施す
$f^*(t)$に逆ルジャンドル変換を行うと
\begin{align}
F(m) &= -\min_t\{f^*(t)-mt\} \\
&= -\min_t\{-mt\},\; t\in [-1,0],\;m\in \mathbb{R}
\end{align}
$t$の値域を改めて
\begin{align}
\boldsymbol{g}(t) &= (-1-t, t)^{\mathrm{T}}\\
\boldsymbol{g}(t) &\le0
\end{align}
と書き直す。
Wolf変換を施す
この最小値問題はWolf変換により
\begin{align}
L(m,t,\boldsymbol{z}) &= -mt +\boldsymbol{z}\cdot \boldsymbol{g}(t) = -mt -z_1(t+1) + z_2t
\end{align}
に対して
\begin{align}
-\max_{t,\boldsymbol{z}} &\{L(m,t,\boldsymbol{z})\}\\
{\rm subject \;to\;} & \nabla_{t} L(m,t,\boldsymbol{z})=0 \Leftrightarrow -m -z_1 + z_2 = 0 \\
&\boldsymbol{z} \ge 0, \; \boldsymbol{g}(t) \le0
\end{align}
を求める問題になる。
ここまでで新たに導入された$t,z_1,z_2$は連続変数ですが、実装時にはバイナリ変数で展開する必要があります。
ラグランジュの未定乗数法を用いる
上記の条件を考慮してラグランジュ未定乗数を用いて改めて$L(m,t,\boldsymbol{z})$を次の様に定義する。
\begin{align}
L(m,t,\boldsymbol{z}) & -mt -z_1(t+1) + z_2t -M(-m -z_1 + z_2)^2
\end{align}
ここで、$-1\le t \le 0, z_1 \ge 0, z_2 \ge 0 $である。
これより
\begin{align}
-\max_{t,\boldsymbol{z}} \{L(m,t,\boldsymbol{z})\} = \min_{t,\boldsymbol{z}} \{-L(m,t,\boldsymbol{z})\}
\end{align}
したがってコスト関数$C(m)$を含めた最小化は次のようになる。
\begin{align}
\min_{m} \biggl\{C(m)-\max_{t,\boldsymbol{z}} \{L(m,t,\boldsymbol{z})\} \biggr\} &= \min_{m} \biggl\{C(m)+\min_{t,\boldsymbol{z}} \{-L(m,t,\boldsymbol{z})\} \biggr\}\\
&= \min_{m,t,\boldsymbol{z}} \biggl\{C(m)-L(m,t,\boldsymbol{z}) \biggr\}
\end{align}
ただし、最小化は$-1\le t \le 0, z_1 \ge 0, z_2 \ge 0 $の範囲で行う。
後は$t, z_1, z_2$をバイナリ変数で展開すれば$L(m,t,\boldsymbol{z})$はReLUのQUBOになっている。
補助変数をバイナリ変数で展開する
連続変数をバイナリ変数で展開する方法にいて説明します。
区間を何等分するかを決める量として
\begin{align}
\delta_i(k) = \frac{2^{k-1}}{2^{d_i} -1}
\end{align}
を定義する。これは考えている区間を$2^{d_i} -1$等分することを意味している。$a(k) = r^{k-1}a_0$に対して$\sum^{n}_{k=1}a(k) = \frac{a_0(1-r^n)}{1-r}$より
\begin{align}
\sum^{d_i}_{k=1} \delta_i(k) &= \sum^{d_i}_{k=1} \frac{2^{k-1}}{(2^{d_i}) -1}\\
&= \frac{1-2^{d_i}}{1-2}\times \frac{1}{2^{d_i}-1} = 1
\end{align}
である。これらを使ってバイナリ変数で連続変数を表すと
\begin{align}
t \rightarrow \dot{t} &= \alpha_t \sum^{d_t}_{k=1}t_k\delta(k)_t + \beta_t\\
\boldsymbol{z} \rightarrow \dot{\boldsymbol{z}} &= \alpha_{\boldsymbol{z}} \sum^{d_{\boldsymbol{z}}}_{k=1}{\boldsymbol{z}}_k\delta_{\boldsymbol{z}}(k) + \beta_{\boldsymbol{z}}
\end{align}
となる。これにより区間$[\beta_i, \alpha_i +\beta_i]$を$1/(2^{d_i} -1)$ずつ増える離散変数として表すことができた。上記の方法では$1/(2^{d_i} -1)$とずつ増えていくという、ある意味キリの悪い数値になってしまっている。これN等分にするためには
\begin{align}
2^M \ge N \ge 2^{M+1}
\end{align}
となる$M$に対して
\begin{align}
\dot{x} = \frac{\alpha_x}{N}\biggl(\sum^{M-1}_{k=0}2^{k-1}x_k + (N + 1 - 2^{M})x_M \biggr) + \beta_x
\end{align}
と変更すれば良い。
一般的なReLU関数に適用する
今までのReLUは原点から傾き1で増加している場合を考えていたが、原点以外で傾きも1以外の場合も考える。$a>0,v>0$として次のReLUを考える。今回は最初のReLUの定義とは逆に右肩上がりの関数を考えます。
\begin{align}
f(m) &= -\min\{0,-a(m-v)\}\\
&= \begin{cases}
0 & (m <v) \\
a(m-v) & (m \ge v)
\end{cases}
\end{align}
基本的にしている計算は前に説明したものと同じです。
まず上の関数にルジャンドル変換を行うと
\begin{align}
f^*(t) = vt, \; t\in [0,a]
\end{align}
となる。$f^*(t)$に逆ルジャンドル変換を行うと
\begin{align}
F(m) = -\min_t\{vt-mt\},\; t\in [0,a],\;m\in \mathbb{R}
\end{align}
制約条件を改めて
\begin{align}
\boldsymbol{g}(t) = (-t, t-a)^{\mathrm{T}}\\
\boldsymbol{g}(t) \le0
\end{align}
と書き直すことができる。この最小値問題はWolf変換により、ラグランジュ関数
\begin{align}
L(m,t,\boldsymbol{z}) &= vt-mt +\boldsymbol{z}\cdot \boldsymbol{g}(t) = -t(m-v) -z_1t + z_2(t-a)
\end{align}
に対して
\begin{align}
-\max_{t,\boldsymbol{z}} &\{L(m,t,\boldsymbol{z})\}\\
{\rm subject \;to\;} & \nabla_t L(m,t,\boldsymbol{z})=0 \Leftrightarrow -m + v -z_1 + z_2 = 0 \\
&\boldsymbol{z} \ge 0, \; \boldsymbol{g}(t) \le0
\end{align}
を求める問題になる。上の制約条件を考慮した$L(m,t,\boldsymbol{z})$は改めて
\begin{align}
L(m,t,\boldsymbol{z}) &= -t(m-v) -z_1t + z_2(t-a) -M(-m + v -z_1 + z_2)^2 \\
-L(m,t,\boldsymbol{z}) &= t(m-v) +z_1t - z_2(t-a) +M(-m + v -z_1 + z_2)^2\\
&\boldsymbol{z}\ge 0,0\le t \le a
\end{align}
と書ける。
L1ノルム関数に適用する
L1ノルム関数のQUBO化を行う。L1L1ノルム関数を次のように定義する。
\begin{align}
f(m) &= |m| \\
&= \begin{cases}
-m & (m <0) \\
m & (m \ge 0)
\end{cases}
\end{align}
この関数にルジャンドル変換を行うと
\begin{align}
f^*(t) = 0, \; t\in [-1,1]
\end{align}
となる。$f^*(t)$に逆ルジャンドル変換を行うと
\begin{align}
F(m) = -\min_t\{-mt\},\; t\in [-1,1],\;m\in \mathbb{R}
\end{align}
制約条件を改めて
\begin{align}
\boldsymbol{g}(t) = (-1-t, t-1)^{\mathrm{T}}\\
\boldsymbol{g}(t) \le0
\end{align}
と書き直すことができる。この最小値問題はWolf変換により、ラグランジュ関数
\begin{align}
L(m,t,\boldsymbol{z}) &= -mt +\boldsymbol{z}\cdot \boldsymbol{g}(t) = -mt -z_1(t+1) + z_2(t-1)
\end{align}
に対して
\begin{align}
-\max_{t,\boldsymbol{z}} &\{L(m,t,\boldsymbol{z})\}\\
{\rm subject \;to\;} & \nabla_t L(m,t,\boldsymbol{z})=0 \Leftrightarrow -m -z_1 + z_2 = 0 \\
&\boldsymbol{z} \ge 0, \; \boldsymbol{g}(t) \le0
\end{align}
を求める問題になる。上の制約条件を考慮した$L(m,t,\boldsymbol{z})$は改めて
\begin{align}
L(m,t,\boldsymbol{z}) &= -mt -z_1(t+1) + z_2(t-1)-M(-m -z_1 + z_2)^2 \\
\boldsymbol{z} \ge 0 &, \;\; -1 \le t \le 1
\end{align}
と書ける。
個人的な解釈・疑問
この論文の定式化は結局、スラック変数を導入してエネルギーの差分を埋めようとしているんじゃないかと解釈しました。
「一般的なReLU関数に適用する」で作ったQUBOを例に弄ってみます。
\begin{align}
-L = (m-v)t + z_1t-z_2(t-a) +M(-m+v-z_1 +z_2)^2
\end{align}
$-m+v-z_1 +z_2=0$の拘束条件が満たされる場合を考え、
\begin{align}
\delta = m-v =z_2 - z_1
\end{align}
とおき、$z_2=z_1 +\delta$を代入すると
\begin{align}
-L = a(z_1 +\delta) = a(z_1 + m-v)
\end{align}
となる。結局$m-v$の差を埋めるようにスラック変数$z_1$を導入していることに対応している(？)。
実装上の課題
ルジャンドル変換を施す時に導入される連続変数$t$、ラグランジュの未定乗数$\boldsymbol{z}$をバイナリ変数で展開する時のビット数の規模が実装上の課題になってくるのかなと思います。
$t$は値域が決まっているので、後はどれだけその範囲を細かく刻むかで使うビット数が決まってきます。
一方で$\boldsymbol{z}$は0以上の実数であることしか分からないないので、どれだけの範囲をカバーできるビットを用意すれば良いか分かりませんよね。
ここはトライアンドエラーで見つけていくしかないんでしょうか。
無理やりReLUを使うよりも二乗ペナルティにスラック変数を導入した方がよっぽど特殊な定式化が必要となる時以外は無難かもしれません。
元々先行研究ではラベルノイズに強い耐性を持つようにするためにReLU関数が導入されていました。
単純なスラック変数の導入ではペナルティの線形増加にはならないので、このような場合だとReLUが必要になってくる気がします。
Wolf変換は他の問題にも使えそうな気がしてきました。