機器學習 高斯混合模型引數估計的EM演算法

2021-07-24 16:07:52 字數 3566 閱讀 8491

看理論之前先來【舉個例子】

對於乙個未知引數的模型,我們觀測他的輸出,得到下圖這樣的直方圖:

我們先假設它是由兩個高斯分布混合疊加而成的,那麼我們該怎麼去得到這兩個高斯分布的引數呢?

em演算法!!

假設觀測資料 y1

,y2,

...,

yn是由高斯混合模型生成的。 p

(y|θ

)=∑k

=1kα

kθ(y

|θk)

其中,θ

= 。表示的是高斯模型的引數,em演算法也正是要用來估計高斯混合模型的這個引數。還是以上面的例子來說,對於我們的觀測資料 yi

,i=1

,2,.

..,n

來說,該資料肯定是由分模型的資料疊加得到的。那麼我們設想 yi

是這樣產生的:

1> 首先依概率 αk

選擇第

k 個高斯模型 ϕ(

y|θk

) ;

2> 然後依第

k 個分模型的概率分布 ϕ(

y|θk

) 生成觀測資料。

這時候觀測資料是已知的,反應觀測資料yj

來自第k 個分模型的資料是未知的,k=

1,2,

...,

k , 以隱變數 γj

k 來表示。

可以得到完全似然函式:

jk和 ∑n

j=1e

γjk 替換,得到q函式。eγ

jk表示分模型

k 對觀測資料yj

的響應度。

迭代m步就是求函式q(

θ,θ(

i)) 對

θ 的極大值,即求新一輪迭代的模型引數: θ

(i+1

)=ar

gmax

θq(θ

,θ(i

))

每一次迭代中引數計算公式表示可得到:

例項**:

import math

import copy

import numpy as np

import matplotlib.pyplot as plt

isdebug = false

# 指定k個高斯分布引數,這裡k=2。2個高斯分布具有相同均方差sigma,均值分別為mu1,mu2。

defini_data

(sigma,mu1,mu2,k,n):

global x

global mu

global expectations

x = np.zeros((1,n))

mu = np.random.random(2)

expectations = np.zeros((n,k))

for i in xrange(0,n):

if np.random.random(1) > 0.5:

x[0,i] = np.random.normal()*sigma + mu1

else:

x[0,i] = np.random.normal()*sigma + mu2

if isdebug:

print

"***********"

print

u"初始觀測資料x:"

print x

# em演算法:步驟1,計算e[zij]

defe_step

(sigma,k,n):

global expectations

global mu

global x

for i in xrange(0,n):

denom = 0

for j in xrange(0,k):

denom += math.exp((-1/(2*(float(sigma**2))))*(float(x[0,i]-mu[j]))**2)

for j in xrange(0,k):

numer = math.exp((-1/(2*(float(sigma**2))))*(float(x[0,i]-mu[j]))**2)

expectations[i,j] = numer / denom

if isdebug:

print

"***********"

print

u"隱藏變數e(z):"

print expectations

# em演算法:步驟2,求最大化e[zij]的引數mu

defm_step

(k,n):

global expectations

global x

for j in xrange(0,k):

numer = 0

denom = 0

for i in xrange(0,n):

numer += expectations[i,j]*x[0,i]

denom +=expectations[i,j]

mu[j] = numer / denom

# 演算法迭代iter_num次,或達到精度epsilon停止迭代

defrun

(sigma,mu1,mu2,k,n,iter_num,epsilon):

ini_data(sigma,mu1,mu2,k,n)

print

u"初始:", mu

for i in range(iter_num):

old_mu = copy.deepcopy(mu)

e_step(sigma,k,n)

m_step(k,n)

print i,mu

if sum(abs(mu-old_mu)) < epsilon:

break

if __name__ == '__main__':

run(6,40,20,2,1000,1000,0.0001)

plt.hist(x[0,:],50)

plt.show()

引數估計 CIR模型的引數估計

cox ingersoll ross cir 模型是量化金融風控中,特別是在利率和信用風險的期限結構模型中經常用到的一種模型。與其他模型如ho lee,vasicek等相比,它的特點是其解總是非負的 如果滿足feller條件則以概率為1為正 並且滿足均值回歸性質。cir 的基本形式是如下的sde 其...

機器學習 單高斯分布引數估計

對於單維高斯分布而言,其概率密度函式可以表示成 p x frac sigma e 其中 u 表示均值,sigma 2 表示方差。對於多維高斯分布而言,其概率密度函式可以表示成 p x frac lvert sigma rvert e x u t sigma x u 其中p表示維度,首先介紹如何根據極...

機器學習之引數估計

引數估計 parameter estimate 就是通過一系列演算法,來求出模型的最優引數。在各個機器學習深度學習的框架裡,都變成了optimizer的活了。其實這個名字很奇怪,但是在比較早的機器學習 裡都是這麼叫的,我們重點來關注下裡面涉及的一些演算法。這裡主要關注的是 二乘是平方的意思,感覺最小...