CS229/EM.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# EM算法\n",
    "\n",
    "## 1.Jensen不等式\n",
    "\n",
    "回顾凸函数的定义，假设$f(X)$是定义在实数集上的函数，如果$f''(x)\\geq 0$对于所有的$x\\in\\mathbb{R}$都成立，则称$f(X)$是凸函数；如果$f''(x)> 0$对于所有的$x\\in\\mathbb{R}$都成立，则称$f(X)$是严格凸函数。当$f(X)$定义在实数向量空间上时，如果其海赛矩阵$H(x)$是半正定的，则称$f(X)$是凸函数；如果$H(x)$是正定矩阵，则称$f(X)$是严格凸函数。\n",
    "接着我们给出Jensen不等式的定义：  \n",
    "如果$f(X)$是凸函数，那么对于任意的$X$都有下式成立：   \n",
    "$$E[f(X)]\\geq f(E[X])$$\n",
    "\n",
    "下图给出了这个定理的直观解释： \n",
    "\n",
    "<img src=http://7xikew.com1.z0.glb.clouddn.com/convex_function.png width=600px>\n",
    "\n",
    "图中，点$X$是一个随机变量，以0.5的概率取值点a，0.5的概率取值点b，那么$E[X]$就落在$[a,b]$的中点。而$E[f(X)]$则位于$y$轴上$f(a)$和$f(b)$的中点。观察可以发现，对于凸函数来说，$E[f(X)]$的值一定在$f(E[X])$的上方。\n",
    "\n",
    "这个定理可以总结为：“对于凸函数而言，函数值的期望总是大于等于期望的函数值”。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.EM算法\n",
    "\n",
    "假设我们有$m$个独立样例构成的训练集$\\{x^{(1)},...,x^{(m)}\\}$，我们希望从数据中学习到模型$p(x)$的参数$\\Theta$，一般都是通过最大化如下的似然函数来找到参数$\\Theta$的最大似然估计：\n",
    "\n",
    "$$\\ell(\\Theta)=\\sum_{i=1}^m \\log p(x^{(i)};\\Theta)$$\n",
    "\n",
    "然而并不是所有情况下都有解析解的。EM算法给出的解决办法是引入隐变量$z$，具体来说就是估计$p(x,z)$的参数，其似然函数为\n",
    "\n",
    "$$\\ell(\\Theta)=\\sum_{i=1}^m \\log p(x^{(i)};\\Theta)=\\sum_{i=1}^m \\log \\sum\\limits_{z^{(i)}} p(x^{(i)},z^{(i)};\\Theta)$$\n",
    "\n",
    "此处，$z$是离散随机变量，如果$z$是连续随机变量，求和符号变更为积分符号。\n",
    "\n",
    "EM算法的原理是：  \n",
    "1）E-step：给出$\\ell(\\Theta)$的下界  \n",
    "2）M-step：优化这个下界  \n",
    "\n",
    "下面我们来详细叙述每个步骤的推导过程：  \n",
    "### E-step\n",
    "对于每个$i$，设$Q_i$是$z$上的分布（$\\sum\\limits_z Q_i(z)=1,Q_i(z)\\geq 0$)，那么  \n",
    "$$\\begin{aligned}\\sum_{i=1}^m \\log p(x^{(i)};\\Theta)&=\\sum_{i=1}^m \\log\\sum\\limits_{z^{(i)}} p(x^{(i)},z^{(i)};\\Theta)\\\\\n",
    "&=\\sum_{i=1}^m \\log\\sum\\limits_{z^{(i)}} Q_i(z^{(i)}) \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})}\\\\&=\\sum_{i=1}^m \\log \\mathbb{E}_{z^{(i)}\\sim Q_i}\\bigg[\\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})}\\bigg]\\\\&\\geq \\sum_{i=1}^m  \\mathbb{E}_{z^{(i)}\\sim Q_i}\\bigg[ \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})}\\bigg]\\\\&=\\sum_{i=1}^m \\sum_{z^{(i)}}Q_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})} \\end{aligned}$$\n",
    "\n",
    "最后一步我们用到了Jensen不等式，同时注意到$\\log x$是凹函数，因此不等式方向与凸函数版本的Jensen不等式相反。\n",
    "\n",
    "由此，我们得到了$\\ell(\\Theta)$的下界，无论我们选择什么样的$Q(z)$，上述不等式都成立。现在我们的问题是，对于特定的$\\theta$选择什么样的$Q$能够使得不等式两端的差距最小，即不等式的等号成立。实际上，只要$\\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})}$是一个与$z^{(i)}$无关的常量，那么在求期望时就可以把$z^{(i)}$消除，因此我们假设：\n",
    "$$\\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})}=c(x^{(i)};\\Theta)$$\n",
    "当我们只考虑$z$时，$c(x^{(i)};\\Theta)$可以视为一个常量。结合$\\sum_{z^{(i)}} Q_i(z^{(i)})=1$的事实，我们知道\n",
    "$$c(x^{(i)};\\Theta)=\\sum_{z^{(i)}} p(x^{(i)},z^{(i)};\\Theta)=p(x^{(i)};\\Theta)$$\n",
    "于是\n",
    "$$Q_i(z^{(i)})=\\frac{p(x^{(i)},z^{(i)};\\Theta)}{p(x^{(i)};\\Theta)}=p(z^{(i)}|x^{(i)};\\Theta)$$\n",
    "这告诉我们，对于特定的$\\Theta$只要$Q_i(z^{(i)})$等于给定$x^{(i)}$下$z^{(i)}$的后验分布即可让不等式的等号成立。\n",
    "\n",
    "\n",
    "### M-step\n",
    "\n",
    "在上一步，对于特定的$\\Theta$的取值，我们通过选择$Q_i(z^{(i)})$为后验分布$p(z^{(i)}|x^{(i)};\\Theta)$，找到了似然函数的下界，即以下不等式的等号成立：\n",
    "$$\\sum_{i=1}^m \\log p(x^{(i)};\\Theta)\\geq \\sum_{i=1}^m \\sum_{z^{(i)}}Q_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})} $$\n",
    "\n",
    "从而得到似然函数的准确估计。接下来我们找到令上式的右端最大的$\\Theta$\n",
    "\n",
    "$$\\Theta:=\\arg\\max_{\\Theta} \\sum_{i=1}^m \\sum_{z^{(i)}}Q_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})} $$\n",
    "\n",
    "由于我们事先不知道$\\Theta$的值，所以后验分布$p(z^{(i)}|x^{(i)};\\Theta)$一开始不会估计的很准确，因此我们通过如下迭代的方式来逐渐更新$\\Theta$的估计\n",
    "\n",
    "Repeat until convergence {  \n",
    "(E-step) For each i,set  \n",
    "$$Q^{(t)}_i(z^{(i)}):=p(z^{(i)}|x^{(i)};\\Theta^{(t)})$$\n",
    "\n",
    "(M-step) Set\n",
    "\n",
    "$$\\Theta^{(t+1)}:=\\arg\\max_{\\Theta} \\sum_{i=1}^m \\sum_{z^{(i)}}Q^{(t)}_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q^{(t)}_i(z^{(i)})}$$\n",
    "\n",
    "$$t:=t+1$$\n",
    "\n",
    "}\n",
    "\n",
    "### 收敛性的证明\n",
    "\n",
    "EM算法是保证能收敛的，我们如何证明这个结论？实际上，我们只要证明似然函数在每次迭代后都是单调递增的，即$\\ell(\\Theta^{(t+1)})\\geq \\ell(\\Theta^{(t)})$。  \n",
    "证明过程实际上是很简单的，关键在于$Q^{(t)}_i(z^{(i)})$的选择，通过$Q^{(t)}_i(z^{(i)}):=p(z^{(i)}|x^{(i)};\\Theta^{(t)})$，我们保证了Jensen不等式的等号成立：\n",
    "$$\\ell(\\Theta^{(t)})=\\sum_{i=1}^m \\sum_{z^{(i)}}Q^{(t)}_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta^{(t)})}{Q^{(t)}_i(z^{(i)})}$$\n",
    "又因为$\\Theta^{(t+1)}=\\arg\\max_{\\Theta} \\sum_{i=1}^m \\sum_{z^{(i)}}Q^{(t)}_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q^{(t)}_i(z^{(i)})}$，我们有\n",
    "\n",
    "$$\\begin{aligned}\\ell(\\Theta^{(t+1)})=\\sum_{i=1}^m \\log p(x^{(i)};\\Theta^{(t+1)})&\\geq \\sum_{i=1}^m \\sum_{z^{(i)}}Q^{(t)}_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta^{(t+1)})}{Q^{(t)}_i(z^{(i)})}\\\\&\\geq \\sum_{i=1}^m \\sum_{z^{(i)}}Q^{(t)}_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta^{(t)})}{Q^{(t)}_i(z^{(i)})} =\\ell(\\Theta^{(t)})\\end{aligned}$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "如果我们定义\n",
    "$$J(Q,\\Theta)=\\sum_{i=1}^m \\sum_{z^{(i)}}Q_i(z^{(i)})  \\log \\frac{p(x^{(i)},z^{(i)};\\Theta)}{Q_i(z^{(i)})} $$\n",
    "\n",
    "从前面的推导我们已经知道$\\ell(\\Theta)\\geq J(Q,\\Theta)$，那么EM算法可以看作是$J(Q,\\Theta)$上的坐标上升算法，E-step首先关于$Q$最大化$J$，在M-step中关于$\\Theta$优化$J$"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}