確率論入門 第15章 §15.5 ★★ 詳細つり合いを満たすように遷移を設計すると、目的の分布を定常分布にもつ連鎖が作れる(MCMC)。対称性を使えば大きな状態空間の連鎖を1次元の連鎖へ還元できる
詳細つり合いと可逆性
定常分布を求めるには πP=π を解けばよいが、状態が多いと連立方程式は大きくなる。次の条件は、より強いが確かめやすい十分条件である。
命題 15.26(詳細つり合い ★★★)
確率分布 π が、すべての i,j について
πipij=πjpji(15.32)
を満たすなら、π は定常分布である。
証明
(15.32) を i について足すと
i∑πipij=πji∑pji=πj
である。
(15.32) は「定常状態で、2つの状態のあいだを行く確率の流れと帰る流れが釣り合う」ことを表す。このとき X0∼π の連鎖は、時間を逆向きにたどっても同じ法則に従う。実際、経路の確率に (15.32) を繰り返し使うと
P(X0=i0,…,Xn=in)=πi0pi0i1⋯pin−1in=πinpinin−1⋯pi1i0
となる。このような連鎖を可逆であるという。
定常分布があっても、詳細つり合いが成り立つとは限らない(確認テスト 15.5-1)。
例
出生死滅連鎖。 {0,1,2,…} 上で隣にしか移らない連鎖
pi,i+1=bi,pi,i−1=di,pii=1−bi−di
を考える。(15.32) は隣どうしの式 πibi=πi+1di+1 だけになり、
πk=π0ℓ=0∏k−1dℓ+1bℓ
と順に決まる。右辺の総和が有限なら、規格化して定常分布が得られる。例題 15.5 のグラフ上のランダムウォークも可逆である。隣どうしなら
πipij=2Edegi⋅degi1=2E1=πjpji
となるからである。
MCMC と Metropolis 法
有限集合 S 上の分布 π から標本を作りたいとする。ただし π は定数倍を除いてしか分からない。
πi=Zwi,Z=i∑wi
重み wi>0 は計算できるが、規格化定数 Z は状態が多すぎて計算できない、という状況である。統計力学の Boltzmann 分布 wi=e−βEi が典型例である。
マルコフ連鎖モンテカルロ法(MCMC)は、π を定常分布にもつ既約な連鎖を作って長く走らせる。そしてエルゴード定理 15.25 により、時間平均から ∑jπjf(j) を推定する。
定理 15.27(Metropolis 法 ★★★)
q=(qij) を既約な遷移確率行列とし、対称 qij=qji であるとする。これを提案の行列と呼ぶ。i=j について
pij=qijmin(1,wiwj),pii=1−j=i∑pij(15.33)
と定める。この P は、πi=wi/Z について詳細つり合いを満たし、既約である。したがって π は P のただ1つの定常分布である。
証明
i=j とし、wj≤wi と仮定しても一般性を失わない。このとき pij=qijwj/wi、pji=qji⋅1 なので
πipij=Zwiqijwiwj=Zwjqij,πjpji=Zwjqji=Zwjqij
で、両者は等しい。
次に、P が確率行列であることを見る。∑j=ipij≤∑j=iqij≤1 なので、pii≥0 である。
重みは正なので、qij>0 なら pij>0 である。よって提案の行列が既約なら P も既約である。命題 15.26 と系 15.20 から結論を得る。
アルゴリズムとしては、次の手順をくり返す。
- 現在の状態 i から、候補 j を提案の行列に従って選ぶ。
- 確率 min(1,wj/wi) で受け入れて、候補に移る。
- 受け入れなければ、今の状態に留まる。
重みの比 wj/wi しか使わないので、Z は要らない。重みが増える提案は必ず受け入れ、減る提案は比の確率でだけ受け入れる。提案が対称でないときは、受け入れ確率を
min(1,wiqijwjqji)
に変えればよい(Metropolis–Hastings 法、確認テスト 15.5-5)。
落とし穴
- 連鎖の標本は独立でない。 時間平均は (15.31) により収束するが、隣り合う標本は相関しているので、同じ個数の独立な標本より推定の誤差は大きい。
- 初期値の影響が残る。 定理 15.21 の収束は漸近的なので、最初のしばらくの値を捨ててから平均をとることが多い。
- 留まる確率 pii>0 があると非周期的になる。 受け入れを拒む確率が正の状態が1つでもあれば、既約な連鎖は非周期的で、各時刻の分布も π に収束する。
例題 15.7 3状態の Metropolis 法 ★★
オリジナルS={1,2,3} 上の分布 πi∝i を目標とし、提案は「今の状態以外の2つから等確率で1つ選ぶ」とする。
- Metropolis 法の遷移確率行列 P を書き、πP=π を確かめよ。
- 連鎖を長く走らせたときの Xk の時間平均の極限を求めよ。
解答
(1) 目標は π=(61,62,63) である。提案は、異なる2状態のあいだで qij=21 である。(15.33) から、遷移確率は次のとおりである。
- 重みが増える向きは pij=21 である。
- 重みが減る向きは pij=21⋅ij である。
P=04161214131212121
対角成分は、各行の和が 1 になるように決めた。
p22=1−41−21=41,p33=1−61−31=21
πP の各成分は
(πP)1(πP)2(πP)3=61⋅0+62⋅41+63⋅61=121+121=61,=61⋅21+62⋅41+63⋅31=121+121+61=31,=61⋅21+62⋅21+63⋅21=21
で、πP=π である。
(2) 連鎖は既約なので、エルゴード定理 15.25 を f(i)=i に使える。時間平均は
i∑πii=61+4+9=37
に収束する。
多次元の状態空間を1次元へ還元する
状態空間が大きくても、知りたい量が状態の一部の情報だけで決まり、その情報がそれ自体でマルコフ連鎖になっていれば、小さな連鎖で計算できる。
定理 15.28(集約できる条件 ★★★)
状態空間の分割 S=⨆aBa をとる。各状態にそのブロックの名前を返す関数を ϕ と書く。つまり i∈Ba なら ϕ(i)=a である。すべての a,b について
j∈Bb∑pij が i∈Ba によらない(15.34)
とき、その値を p^ab と書く。このとき、もとの連鎖がどの初期分布から出発しても、Yn=ϕ(Xn) は遷移確率 p^ab のマルコフ連鎖である。
証明
事象 {Y0=a0,…,Yn=an} は、事象 {X0∈Ba0,…,Xn∈Ban} と同じである。これを、時刻 n の値 Xn=i で排反に分ける。マルコフ性から、各部分の上で
P(Xn+1∈Bb∣X0∈Ba0,…,Xn−1∈Ban−1,Xn=i)=j∈Bb∑pij=p^anb
である。右辺は i によらないので、部分ごとの確率を足し合わせても
P(Yn+1=b∣Y0=a0,…,Yn=an)=p^anb
となる。
(15.34) は、同じブロックに属する状態はどれも、次にどのブロックに移るかについて区別がつかないという条件である。状態空間の対称性からこの条件が出ることが多い。
例題 15.8 Ehrenfest の壺 ★★
オリジナル2つの壺 A、B に、区別できる N 個の玉が入っている。各時刻に N 個の玉から1個を等確率で選び、もう一方の壺へ移す。
- 全体の状態を、各玉がどちらの壺にあるかを表す x∈{0,1}N とする。この 2N 状態の連鎖が、壺 A の玉の数 K=∑ℓxℓ の連鎖に集約できることを示せ。また K の遷移確率を求めよ。
- K の定常分布を詳細つり合いから求めよ。
- N=4 とする。壺 A が空の状態から出発して、再び空になるまでの平均時間 m0 を求めよ。また m2 を求めよ。
解答
(1) 全体の連鎖は、x の成分を1つ選んで反転させる連鎖である。A に k 個の玉がある状態から1個を選ぶ。
- A にある k 個のどれかを選べば、A の玉は1個減る。
- B にある N−k 個のどれかを選べば、A の玉は1個増える。
よって
y:K(y)=k−1∑pxy=Nk,y:K(y)=k+1∑pxy=NN−k
で、どちらも K(x)=k にしかよらない。条件 (15.34) が成り立つので、定理 15.28 から K は次の遷移確率の連鎖である。
p^k,k−1=Nk,p^k,k+1=NN−k
(2) 出生死滅連鎖なので、詳細つり合い πkNN−k=πk+1Nk+1 を使う。これから
πkπk+1=k+1N−k=(kN)(k+1N)
となり、πk=(kN)2−N である。これは二項分布 Bin(N,21) で、各玉が独立に確率 21 でどちらかの壺にある状態に対応する。
(3) K の連鎖は既約で有限なので、定理 15.19 から
mk=πk1=(kN)2N
である。N=4 では
m0=16,m2=616=38
である。玉を増やすと m0=2N は急速に増える。100 個なら 2100≈1.3×1030 である。片方の壺が空の状態は再帰的なので必ず戻ってくるが、玉の数が多いと実際には観測できないほど先になる。
吸収を問う問題と長時間を問う問題の見分け
- 「いつか〜に達する確率」「〜より先に」「達するまでの平均時間」を問うとき。 吸収状態か到達先の集合 A を決め、一次ステップ解析(定理 15.10・15.11)を使う。有限なら基本行列 N=(I−Q)−1 でまとめて計算できる。
- 「長い時間のうち〜にいる割合」「十分時間がたったときの確率」「平均再帰時間」を問うとき。 既約で正再帰的かを確かめ、πP=π(または詳細つり合い)を解く。各時刻の分布の極限には非周期性を確かめる。
- 閉じたクラスが複数あり、過渡的な状態もあるとき。 2つを組み合わせる。閉じたクラス C が非周期的・正再帰的で、定常分布 πC をもつとする。i から出発してそのクラスの状態 j にいる確率の極限は
n→∞limpij(n)=Pi(τC<∞)πjC
である。第1因子は一次ステップ解析で求め、第2因子は定常分布で求める。クラスに入ったあとは、強マルコフ性でクラスの中の連鎖になるからである。
確認テスト
オリジナル
- 15.5-1 S={1,2,3} で p12=p23=p31=1 の連鎖を考える。この連鎖は一様分布を定常分布にもつが、詳細つり合いを満たさないことを示せ。
- 15.5-2 S={0,1,2,3} 上の分布 πk∝2−k を目標とする。提案は左右の隣を確率 21 ずつ選び、範囲の外を選んだら必ず拒んでその場に留まるとする。Metropolis 法の遷移確率をすべて求めて、詳細つり合いを確かめよ。また Xk の時間平均の極限を求めよ。
- 15.5-3 例題 15.8 で N=6 とする。π3 と m3 を求めよ。また m0 を求めよ。
- 15.5-4 円周 Z5={0,1,2,3,4} の上を2つの粒子が独立に動く。各時刻にそれぞれが、確率 21 ずつで左右に1動く(法5)。差 Dn=Xn−Yn(mod5) がマルコフ連鎖であることを示し、遷移確率を求めよ。また D0=1 から、2つの粒子が初めて同じ点に来るまでの時間の期待値を求めよ。
- 15.5-5 提案 q が対称とは限らない場合を考える。i=j かつ qij>0 のとき
pij=qijmin(1,wiqijwjqji)
と定めると、詳細つり合いが成り立つことを示せ。
AI(Claude)が書いた教材を、木下が勉強しながら公開しています。式と数値は検算していますが、誤りが残っているかもしれません。誤りのご指摘・削除のご依頼をいただければ直します。 連絡先:support@tecktom.com