読者です 読者をやめる 読者になる 読者になる

最適化:非線形計画+組み合わせ最適化のまとめのメモ

はじめに

{
 目的関数:  min   f(x) \\
 制約条件:  c_i(x) \\
 x はすべてn次元ベクトル
}
としてほとんどの問題を定式化する。細かい記号の定義や意味は略してます。ところどことミスがあるかも。参考文献は文末、以下は教科書。

[0] 福島雅夫:『数理計画入門』, 朝倉書店, 1996.

非線形計画

問題の定式化

{
 目的関数:  min   f(x) \\
 制約条件:  c_i(x) \\
 x はすべてn次元ベクトル
}
として、各問題を定式化。

凸関数と凸領域

凸関数

{
 x,y \in R^n(もしくは実行可能領域)   0 \leq \alpha \leq 1 \\
 -> f(\alpha x + (1 - \alpha)y ) \leq \alpha f(x) + (1 - \alpha)f(y)
}
を満たす関数。

凸領域

{
x,y \in R^n(もしくは実行可能領域)    0 \leq \alpha \leq 1 \\
 -> \alpha x + (1 - \alpha) y \in R^n(もしくは実行可能領域)
}
となる領域。

大域最適解と局所最適解

以下を参照。

大域的最適解とは - OR事典 Weblio辞書


局所的最適解とは - OR事典 Weblio辞書


大域的収束性 - 数理計画用語集

局所的収束性 - 数理計画用語集


ここで「目的関数が凸関数で実行可能領域が凸領域なら、局所最適解と大域最適解は一致する」と言える。

証明:
 背理法、一致しないと仮定して局所最適解{\vec{x} }と大域最適解{ \vec{y}}があるとする。もちろん{ \vec{x} \vec{y} }。ここで凸集合であることから{\vec{x},\vec{y} \in S , 0 \leq \alpha \leq 1 \to \alpha \vec{x} + (1 - \alpha)\vec{y} \in S} が成り立つ。
 一方目的関数fは凸関数なのだから
{ \vec{x},\vec{y} \in R^n , 0 \leq \alpha \leq 1 \\
 \to f(\alpha \vec{x} + (1 - \alpha)\vec{y}) \leq \alpha f(\vec{x}) + (1 - \alpha)f(\vec{y})
  }
が成り立つ。すなわち{ f(\alpha \vec{x} + (1 - \alpha)\vec{y})}は大域最適解よりも小さい値になる、これはおかしい■

勾配とヘッセ行列

{ \displaystyle
\Delta f (\vec{x})= \nabla^2 f (\vec{x})}がヘッセ行列で勾配は{ \nabla f(\vec{x})}と記述する。詳しい定義は*1参照。勾配はその点での関数が最も増加する方向へのベクトル(つまり"等高線"に垂直な方向)であることを後述の最急降下法で利用する。

固有値の符号による行列の種類分け

  • 正定値行列
  • 半正定値行列

凸関数とヘッセ行列の関係

 制約なし問題で、目的関数か凸関数ならば「最適解{x^*}では{\nabla^2 f(x^*)}は半正定値行列である。」ことが成り立つ。これは二次関数での"傾き0の点が局所最適解の候補になりうる"ことに対応。
証明:
 背理法で示す。「最適解{x^*}では{\nabla^2 f(x^*)}は半正定値行列である。」が成り立たないと仮定する。
ここで、目的関数fを任意の点{x'}テイラー展開すると、

{
f(x) \\
= f(x') + \sum_{k=1}^n \frac{d f(x')}{d x_i} (x_i - x'_i) + \frac{1}{2} \sum_{k=1}^n \sum_{l=1}^n \frac{d f(x')}{d x_k} \frac{d f(x')}{d x_l} (x_k - x'_k)(x_l - x'_l) + (三次以上の項)\\
= f(x') + \nabla f(x')^T (x - x') + \frac{1}{2} (x - x')^T \nabla^2 f(x') (x - x') + (三次以上の項)\\
}

と記述できる。今、任意の点{x}は最適解{x^*}と適当な実数{\alpha}とベクトル{d}を用いて{x = x^* + \alpha d}とかける。今点{x}を最適解{x^*}にてテイラー展開すると、上の式に代入して

{
f(x^* + \alpha d) \\
= f(x^*) + \nabla f(x^*)^T (x^* + \alpha d - x^*) + \frac{1}{2} (x^* + \alpha d - x^*)^T \nabla^2 f(x^*) (x^* + \alpha d - x^*) + (三次以上の項)\\
= f(x^*) + \frac{1}{2} \alpha^2 d^T \nabla^2 f(x^*) d + O(\alpha^3)\\
}

とかける。「最適解{x^*}では{\nabla^2 f(x^*)}は半正定値行列である。」が成り立たないと仮定しているのだから、
{\alpha^2 d^T \nabla^2 f(x^*) d < 0}となるようなベクトル{d}が存在するはずである。ここで上のテイラー展開の結果でαの値を十分小さくすると{O(\alpha^3)}の項を無視できる。
すると{\alpha^2 d^T \nabla^2 f(x^*) d < 0}なのだから、

{
f(x^* + \alpha d) = f(x^*) + \frac{1}{2} \alpha^2 d^T \nabla^2 f(x^*) d < f(x^*)
}

であるがこれは最適解が{x^*}であることに矛盾。■

最急降下法ニュートン法準ニュートン法

過去記事参照。

最適化:非線形計画について(最急降下、ニュートン、KKT条件)とNP困難問題に対する動的計画法 - 雑なメモ


もしくは*2など。

アルミホのルールとウルフのルール

文献J. Nocedal and S. J. Wright, Numerical optimization. Springer Series in Operations
Research. Springer-Verlag, 1999.

に説明があります。省略。

ルーシュ・キューン・タッカー条件(KKT条件)*3

カルーシュ・キューン・タッカー条件 - Wikipedia


KKT 条件 - 数理計画用語集

上記のページに詳しく記述されている(ただし出典はこの記事の参考文献におなじらしい)。

問題の定式化

{ 目的関数: f(x) -> min }
{ c_1(x) = 0 , \dots , c_m(x) = 0}
{ d_1(x) \leq 0, \dots , d_n(x) \leq 0}

制約なし問題では、{\nabla f(x) = 0}が最適であるための必要条件のひとつだが、制約ありの場合は{\nabla f(x) = 0}ではなくても実行可能領域の境界上などで最適解がみつかる場合が多いから{\nabla f(x)}という情報だけでは最適解は見つからない。

そこで、最適解である点がいくつもの制約を満たすために制約を満たすような方向へ”引っ張られて釣り合った結果”の点であると考えたら、

{
 \nabla f(x^*) +( u_1 \nabla c_1(x^*) + \dots + u_{m} \nabla c_{m} (x^*) ) + ( u_{m+1} \nabla d_1(x^*) + \dots + u_{m} \nabla d_{m+n} (x^*) ) = \vec{0} 
}

となるような非負の実数{u_1,\dots,u_{(m+n)}}が存在すると言える(証明は参考文献に書いてないし実際学部の範囲ではなさそうなので考え方だけ理解)。{\nabla}をつけるのは勾配は関数が最も増加する方向の方向ベクトルなのだから、最適解において各関数(制約)がもっとも大きく引っ張る方向と個人的に解釈、釣り合うのだからその和は零ベクトルで表される。

そして、有効制約でない制約(最適解である{x^*}を選んだ時に{c_k (x^*) = 0}とならない制約)はこの文脈で言えば最適解の点を制約がなんとか満たすように"引っ張っていない(まだいくらか動ける余裕がある)"のだから

{
 \nabla f(x^*) +( u_1 \nabla c_1(x^*) + \dots + u_{m} \nabla c_{m} (x^*) ) + ( u_{m+1} \nabla d_1(x^*) + \dots + u_{m} \nabla d_{m+n} (x^*) ) = \vec{0} 
}

にて最適解を移動させるための役割を持たないはずだと考えると「有効制約でないものに対応する係数uは0である」と考えられる。参考文献にはもっとしっかりベクトルの一次従属などを用いて説明があった。

そしてもしも制約{d_1(x),\dots,d_n(x)}が存在しないとすると、この"引っぱりあい"を示す式は

{
 \nabla f(x^*) +( u_1 \nabla c_1(x^*) + \dots + u_{m} \nabla c_{m} (x^*) )  = \vec{0} 
}

と書くことができる。これはまさにラグランジュ未定乗数法の式そのもので実際この最適化の文脈でもuをラグランジュ係数と呼ぶらしい。

ラグランジェの未定乗数法 [物理のかぎしっぽ]
ラグランジュの未定乗数法 - Wikipedia

KKT条件の定式化

 以上の考えをもって、KKT条件を以下に記述(参考文献[0]p121)。

KKT条件
{ \displaystyle
  \left\{
    \begin{array}{l}
     \nabla f(x^*) +( u_1 \nabla c_1(x^*) + \dots + u_{m} \nabla c_{m} (x^*) ) + ( u_{m+1} \nabla d_1(x^*) + \dots + u_{m} \nabla d_{m+n} (x^*) ) = \vec{0} \\
\\
     u_k = \left\{
    \begin{array}{lll}
      u_k \geq 0 & if & c_k({x^*}) = 0 \\
      0 & if & d_k({x^*}) < 0 \\
      u_k \geq 0 & if & d_k({x^*}) \leq 0 \\
    \end{array}
  \right.
    \end{array}
  \right.
}

最適解であるための必要・十分条件(参考文献[0]p123)

 KKT条件を満たすすべての点が最適解であるとは限らないのは適当に制約と目的関数を作ってためしてみると分かる。なのでさらに目的関数と制約が凸関数であるという制約を加えることで局所最適解が大域最適解と一致するようにするとKKT条件を満たす点が最適解であると言えるようになる。さらに以下の式Lを考える。

{

L = f(x) + \sum_{k=1}^m u_k c_k (x) + \sum_{k=m+1}^n u_k d_k (x)

}

このLの勾配は

{

\nabla L = \nabla f(x) + \sum_{k=1}^m u_k \nabla c_k (x) + \sum_{k=m+1}^n u_k \nabla d_k (x)

}

であり、まさにKKT条件の”引っ張り合い”の式であって

{

\nabla^2 L = \nabla^2 f(x) + \sum_{k=1}^m u_k \nabla^2 c_k (x) + \sum_{k=m+1}^n u_k \nabla^2 d_k (x)

}

がLのヘッセ行列。そして有効制約の勾配ベクトルと直行するベクトル{  y }{ y^T \nabla^2 L y \geq 0}が成り立つことが知られている。この式は以前の制約なしの凸計画問題の「ヘッセ行列が半正定値行列である」ことと対応しているとみると、見つけた点が最適解であるための十分条件は凸計画問題の「ヘッセ行列が正定値行列である」ことに対応して「{ y^T \nabla^2 L y > 0}」と見当がつく。

追記:KKT条件を満たす点が最適解と保証されるのは

  • 目的関数が凸関数
  • 等式の制約を与える(***=0)制約はすべて一次式
  • 不等式制約を与える(***≦0)制約はすべて凸関数

の場合です。

組み合わせ最適化

問題の定式化

{
 目的関数:  min   f(x) \\
 制約条件:  c_i(x) \\
 x はすべてn次元ベクトル
}
として、各問題を定式化。

ナップザック問題

{
目的関数 :max   \sum_{i = 1}^n z_i a_i  \\
制約条件 :z_i   \in \{0,1 \}   \sum_{i=1}^n z_i c_i \leq b 
}
※記号の意味は後述

巡回セールスマン問題*4

{
d_{ij}:= 地点iと地点jの距離 \\
z_{ij}:= 地点iから地点jへと通る場合は1、通らないなら0 \\
}
とする。巡回セールスマンは同じ場所を通らないし、すべての地点に一回だけ行くことを考慮して制約条件には
{
 z_{ij} \in \{0,1\} \\
 \sum_{i=1}^{町の数} \sum_{j=1}^{町の数} z_{ij} = n - 1 = 町の数 - 1 \\
 \sum_{i=1}^{町の数} z_{ij} = 1   (一回だけしか通らないから) \\
 \sum_{j=1}^{町の数} z_{ij} = 1   (一回だけしか通らないから) \\
 \sum_{i \in S} \sum_{j \notin S} x_{ij} \geq 1 \\
※Sは空でない任意の地点の集合
}
が記述できる。最後の条件がないと、ひとつのグラフ内部でいくつもの循環がある場合でも制約条件を満たしてしまうが巡回セールスマン問題はあくまで「全ての地点を一回だけ通って最初の都市に戻るルート」であることを考慮する。

 目的関数は総距離を最小化したいのだから
{
目的関数: min \sum_i \sum_j z_{ij} d_{ij}
}
と書ける.*5

英文の単語の最適配置問題

 英語の文を印刷するときに、単語の途中で改行させずになおかつ単語間のスペース幅の総和をなるべく少なくしたい。
 今、印刷する文章の一行の幅を{L}として、単語{w_i}{l_i}の幅をもっている。このとき、ある行にあるスペースの幅の総和は

{
 s(i,j) := \frac{ L   -   \sum_{i}^j l_i }{ max\{j - i , 1\}}
}

と書ける。そしてスペースがあるほどコストがかかるとして、

{
 c(i,j) = 
  \left\{
    \begin{array}{l}
      \infty   if   s(i,j)\leq 0 \\
      (s - s(i,j))^2   if   s(i,j)\ge 0\\
      max\{ s - s(i,j),0\}^2 if   s(i,j)\ge 0 かつ j = n 
    \end{array}
  \right. \\
s := ideal \ space \ width
}
最後の条件は文章の最後の行はスペースがたくさんあっても無視するのならば必要な条件となる。動的計画法を適用するためには漸化式を記述して、なにを「メモ化」するか決める必要あり、後述。

双対

問題
{
 目的関数:  min   f(x) \\
 制約条件:  c_i(x) \\
}
について、制約条件をまとめて一つの行列として表し目的関数もベクトルの内積に書き換えたとして、
{
 目的関数:  min   c^Tx \\
 制約条件:  Ax \geq b 
}
と定義する。この問題の”双対”にあたる問題は
{
 目的関数:  max   b^Ty \\
 制約条件:  A^Ty \leq c 
}
と書ける。*6

そして、弱双対定理によると{c^Tx - b^T y \geq 0}.等号が成り立つのは二つの問題に最適解があって、ともに最適解の場合のみ。制約条件を満たすベクトルはすべてこの不等号を満たす。なぜなら、任意の制約条件を満たすベクトルに対して{ y^Tb \leq y^TAx \leq c^Tx  }が成り立つから。
 そして、このような二つの問題を考えるとき、その二つの最適解が存在する、もしくは目的関数が無限に大きく(小さく)なり一方は最適解を持たない、共に解を持たない場合に限られる。弱双対定理から、最適解の下限(上限)を得ることができる。

ナップザック問題での分枝限定法と動的計画法

はじめにナップザック問題(0-1ナップサック問題)*7を定式化。

{
目的関数 :max   \sum_{i = 1}^n z_i a_i  \\
制約条件 :z_i   \in \{0,1 \}   \sum_{i=1}^n z_i c_i \leq b 
}

zはナップザックにアイテムがある時に1、ないときは0に対応した変数。
{a_i:item_iの価値 \ \ b:ナップザックのサイズ \ \ c_i:item_iの大きさ}
に対応する。

貪欲法*8

 最適解だと思うものを優先して(先のことを考慮せずに)選ぶ手法。ナップザック問題ならば例えば「価値の大きいものを優先して先に選ぶ」など。

いくつかのアルゴリズムは貪欲法を基本戦略としているものの、厳密解が求まることが証明されている。
ダイクストラ
プリム法
クラスカル*9

分枝限定法

 「解いても現在得られると分かっている最適解よりも”よい”最適解をえることができないと保証された部分問題が検出されたとき、その部分問題は解かない」ことで計算時間を短縮する。

  • はじめに、アイテムを基準に沿って並べ替え。

 アイテムを単位大きさあたりの価値が最大になるように、すなわち
{
\ \frac{a_1}{c_1} \ \leq \frac{a_2}{c_2} \ \leq \frac{a_3}{c_3} \ \leq \frac{a_4}{c_4} 
}
となるようにする。

  • 探索木を作る。

 {z_i \in \{0,1\}}に基づいた二分木を作成する。木の根のはじめの分岐はもっとも単位重さあたりの価値があるアイテム1({価値a_i, 大きさc_i})として{x = 1}ならばナップザックにアイテムを入れたことを表す。

  • 部分問題を解くべきかの基準を計算する。

 「現在得られると分かっている最適解よりも”よい”最適解をえることができないと保証された部分問題」を判定するための値を計算する。ナップザック問題の例だと「この部分問題が与えうる最適解の上限」と「この部分問題が与えうる最適解の下限」のふたつが部分問題をさらに解く前に、つまり探索木のそれぞれのノードに達したときに計算される必要がある。
 「この部分問題が与えうる最適解の下限{lower-bound}」:とりあえず何でもいいから鞄に入ればそれが解になりうるのだからここでは【単位大きさあたりの価値に沿って並んだアイテムを(単位大きさあたりの価値が大きい)順にいれていってナップザックをはみ出ない最大量の価値】をその値にする。その時、最後に入れたアイテムを{item_j}とする。
 「この部分問題が与えうる最適解の上限{upper-bound}」:【単位大きさあたりの価値に沿って並んだアイテムを(単位大きさあたりの価値が大きい)順にいれていってナップザックをはみ出ない最大量の価値】+ {item_{j+1}}の価値 とする。もともとアイテムは単位あたりの価値が大きい順に並んでいるのだから、どうやってもナップザックにはいる量の範囲内でこの値は超えることができない。

  • 解く必要のない部分問題を決定する

 ノードに幾たびに上記の{upper-bound}{lower-bound}を計算する。ここで明らかに

{既存のあるlower-bound(lb) \geq 今いるノードのupper-bound}

ならば、その部分問題を解いたとして最適解の候補{既存のあるlower-bound(lb)}が解として得られることが分かっているのだからその部分問題を解く意味がない。

 結局以下のようにまとめられる。

  1. はじめに、なにか基準(講義では単位あたりの価値が大の昇順)を定義して、それに基づいてアイテムを並び替える。
    1. 基準を正確に決めないと「{upper bound \leq 既知のlower bound}だからそのノードの子孫を切り捨ててもいい」と言えない..!
  2. その並びに基づいてまだいれるかいれないかを決めていない({z_kの値が定まっていない})アイテムを入れていく。
  3. ナップザックにこれ以上アイテムが入らないとなったとき...以下の値をノードに保存
    1. lower bound := その時の目的関数の値
    2. upper bound := ナップザックをはみ出てもいいので順番に並んだアイテムの次のアイテムを入れた場合の目的関数の値
  4. ノードを移動して2.3.を実効する
  5. そのとき、求められた{upper bound \leq 既知のlower bound}ならばそのノードの子孫は計算せずに切り捨てる
  6. 以上を繰り返し、もっともよい{lower bound}が解となる
木の探索の順番
  • 最良優先探索
    • 最も良い解が得られそうな部分問題のあるノードへと優先的に進む。
    • 適切な基準で部分問題を選べば計算時間は少なくすむ。
  • 深さ優先探索
    • 部分問題の根までの深さがより長い部分問題のあるノードへと優先的に進む。
    • メモリ消費が少なく実装しやすいが、深さの限界が保証されていない場合や計算がたくさん必要な場合あり。

動的計画法

 「繰り返し出現する部分問題の解を記憶することで、部分問題は一回だけしか計算しない」ことで計算時間を短縮する。

対象となる問題を複数の部分問題に分割し、部分問題の計算結果を記録しながら解いていく手法で..(中略)..下記2条件を満たすアルゴリズム

  1. 分割統治法:部分問題を解き、その結果を利用して、問題全体を解く
  2. メモ化:部分問題の計算結果を再利用する

である。
*10

 部分問題を解く、という記述から単純に再帰的にプログラムを書くと計算時間はものすごく大きくなる(NP完全の部分問題もやはりNP完全)。だから、「メモ化」を利用してなんども同じ問題を計算することを回避する。

ナップザック問題の場合の動的計画法概略
  1. アイテムをなにか基準(講義では単位あたりの価値が大の昇順)を定義して、それに基づいてアイテムを並び替える。(さっきと違って、別に並び替えなくてもよい)。
  2. アイテムを縦軸(総アイテム数n)、ナップザックのサイズを横軸(ナップザックのサイズはb)としてn*bの表をつくる。
  3. 以降は関連資料のスライドシェア参照。
ナップザック問題への動的計画法

 いま、{item_i}{価値a_i, コストc_i})を入れることを考えるとして、

{
 f^*(i,j) =   \left\{
    \begin{array}{l}
       0    if \ a_i \ge j \\
       max\{  f^*(i - 1,j) ,  f^*(i-1,j-c_i) + a_i\}   otherwise
    \end{array}
  \right. \\
}
となる。計算時間はメモ化したものは{O(1)}で参照できる、ナップザックのサイズの回数だけアイテムが入るかどうかを試す、アイテムの数だけアイテムが入るかを試すことを考えて{O(1*n*b)}と予想できる。

英文の単語の最適配置問題への動的計画法

{
 s(i,j) := \frac{ L  -  \sum_{i}^j l_i }{ max\{j - i , 1\}}
}
{
 c(i,j) = 
  \left\{
    \begin{array}{l}
      \infty   if  s(i,j)\leq 0 \\
      (s - s(i,j))^2    if  s(i,j)\ge 0\\
      max\{ s - s(i,j),0\}^2 if   s(i,j)\ge 0 かつ j = n  
    \end{array}
  \right. \\
s := ideal \ space \ width
}

と問題を定式化したとして、さらに
{
 f^*(i) =   \left\{
    \begin{array}{l}
       0 \ \ \ if \ i = n + 1\\
       min_{i\leq j \leq n} \{ c(i,j) + f^*(j+1)\} \ \ otherwise
    \end{array}
  \right. \\
}
と記述。つまりこれは単語{w_i, \dots, w_n}で単語の配置を考える時の最適解で、{f^*(i)}をメモ化してそれを参照していくことで計算時間を削減する。

緩和問題

 { z_i \in \{0,1\}}の条件について、これを{0\leq z_i \leq 1}とした場合は目的関数、制約条件はともに連続な関数になるのだから線形計画問題として解くことができる。そうすると、Leonid Khachiyanによって線形計画問題を最悪の場合でも多項式時間で解くアルゴリズム*11が提案されているのだから、この問題は多項式時間でとける保証がある。ただしこのアルゴリズムは実用に向かない。

 問題の条件を緩和したのだから、実行可能解は”元の問題の実行可能解の領域を含んだ状態で”もとの領域より大きい。だからこの緩和した問題の最適解はもとの問題の最適解よりも「解を求める計算時間が少なく」「もとの最適解よりも”よい”値を最適解にもつ」ことがわかりこれは{upper-bound}の計算に用いることができる。

ナップザック問題の緩和問題

元の問題を
{
目的関数 :max   \sum_{i = 1}^n z_i a_i  \\
制約条件 :z_i   \in \{0,1 \}   \sum_{i=1}^n z_i c_i \leq b 
}
としたとき、これに対応する緩和問題(の一例)は
{
目的関数 :max \sum_{i = 1}^n z_i a_i  \\
制約条件 :z_i  \in [0,1 ]  \sum_{i=1}^n z_i c_i \leq b 
}
となる。この場合、{upper-bound}を計算したときの「はみ出た部分」の評価をせずに{upper-bound}を定めることができるからその分だけ削除できる部分問題が増加しうる。

緩和問題と双対問題と元の問題の関係

 一番大切なところなきがするけど省略(時間的に)。けど緩和(制約条件の内で変数を離散に限定するための条件を変数が連続になりなおかつ元の離散値をすべて含むような領域に変更)したら実行可能領域が増えるのだから最適解が良くなること、線形計画問題になるのだから多項式時間で解が求まることが分かる。

関連資料および参考文献

関連資料

分枝限定

組合せ最適化

参考文献

[0] 福島雅夫:『数理計画入門』, 朝倉書店, 1996.
[1] 酒匂貴市:最適化 平成20年1 月 3 日Project AFS