EMアルゴリズムとは †確率モデルのパラメータを最尤法に基づいて推定する手法のひとつであり、 観測不可能な潜在変数に確率モデルが依存する場合に用いられる。
例:混合正規分布でy=ax+uのモデルの場合の最尤法 †同一正規分布の場合観測値yをm次元ベクトル、入力(隠れ状態) xをn(n>m)次元ベクトルとする。データセット{xt, yt t= 1,..., T} から、パラメータa を求める問題を考える。ただしu は分平均0、分散σ^2Iのガウスノイズである。このとき、y は中心位置ax、共分散σ^2Iの正規分布N(y|ax,σ^2) にしたがう。I:単位行列とする。 P(y|x,a)=N(y|ax,σ^2)=SQRT(2πσ^2)^(-m/2)EXP[-||y-ax||^2/)2σ^2)] 中心位置が入力x に依存しているため上式で与えられるy の分布はx の条件付きである。対数尤度を取ると LogP(y|x,a)=(-1/2σ^2)||y-ax||^2 + aに無関係な定数 で表わされる。 観測データy1;t=(y1,y2,...yt) x1;t=(x1,x2,...,xt)の生起確率は、独立なガウス分布に従うことから、y1;tの正規確率は上記の確率密度関数ので表わされるので、対数尤度関数は和の形になり P(y1,y2,...,yt|a)=P(y1|a)P(y2|a)....P(yt|a) であることから、観測データ所与のもとでの尤度関数は L(a)=(-1/2σ^2) Σ||yt-axt||^2 + aに無関係な定数 このことから、L(a) をa について最大化することと、二乗誤差E=Σ||yt-axt||^2を最小化することは等価である。したがってこの場合の最尤推定法は最小二乗法と等価である。 混合分布とは †さて、ここで、M個の正規分布、P(X|i)=N(x|μi,σi^2),i=1,Mがあり、この混合分布を示す。 混合モデルとは、M個の分布の重ね合わせである。 P(X|θ)=ΣP(x,i|θ)=ΣP(X|i,θ)P(i|θ) i=1~Mの合計
対数尤度関数と最尤法 †データセットx1;t=(x1,x2,...,xt)が、観測された場合、 尤度は P(x1;t|θ)=P(xt,x1;t-1|θ)=P(xt,x1;t-1|θ) =P(xt|x1;t-1,θ)P(x1;t-1|θ) =P(xt|x1;t-1,θ)......P(x3|x1;2,θ)P(x2|x1;1,θ) =πP(xt|x1:t-i,θ) i=1~tの積 =π {πP(k|x1:t-i,θ)P(xt,k|x1:t-i,θ)} であるので、その対数尤度は L(θ)=LogP(x1;t|θ)=Σlog P(xt|x1;t-i,θ) = ΣLog {Σp(k|x1:t-i,θ)P(xt,k|x1:t-i,θ)} 内側のΣは隠れ変数k=1~Mについて、外側のΣはt=1~tについての合計である。 である。 これをを最大化するθを求めるのが最尤法である。
逐次表示でない場合、独立性を仮定してP(x1;t|θ)=P(Xt|θ)P(xt-1|θ)...P(x2|θ)p(x1|θ)として、L(θ)=ΣP(x1;t|θ)=Σlog Σp(xt,i|θ) から必要条件である停留条件はθで偏微分して0とおいて次のように計算できる。 [ΣΣ∂p(xt,i|θ)/∂θ]/Σp(xt,j|θ)=0 ΣΣ{p(i|Xt,θ)∂p(xt,i|θ)/∂θ}=0 a式 ただし、p(i|Xt,θ)=P(xt,i|θ)/Σp(xt,j|θ) b式 p(i|Xt,θ)は、データx1;t を観測した際の、隠れ変数i の事後確率(posterior) と呼ばれる。 それぞれのパラメータについて、最尤推定を行ってみる。Xベクトルの次数をdとする。
EMアルゴリズム †一般に、上記の方程式は、非線形であり、EMアルゴリズムでは、式の形に注目して、以下のような繰り返しアルゴリズムを考える。
前記の、パラメータθ=(γ1,...,γm;μ1,...,μmσ1,...,σm)の正規混合分布の場合を計算してみる。但し、Xtベクトルの次元をd次元とする。 P(x,i|θ)=vi・N(μi,σi^2)=vi・SQRT(2πσi^2)^(-d/2)EXP[-(x-μi)^2/2σ^2] P(x|θ)=Σp(x,i|θ)=ΣP(x|i,θ)p(i|θ)=Σvi・N(x|μi,σi^2) Vi=EXP(γi) である。
対数尤度を安定的に減少させるEMアルゴリズム †EとMのステップを1回繰り返すことで、対数尤度が増加する。最初の節で述べたように対数尤度関数は、マイナスの誤差二乗和に比例する形であるので、誤差二乗和が安定的に減少しているかチェックできる。そこで E=-L(θ)とおいて、パラメータの更新によって、これがどの程度減少するかを調べてみよう。 Et+1-Et=L(θt)-L(θt+1) =-ΣLogp(xt+1|θt+1)+ΣLogp(xt|θt) =-ΣLog{p(xt+1|θt+1)/p(xt|θt)} =-ΣLog{Σ[P(j|θt+1)p(xt+1|j,θt+1)P(j|xt,θt)]/[P(xt|θt)p(j|xt,θt)]} Jensenの不等式をここで使ってみる。 λj>o,Σλj=1 ---->log[Σλjxj] > Σλjxj より Et+1-Et<-ΣΣP(j|xt,θt)log{P(j|θt+1)p(xt+1|j,θt+1)/[P(xt|θt)p(j|xt,θt)]}=Q となるので、Qを最小化するところの、Mステップによって、Qの値に応じて、誤差二乗和が減少することがわかる。EMアルゴリズムは、安定的に誤差二乗和を減少させる最適化アルゴリズムになっていることが理解できる。 例:隠れマルコフモデル †隠れマルコフモデルは,観察データの列y1:T={y1,y2,...,yT} をモデル化するために用いる.なお,各時刻t で観察されるデータyt は,全部でM 種類あるシンボルw1,...,wM のうちのひとつである.つまり,モデル化したい観察データは,長さT のシンボル列である。 簡単のため、2種類のシンボル、(降雨あり、降雨なし)を考える。この場合、観測データーは、毎日の降雨の有無の観測結果である。そして、隠れ状態としては、(晴、曇り、雨)を考える。
初期の天気の隠れ状態の確率をP(s1)で表わし、St-1からStに推移する確率が、1次のマルコフ連鎖で表わされとしよう。この時、 Prob(St=ej|St-1=ei)=aij Σaij=1 i=1,nの式 Σはjに関する合計。 である。これらの確率をまとめてA = {aij} 行列であらわす。 また、隠れ状態がeiのとき、シンボルwjが表われる確率を、cijと書こう。これらの確率をまとめてC = {cij} 行列であらわす。 Prob{yt=wj|St=ei}=cij Σcij=1 j=1,2の式 Σはiに関する合計 そして、最初の隠れ状態がek となる確率P(s1 = ek) を,πk と書き、これを初期状態π0=(π1,π2,π3)と呼ぶ。 仮定より、長さT の完全データ(つまり,観察されたデータも観察されていないデータも含めたすべてのデータ)の尤度は,次のように書くことができる. P(s1:T , y1:T ) = P(s1)P(y1|s1)......P(yT|sT) 以上すべては,隠れマルコフモデルのパラメータである.つまり,θ = (A; C; π0) となる。 パラメータの数は、K^2 + KM + K 個(k=3,M=2なので合計)のパラメータの値を定める必要がある 参考 †
|