在对MLlib中的lda进行源码分析之前，我们先熟悉下blei的LDA paper以及相关的概念。
介绍
首先，先做下约定。
- 一个词(word)是离散数据的基本单位，被定义成词汇表中的一个item，由{1,…,V}进行索引。我们将词表示成one-hot格式。使用上标来表示components，在词汇表中第v个词被表示成一个V维的向量w，其中$ w^v=1, w^u=0(u\neq{v}) $。
- 一个文档(document)是一个N个词组成的序列，表示成W=(w1,w2,…,wN)，其中wn是序列上的第n个词。
- 一个语料(corpus)是一个关于M个文档的集合，被表示成D={W1,W2,…,Wm}.
(blei论文中用黑体的w表示文档，为便于区分，此处用大写W)。
我们希望找到一个关于该语料的概率模型，它不仅可以将高概率分配给该语料的组成文档，还可以将高概率分配给其它“相似”文档。
LDA
LDA是一个关于语料的生成概率模型。基本思想是，文档可以被表示成在隐主题(latent topics)上的随机混合，其中每个主题(topic)都以一个在词(words)上的分布进行表示。
对于在语料D中的每个文档W，LDA假设如下的生成过程：
- 1.选择 N ~ Poisson(ξ)
- 2.选择 θ ~ Dir(α)
- 3.对于每个文档（N个词汇wn）:
- (a) 选择一个主题$z_n$ ~ Multinomial(θ)
- (b) 以主题$z_n$为条件，使用多项分布条件概率（multinomial probability）：$ P(w_n \mid z_n, \beta) $选中一个词$w_n$
在该基础模型上，做了一些简化假设。
- 首先，Dirichlet分布的维度k(主题变量z的维度)，假设是已知并且确定的。
- 第二，词概率由一个k x V的矩阵β进行参数化，其中$ \beta_{ij}=p(w^j=1 \mid z^i=1) $, 目前当成是一个待估计的确定量。
- 最后，Poisson猜想是不严格的，可以使用更多的实际文档长度分布。注意N对于所有其它数据生成变量(θ和z)是独立的。它是个辅助变量，我们通常忽略它的随机性。
一个k维的Dirichlet随机变量θ，它的取值为(k-1)-simplex，(一个k维向量θ，它在泛化三角形(k-1)-simplex之内，其中：$ \theta_i \geq 0, \sum_{i=1}^k\theta_i=1 $) ，在该simplex上具有以下的概率密度：
\(p(\theta|\alpha)=\frac{\Gamma(\sum_{i=1}^{k}\alpha_{i})}{\prod_{i=1}^{k}\Gamma(\alpha_{i})}{\theta_{1}}^{\alpha_{1}-1}...{\theta_{k}}^{\alpha_{k}-1}\) ……(1)
其中，参数α是一个k维的vector，相应的components上: αi > 0, 其中 Γ(x)为Gamma函数。Dirichlet是在simplex上的一个合适分布——它在指数族内，具有有限维的充分统计量（finite dimensional sufficient statistics），与multinomial分布共轭。在paper第5部分，这些属性将有助于LDA的inference和parameter estimation算法。
给定参数 α 和 β，我们可以给出关于一个主题混合θ，一个关于N个主题的z集合(每个词都有一个主题)，一个N个词W的的联合分布：
\(p(\theta,z,W|\alpha,\beta)=p(\theta|\alpha) \prod_{n=1}^{N}p(z_n|\theta)p(w_n|z_n,\beta)\) ……(2)
其中$p(z_n \mid \theta) $可以简单认为：θi对于唯一的i，有$ {z_{n}}^{i}=1 $。在θ上积分(θ上连续)，并在z上进行求和（z上离散），我们可以得到一个关于文档的边缘分布（marginal distribution）：
\(p(W|\alpha,\beta)=\int p(\theta|\alpha)(\prod_{n=1}^{N}\sum_{z_n}p(z_n|\theta)p(w_n|z_n,\beta))d\theta\) ……(3)
最后，将所有单个文档（document）的边缘分布进行连乘运算，我得到整个语料（corpus）的概率：
\[p(D|\alpha,\beta)=\prod_{d=1}^{M}\int p(\theta_{d}|\alpha)(\prod_{n=1}^{N}\sum_{z_{dn}}p(z_n|\theta)p(w_n|z_n,\beta))d\theta_{d}\]LDA可以表示成图1所示的概率图模型。有三个级别的LDA表述。
- 语料级参数：α 和 β是语料级参数（corpus-level parameters），假设在生成一个语料的过程中只抽样一次。
- 文档级变量：θd是文档级变量（document-level variable），每个文档抽样一次。
- 词级别变量：$z_{dn}$和$w_{dn}$是词级别变量（word-level variables），对于每个文档中的每个词抽样一次。
将LDA与一个简单的Dirichlet-multinomial clustering模型相区分很重要。一个经典的clustering模型将涉及到一个两层模型(two-level)，对于一个语料只抽样一次Dirichlet；对于语料中的每个文档，只选一次multinomial clustering变量；对于在基于该cluster变量条件的文档上只选择一个词集合。有了多个clustering models后，这样的模型会限制一个文档只与单个主题相关联。而LDA涉及三层(three-level)，尤其是主题节点(topic node)在单个文档中被重复抽样。在该模型下，文档可以与多个主题相关联。
图一：LDA的图形化模型表示。”plates”表示可重复。outer plate表示文档，inner plate表示在单个文档中关于主题和词汇的可重复的选择。
通常在贝叶斯统计建模（Bayesian statistical modeling）中学到与图1相似的结构。它们被称为层级模型（hierarchical models），或更精确地称为：条件独立层级模型（con-ditionally independent hierarchical models）(Kass and Steffey, 1989).这样的模型经常作为参数期望贝叶斯模型（parametric empirical Bayes models）被提到，它也被用于参数估计中。在第5部分，我们将采用经验贝叶斯方法（empirical Bayes approach），来估计在LDA的简单实现中的α 和 β, 我们也会考虑更完整的Bayesian方法。
2.1 LDA与可交换性
随机变量 {z1,…,zN}的一个有限集合，如果联合分布在置换(permutation)后不变，我们称之为是可交换的（exchangeable）。如果π是从整数1到N的一个置换(permutation)：
\[p(z_1,...,z_N)=p(z_{\pi(1)},...,z_{\pi(N)})\]一个无限的随机变量序列，如果每个有限的子序列是可交换的，那么就是无限可交换的（infinitely exchangeable）。
De Finetti’s representation理论声称：一个无限可交换的随机变量序列的联合概率，如果一个随机参数从这些分布中抽取，那么问题中的随机变量，在其于该参数条件下是独立同分布的（iid: independent and identically distributed）。
在LDA中，我们假设，词由主题（确定的条件分布）生成，这些主题在一个文档中是无限可交换的(infinitely exchangeable)。由de Finetti的理论，一个词和主题的序列的概率，必须具有以下形式：
\[p(W,z)=\int p(\theta) (\prod_{n=1}^{N} p(z_n|\theta) p(w_n|z_n)) d\theta\]其中，θ是在这些主题上的一个服从多项分布的随机变量。我们可以在等式(3)中获取LDA分布，通过对主题变量边缘化，并让θ服从一个Dirichlet分布。
2.2 其它模型
unigram模型
在unigram模型中，每个文档中的词都从一个单独的multinomial分布独立抽取得到：
\[p(W)=\prod_{n=1}^{N} p(w_n)\]Mixture of unigrams
如果unigram模型和一个离散随机主题变量z混合，我们就得到了一个mixture of unigrams model。在该混合模型下，每个文档的生成：首先选中一个主题z，然后以条件多项概率$p(w \mid z)$独立生成N个词。一个文档的概率为：
\[p(W)=\sum_{z} p(z) \prod_{n=1}^{N} p(w_n|z)\]当从一个语料估计时，词分布可以看成是：假设每个文档只有一个主题，在该假设下的主题表示。该假设对于大语料集来说是一个很受限的模型。
相反的，LDA允许文档展现多个主题。这会额外带来另一个参数形销：有k-1参数与mixture of unigrams中的p(z)有关，而LDA模型中有k个参数与$ p(\theta \mid \alpha)$有关。
pLSI
pLSI是另一个广泛使用的文档模型(Hofmann 1999)。pLSI模型认为，一个文档d，词wn，对于一个未观察到的主题z来说是条件独立的：
\[p(d,w_n)=p(d)\sum_z p(w_n|z) p(z|d)\]pLSI模型尝试放宽maxture of unigrams模型中作出的简化猜想：每个文档只由一个主题生成。在某种意义上，它会捕获概率：一个文档可以包含多个主题，因为p(z | d)可看成是对于一个特定文档d的主题权重的混合。然后，注意，d是一个dummy index，指向在训练集中文档列表。这样，d是一个多项分布随机变量，值尽可能与训练文档一样多，模型学到的主题混合 p(z | d)只对应于被训练的文档。出于该原因，pLSI并不是一种定义良好的（well-defined）文档生成模型；天然不支持使用它来分配概率给一个之前未见过的文档。
图3: 不同模型的图表示
2.3 几何学解释
LDA与其它隐主题模型不同，可以通过隐空间(latent space)的几何表示来解释。来看下在每种模型下，一个文档在几何上是如何表示的。
上面所描述的所有4种模型（unigram, mixture of unigrams, pLSI, LDA），都是在词的空间分布上操作。每个这样的分布，可以看成是在(V-1)-simplex上的一个点，我们称之为word simplex。simplex是泛化三角形，关于simplex，参见simplex的wiki
unigram模型是在word simplex上找到单个点，假定语料中的所有词都来自该分布。隐变量模型（letent variable models）会考虑在word simplex上的k个点，并基于这些点形成一个子simplex(sub-simplex)，我们称之为topic simplex。注意，在topic simplex上的任意点，在word simplex也是一个点。使用topic simplex的不同隐变量模型，以不同的方式生成一个文档。
- 1-gram混合模型（mixture of unigram model）：假设对于每个文档，在word simplex上的k个点（也就是说：topic simplex的其中一个角落）被随机选中，该文档的所有词汇都从这些点的分布上抽取得到。
- pLSI模型：假定一个训练文档的每个词都来自一个随机选中的topic。这些主题(topics)本身从一个在这些主题上的指定文档分布上抽取到，例如，在topic simplex上的一个点。对于每个文档都有一个这样的分布；训练文档的集合定义了在topic simplex上的一个经验分布。
- LDA: 假定，已观察到和未观察到的文档上的每个词，都由一个随机选中的topic生成，这个topic从一个随机选中的参数的分布上抽取。该参数从来自topic simplex的一个平滑分布上的每个文档中抽样得到。
图4: 嵌在包含三个词的word simplex上的三个主题的topic simplex。word simplex的角(corners)对应于三个分布，每个词各自都具有一个概率分布。topic simplex的三个顶点对应于在词上的三个不同分布。The mixture of unigrams模型会将每个文档放到topic simplex上其中一个角(corners)上。pLSI模型会引入根据x的topic simplex的一个经验分布。LDA则在由等高线表示的topic simplex上放置一个平滑分布。
3.推断与参与估计
3.1 inference
为了使用LDA，我们需要解决的核心推断问题是，对于一个给定文档，计算这些隐变量的后验分布(主题分布)：
\[p(\theta,z|W,\alpha,\beta)=\frac{p(\theta,z,W|\alpha,\beta)}{p(W|\alpha,\beta)}\]不幸的是，该分布很难计算。为了归一化该分布，我们对隐变量边缘化（marginalize），并将等式(3)表示成模型参数的形式：
\[p(W|\alpha,\beta)=\frac{\Gamma(\sum_i\alpha_i)}{\prod_{i}\Gamma(\alpha_i)} \int (\prod_{i=1}^{k}\theta_{i}^{\alpha_i-1}) (\prod_{n=1}^{N}\sum_{i=1}^{k}\prod_{j=1}^{V}(\theta_i\beta_{ij})^{w_n^j}) d\theta\]该函数很难解，因为θ 和 β在隐主题的求和上相耦合(Dickey, 1983)。Dickey展示了该函数是一个在Dirichlet分布(可以表示成一个超几何函数)的特定扩展下的期望。它被用在一个Beyesian上下文中…
尽管后验分布很难精准推断（inference），有许多近似推断算法可以用于LDA，包括Laplace近似，变分近似（variational approximation），马尔可夫链蒙特卡罗法MCMC（jordan,1999）。本节我们描述了一种在LDA中简单的凸变分推断算法（convexity-based variational algorithm for inference）。其它方法在paper 第8部分讨论。
3.2 变分推断
凸变分推断算法的基本思想是，充分利用Jensen不等式来获取一个在log likelihood上的可调下界（jordan et al.,1999）。本质上，可以考虑一个关于下界的家族(a family of lower bounds)，由一个变分参数（variational parameters）集合进行索引。变分参数由一个优化过程进行选择，它会尝试找到最紧可能的下界（tightest possible lower bound）。
获取一个可控下界家族的一个简单方法是，考虑将原始的图示模型进行简单修改，移去一些边（edge）和节点（node）。将图1所示的LDA模型进行修改， θ 和 β之间的耦合，由于边θ, z, 和 W之间存在边(edge)。通过抛弃这些边以及W节点，生成一个更简单的图示模型，它有自由的变分参数，我们可以在隐变量上获取一个分布族。该分布族由以下的变分分布组成：
\(q(\theta,z|\gamma,\phi)=q(\theta|\gamma)\prod_{n=1}^{N}q(z_n|\phi_n)\) ……(4)
其中，Dirichlet分布参数γ和multinomial分布参数(φ1 , . . . , φN)，是自由变分参数。
图4: (左)LDA的图示模型 (右)采用变分分布来近似LDA后验的图示模型
在指定了一个简化版本的概率分布族后，下一步是建立一个优化问题，来确定变分参数γ 和 φ的值。正如在附录A中展示的，找到一个在log likelihood上的紧下限，可以直接翻译成下面的优化问题：
\((\gamma^{*},\phi^{*})=arg min_{\gamma,\phi}^{} D(q(\theta,z|\gamma,\phi)||p(\theta,z|W,\alpha,\beta))\)……(5)
变分参数的最优值，可以通过对变分分布$ q(\theta,z \mid \gamma,\phi) $和真实后验概率$ p(\theta,z \mid W,\alpha,\beta) $间的KL散度进行最小化得到。该最小化可以通过一个迭代型的定点方法(fixed-point method)完成。特别的，我们在附录A.3中展示了：通过计算KL散度的导数，将它们等于0, 可以得到以下更新等式的pair：
\(\phi_{ni} \propto \beta_{iw_n}exp\{E_q[log(\theta_i)|\gamma]\}\)……(6)
\(\gamma_{i}=\alpha_{i}+\sum_{n=1}^N\phi_{ni}\) ……(7)
在附录A.1中，多项分布的期望更新可以以如下方式计算：
\(E_q[log(\theta_i)|\gamma]=\Phi(\gamma_i)-\Phi(\sum_{j=1}^{k}\gamma_{j})\) ……(8)
其中，Ψ 是logΓ函数通过Taylor展开式的首个导数项。
等式(6)和(7)具有一个吸引人的直觉解释。Dirichlet分布更新(Dirichlet update)是一个在给定变分分布($ E[z_n \mid \phi_n] $>)下给定的期望观察的后验Dirichlet。多项分布更新（multinomial update）类似于使用贝叶斯理论，$ p(z_n \mid w_n) \propto p(w_n \mid z_n) $，其中p(zn)是在变分分布下的期望值的log的指数近似。
注意，变分分布实际上是一个条件分布，由一个关于W的函数区分。因为在等式(5)中的优化问题，由确定的W所管理，因而，最优解($ \gamma^{\ast},\phi^{\ast} $)是一个关于W的函数。我们可以将产生的变分分布写成：$ q(\theta,z \mid \gamma^{\ast}(W),\phi^{\ast}(W)) $，其中，我们已经在W上显式地作出了独立性。这样，变分分布可以看成是一个对后验分布$ p(\theta,z \mid W,\alpha,\beta)$的近似。
在文本语言中，最优参数（$ \gamma^{\ast}(W), \phi^{\ast}(W) $）是由文档决定的(document-specific)。特别的，我们可以将Dirichlet参数$ \gamma^{\ast}(W) $看成是在topic simplex中提供一个文档表示。
图6: LDA的variational inference算法
图6展示了变分推断过程的总结，对于 γ 和 φn具有合适的起始点。从该伪代码我们可以知道：LDA的变分推断算法，每次迭代需要O((N+1)k)次操作。经验上，我们可以找到对于单个文档的迭代数，与文档中词的数目相近似。因而总的操作数与$ N^2k $相近。
3.3 参数估计
本节描述了一种在LDA模型中用于参数估计的经验贝叶斯方法(empirical Bayes method)。特别的，给定一个文档集的语料库：D={W1,W2,…WM}。我们希望找到α 和 β，来最大化整个语料数据的(marginal)log likelihood:
\[l(\alpha,\beta)=\sum_{d=1}^{M}logp(W_d|\alpha,\beta)\]如上所述，p(W | α,β)的计算很困难。然而，变分推断提供给我们一个在log likelihood上的下界，我们可以分别根据α和β进行最大化。我们找到了近似经验贝叶斯估计，通过一种交替变分EM过程(alternating variational EM procedure)，来分别对变分参数γ 和 φ，最大化下界。接着，变分参数确定下来后，再分别根据模型参数α 和 β，最大化下界。
我们在附录A.4中提供了一个用于LDA的变分EM算法的导数的详细推导。推导以下面的迭代形式描述：
- 1.(E-step): 对于每个文档，找到变分参数{$ \gamma_d^{\ast},\phi_d^{\ast}: d \in D $}。这个在上一节已描述。
- 2.(M-step): 分别根据模型参数α 和 β，最大化在log likelihood上产生的下界。这相当于，在E-step中计算出的合适的后验概率下，为每个文档使用期望的足够统计量，找到极大似然估计。
这两个step会一直循环，直到log likelihood的下界收敛为止。
在附录A.4,我们展示了对于条件多项分布参数β进行M-step更新（M-step update），可以写成：
\(\beta_{ij} \propto \sum_{d=1}^{M}\sum_{n=1}^{N_d}\phi_{dni}^{*}w_{dn}^j\) ……(9)
我们会进一步展示Dirichlet参数α的M-step更新，它可以使用一种有效在线性时间上增长的Newton-Raphson方法，其中的Hessian是可逆的。
3.4 Smoothing
对于许多文档的大语料，具有很大的词汇表size，通常会有严重的sparsity问题。对于一个新文档，它包含的词汇很可能没有出现在训练语料中。对于这些词来说，多项分布参数（multinomial parameters）的最大似然估计会分配概率为0，因而，这些新文档的概率为0。应付该问题的方法是，对多项分布参数进行”平滑（smooth）”，对于所有词汇表item都分配一个正的概率，不管它们是否出现在训练集当中。最常用的方法是Laplace smoothing；这本质上会产生在多项分布参数上的均匀Dirichlet先验分布下，对后验分布求平均。
不幸的是，在混合模型设置中，简单的Laplace smoothing不再适合最大化后验方法。事实上，在多项分布参数上放置一个Dirichlet先验，我们可以在混合模型设置上得到一个难解的后验，这与在基本的LDA模型上得到一个难解的后验的原因基本类似。我们提出解决该问题的方法是，简单地应用变分推断方法到扩展模型上，它包括了在多项分布参数上的Dirichlet smoothing。
图7:带平滑的LDA图示模型
在LDA设置中，我们可以获得如图7所示的扩展模型。我们将β看成是一个k x V的随机矩阵（每行表示一个mixture component），其中我们假设每行都是从一个可交换的Dirichlet分布上独立抽取的。我们现在扩展我们的inference过程，将βi看成是随机变量，它天然具有一个基于训练语料条件的后验分布。此时我们来看下这个更完整的LDA Bayesian方法。
将一个变分过程看成是一个Bayesian inference，它在随机变量β, θ 和 z上放置一个独立分布(Attias,2000):
\[q(\beta_{1:k},z_{1:M},\theta_{1:M}|\lambda,\phi,\gamma)=\prod_{i=1}^{k}Dir(\beta_i|\lambda_i) \prod_{d=1}^M q_d(\theta_d,z_d|\phi_d,\gamma_d)\]其中 $ q_d(\theta,z \mid \phi,\gamma) $是在等式(4)中定义的变分分布。很容易验证，产生的变分推断过程会产生等式(6)和(7)来分别作为对变分参数φ 和 γ更新等式，同时也会对一个新的变分参数λ做一个额外更新：
\[\lambda_{ij}=\eta + \sum_{d=1}^{M} \sum_{n=1}^{N_d} \phi_{dni}^{*} w_{dn}^j\]迭代这些等式直到收敛，会产生一个在β, θ 和 z上的一个合适后验分布。
现在我们将超参数η用于可交换的Dirichlet，如同前面提到的超参数α。设置这些超参数的方法是：我们使用变分EM来寻找基于边缘似然上的这些参数的最大似然估计。这些过程在附录A.4描述。
4.online算法
对于主题模型，后验概率是很难计算的，可以通过近似后验推断来计算。当前的近似后验推导算法分为两大类：抽样方法（sampling approaches）和优化解法(optimization approaches)。抽样方法有MCMC，优化解法有VB（Variational Bayes）。VB在经验上比MCMC更快，与MCMC准确率相当，对于大数据集来说，VB是更吸引人。
但是，使用VB在大数据上计算很困难。标准的“batch”VB算法，会迭代分析每个观察（observation），并迭代更新数据集范围的变分参数。batch算法的每次迭代开销，对于非常大的数据集来说不切实际。在主题建模（topic modeling）应用中，这一点特别重要。主题建模会对不能手工标注的大型文档集合进行归纳隐式结构。
为此，blei等提出了一种online VB算法，这种算法基于online随机优化(online stochastic optimization)，在大数据集上，它比batch算法更快地产生好的参数估计。Online LDA可以方便地分析大量文档集合，不需要本地存储，或者收集文档，每个文档可以以流式(streaming)到达，看完后抛弃。
以下的等式重新标号。
假设有K个主题。每个主题定义了一个在词汇表上的多项分布，假设这些主题从一个Dirichlet上抽取得到: βk ∼ Dirichlet(η)。给定主题后，LDA为每个文档d假设了以下的生成过程。
- 首先，从主题$ \theta_d \sim Dirichlet(\alpha) $上抽取一个分布；
- 接着，对于在文档中的每个词i，从主题权重$ z_{di} \sim \theta_d $ 上抽取一个主题索引$ z_{di} \in {1, . . . , K }$，从选中的主题上抽取观察词$w_{di}$，$w_{di} \sim \beta_{z_{di}}$。出于简洁性，我们假设在θ 和 β上具有对称先验，但该假设很容易放宽。
注意，如果我们在主词分配z上进行求和，接着，我们得到：$p(w_{di} \mid \theta_d,\beta)=\sum_k\theta_{dk}\beta_{kw}$。这会导致LDA的”multinomial PCA”解释。我们可以将LDA看成是将一个词频n的矩阵(其中$n_{dw}$是词w出现在文档d中的次数)，概率分解(probabilistic factorization)成一个关于主题权重θ的矩阵，以及一个主题字典β。我们的工作可以看成是在线矩阵分解技术的一个扩展，它会最优化squared error来得到更泛化的概率公式。
- 主题(the topics): β
- 主题比例(topic proportions): θ
- 主题分配(topic assignments): z
4.1 Batch算法
在VB推断中，真实的后验由一个更简单的分布q(z,θ,β)来近似，它由一个自由参数集进行索引。这些参数用于优化：最大化置信下界(ELBO:Evidence Lower BOund)
\(logp(W|\alpha,\eta) \geq \mathcal{L} (W,\phi,\gamma,\lambda) \triangleq E_q[logp(W,z,\theta,\beta|\alpha,\eta)] - E_q[logq(z,\theta,\beta)]\) …… (1)
最大化ELBO，相当于最小化q(z,θ,β)与后验$p(z,\theta,\beta \mid W,\alpha,\eta)$间的KL散度。我们选择一个完整的因子分布q（fully factorized distribution）:
\[q(z_{di}=k)=\phi_{dw_{di}k}\] \[q(\theta_d)=Dirichlet(\theta_d;\gamma_d)\] \[q(\beta_k)=Dirichlet(\beta_k;\lambda_k)\]……(2)
每个词的主题分配z的后验，可以由φ参数化，在每个文档上的主题权重θ可以由γ进行参数化，在主题上的后验β可以由λ进行参数化。为便于记忆，我们将λ作为”主题(the topics)”。等式1分解为：
\[\mathcal{L} (W,\phi,\gamma,\lambda) = \sum_{d} \{ E_q[log p(w_d | \theta_d,z_d,\beta)] \\ + E_q [log p(z_d | \theta_d)] - E [logq(z_d)] + E [logp(\theta_d|\alpha)] - E [log q(\theta_d)] \\ + ( E [logp(\beta |\eta)] - E [log q(\beta)])/D \}\]……(3)
注意，我们在对文档求和中引入了每个语料的项，可以通过文档D的数目进行划分。该步可以帮助我们对online VB算法求导。
接着，我们将上述期望展开成变分参数的函数。变分目标函数只依赖于$ n_dw $, 词w在文档中的出现次数。当使用VB时，文档可以由词数（word counts）来进行归纳：
\[\mathcal{L} = \sum_d \sum_w n_{dw} \sum_k \phi_{dwk}( E_q[log\theta_{dk}] + E[log\beta_{kw}] - log\phi_{dwk}) \\ - log\Gamma(\sum_k\gamma_{dk}) + \sum_k(\alpha-\gamma_{dk})E[log\theta_{dk}] + log\Gamma(\gamma_{dk}) \\+ (\sum_k -log\Gamma(\sum_w\lambda_{kw}) + \sum_w(\eta-\lambda_{kw})E[log\beta_{kw}] + log\Gamma(\lambda_{kw}))/D \\ + log\Gamma(K\alpha) - K log\Gamma(\alpha) + (log\Gamma(\mathcal{V}\eta)-\mathcal{V}log\Gamma(\eta))/D \\ \triangleq \sum_d \ell(n_d,\phi_d,\gamma_d,\lambda)\]……(4)
其中V是词汇表的size，D是是文档数目。$ \ell(n_d,\phi_d,\gamma_d,\lambda) $ 表示文档d对ELBO的贡献度。
L可以使用在变分参数φ, γ, λ 上的坐标上升法进行优化：
\[\phi_{dwk} \propto exp \{ E_q[log\theta_{dk}] + E_q[log\beta_{kw}]\}\] \[\gamma_{dk}=\alpha + \sum_{w}n_{dw}\phi_{dwk}\] \[\lambda_{kw}=\eta+\sum_d n_{dw} \phi_{dwk}\]……(5)
在q分布下logθ和logβ的期望为：
\[E[log\theta_{dk}]=\Phi(\gamma_{dk})- \Phi(\sum_{i=1}^K \gamma_{di})\] \[E[log\beta_{kw}]=\Phi(\gamma_{kw})- \Phi(\sum_{i=1}^V \gamma_{ki})\]……(6)
其中Ψ表示digamma函数。
等式5的更新可以保证ELBO收敛到一个稳定点。通过EM算法（Expectation-Maximization）来完成，我们可以将这些更新划分到”E-step”：它会保持λ固定，迭代更新γ 和 φ 直到收敛，接着会进行”M-step”：由给定φ更新λ。实例上，如果在每个E-step上重新初始化γ 和 φ，该算法会收敛到一个更好的解。算法1就是batch VB for LDA:
Online算法
算法1具有常数量的内存开销，经验上，会比batch collpased Gibbs sampling收敛更快。然而，它仍需要在每次迭代时，将整个语料库传进去。因而，如果对于非常大的数据集会很慢，天然不适合不断到来的新数据。我们提出了一种在线变分推断算法来拟合λ，该参数是在主题分布β上的变分后验。我们的算法几科与batch VB算法一样简单，但对于大数据集收敛更快。
要想让主题参数λ的设置更好，就要让在算法1中的E-step中拟合得到每个文档的变分参数γ和φ之后，衡量ELBO的L要尽可能地高。将γ(nd, λ)和φ(nd, λ)设置为由E-step产生的γd 和 φd。我们的学习目标是设置λ，使下面目标最大化：
\(\mathcal{L}(n,\lambda) \triangleq \sum_d \ell (n_d,\gamma(n_d,\lambda), \phi(n_d,\lambda), \lambda)\) ……(7)
$ \ell(n_d,\gamma_d,\phi_d,\lambda) $是第d个文档对等式(4)的变分上界的贡献度。这类似于最小平方矩阵分解的目标函数，尽管LDA的ELBO比简单的squared loss function更不方便。
Online VB for LDA在算法2中描述。
第t个词频向量$n_t$被观察看，保持λ固定，我们执行一个E-step来找到局部最优解$\gamma_{t}$和$\phi_{t}$。接着，我们计算$\hat{\lambda}$，假如整个语料由单个文档$n_t$重复D次组成，那么λ的设置就是最优的(给定φt)。D是提供给算法的唯一文档数目，比如：一个语料的size。（在真实online的例子中：D->∞，对应于β的经验Bayes估计）。我们接着更新λ，使用一个λ的先前值和$\hat{\lambda}$进行加权平均得到。$\hat{\lambda}$的权重为：$\rho_t \triangleq (\tau_0 + t)^{-\kappa}$，其中$\kappa \in (0.5,1]$控制着$\hat{\lambda}$的旧值被遗忘的比率，τ0 ≥ 0会减缓算法的早期迭代。条件κ ∈ (0.5, 1]需要保证收敛。我们会在下一节展示，online LDA对应于在变分目标函数L上的一个随机自然梯度算法（stochastic natural gradient algorithm）。
该算法很像paper[16]中提出在online VB，在模型上有隐数据——最重要的区别是，我们使用一个近似的E-step来优化γt 和 φt, 因为我们不能精确计算条件分布$ p(z_t,\theta_t \mid \beta,n_t,\alpha) $。
Mini-batches: 在随机学习中的常用技术是，在每次更新时考虑多个观察（observations）以降噪。在online LDA中，这意味着计算$\hat{\lambda}$使用S>1的观察：
\(\hat{\lambda}_{kw}= \eta + \frac{D}{S} \sum_{S} n_{tsk} \phi_{tskw}\) ……(8)
其中$n_{ts}$是在mini-batch t中的第s个文档。该文档的变分参数$\phi_{ts}$和$\gamma_{ts}$是使用普通E-step进行拟合。注意，当S=D 和 κ = 0时，online VB就成了batch VB。
超参数估计：在batch 变分LDA中，超参数α 和 η的点估计，可以拟合给出的γ 和 λ，它使用线性时间的Newton-Raphson法。我们同样将α 和 η的更新并入到online LDA:
\[\alpha \leftarrow \alpha - \rho_{t} \hat{\alpha}(\gamma_t)\] \[\eta \leftarrow \rho_{t} \hat{\eta}(\lambda)\]……(9)
其中，$\hat{\alpha}(\gamma_t)$是Hessian逆矩阵乘以梯度$\triangledown_{\alpha}l(n_t,\gamma_t,\phi_t,\lambda$，$ \hat{\eta}(\lambda))$是Hessian的逆矩阵乘以梯度$\triangledown_{\eta}\mathcal{L}$，$\rho_{t} \triangleq (\tau_0 + t)^{-\kappa}$。
5.3 收敛分析
此处不详说，见paper。