X - 九州大学

九州大学 工学部地球環境工学科
船舶海洋システム工学コース
システム設計工学 （担当：木村）
「数値最適化」
2006.07.20
木曜３限（13：00～14：30)
場所：船１講義室
1
【最適化とは？】
「最も良い状態になるように、手段や方法を考えて行動すること」
最適化以前に、
その目的である
「評価」を定義
することが重要
解候補－－－評価値
例） 製品のデザイン－－－売上げ高
船の形状－－－推進抵抗
モデル化に関連
一般的な最適化のプロセス：
新しい
解候補の生成
解候補の評価
以前の解候補の評価値を参考にして、より有望な解候補を生成することで
解候補を改善していく
2
【最適化とは？】
１）連続関数最適化問題
例： 石油を汲み上げるための油井を掘る
油を含む地層が最も地上に近い部分を掘ると利益が最大になる
予算の都合で、試掘の回数が制限される
→ 限られた試掘回数（６回）の範囲内で、最も地上に近い油層を探索する
解候補1
解候補5
解候補6
解候補2
解候補4
解候補3
x
真の最適解
利益関数 f(x)
探索で発見
した最良解
実際のxはn次元空間上を探索
3
【制約のない関数最適化】 1次元の最小化：ラインサーチ
黄金分割法
f(x)
１次元空間上
（直線）を探索
a
１） a < b < c の3点をとる
b
c
x
f(b) が最小なら、区間 [a,c] に必ず極小がある。
２） [b, c] の区間で新しい点 x をとり、f(x) を計算する。
３） f(b) < f(x) なら、 a < b < x に狭められた範囲に極小がある
f(x) < f(b) なら、 b < x < c に狭められた範囲に極小がある
4
【制約のない関数最適化】 1次元の最小化：ラインサーチ
黄金分割法
f(x)
１次元空間上
（直線）を探索
a
１） a < b < c の3点をとる
b
c
x
f(b) が最小なら、区間 [a,c] に必ず極小がある。
２） [b, c] の区間で新しい点 x をとり、f(x) を計算する。
３） f(b) < f(x) なら、 a < b < x に狭められた範囲に極小がある
f(x) < f(b) なら、 b < x < c に狭められた範囲に極小がある
5
【制約のない関数最適化】 1次元の最小化：ラインサーチ
黄金分割法
f(x)
１次元空間上
（直線）を探索
a
１） a < b < c の3点をとる
b
x
c
x
f(b) が最小なら、区間 [a,c] に必ず極小がある。
２） [b, c] の区間で新しい点 x をとり、f(x) を計算する。
３） f(b) < f(x) なら、 a < b < x に狭められた範囲に極小がある→ a,b,x をa,b,c に置換
6
f(x) < f(b) なら、 b < x < c に狭められた範囲に極小がある
【制約のない関数最適化】 1次元の最小化：ラインサーチ
黄金分割法
f(x)
１次元空間上
（直線）を探索
a
１） a < b < c の3点をとる
b
x
c
x
f(b) が最小なら、区間 [a,c] に必ず極小がある。
２） [b, c] の区間で新しい点 x をとり、f(x) を計算する。
３） f(b) < f(x) なら、 a < b < x に狭められた範囲に極小がある
f(x) < f(b) なら、 b < x < c に狭められた範囲に極小がある→ b,x,c をa,b,c に置換
【制約のない関数最適化】 1次元の最小化：ラインサーチ
黄金分割法
f(x)
a
b や x をどのように決めるか？
１次元空間上
（直線）を探索
z
w
x
b
c
x
1
次の新しい区間長は w+z または 1-w 。これらを等しくすると z = 1 - 2w
区間[b,c]に対するxの位置関係は、区間[a,c]に対するbの位置関係に等しい
1 : w = (1 – w) : z
黄金比
これらの連立方程式を解くと、
w 2 − 3w + 1 = 0
w=
3− 5
2
8
【制約のない関数最適化】 n次元関数の最小化
探索点での関数の値（と勾配の大きさ）に基づいて次の探索点を順次決定する
手順により、あたかも山の頂上を目指すように、関数の極大値（あるいは極小値）を
すみやかに発見する手法
１）勾配情報を用いる方法：
・最大勾配法（最急降下法）
勾配方向へ少しずつ解を改善
●ラインサーチと方向転換を繰り返す方法
・最適勾配法：勾配方向にラインサーチを行い、その極小点で再び勾配方向へラインサーチ
だいたいうまく行くが、２次関数において極めて効率の悪い場合がある
・共役勾配法：方向転換するとき、その点の
勾配方向だけでなく、今まで
下ってきた方向も考える
（関数の2階偏導関数を考慮に入れることと等価）
対象とする関数がn変数の２次形式である場合、
ラインサーチをn回繰返すことで最大値が求まることが保障される
振動が発生
２）勾配情報を用いない方法： 滑降シンプレックス法
9
たいした計算を要しないで手っ取り早く動くアルゴリズムが欲しい場合に最善の方法
【復習】 関数の勾配とは？
2変数関数
f ( x1 , x2 )
勾配ベクトル
パラメータベクトル
x = ( x1 , x2 )
⎛ ∂
⎞
∂
g = ∇f ( x1 , x2 ) = ⎜⎜
f ( x1 , x2 ),
f ( x1 , x2 ) ⎟⎟
∂x2
⎝ ∂x1
⎠
x2
∂
f ( x1 , x2 )
∂x2
( x1 , x2 )
勾配ベクトルg
∂
f ( x1 , x2 )
∂x1
x1
10
【復習】 「2次形式」と「共役」な方向とは？
f ( x1 , x2 , L , xn ) = X t AX + B t X + C
ただし
⎡ a11
⎡ x1 ⎤
⎢a
⎢x ⎥
X = ⎢ 2 ⎥, A = ⎢ 21
⎢ M
⎢M⎥
⎢
⎢ ⎥
x
⎣ an1
⎣ n⎦
a12
a22
M
an 2
2次形式
L a1n ⎤
⎡ b1 ⎤
⎢b ⎥
L a2 n ⎥⎥
, B = ⎢ 2⎥
⎢M⎥
O M ⎥
⎢ ⎥
⎥
L ann ⎦
⎣bn ⎦
N×N 対称行列 A に対して、２つの方向を表すベクトル p と q が
p t Aq = 0 を満たすとき、 p と q は A に対して共役であるという
参考：A を単位行列とすると、上の式は直交条件となり、共役性は直交性の拡張概念
11
【制約のない関数最適化】 共役勾配法
xi
i 番目の点の座標ベクトル
gi
点 xi における勾配ベクトル
pi
点 xi からラインサーチを
行う方向ベクトル
x0
まず最初は勾配の反対方向へラインサーチ
g i +1
p0 = − g 0
xi
の次の点 xi +1 を、ラインサーチで見つけた最小点とする
xi +1 = xi + α i pi
pi
pi +1
xi +1
xi
次の点 xi +1での新しい探索方向 pi +1は、前の探索方向 pi と共役であるようにする。
そのため、点 xi +1 での勾配と前の探索方向 pi とを結合する：
pi +1 = − g i +1 + β i +1 pi
ただし、重み係数 β i +1は次の式のどちらかで与える：
2
g i +1
g iT+1 g i +1
β i +1 = T
=
2
gi gi
gi
（Fletcher-Reeves法）
（Polak-Ribiere法）
( g i +1 − g i )T g i +1
β i +1 =
g iT g i
【共役勾配法が2次形式関数でｎ回のラインサーチで最適解を見つける理由】
次式を満たす n 個の
⎧const.
Z i AZ j = ⎨
⎩ 0
t
n 次元ベクトル Z1 , Z 2 , L Z n を考える：
for i = j
Z1 , Z 2 , L Z n
for i ≠ j
X を X = ∑ j =1α j Z j
は共役ベクトル集合
n
次に、ベクトル
と表すと、2次形式の評価式は
f ( x1 , x2 , L , xn ) = X t AX + B t X + C
n
⎞ ⎛ n
⎛ n n
⎞
2
t
t
t
= ⎜⎜ ∑∑ α iα j Z i AZ j ⎟⎟ + ⎜ ∑ α i B Z i ⎟ + C = ∑ α i Z i AZ i + α i B t Z i + C
i =1
⎠
⎠ ⎝ i =1
⎝ i =1 j =1
ここで
g i (α i ) = α i Z i AZ i + α i B t Z i
2
t
n
f ( x1 , x2 , L , xn ) = ∑ g i (α i ) + C
i =1
探索開始点から順次 Z1 , Z 2 , L Z n 方向
の関数の断面について α1 , α 2 Lα n の最小化
を行えばn回のラインサーチで最適解を得る
とおくと、上の式は以下のようになる：
これより、関数 f はそれぞれ α i
に対する n 個の関数の代数和
よって g1 (α1 ), g 2 (α 2 ), L g n (α n )
のそれぞれについて最小値を求める
13
【制約のない関数最適化】 滑降シンプレックス法
関数の勾配を計算せずに最小値を得る
N次元変数関数に対してN+1コ以上のシンプレックス（探索点） x j について以下を定義
1) 目的関数値が最大の点 xh
← 関数最小化なので、最悪の点
2) 目的関数値が2番目に大きい点 xs
3) 目的関数値が最小の点 xl
← 関数最小化なので、最良点
4) $i = h$なる点を除いた全ての x j の重心 x g
xh を以下の xr で置き換える：
xr = (1 + α ) x g − αxh ，ただし α > 0 は反射係数
(操作2：拡大) x g − xr 方向に沿って xr を以下の xe に置き換える：
γ > 0 は拡大係数
，ただし
xe = γxr + (1 − γ ) x g
(操作3：縮小) xh を以下の xc で置き換える：
xc = β xh + (1 − β ) x g ，ただし 0 < β < 1 は縮小係数
(操作4：収縮) シンプレックス全体を xl の方向へ半分に縮小する
xi = 0.5( xl + xi )
，ただし i = 1, L n + 1
(操作1：反射)
これらの操作を、次に述べる手順で組合せ，シンプレックスを更新する。
基本的に，最悪点をシンプレックス重心の逆側へ移動して関数値を小さくしていく。
α = 1, β = 0.5, γ = 2 が経験的に良い
14
滑降シンプレックス法の動作
最悪点
xh
〔スタート〕
・シンプレックスを構成する探索点の中で、
評価の最悪点、２番目に悪い点、最良点
シンプレックス重心座標を求める
xl
最良点
xg
2番目に
悪い点
xs
最悪点を除いた
全シンプレックスの重心
15
滑降シンプレックス法の動作
〔まずは反射操作〕
最悪点
xh
xl
最良点
xg
2番目に
悪い点
xs
〔反射〕
xr
16
滑降シンプレックス法の動作
最悪点
xh
xg
2番目に
悪い点
〔反射の結果が最良点より
良ければ、拡大操作〕
xl
最良点
xs
xr
xe
〔拡大〕
17
滑降シンプレックス法の動作
最悪点
xh
・拡大の結果が反射の結果
より悪ければ、最悪点を
反射点と置き換えて最初へ
xg
2番目に
悪い点
・拡大の結果が反射の結果
より良ければ、最悪点を
拡大点と置き換えて最初へ
xl
最良点
xs
xr
xe
〔拡大〕
18
滑降シンプレックス法の動作
最悪点
xh
xc
〔縮小〕
xg
2番目に
悪い点
〔反射の結果が最悪点より
悪ければ、縮小操作〕
xl
最良点
xs
xr
xe
19
滑降シンプレックス法の動作
最悪点
xh
xc
・さもなければ〔収縮〕操作
新たに生成したシンプレックス
〔縮小〕 の関数値を評価して最初へ
xg
2番目に
悪い点
・縮小の結果が最悪点より良
ければ、最悪点を縮小点と
置き換えて最初へ
xl
最良点
xs
xr
xe
20
滑降シンプレックス法の動作
最悪点
xh
〔収縮〕
2番目に
悪い点
xl
最良点
xs
21
シンプレックス全体を最良点の方向へ移動 （最良点は動かない）
滑降シンプレックス法の動作
最悪点
xh
〔反射の結果が最悪点よりまし
だが２番目に悪い点より悪い
場合、最悪点を反射点に置換
えて縮小操作〕
〔縮小〕
xl
最良点
xg
2番目に
悪い点
xs
xc
xr
xe
22
滑降シンプレックス法の動作
最悪点
xh
・縮小の結果が最悪点より良
ければ、最悪点を縮小点と
置き換えて最初へ
・さもなければ〔収縮〕操作
新たに生成したシンプレックス
〔縮小〕 の関数値を評価して最初へ
xl
最良点
xg
2番目に
悪い点
xs
xc
xr
xe
23
滑降シンプレックス法の動作
最悪点
〔反射の結果が２番目に悪い
点よりましだが最良点ほどで
はない場合、最悪点を反射点
に置き換えて最初へもどる〕
xh
xl
最良点
xg
2番目に
悪い点
xs
〔反射〕
xr
24
【ランダムサーチ・網羅探索による最適化】
モンテカルロ法とも呼ばれる
全数探索と同義な場合と区別する場合がある
・最も単純な方法で、どんな問題に対しても適用可能
・時間さえかければ、理屈の上ではいつかは最適解を発見できる
・今まで生成した解候補と評価値を参考にしないので、探索効率は悪い
探索性能の下限としての比較対象
ランダムサーチによる最適化のプロセス：
新しい解候補を
ランダムに生成
解候補の評価
それまで見つけた最良の解候補よりも高い評価値だった場合、
最良解候補を置き換える
25
【焼きなまし法(simulated annealing)による最適化】
統計物理学からのアナロジー：
・物質の分子結晶はいろいろな状態を遷移
・エネルギーの低い状態へ遷移しやすい傾向
・温度を下げるとエネルギー極小状態へ到達
・温度を一気に下げると局所的極小状態へ
・温度をゆっくり下げると大域的極小状態へ
「状態」
→ 解候補
「エネルギー」 → 評価値
エネルギー最小状態＝最適解
とした最適化のプロセス
システムがある状態 xi になっている確率
自由エネルギー
exp( − f ( xi T ))
P( xi | T ) =
∑ exp(− f ( x) T )
x∈S
ギブス分布（ボルツマン分布）
Tは温度パラメータ
f ( x)
状態空間
x
エネルギーの低い状態ほど確率が高い
特徴：
１）時間をかければかけるほど解候補が改善され、
十分時間をかければ最適解の発見が保障される
２）連続関数最適化，および組合せ問題最適化の両方に適用できる
26
【参考】
大腸菌の走性
2～3秒間の直線的運動と、
0.1秒ほどの方向転換の繰り返し
化学的濃度が低い（要するにエサが少ない）場合には活発に動き回り、
濃度の高い領域ではあまり動かないようにすることでエサのある場所へ集まる
焼きなまし法によく似た探索をしている
27
【焼きなまし法(simulated annealing)による最適化】
処理手順
１） 初期解候補 x を生成
２） 何らかの確率分布によって解候補 x から次の解候補
（ランダムに x の近傍を選ぶことが多い）
３） 次の確率に従って乱数Zを生成：
ただしZは0または1の値をとる
ただし温度関数
４）
Z =1
のとき
T (t ) =
k
ln(t + 2)
x' を生成
1
P( Z = 1) =
1 + exp(( f ( x' ) − f ( x)) T (t )
x を x' に置換える、
５）時刻 t を1進めて２）から繰返す
これでは収束が遅すぎて非実用的
k/t の関数にすることもある
28
「近傍」の解候補の生成方法
１） 数値最適化問題の場合： 解候補 x は連続値パラメータなので、
中心値 x , 分散σ 2 の正規分布によって生成： x ' = N ( x, σ 2 )
分散のパラメータは設計者が問題に応じて適当に与える
時間とともに小さくしていく場合もある
x' x
２） 組合せ最適化の場合：
まず解候補の「近傍」を定義する必要がある
→ 問題毎に個別に設定： 実行可能な表現になるよう工夫を要する
例１） ＴＳＰ
順列中の都市を入れ替える
順列中の部分順列を移動するなど
例２） ＳＡＴ
論理変数の論理を反転させる
ＡＢＣＤＥＦ → ＡＥＣＤＢＦ
ＡＢＣＤＥＦ → ＡＣＤＥＢＦ
１０１１１０１ → １００１１０１
29
【制約つき関数最適化】
目的関数
f (x)
を制約条件
ただし
非線形計画法（nonlinear programming）
x = [x1 , x2 , L xn ]
hi ( x) = 0
ただし
i = 1, 2, L m
（等式制約条件）
g j ( x) ≥ 0
ただし
j = 1, 2, L l
（不等式制約条件）
のもとで最小化する問題。
等式制約条件だけの場合：まずラグランジュの未定乗数法を試みる
ラグランジュ関数
この関数が
L(x, λ ) = f ( x) + ∑i =1 λi hi (x)
x* , λ *
m
を導入する。
において局所的最小点となるための必要条件が
∂
∂
∂
m
L(x* , λ * ) =
f (x* ) + ∑i =1 λ*i
hi (x* ) = 0
∂ xj
∂ xj
∂ xj
∂
∂ λj
L( x* , λ * ) = hi (x* ) = 0
であることを利用し、この連立方程式を解いて
30
最小点 x * を求める
【制約つき関数最適化】
目的関数
を制約条件
f ( x)
ただし
非線形計画法（nonlinear programming）
x = [x1 , x2 , L xn ]
hi ( x) = 0
ただし
i = 1, 2, L m
（等式制約条件）
g j ( x) ≥ 0
ただし
j = 1, 2, L l
（不等式制約条件）
のもとで最小化する問題。
不等式制約条件の扱い：ペナルティ関数法
ペナルティ関数
p ( x)
f ( x) + p ( x)
を導入し、目的関数を
に変換して、
制約のない関数最適化問題として解く
ペナルティ関数の代表例： SUMT (sequential unconstrained minimization technique)
1
1
p (x) = α t ∑i =1
+
g i ( x)
αt
l
∑ j =1 h j (x)
m
2
ただし α t は最初大きな値に設定し、探索とともに小さな値にしていく
31