Transcript JNNS-DEX-SMI-玉川 公開講座 「交換モンテカルロ法と

JNNS-DEX-SMI-玉川 公開講座
「交換モンテカルロ法とその応用」
東京工業大学大学院 総合理工学研究科
知能システム科学専攻 博士課程2年
永田賢二
交換モンテカルロ法とは？
•マルコフ連鎖モンテカルロ法（MCMC法）のひとつ
•乱数を用いて、確率分布を再現するための一群の手法
•正規分布など、性質のわかっている分布だけでなく、離散・連続を問わず、
様々な分布に適用できる。汎用性が高い。
•交換モンテカルロ法は、従来のMCMC法の改良アルゴリズム[Hukushima,96]
＜MCMC法の主な目的＞
•サンプリング
確率分布
p(w)
•期待値計算
1 m
g (w(t ) )   g (w) p(w)dw

m t 1
MCMC法の応用例
＜ベイズ統計＞
p( x | w)：確率モデル
データ
p(w) ：パラメータの事前分布
X が与えられたもとでのパラメータ w の条件付き確率（事後分布）
p( X | w) p(w)
p(w | X ) 
Z   p( X | w) p(w)dw
Z
＜統計物理＞
ギブス分布・カノニカル分布における期待値計算
exp(E(w))
p(w) 
Z ( )   exp(E(w))dw
Z ( )
E(w)：エネルギー関数

：温度の逆数（逆温度）
Outline

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析
Outline

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析
マルコフ連鎖モンテカルロ法
d次元空間上の点 w が従う確率分布 p(w) が与えられているとする。
また、点
w の各成分は連続値をとるものとする。
１．確率分布 p(w) に従う点をサンプリング
w  w
(t )
(1)
, w(2) ,w(m) 
２．wの関数 g(w) の確率分布 p(w) についての期待値の計算
1 m
E[g (w)]   g (w) p(w)dw   g (w(t ) )
m t 1
メトロポリス法
（例）以下の確率分布 p(w) に従うサンプル生成
３．
２．
現在の点
密度の比較により、次の点
wを生成する。
w(t 1) を決める。
１． w
の初期値w(tw
) から、以下の式で候補
( 0) を設定する。
)p(w
w(t))  
(t )()w
p(w(t ) p) (wp
rで:w平均0の一様乱数、正規乱数など
w(t w
 w r  p(w)

1
)
(t 1)
確率1 r で w(t 1)  w(t )
p(w(t ) )
確率
pp(w
(w
(11
))
p(wr)
p(w(0) )
0
w(0)w(2) w(1)
w
メトロポリス法
（例）以下の確率分布 p(w) に従うサンプル生成
３． 密度の比較により、次の点 w(t 1) を決める。
p(w(t ) )  p(w)
r で w(t 1)  w
確率1 r で w(t 1)  w(t )
確率
p(w)
r
p(w(t ) )
1
r
0
w(0)w w(1)
w( 2)
w
メトロポリス法のアルゴリズム
＜確率分布 p(w) に従うサンプリング・アルゴリズム＞
１．
wの初期値w(0) を適当に設定する。
２．現在の点 w(t ) から、以下の式で候補
w  w(t )  
w
を生成する。
 : 平均0の一様乱数、正規乱数など
３．密度の比較により、次の状態 w(t 1) を決める。
p(w(t ) )  p(w)
w(t 1)  w
p(w(t ) )  p(w)
確率
４．ステップ２に戻り、繰り返す。
r で w(t 1)  w
確率1 r で w(t 1)  w(t )
p(w)
r
p(w(t ) )
ステップサイズ
ステップサイズ：候補を選ぶ際の範囲の大きさ
（例）以下の２次元の確率分布からのサンプリング
・大きすぎると・・・
・大きすぎると、ほとんどの候補が採択されなくなる。
・小さすぎると、一回の更新が少ないため、遠くに行きにくい。
・小さすぎると・・・
ステップサイズ
(例) 右の目標分布から1000個のサンプルを生成
 :２次元の一様分布からランダムに選ぶ。
ステップサイズ：0.05
ステップサイズ：0.5
ステップサイズ:5
・実際には、採択される割合が40%～60%程度になるように設定
・要素ごとに更新するのも一つの手。
「メトロポリス法」のまとめ

確率的に「候補」を選んで、それを採用するかどうかを、確率
的に決定する。

目標分布の情報は、密度の比のみしか用いないため、密度
さえ計算できれば、どんな分布にも適用できる。規格化され
ていなくても大丈夫。

ステップサイズの設定は、アルゴリズムの効率アップのため
に、重要。大きすぎず、小さすぎず。要素ごとの更新を考えて
もいいかも。
Outline

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析
マルコフ連鎖
マルコフ連鎖：直前の点 w(t 1) のみに依存して、次の点
w(t ) を決定する。
w(t 1)  w(t )  w(t 1) 
・遷移確率
 (w  w)： 点 w から点 wに移る確率
（メトロポリス法の場合）
Case1:
 :[-D,D]の範囲の一様分布からランダムに選ぶ。
 (w  w)  0
p(w)
1
Case2:  (w  w) 
2D
1 p(w)

Case3:  (w  w ) 
2D p(w)
w w w
w
マルコフ連鎖の原理
＜遷移確率  (w  w) が満たすべき条件＞
１．詳細つりあい条件
p(w) (w  w)  p(w) (w  w)
確率分布 p(w) に従う点がたくさんある状況を考える。
（左辺）： w から
w に移る個数
（右辺）： wから
w
に移る個数
任意の２つの位置での
「流入」と「流出」がつりあっている。
「確率分布 p(w) を不変にする」
それぞれの点を更新
 (w  w)
p(w)
p(w)
 (w  w)
w
w
マルコフ連鎖の原理
＜遷移確率  (w  w) が満たすべき条件＞
２．エルゴード性
任意の2つの点 w と w の間の遷移確率がゼロでないか、
有限個のゼロでない遷移確率の積で表すことができる。
・何回かの更新で、どこへでも
到達することが可能である。
・どんな初期値から始めても
唯一の分布に収束する。
w
w
メトロポリス法における詳細つりあい条件
p(w) (w  w)  p(w) (w  w)
（先のメトロポリス法の場合）
 :[-D,D]の範囲の一様乱数
１．w  w
p(w)
 D の場合
 (w  w)   (w  w)  0
２．
p(w)  p(w) の場合
1
2D
1 p(w)
 (w  w) 
2D p(w)
 (w  w) 
w w
 (w  w) p(w)

 (w  w) p(w)
w
MCMC法のアルゴリズム

遷移確率の満たすべき条件は緩くて、一意に決定できない。


詳細つりあい条件
エルゴード性
MCMC法のアルゴリズムは、たくさん存在する。
（例）






メトロポリス法
メトロポリス・ヘイスティングス法
ギブスサンプラー、熱浴法
独立サンプラー
ハミルトニアン・モンテカルロ法
「マルコフ連鎖の原理」のまとめ

マルコフ連鎖



直前の状態にのみ依存して、次の状態が決まる系列
遷移確率で特徴づけられる。
マルコフ連鎖の基本原理

詳細つりあい条件


エルゴード性


任意の２つの位置での「流入」と「流出」がつりあっている。
有限回のステップで、任意の２点間を行き来できる。
条件は緩く、いろいろなアルゴリズムが存在する。
Outline

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析
遅い緩和の問題
メトロポリス法の基本は、「少し変えて、選ぶかどうかを確率的に決める。」
ある確率分布に対しては、ものすごく効率が悪くなってしまう。
（例１）密度の高い領域が、いくつもあり、互いに離れている場合
（多峰性のある確率分布）
・ある領域から、他の領域に到達するには、
密度の低い領域を通る必要がある。
⇒サンプリング効率の悪化
遅い緩和の問題
（例１）基底状態が一点ではなく、次元を持った集合になっている。
（ベイズ学習でみられる問題）
exp(E(w))
p(w) 
Z ( )
基底状態：エネルギー関数 E(w) を最小にする点 w のこと
＜一点の例＞
＜集合の例＞
拡張アンサンブル法
確率分布によっては、
遷移確率が著しく小さくなる。
サンプリング精度が悪くなり、
期待値の評価に影響を与えてしまう。
＜実質的なエルゴード性の破れ・遅い緩和の問題＞
＜拡張アンサンブル法＞ 上記の問題を解決する一群の手法
確率分布を拡張したり、混合したものを考える。
・マルチカノニカル法
・シミュレーテッド・テンパリング法
・交換モンテカルロ法 [Hukushima-Nemoto,96]
交換モンテカルロ法のアイデア（温度の導入）
・確率分布 p(w) に対して、 p (w) 
ギブス分布の場合：
＜高温状態＞
p(w)
p (w)  exp(E(w))
＜低温状態＞
・エネルギーの低い点は探せない
・大域的に行き渡れる
・常にエネルギーの低い点にいる
・局所領域に留まりやすい
サンプリング中に
温度を上げ下げする。
＜問題＞
・温度を上げ下げする過程で、詳細つりあい条件を破ることになるので、
目標分布からのサンプリングの保証がなくなる！
交換モンテカルロ法[Hukushima,96]
1
目標分布： p(w)  exp E(w)
Z
{w}  {w1,, wL}
拡張された確率分布 p{w}
L
L
1
p{w}   pl (wl )  
exp l E(wl ) （ 
l 1
l 1 Z (l )
: 逆温度）
＜アルゴリズム＞
１．（通常の更新）それぞれの確率分布について、状態の更新
wl  wl
２．隣り合った分布間で、状態の交換を行う。
wl , wl 1 wl 1, wl 
交換モンテカルロ法の詳細つりあい条件
＜詳細つりあい条件＞
p{, wl , wl 1,} {, wl , wl 1,}  {, wl 1, wl ,}
 p{, wl 1, wl ,} {, wl 1, wl ,}  {, wl , wl 1,}
 {wl , wl 1}  {wl 1, wl } pl (wl 1 ) pl 1 (wl )

 {wl 1, wl }  {wl , wl 1} pl (wl ) pl 1 (wl 1 )
 exp (l 1  l )(E(wl )  E(wl 1 ))  r
r
<Case2> r
<Case1>
:小
:大
<Case3> r  1
交換の採択確率
pl (wl 1 ) pl 1 (wl )
r
 exp (l 1  l )(E(wl )  E(wl 1 ))
pl (wl ) pl 1 (wl 1 )
1.メトロポリス型
u1  min(1, r)
(1)
交換前
交換後
必ず交換する。
(2)
交換前
交換後
確率
2.熱浴型
u1  1
r
u1  r
で交換する。
r
u2 
1 r
（交換後）
0
（交換前）
u2
1
交換モンテカルロ法のイメージ
<メトロポリス法>
p(w)
<交換モンテカルロ法>
p1(w1)
p2 (w2 )
p3 (w3 )
p4 (w4 )
１．（通常の更新）
メトロポリス法により、状態の更新
２．（状態の交換）
隣り合った分布間で、状態の交換
wl , wl 1 wl 1, wl 
＜交換の採択確率＞
u  min(1, r )
r  exp (l 1  l )(E(wl )  E(wl 1 ))
交換モンテカルロ法の挙動（イメージ）
＜前に出した例では・・・＞
高温
低温
・低温での点が、高温に移ることで、大域的なサンプリングが可能に。
・詳細つりあいを満たしているので、サンプリングの保証つき。

交換モンテカルロ法の
実験結果の例
右の確率分布から10000個のサンプルを生成
メトロポリス法
交換モンテカルロ法
ベイズ学習での交換モンテカルロ法
＜混合正規分布モデルにおけるベイズ学習＞
推定
<学習データ>
1000個の3次元データを生成
<正規分布の数>
データ生成: 4個
学習モデル: 10個
汎化誤差：真の構造と予測結果の相違
アルゴリズム
汎化誤差
Gibbs sampler 0.011188 0.003249
0.009809 0.002989
0.010500
理論値(上限)
交換法
その他の応用例

ポリマーの構造推定
 A. Irb
ack, E.Sandelin, J. Chem. Phys.110 (1999).

タンパク質の立体構造推定
 A. Mitsutake,Y. Sugita, Y. Okamoto,Biopolymer60 (2001).

スピングラス・シミュレーション
 K. Hukushima,K. Nemoto,J. Phys. S
oc. Jpn.65 (1996).

組み合わせ最適化問題
wski, Int. J. Mod. Phys.C 9 (1998).
 K. Pinn, C. Wieczerko
 K. Hukushima,Comp. Phys
. Comm.147 (2002).
「交換モンテカルロ法」のまとめ

遅い緩和の問題



密度の高い領域が複数存在し、互いに離れている場合
エネルギーの基底状態が、次元をもった集合になっている場合
交換モンテカルロ法の原理


温度パラメータを導入することで、大域的なサンプリングが可能に
同時分布の詳細つりあい条件を考える
Outline

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析
交換法が有効に働くには・・・
＜交換が、ある程度の確率で行われている。＞
交換が行われないと、通常のメトロポリス法を行っているのと同じ。
r  exp (l 1  l )(E(wl )  E(wl 1 ))
＜ E(w) のヒストグラム＞
低温
温度パラメータの値によって
交換の頻度が決まる。
高温
＜温度パラメータの設定＞
・各  の間隔は？
・温度パラメータの総数は？
温度パラメータ設定の際の基準例
＜平均交換率＞
r  exp (l 1  l )(E(wl )  E(wl 1 ))
各温度間で、交換が行われた頻度
1
2
3
4

2
3
1
2
2
3
・もし、各温度で定常分布に収束していると、
平均交換率は、２つの温度パラメータによって定まる関数。
平均交換率の理論解析
 がある程度大きい状況
＜低温同士での平均交換率＞


w R 


1  2 
1



p
w

exp  1Eˆ (w)  (w)
 1
Z (1 )


 p2 w  1 exp  2 Eˆ (w)  (w)

Z ( 2 )

Eˆ (w) ：
(w) ：
d
ˆ (w0 ) をもつ関数
ある点 w0 において、最小値 E
任意の確率分布
平均交換率の理論解析
＜交換に関する採択確率＞
u1  min(1, r) for Metropolis type
r
u2 
for heat bath type
1 r
r  exp  (2  1)(Eˆ (w1)  Eˆ (w2 ))

＜平均交換率＞
Metropolis type: J1   u1 p1 (w1 ) p2 (w2 )dw1dw2
heat bath type: J 2   u2 p1 (w1 ) p2 (w2 )dw1dw2

メトロポリス型における平均交換率
＜定理１＞[Nagata, 2008]
平均交換率
J1 は、 1, 2   において以下の式に収束する。
 2  1  22   2  1 


J1  1 
A ,
2
1    
1 

1
s  1
A, a  
ds
0 1  a  s 2

熱浴型における平均交換率
＜定理２＞[Nagata, 2008]
平均交換率 J 2 は、 1, 2
  において以下の式に収束する。


1   2  1 
1

J 2  1  1 
2
2 

1
  

 2  1 

B ,

1 




1

B, a   ds1  ds2 tanh as2  s1  exp s1 exp 1  as2 s1 1s2 1
0
0
2

平均交換率の挙動
 2  1  22   2  1 


Metropolis type: J1  1 
A ,
2
1    
1 





 2  1 
1
2  1
1




heat bath type: J 2  1  1 
B

,

2


2 


1
1
   




平均交換率は、
平
均
交
換
率
2  1 / 1 と 
の関数

  5.0
2  1 / 1

 2 1
1

 2.0
平均交換率と温度パラメータ
平均交換率は、
各温度で
2  1
1
1 2
l 1  l
l
3   2
2
2  1 / 1 の関数
が一定ならば、平均交換率は同じ値になる。
 4  3
3
3
5   4
4
5
4
このとき、温度パラメータの値は、 l
 1  1  
指数的に区切れば、平均交換率が一定の値になる。
l 1

 ってなに？
「自由エネルギー」・「周辺対数尤度」・「確率的複雑さ」


F ( )   log Z ( )   log  exp  Eˆ (w)  (w)dw
  log   o(log  ).
・
  
 の値は、主に Eˆ (w)の性質によって定まる。
＜正則なケース＞
ˆ (w)のヘッセ行列が正定値
・任意の wにおいて E
  2 Eˆ (w) 

ヘッセ行列： 
 wi w j 



d
2
 ってなに？
＜特異なケース＞
・ヘッセ行列が縮退する w が存在する場合
d

2
右の例では、

1
(d  2)
2
・特異なケースでの  の解析法
Im z
 (z)の極を調べる。


z
ˆ
 ( z)   E(w) (w)dw
代数幾何学の手法である
特異点解消を行えば求められる。
 0
Re z
ベイズ学習での
 の性質
＜汎化誤差[Watanabe, 2001]＞
 

G X n  KL q( x) || p( x | X n )

：真の構造と予測分布のカルバック距離
予測の精度を示す尺度の一つ
漸近形： G (n)  E
 1
GX   n  o n  n  
Xn
 
様々な学習モデルにおいて、
n
 の値を求める研究がなされている。
＜厳密解＞
・ニューラルネットワーク
・隠れマルコフモデル
・縮小ランク回帰モデル
・混合二項分布モデル
＜上限値＞
・一般混合分布モデル
・確率文脈自由文法
・ベイジアンネットワーク
平均交換率の理論値の検証
＜縮小ランク回帰モデル（線形ニューラルネットワーク）のベイズ学習＞
y(1) y(2)

パラメータ
w  {ai,h},{bh, j }
パラメータの次元： d

N
H
y(N )
A
 (M  N ) H
 N  20, H  12 のとき
（真の構造は H  4 ）
0
学習モデルが M
B
M

x(1) x( 2)
x(M )
y  ABx  (noise)

温度パラメータ設定に関する研究

平均交換率を均等にするよう設定



関数 E(w)などの挙動をもとにした繰り返しアルゴリズム



K. Y.Sanbonmatsu, A. E. Garcia, Proteins,46 (2002).
A. Schug et al, Proteins,57 (2004).
最適な平均交換率の値




Y.Sugita, Y. Okamoto,Chem. Phys. Lett.,314 (1999).
D. A. Kofke, J. Chem. Phys.,117 (2002).
20%～25%くらいが最適らしい。
N. Rathore et al, J.Chem.Phys.,122 (2005).
A. Kone,D. A. Kofke, J . Chem.Phys.,122 (2005).
温度の端から端まで動くための時間の最小化

H. G. Katzgraberet al, J . Stat.Mech.,(2006).
全体のまとめ

マルコフ連鎖モンテカルロ法



交換モンテカルロ法



メトロポリス法
マルコフ連鎖の原理
遅い緩和の問題
交換モンテカルロ法の原理
交換モンテカルロ法の設計に関する理論


温度パラメータの設定
平均交換率の理論解析