確率的グラフィカルモデルの変分近似法佐野 崇
東洋⼤学 情報連携学部

⽬次
1. 確率的グラフィカルモデルとは
2. 確率的グラフィカルモデルの推論
3. 確率的グラフィカルモデルの学習

Copyright © 2022 by INIAD

今⽇の講義の⽬標
l 物理屋向けに、確率的グラフィカルモデルの⾯⽩いところ(統計⼒学が応⽤されているところ)を概観して紹介する
l ぜひとも物理屋に解決してもらいたい問題も紹介
l 逆に、この⽅法で物理の研究が進むのでは? という視点で機械学習を紹介するものではない

Copyright © 2022 by INIAD

確率を⽤いたデータのモデル
l 収集したデータを説明するモデルがほしい
l 多くの場合、データは定まった値を持たず、確率的にゆらぐ■ 例:理論的には定まった値が得られるはずでも、⼀般には測定誤差が乗る

l 測定すべき量を確率変数ととらえ、測定されたデータを確率変数の実現値と考える
l データから確率変数の従う確率分布を推定(学習)したい
Copyright © 2022 by INIAD

l 例: ある物理量(⻑さ)を知りたい。しかし、測定誤差が存在する。誤差は正規分布(ガウス分布)すると仮定できる
𝑋−𝜇 !
!
𝑃 𝑋 =
exp −
≡
𝑁(𝜇,
𝜎
)
!
2𝜎
2𝜋𝜎 !
1

• 未知のパラメータは 平均𝜇 と分散 𝜎 !
• 複数回の測定値 𝑥" , 𝑥! , … , 𝑥# が得られたとき、ここから平均𝜇 と分散 𝜎 ! を推定する• 測定値が独⽴に、同じ正規分布に従うのであれば、それらが同時に測定される確率は#𝑃 𝑥! , 𝑥" , … , 𝑥# = & 𝑃(𝑋 = 𝑥$ )$%!
この尤度関数を最⼤化するようにパラメータを決定することを最尤推定という

Copyright © 2022 by INIAD

l 尤度関数の代わりに対数尤度を最⼤化しても同じ#
log 𝑃 𝑥! , 𝑥" , … , 𝑥# = - ln 𝑃(𝑋 = 𝑥$ )$%!
平均と分散について1階微分をゼロと置くと次を得る#

1
𝜇 = 5 𝑥$
𝑛
$%"

#

1
𝜎 = 5(𝑥$ −𝜇)!
𝑛
!

$%"

最尤推定によって、パラメータを求めることができた
これは(広い意味で)学習である(よく⾒ると不変推定量ではない)
Copyright © 2022 by INIAD

多変数の場合
l ⼀般には、多数の測定値が同時に得られる。それらを確率変数の実現値だと思うと、多数の確率変数を持つ同時確率分布が背後にある、と考えることができる
l N個の変数を考えるなら 𝑃(𝑋! , 𝑋" , … , 𝑋# )l 仮に各変数が2値(0か1)を取りうるとしても、この確率分##
布は 2 個の確率値を保存しておく必要がある( 𝒪(2 ) のパラメータが必要)
l ⼀般的には⾮現実的
Copyright © 2022 by INIAD

因⼦化の仮定l 多くの場合は、お互いに関係のない(独⽴な)確率変数がある
l それらの同時確率分布は、単独の確率分布の単純な積l 例: 2個の確率変数 𝑋6 , 𝑋7 が独⽴であれば𝑃 𝑋6 , 𝑋7 = 𝑃 𝑋6 𝑃(𝑋7 )Copyright © 2022 by INIAD

因⼦化の図形的表現
l 例えば4個の確率変数 𝑋6 , 𝑋7 , 𝑋8 , 𝑋9 があって、独⽴でない組み
合わせは𝑋6 , 𝑋7 、 𝑋7 , 𝑋8 、 𝑋8 , 𝑋9 のみであるとする。同時確率
分布は
𝑃 𝑋6 , 𝑋7 , 𝑋8 , 𝑋9 = 𝑃 𝑋6 , 𝑋7 𝑃(𝑋7 , 𝑋8 )𝑃(𝑋8 , 𝑋9 )と因⼦化できる。このとき、確率変数を頂点(ノード)とし、独⽴でない組み合
わせの間に辺(エッジ)を書き⼊れてグラフとして表現すると...

𝑋!

𝑋"

𝑋#

𝑋$

同時確率分布がどのように因⼦化されるかをグラフとして表したものを確率的グラ
フィカルモデルという。これは後述するマルコフ確率場の⼀例
Copyright © 2022 by INIAD

ベイジアンネットワーク

l いろいろな確率変数が実現値を取るときに、背後の現象を反映していることがあるl 例: 以下の4つの2値(真または偽)確率変数を考える■ X! : 天気が⾬である■ X" : 気温が⾼い■ X# : スプリンクラーが作動している■ X$ : 芝⽣が濡れている

l この4つには因果関係があるはず■ 気温が⾼いとスプリンクラーが作動する(X# はX" に依存)■ スプリンクラーが作動するか⾬が降ると芝⽣が濡れる(X$ は X# とX! に依存)

l この因果関係を図に書き表してみる

Copyright © 2022 by INIAD

ベイジアンネットワーク

⾬

𝑋!
𝑋$

𝑋"

⾼気温

𝑋#

スプリンクラー

濡れた芝⽣

𝑃 𝑋! , 𝑋" , 𝑋& , 𝑋' = 𝑃 𝑋' 𝑋! , 𝑋& 𝑃 𝑋& 𝑋" 𝑃(𝑋! )𝑃(𝑋" )このような確率変数の因果関係(条件付き確率による同時確率分布の因⼦化)を図形的に表したものをベイジアンネットワークというCopyright © 2022 by INIAD

ベイジアンネットワークの⾒⽅
l ある変数に⽮印を向けている変数を親変数というl ⼦変数の確率は、親変数集合の条件付き確率になる
⾬

𝑋!

𝑋#
𝑋$

スプリンクラー

濡れた芝⽣

Copyright © 2022 by INIAD

𝑃(𝑋9 |𝑋6 , 𝑋8 )l すべての親⼦関係を条件付き確率に変換
l 親のない変数には事前分布を割り当てる
l それらの積が同時確率分布
𝑃(𝑋8 |𝑋7 )𝑋"
𝑋!

𝑃(𝑋9 |𝑋6 , 𝑋8 )𝑋#
𝑋$

𝑃 𝑋! , 𝑋" , 𝑋# , 𝑋$ = 𝑃 𝑋$ 𝑋! , 𝑋# 𝑃 𝑋# 𝑋" 𝑃(𝑋! )𝑃(𝑋" )Copyright © 2022 by INIADASIA ネットワーク• 肺疾患に関するベイジアンネット
ワーク• 例えば dyspnea(呼吸困難)、XRay
result(所⾒)のあるなしといった観測を与えたときに、肺がんである確率
がどのくらいあるかを推論できる

Copyright © 2022 by INIAD

無向グラフ: マルコフ確率場

l 辺で結ばれた変数は因⼦化できず、同じ関数(ポテンシャル関数)に所属するとする
l 同じポテンシャルに所属する変数同⼠を向きのない辺で結ぶ
𝜓"! (𝑋" , 𝑋! ) 𝜓!& (𝑋! , 𝑋& ) 𝜓&' (𝑋& , 𝑋' )

𝑋!

𝑋"

𝑋#

𝑃 𝑋! , 𝑋" , 𝑋# , 𝑋$
1
= 𝜓!" 𝑋! , 𝑋" 𝜓"# 𝑋" , 𝑋# 𝜓#$ (𝑋# , 𝑋$ )𝑍

𝑋$

Zは確率を規格化するための分配関数
𝑍 = 5 𝜓"! 𝑋" , 𝑋! 𝜓!& 𝑋! , 𝑋& 𝜓&' (𝑋& , 𝑋' )𝑿
Copyright © 2022 by INIAD

マルコフ確率場の例:ボルツマンマシン
l ホップフィールドネットワーク(連想記憶モデル)の確率版
l 例えば2次元の画像画像の各ピクセルが2値確率変数であり、その確率が隣接するピクセルの状態で決まるとする
𝑋!!

𝑋!"

𝑋!#

𝑋"!

𝑋""

𝑋"#

𝑋#!

𝑋#"

𝑋##

MNIST hand-written digits...
Copyright © 2022 by INIAD

...

𝑋!!

𝑋!"

1
𝑃 𝑿 = = 𝜓$* (𝑋$ , 𝑋* ) = 𝜓$ (𝑋$ )
𝑍

𝑋!#

𝑋"!

𝑋""

𝑋"#

𝑋#!

𝑋#"

𝑋##

...

1
= exp
𝑍

$,* ∈,$∈-

5 𝑊$* 𝑋$ 𝑋* + 5 𝑏$ 𝑋$$,* ∈,

1≡ exp −𝐸(𝑿)𝑍$∈-

𝑿 = {𝑋1 , … , 𝑋𝑁 } :2値確率変数(0または1)
𝑊𝑖𝑗 , 𝑏𝑖 :パラメータ (重みとバイアス)
𝐸: 辺(edge)集合、𝑉:頂点(vertex)集合

l 画像の⼀部を復元したり、ノイズの除去が可能
l Ising (spin glass)モデルと等価
Copyright © 2022 by INIAD

因⼦グラフl より⼀般には、3つ以上の確率変数が同じポテンシャル関数に所属しても良いl その場合は変数頂点と因⼦頂点を持つ2部グラフとして表現できる(因⼦グラフ)

𝑋!

𝑋#
𝜓"

𝜓&

𝑋"

𝜓!

1
𝑃 𝑿 = + 𝜓% (𝑿% )
𝑍

𝑋$

%∈'

1
= 𝜓! (𝑋! , 𝑋" , 𝑋# )𝜓" (𝑋# , 𝑋$ )𝜓# (𝑋! )𝑍
𝐹: 因⼦の集合
𝑿𝑓 : 因⼦fに属する変数の集合• マルコフ確率場やベイジアンネットワークは因⼦グラフの特別な場合• 因⼦グラフで⼀般的な理論を作り、そこからマルコフ確率場やベイジアン
ネットワークに適⽤することもよく⾏われる
Copyright © 2022 by INIAD

統計⼒学との対応
l 因⼦グラフの因⼦を指数の肩にあげてやると...
1
1
𝑃 𝑿 = = 𝜓. (𝑿. ) = exp −𝐸 𝑿
𝑍𝑍.∈/

ここで 𝐸 𝑿 = ∑.∈/ 𝐸. (𝑿. ), 𝐸. 𝑿. = −log 𝜓. (𝑿. )• すなわち、状態𝑋に対してエネルギー(ハミルトニアン)が𝐸 𝑿 と与えれられたときに、 𝑋はボルツマン分布に従う(温度=1)• すぐに分かるように、この対応から統計⼒学的⼿法をグラフィカルモ
デルに使うことができる(同じ困難を共有しているとも⾔える)

Copyright © 2022 by INIAD

結局確率的グラフィカルモデルとは?
l 確率分布に何らかの因⼦化を仮定すると、それをグラフとして表現することができる
l 極論すれば、すべての確率分布(確率モデル) は確率的グラフィ
カルモデル、とも⾔える
l 因⼦化を仮定してから、それをグラフに表現して考察したり、逆にグラフを作ってから、対応する確率分布を書き下したり、といった使い⽅をする
𝑋!

𝑃 𝑋" , 𝑋! , … , 𝑋0
1
= 𝜓(𝑋" , 𝑋! , … , 𝑋0 )
𝑍

...

𝜓

𝑋"
Copyright © 2022 by INIAD

𝑋#

確率的グラフィカルモデルの推論

Copyright © 2022 by INIAD

確率的グラフィカルモデルの推論

l 確率的グラフィカルモデルが与えられたとする(条件付き確率分布やポテン
シャル関数が定まっているとする)
l このとき、しばしば「最も⾼い確率値を持つ確率変数の実現値」に興味がある

𝑋 = argmax𝑿 𝑃(𝑿)

• 実際上は、いくつかの確率変数が観測されたときに、観測されていない変数について推定を⾏いたい• 呼吸困難の所⾒があり、アジアへの訪問歴のある⼈は肺がんであるか?
• 画像の半分が⽋落したときに、残りの部分から復元を⾏う

5)
𝑋 = argmax𝑿 𝑃(𝑿|𝑿• MAP(Maximum a posteriori:最⼤事後確率)推定という• 素朴にやると、2値変数がN個あるときに𝒪(20 )の確率値の⽐較が必要であり、困難
Copyright © 2022 by INIAD

肺がんである(L)といえるか? あるいはその確率は?

X線所⾒(X)、アジアの渡航歴(A) が分かったとき

ASIAネットワーク(Yedidia 2001より)

Copyright © 2022 by INIAD

周辺確率
l MAP推定の他にも、注⽬する変数以外の変数を積分してしまった周辺確率に興味があることもある■ 例: 呼吸困難の所⾒があり、アジアへの訪問歴のある⼈が肺がんである確率は?

𝑃" (𝑋" ) = Σ𝑿 #{%$ } 𝑃(𝑿)
• こちらも素朴にやると、𝒪(20 )の和が必要• MAP推定とよく似ている(maxと和の置き換え)が、周辺確率のargmaxを取ってもMAPが得られるとは限らない

Copyright © 2022 by INIAD

ベイジアンネットワークの例

⾬

𝑋!
𝑋$

𝑋"

⾼気温

𝑋#

スプリンクラー

濡れた芝⽣

𝑃 𝑋' , 𝑋( , 𝑋) , 𝑋*
= 𝑃 𝑋* 𝑋' , 𝑋) 𝑃 𝑋) 𝑋( 𝑃(𝑋' )𝑃(𝑋( )• 実際の推論はいくつかの変数が観測されたときに、観測されていない変数の値や確率に興味がある• 芝⽣が濡れている(𝑋' = 1)ときに、⾬である確率𝑃(𝑋" = 1)は?
• 以下では簡単のため、観測された変数は陽には現れないが、観測があっても理論は同じ
Copyright © 2022 by INIAD

統計⼒学との関係
l 統計⼒学(や場の理論)で問題になるのは、ある変数(例えばスピン)が、与えられた確率分布のもとでどのような期待値を持つかだった
6< 𝑠G >= Σ𝒔 𝑠G 𝑃 𝒔 = ∑𝒔 exp[−𝐸 𝒔 ]I• これは周辺確率を求める問題とほぼ同じ• よく知られているように、⼀般にはこの 𝒪(20 )の和を取ることはできない• そのため、物理屋は⾊々な近似法を考案してきた• サンプリング法(モンテカルロ法)や変分近似法(平均場近似)はどちらも統計的推論に応⽤されている• 今回は特に、変分近似法について⾒ていく
Copyright © 2022 by INIAD

変分近似法l ⼿に負えない確率分布𝑃(𝑋)に対して、より扱いやすい確率分布𝑄(𝑋)を導⼊する
l 𝑄(𝑋)をできるだけ𝑃(𝑋)に近づける
l 2つの確率分布の近さとして、 Kullback-Leibler divergence (KL情報量)を⽤いる
𝑄(𝑿)
KL(𝑄|𝑃) ≡ 5 𝑄(𝑿) log
𝑃(𝑿)1

KL情報量は距離ではないが、距離に似た性質をもつ• KL 𝑄 𝑃 ≥ 0• KL 𝑄 𝑃 = 0 ↔ 𝑄 𝑿 = 𝑃 𝑋 ∀ 𝑿
ただし⼀般に KL 𝑄 𝑃 ≠ KL 𝑃 𝑄※ 以下では主に離散確率変数を扱うが、和を積分に置き換えれば連続確率変数でもほぼ同様のことが成⽴
Copyright © 2022 by INIAD

変分分布
l 𝑄を適当にパラメータ付けして、KL情報量を最⼩化するようにパラメータを求めたいl KL情報量をパラメータで変分(微分)して最⼩値を探すので、変分近似法という。𝑄を変分分布というl よく使われる𝑄の形は、次の形の完全に因⼦化されたもの*

𝑄(𝑿) = + 𝑞( (𝑋( )()!

• これを平均場近似という• KL情報量が最⼩になるように、パラメータ𝑞( を決定する。 𝑞( の満たす⽅程式は𝑞$ (𝑋$ ) ∝ exp[ 5 = 𝑞* 𝑋* ln 𝑃] 𝑋! 以外でln 𝑃の期待値を取ると、ln 𝑞! に⼀致する
𝑿21! *3$
Copyright © 2022 by INIAD

平均場近似?
l 統計⼒学を勉強すると平均場近似が出てくるが、それとはかなり違うように⾒えるので、チェックしてみる
l 簡単のためIsingモデルを考えることにして、その教科書
的な平均場近似を考察する
l 相互作⽤と外場が局所的に異なるspin glassモデルを考えてみる
0

𝐸 = − 5 𝐽$* 𝑠$ 𝑠* − 5 ℎ$ 𝑠$
4$*5

$%"

< 𝑖𝑗 > は接続のあるスピンの組み合わせの集合

Copyright © 2022 by INIAD

l よくある説明は、スピン変数を、平均値𝑚! とその周りのゆらぎに分解し、ゆらぎを⼩さいと考える
𝑠! = 𝑚! + 𝛿𝑠!
𝐸 = − 5 𝐽$* 𝑚$ + 𝛿𝑠$ 𝑚* + 𝛿𝑠* − 5 ℎ$ 𝑠$
4$*5$≅ − 5 𝐽$* 𝑚$ 𝑚* + 𝛿𝑠$ 𝑚* + 𝑚$ 𝛿𝑠* − 5 ℎ$ 𝑠$
4$*5$= − 5 𝐽$* 𝑠$ 𝑚* + 𝑚$ 𝑠* − 𝑚$ 𝑚* − 5 ℎ$ 𝑠$ ≡ 𝐸6/
4$*5$

近似ハミルトニアンは 𝑠$ の1次のみになった。つまり($)
𝑃(𝒔) ∝ exp −𝐸6/ = = exp[−𝐸6/ (𝑠$ )]$と全体の確率が1体の確率の積の形で書ける à Qの因⼦化と同じ!
Copyright © 2022 by INIAD

𝛿𝑠の2次を無視

l 次に、KL情報量の最⼩化の意味を考えてみる。⼀般の確率分布Qに対して
𝑄 𝑿
KL 𝑄 𝑃 = 5 𝑄 𝑿 ln
= 5 𝑄 𝑿 (ln 𝑄 𝑿 − ln 𝑃(𝑿))
𝑃 𝑿
𝑿

𝑿

"

ここで、𝑃 = exp(−𝐸)であることを使うと9KL 𝑄 𝑃 = 5 𝑄 𝑿 ln 𝑄 𝑿 + 5 𝑄 𝑿 (𝐸 + ln 𝑍)
𝑿

𝑿= −𝐻 𝑄 + 5 𝑄 𝑿 𝐸 𝑿 + ln 𝑍
𝑿

ただし 𝐻 𝑄 = − ∑𝑿 𝑄(𝑿) ln 𝑄(𝑿) : エントロピー関数
Copyright © 2022 by INIAD

まとめるとKL 𝑄 𝑃 = −𝐻 𝑄 + 5 𝑄 𝑿 𝐸 𝑿 + ln 𝑍
𝑿

=𝐺 𝑄 −𝐹

𝐹 = − ln 𝑍: (真の)⾃由エネルギー𝐺 = ∑𝑿 𝑄 𝑿 𝐸 𝑿 − 𝐻 𝑄 : 変分⾃由エネルギー
すなわち、KL情報量の最⼩化とは、変分分布𝑄を動かして、変分⾃由エネル
ギー𝐺を最⼩化し、できるだけ真の⾃由エネルギーに近づけようすることだった• KL情報量が⾮負であることから、𝐺 ≥ 𝐹, つまり変分⾃由エネルギーは真の⾃由エネルギーの上界であることも分かる• 等号が成⽴する𝐺 = 𝐹のは𝑄が𝑃に⼀致するときであるが、Qの形に制限を置いてしまうと、⼀般に等号は成⽴しない
Copyright © 2022 by INIAD

l Isingモデルの場合に、この平均場近似が、self-consistency ⽅程式を導くことをチェックする
l 変分分布を次の形に仮定
0

0

$%"

$%"

1
𝑄(𝒔) = = 𝑞$ (𝑠$ ) = = (1 + 𝑠$ 𝑚$ )2
このとき変分⾃由エネルギー𝐺は
𝐺 = − ; 𝐽!$ 𝑚! 𝑚$ − ; ℎ! 𝑚! − ; 𝐻[𝑞! ]
#!$%

!

!

𝑚! について微分して0とおくと
1 1 + 𝑚!
− ; 𝐽!$ 𝑚$ − ℎ! + ln
=0
2 1 − 𝑚!
$∈' !

Copyright © 2022 by INIAD

𝑁 𝑖 : iと接続のある変数の集合

l 逆双曲線関数 tanh

a6

6
6bc
𝑥 = lnであるから
76ac

0 = − ; 𝐽!$ 𝑚$ − ℎ! + tanh() 𝑚!
$∈' !

𝑚$ = tanh( 5 𝐽$* 𝑚* + ℎ$ )*∈0 $特に、相互作⽤や外場が⼀様な場合(𝐽!$ = 𝐽, ℎ! = ℎ)、期待値も⼀様 𝑚! = 𝑚 であるから

𝑚 = tanh(𝑧𝐽𝑚 + ℎ)𝑧 = ∑𝑗∈𝑁 𝑖 1: 変数 𝑠! と相互作⽤を持つ変数の数

有名なself-consistency ⽅程式が再現できた!

Copyright © 2022 by INIAD

平均場近似を解くアルゴリズム
l 与えられた確率分布Pに対して、
𝑞! (𝑋! ) ∝ exp[ ; L 𝑞$ 𝑋$ ln 𝑃]
𝑿(-! $.!

が成⽴するような𝑞! (𝑋! ) (𝑖 = 1, … , 𝑁)を⾒つける• それには、適当な値でln 𝑞! (𝑋! )を初期化してln 𝑞! (𝑋! ) ←; L 𝑞$ 𝑋$ ln 𝑃
𝑿(-! $.!

のように反覆更新していくのが簡便(規格化は適当に⾏う)• 変分⾃由エネルギー𝐺は𝑞! (𝑋! )について凸であるため、この更新は収束性が⾮常に良い• なお、平均場近似では、周辺確率を最⼤化する実現値は、この近似でのMAP推定値でもある
Copyright © 2022 by INIAD

グラフ構造を利⽤する
l 平均場近似はどのようなモデルでも⽤いることができるが、もっと積極的にグラフ構造を⽤いた近似⽅法がないだろうか
l 例えば、以下のようなマルコフ確率場を考えてみる
𝑋!

𝑋"

𝑋#

𝑋$

...

𝑋*

𝑃 𝑋 ∝ 𝜓!" 𝑋! , 𝑋" 𝜓"# 𝑋" , 𝑋# … 𝜓*,!* (𝑋*,! , 𝑋* )いま、周辺確率𝑃(𝑋/ )を求めたいとする。
闇雲に計算すると、N-1個の確率変数について和を取ることになり、𝑂(2' )の計算が必要になる
Copyright © 2022 by INIAD

l ⼀⽅で、グラフ構造にそって端から計算を⾏うとどうなるか?
l まず𝑋! について周辺化
𝑃 𝑋! , … , 𝑋0 ∝ Σ1" 𝜓"! 𝑋" , 𝑋! 𝜓!& 𝑋! , 𝑋& … 𝜓02"0 𝑋02" , 𝑋0
= 𝜇: (𝑋! )𝜓!& 𝑋! , 𝑋& … 𝜓02"0 𝑋02" , 𝑋0
ここで 𝜇𝛼 (𝑋2 ) = Σ𝑋1 𝜓12 𝑋1 , 𝑋2 で、𝑂(2)で計算できる• 次に𝑋4 で周辺化

𝑃 𝑋& , … , 𝑋0 ∝ Σ1# 𝜇: (𝑋! )𝜓!& 𝑋! , 𝑋& … 𝜓02"0 𝑋02" , 𝑋0
= 𝜇: (𝑋& )𝜓&' 𝑋& , 𝑋' … 𝜓02"0 𝑋02" , 𝑋0

ここで 𝜇𝛼 (𝑋5 ) = Σ𝑋# 𝜇6 (𝑋4 )𝜓45 𝑋4 , 𝑋5 である。計算量は𝑂(24 ) (2x2⾏列とベクトルの積)• これを繰り返すと𝜇𝛼 (𝑋𝑛 ) = Σ𝑋𝑛−1 𝜇𝛼 (𝑋𝑛−1 )𝜓𝑛−1,𝑛 𝑋𝑛−1 , 𝑋𝑛 で𝜇: (𝑋# )を順に計算していき
𝑃 𝑋𝑛 , … , 𝑋𝑁 ∝ 𝜇𝛼 (𝑋𝑛 )𝜓𝑛,𝑛+1 𝑋𝑛 , 𝑋𝑛+1 … 𝜓𝑁−1,𝑁 𝑋𝑁−1 , 𝑋𝑁 を得る
Copyright © 2022 by INIAD

l 𝜇" は、グラフの左側からの「メッセージ」と解釈できる
𝑋!

𝑋#

𝑋"
𝜇" 𝑋#

𝜇" 𝑋$

𝑋-

...
𝜇" 𝑋&

...

𝜇" 𝑋%

変数がメッセージを受け取ると、ポテンシャル関数をかけ、左側の変数について和をとって、⾃分のメッセージとして次の変数に送り出す

𝜇𝛼 (𝑋𝑛 ) = Σ𝑋𝑛−1 𝜇𝛼 (𝑋𝑛−1 )𝜓𝑛−1,𝑛 𝑋𝑛−1 , 𝑋𝑛
メッセージを使って確率が表せる

𝑃 𝑋𝑛 , … , 𝑋𝑁 ∝ 𝜇𝛼 (𝑋𝑛 )𝜓𝑛,𝑛+1 𝑋𝑛 , 𝑋𝑛+1 … 𝜓𝑁−1,𝑁 𝑋𝑁−1 , 𝑋𝑁

Copyright © 2022 by INIAD

l 同じことを、今度は𝑋# から逆向きに計算していく
𝑋!

𝑋#

𝑋"
𝜇" 𝑋#

𝜇" 𝑋$

𝑋-

...
𝜇" 𝑋&

𝜇" 𝑋%

𝑋*,!

...

𝜇' 𝑋%

𝜇' 𝑋()#

メッセージの計算式:

𝜇B (𝑋𝑛 ) = Σ𝑋𝑛'1 𝜇B (𝑋𝑛C1 )𝜓𝑛,𝑛C" 𝑋𝑛 , 𝑋𝑛C"
最終的に求めたい周辺確率は

𝑃 𝑋𝑛 ∝ 𝜇𝛼 (𝑋𝑛 )𝜇B (𝑋# )メッセージをすべて計算するのにかかる計算量は𝑂 𝑁 24 、それらの積を取ることですべての変数の周辺確率を求められる!
確率伝播法という
Copyright © 2022 by INIAD

𝜇' 𝑋()*

𝑋*

l もしグラフに枝分かれがあっても、メッセージの集約を⾏え新しいメッセージを計算できる
...

𝑋-

𝑋(

...
𝜇" 𝑋%,*

𝜇" 𝑋%

𝜇+ 𝑋%

𝑋.
Σ-$ 𝜇6 𝑋/ 𝜇7 𝑋/ 𝜓/,/9) 𝑋/ , 𝑋/9) 𝜓/9),/94 𝑋/9) , 𝑋/94 …
= 𝜇6 (𝑋/9) )𝜓/9),/94 𝑋/9) , 𝑋/94 …

𝜇𝛼 𝑋𝑛C" = Σ1( 𝜇: 𝑋# 𝜇D 𝑋# 𝜓#,#C" 𝑋# , 𝑋#C"
Copyright © 2022 by INIAD

l 同じアイデアは、もっと複雑なグラフにも適⽤できる
l ⼀般に因⼦グラフの場合を紹介する
l 𝐸を因⼦の集合、𝑉を変数の集合とし、a,b,c...を因⼦の添字、i,j,k...を変数の添字に使う
l メッセージには、因⼦から変数へのメッセージ 𝑀!" と変数から因⼦へのメッセージ 𝑀"! がある。

𝑋!
𝑀*. (𝑋*)...

𝑀.- (𝑋- )𝑋(

𝑋*, … , 𝑋- , … , 𝑋(! が所属する因⼦𝜓. から変数𝑋- へ送られるメッセージ𝑀.- (𝑋- )は、変数から因⼦へ送られるメッセージを⽤いて

𝜓E (𝑋" , … , 𝑋$ , … , 𝑋0) )𝑀(! . (𝑋(! )𝑋*#

𝑀:! 𝑋! =;{-% ,…,-&' }(-!

Copyright © 2022 by INIAD

𝜓: (𝑋) , … , 𝑋! , … , 𝑋'' ) L 𝑀$: (𝑋$ )$.!

l 変数𝑋! から因⼦𝜓$ のメッセージは、 𝜓$ 以外の因⼦からやってきたメッセージの積𝜓"𝑀!: 𝑋! = L 𝑀>! (𝑋! )

...
...

𝑀*- (𝑋- )𝑀-. (𝑋- )𝑋(

𝑀(" - (𝑋- )

>.:

𝜓0!

𝜓E

Copyright © 2022 by INIAD

確率伝播法
l 最終的に、メッセージを⽤いて、周辺確率は次のように計算できる
𝑃 𝑋! ∝ L 𝑀:! (𝑋! ):∈'(!)𝑃 𝑋) , … 𝑋'' ∝ 𝜓: (𝑋) , … , 𝑋'' ) L 𝑀!: (𝑋! )!∈'(:)

このような計算⽅法を確率伝播法(Belief Propagation, BP)や積和アルゴリズムと呼ぶ
特に1988年にPearlによってベイジアンネットワークの推論⽅法として導⼊され有名になった
命題: 確率伝播法は、グラフ構造が⽊のとき(閉路を含まないとき)、端のノードから順
番にメッセージを計算していくことで、厳密な周辺確率を求めることができる
隠れマルコフモデルの推論アルゴリズムなども、実は確率伝播法と同等
Copyright © 2022 by INIAD

ベイジアンネットワークの例

⾬

𝑋!

メッセージ
(すべての辺を流れる)

𝑋$

𝑋"

⾼気温

𝑋#

スプリンクラー

濡れた芝⽣

𝑃 𝑋' , 𝑋( , 𝑋) , 𝑋*
= 𝑃 𝑋* 𝑋' , 𝑋) 𝑃 𝑋) 𝑋( 𝑃(𝑋' )𝑃(𝑋( )𝑋4 = 1(濡れている) :値が確定• ベイジアンネットワークでも確率伝播法は(もちろん)実⾏できる• 値が確定しているノードの情報がメッセージによって伝わり、その他の確率変数の周辺確率が求まる• メッセージは⽮印の順⽅向にも逆⽅向にも流れる

Copyright © 2022 by INIAD

ループあり確率伝播法
l 確率伝播法で⾏う計算は局所的なものなので、グラフにループがあっても実⾏は可能
l もちろん収束するとは限らない
l 収束すれば近似解が得られると期待できるループあり確率伝播法 (Loopy Belief Propagation)
1. メッセージを適当に初期化する
2. 右の式でメッセージを収束するまで更新
3. 収束したメッセージを⽤いて近似周辺確率を計算

𝑀.- 𝑋- =2{0# ,…,0$! })0"

𝑀-. 𝑋- = 3 𝑀6- (𝑋- )𝑃 𝑋- ∝ 3 𝑀.- (𝑋- )65.

.∈((-)𝑃 𝑋*, … 𝑋(! ∝ 𝜓. (𝑋*, … , 𝑋(! ) 3 𝑀-. (𝑋- )-∈((.)
Copyright © 2022 by INIAD

𝜓. (𝑋*, … , 𝑋- , … , 𝑋(! ) 3 𝑀4. (𝑋4 )45-

X線所⾒(X)、アジアの渡航歴(A) が分かったとき

肺がんである(L)といえるか? あるいはそ
の確率は?

メッセージ
(すべての辺を流れる)

ここにループASIAネットワーク(Yedidia 2001より)

Copyright © 2022 by INIAD

l ループあり確率伝播法は、ループが少ないグラフィカルモデル
では良い近似を与える
l 例えばTurbo codeというエラー訂正符号化アルゴリズムは、実はループあり確率伝播法である

McEliece+ 1998 より
Copyright © 2022 by INIAD

l ループあり確率伝播法はヒューリスティックに導⼊されたが、すこし時代が下ると、理論的背景が⾊々と明らかにされた
l 驚くべきことに、ループあり確率伝播法とは、統計⼒学で有名なBethe近似と同じものであった
l Bethe近似は平均場近似の改良であり、ゆらぎの⾼次の項を取り込む

Copyright © 2022 by INIAD

確率伝播法の理論

Yedidia+ 2001

l 表記を簡単にするため、マルコフ確率場の場合を考える
l 確率伝播法は、もし収束すれば、以下のBethe⾃由エネルギーの極値になる𝑄を与える
𝐸: 辺の集合
𝐹FGHIG = 5 𝑄 𝑿 𝐸 𝑿 − 𝐻 𝑄
𝑿

ただし

𝑄(𝑿) =∏(!$)∈A 𝑞!$ (𝑋! , 𝑋$ )∏!∈B 𝑞! 𝑋! C!()

𝑉: 変数の集合
𝑑- : 𝑋- に隣接する変数の数

𝑞- 𝑋- = Σ0% 𝑞-4 (𝑋- , 𝑋4 )• 実は、これはBethe近似の⽬的関数に他ならない• 平均場近似と⽐較すると、変分分布𝑄(𝑿)の因⼦化が異なる。Bethe
近似では隣接する変数の相関も取り込み、最適化を⾏う• この因⼦化のせいで、Bethe⾃由エネルギーはもはや凸ではない• Bethe近似はBethe格⼦(⽊構造)で厳密になるが、このことは確率伝播法が⽊構造のグラフ上で厳密であることと対応している

cf. 平均場:𝑄(𝑿) = ∏𝑁𝑖=1 𝑞𝑖 (𝑋𝑖 )https://en.wikipedia.org/wiki/Bethe_lattice
Copyright © 2022 by INIAD

平均場近似と確率伝播法

l 平均場近似も確率伝播法(Bethe近似)も、どちらも変分分布を仮定して、真の確率分布とのKL情報量が⼩さくなるように(=変分⾃由エネ
ルギー𝐺が⼩さくなるように)、変分パラメータ𝑞% , 𝑞%& を決定している
𝑁

𝐺[𝑄] = 5 𝑄 𝑿 𝐸 𝑿 − 𝐻 𝑄

𝑄(𝑿) = 3 𝑞𝑖 (𝑋𝑖 )𝑖=1

𝑿

𝑄(𝑿) =∏(!$)∈A 𝑞!$ (𝑋! , 𝑋$ )∏!∈B 𝑞! 𝑋! C!()

平均場
確率伝播法

Welling & Teh 2003

• どちらも変分近似法の⼀種として統⼀的に理解できた• すると、いろいろな変分分布Qを⽤いることで、より良い近似ができそうである• 実際、グラフ構造をクラスタに分けたKikuchi近似を元にした変分分布や、もっと⾃由に設定された⼀般化確率伝播法も提案されている
Yedidia+ 2000, 2001

Copyright © 2022 by INIAD

確率伝播法の実際
l Bethe近似が平均場近似の改善であることから、常に確率伝播法を使えばよいのでは? と思われるかもしれない
l しかしながら、少しループの多いグラフィカルモデルでは確率伝播法は収束性が⾮常に悪く、また収束しても精度の悪い解しか得られないことが多い■ 経験的には、変数が数百くらいだと厳しい(グラフ構造にもよるが)l 平均場近似は⾃由エネルギーが凸であるため、収束性は⾮常に良い。そのために、複雑なモデルでは平均場近似が⾮常に多く⽤いられているのが現状である

Copyright © 2022 by INIAD

確率的グラフィカルモデルの学習

Copyright © 2022 by INIAD

グラフィカルモデルの学習とはl 推論とは、グラフィカルモデルのパラメータ(ポテンシャル関数)が分かっているときに、確率変数の周辺確率やMAPを推定することだった
l 学習はその逆で、データとして確率変数の実現値(観測値)が与えられたときに、グラフィカルモデルのパラメータを求めることである
l Isingモデルだと、パラメータ𝐽%& , ℎ% が与えられたときに期待値< 𝑠% >, < 𝑠% 𝑠& >を求めるのが推論であったが、逆にデータから計算される
統計量 < 𝑠% >'()( , < 𝑠% 𝑠& >'()( を⽤いて𝐽%& , ℎ% を求めるのが学習である
l その意味で、学習のことを逆問題と呼ぶこともある

Copyright © 2022 by INIAD

すべてのデータが与えられた場合
l 確率的グラフィカルモデルのすべての確率変数について、実現値の観測データが与えられる場合を考える
l 尤度関数を最⼤化することで、パラメータを学習(推定)することにしよう(最尤推定)@

ln 3 𝑃 𝑿 𝒅>?*@= 2 ln 𝑃(𝑿 𝒅 )>?*𝑿 𝒅 : 確率変数の観測値。𝑑 = 1, … , 𝐷• マルコフ確率場や因⼦グラフには分配関数𝑍が含まれる• 分配関数の計算には指数的なコストがかかるため、尤度関数や、そのパラメータ微分を求めることは容易ではない
𝑍 = 2 3 𝜓. (𝑿. )• ここで変分近似法を⽤いることができる
𝑿 .
Copyright © 2022 by INIAD

Isingモデルの学習
l Isingモデルについて、対数尤度関数をパラメータ𝐽!% 微分
してみる
@

@

𝜕 1
𝜕 1
𝒅
2 ln 𝑃(𝒔 ) =2 𝐸 𝒔𝒅𝜕𝐽-4 𝐷𝜕𝐽-4 𝐷>?*@− ln 𝑍>?*1𝜕= 2 𝑠- 𝑠4 −ln 𝑍 =< 𝑠- 𝑠4 >CDED−< 𝑠- 𝑠4 >
𝐷
𝜕𝐽-4
>?*

この微分=0と置くと< 𝑠- 𝑠4 >CDED=< 𝑠- 𝑠4 >ℎ- 微分からも同様に< 𝑠- >CDED=< 𝑠- >moment matching 条件• 左辺はデータを集計した平均値、右辺はIsingモデルで計算した期待値である• 左辺は容易に計算できるが、右辺は⼀般には計算が⼤変である• しかし、推論で⽤いた変分近似法を⽤いれば近似的に計算可能である
Copyright © 2022 by INIAD

平均場近似を⽤いた学習< 𝑠- 𝑠4 >CDED=< 𝑠- 𝑠4 >< 𝑠- >CDED=< 𝑠- >moment matching 条件

l アイデア: 右辺を平均場近似によってパラメータの関数として求めてしまい、これらを連⽴してパラメータについて解いてしまえば近似学習が可能である
l しかし、相関< 𝑠" 𝑠$ > は平均場近似で計算できるのか? 𝑚" 𝑚$ のようになってしまわないか?
l 実は、線形応答関係式を⽤いると、平均場近似でも< 𝑠" 𝑠$ > が計算できる
𝜕#𝜕𝜕
𝜕
ln 𝑍 =
ln 𝑍 −
ln 𝑍
ln 𝑍
𝜕ℎ- 𝜕ℎ4
𝜕𝐽-4
𝜕ℎ𝜕ℎ4=< 𝑠- 𝑠4 > − < 𝑠- >< 𝑠4 >Kappen&Rodriguez 1997, Takana 1998

• この右辺は、moment matching 条件を⽤いて観測値の統計量と⼀致させる• 左辺は、いわゆる感受率𝜕𝑚4 /𝜕ℎ- であり、self-consistency ⽅程式を微分すると求められる

𝐽$* =

J!*
"2K!#

− [𝐶

2"

]$*

𝐶 =< 𝑠- 𝑠4 >CDED−< 𝑠- >CDED< 𝑠4 >CDED
𝑚- =< 𝑠- >CDED
Copyright © 2022 by INIAD

Bethe近似を⽤いた学習l どんな近似⽅法でも、モデルで計算された相関(期待値)と、データの統計量を⼀致させる(moment matching)ことで学習が可能であるl Bethe近似(確率伝播法)を⽤いた学習法も、もちろん⾊々提案されていたl 2012年、 Bethe近似を⽤いたIsingモデルの学習について、⾮常に重要な論⽂が出る: Ricci-Tersenghi 2012

■ Bethe近似のもとで線形応答関係式とmoment matching条件を連⽴させたとき、その解として最適なパラメータ𝐽!" が閉形式で求まる
Copyright © 2022 by INIAD

Bethe近似によるIsingモデルの学習⽅程式
𝐶 =< 𝑠- 𝑠4 >CDED−< 𝑠- >CDED< 𝑠4 >CDED
𝑚- =< 𝑠- >CDEDl Jについての4次⽅程式であるが、閉形式で求まる事が分かった
l グラフ構造にもよるが、単純な平均場と⽐較して精度が良かった
l ⼀⽅、平⽅根の中が負になってしまう場合があるCopyright © 2022 by INIAD

Ricci-Tersenghi 2012

隠れ変数がある場合
l すべての変数の観測値が与えられておらず、観測されない
隠れ(潜在)変数がある場合は多い■ 可視変数集合を𝑿、隠れ変数集合を𝒁と書く

l 隠れ変数がある場合、それについて周辺化を⾏い、最尤推定を⾏う。尤度関数Lは
E

E

CD)

CD)

1
1
𝒅
𝐿 = ; ln 𝑃(𝑿 ) = ; ln ; 𝑃(𝑿 𝒅 , 𝒁)𝐷
𝐷𝒁

• 分配関数の時と同様に、確率変数Zについての和に計算コストがかかる• これについても、変分近似法による取り扱いがある
Copyright © 2022 by INIAD

変分ベイズ法l 尤度関数を次のように変形してみる
𝑃 𝑿 𝒅 , 𝒁 = 𝑃 𝒁|𝑿 𝒅 𝑃 𝑿 𝒅

ln 𝑃 𝑿 𝒅= ln 𝑃 𝑿 𝒅 , 𝒁 − ln 𝑃 𝒁|𝑿 𝒅

両辺に適当な確率分布𝑄(𝒁|𝑿(𝒅) )をかけ、𝒁について和を取ると、左辺はそのままなので、
ln 𝑃 𝑿 𝒅= 2 𝑄(𝒁|𝑿(𝒅) )[ln 𝑃 𝑿 𝒅 , 𝒁 − ln 𝑃 𝒁|𝑿 𝒅 ]
𝒁

𝒅 ,𝒁
𝒅
𝑃
𝑿
𝑃
𝒁|𝑿
= 2 𝑄(𝒁|𝑿(𝒅) )[ln
− ln
]
(𝒅)
(𝒅)
𝑄(𝒁|𝑿 )
𝑄(𝒁|𝑿 )
𝒁

= 2 𝑄(𝒁|𝑿 𝒅 )[ln𝒁

𝑃 𝑿 𝒅 ,𝒁
𝑄 𝒁𝑿𝒅

] + KL(𝑄|𝑃)

KL(𝑄|𝑃)= ∑𝒁 𝑄(𝒁|𝑿

𝒅

) ln

H(𝒁|𝑿(𝒅) )J 𝒁|𝑿 𝒅

KL情報量は⾮負であるから
ln 𝑃 𝑿 𝒅≥ 2 𝑄(𝒁|𝑿 𝒅 )[ln𝒁

𝑃 𝑿 𝒅 ,𝒁𝑄 𝒁𝑿𝒅]

が常に成⽴。等号は𝑄 𝒁 𝑿 𝒅Copyright © 2022 by INIAD= 𝑃 𝒁|𝑿 𝒅 のときのみ

ln 𝑃 𝑿 𝒅

≥ 2 𝑄(𝒁|𝑿 𝒅 )[ln𝒁

𝑃 𝑿 𝒅 ,𝒁
𝑄 𝒁𝑿𝒅

] = 2 𝑄 𝒁 𝑿 𝒅 ln 𝑃 𝑿 𝒅 , 𝒁 − 2 𝑄 𝒁 𝑿 𝒅 ln 𝑄 𝒁 𝑿 𝒅𝒁𝒁l 右辺をEvidence Lower BOund (ELBO)と呼ぶ
l ELBOはどんな𝑄についても、対数尤度の下界になっている
l 適当な𝑄を選べば、ELBOの和は評価できる場合がある(サンプリングを⾏うなど)
l そこで、𝑄を適当にパラメータ付けし、𝑃のパラメータと同時に、ELBOを最⼤化するように勾配法(1階微分を使う)で学習を⾏う
l このような⽅法を変分ベイズ法という
l 変分ベイズ法で学習を⾏うと、 𝑄 𝒁 𝑿 𝒅 によって可視変数から隠れ変数を推定することもできる
多くの場合、 𝑄に平均場近似を⽤いる
𝑄 𝒁𝑿𝒅= 3 𝑞- 𝑍𝒊 𝑿 𝒅

𝑍"

𝑍!

...
𝑄 𝒁𝑿

𝑃 𝑿𝒁

-

𝑋!

𝑋"

...

トピックモデルとして有名なLatent Dirichlet allocation はベイジアンネットワークを平均場近似を⽤いた変分ベイズ法で学習する
Blei+ 2003
Copyright © 2022 by INIAD

Stochastic Gradient Variational Bayesl 変分分布は可視変数の関数として複雑なものでも良い
l 複雑で⾮線形な関数近似器としてニューラルネットワークがある
l 確率的グラフィカルモデルPと変分分布Qのパラメータを
ニューラルネットワークの出⼒で表し、ELBOを最⼤化するようにニューラルネットワークを訓練する⽅法がある■ 変分オートエンコーダ(VAE)

Kingma & Weilling 2014l この学習アルゴリズムはStochastic Gradient Variational
Bayes と名付けられた
Copyright © 2022 by INIAD

変分オートエンコーダ
l あまり変分オートエンコーダが確率的グラフィカルモデルであると認識されないが、元論⽂には書いてある
l 連続な隠れ確率変数への対応や、サンプリングを⾏った場合でもbackpropagationを⾏えるような⼯夫など、⾮常に⾰新的な論⽂だった
l 変分オートエンコーダでも変分分布は平均場近似を⽤いている
Kingma & Weilling 2014

𝑍"

𝑍!
𝑋!

𝑋"

確率変数を少しずつ変えると、⽣成される可視変数も連続的に
変わる

...

𝑋*

可視変数は画素
Copyright © 2022 by INIAD

変分オートエンコーダの拡張
l 変分オートエンコーダが確率的グラフィカルモデルであることを意識すれば、その拡張は⾊々考えられる
l 例えば条件付き変分オートエンコーダは、隠れ変数ZとラベルYが協⼒して可視変数Xを⽣成する、というベイジアンネットワークとして理解できる
l Jointly Multimodal VAE(JMVAE)は、隠れ変数𝑍から異なるモダリティ(画像とテキスト)が⽣成されるベイジアンネットワークを⽤いている

𝒁𝑌
Suzuki+ 2017

𝑿
Kingma+ 2014
Copyright © 2022 by INIAD

平均場近似で⼗分なのか?
l ほとんどの変分ベイズ法による学習は、変分分布を平均場近似の形にしている
l 例えばBethe近似の形にしても良いはずである
l そのような研究もなくはないが、あまり流⾏っていない
l 平均場近似の簡単さ、収束性の良さのためと思われる。今後改善された⼿法が現れるだろうか?

Copyright © 2022 by INIAD

まとめ
l 確率的グラフィカルモデルは視覚的で分かりやすい⽣成モ
デルl しかし、推論にも学習にもコストがかかる
l 変分近似法として、平均場近似やBethe近似といった統計⼒学に現れた⼿法は強⼒であり、広く使われている
l ⼀⽅、このような近似の改善法は、存在はするがあまり応⽤されていない
l 今後より良い⽅法が⽣まれるか? 理論的な研究に期待
Copyright © 2022 by INIAD

参考⽂献(教科書)l Bishop, 「パターン認識と機械学習」
l ⻄森, 「相転移・臨界現象の統計物理学」
l Koller & Freidman, "Probabilistic Graphical Models"

Copyright © 2022 by INIAD