確率モデルによる 画像処理技術入門

物理フラクチュオマティクス論
Physical Fluctuomatics
応用確率過程論
Applied Stochastic Process
第8回 マルコフ連鎖モンテカルロ法
8th Markov chain Monte Carlo methods
東北大学 大学院情報科学研究科 応用情報科学専攻
田中 和之(Kazuyuki Tanaka)
[email protected]
http://www.smapip.is.tohoku.ac.jp/~kazu/
物理フラクチュオマティクス論(東北大)
1
確率的情報処理における計算の壁
不確実性を確率・統計を用いて表現することの代償
起こりやすいことも起こりづらいこともまじめ
に考慮して計算
計算量的困難
物理フラクチュオマティクス論(東北大)
2
計算困難のポイントは何か
a  0;
for(x1  0, 1){
2L 通りの和が計算できるか?
    f x , x ,, x 
x1 0,1 x2 0,1
xL 0,1
1
2
for(x2  0,1){

for(xL  0,1){
L
a  a  f  x1 , x2 , , xL ;
このプログラムでは
L=10個のノードで1秒かかるとしたら
L=20個で約17分,
L=30個で約12日,
L=40個で約34年かかる.
厳密に計算するのは一部の特殊な例を
除いて難しい.
マルコフ連鎖モンテカルロ法
確率伝搬法
}

}
}
L 重ループ
今回
次回
物理フラクチュオマティクス論(東北大)
3
乱数による円周率の計算(モンテカルロ法)
1
n0 m0
nn+1
-1
区間[-1,1]において2個の
一様乱数a,bを発生
a2+b2≤1
No
1
-1
緑の四角の枠の中に
ランダムに点をプロットし,
赤い円の内部の点の個数
をカウントする.
Yes
mm+1
0
円周率
S  R 2
R 1
物理フラクチュオマティクス論(東北大)
 S
4m
n
試行回数 n が大きいほど
精度が上がる
4
大数の法則
X1,X2,...,Xn は平均 , 分散 s2 の互いに独立な同一の確
率変数であるとき
Yn 
1
( X 1  X 2    X n )   (n   )
n
中心極限定理
X1,X2,...,Xn は平均 , 分散 s2 の互いに独立な同一の確率
変数であるとき
Yn 
1
( X1  X 2    X n )
n
は n が大きいとき平均 , 分散 s2/n の正規分布に従う.
物理フラクチュオマティクス論(東北大)
5
乱数による円周率の計算(モンテカルロ法)
1
1個1個の点をランダムにプロットすると
いう試行を行った時の x 座標,y 座標
の平均 0, 分散は 1/2.
ランダムに n 個プロットした標本
平均は平均が 0, 分散が 1/2n
(中心極限定理)なので,その確
率密度関数の幅は 1/n1/2 で減少
する.
円周率の推定誤差は O(1/n1/2)
-1
0
1
-1
緑の四角の枠の中に
ランダムに点を n 個プロットし,
赤い円の内部の点の個数 m を
カウントする.
円周率
S  R 2
R 1
物理フラクチュオマティクス論(東北大)
 S
4m
n
試行回数 n が大きいほど
精度が上がる
6
積分の計算(モンテカルロ法)
1
n0 m0
I 
1 1

1 1
f ( x, y)dxdy
nn+1
-1
区間[-1,1]において2個の
一様乱数a,bを発生
4m
n
1
-1
緑の四角の枠の中に
ランダムに点をプロットし,
その点の被積分関数の
値を求める操作を繰り返す.
mm + f(a,b)
I
0
試行回数 n が大きいほど
精度が上がる
物理フラクチュオマティクス論(東北大)
推定誤差は O(1/n1/2)
積分が高次元でも同様
7
確率推論でよく必要となる計算
周辺確率
P1 x1  
    Px , x ,, x 
x2 0,1 x3 0,1
P12 x1 , x2  
P13 x1 , x3  
xL 0,1
1
2
L
    Px , x , x ,, x 
x3 0,1 x4 0,1
xL 0,1
1
2
3
L
    Px , x , x ,, x 
x2 0,1 x4 0,1
xL 0,1
1
物理フラクチュオマティクス論(東北大)
2
3
L
8
計算のポイントは
与えられた確率分布
P x   Px1 , x2 ,, xL 
に従うN個の独立なサンプル x を生成するアル
ゴリズムが作れるか?
但し,1個のサンプルを生成するための計算量は
L の多項式オーダーでなければならない.
物理フラクチュオマティクス論(東北大)
9
簡単な確率過程:マルコフ連鎖
推移確率 w(x|y)≥0 (x,y=0,1) が与えられたとき
 wz x   1
任意の P0(x) を初期値として
Pt ( x) 
1
 w( x | z ) Pt 1 ( z )
z 0
( x  0,1; t  1,2, )
z 0
 Pt (0)   w(0 | 0) w(0 | 1)  Pt 1 (0) 

  


 Pt (1)   w(1 | 0) w(1 | 1)  Pt 1 (1) 
推移確率行列
物理フラクチュオマティクス論(東北大)
10
簡単な確率過程:マルコフ連鎖
Pt ( x) 

1
 w( x | z[t  1])Pt 1 ( z[t  1])
z[t 1]  0
1

1

z[t 1]  0 z[t  2]  0
1

 w( x | z[t  1])w( z[t  1] | z[t  2]) w( z[1] | z[0])P0 ( z[0])
z[t 1]  0
 Pt (0) 
 P0 (0) 
t

  W 

 Pt (1) 
 P0 (1) 
 Pt (0) 
 1 0  1  P0 (0) 

  U 

t U 
0  
 Pt (1) 
 P0 (1) 
遷移確率行列は
 1 0  1
U
W  U 
(   1)
0  
と対角化可能
 Pt (0) 
 P(0) 
 1 0  1  P0 (0) 
  U 


  lim 
U 
 P(1)  t   Pt (1) 
 0 0
 P0 (1) 
極限分布
物理フラクチュオマティクス論(東北大)
11
簡単な確率過程:マルコフ連鎖
Pt ( x) 
P( x) 
1
 w( x | z[t  1])Pt 1 ( z[t  1])
z[t 1]  0
1
 w( x | z[t  1])P( z[t  1])
z[t 1]  0
 P(0) 
 P(0) 

  W 

 P(1) 
 P(1) 
定常分布または平衡分布
 Pt (0) 
 P(0) 
 が存在すれば定常分布は存在する


 lim 
極限分布 

 P(1)  t   Pt (1) 
定常分布が存在しても極限分布が存在するとは限らない.
0 1
 P(0)  1 / 2 
Example

W  

  
 は定常分布
1
0


 P(1)  1 / 2 
物理フラクチュオマティクス論(東北大)
12
マルコフ連鎖の定常分布と詳細釣り合い
Pt x  
1
 wx y Pt 1 ( y)
ただし
y 0
1
 wy x   1
y 0
P1(x), P2(x), P3(x),…: マルコフ連鎖
wy xP( x)  wx y P( y)
詳細釣り合い
を満たすように推移確率 w(x|y) を選ぶと
その定常分布は P(x) になる.
P x  
1
 wx y P y 
マルコフ連鎖の定常分布
y 0
物理フラクチュオマティクス論(東北大)
13
マルコフ連鎖モンテカルロ法

確率分布P( x )  P( x1 , x2 ,, x L ) が与えられたとき

任意の P0 ( x ) を初期値として



Pt ( x )   wx y Pt 1 ( y ) (t  1,2,3,)

y


を繰り返し計算して, lim Pt ( x )  P( x ) となるように

wx y 
t  
を選ぶにはどうしたらよいか?
物理フラクチュオマティクス論(東北大)
14
マルコフ連鎖モンテカルロ法



Pt x    wx y Pt 1 ( y )

y

ただし  w y x   1

y
P1(x), P2(x), P3(x),…: マルコフ連鎖
 
 
詳細釣り合い wy x Px   wx y P y 

wx y  を選ぶとその推移確率
を満たすように推移確率

の定常分布は P (x ) になる.


 

w
x
y の
推移確率
P x    wx y P y 

マルコフ連鎖の定常分布
y
極限分布がただ一つ存在するとしたら定常分布のひとつが極限分布である.
物理フラクチュオマティクス論(東北大)
15
マルコフ連鎖モンテカルロ法

xt  1
 
wxt  xt  1

x t 
無駄弾

x[1]

x[t  1]


x[2]

x[t  2]





x[t ]

x[2t ]




x[(N  1)t  1] x[(N  2)t  1]  x[ Nt  1]
精度は
1/2
O(1/t )

P (x ) からのサンプリング
ランダムに生成
τ が十分大きい
ときサンプルx[t],
x[2t], x[3t], …,
x[Nt] は独立
独立とみ
なせるか
t をどれだけ大き
くとるべきか→緩
和時間
とみなすことができる
物理フラクチュオマティクス論(東北大)
16
マルコフ連鎖モンテカルロ法

x[1]

x[2]

x[t  1]


x[t  2]


x[t ]

x[2t ]


x[ Nt ]





x[(N  1)t  1] x[(N  1)t  2] 
P1  X 1  
    P X 1 , X 2 , X 3 ,, X L 
X2 X3
XL
頻度
N
1

 x1 nt , X 1

N n 1
Xi
ヒストグラムから周辺確率も計算できる
物理フラクチュオマティクス論(東北大)
17
マルコフ連鎖モンテカルロ法

P( x)  P( x1 , x2 ,, xL ) 

ij
( xi , x j )
{i , j}E
Px1 , x2 ,, xi 1 , xi , xi 1 ,, xL 
Pxi x1 , x2 ,, xi 1 , xi 1 ,, xL  
Px1 , x2 ,, xi 1 , xi 1 ,, xL 
Px1 , x2 ,, xi 1 , xi 1 ,, xL    Px1 , x2 ,, xi 1 , zi , xi 1 ,, xL 
zi
V:すべてのノードの集合
E:すべてのリンクの集合
物理フラクチュオマティクス論(東北大)
18
マルコフ連鎖モンテカルロ法

P( x)  P( x1 , x2 ,, xL ) 

{i , j}
( xi , x j )
{i , j}E
Pxi x1 , x2 , , xi 1 , xi 1 , , xL
  x , x 

 P x x
  z , x 
{i , j }
i
ji
i
{i , j }
zi
(V , E )
リンクの集合
i
j

j  i
j
ji
マルコフ確率場
E:すべての
j
(V , E )
i :ノード i とリンクで結ば
れたすべてノードの集合
物理フラクチュオマティクス論(東北大)
19
マルコフ連鎖モンテカルロ法

P( x)  P( x1 , x2 ,, xL ) 

{i , j}
( xi , x j )
{i , j}E
wx1, x2 , , xL x1 , x2 , , xL
xk  xk k  i 
(V , E )
 x, x 




      ( x , x ) 

   z , x 
{i , j }
i
j
ji
iV
k
k
kV / i
{i , j }
zi
i
j
ji
  


wx' x Px   wx x' Px' 

xt  1
 
wxt  xt  1
物理フラクチュオマティクス論(東北大)

x t 
20
マルコフ連鎖モンテカルロ法
 
wxt  xt  1


x t 
x[t  1]
True
xi = ○ or
False
●

x

1回の更新にかか
る計算量は変数の
次元 L に対して
O(1)
x’
物理フラクチュオマティクス論(東北大)
21
マルコフ連鎖モンテカルロ法
によるサンプリング
p
p



p が小さい
p が大きい
無秩序状態
秩序状態
マルコフ連鎖モン
テカルロ法による
サンプリング
More is different.
ある p の値の付近で
ゆらぎが大きくなる.
量が増えれば質が変わる
物理フラクチュオマティクス論(東北大)
22
今回の参考文献
マルコフ連鎖モンテカルロ法の詳細
伊庭幸人, 種村正美: 統計科学のフロンティア/計算
統計 II ---マルコフ連鎖モンテカルロ法とその周辺,
岩波書店, 2005.
マルコフ過程
伏見正則, 確率と確率過程,第7章, 講談社, 1987.
物理フラクチュオマティクス論(東北大)
23
本日のまとめ
モンテカルロ法による円周率の計算
大数の法則と中心極限定理
マルコフ連鎖モンテカルロ法
次回
第9回 確率伝搬法
第10回 確率的画像処理と確率伝搬法
第11回 確率推論におけるベイジアンネットと確率伝搬法
物理フラクチュオマティクス論(東北大)
24
演習問題8-1
詳細釣り合い
および
wx' x P x   wx x'P x'
 wx x   1 を満たすとき
P x    w x x' P x' 
x
x'
が成り立つことを示せ.
物理フラクチュオマティクス論(東北大)
25
演習問題8-2
1 2 1

 により与えられるとき,
W

推移確率行列が
3  1 2
 Pt (0) 
 P(0) 
 を求めよ.
  lim 
その極限分布 
P
(
1
)
P
(
1
)

 t   t 
 Pt (0) 
 Pt 1 (0) 



W
但し, P (1)   P (1) 
 t 
 t 1 
である.
物理フラクチュオマティクス論(東北大)
26
演習問題8-3
以下の図で与えられる正方格子 L=Lx×Ly によるグラフのすべてのリンクからな
る集合を B として各ノードの確率変数 Fi が 0,1,2,…,Q-1 の整数値をとり,結合
確率が
PrF1 , F2 ,, FL  
1
Z Prior
 1
2
exp    Fi  Fj  
 2 {i , j}E

により与えられる無向グラフ上の確率モデルから互いに独立な状態 (f1,f2,…,fL)
をランダムに N 個生成するプログラムをマルコフ連鎖モンテカルロ法により作成
し,様々のα の値に対して数値実験を実行せよ.
Q=256, α=0.0005 に
対して生成されたサン
プルの例
Ly
L=Lx×Ly
Lx
Q=2, α=2 に対して生
成されたサンプルの
例
物理フラクチュオマティクス論(東北大)
27