機械学習分科会、第八章の四,グラフィカルモデルによる推論

要約

木DPをすればよいです。

以下では、各変数(ノード)はK種類の離散状態をとるものとし、ノード(変数)はN個あるものとします。
あと、以下のコードは、pythonっぽい疑似コードもどきです。

ある変数xが状態aを取る確率p(x=a)を求める

モデルが一本鎖のとき

#p(xk=a)の値を求める
#p(x1,x2, ... ,xN) = phi[1,2](x1,x2) * phi[2,3](x2,x3) ... phi[N-1,N](xN-1,xN)
ans = 0
for x1 in xrange(0,K):
  ....
    for xk-1 in xrange(0,K):
      for xk+1 in xrange(0,K):
        ....
          for xN in xrange(0,K):
            ans += phi[1,2](x1,x2) * phi[2,3](x2,x3) ... phi[N-1,N](xN-1,xN)

ですが、なんやかんやすると、

arpha = 0
for x1 in xrange(0,K):
  ....
    for xk-1 in xrange(0,K):
      arpha += phi[1,2](x1,x2) * ... * phi[k-1,k](xk-1,xk=a)

beta = 0
for xk+1 in xrange(0,K):
  ....
    for xN in xrange(0,K):
      beta = phi[k,k+1](xk=a,xk+1) * ... * phi[N-1,N](xN-1,xN)

ans += arpha * beta

となって、さらにdpして、

for x0 in xrange(0,K):
	arpha[0][x0] = 1
for i in xrange(1,k+1):
	for xi in xrange(0,K):
		arpha[i][xi] = 0
		for xi-1 in xrange(0,K):
			arpha[i][xi] += arpha[i-1][xi-1] * phi[i-1,i](xi-1,xi)

for xN+1 in xrange(0,K):
	beta[N+1][xN+1] = 1
for i in xrange(k,N).reverse:
	for xi in xrange(0,K):
		beta[i][xi] = 0
		for xi+1 in xrange(0,K):
			beta[i][xi] += beta[i+1][xi+1] * phi[i,i+1](xi,xi+1)

ans = arpha[k][xk=a] * beta[k][xk=a]

となる。

モデルが木のとき

黒頂点(■)と白頂点(○)についてのdpみたいにすればいけます。
確率を求めたいノードを根にしてやって木を吊るします。
で、各白頂点は、子の黒頂点を、各黒頂点は、子の白頂点を持ってるものとします。
このとき、疑似コードはこんな感じ。

def root(x=a):
	res = 1
	for to in rootnode:
		 res *= kuro(no=to,state=a)
	return res

def kuro(no,state):
	#no番目の黒ノードの親白ノードが状態stateを取るときの、no番目の黒ノードを根とする部分木の確率の総和
	res = 0
	tos = kuro[no]
	ls = len(tos)
	for t1 in xrange(0,K):
		... 
			for tls in xrange(0,K):
				res = kuro_function[no](state,t1,t2, ... tls) * siro(no=tos[1],state=t1) *  ... * siro(no=tos[ls],state=tls)

def siro(state,no):
	#no番目の白ノードが状態stateを取るときの、no番目の白ノードを根とする部分木の確率の総和
	res = 1
	for to in rootnode:
		 res *= kuro(no=to,state=state)
	return res

def leafsiro(state,no):
	return 1


def leafkuro(state,no):
	return kuro_function[no]()

で、計算量ですが、メモ化なりなんなりすると、状態はK*Nパターンしか起こらないので、あとは遷移に掛かる時間に寄りますね。
具体的には、たくさん(dこ)子を持つ黒ノードがあった時に、K^dかかってしまうわけです。
で、因子グラフが木構造をもとにしてできていれば、dがたかだか2で抑えられて幸せという寸法。
(逆に、全結合みたいなやつだと、d=NとなってK^Nかかるのでえらいことになる)

最大の確率をとる変数ベクトル(x1,...,xN)とその時の確率p(x1,...,xN)を求める(max-sumというやつ)

これは、dp->経路復元みたいなのをしてやればよい。

モデルが一本鎖のとき

問題は、
『最大の確率 p(x1,...,xN) を求めよ』となり、
『最大の確率 phi[1,2](x1,x2) * ... * phi[N-1,N](xN-1,xN) を求めよ』となり、logは単調増加なので、logとってもよく、
『最大の log(phi[1,2](x1,x2)) + ... + log(phi[N-1,N](xN-1,xN)) を求めよ』となる。
で、各log(phi[a,a+1](b,c)) を、【頂点 a[b] から頂点 a+1[c] 間の辺の重み】みたいにみると、
『頂点1から頂点Nまでの経路のうち、重み最大のものの重みと、それを与える経路を求めよ』
みたいな感じになります。
で、これはもうdpしてからの経路復元で解けますね。

モデルが木のとき

同様のアナロジーをすると、
『木の各頂点に適切な状態を割り振ることにより、木全体の辺の重みを最大化しろ』
となります。
で、これは、
kuro[x][s] .. 頂点xが状態sを取るときのxからの部分木の辺の重み和の最大値
みたいなdpをすれば、解けますね。