原論文和譯及び解說 Geman-Gemanの定理 (シュミレーテッドアニーリングの収束定理)

原論文:http://image.diku.dk/imagecanon/material/GemanPAMI84.pdf

まずは論文概要や導入の概略を見てみましょう。原論文そのものは、劣化画像の復元のための手段としてシュミレーテッドアニーリングを捉えていることを念頭に置いておくと理解しやすいでしょう。原論文に数式などの不備や書き間違いが多少あるので、そこは修正しています。また、シュミレーテッドアニーリングそのものの理解にあまり必要ないと思われる章は飛ばしています。補足は、私が書いたものです。

 

1)直近の「近傍」のピクセルと境界要素の現在の値に基づいて、画像に局所的な変更が加えられます。この変更はランダムであり、ローカル条件付き確率分布からサンプリングすることにより生成されます。

2)局所条件付き分布は、「温度」と呼ばれるグローバル制御パラメータ\(T\)に依存します。低温では、局所条件付き分布は目的関数を増加させる状態に集中しますが、高温では分布は本質的に均一である。極限を考えると、\(T=0\)及び\(T=\infty\)は、それぞれ貪欲アルゴリズム(勾配上昇など)及び無向(つまり「純粋にランダム」)な変更に対応します。(高温は、隣接するピクセル間の疎結合と画像の混沌とした外観を引き起こします。低温では、結合はより緊密になり、画像はより規則的に見えます。)


3)我々の画像復元法は、確率的変化の多くが起こる高温で開始することにより、局所的な最大値を回避し、実際に目的関数を減少させます。緩和が進むと温度は徐々に低下し、過程は繰り返し改善されていくように振る舞います。(「焼き鈍し」をシュミレートしたこの温度の段階的低下は、特定の化学システムを低エネルギーの非常に規則的な状態に遷移させる過程です。)


我々の「アニーリング定理」は、事後分布のグローバルな最大値への収束を保証する温度を下げるためのスケジュールを規定しています。
このスケジュールは実際に適用するには遅すぎ、温度と時間の依存関係の関数形式を選択する際のガイドとしてのみ使用します。


また、私自身は引用論文は読んでいませんが、シュミレーテッドアニーリングの考え方の源流がどこにあるかも見て取れます。


温度と焼き鈍しのシュミレートを導入する発想は、Cerny[8]及び Kirkpatrick et al.[40]に依る。彼らはどちらも最適化問題に用いていて、Kirkpatrickはコンピュータ設計にも適用している。

[8] V. Cerny, “A thermodynamical approach to the traveling salesman problem: an efficient simulation algorithm,” preprint, Inst.
Phys. & Biophys., Comenius Univ., Bratislava, 1982.

[40] S. Kirkpatrick, C. D. Gellatt, Jr., and M. P. Vecchi, “Optimization by simulated annealing,” IBM Thomas J. Watson Research
Center, Yorktown Heights, NY, 1982.


いよいよ、定理の主張を見ていきます。

\(T\)が減少すると、事後分布からのサンプルは最小エネルギー構成に向かって強制されていきます。これらは、分布のモードに対応します。定理Bでこれを厳密にし、私たちの知る限り、この性質の最初の理論的結果です。大まかに言えば、\(k\)番目のサイト置換の実行に使用される温度\(T(k)\)(つまり、反復スキームの\(k\)番目のイメージ)が次の条件を満たす場合、

\[T(k)\geq \frac{c}{\log(1+k)}\]

すべての\(k\)について、\(c\)は\(k\)に依存しない定数であり、確率が1に収束する(\(k\rightarrow \infty\)とした)場合、アルゴリズムによって生成される構成は最小エネルギーの構成になります。
別の言い方をすれば、アルゴリズムは、最小エネルギー構成で分布が均一な尺度に収束するマルコフ連鎖を生成します。(確率1の収束は、一般に不可能であることを強調する必要があります。)

また、実際の収束率に関する定数\(c\)の性質についても議論します。

基本的に我々は、対数レートが最適であると考えます。


では定理の証明に関する話題に入っていきましょう。まずは数学的準備です。

3.グラフと近傍

\(S=\{s_1,s_2,…,s_N\}\)をサイトの集合とし、\(G=\{G_s,s \in S\}\)を\(S\)の近傍、つまり\(s \notin G_s,かつs\in G_r \Leftrightarrow r \in G_s\)を満たすような\(S\)の部分集合のそれぞれの集まりと定義する。

明らかに、\(G_s\)は\(s\)の近傍であり、\(\{S,G\}\)は普通の意味でのグラフを成す。

部分集合\(c  \subseteq S\)は、\(c\)のそれぞれのサイトの全てのペアが隣接している場合、クリークであり、\(C\)をクリーク\(c\)全ての集合とする。

 

4.マルコフ確率場とギブス分布

 

3章と同じように、\(\{S,G\}\)は任意のグラフを表す。\(X=\{X_s,s \in S\}\)を\(S\)を添え字集合としたランダムな値の族とする。

簡単のために、全ての\(s\)に対して\(X_s\in A\)なる\(A\equiv \{0,1,2,…,L-1\}\)のような共通の状態空間を想定する。

(補足:普通のシュミレーテッドアニーリングにおける、各サイトの取りうる値\(\{0,1\}\)に対応します。(quboの場合。量子ビットでは\(\{-1,1\}\)))

 

\(\Omega\)を全ての可能な構成の集合とする。

\[\Omega =\{ \omega = (x_{s_1},..,x_{s_N}):x_{s_i} \in A ,1 \leq i \leq N\}\]

(補足:普通のシュミレーテッドアニーリングでは、サイト数\(N\)のとき\(\Omega=\{(0,0,..,0),…(1,1,..,1)\}で、|\Omega|=2^N\)になりますね)

 

\(X\)は以下を満たすとき、\(G\)に関するMRF(マルコフ確率場)である。

全ての\(\omega \in \Omega\)に関して、\(P(X=\omega)>0\)

全ての\(S^{‘},S^{”} \subset S(r \neq s)\)と\((x_{s_1},..,x_{s_N})\in \Omega\)に関して、グラフ\(\{S,G\}\)上で\(B\)が\(S^{‘},S^{”}\)を分断している、つまり、\(S^{‘}\)から\(^{”}\)に行く全ての経路が\(B\)を通るとすると、

\(X(B)\)が与えられた時、\(X(S^{‘})\)と\(X(S^{”})\)は条件付き独立になる。

(補足:大域マルコフ性。上のは書き直してありますが、論文を見るときは論文の説明がまずいので、ウィキペディアの方を参考にして下さい。)

 

\(\{S,G\}\)に関するギブス分布は、以下の表式を持つ\(\Omega\)上の確率測度である。

\[\pi(\omega)=\frac{1}{Z}e^{-U(\omega)/T}\]

\(T\)は定数であり、\(U\)はエネルギー関数である。\(Z\)は規格化定数であり、

\[Z=\sum_{\omega}e^{-U(\omega)/T}\]

 分配関数と呼ばれる。

 

11.ギブスサンプラー:概要

 

一般的な計算上の問題は、

A)分布\(\pi\)からサンプルする

B)\(\Omega\)上で\(U\)を最小にする

C)期待される値を計算する

である。

 

時刻\(t\)で、\(X_{s_t}\)は\(A\)のいずれかの値である。各更新では、一つのサイトのみが値を変えることができる。つまり、\(X(t-1)\)と\(X(t)\)では一か所のサイトしか変化していない。\(n_1,n_2,..\)を、どのサイトが変化するために「訪問」されたかを示す数列とする。つまり、\(n_t \in S\)であり、\(i \neq n_t\)の時\(X_{s_i}(t)=X_{s_i}(t-1)\)である。

 

12.ギブスサンプラー:数学的基礎

 

定理A(省略)

一様マルコフ連鎖の収束性について述べています。ペロン・フロベニウスの定理の証明をもう一回やっているだけです。別の文献がいくらでもあります。

 

定理Bは、「アニーリングスケジュール」または温度低下率であり、システムを強制的に最低エネルギー状態にする。

\(\pi\)の\(T\)依存性を\(\pi_T\)と書き、\(T(t)\)をステージ\(t\)での温度とする。

最小エネルギー状態\(\Omega_0=\{\omega \in \Omega:U(\omega)=\min _{\eta}U(\eta)\}\)とし、

\(\pi_0\)を\(\Omega_0\)上の一様分布とする。

(補足:最小エネルギー状態が\(n\)個あるとき、最小エネルギー状態の取りうる確率が\(1/n\),それ以外が0になるような確率測度)

最後に、

\(U^*=\max_{\omega}U(\omega)\)

\(U_*=\min_{\omega}U(\omega)\)

\(\Delta=U^*-U_*\)

と定義する。

 

定理B(アニーリング)

全ての\(t=0,1,2…\)で\(S \subset \{n_{t+1},n_{t+2},…,n_{t+\tau}\}\)となる\(\tau>N\)の存在を仮定する。

(補足:そういう風にマルコフ連鎖を構成すればいいだけなので、容易に前提条件は満たせます。例えばそれぞれのサイトを順番に回るなど。)

 

\(T(t)\)を、以下を満たすような温度の減少列とする。

a)\(t \rightarrow \infty\)で\(T(t)\rightarrow 0\)

b)\(t_0\geq 2\)なる適当な整数に関して、全ての\(t \geq t_0\)で\(T(t)\geq N\Delta/\log t\)

を満たすならば、

任意の初期分布\(\eta\in \Omega\)と全ての\(\omega \in \Omega\)に関して、

\(\lim_{t\rightarrow \infty }P(X(t)=\omega|X(0)=\eta)=\pi_0(\omega)\)

 

定理C(エルゴード性)

省略。論文を参照してください。

 

定理B-証明

 

まずは二つの補題を示す。

(補足:非一様マルコフ連鎖の収束性について述べています。シュミレーテッドアニーリングや経路積分モンテカルロ法による量子アニーリング、Green関数モンテカルロ法による量子アニーリングなど、古典的に量子アニーリングをシュミレートする、本質的に非一様マルコフ連鎖と関連するアルゴリズムの全ての収束性について関係があります。量子アニーリングの収束定理参照。)

 

補題2

全ての\(t_0=0,1,2,…,\)に対して、

\(\lim_{t\rightarrow \infty}\sup _{\omega,\eta^\prime,\eta^{\prime \prime}}|P(X(t)=\omega|X(t_0)=\eta^\prime )-P(X(t)=\omega|X(t_0)=\eta^{\prime \prime}))|=0\)

 (補足:補題2と補題3を証明すれば、定理Bが証明されます。補題2の必要条件として定理Bの条件bが得られます。補題3の必要条件として、条件aが使われます。)

 

補題3

||は\(L_1\)距離を表す。

\(\lim_{t_0\rightarrow \infty}\sup_{t\geq t_0}||P(t|t_0,\eta)-\pi_0||=0\)

 

まず、これらの補題を如何に定理Bに適用するかを示す。\(\sup\)は\(\omega,\eta\)についてとる。

\(\limsup_{t\rightarrow \infty}||P(X(t)=\omega|X(0)=\eta)-\pi_0(\omega)||\)

\(=\limsup _{t_0 \rightarrow \infty }\limsup_{t\rightarrow \infty,t\geq t_0}\)

\(||\sum_{\eta^\prime} P(X(t)=\omega|X(t_0)=\eta) \cdot P(X(t_0)=\eta|X(0)=\omega)-\pi_0||\)

\(\leq \limsup _{t_0 \rightarrow \infty }\limsup_{t\rightarrow \infty,t\geq t_0}\)

\(||\sum_{\eta^\prime} P(X(t)=\omega|X(t_0)=\eta^\prime) \cdot P(X(t_0)=\eta^\prime|X(0)=\omega)-P(X(t)=\omega)|t_0,\pi_0)||\)

\(+||P(X(t)=\omega|X(t_0)=\pi_0)-\pi_0||\)

第二項は、補題3から0になる。第一項について、

 \(||\sum_{\eta^\prime} P(X(t)=\omega|X(t_0)=\eta) \cdot P(X(t_0)=\eta|X(0)=\omega))-P(X(t)=\omega)|t_0,\pi_0)||\)

\(=\sum_\omega \sup_{\eta^{\prime \prime}}|\sum_{\eta^\prime}(P(X(t)=\omega|X(t_0)=\eta^\prime) -P(X(t)=\omega|X(t_0)=\eta^{\prime\prime}) )\)

\(\times (P(X(t)=\eta^\prime|X(t_0)=\eta) -\pi_0(\eta^\prime))|\)

\(\leq 2\sum_\omega \sup_{\eta^\prime,\eta^{\prime\prime}}|P(X(t)=\omega|X(t_0)=\eta^\prime) -P(X(t)=\omega|X(t_0)=\eta^{\prime\prime}) |\)

よって、第一項も補題2から0になるから、補題2,3が成り立てば定理Bを示せることはわかった。

 

補題2の証明

\(\delta (t)=\inf_{1\leq i \leq N,(x_{s_1},…,x_{s_N})\in \Omega} \pi_{T(t)}(X_{s_i}=x_{s_i}|X_{s_j}=x_{s_j},i \neq j)\)

と定義する。

\(\delta (t)\geq \frac{e^{-U^*/T(t)}}{{e^{-U_*/T(t)}}}=e^{-\Delta/T(t)}\)

であることは容易に示せる。

(補足:ベイズの定理を用いる。論文では分母に\(L\)が来ていますが、必要ありません。どちらにせよ結果に影響はありません。)

 

今、\(k\)を固定して考え、\(m_i\)を以下のように定義する。

\(m_i=\sup \{t:t \leq T_k ,n_t=s_i\}\ ,\ 1\leq i \leq N\)

(補足:\(k\)は直観的に言えば、マルコフ連鎖で全てのサイトを今まで何周したかを表すものであることを思い出してください。つまり\(m_i\)は、\(k\)周目で”訪問”するサイト\(i\)の最後に”訪問”する時刻を表しています。)

 また、\(m_1>m_2>m_3>\cdots >m_N\)と仮定する。

(補足:そうなるように最初からマルコフ連鎖を構成すればよいだけです。逆に言えば、シュミレーテッドアニーリングのマルコフ連鎖を構成するときは、各サイトを順番に周って値を更新する必要があることを言っています。)

このとき、

\(P(X(T_k)=\omega|X(t_{k-1})=\omega^\prime)=P(X_{s_1}(m_1)=x_{s_1},\cdots X_{s_N}(m_N)=x_{s_N}|X(t_{k-1})=\omega^\prime)\)

\(=\prod _{j=1}^N P(X_{s_j}(m_j)=x_{s_j}|X_{s_{j+1}}(m_{j+1})=x_{s_{j+1}},\cdots ,X_{s_N}(m_N)=x_{s_N},X(T_{k-1})=\omega^\prime)\)(条件付き分布の定義から)

\(\geq \prod_{j=1}^N \delta(m_j)\)(\(\delta\) の定義から)

\(\geq  \prod_{j=1}^N e^{-\Delta/T(m_j)}\)

\(\geq \exp\{-\Delta N/(t_0+k\tau)\}\)(\(m_j \leq T_k=t_0+k\tau ,j=1,2,…,N\)であり、かつ\(Tは\)減少関数であるから)

\(\geq (t_0+k\tau)^{-1}\)(定理の仮定Bより。逆に言えば、以下の議論を追えば、このアニーリングスケジュールがシュミレーテッドアニーリングの収束の必要条件になることがわかります。)

 

\(t>t_0\)なる\(t\)に対して、\(K(t)\)を

\(K(t)=\sup\{k:T_k<t\}\)と定義する。

(補足:時刻\(t\)ですでに各サイトを何周回ったか)

 

\(\sup _{\omega,\eta^\prime,\eta^{\prime\prime}}|P(X(T_k)=\omega|X(t_0)=\eta^\prime)-P(X(T_k)=\omega|X(t_0)=\eta^{\prime\prime})|\)

\(\leq \prod_{k=1}^{K(t)} (1-(t_0+k\tau)^{-1})\)

 

\(t\rightarrow \infty \Rightarrow K(t) \rightarrow \infty\)であることと、

\(\prod_{k=1}^m (1-(t_0+k\tau)^{-1})\)

\(=\prod_{k=1}^m \exp \log (1-(t_0+k\tau)^{-1})\)

\(=\exp(\sum_{k=1}^m \log (1-(t_0+k\tau)^{-1}))\)

\(\leq \exp(\sum_{k=1}^m (t_0+k\tau)^{-1})\)(\(x>0\)で\(\log (1-x)<-x\)であることから)

\(\rightarrow 0(m \rightarrow \infty)\)

 

これより、補題2は示された。

 

補題3の証明

\(P(t,\omega|t_0,\pi_0)\)

\(=\sum_\eta P(t,\omega|X(t_0)=\eta)\pi_0 (\eta)\)

(補足:初期状態が\(\pi\)によって決定された後、時刻\(t\)で\(\omega\)である確率。証明を構成するためにこのような確率を考えますが、結局初期状態が何であっても収束することが示せればよいので、初期状態の選び方そのものは特殊な状況を考えて問題ないということを言っています。)

 

 \(P(t,\omega|t_0,\pi_0)\)を以下、\(P_{t_0,t}\)と略記する。

始めに、\(t>t_0\geq 0\)に対し、以下を示す。

\(||P_{t_0,t}-\pi_{T(t)}|| \leq ||P_{t_0,t-1}-\pi_{T(t)}||\)

 

簡単のために\(n_t=s_1\)として、

\(||P_{t_0,t}-\pi_{T(t)}||\)

\(=\sum _{(x_{s_1},…,x_{s_N})} |\pi_{T(t)}(X_{s_1}=x_{s_1}|X_s=x_s,s \neq s_1) \cdot P_{t_0,t-1}(X_s=x_s,s \neq s_1)-\pi_{T(t)} (X_s=x_s,s \in S)|\)

\(=\sum _{(x_{s_2},…,x_{s_N})}|P_{t_0,t-1}(X_s=x_s,s \neq s_1)-\pi_{T(t)} (X_s=x_s,s \in S)|\)

\(\leq \sum_{(x_{s_1},…,x_{s_N}) \in \Omega} |P_{t_0,t-1}(X_s=x_s,s \neq s_1)-\pi_{T(t)} (X_s=x_s,s \in S)|\)

\(=||P_{t_0,t-1}-\pi_{T(t)}||\)

 

次に、\(t \rightarrow \infty\)で\(||\pi_0-\pi_{T(t)}||\rightarrow 0\)を示す。

\(\pi_{T(t)} (\omega)=\frac{\exp(-U(\omega)/T(t))}{\sum_{\omega^\prime \in \Omega_0}\exp(-U(\omega)/T(t)+\sum_{\omega^\prime \in \Omega \backslash \Omega_0}\exp(-U(\omega)/T(t)}\)

\(=\frac{\exp(-(U(\omega)-U_*)/T(t))}{|\Omega_0|+\sum_{\omega^\prime \in \Omega \backslash \Omega_0}\exp(-(U(\omega)-U_*)/T(t)}\)

よって、\(t \rightarrow \infty\)の極限で、\(\pi_{T(t)} (\omega)\)は\(\Omega_0\)上の一様分布に近づく。

 

(原論文のこれ以降の説明は、ある一定より先は単調減少して0に近づく数列の和が有界であるなどと主張している(簡単な反例が挙げられる)ことや、\(\pi_T(t)\)の時間依存性についての議論がデタラメなので、修正しています。)

 

\(t>t_0\geq 0\)において、

\(||P_{t_0,t}-\pi_0||\leq||P_{t_0,t}-\pi_{T(t)}||+||\pi_{T(t)}-\pi_0|| \)

\(\leq ||P_{t_0,t_0}-\pi_T(t)||+||\pi_{T(t)}-\pi_0||\)(\(||P_{t_0,t}-\pi_{T(t)}|| \leq ||P_{t_0,t-1}-\pi_{T(t)}||\)であるから)

\(=||P_{t_0,t_0}-\pi_{T(t_0)}||+||\pi_{T(t_0)}-\pi_{T(t)}||+||\pi_{T(t)}-\pi_0||\)

 ここで、\(P_{t_0,t_0}=\pi_0\)であり、\(t\rightarrow \infty\)で\(||\pi_{T(t)}-\pi_0||\rightarrow 0\)であるから、上式の\(\lim_{t_0\rightarrow \infty}\sup_{t\geq t_0}\)をとれば良い。