PythonでやるML Cycle6: SVM(2)-制約あり最適化問題とその解き方

制約ありの最適化問題について。ラグランジュの未定乗数法、KKT条件、内点法、双対問題、双対定理など。

以下、考える関数はすべて狭義凸関数とします。

# ゴール

* 等式制約の制約あり最適化問題を解く
* 一般的な制約あり最適化問題を解くための条件(KKT条件)を理解する
* 制約あり最適化問題を内点法で解く
* 主問題と双対問題の解が一致することを実感する(双対定理を理解する)

# 制約あり問題

# 制約あり問題の概観

[前回](https://qiita.com/kenji-kondo/items/77749b462d033e1a3189)のつづきです。前回は**制約なし**最適化問題について考え、その解き方(ニュートン法)について説明しました。今回は、以下のような**制約あり**最適化問題を考えます。まずはざっくり概観してみましょう。

```math
\begin{align}
\min_{\vec{x}}\  & f(\vec{x}) \\
\text{s.t.}\  & g_i(\vec{x}) \leq 0, \ i = 1, 2, 3, \cdots, M
\end{align}  
```

まず、この制約あり問題を**素朴に**解くことが、制約なし問題を解くことに比べていかに難しいか、その気持ちを説明します。

前回も紹介した簡単な例で説明します。以下の制約あり最適化問題を考えてます。

```math
\begin{align}
\min_{x} \ & x^2 \\
\text{s.t.} \ & -(x-2) \leq 0
\end{align}
```

頭の中で考えてみると、答えは2です。それでは、この問題をどのように解いたでしょうか。

$y=x^2$のグラフを頭に描いたのではないでしょうか。そして、その頂点の座標が$(0,0)$であることに注意して、また制約条件$x \geq 2$であることから頂点が最小値にならないことがわかって、結果$x=2$となる場所が$x^2$の最小値となる、という風に考えたのではないでしょうか。

次に、制約条件を少しだけ変えた問題を考えます。

```math
\begin{align}
\min_{x} \ & x^2 \\
\text{s.t.} \ & -(x+2) \leq 0
\end{align}
```

この場合、制約条件が$x \geq -2$となるので、制約条件の中に$x^2$の頂点の座標が含まれます。従って、この問題の解は頂点$x=0$となります。

このように、制約あり問題を素朴に解こうとした場合、考えている関数$f(\vec{x})$や制約条件の式の**具体的な形を考慮する**必要があります。

上の例にあげた簡単な1変数の関数の問題なら、まだグラフのイメージを描けるのでいいでしょう。では次のような問題はどうでしょう。

```math
\begin{align}
\min_{x,y} \ & 2 x + 3y -1 \\
\text{s.t.} \ & x^2 + y^2 -3 \leq 0 \\
& -x \leq 0 \\
& -y \leq 0
\end{align}
```

また次のような問題はどうでしょう。

```math
\begin{align}
\min_{\vec{w}, b, \vec{\xi}} \ &||\vec{w}|| + C \vec{1}^T \vec{\xi} \\
\text{s.t.} \ & -(y_i (\vec{w}^T \vec{x_i} + b) - 1 + \xi_i)) \leq 0, \ i = 1, 2, \cdots, N \\ 
& - \xi_i \leq 0, \ i = 1, 2, \cdots, N 
\end{align}
```

解く気が失せます。具体的なグラフのイメージをするなんてできません。(実は↑はSVMの最適化問題です)

一方、制約なし問題の解法はどうだったでしょうか。制約なし最適化問題は以下のようなものでした：

```math
\min_{\vec{x}} \ f(\vec{x})
```

これの解は、以下の(連立)方程式を解くことで得られました。

```math
\nabla f(\vec{x}) = \vec{0}
```

そして上の連立方程式を解くにはニュートン法を使えば近似解を得られました。

ここで注目したいのは、「$f(\vec{x})$がどんな形であっても、その勾配ベクトルが0ベクトルとなる式を解くだけで、解を得られる。出てくる連立方程式はニュートン法を使って解けばよい」なる**汎用的な解法が存在する**というところです。

となると、制約あり問題にも、制約なし問題のような、汎用的な解法があることを期待します。そんなものがあるのか？---あります。それが**ラグランジュの未定乗数法**、また**KKT条件**、**内点法**と呼ばれるものです。これからこの3つの解法について説明します。

また、制約あり問題には特有の面白い現象があります。

やはり以下の簡単な問題を考えます。

```math
\begin{align}
\min_{x} \ & x^2 \\
\text{s.t.} \ & -(x-2) \leq 0
\end{align}
```

唐突ですが、以下の問題を考えます。

```math
\begin{align}
\max_{\lambda} \ & -\frac{\lambda^2}{4} + 2 \lambda -2 \\
\text{s.t.} \ & -(\lambda - 4) \leq 0
\end{align}
```

二つ目の問題が$\max$に関する最適化問題になっていることに注意してください。

実は、この2つ目の問題を解けば、1つ目の問題を解いたことになります。つまり2つ目の問題は、1つ目の問題と双子のような関係にあります。1つ目の問題を**主問題**、二つ目の問題を主問題に対する**双対(そうつい)問題**と呼びます。そして、双対問題を解けば主問題を解いたことになる、またその逆も成立する、という関係を**双対定理**と呼びます。

双対問題なんていう、かえって複雑になった問題をなぜ考えるのか？と思われるかもしれません。一番の理由は、SVMの双対問題を拡張することによって、線形分離不可能なデータの分類をうまくやってくれるようになるからです。これは後で説明します。

以上、今回お話する概観でした。まずはラグランジュの未定乗数法から考えていきます。

# 制約あり問題-等式制約の場合

## 説明

以下の制約あり問題を考えます。

```math
\begin{align}
\min_{\vec{x}} \ & f(\vec{x}) \\
\text{s.t.} \ & h_i(\vec{x}) = 0, \ i = 1, 2, \cdots, M
\end{align}
```

制約条件が不等式ではなく等式になっています。これを**等式制約**と呼びます。

等式制約の最適化問題の解き方は以下です。ラグランジュの未定乗数法と呼ばれています。

1. $\vec{\lambda} = (\lambda_1, \lambda_2, \cdots, \lambda_M)^T$, $\lambda_i \geq 0, \ i=1,2,\cdots,M$として、$L(\vec{x}, \vec{\lambda}) = f(\vec{x}) + \sum_{i=1}^{M} \lambda_i h_i (\vec{x})$ という関数$L$を考える。これを**ラグランジュ関数**と呼ぶ。
2. $\nabla_{\vec{x}} L (= \nabla_{\vec{x}} f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla_{\vec{x}} h_i  (\vec{x})) = \vec{0}$と$h_i(\vec{x}) = 0, \ i = 1, 2, \cdots, M$の連立方程式を解いて、$\vec{x}$(と$\vec{\lambda}$)を求める。


さて、なぜこのラグランジュ関数の勾配$=\vec{0}$を考える必要があるのか。それを説明します。以下は直観的な理解を優先しているため、厳密でないところがあります。

すべての$i$に対して$h_i(\vec{x}) = 0$を満たしつつ、$f(\vec{x})$を最小にする変数を$\vec{y}$とおきます。つまり、$h_i(\vec{x}) = 0$を満たすあらゆる$\vec{x}$に対して、

```math
f(\vec{x}) > f(\vec{y})
```

が成立します。

さらにラグランジュ関数の勾配が$\vec{0}$でないと仮定しましょう。すなわち、あらゆる$\vec{\lambda} \geq \vec{0}$(これは$\vec{\lambda}$のすべての成分が0以上、という意味です)に対し、

```math
\nabla f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla h_i(\vec{x}) \neq \vec{0}
```

とします。

こらから先の展開を説明しておくと、これらを仮定すると、あるところで矛盾が生じます。矛盾が生じてしまったということは、仮定が間違っていた、つまり本当はある$\vec{\lambda}$で$\nabla f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla h_i(\vec{x}) = \vec{0}$が成立していたはずだ、という展開です。(いわゆる背理法です。)

では矛盾を導きましょう。


簡単のため、$H(\vec{x}) = (h_1(\vec{x}), h_2(\vec{x}), \cdots, h_M(\vec{x}))^T$と表記します。まず$H(\vec{y} + \vec{\epsilon}) = \vec{0}$となる$\vec{\epsilon}$を考えます。$||\vec{\epsilon}||$がめちゃくちゃ小さい状態を考えると、$H(\vec{y} + \vec{\epsilon})$を1次の項までテイラー展開すると、

```math
H(\vec{y} + \vec{\epsilon}) = H(\vec{y}) + \partial H(\vec{y}) \vec{\epsilon}
```

となります(前回のも同じような式展開をしました)。$H(\vec{y} + \vec{\epsilon}) = \vec{0}$, $H(\vec{y}) = \vec{0}$であることを考慮すると以下が成立します。

```math
\partial H(\vec{y}) \vec{\epsilon} = \vec{0}.
```

つまり、$\vec{\epsilon}$はすべての$\nabla h_i(\vec{y})$と直交しています。

ここで仮定を思い出します。$\nabla f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla h_i(\vec{x}) \neq \vec{0}$ということは、$\nabla f(\vec{y})$は$h_i(\vec{y})$で張られる空間に含まれていません。このことから$\vec{\epsilon}$として、$\nabla f(\vec{y})$の成分を含む、以下のようなものを選んでみます(必ずしも選べるとは限りません)。

```math
\vec{\epsilon} = \alpha \nabla f(\vec{y}) + \beta \nabla f(\vec{y})^{\perp}, \ \alpha > 0.
```

$\alpha > 0$に注意してください。ちなみに$f(\vec{y})^{\perp}$は$\nabla f(\vec{y})$と直交したベクトルの一つで、$\nabla f(\vec{y})^T \nabla f(\vec{y})^{\perp} = 0$が成立するものです。

ここで$f(\vec{y} + \vec{\epsilon})$のテイラー展開を考えます。

```math
\begin{align}
f(\vec{y} + \vec{\epsilon}) &= f(\vec{y}) + \nabla f(\vec{y})^T \vec{\epsilon} \\
&= f(\vec{y}) + \alpha || \nabla f(\vec{y}) ||^2 .
\end{align}
```

$\alpha > 0$, $|| \nabla f(\vec{y}) || > 0$であることに注意すれば、

```math
f(\vec{y} + \vec{\epsilon}) < f(\vec{y})
```

となってしまうので、これは最初に置いた仮定「$\vec{x} = \vec{y}$で$f(\vec{x})$は最小となる」に矛盾します。(ちなみに$||\nabla f(\vec{y})||^2 = 0$となる場合は考える必要がありません。なぜならどんな$\vec{\lambda}$に対しても$\nabla f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla h_i(\vec{x}) \neq \vec{0}$なので、$\vec{\lambda} = \vec{0}$である場合を考えても$\nabla f(\vec{y}) \neq \vec{0}$だからです。)

以上が等式が成立する理由です。まとめると、「もしも等式が成立しないならば、少し動かしたらさらに$f(\vec{x})$が小さくなるようなものを見つけられてしまう」ということです。多変数の場合でやったので少々煩雑な説明になった気もします。$M=1$の場合でこの説明をトレースすると気持ちがわかるかもしれません。

ところで、$\vec{\lambda} \geq \vec{0}$の仮定はどういう意味を持つでしょう。それはこのページの一番最初に置いた仮定「考えるあらゆる関数は狭義凸関数である」から来ています。$\vec{\lambda}$がすべて正であることにより、ラグランジュ関数もまた凸性を持ちます。ちょっとテクニカルですが、[このページ](https://qiita.com/kenji-kondo/items/d9aae38b461d3c340086)にまとめました。


## 確認

1. ラグランジュ未定乗数法の上の説明を、2次元空間上で$M=1$と置き換えて説明できるようになってください。
2. 以下の最適化問題をラグランジュ未定乗数を使って解いてください。

```math
\begin{align}
\min_{x, y} \ & x^2 + y^2 \\
\text{s.t.} \ & x^2 + y^2 -2x -2y + 1 = 0
\end{align}
```

# 制約あり問題-不等式制約

それでは等式制約ではない不等式制約の最適化問題について考えます。以下のようなものです。

```math
\begin{align}
\min_{\vec{x}}\ & f(\vec{x}) \\
\text{s.t.} \ & g_i(\vec{x}) \leq 0, \ i = 1, 2, \cdots, M.
\end{align}
```

これは以下の式をすべて満たす$\vec{x}$が解になります。これは**KKT(カルーシュ・クーン・タッカー)条件**と呼ばれます。式がたくさんありますが、ほとんどラグランジュの未定乗数法と同じです。一つだけ条件が追加されました。

```math
\begin{align}
(1)\  & \nabla_{\vec{x}} L(\vec{x}, \vec{\lambda}) = \nabla f(\vec{x}) + \sum_{i=1}^{M} \lambda_i \nabla g_i(\vec{x}) = \vec{0} \\
(2)\ & g_i(\vec{x}) \leq 0, \ i = 1,2,\cdots,M \\
(3)\ & \vec{\lambda} \geq \vec{0} \\ 
(4)\ & \lambda_i g_i(\vec{x}) = 0, \ i = 1,2,\cdots,M.
\end{align}
```

(1)-(3)は、ラグランジュの未定乗数法で考えた条件と同じです。(4)の式が新しく追加されました。この条件を**相補性条件**と呼びます。

相補性条件について説明します。この問題の解を$\vec{y}$と置いて、かつある$i$で$g_i(\vec{x}) < 0$が成立したとしましょう。この相補性条件に当てはめると、$g_i(\vec{x}) < 0$と真に$g_i(\vec{x})$が0よりも小さくなっているので、相補性条件を成立させるためには$\lambda_i = 0$としなければなりません。これはなんのためにあるかというと、この相補性条件によって、無駄な制約を取り除いてくれる、というかんじです。例えば制約条件が$-(x+2)\leq 0$であるとき、つまり$-2 \leq x$であるときの$x^2$の最小値を考えると、解は0になります。この問題は、この制約があってもなくても同じ解を得ますよね。この問題のラグランジュ関数を考えると、

```math
L(x, \lambda) = x^2 -\lambda (x+2)
```

となりますが、相補性条件

```math
\lambda \cdot (- (x+2)) = 0
```

より、これの解を$y$とすれば、$-(y+2) < 0$が成立するので相補性条件により$\lambda = 0$となります。すなわちこのラグランジュ関数は

```math
L(x, \lambda) = x^2
```

となり、制約なし問題へと変わります。相補性条件が、考えなくてもよい制約を省いてくれました。

以下の確認を解いてみましょう。

## 確認

1. 以下の最適化問題をKKT条件を解くことで解いてください。(ヒント: 2番目の制約条件が0よりも小さくなることを先に示して相補性条件を使うと簡単にできます。)

```math
\begin{align}
\min_{x,y}\ & x^2 + y^2 \\
\text{s.t.}\ & (x-1)^2 + (y-1)^2 - 1 \leq 0 \\
& x^2 + y^2 -4 \leq 0
\end{align}
```

# 内点法

ここからは、実装に向けて、具体的に制約あり問題を解くアルゴリズムを考えます。具体的には**内点法**について説明します。

内点法は、簡単に言えば、不等式制約の制約あり問題を等式制約の問題に置き換えるアルゴリズムです。

具体的に見ていきます。

以下のような最適化問題を考えます。もう覚えましたかね。

```math
\begin{align}
\min_{\vec{x}} \ & f(\vec{x}) \\
\text{s.t.} \ & g_i(\vec{x}) \leq 0, \ i = 1,2,\cdots,M
\end{align}
```

内点法は以下のように問題を変更します。

$\mu > 0$を決めて、$\vec{s} = (s_1, s_2, \cdots, s_M)^T > \vec{0}$という変数を導入して、

```math
\begin{align}
\min_{\vec{x}, \vec{s}} \ & f(\vec{x}) - \mu \sum_{i=1}^M \log s_i \\
\text{s.t.} \ & g_i(\vec{x}) + s_i = 
0, \ i=1,2,\cdots,M \\
& s_i > 0, \ i=1,2,\cdots, M
\end{align}
```

何をしたかというと、$s_i > 0, \ i=1,2,\cdots,M$という変数を使って、不等式制約を等式制約に無理やり変えました。これで何が嬉しいかというと、等式制約でみたラグランジュの未定乗数法がそのまま使えるというわけです。思い出してみれば、等式制約の問題は、ただの連立方程式になりました。ということは、前回やったニュートン法がそのまま使えます。また$\mu > 0$の条件ですが、これは予め定めておき、問題の中では定数として扱います。$\mu$が0に近づくほど元の問題に近づくことを念頭に置いておいてください。

この問題のラグランジュ関数は以下のようになります：

```math
L(\vec{x}, \vec{s}, \vec{\lambda}) = f(\vec{x}) - \mu \sum_{i=1}^{M} \log s_i + \sum_{i=1}^{M} \lambda_i g_i(\vec{x})
```


アルゴリズムは以下です。

1. $\mu > 0$を適当に決める。
2. 内点法のために置き換えた上の最適化問題をとく。($\nabla_{\vec{x}, \vec{s}} L = \vec{0}$と、等式制約$g_i(\vec{x}) = 0$の連立方程式で$\vec{x}, \vec{s}, \vec{\lambda}$をそれぞれ求める。)
3. $\mu$をちょっとだけ0に近づける。(例えば$\mu$を1/2倍するなどして)
4. 再度2を解く。
5. $\mu$がある程度小さくなるまで(例えば0.001を下回るまで)3.と4.を繰り返す。

## 確認

1. 以下を内点法で解いてください。前に解いた問題と近い解がでれば成功です。

```math
\begin{align}
\min_{x,y}\ & x^2 + y^2 \\
\text{s.t.}\ & (x-1)^2 + (y-1)^2 - 1 \leq 0 
\end{align}
```

# 双対問題

書きかけ