下面这种分布,它并不属于任何一种标准分布:
image.Png
它看起来像多个高斯模型,能用多个高斯模型来表达这个分布呢?这样做有一个问题就是,我们不知道哪个数据属于第一个高斯模型,哪个数据属于第二个高斯模型。

这个问题可以从产生数据的角度思考,假设现在有这样K个高斯分布:
image.Png
现在要根据这K个高斯分布产生数据:

  • 依照概率a_k随机选择一个高斯分布\phi
  • 依照这个高斯分布随机产生数据

重复这个操作N次,于是产生了数据y_ { j } , \quad j = 1,2 , \cdots ,  N.

用数学公式描述上述的事情:

P ( y ,  z | \theta ) =\prod _ { j = 1 } ^ { N }\prod _ { k = 1 } ^ { K } [a_k\phi ( y_i | \theta _ { k } )]^{ z_{jk}}

其中z为隐变量,z_{jk} =1代表第j个数据选择了第k个模型.

这种逆向的思想非常重要,我们直觉上习惯根据数据去区分是属于哪个分布的,但是这很难做到.先假设模型存在,然后根据这些模型产生的效果去确定模型的参数,这就是极大似然法的思想.

因为:

P (  z ) = \prod _ { k = 1 } ^ { K } P \left(  z _ { k } = 1 \right) ^ {  z _ { k } } = \prod _ { k = 1 } ^ { K } \alpha _ { k } ^ {  z _ { k } }

P ( y | z ) = \prod _ { k = 1 } ^ { K } P ( y |  z _ { k } = 1 ) ^ {  z _ { k } } = \prod _ { k = 1 } ^ { K }\phi ( y | \theta _ { k } )^ {  z _ { k } }

所以有:

P ( y | \theta ) =\sum _ { k = 1 } ^ { K } P ( y |  z _ { k } = 1 ) P \left(  z _ { k } = 1 \right)=\sum _ { k = 1 } ^ { K } \alpha _ { k } \phi ( y | \theta _ { k } )

这个模型称为高斯混合模型. 其中\theta = \left\{\alpha _ { k },\mu _ { k } , \Sigma _ { k } \right\} _ { k = 1 } ^ { K }

使用EM算法对高斯混合模型进行参数估计

高斯混合模型的对数似然函数:

\log P(y|\theta)= \sum _ { i = 1 } ^ { N } \log P \left( y _ { i } | \theta \right)=\sum _ { i = 1 } ^ { N } \log\sum _ { k = 1 } ^ { K } \alpha _ { k } P ( y_i |\theta_ {k})

这个似然函数很复杂,无法直接求得其极值,事实证明使用梯度下降来求解模型参数的效果也并不好.对于高斯混合模型通常使用EM算法求解.

对于EM算法,E步计算Q函数:

Q(\theta,\theta^{(i)}) = E_Z[\log P(Y,Z|\theta)|Y,\theta^{(i)}]

M步极大化Q函数求得模型参数.

首先计算完全数据的似然函数:

\begin{aligned} P ( y ,  z | \theta ) & = \prod _ { j = 1 } ^ { N }\prod _ { k = 1 } ^ { K } [a_k\phi ( y_i | \theta _ { k } )]^{ z_{jk}} \\ 
& = \prod _ { k = 1 } ^ { K } \alpha _ { k } ^ { n _ { k } } \prod _ { j = 1 } ^ { N } \left[ \phi \left( y _ { j } | \theta _ { k } \right) \right] ^ {  z _ { jk } }\\
&= \prod _ { k = 1 } ^ { K } \alpha _ { k } ^ { n _ { k } } \prod _ { j = 1 } ^ { N } \left[ \frac { 1 } { \sqrt { 2 \pi } \sigma _ { k } } \exp \left( - \frac { \left( y _ { j } - \mu _ { k } \right) ^ { 2 } } { 2 \sigma _ { k } ^ { 2 } } \right) \right] ^ {  z _ { jk } }\end{aligned}

其中n _ { k } = \sum _ { j = 1 } ^ { N }  z _ { j k },\sum _ { k = 1 } ^ { K } n _ { k } = N

完全数据的对数似然函数为:

\log P ( y ,  z | \theta ) = \sum _ { k = 1 } ^ { K } n _ { k } \log \alpha _ { k } + \sum _ { j = 1 } ^ { N }  z _ { j k } \left[ \log \left( \frac { 1 } { \sqrt { 2 \pi } } \right) - \log \sigma _ { k } - \frac { 1 } { 2 \sigma _ { k } ^ { 2 } } \left( y _ { j } - \mu _ { k } \right) ^ { 2 } \right]

EM算法的E步:确定Q函数

\begin{aligned} Q \left( \theta , \theta ^ { ( i ) } \right) & = E [ \log P ( y ,  z | \theta ) | y , \theta ^ { ( i ) } ] \\ & = E \left\{ \sum _ { k = 1 } ^ { K } n _ { k } \log \alpha _ { k } + \sum _ { j = 1 } ^ { N }  z _ { j k } \left[ \log \left( \frac { 1 } { \sqrt { 2 \pi } } \right) - \log \sigma _ { k } - \frac { 1 } { 2 \sigma _ { k } ^ { 2 } } \left( y _ { j } - \mu _ { k } \right) ^ { 2 } \right] \right\} \\ & = \sum _ { k = 1 } ^ { K } \left\{ \sum _ { j = 1 } ^ { N } \left( E  z _ { j k } \right) \log \alpha _ { k } + \sum _ { j = 1 } ^ { N } \left( E  z _ { j k } \right) \left[ \log \left( \frac { 1 } { \sqrt { 2 \pi } } \right) - \log \sigma _ { k } - \frac { 1 } { 2 \sigma _ { k } ^ { 2 } } \left( y _ { j } - \mu _ { k } \right) ^ { 2 } \right] \right\} \end{aligned}

这里对E展开的时候将z _ { j k }替换为E z _ { j k },因为整个式子展开相当于K\times N\beta_ { j k }z _ { j k }的连加,所以可以这样替换,但是若展开后不是连加式则不能这样替换.

这里需要计算\mathrm { E } \left(  z _ { j \mathrm { k } } | \mathrm { y } , \theta \right),记为\hat {  z } _ { j k }

\begin{aligned} \hat {  z } _ { j k } & = E \left(  z _ { j k } | y , \theta \right) = P \left(  z _ { j k } = 1 | y , \theta \right) \\ & = \frac { P \left(  z _ { j k } = 1 , y _ { j } | \theta \right) } { \sum _ { k = 1 } ^ { K } P \left(  z _ { j k } = 1 , y _ { j } | \theta \right) } \\ & = \frac { P \left( y _ { j } |  z _ { j k } = 1 , \theta \right) P \left(  z _ { j k } = 1 | \theta \right) } { \sum _ { k = 1 } ^ { K } P \left( y _ { j } |  z _ { j k } = 1 , \theta \right) P \left(  z _ { j k } = 1 | \theta \right) } \\ & = \frac { \alpha _ { k } \phi \left( y _ { j } | \theta _ { k } \right) } { \sum _ { k = 1 } ^ { K } \alpha _ { k } \phi \left( y _ { j } | \theta _ { k } \right) } , \quad j = 1,2 , \cdots , N ; \quad k = 1,2 , \cdots , K \end{aligned}

\hat {  z } _ { j k }是在当前模型参数下第j个观测数据来自第k个分模型的概率,称为分模线k对观测数据y_j的响应度.

\hat {  z } _ { j k } = E  z _ { j k }n _ { k } = \sum _ { j = 1 } ^ { N } E  z _ { j k }带入得:

Q \left( \theta , \theta ^ { ( i ) } \right) = \sum _ { k = 1 } ^ { K } n _ { k } \log \alpha _ { k } + \sum _ { k = 1 } ^ { N } \hat {  z } _ { j k } \left[ \log \left( \frac { 1 } { \sqrt { 2 \pi } } \right) - \log \sigma _ { k } - \frac { 1 } { 2 \sigma _ { k } ^ { 2 } } \left( y _ { j } - \mu _ { k } \right) ^ { 2 } \right]

EM算法的M步

M步需要极大化Q函数求出\theta:

\theta ^ { ( i + 1 ) } = \arg \max _ { \theta } Q \left( \theta , \theta ^ { ( i ) } \right)

分别对\boldsymbol { u } _ { \mathrm { k } } , \sigma _ { k } ^ { 2 },\hat { \alpha } _ { k }求偏导,即可求得:

\hat { \mu } _ { k } = \frac { \sum _ { i = 1 } ^ { N } \hat {  z } _ { \mu } y _ { j } } { \sum _ { j = 1 } ^ { N } \hat {  z } _ { k } } , \quad k = 1,2 , \cdots , K

\hat { \sigma } _ { k } ^ { 2 } = \frac { \sum _ { j = 1 } ^ { N } \hat {  z } _ { j k } \left( y _ { j } - \mu _ { k } \right) ^ { 2 } } { \sum _ { j = 1 } ^ { N } \hat {  z } _ { j k } } , \quad k = 1,2 , \cdots , K

\hat { \alpha } _ { k } = \frac { n _ { k } } { N } = \frac { \sum _ { j = 1 } ^ { N } \hat {  z } _ { j k } } { N } , \quad k = 1,2 , \cdots , K


编程实现:
https://www.kaggle.com/swimmingwhale/gaussian-mixture-model


参考:
李航<统计学习方法>

posted @ 2018-08-16 18:26:58
评论加载中...

发表评论