简单易学的机器学习算法——EM算法

一、机器学习中的參数预计问题

    在前面的博文中，如“简单易学的机器学习算法——Logistic回归”中，採用了极大似然函数对其模型中的參数进行预计，简单来讲即对于一系列样本，Logistic回归问题属于监督型学习问题，样本中含有训练的特征以及标签。在Logistic回归的參数求解中。通过构造样本属于类别和类别的概率：

这样便能得到Logistic回归的属于不同类别的概率函数：

此时，使用极大似然预计便可以预计出模型中的參数。

可是。假设此时的标签是未知的。称为隐变量，如无监督的学习问题，典型的如K-Means聚类算法，此时不能直接通过极大似然预计预计出模型中的參数。

二、EM算法简单介绍

    在上述存在隐变量的问题中，不能直接通过极大似然预计求出模型中的參数，EM算法是一种解决存在隐含变量优化问题的有效方法。

EM算法是期望极大(Expectation Maximization)算法的简称，EM算法是一种迭代型的算法。在每一次的迭代过程中。主要分为两步：即求期望(Expectation)步骤和最大化(Maximization)步骤。

三、EM算法推导的准备

1、凸函数

    设是定义在实数域上的函数，假设对于随意的实数。都有

那么是凸函数。若不是单个实数，而是由实数组成的向量，此时。假设函数的Hesse矩阵是半正定的，即

那么是凸函数。特别地。假设或者。那么称为严格凸函数。

2、Jensen不等式

    假设函数是凸函数，是随机变量，那么

特别地，假设函数是严格凸函数，那么当且仅当

即随机变量是常量。

(图片来自參考文章1)

注：若函数是凹函数。上述的符号相反。

3、数学期望

3.1随机变量的期望

   设离散型随机变量的概率分布为：

当中。，假设绝对收敛，则称为的数学期望，记为，即：

   若连续型随机变量的概率密度函数为。则数学期望为：

3.2随机变量函数的数学期望

   设是随机变量的函数。即，若是离散型随机变量，概率分布为：

则：

   若是连续型随机变量，概率密度函数为。则

四、EM算法的求解过程

    如果表示观測变量，表示潜变量，则此时即为全然数据，的似然函数为，当中，为须要预计的參数，那么对于全然数据，的似然函数为。

    构建好似然函数，对于给定的观測数据，为了预计參数，我们能够使用极大似然预计的方法对其进行预计。由于变量是未知的。我们仅仅能对的似然函数为进行极大似然预计，即须要极大化：

上述式子中无法直接对求极大值，由于在函数中存在隐变量，即未知变量。若此时，我们可以确定隐变量的值，便可以求出的极大值，可以用过不断的改动隐变量的值，得到新的的极大值。这便是EM算法的思路。通过迭代的方式求出參数。

    首先我们须要对參数赋初值，进行迭代运算，如果第次迭代后參数的值为，此时的log似然函数为，即：

在上式中，第二行到第三行使用到了Jensen不等式，因为log函数是凹函数，由Jensen不等式得到：

而

表示的是的期望。当中，表示的是隐变量满足的某种分布。这样，上式的值取决于和的概率。

在迭代的过程中。调整这两个概率，使得下界不断的上升，这样就能求得的极大值。

注意，当等式成立时。说明此时已经等价于。由Jensen不等式可知，等式成立的条件是随机变量是常数，即：

已知：

所以：

则：

至此，我们得出了隐变量满足的分布的形式。这就是EM算法中的E步。

在确定了后，调整參数使得取得极大。这便是M步。EM算法的步骤为：

1. 初始化參数。開始迭代；
2. E步：如果为第次迭代參数的预计值，则在第次迭代中。计算：
3. M步：求使极大化的，确定第次的參数的预计值：

五、EM算法的收敛性保证

迭代的过程是否能保证最后找到的就是最大的似然函数值呢？即须要证明在整个迭代的过程中，极大似然预计是单调添加的。假定和是EM算法的第次和第次迭代后的结果，选定，进行迭代：

1. E步：
2. M步：

固定。将看成变量：

上式中，第一个大于等于是由于：

六、利用EM算法參数求解实例

    如果有有一批数据各自是由两个正态分布：

产生，当中，和未知，。可是不知道详细的是第产生，即能够使用和表示。

这是一个典型的涉及到隐藏变量的样例，隐藏变量为和。

能够使用EM算法对參数进行预计。

1. 首先是初始化和；
2. E步：，即求数据是由第个分布产生的概率：
3. M步：，即计算最大的期望值。

   然而我们要求的參数是均值，能够通过例如以下的方式预计：

Python代码

#coding:UTF-8
'''
Created on 2015年6月7日

@author: zhaozhiyong
'''
from __future__ import division
from numpy import *
import math as mt
#首先生成一些用于測试的样本
#指定两个高斯分布的參数，这两个高斯分布的方差同样
sigma = 6
miu_1 = 40
miu_2 = 20

#随机均匀选择两个高斯分布，用于生成样本值
N = 1000
X = zeros((1, N))
for i in xrange(N):
    if random.random() > 0.5:#使用的是numpy模块中的random
        X[0, i] = random.randn() * sigma + miu_1
    else:
        X[0, i] = random.randn() * sigma + miu_2

#上述步骤已经生成样本
#对生成的样本，使用EM算法计算其均值miu

#取miu的初始值
k = 2
miu = random.random((1, k))
#miu = mat([40.0, 20.0])
Expectations = zeros((N, k))

for step in xrange(1000):#设置迭代次数
    #步骤1。计算期望
    for i in xrange(N):
        #计算分母
        denominator = 0
        for j in xrange(k):
            denominator = denominator + mt.exp(-1 / (2 * sigma ** 2) * (X[0, i] - miu[0, j]) ** 2)
        
        #计算分子
        for j in xrange(k):
            numerator = mt.exp(-1 / (2 * sigma ** 2) * (X[0, i] - miu[0, j]) ** 2)
            Expectations[i, j] = numerator / denominator
    
    #步骤2。求期望的最大
    #oldMiu = miu
    oldMiu = zeros((1, k))
    for j in xrange(k):
        oldMiu[0, j] = miu[0, j]
        numerator = 0
        denominator = 0
        for i in xrange(N):
            numerator = numerator + Expectations[i, j] * X[0, i]
            denominator = denominator + Expectations[i, j]
        miu[0, j] = numerator / denominator
        
    
    #推断是否满足要求
    epsilon = 0.0001
    if sum(abs(miu - oldMiu)) < epsilon:
        break
    
    print step
    print miu
    
print miu

终于结果

[[ 40.49487592  19.96497512]]

參考文章：

1、(EM算法)The EM Algorithm (http://www.cnblogs.com/jerrylead/archive/2011/04/06/2006936.html)

2、数学期望(http://wenku.baidu.com/view/915a9c1ec5da50e2524d7f08.html?re=view)