The 19th Annual Conference of the Japanese Society for Artificial Intelligence, 2005

2F2-05

自己組織化サンプリングの提案とその最適化手法-多峰性関数への適用Self-Organaized Sampling as an Optimization Method:
An Approch towards Multimodal Functions
比護 貴之

高玉 圭樹

Takayuki Higo

Keiki Takadama

東京工業大学 大学院総合理工学研究科 知能システム科学専攻
Department of Computational Intelligence and System Science
Interdisciplinary Graduate School of Science and Engineering
Tokyo Institute of Technology
This paper improves Self-Organized Sampling (SOS) method, which is our proposed optimization method for
multimodal functions, by using new parameter instead of the inverse temperature parameter on Bolztman distribution that causes the problem of an appropriate parameter setting. Through intensive simulations , we have foundan effectiveness of the proposed parameter. Concretely, SOS with the proposed parameter can optimize both theoriginal objective function and the linearly scaled function using the same parameter value, while SOS with theinverse temperature parameter cannot always optimize both functions using the same parameter value.
Keywords: global optimization, multimodal function, bolzman distribution, annealing, inverse temperature

1.

はじめに案し、それを SOS に組み込む方法を説明する。4 章では実験を行い、5 章で考察を行う。そして、6 章で本研究をまとめる。

近年最適化問題の中でも特に多峰性関数最適化が盛んに研究されている。多峰性関数最適化問題に対する有効なアプロー
チとして、統計量の遺伝に着目した遺伝的アルゴリズム (Genetic Algorithms: GA)[Kita 1999] と、GA と関連性の強い分布推定アルゴリズム (Estimation of Distribution Algorithms:EDA)[Larranaga 2001] が挙げられる [佐久間 2003]。統計量の遺伝に基づく GA では、理想的な交叉アルゴリズムであるSPX[樋口 2001] が提案され有効性が確認されているが、多峰性関数最適化においては問題が残る。一方、EDA では探索点の背景にある確率分布を推定し、その分布に従うサンプルを発生させるという設計指針に基づいたアルゴリズムが研究されている。
これらに対して、著者らの提案する自己組織化サンプリング(Self-Organized Samling: SOS) は、サンプルの期待評価値とエントロピーを最大化する点で最適なアルゴリズムであり、一部の多峰性関数最適化問題に対して有効性が確認されつつある[比護 2004]。しかし、この手法には SOS のパラメータである逆温度を適切に設定することが難しいという欠点がある。逆温度は探索点の収束度合を決めるパラメータであるため、最適化効率に大きな影響をあたえる。しかしながら、ある逆温度における探索点の収束度合は目的関数に依存する上、見積もることが困難であるため、試行錯誤を要する。
そこで、本研究では逆温度の問題を解決するために、逆温度の代替パラメータと、それを用いた計算アルゴリズムを提案する。また、この改良により SOS が目的関数の線形スケーリングに対して不変なアルゴリズムとなることを示す。
以下 2 章では、SOS とそれを用いた最適化手法及びその問題点についてまとめる。3 章では逆温度の代替パラメータを提2.

自己組織化サンプリング (SOS)

2.1

SOS の特徴

SOS の特徴としては次の３つが挙げられる。
• 目的関数のボルツマン分布に従うサンプリング。
• 目的関数を階段関数で近似する。
• 確率密度和を一定にするような分割を行う。
これらについて以下で簡単な説明を行う。まず、ボルツマン分布は次のような分布である。

exp(−f (x)β)Z(1)

exp(−f (x)β)dx(2)p(x) =
Z=Z

f(x) はエネルギー関数、β は逆温度、Z は分配関数と呼ばれる。
R
ボルツマン分布の特徴はエネルギー関数の期待値 f (x)p(x)dx
Rを最小化し、エントロピー p(x) log p(x)dx を最大化する分布であり、そのトレードオフは逆温度で決まる [深尾 1987]。エネルギー関数に最小化問題の目的関数を用いる場合（最大化問題の場合は目的関数を-1 倍する）、ボルツマン分布の特徴は最適化において良い性質と考えられる。
一般的にボルツマン分布を計算することは困難であり、既存手法にマルコフ連鎖モンテカルロ法 [伊庭 2003] が挙げられるがいくつかの問題がある。そこで、SOS は領域を単体に分割し、単体内部では確率密度が一定となるようなモデルで近似する。このモデルは、ボルツマン分布のエネルギー関数を階段関数とすることで構成できる。
近似においてはその精度が問題となる。そこで、それぞれの単体領域における確率密度の和を一定化するように、単体へ連 絡 先: 比 護 貴 之 , 東 京 工 業 大 学 大 学 院 総 合 理 工 学 研
究 科 知 能 シ ス テ ム 科 学 専 攻, 〒 226-8502 神 奈 川
県横浜市緑区長津田町 4259 番地, Tel:045-924-5216,
higo@cas.dis.titech.ac.jp1

The 19th Annual Conference of the Japanese Society for Artificial Intelligence, 2005の分割を変えることで、確率密度の高い領域を正確に近似し、確率密度の低い領域を粗く近似する。

2.2

2.3

アルゴリズムまず、初期化について説明する。SOS は線形制約を想定しており、それは凸多面体となる。そこで、凸多面体の頂点に疑似探索点を配置し、内部に探索点を適当に配置する。探索点と疑似探索点の違いは、探索点が最適化の過程でその位置を修正するのに対し、疑似探索点は位置を変えない。そして、探索点と疑似探索点に対してドローネ三角形分割 [今井 2003][譚 2001]を適用する。ドローネ三角形分割により探索点及び疑似探索点の間に枝が定義され、探索領域が単体に分割される。単体というのは d 次元の空間において d+1 個の頂点を持つ多面体を指す。
次に、最適化過程について説明する。最適化過程ではランダムに探索点を選び、それを適切な位置、つまり全ての単体の確率密度和が一定となるような位置へ修正することを繰り返す。
ある探索点を適切な位置へ修正するためには、その探索点に隣接する単体から M 個のサンプルを発生させ、そのサンプル情報を用いて次のように更新を行う：(1) 探索点に隣接する単体
si の重心 c(si ) を求める。重心とは頂点の座標を平均したものを指す。(2) サンプルから単体 si の平均確率密度 P (si ) を求める。平均確率密度は次のように計算する。

X 1βt = β0 rt

j=1Mexp(−f (xij )β)3.

X

si ∈Nxp(si )
c(si )p(sk )s ∈NxP

接近率の提案

逆温度の問題を解決するために、次のような代替パラメー
タ γ を提案する。以下では目的関数 f (x) の最大化問題を想定して話を進める。最大化問題ではボルツマン分布の分子がexp(f (x)β) となることに注意されたい。

(3)

ただし、xij は単体 si の内部から一様に発生させた M 個のサンプルであり、f (x) は最小化問題の目的関数とする。単体から一様にサンプリングする方法は [樋口 2001] を参照されたい。
(3) 探索点を次のように更新する。

x=(5)

ただし、t は時間（関数評価回数）を表す。この方法の有効性は経験的に知られているが、初期値 β0 と変化率 r をどのように設定すべきかという問題がある。直感的には、ボルツマン分布が一様分布に近い逆温度を初期値とし、ボルツマン分布がある一点のみで確率 1 である分布とみなせるような逆温度を終了値とし、そこから変化率 r を見積もればよいのだが、初期値と終了値は目的関数に依存している上に、見積もるのが困難である。また、初期値と終了値の範囲は 0 から無限大の範囲にあり、あまりにも大きい逆温度を計算機で扱うのは難しい。
このような問題が SOS に内在している。

Mp(si ) =アニーリングアニーリングとは、逆温度の低い状態から高い状態へゆっくり変化させる操作を意味し、SOS では局所解へ陥ることを防ぐ役割をする。逆温度のアニーリングは焼き鈍し法 (SimulatedAnnealing: SA)[白石 2002] でも登場する。SA では指数的に逆温度を変えるアニーリングがよく利用され、それは次のように書ける。

γ==exp(f (xu )β)exp(f (xu )β) + exp(f (xv )β)1
1 + exp(δβ)

(6)
(7)

ただし、xu , xv は適当に決めるとし、δ = f (xv ) − f (xu ) とする。γ はある２点の評価値の比を表している。これは SOS のアルゴリズムにおける探索点の更新において、ある探索点を更新する際に隣接する単体の重心への近付き方の比と解釈できる。その様子を図 2 に示す。以下ではこのパラメータ γ を接近率と呼ぶことにする。

(4)k

ただし、Nx は探索点 x に隣接する単体の集合である。
最後に、アルゴリズムの手順を図 1 に示す。

• 初期化
1. 線形制約が作る凸多面体の頂点に疑似探索点を配置する。
2. 凸多面体内部に探索点を配置する。
3. 探索領域をドローネ三角形分割する。
図 2: 接近率の解釈

4. xi の隣接単体集合 Nxi を求める。
• 最適化

逆温度は式 (7) から接近率を用いて次のように書ける。

1. 探索点 xp をランダムに選ぶ。

β

2. 隣接単体集合 Nxp に含まれる単体の内部から一様に M 個サンプルを発生させる=

11log ( − 1)
δγ(8)f (xv ) < f (xu ) を仮定すると、δ < 0 となり、この時 γ = 12 でβ = 0、γ = 1 で β = ∞ になることがわかる。従って、γ を 12から１まで動かすことでアニーリングが行える。また、γ = 1の時 β = ∞ となることから探索点が収束することは明らかである。つまり、γ = 1 のときは、少しでも評価値に差があれば探索点は評価値の高い方へ完全に移動するため、収束することが保証される。ただし、(準) 最適解が二つありその評価値が全く同じ場合は例外である。

3. 探索点 xp を式 (4) に従って更新する。
4. １に戻る
図 1: SOS アルゴリズム

2

The 19th Annual Conference of the Japanese Society for Artificial Intelligence, 2005
次にこれをアルゴリズムに組み込む方法を説明する。逆温度β が必要となるのは式 (4) の計算である。式 (4) は各単体の重心ベクトル c(si ) を単体の平均確率密度 p(si ) で平均することを意味しており、p(s1 ) : p(x2 ) : . . . : p(xn ) の比がわかればよ4.2p(s )い。つまり、各 i,j に対して p(sji ) を計算できればよい。これは次のように計算できる。

p(sj )
p(si )==

実験 1 逆温度を 10−2 から 1010 まで式 (5) で変化させる。

exp(f¯j β)
exp(f¯i β)exp((f¯j − f¯i )β)(9)

実験 2 目的関数を 10−6 倍した関数 g(x) = 10−6 f (x)に対して、逆温度を 10−2 から 1010 まで式 (5)で変化させる。

(10)=1f¯j − f¯iexp(log( − 1))δγ=(¯設定

先ほど定義した関数に対して、探索点数 15、サンプル数 1の SOS を適用し 6 万回のサンプリングを行う。ただし、サンプル数というのは単体内部から一度に発生させるサンプルの数である。この条件のもと以下で示す 4 つの実験を行う。

(11)

実験 3 接近率を 0.5 から 1 まで式 (16) で変化させる。

¯fj −fi1− 1) δγ(12)実験 4 目的関数を 10−6 倍した関数 g(x) = 10−6 f (x)に対して、接近率を 0.5 から 1 まで式 (16) で変化させる。

ただし、f¯ はサンプルから求まる単体の平均評価値である。また、あらかじめ p(si ) を大きい順に並び替え、最も大きいものを 1 と置き、大きい順に比を求めていけば値は発散しない。たf¯ −f¯だし、δ が f¯j − f¯i と比べて小さくなり過ぎると j i が無限γt = 0.5 +δ大になるので注意されたい。
最後に、便宜上 δ = f (xu ) − f (xv ) と置き直せば (符号を反転させれば) 式 (12) は次のように書ける。
f¯i −f¯jp(sj )1= ( − 1) δp(si )
γ

実験

4.1

問題

(13)

δt+1 = δt + α((f¯f st − f¯sec ) − δt )f (x) = 0.5(cos((x − c) ∗ 400) − 1) + 10g(x)c1−x1−c(14)x<cc<x(15)ただし、c = 0.3 で、0 < x < 1 とする。この関数は図 3 のようになる。本実験ではこの関数の最大化問題を考える。
10x8f(x)6420-200.10.20.30.40.5x0.60.70.80.9

結果

実験 1、2 の結果を図 4、5 に示す。図の縦軸は探索空間、横軸は逆温度を表し、最適化過程でサンプリングした点をプロットする。実験 1 では、逆温度が 10−2 から既に収束が始まり、102 でほぼ収束しているのがわかる。一方、実験２では逆温度が 104 辺りから収束が始まり、107 でほぼ収束しているのがわかる。実験 1 と 2 で収束の始まりと終わりの逆温度が異なるのは、目的関数を 10−6 倍しているためである。従って、逆温度の初期値や終了値は関数のスケールに依存すると言える。

次のような関数を用いる。

g(x) =(17)ただし、α = 0.05 とする。

4.3

½ x(16)

接近率における δ の値は、探索点を更新する際に平均確率密度の 1 番大きい単体と、2 番目の単体の平均評価値 f¯f st , f¯secを使って次のように更新する。

以後 δ はこの定義を用いることにする。

4.

0.5t60000

110.90.80.70.60.50.40.30.20.1
010-2 100

102 104 106 108 1010
inverse temperature

図 4: 逆温度を用いた場合。
次に実験 3、4 の結果を図 6、7 に示す。横軸が接近率になっていることに注意されたい。実験 3、4 を比較すると結果に大きな差は見られないので、目的関数のスケールに依存しないと予想される。接近率が 0.5 の付近では一様な探索をしており、

図 3: テスト用の関数3

10.90.8
0.70.60.50.40.30.20.1
010-2 100

10.90.80.70.60.50.40.30.20.10xx

The 19th Annual Conference of the Japanese Society for Artificial Intelligence, 2005

102 104 106 108 1010
inverse temperature

0.50.60.70.80.91approach rate
図 7: 目的関数を 10−6 倍し、接近率を用いた場合。

接近率を 1 に近づけるに従い徐々に x=0.3 付近に収束していることがわかる。これより、目的関数によらず接近率を 0.5 から 1 へゆっくり変化させることで、アニーリングが実現でき、接近率は逆温度に比べて利便性が高いことがわかる。

[Wolpert 1995] のような結論も出ており慎重な議論が必要である。

x

図 5: 目的関数を 10−6 倍し、逆温度を用いた場合。

10.90.80.70.60.50.40.30.20.106.

本論文では SOS で扱うボルツマン分布に表れる逆温度が応用上問題となることに着目し、代替パラメータである接近率を提案した。そして、逆温度と接近率の比較実験によりその有効性を示すとともに、接近率を用いた SOS が目的関数の線形スケーリングに対して不変であることを示した。
今後の課題としては、初期化における三角形分割の計算量が次元の指数オーダーとなるため、これに代わるアルゴリズムを提案することが先決である。また、接近率のアニーリングを調整するアルゴリズムの提案も課題の一つである。

謝辞
0.5

0.6

0.7

0.8

0.91

本研究の一部は文部科学省の科学研究費補助金 (基盤研究(B)，課題番号 17360424) の支援によって行われた．

approach rate

参考文献

図 6: 接近率を用いた場合。

5.

考察

5.1

線形スケーリングに対する不変性[深尾 1987] 深尾 毅: 分散熱力学システム論, 昭晃堂 (1987).
[比護 2004] 比護 貴之, 高玉 圭樹: “適応的探索領域に基づく最適化手法の提案：多峰性関数へのアプローチ”, 計測自動制御学会, システム・情報部門 学術講演会 2004, 95/100 (2004).
[樋口 2001] 樋口 隆英, 筒井 茂義, 山村 雅幸: “実数値 GA におけるシンプレクス交叉の提案”, 人工知能学会論文誌, 16, 1Q, 147/155(2001).
[伊庭 2003] 伊庭 幸人: ベイズ統計と統計物理, 岩波書店 (2003).
[今井 2003] 今井 浩, 今井 圭子: 計算幾何学, 共立出版 (2003).
[Kita 1999] Kita, H. and Yamamura, M.: “A Functional Specialization Hypothesis for Designing Genetic Algorithms”,
Proc.1999 IEEE International Conference on Systems,
Man, and Cybernetics, 579/584 (1999).
[Larranaga 2001] Larranaga P., Lozano,J.A.: Estimation of Distribution Algorithms: A New Tool for Evolutionary Optimization. Kluewr Academic Publishers, Boston (2001).
[佐久間 2003] 佐久間 淳, 小林 重信: “確率分布推定に基づく実数値GA の新展開”, 人工知能学会誌, Vol.18, No.5, 479/486 (2003).
[白石 2002] 白石 洋一: 組合せ最適化アルゴリズムの最新手法, 丸善株式会社 (2002).
[譚 2001] 譚 学厚, 平田 富夫: 計算幾何学入門, 森北出版 (2001).
[Wolpert 1995] Wolpert, D.H. and Macready, W.G.: “No Free
Lunch Theorems for Search”, Technical Report SFI-TR-9502-010, Santa Fe, NM (1995).

実験より、接近率を用いた SOS は目的関数の線形スケーリング、すなわち目的関数 f (x) を af (x) + b と置き直しても、もとの目的関数に適用した場合と同じように動作すると予想される。これは式 (13) に、δ = f (xu ) − f (xv ) を代入し、関数f (x) を af (x) + b と置き換えることで値が変わらないことかPtら証明できる。より一般的に δt =w (f (xut ) − f (xvt ))i=1 tPを考えよう。ただし、 t wt = 1 である。これは、アルゴリ
ズムの過程で得られた xut 、xvt に対して、δ を計算し、w で加重平均したものである。これを式 (13) に代入し目的関数を置き換えれば同様に値が変わらないことがわかる。そして、本実験で用いた δ の更新式 (17) もこれに含まれる。

5.2

おわりに問題点

本論文の実験では接近率を線形的に変化させたが、この方法が最適である保証はない。従って、どのように接近率を調整するかは今後の課題となるが、No Free Lunch Theorem

4