NLP/latent-dirichlet-allocation-2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# LDA算法系列——Gibbs采样"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## LDA简要回顾\n",
    "\n",
    "这一节我们将用Gibbs采样来估计LDA模型的参数。回顾上一节，LDA算法的工作流程如下图所示：\n",
    "\n",
    "![](http://7xikew.com1.z0.glb.clouddn.com/topic_model8.png)\n",
    "\n",
    "LDA产生语料集的过程如下所示：\n",
    "\n",
    ">for $k$ in $1,2,...,K$:  \n",
    "&emsp;$\\boldsymbol{\\beta}_k\\sim DIR(\\boldsymbol{\\eta})$  \n",
    "for $d$ in $1,2,...,D$:  \n",
    "&emsp;$\\boldsymbol{\\theta}_d \\sim DIR(\\boldsymbol{\\alpha})$  \n",
    "&emsp;for $n$ in $1,2,...,N_d$:  \n",
    "&emsp;&emsp;$Z_{d,n}\\sim Multi(\\boldsymbol{\\theta}_d)$  \n",
    "&emsp;&emsp;$W_{d,n}\\sim Multi(\\boldsymbol{\\beta}_{Z_{d,n}})$\n",
    "\n",
    "\n",
    "其中$N_d$表示第$d$篇档的单词数，$\\boldsymbol{\\beta}_k\\in\\mathbb{R}^V$表示第$k$个主题下词的分布，$\\boldsymbol{\\theta}_d\\in\\mathbb{R}^K$表示第$d$篇文档的主题分布,$W_{d,n}$表示第$d$篇文档的第$n$个词，$Z_{d,n}$表示$W_{d,n}$对应的主题，$\\boldsymbol{\\alpha}=(\\alpha_1,...,\\alpha_K)^\\top$是$\\boldsymbol{\\theta}_d$的先验分布的超参数，$\\boldsymbol{\\eta}=(\\eta_1,...,\\eta_V)^\\top$是$\\boldsymbol{\\beta}_k$的先验分布超参数（一般我们会做一个简化，认为各个维度的超参数均相等，即$\\alpha_1=...=\\alpha_K=\\alpha$，$\\eta_1=...=\\eta_V=\\eta$）\n",
    "\n",
    "我们的目标是估计LDA模型的三组参数：\n",
    "\n",
    "* 主题-词分布矩阵：$\\boldsymbol{B}=\\{\\boldsymbol{\\beta}_k\\}_{k=1}^K$\n",
    "* 文档-主题分布矩阵：$\\boldsymbol{\\Theta}=\\{\\boldsymbol{\\theta}_d\\}_{d=1}^D$\n",
    "* 隐变量矩阵：$\\boldsymbol{Z}=\\{\\boldsymbol{z}_{d}\\}_{d=1}^D$\n",
    "\n",
    "\n",
    "更具体来说，就是根据给定的$\\boldsymbol{W}, \\boldsymbol{\\alpha}, \\boldsymbol{\\eta}$求解$\\boldsymbol{Z},\\boldsymbol{\\beta},\\boldsymbol{\\theta}$的后验分布"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 直接求解后验的可行性分析\n",
    "\n",
    "首先我们来求解隐变量$\\boldsymbol{z}$的后验分布，这里我们假定主题-词分布$\\boldsymbol{\\beta}$和文档-主题分布$\\boldsymbol{\\theta}$已知\n",
    "\n",
    "\n",
    "考察$\\boldsymbol{z}$的后验概率分布\n",
    "$$p(\\boldsymbol{z}|\\boldsymbol{w})=\\frac{p(\\boldsymbol{w},\\boldsymbol{z})}{p(\\boldsymbol{w})}$$\n",
    "\n",
    "这里，我特别用粗体的$\\boldsymbol{w}$和$\\boldsymbol{z}$来分别表示语料库所有的单词构成的向量以及其对应的主题构成的向量\n",
    "\n",
    "我们来考察后验分布的分母$p(\\boldsymbol{w})$，根据全概率公式，可以将其写为\n",
    "$$p(\\boldsymbol{w})=\\sum_{z}p(\\boldsymbol{w},\\boldsymbol{z})$$\n",
    "\n",
    "具体地，我们来看其中的一个单词$w_j$的概率\n",
    "$$p(w_j)=\\sum_{k=1}^K p(w_j,z_j=k)=\\sum_{k=1}^K p(w_j|z_j=k)p(z_j=k)$$\n",
    "上一篇文章提到，LDA是一种词袋模型，因此词的产生过程是独立的，那么\n",
    "$$p(\\boldsymbol{w})=\\prod_{j=1}^{N_{w}} p(w_j)=\\prod_{j=1}^{N_{w}} \\sum_{k=1}^K p(w_j|z_j=k)p(z_j=k)$$\n",
    "其中，$N_{w}$是文档集的总单词数。我们看到，计算$p(\\boldsymbol{w})$的时间复杂度是$O(K^{N_w})$，假设有200个主题，10万个词，那么时间复杂度就是$O(200^{100000})$，因此直接求解后验是不可行的。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Collapsed Gibbs Sampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "针对此问题，不同的研究者提出了不同的近似方法来估计后验分布，其中最主要的流派有2种：  \n",
    "1）变分贝叶斯  \n",
    "2）MCMC  \n",
    "\n",
    "本文要介绍的Collapsed Gibbs Sampling（坍缩吉布斯采样）算法就属于MCMC算法的一个变种，关于MCMC及Gibbs sampling算法的原理，请参考[这篇笔记](../YiDaXu-ML/sampling-methods-part2.ipynb)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "我们先来考察生成语料库的概率，即$\\boldsymbol{w},\\boldsymbol{z}$的联合概率分布$p(\\boldsymbol{w},\\boldsymbol{z}|\\boldsymbol{\\alpha},\\boldsymbol{\\eta})$\n",
    "\n",
    "$$\\begin{aligned}\n",
    "p(\\boldsymbol{w},\\boldsymbol{z}|\\boldsymbol{\\alpha},\\boldsymbol{\\eta})&=p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{\\eta})p(\\boldsymbol{z}|\\boldsymbol{\\alpha})\\\\&=\n",
    "\\int p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})p(\\boldsymbol{B}|\\boldsymbol{\\eta})d\\boldsymbol{B} \\cdot \\int p(\\boldsymbol{z}|\\boldsymbol{\\Theta})p(\\boldsymbol{\\Theta}|\\boldsymbol{\\alpha})d\\boldsymbol{\\Theta}\n",
    "\\end{aligned}\n",
    "$$\n",
    "\n",
    "\n",
    "这里，有些读者可能会疑惑\n",
    "$$p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{\\eta})=\\int p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})p(\\boldsymbol{B}|\\boldsymbol{\\eta})d\\boldsymbol{B}$$\n",
    "\n",
    "是如何得到的，这里简要给出证明:\n",
    "\n",
    "根据全概率公式有\n",
    "\n",
    "$$ p(\\boldsymbol{w}|\\boldsymbol{z}) = \\int p(\\boldsymbol{w},\\boldsymbol{B}|\\boldsymbol{z})d\\boldsymbol{B} $$\n",
    "\n",
    "注：由于$\\eta$是常数，故可以省略\n",
    "\n",
    "根据条件概率公式有\n",
    "\n",
    "$$p(\\boldsymbol{w},\\boldsymbol{B}|\\boldsymbol{z})=\\frac{p(\\boldsymbol{w},\\boldsymbol{B},\\boldsymbol{z})}{p(\\boldsymbol{z})}$$\n",
    "\n",
    "回顾开篇的概率图模型，我们知道$W$与$Z,\\boldsymbol{B}$相关，因此希望将上式表示成关于条件概率$p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})$的形式，那么有\n",
    "\n",
    "$$p(\\boldsymbol{w},\\boldsymbol{B}|\\boldsymbol{z})=\\frac{p(\\boldsymbol{w},\\boldsymbol{B},\\boldsymbol{z})}{p(\\boldsymbol{z})}=\\frac{p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})p(\\boldsymbol{z},\\boldsymbol{B})}{p(\\boldsymbol{z})}=p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})p(\\boldsymbol{B})$$\n",
    "\n",
    "注意到以上推导用到了$Z$与$\\boldsymbol{B}$相互独立的性质：\n",
    "$$p(\\boldsymbol{z},\\boldsymbol{B})=p(\\boldsymbol{z})p(\\boldsymbol{B})$$\n",
    "\n",
    "综合上述推导过程不难得出结论"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "根据LDA生成文档的流程可知，$W$服从多项式分布$Multi(\\boldsymbol{\\beta}_{Z_{d,n}})$，其参数由$Z$和$\\boldsymbol{\\beta}$共同决定，因此\n",
    "\n",
    "$$p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})=\\prod_{d=1}^D\\prod_{n=1}^N p(w_{d,n}|\\boldsymbol{\\beta}_{z_{d,n}})$$\n",
    "\n",
    "又由于$\\boldsymbol{\\beta}_k$相互独立，于是\n",
    "$$p(\\boldsymbol{B}|\\boldsymbol{\\eta})=\\prod_{k=1}^K p(\\boldsymbol{\\beta}_k|\\boldsymbol{\\eta})$$\n",
    "\n",
    "同理可得\n",
    "\n",
    "$$ p(\\boldsymbol{z}|\\boldsymbol{\\Theta}) = \\prod_{d=1}^D\\prod_{n=1}^N p(z_{d,n}|\\boldsymbol{\\theta}_d)$$\n",
    "\n",
    "$$ p(\\boldsymbol{\\Theta}|\\boldsymbol{\\alpha}) = \\prod_{d=1}^D p(\\boldsymbol{\\theta}_d|\\boldsymbol{\\alpha}) $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "综上\n",
    "\\begin{equation}\\begin{aligned}p(\\boldsymbol{w},\\boldsymbol{z}|\\boldsymbol{\\alpha},\\boldsymbol{\\eta})&=\\int p(\\boldsymbol{w}|\\boldsymbol{z},\\boldsymbol{B})p(\\boldsymbol{B}|\\eta)d\\boldsymbol{B} \\cdot \\int p(\\boldsymbol{z}|\\boldsymbol{\\Theta})p(\\boldsymbol{\\Theta}|\\boldsymbol{\\alpha})d\\boldsymbol{\\Theta}\\\\&=\\int \\prod_{d=1}^D\\prod_{n=1}^N p(w_{d,n}|\\boldsymbol{\\beta}_{z_{d,n}})\\cdot \\prod_{k=1}^K p(\\boldsymbol{\\beta}_k|\\boldsymbol{\\eta})d\\boldsymbol{B}\\\\&\\cdot\\int\\prod_{d=1}^D\\big(\\prod_{n=1}^N p(z_{d,n}|\\boldsymbol{\\theta}_d)p(\\boldsymbol{\\theta}_d|\\boldsymbol{\\alpha})\\big)d\\boldsymbol{\\Theta}\n",
    "\\end{aligned}\\end{equation}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "由于$w_{d,n}\\sim Multi(\\beta_{z_{d,n}}) $，设$z_{d,n}=k$，那么\n",
    "\n",
    "$$p(w_{d,n}=v|\\boldsymbol{\\beta}_{k})=\\boldsymbol{\\beta}_{k,v}$$\n",
    "\n",
    "因此\n",
    "\n",
    "\\begin{equation}\\prod_{d=1}^D\\prod_{n=1}^N p(w_{d,n}|\\boldsymbol{\\beta}_{z_{d,n}})=\\prod_{k=1}^K\\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{C_{k,v}}\\end{equation}\n",
    "\n",
    "其中$C_{k,v}$表示词$w_v$被指派到主题$k$的次数  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "由于$\\boldsymbol{\\beta}_k\\sim DIR(\\boldsymbol{\\eta})$,因此\n",
    "\n",
    "\\begin{equation}\\prod_{k=1}^K p(\\boldsymbol{\\beta}_k|\\boldsymbol{\\eta})=\\prod_{k=1}^K\\frac{\\Gamma(\\sum_{v=1}^V\\eta_v)}{\\prod_{v=1}^V\\Gamma(\\eta_v)}\\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{\\eta_v-1}=\\frac{\\Gamma(V\\eta)}{\\big(\\Gamma(\\eta)\\big)^V}\\prod_{k=1}^K \\prod_{v=1}^V \\beta_{k,v}^{\\eta-1}\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "由于$z_{d,n}$服从多项式分布$z_{d,n}\\sim Multi(\\boldsymbol{\\theta}_d)$，故\n",
    "\n",
    "$$p(z_{d,n}=k|\\boldsymbol{\\theta}_d)=\\boldsymbol{\\theta}_{d,k}$$\n",
    "\n",
    "由$\\boldsymbol{\\theta}_d\\sim DIR(\\alpha_1,...,\\alpha_K)$得\n",
    "\n",
    "$$\\prod_{d=1}^D p(\\boldsymbol{\\theta}_d|\\alpha) = \\prod_{d=1}^D \\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\prod_{k=1}^K\\boldsymbol{\\theta}_{d,k}^{\\alpha_k-1}\\big)$$\n",
    "\n",
    "于是\n",
    "\n",
    "\\begin{equation}\\prod_{d=1}^D\\big(\\prod_{n=1}^N p(z_{d,n}|\\boldsymbol{\\theta}_d)p(\\boldsymbol{\\theta}_d|\\boldsymbol{\\alpha})\\big) = \\prod_{d=1}^D\\prod_{k=1}^K \\boldsymbol{\\theta}_{d,k}^{C_{d,k}}\\cdot \\prod_{d=1}^D \\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\prod_{k=1}^K\\boldsymbol{\\theta}_{d,k}^{\\alpha_k-1}\\big)\\end{equation}\n",
    "\n",
    "其中$C_{d,k}$表示第$d$篇文档中被指派到主题$k$的单词数，$C_{d,k}=\\sum_{n=1}^N \\mathbb{1}(z_{d,n}=k)$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "将$(2)$,$(3)$,$(4)$代入$(1)$得：\n",
    "\n",
    "$$\\begin{equation}\\begin{aligned}p(\\boldsymbol{w},\\boldsymbol{z}|\\alpha,\\eta)&=\n",
    "\\int\\prod_{k=1}^K\\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{C_{k,v}}\\prod_{k=1}^K \\big( \\frac{\\Gamma(\\sum_{v=1}^V \\eta_v)}{\\prod_{v=1}^V \\Gamma(\\eta_v)}\\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{\\eta_v-1} \\big)d\\boldsymbol{\\beta}\\\\&\n",
    "\\cdot \\int \\prod_{d=1}^D\\prod_{k=1}^K \\boldsymbol{\\theta}_{d,k}^{C_{d,k}}\\prod_{d=1}^D\\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\prod_{k=1}^K\\boldsymbol{\\theta}_{d,k}^{\\alpha_k-1}\\big) d\\boldsymbol{\\theta} \\\\&=\n",
    "\\big(\\frac{\\Gamma(\\sum_{v=1}^V \\eta_v)}{\\prod_{v=1}^V \\Gamma(\\eta_v)})^K\\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\big)^D\\\\&\\cdot\n",
    "\\int \\prod_{k=1}^K\\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{C_{k,v}+\\eta_v-1}d\\boldsymbol{\\beta}\\int\\prod_{d=1}^D\\prod_{k=1}^K\\boldsymbol{\\theta}_{d,k}^{C_{d,k}+\\alpha_k-1}d\\boldsymbol{\\theta}\\\\&=\\big(\\frac{\\Gamma(\\sum_{v=1}^V \\eta_v)}{\\prod_{v=1}^V \\Gamma(\\eta_v)}\\big)^K\\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\big)^D\\\\&\\cdot \n",
    "\\prod_{k=1}^K\\int \\prod_{v=1}^V \\boldsymbol{\\beta}_{k,v}^{C_{k,v}+\\eta_v-1}d\\boldsymbol{\\beta}\\cdot \\prod_{d=1}^D\\int\\prod_{k=1}^K\\boldsymbol{\\theta}_{d,k}^{C_{d,k}+\\alpha_k-1}d\\boldsymbol{\\theta} \\\\&=\n",
    "\\big(\\frac{\\Gamma(\\sum_{v=1}^V \\eta_v)}{\\prod_{v=1}^V \\Gamma(\\eta_v)})^K\\big(\\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}\\big)^D\\\\&\\cdot\n",
    "\\prod_{k=1}^K \\frac{\\prod_{v=1}^V \\Gamma(C_{k,v}+\\eta_v)}{\\Gamma(\\sum_{v=1}^V (C_{k,v}+\\eta_v))}\\cdot \\prod_{d=1}^D\\frac{\\prod_{k=1}^K \\Gamma(C_{d,k}+\\alpha_k)}{\\Gamma(\\sum_{k=1}^K (C_{d,k}+\\alpha_k))}\\end{aligned}\\end{equation} $$\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(5)式的结果太冗长，我们可以做个简化，令\n",
    "\n",
    "$$\\Delta(\\boldsymbol{\\alpha})=\\frac{\\prod_{k=1}^K\\Gamma(\\alpha_k)}{\\Gamma(\\sum_{k=1}^K\\alpha_k)}$$\n",
    "\n",
    "则\n",
    "\n",
    "$$\\frac{\\Gamma(\\sum_{v=1}^V \\eta_v)}{\\prod_{v=1}^V \\Gamma(\\eta_v)}= \\frac{1}{\\Delta(\\boldsymbol{\\eta})}\\\\ \\frac{\\Gamma(\\sum_{k=1}^K \\alpha_k)}{\\prod_{k=1}^K \\Gamma(\\alpha_k)}= \\frac{1}{\\Delta(\\boldsymbol{\\alpha})}\\\\\\frac{\\prod_{v=1}^V \\Gamma(C_{k,v}+\\eta_v)}{\\Gamma(\\sum_{v=1}^V (C_{k,v}+\\eta_v))}=\\Delta(\\boldsymbol{C}_{k,\\cdot}+\\boldsymbol{\\eta})\\\\\\frac{\\prod_{k=1}^K \\Gamma(C_{d,k}+\\alpha_k)}{\\Gamma(\\sum_{k=1}^K (C_{d,k}+\\alpha_k))}=\\Delta(\\boldsymbol{C}_{d,\\cdot}+\\boldsymbol{\\alpha})$$\n",
    "\n",
    "于是\n",
    "\n",
    "$$p(\\boldsymbol{w},\\boldsymbol{z}|\\boldsymbol{\\alpha},\\boldsymbol{\\eta})= \\prod_{k=1}^K \\frac{\\Delta(\\boldsymbol{C}_{k,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{\\eta})}\\cdot\n",
    "\\prod_{d=1}^D\\frac{\\Delta(\\boldsymbol{C}_{d,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{\\alpha})}$$\n",
    "\n",
    "其中\n",
    "* $\\boldsymbol{C}_{k,\\cdot}=(C_{k,1},...,C_{k,v})^\\top$ \n",
    "* $\\boldsymbol{C}_{d,\\cdot}=(C_{d,1},...,C_{d,k})^\\top$\n",
    "* $\\Delta(\\cdot)$是多元Beta函数\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来我们用Collapsed Gibbs采样算法来推导$z_{m,n}$的参数估计。记$\\boldsymbol{z}_{\\neg(m,n)}$表示$\\boldsymbol{z}$去掉$z_{m,n}$后剩余隐变量构成的集合，并假设$w_{m,n}=j$\n",
    "\n",
    "\n",
    "$$\\begin{equation}\\begin{aligned}p(z_{m,n}=i|\\boldsymbol{w},\\boldsymbol{z}_{\\neg(m,n)})&=\n",
    "\\frac{p(\\boldsymbol{z},\\boldsymbol{w})}{p(\\boldsymbol{w},\\boldsymbol{z}_{\\neg(m,n)})}\\\\&=\n",
    "\\bigg(\\prod_{k=1}^K \\frac{\\Delta(\\boldsymbol{C}_{k,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{\\eta})}\\cdot\n",
    "\\prod_{d=1}^D\\frac{\\Delta(\\boldsymbol{C}_{d,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{\\alpha})}\\bigg) \\bigg/ \\bigg(\\prod_{k=1}^K \\frac{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{k,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{\\eta})}\\cdot\n",
    "\\prod_{d=1}^D\\frac{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{d,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{\\alpha})}\\bigg)\\\\&=\\prod_{k=1}^K\\frac{\\Delta(\\boldsymbol{C}_{k,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{k,\\cdot}+\\boldsymbol{\\eta})}\\cdot \\prod_{d=1}^D\\frac{\\Delta(\\boldsymbol{C}_{d,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{d,\\cdot}+\\boldsymbol{\\alpha})}\\\\&=\n",
    "\\frac{\\Delta(\\boldsymbol{C}_{i,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{i,\\cdot}+\\boldsymbol{\\eta})}\\cdot \\frac{\\Delta(\\boldsymbol{C}_{m,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{m,\\cdot}+\\boldsymbol{\\alpha})}\\end{aligned}\\end{equation}$$\n",
    "\n",
    "对上面公式的推导过程做一下说明：\n",
    "\n",
    "* 当$k\\neq i$时，去掉$z_{m,n}$不会对除去主题$i$外的其他主题造成影响，因此\n",
    "$$\\Delta(\\boldsymbol{C}_{k,\\cdot}+\\boldsymbol{\\eta})=\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{k,\\cdot}+\\boldsymbol{\\eta}), \\quad   (k\\neq i)$$\n",
    "\n",
    "* 同理，当$d\\neq m$时，去掉$z_{m,n}$不会对第$m$篇文章外的其他文章造成影响，于是\n",
    "$$\\Delta(\\boldsymbol{C}_{d,\\cdot}+\\boldsymbol{\\alpha})=\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{d,\\cdot}+\\boldsymbol{\\alpha}), \\quad (d\\neq m)$$\n",
    "\n",
    "因此，当$k\\neq i$或$d\\neq m$时，分子分母均相等，于是连乘符号$\\prod$全部抵消\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "接下来考察式(6)中第一项的分子：\n",
    "$$\\Delta(\\boldsymbol{C}_{i,\\cdot}+\\boldsymbol{\\eta})=\\frac{\\prod_{v=1}^V \\Gamma(C_{i,v}+\\eta_v)}{\\Gamma(\\sum_{v=1}^V (C_{i,v}+\\eta_v))}=\\frac{\\Gamma(C_{i,1}+\\eta_1)\\Gamma(C_{i,2}+\\eta_2)...\\Gamma(C_{i,V}+\\eta_{V})}{\\Gamma\\big(\\sum_{v=1}^V (C_{i,v}+\\eta_v)\\big)}$$\n",
    "\n",
    "第一项的分母：\n",
    "$$\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{i,\\cdot}+\\boldsymbol{\\eta})=\\frac{\\prod_{v=1}^V \\Gamma(C_{i,v}^{\\neg (m,n)}+\\eta_v)}{\\Gamma(\\sum_{v=1}^V (C_{i,v}^{\\neg (m,n)}+\\eta_v))}=\\frac{\\Gamma(C_{i,1}+\\eta_1)\\Gamma(C_{i,2}+\\eta_2)...\\Gamma(C_{i,j-1}+\\eta_{j-1})\\Gamma(C_{i,j}-1+\\eta_{j})\\Gamma(C_{i,j+1}+\\eta_{j+1})...\\Gamma(C_{i,V}+\\eta_{V})}{\\Gamma\\big(\\sum_{v=1}^V (C_{i,v}+\\eta_v)-1\\big)}$$\n",
    "\n",
    "\n",
    "\n",
    "我们来考察$C_{i,v}^{\\neg (m,n)}$\n",
    "\n",
    "* 当$v\\neq j$时，$C_{i,v}^{\\neg (m,n)}=C_{i,v}$\n",
    "* 当$v=j$时，$C_{i,v}^{\\neg (m,n)}=C_{i,j}-1$\n",
    "\n",
    "于是\n",
    "$$\\begin{aligned}\\frac{\\Delta(\\boldsymbol{C}_{i,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{i,\\cdot}+\\boldsymbol{\\eta})}&=\n",
    "\\frac{\\Gamma(C_{i,1}+\\eta_1)\\Gamma(C_{i,2}+\\eta_2)...\\Gamma(C_{i,V}+\\eta_{V})}{\\Gamma(C_{i,1}+\\eta_1)\\Gamma(C_{i,2}+\\eta_2)...\\Gamma(C_{i,j-1}+\\eta_{j-1})\\Gamma(C_{i,j}-1+\\eta_{j})\\Gamma(C_{i,j+1}+\\eta_{j+1})...\\Gamma(C_{i,V}+\\eta_{V})}\\\\&\\cdot\n",
    "\\frac{\\Gamma\\big(\\sum_{v=1}^V (C_{i,v}+\\eta_v)-1\\big)}{\\Gamma\\big(\\sum_{v=1}^V (C_{i,v}+\\eta_v)\\big)}\\\\&=\\frac{\\Gamma(C_{i,j}+\\eta_j)}{\\Gamma(C_{i,j}-1+\\eta_{j})}\\cdot \\frac{\\Gamma(\\sum_{v=1}^V (C_{i,v}+\\eta_v)-1)}{\\Gamma(\\sum_{v=1}^V (C_{i,v}+\\eta_v))}\\\\&=\\frac{C_{i,j}+\\eta_{j}-1}{\\sum_{v=1}^V (C_{i,v}+\\eta_v)-1}\\end{aligned}$$\n",
    "\n",
    "论文中一般习惯将式中的\"-1\"隐去，表示为\n",
    "\n",
    "$$\\frac{\\Delta(\\boldsymbol{C}_{i,\\cdot}+\\boldsymbol{\\eta})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{i,\\cdot}+\\boldsymbol{\\eta})}=\\frac{C^{\\neg(m,n)}_{i,j}+\\eta_{j}}{\\sum_{v=1}^V (C^{\\neg(m,n)}_{i,v}+\\eta_v)}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "考察第二项的分子：\n",
    "$$\\Delta(\\boldsymbol{C}_{m,\\cdot}+\\boldsymbol{\\alpha})=\\frac{\\prod_{k=1}^K \\Gamma(C_{m,k}+\\alpha_k)}{\\Gamma(\\sum_{k=1}^K(C_{m,k}+\\alpha_k))}$$\n",
    "第二项的分母：\n",
    "$$\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{m,\\cdot}+\\boldsymbol{\\alpha})=\\frac{\\prod_{k=1}^K \\Gamma(C_{m,k}^{\\neg(m,n)}+\\alpha_k)}{\\Gamma(\\sum_{k=1}^K(C_{m,k}^{\\neg(m,n)}+\\alpha_k))}$$\n",
    "\n",
    "值得指出的是上式中的$C_{m,k}^{\\neg(m,n)}$，在$k$不同取值时的变化：\n",
    "\n",
    "* 当$k\\neq i$时，$C_{m,k}^{\\neg(m,n)}=C_{m,k}$\n",
    "* 当$k=i$时，$C_{m,k}^{\\neg(m,n)}=C_{m,i}-1$\n",
    "\n",
    "类似地\n",
    "\n",
    "$$\\frac{\\Delta(\\boldsymbol{C}_{m,\\cdot}+\\boldsymbol{\\alpha})}{\\Delta(\\boldsymbol{C}^{\\neg(m,n)}_{m,\\cdot}+\\boldsymbol{\\alpha})}=\\prod_{k=1}^K \\frac{\\Gamma(C_{m,k}+\\alpha_k)}{\\Gamma(C_{m,k}^{\\neg(m,n)}+\\alpha_k)}\\cdot \\frac{\\Gamma(\\sum_{k=1}^K(C_{m,k}^{\\neg(m,n)}+\\alpha_k))}{\\Gamma(\\sum_{k=1}^K C_{d,k}+\\alpha_k)}=\\frac{C_{m,i}+\\alpha_i-1}{\\sum_{k=1}^K(C_{m,k}+\\alpha_k)-1}=\\frac{C^{\\neg(m,n)}_{m,i}+\\alpha_i}{\\sum_{k=1}^K(C^{\\neg(m,n)}_{m,k}+\\alpha_k)}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "于是\n",
    "$$p(z_{m,n}|\\boldsymbol{w},\\boldsymbol{z}_{\\neg(m,n)})=\\frac{C^{\\neg(m,n)}_{i,j}+\\eta_{j}}{\\sum_{v=1}^V (C^{\\neg(m,n)}_{i,v}+\\eta_v)}\\cdot \\frac{C^{\\neg(m,n)}_{m,i}+\\alpha_i}{\\sum_{k=1}^K(C^{\\neg(m,n)}_{m,k}+\\alpha_k)}$$\n",
    "\n",
    "如果我们引入对称超参数，即$\\alpha_1=\\alpha_2=...=\\alpha_K=\\alpha$，$\\eta_1=\\eta_2=...=\\eta_V=\\eta$，那么上式可写为\n",
    "\n",
    "$$\\begin{aligned}p(z_{m,n}=i|\\boldsymbol{w},\\boldsymbol{z}_{\\neg(m,n)})&=\\frac{C^{\\neg(m,n)}_{i,j}+\\eta}{\\sum_{v=1}^V C^{\\neg(m,n)}_{i,v}+ V\\eta}\\cdot \\frac{C^{\\neg(m,n)}_{m,i}+\\alpha }{\\sum_{k=1}^K C^{\\neg(m,n)}_{m,k}+K\\alpha}\\\\ &=\\frac{C^{\\neg(m,n)}_{i,j}+\\eta}{C^{\\neg(m,n)}_i+ V\\eta}\\cdot \\frac{C^{\\neg(m,n)}_{m,i}+\\alpha }{C^{\\neg(m,n)}_m+K\\alpha}\\end{aligned}$$\n",
    "\n",
    "其中$C_i=\\sum_{v=1}^V C_{i,v}$表示被指派到主题$i$的token数，$C_m=\\sum_{k=1}^K C_{m,k}$表示文档$m$的token数\n",
    "\n",
    "注意到对于同一篇文章来说，$C^{\\neg(m,n)}_m+K\\alpha$是一个常数，因此在对$z_{m,n}$采样时可以将其忽略从而节省计算量，故进一步简化为\n",
    "$$p(z_{m,n}=i|\\boldsymbol{w},\\boldsymbol{z}_{\\neg(m,n)})\\propto \\frac{C^{\\neg(m,n)}_{i,j}+\\eta}{C^{\\neg(m,n)}_i+ V\\eta}\\cdot (C^{\\neg(m,n)}_{m,i}+\\alpha )$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## $\\boldsymbol{\\theta}_d$、$\\boldsymbol{\\beta}_k$的参数估计\n",
    "\n",
    "上一节中我们已经得出$z_{m,n}$的采样公式，所以在多次循环采样后，我们就能得到文档集中的每个词的主题赋值，接着我们就可以获得$\\boldsymbol{\\theta}_d$、$\\boldsymbol{\\beta}_k$的参数估计\n",
    "\n",
    "在上一篇文章，根据共轭先验的原理我们知道$\\boldsymbol{\\theta}_d,\\boldsymbol{\\beta}_k$的后验服从Dirichlet分布，那么可以根据Dirichlet分布的性质计算得出它们的参数估计：\n",
    "\n",
    "\n",
    "$$\\begin{aligned}\\boldsymbol{\\theta}^*_{d,k}&=\\mathbb{E}[\\boldsymbol{\\theta}_{d,k}]\\\\&=\\frac{C_{d,k}+\\alpha_k}{\\sum_{k=1}^K (C_{d,k}+\\alpha_k)}\\\\&\\text{(引入对称超参数)}\\\\&=\\frac{C_{d,k}+\\alpha}{\\sum_{k=1}^K C_{d,k}+K\\alpha}\\end{aligned}$$\n",
    "\n",
    "\n",
    "\n",
    "$$\\begin{aligned}\\boldsymbol{\\beta}^*_{k,v}&=\\mathbb{E}[\\boldsymbol{\\beta}_{k,v}]\\\\&=\\frac{C_{k,v}+\\eta_v}{\\sum_{v=1}^V( C_{k,v}+\\eta_v)}\\\\&\\text{(引入对称超参数)}\\\\&=\\frac{C_{k,v}+\\eta}{\\sum_{v=1}^V C_{k,v}+V\\eta}\\end{aligned}$$\n",
    "\n",
    "当然它们的共轭参数估计也可以根据Gibbs采样求得，但是这么做显然没有意义。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 参考\n",
    "1.LDA数学八卦  \n",
    "2.LDA漫游指南"
   ]
  }
 ],
 "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.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}