如何感性地理解EM算法?

作者: 工程师milter | 来源:发表于2017-03-30 19:15 被阅读11815次

    如果使用基于最大似然估计的模型,模型中存在隐变量,就要用EM算法做参数估计。个人认为,理解EM算法背后的idea,远比看懂它的数学推导重要。idea会让你有一个直观的感受,从而明白算法的合理性,数学推导只是将这种合理性用更加严谨的语言表达出来而已。打个比方,一个梨很甜,用数学的语言可以表述为糖分含量90%,但只有亲自咬一口,你才能真正感觉到这个梨有多甜,也才能真正理解数学上的90%的糖分究竟是怎么样的。如果EM是个梨,本文的目的就是带领大家咬一口。

    001、一个非常简单的例子

    假设现在有两枚硬币1和2,,随机抛掷后正面朝上概率分别为P1,P2。为了估计这两个概率,做实验,每次取一枚硬币,连掷5下,记录下结果,如下:

    硬币 结果 统计
    1 正正反正反 3正-2反
    2 反反正正反 2正-3反
    1 正反反反反 1正-4反
    2 正反反正正 3正-2反
    1 反正正反反 2正-3反

    可以很容易地估计出P1和P2,如下:

    P1 = (3+1+2)/ 15 = 0.4
    P2= (2+3)/10 = 0.5

    到这里,一切似乎很美好,下面我们加大难度。

    010、加入隐变量z

    还是上面的问题,现在我们抹去每轮投掷时使用的硬币标记,如下:

    硬币 结果 统计
    Unknown 正正反正反 3正-2反
    Unknown 反反正正反 2正-3反
    Unknown 正反反反反 1正-4反
    Unknown 正反反正正 3正-2反
    Unknown 反正正反反 2正-3反

    好了,现在我们的目标没变,还是估计P1和P2,要怎么做呢?

    显然,此时我们多了一个隐变量z,可以把它认为是一个5维的向量(z1,z2,z3,z4,z5),代表每次投掷时所使用的硬币,比如z1,就代表第一轮投掷时使用的硬币是1还是2。但是,这个变量z不知道,就无法去估计P1和P2,所以,我们必须先估计出z,然后才能进一步估计P1和P2。

    但要估计z,我们又得知道P1和P2,这样我们才能用最大似然概率法则去估计z,这不是鸡生蛋和蛋生鸡的问题吗,如何破?

    答案就是先随机初始化一个P1和P2,用它来估计z,然后基于z,还是按照最大似然概率法则去估计新的P1和P2,如果新的P1和P2和我们初始化的P1和P2一样,请问这说明了什么?(此处思考1分钟)

    这说明我们初始化的P1和P2是一个相当靠谱的估计!

    就是说,我们初始化的P1和P2,按照最大似然概率就可以估计出z,然后基于z,按照最大似然概率可以反过来估计出P1和P2,当与我们初始化的P1和P2一样时,说明是P1和P2很有可能就是真实的值。这里面包含了两个交互的最大似然估计。

    如果新估计出来的P1和P2和我们初始化的值差别很大,怎么办呢?就是继续用新的P1和P2迭代,直至收敛。

    这就是下面的EM初级版。

    011、EM初级版

    我们不妨这样,先随便给P1和P2赋一个值,比如:
    P1 = 0.2
    P2 = 0.7

    然后,我们看看第一轮抛掷最可能是哪个硬币。
    如果是硬币1,得出3正2反的概率为 0.2*0.2*0.2*0.8*0.8 = 0.00512
    如果是硬币2,得出3正2反的概率为0.7*0.7*0.7*0.3*0.3=0.03087
    然后依次求出其他4轮中的相应概率。做成表格如下:

    轮数 若是硬币1 若是硬币2
    1 0.00512 0.03087
    2 0.02048 0.01323
    3 0.08192 0.00567
    4 0.00512 0.03087
    5 0.02048 0.01323

    按照最大似然法则:
    第1轮中最有可能的是硬币2
    第2轮中最有可能的是硬币1
    第3轮中最有可能的是硬币1
    第4轮中最有可能的是硬币2
    第5轮中最有可能的是硬币1

    我们就把上面的值作为z的估计值。然后按照最大似然概率法则来估计新的P1和P2。

    P1 = (2+1+2)/15 = 0.33
    P2=(3+3)/10 = 0.6

    设想我们是全知的神,知道每轮抛掷时的硬币就是如本文第001部分标示的那样,那么,P1和P2的最大似然估计就是0.4和0.5(下文中将这两个值称为P1和P2的真实值)。那么对比下我们初始化的P1和P2和新估计出的P1和P2:

    初始化的P1 估计出的P1 真实的P1 初始化的P2 估计出的P2 真实的P2
    0.2 0.33 0.4 0.7 0.6 0.5

    看到没?我们估计的P1和P2相比于它们的初始值,更接近它们的真实值了!

    可以期待,我们继续按照上面的思路,用估计出的P1和P2再来估计z,再用z来估计新的P1和P2,反复迭代下去,就可以最终得到P1 = 0.4,P2=0.5,此时无论怎样迭代,P1和P2的值都会保持0.4和0.5不变,于是乎,我们就找到了P1和P2的最大似然估计。

    这里有两个问题:
    1、新估计出的P1和P2一定会更接近真实的P1和P2?
    答案是:没错,一定会更接近真实的P1和P2,数学可以证明,但这超出了本文的主题,请参阅其他书籍或文章。
    2、迭代一定会收敛到真实的P1和P2吗?
    答案是:不一定,取决于P1和P2的初始化值,上面我们之所以能收敛到P1和P2,是因为我们幸运地找到了好的初始化值。

    100、EM进阶版

    下面,我们思考下,上面的方法还有没有改进的余地?

    我们是用最大似然概率法则估计出的z值,然后再用z值按照最大似然概率法则估计新的P1和P2。也就是说,我们使用了一个最可能的z值,而不是所有可能的z值。

    如果考虑所有可能的z值,对每一个z值都估计出一个新的P1和P2,将每一个z值概率大小作为权重,将所有新的P1和P2分别加权相加,这样的P1和P2应该会更好一些。

    所有的z值有多少个呢?显然,有2^5=32种,需要我们进行32次估值??

    不需要,我们可以用期望来简化运算。

    轮数 若是硬币1 若是硬币2
    1 0.00512 0.03087
    2 0.02048 0.01323
    3 0.08192 0.00567
    4 0.00512 0.03087
    5 0.02048 0.01323

    利用上面这个表,我们可以算出每轮抛掷中使用硬币1或者使用硬币2的概率。比如第1轮,使用硬币1的概率是:
    0.00512/(0.00512+0.03087)=0.14
    使用硬币2的概率是1-0.14=0.86
    依次可以算出其他4轮的概率,如下:

    轮数 z_i=硬币1 z_i=硬币2
    1 0.14 0.86
    2 0.61 0.39
    3 0.94 0.06
    4 0.14 0.86
    5 0.61 0.39

    上表中的右两列表示期望值。看第一行,0.86表示,从期望的角度看,这轮抛掷使用硬币2的概率是0.86。相比于前面的方法,我们按照最大似然概率,直接将第1轮估计为用的硬币2,此时的我们更加谨慎,我们只说,有0.14的概率是硬币1,有0.86的概率是硬币2,不再是非此即彼。这样我们在估计P1或者P2时,就可以用上全部的数据,而不是部分的数据,显然这样会更好一些。

    这一步,我们实际上是估计出了z的概率分布,这步被称作E步。

    结合下表:

    硬币 结果 统计
    Unknown 正正反正反 3正-2反
    Unknown 反反正正反 2正-3反
    Unknown 正反反反反 1正-4反
    Unknown 正反反正正 3正-2反
    Unknown 反正正反反 2正-3反

    我们按照期望最大似然概率的法则来估计新的P1和P2:

    以P1估计为例,第1轮的3正2反相当于
    0.14*3=0.42正
    0.14*2=0.28反
    依次算出其他四轮,列表如下:

    轮数 正面 反面
    1 0.42 0.28
    2 1.22 1.83
    3 0.94 3.76
    4 0.42 0.28
    5 1.22 1.83
    总计 4.22 7.98

    P1=4.22/(4.22+7.98)=0.35

    可以看到,改变了z值的估计方法后,新估计出的P1要更加接近0.4。原因就是我们使用了所有抛掷的数据,而不是之前只使用了部分的数据。

    这步中,我们根据E步中求出的z的概率分布,依据最大似然概率法则去估计P1和P2,被称作M步。

    101、总结

    以上,我们用一个实际的小例子,来实际演示了EM算法背后的idea,共性存于个性之中,通过这个例子,我们可以对EM算法究竟在干什么有一个深刻感性的认识,掌握EM算法的思想精髓。

    Reference:
    http://pan.baidu.com/s/1i4NfvP7

    相关文章

      网友评论

      • 豆子1980:怒赞!
      • e35b0d2923fa:受益良多,但是请换P1=0.6,P2=0.3试试
        e35b0d2923fa:@909c00fa083f 如果是这样,那我觉得EM根本没有实用性?
        909c00fa083f:@剑人一箩筐 2、迭代一定会收敛到真实的P1和P2吗?
        答案是:不一定,取决于P1和P2的初始化值,上面我们之所以能收敛到P1和P2,是因为我们幸运地找到了好的初始化值。

        按照作者的说法,EM的结果很受初始值影响。
        e35b0d2923fa:采用上述初始化会导致相反结论,怎么解决?
      • 榴莲love臭豆腐:谢谢作者,写的通俗易懂,辛苦了。我有个很傻瓜式的问题,可能自己基础不好吧,就是为什么我们假设的初始概率利用最大似然估计计算出来哪个是第一个硬币,哪个是第二个硬币后,为何在利用这个信息反推两个硬币正面概率的值会和我们初始的P不一样呢,用一次极大似然估计就估计哪个硬币是哪个,我在用极大似然估计往回估计概率缺不是原来的p了呢?
        e35b0d2923fa:@剑人一箩筐 采用上述初始化会导致相反结论,怎么解决?
        e35b0d2923fa:受益良多,但是请换P1=0.6,P2=0.3试试
      • 草里有只羊:爱你,写得真棒!! 通俗易懂
        e35b0d2923fa:@剑人一箩筐 采用上述初始化会导致相反结论,怎么解决?
        e35b0d2923fa:受益良多,但是请换P1=0.6,P2=0.3试试
      • 30fbdecfc3ba:我这几天把鲍姆-韦尔奇算法用python重新实现了一次,我才真正了解了为什么P1和P2最终都收敛到了0.44,因为EM算法不能无限制地计算下去,到了一个人为限制的条件就需要停止,这时候我们就认为他收敛了到了一个相对较好的值了。否则无限制地计算下去最终的极值都是完全错误的。具体代码可以参考我的CSDN博客,地址是:
        https://blog.csdn.net/bian_h_f612701198412/article/details/81106807。我用了别人的一个初始值,然后自己写代码实现的。所以EM算法到了一定的限定条件以后,我们就需要停止计算了。
        kQYsnJc:什么鬼。。。。。em算法会收敛的。。。。只是有可能不是收敛到最优值罢了。。。。
        30fbdecfc3ba:@milter :smile: 以前我以为只要EM算法一直计算下去,就会收敛到一个正确的值,但是经过实践以后,原来需要我们设置停止计算的值,那我们就认为他收敛了,否则就到了一个完全错误的值,其实也不一定就完全错误,而是它不是我们希望看到的值。
        工程师milter:@边宏飞 赞!!!
      • 30fbdecfc3ba:根据楼主的数据,P1,P2最终都收敛到了0.44。所以我觉得这些数据一定要做置信,或者假设检验,才能确定楼主提供的数据是否可信。
        工程师milter:@边宏飞 em算法本来就不保证能逼近全局最优的。你可以换换初始化值
        30fbdecfc3ba:@milter 我把代码写出来以后,不管怎么逼近,都不会是真正的结果,不知道为什么?最后P1和P2都是0.44左右。
        工程师milter:@边宏飞 你说的假设检验具体怎么做呢?
      • 30fbdecfc3ba:文章写得真好,最害怕类似于李航的书,一上来就给我摆上一桌子公式,还有一瓶烧刀子。2分钟我就趴下来了。楼主在EM初级版,那里计算抛掷5次,2正3反的概率的时候我觉得应该用二项分布去计算概率吧?不应该0.2*0.2*0.8*0.8*0.8。
        工程师milter:@边宏飞 这里计算的是五次投掷情况,前2次为正,后三次为反。有顺序的,所以不用二项分布
      • 371f1e38c88b:为什么经典的书(比如李航统计学习方法)就不能说的像博主一样通俗易懂点呢:joy:
        工程师milter:@微积分_a31c 因为大牛不屑😉
      • 4e5897781711:看了一整天EM,这个感觉终于能看懂了
        工程师milter:@Deadlock_83fb 那这篇文章就算实现了它的价值了😉
      • DaydayHoliday:大大的赞。:+1: :grin: :grin: :grin:
      • 1a377bbe18b1:你好,我复制了@冯顿 的问题,刚好也对这个有疑问,能否给点简单的解说呢?
        有个疑问:
        “以P1估计为例,第1轮的3正2反相当于
        0.14*3=0.42正
        0.14*2=0.28反”
        这一部分,是后验概率P(w=1|x,Θ)=0.14,为什么要用0.14乘以第一轮中反面的出现的两次呢?
        王简之:P1: 硬币1抛正面的概率
        如何求呢?
        P1 = (用硬币1抛正面的次数)/(用硬币1抛的所有次数)

        概念的阻碍是次数是整数,但这里引入了概率,所以有了0.14
        原来:要么是硬币1,要么是硬币2, 如果第一轮确定用硬币1,则次数正面次数是3,反面是2
        现在:用硬币1的概率是0.14, 所以第一轮用硬币1抛正面的次数是0.14 * 3, 反面是0.14 * 2
        现在是对「次数」 乘了 「用硬币1的概率」

        但本质上,还是用下面的公式在计算:
        P1 = (用硬币1抛正面的次数)/(用硬币1抛的所有次数)
      • 25d5987aa6e2:有个疑问:
        “以P1估计为例,第1轮的3正2反相当于
        0.14*3=0.42正
        0.14*2=0.28反”
        这一部分,是后验概率P(w=1|x,Θ)=0.14,为什么要用0.14乘以第一轮中反面的出现的两次呢?
        1a377bbe18b1:@冯顿 你好,刚好这一块我也有疑问,可否再详细一点说说?谢谢
        25d5987aa6e2:这一块我弄明白了,这个结果其实是类似于最大似然估计的,Mstep需要最大化,求导以后就可以得
      • 云端漫步者2011:A good article. And thanks for the attached paper!
      • f895a9d03245:看懂了,大赞
        工程师milter:@Hydra_5572 那就拜托多帮帮推广下哈,让更多的人受益
        Hydra_5572:终于看到篇 能看懂的了 :flushed:
        工程师milter:@canopythonic 那就多多推广一下我的号吧,谢谢!我会持续写的
      • f895a9d03245:终于看明白了~感性认识是理解公式推导的前提。大赞大赞
      • easytoremember:用猜测的结果去与已经得到的结果去对比,感觉怪怪的,那EM算法的实质作用是用来干嘛的
        扇子和杯子:真的赞:+1:
        陶大明:EM算法实质就是在包含隐变量的条件下求各个参数
      • 9357f4d57301:写的真好!没有笨学生,只有不愿用心教的老师。
        工程师milter: @Firefly_T 期待一起交流,共同进步
      • 5c5d44761f80:楼主,在计算EM初级版的时候,使用你提供的P1,P2设定为何最终只能收敛到p1=0.33,P2=0.6之后将无法再逼近真实值了,是我计算错了吗? 还希望楼主能抽空简答一下.
        工程师milter:@tompusheng 这个只能根据具体的问题来决定,一般来说,当反复迭代逼近的距离很小时,就可以停下来。不知道你用不用xgboost,里面有一个参数early_stop_round就是类似的意思。需要说明的是,即使是升级版,也不能保证逼近全局最优解
        5c5d44761f80:@milter 楼主,对于升级版,我有个疑问就是,我们不断迭代p1,p2时,那什么时候来判断出它停下来呢?
        工程师milter:@tompusheng 初级版是有局限性的,可以对比升级版
      • 7dea46513000:谢谢分享,收益颇丰!
      • 小亮_85e5:很棒!!!
      • 51900d543390:0.14*3=0.42正
        0.14*2=0.28反
        这是什么意思
        928643364fb0:我也是这一块没看懂
        坑底Z蛙:@milter 你这个意思就是0.14代表用第一枚硬币投掷为正面朝上的概率?
        工程师milter: @睡觉说梦话XMFLSer 0.14表示的是第一轮投掷用的是硬币1的概率,3表示第一轮投掷有三次为正,那么,可以认为第一次投掷中用硬币1投掷得出正面的“次数”为0.14*3=0.42次
      • qiubo_ml:撇开公式推导这一块 EM算法概念解释得通俗易懂 非常好的一篇文章
        工程师milter: @Namely_20b9 感谢夸奖,还请多指点
      • 6935ae3905c3:终于搜到一篇能能看懂的文章了
        工程师milter: @睡在半山腰 哈哈,多多交流!
      • yzqlambda:011、EM初级版
        这一节的例子里面是不是第一轮过后就收敛了?
        工程师milter: @yzqlambda 我都是手动计算的,不知道再算一轮会如何
      • b19707134332:蛮不错的文章,没有数学公式。
        只是需要一开始就说明一下Z的意义是正反。
        工程师milter: @Phigo 硬币1和2的flag
        9f5dd76c1445:Z是正反还是硬币1和2的flag?
        工程师milter: @TensorFlow教室 感谢指点!
      • dc0b4fa52486:最近在学pgm:fearful::fearful:各种蒙蔽
        dc0b4fa52486:嗯嗯!:heart_eyes::heart_eyes:
        工程师milter: @无语_00d2 找点适合自己的入门资料
      • CoderJmd:怎么一点都看不懂,。
        工程师milter: @CoderJmd 可能是我写得不够好,文末的参考资料很不错,可以看看

      本文标题:如何感性地理解EM算法?

      本文链接:https://www.haomeiwen.com/subject/fcguottx.html