第20章 潜在狄利克雷分配

习题20.1

  推导狄利克雷分布数学期望公式。

解答：

解答思路：

1. 给出狄利克雷分布
2. 写出狄利克雷分布的数学期望计算公式
3. 利用概率分布的归一化性质，推导规范化因子的形式
4. 将规范化因子代入原式，分别计算多元随机变量每个分量的期望

解答步骤：

第1步：狄利克雷分布

  根据书中第20章的定义20.2的狄利克雷分布的定义：

  定义20.2（狄利克雷分布） 若多元连续随机变量$\theta =({\theta }_1,{\theta }_2,\cdots ,{\theta }_k)$的概率密度函数为

$$\begin{array}{}\text{(20.2)}& p(\theta |\alpha )=\frac{\mathrm{\Gamma }\left({\displaystyle \sum _{i=1}^k{\alpha }_i}\right)}{{\displaystyle \prod _{i=1}^k\mathrm{\Gamma }({\alpha }_i)}}\prod _{i=1}^k{\theta }_i^{{\alpha }_i-1}\end{array}$$

其中${\displaystyle \sum _{i=1}^k{\theta }_i=1,{\theta }_i⩾0,\alpha =({\alpha }_1,{\alpha }_2,\cdots ,{\alpha }_k),{\alpha }_i>0,i=1,2,\cdots ,k}$，则称随机变量$\theta$服从参数为$\alpha$的狄利克雷分布，记作$\theta \sim \text{Dir}(\alpha )$。

  令

$$\begin{array}{}\text{(20.3)}& \text{B}(\alpha )=\frac{{\displaystyle \prod _{i=1}^k\mathrm{\Gamma }({\alpha }_i)}}{\mathrm{\Gamma }\left({\displaystyle \sum _{i=1}^k{\alpha }_i}\right)}\end{array}$$

则狄利克雷分布的密度函数可以写成

$$\begin{array}{}\text{(20.4)}& p(\theta |\alpha )=\frac{1}{\text{B}(\alpha )}\prod _{i=1}^k{\theta }_i^{{\alpha }_i-1}\end{array}$$

第2步：写出狄利克雷分布的数学期望计算公式

  已知概率分布$p(x)$数学期望的计算公式，即函数$x$在概率分布$p(x)$下的期望：${\displaystyle E[p(x)]=\int x\cdot p(x)\text{d}x}$。

  可得，狄利克雷分布的数学期望

$$E[p(\theta |\alpha )]=\int \theta \cdot p(\theta |\alpha )\text{d}\theta =\int \theta \cdot \frac{1}{\text{B}(\alpha )}\prod _{i=1}^k{\theta }_i^{{\alpha }_i-1}\text{d}\theta$$

其中，多元连续随机变量$\theta =({\theta }_1,{\theta }_2,\dots ,{\theta }_k),\alpha =({\alpha }_1,{\alpha }_2,\dots ,{\alpha }_k)$

  分别计算 $\theta$ 在每个分量上的期望，即

$$E[p(\theta |\alpha )]=(E[p({\theta }_1|\alpha )],E[p({\theta }_2|\alpha )],\cdots ,E[p({\theta }_k|\alpha )])$$

第3步：利用概率分布的归一化性质，推导规范化因子的形式

  对于每个分量${\theta }_i$的期望：

$$\begin{array}{}\text{(1)}& \begin{array}{rl}E[p({\theta }_i|\alpha )]& =\int {\theta }_i\cdot \frac{1}{\text{B}(\alpha )}\prod _{j=1}^k{\theta }_j^{{\alpha }_j-1}\text{d}{\theta }_i\\ & =\frac{1}{\text{B}(\alpha )}\int {\theta }_1^{{\alpha }_1-1}{\theta }_2^{{\alpha }_2-1}\cdots {\theta }_i^{{\alpha }_i}\cdots {\theta }_k^{{\alpha }_k-1}\text{d}{\theta }_i\end{array}\end{array}$$

  记上式右边的积分部可对应如下狄利克雷分布的密度函数

$$\begin{array}{rl}p({\theta }_i|{\alpha }^{\prime })& =\frac{1}{\text{B}({\alpha }^{\prime })}{\theta }_1^{{\alpha }_1-1}{\theta }_2^{{\alpha }_2-1}\dots {\theta }_i^{{\alpha }_i}\dots {\theta }_k^{{\alpha }_k-1}\\ & =\frac{\mathrm{\Gamma }({\alpha }_1+{\alpha }_2+\cdots +{\alpha }_i+1+\cdots +{\alpha }_k)}{\mathrm{\Gamma }({\alpha }_1)\mathrm{\Gamma }({\alpha }_2)\dots \mathrm{\Gamma }({\alpha }_i+1)\dots \mathrm{\Gamma }({\alpha }_k)}{\theta }_1^{{\alpha }_1-1}{\theta }_2^{{\alpha }_2-1}\dots {\theta }_i^{{\alpha }_i}\dots {\theta }_k^{{\alpha }_k-1}\end{array}$$

其中${\alpha }^{\prime }=({\alpha }_1,{\alpha }_2,\cdots ,{\alpha }_i+1,\cdots ,{\alpha }_k)$

  根据归一化条件：

$$\int p({\theta }_i|{\alpha }^{\prime })\text{d}\theta =1$$

  可得：

$$\begin{array}{}\text{(2)}& \text{B}({\alpha }^{\prime })=\int {\theta }_1^{{\alpha }_1-1}{\theta }_2^{{\alpha }_2-1}\cdots {\theta }_i^{{\alpha }_i}\cdots {\theta }_k^{{\alpha }_k-1}\text{d}\theta \end{array}$$

第3步：将规范化因子代入原式，分别计算多元随机变量每个分量的期望。

  将公式（2）代入公式（1），可得分量${\theta }_i$的期望：

$$\begin{array}{rl}E[p({\theta }_i|\alpha )]& =\frac{\text{B}({\alpha }^{\prime })}{\text{B}(\alpha )}\\ & =\frac{\mathrm{\Gamma }({\alpha }_1)\mathrm{\Gamma }({\alpha }_2)\cdots \mathrm{\Gamma }({\alpha }_i+1)\cdots \mathrm{\Gamma }({\alpha }_k)}{\mathrm{\Gamma }({\alpha }_1+{\alpha }_2+\cdots +{\alpha }_i+1+\cdots +{\alpha }_k)}\cdot \frac{\mathrm{\Gamma }({\alpha }_1+{\alpha }_2+\cdots +{\alpha }_i+\cdots +{\alpha }_k)}{\mathrm{\Gamma }({\alpha }_1)\mathrm{\Gamma }({\alpha }_2)\cdots \mathrm{\Gamma }({\alpha }_i)\cdots \mathrm{\Gamma }({\alpha }_k)}\\ & =\frac{{\displaystyle {\alpha }_i\prod _{j=1}^k\mathrm{\Gamma }({\alpha }_j)}}{{\displaystyle \sum _{j=1}^k{\alpha }_j\mathrm{\Gamma }\left(\sum _{j=1}^k{\alpha }_j\right)}}\cdot \frac{{\displaystyle \mathrm{\Gamma }\left(\sum _{j=1}^k{\alpha }_j\right)}}{{\displaystyle \prod _{j=1}^k\mathrm{\Gamma }({\alpha }_j)}}\\ & =\frac{{\alpha }_i}{{\displaystyle \sum _{j=1}^k{\alpha }_j}}\end{array}$$

  因此，狄利克雷分布的数学期望：

$$\begin{array}{rl}E[p(\theta |\alpha )]& =(E[p({\theta }_1|\alpha )],E[p({\theta }_2|\alpha )],\cdots ,E[p({\theta }_k|\alpha )])\\ & =(\frac{{\alpha }_1}{{\displaystyle \sum _{i=1}^k{\alpha }_i}},\frac{{\alpha }_2}{{\displaystyle \sum _{i=1}^k{\alpha }_i}},\cdots ,\frac{{\alpha }_k}{{\displaystyle \sum _{i=1}^k{\alpha }_i}})\end{array}$$

习题20.2

  针对17.2.2节的文本例子，使用LDA模型进行话题分析。

解答：

解答思路：

1. 给出LDA吉布斯抽样算法
2. 自编程实现基于吉布斯抽样算法的LDA模型
3. 针对17.2.2节的文本例子，使用LDA模型进行话题分析

解答步骤：

第1步：LDA吉布斯抽样算法

  根据书中第20章的算法20.2的LDA吉布斯抽样算法：

算法20.2（LDA吉布斯抽样算法）
输入：文本的单词序列$w=\{w_1,w_2,\cdots ,w_m,\cdots ,w_M\},w_m=(w_{m1},w_{m2},\cdots ,w_{mn},\cdots ,w_{m{N}_m})$；
输出：文本的话题序列$z=\{z_1,z_2,\cdots ,z_m,\cdots ,z_M\},z_m=(z_{m1},z_{m2},\cdots ,z_{mn},\cdots ,z_{m{N}_m})$的后验概率分布$p(z|w,\alpha ,\beta )$的样本计数，模型的参数$\phi$和$\theta$的估计值；
参数：超参数$\alpha$和$\beta$，话题个数$K$。
（1）设所有计数矩阵的元素$n_{mk},n_{kv}$，计数向量的元素$n_m,n_k$初值为0；
（2）对所有文本$w_m,m=1,2,\cdots ,M$，对第$m$个文本中的所有单词$w_{mn},n=1,2,\cdots ,N_m$，抽样话题$z_{mn}=z_k\sim \text{Mult}\left({\displaystyle \frac{1}{K}}\right)$；增加文本-话题计数$n_{mk}=n_{mk}+1$，增加文本-话题和计数$n_m=n_m+1$，增加话题-单词计数$n_{kv}=n_{kv}+1$，增加话题-单词和计数$n_k=n_k+1$；
（3）循环执行以下操作，直到进入燃烧期。对所有文本$w_m,m=1,2,\cdots ,M$，对第$m$个文本中的所有单词$w_{mn},n=1,2,\cdots ,N_m$；
    （a）当前的单词$w_{mn}$是第$v$个单词，话题指派$z_{mn}$是第$k$个话题；减少计数$n_{mk}=n_{mk}-1,n_m=n_m-1,n_{kv}=n_{kv}-1,n_k=n_k-1$；
    （b）按照满条件分布进行抽样

$$p(z_i|{\mathbf{\text{z}}}_{-i},\mathbf{\text{w}},\alpha ,\beta )\propto \frac{n_{kv}+{\beta }_v}{{\displaystyle \sum _{v=1}^V(n_{k}v+{\beta }_v)}}\cdot \frac{n_{mk}+{\alpha }_k}{{\displaystyle \sum _{k=1}^K(n_{mk}+{\alpha }_k)}}$$

得到新的第$k^{\prime }$个话题，分配给$z_{mn}$；
    （c）增加计数$n_{m{k}^{\prime }}=n_{m{k}^{\prime }}+1,n_m=n_m+1,n_{k^{\prime }v}=n_{k^{\prime }v}+1,n_{k^{\prime }}=n_{k^{\prime }}+1$；
    （d）得到更新的两个计数矩阵$N_{K\times V}=[n_{kv}]$和$N_{M\times K}=[n_{mk}]$，表示后验概率分布$p(\mathbf{\text{z}}|\mathbf{\text{w}},\alpha ,\beta )$的样本计数；
（4）利用得到的样本计数，计算模型参数

$$\begin{gathered} {\theta }_{mk}=\frac{n_{mk}+{\alpha }_k}{{\displaystyle \sum _{k=1}^K(n_{mk}+{\alpha }_k)}} \\ {\phi }_{kv}=\frac{n_{kv}+{\beta }_v}{{\displaystyle \sum _{v=1}^V(n_{kv}+{\beta }_v)}} \end{gathered}$$

第2步：自编程实现基于吉布斯抽样算法的LDA模型

import numpy as np


class GibbsSamplingLDA:
    def __init__(self, iter_max=1000):
        self.iter_max = iter_max
        self.weights_ = []

    def fit(self, words, K):
        """
        :param words: 单词-文本矩阵
        :param K: 话题个数
        :return: 文本话题序列z
        """
        # M, Nm分别为文本个数和单词个数
        words = words.T
        M, Nm = words.shape

        # 初始化超参数alpha, beta，其中alpha为文本的话题分布相关参数
        # beta为话题的单词分布相关参数
        alpha = np.array([1 / K] * K)
        beta = np.array([1 / Nm] * Nm)

        # 初始化参数theta, varphi，其中theta为文本关于话题的多项分布参数，
        # varphi为话题关于单词的多项分布参数
        theta = np.zeros([M, K])
        varphi = np.zeros([K, Nm])

        # 输出文本的话题序列z
        z = np.zeros(words.shape, dtype='int')

        # (1)设所有计数矩阵的元素n_mk、n_kv，计数向量的元素n_m、n_k初值为 0
        n_mk = np.zeros([M, K])
        n_kv = np.zeros([K, Nm])
        n_m = np.zeros(M)
        n_k = np.zeros(K)

        # (2)对所有M个文本中的所有单词进行循环
        for m in range(M):
            for v in range(Nm):
                # 如果单词v存在于文本m
                if words[m, v] != 0:
                    # (2.a)抽样话题
                    z[m, v] = np.random.choice(list(range(K)))
                    # 增加文本-话题计数
                    n_mk[m, z[m, v]] += 1
                    # 增加文本-话题和计数
                    n_m[m] += 1
                    # 增加话题-单词计数
                    n_kv[z[m, v], v] += 1
                    # 增加话题-单词和计数
                    n_k[z[m, v]] += 1

        # (3)对所有M个文本中的所有单词进行循环，直到进入燃烧期
        zi = 0
        for i in range(self.iter_max):
            for m in range(M):
                for v in range(Nm):
                    # (3.a)如果单词v存在于文本m，那么当前单词是第v个单词，
                    # 话题指派z_mv是第k个话题
                    if words[m, v] != 0:
                        # 减少计数
                        n_mk[m, z[m, v]] -= 1
                        n_m[m] -= 1
                        n_kv[z[m, v], v] -= 1
                        n_k[z[m, v]] -= 1

                        # (3.b)按照满条件分布进行抽样
                        max_zi_value, max_zi_index = -float('inf'), z[m, v]
                        for k in range(K):
                            zi = ((n_kv[k, v] + beta[v]) / (n_kv[k, :].sum() + beta.sum())) * \
                                 ((n_mk[m, k] + alpha[k]) / (n_mk[m, :].sum() + alpha.sum()))

                        # 得到新的第 k‘个话题，分配给 z_mv
                        if max_zi_value < zi:
                            max_zi_value, max_zi_index = zi, k
                            z[m, v] = max_zi_index

                        # (3.c) (3.d)增加计数并得到两个更新的计数矩阵的n_kv和n_mk
                        n_mk[m, z[m, v]] += 1
                        n_m[m] += 1
                        n_kv[z[m, v], v] += 1
                        n_k[z[m, v]] += 1

        # (4)利用得到的样本计数，计算模型参数
        for m in range(M):
            for k in range(K):
                theta[m, k] = (n_mk[m, k] + alpha[k]) / (n_mk[m, :].sum() + alpha.sum())

        for k in range(K):
            for v in range(Nm):
                varphi[k, v] = (n_kv[k, v] + beta[v]) / (n_kv[k, :].sum() + beta.sum())

        self.weights_ = [varphi, theta]
        return z.T, n_kv, n_mk

第3步：针对17.2.2节的文本例子，使用LDA模型进行话题分析

  17.2.2节的文本例子中的数据如下：

  使用LDA模型进行话题分析：

gibbs_sampling_lda = GibbsSamplingLDA(iter_max=1000)

# 输入文本-单词矩阵，共有9个文本，11个单词
words = np.array([[0, 0, 1, 1, 0, 0, 0, 0, 0],
                  [0, 0, 0, 0, 0, 1, 0, 0, 1],
                  [0, 1, 0, 0, 0, 0, 0, 1, 0],
                  [0, 0, 0, 0, 0, 0, 1, 0, 1],
                  [1, 0, 0, 0, 0, 1, 0, 0, 0],
                  [1, 1, 1, 1, 1, 1, 1, 1, 1],
                  [1, 0, 1, 0, 0, 0, 0, 0, 0],
                  [0, 0, 0, 0, 0, 0, 1, 0, 1],
                  [0, 0, 0, 0, 0, 2, 0, 0, 1],
                  [1, 0, 1, 0, 0, 0, 0, 1, 0],
                  [0, 0, 0, 1, 1, 0, 0, 0, 0]])

K = 3  # 假设话题数量为3

# 设置精度为3
np.set_printoptions(precision=3, suppress=True)

z, n_kv, n_mk = gibbs_sampling_lda.fit(words, K)
varphi = gibbs_sampling_lda.weights_[0]
theta = gibbs_sampling_lda.weights_[1]

print("文本的话题序列z：")
print(z)
print("样本的计数矩阵N_KV：")
print(n_kv)
print("样本的计数矩阵N_MK：")
print(n_mk)
print("模型参数varphi：")
print(varphi)
print("模型参数theta：")
print(theta)

文本的话题序列z：
[[0 0 2 2 0 0 0 0 0]
 [0 0 0 0 0 2 0 0 2]
 [0 2 0 0 0 0 0 2 0]
 [0 0 0 0 0 0 2 0 2]
 [2 0 0 0 0 2 0 0 0]
 [2 2 2 2 2 2 2 2 2]
 [2 0 2 0 0 0 0 0 0]
 [0 0 0 0 0 0 2 0 2]
 [0 0 0 0 0 2 0 0 2]
 [2 0 2 0 0 0 0 2 0]
 [0 0 0 2 2 0 0 0 0]]
样本的计数矩阵N_KV：
[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
 [2. 2. 2. 2. 2. 9. 2. 2. 2. 3. 2.]]
样本的计数矩阵N_MK：
[[0. 0. 4.]
 [0. 0. 2.]
 [0. 0. 4.]
 [0. 0. 3.]
 [0. 0. 2.]
 [0. 0. 4.]
 [0. 0. 3.]
 [0. 0. 3.]
 [0. 0. 5.]]
模型参数varphi：
[[0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091]
 [0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091 0.091]
 [0.067 0.067 0.067 0.067 0.067 0.293 0.067 0.067 0.067 0.1   0.067]]
模型参数theta：
[[0.067 0.067 0.867]
 [0.111 0.111 0.778]
 [0.067 0.067 0.867]
 [0.083 0.083 0.833]
 [0.111 0.111 0.778]
 [0.067 0.067 0.867]
 [0.083 0.083 0.833]
 [0.083 0.083 0.833]
 [0.056 0.056 0.889]]

习题20.3

  找出LDA的吉布斯抽样算法、变分EM算法中利用狄利克雷分布的部分，思考LDA中使用狄利克雷分布的重要性。

解答：

解答思路：

1. 写出LDA的吉布斯抽样算法中利用狄利克雷分布的部分
2. 写出变分EM算法中利用狄利克雷分布的部分
3. 写出LDA中使用狄利克雷分布的重要性

解答步骤：

  在话题模型中，假设话题是单词的多项分布，文本是话题的多项分布，而狄利克雷分布是多项分布的先验，在进行贝叶斯学习时，需要使用狄利克雷分布作为先验分布，因此要先找出LDA的吉布斯抽样算法和变分EM算法中的多项分布，再找出对应的狄利克雷分布。

第1步：LDA的吉布斯抽样算法

  根据书中第20.3节的LDA模型基本想法：

  LDA模型的学习通常采用收缩的吉布斯抽样方法，基本想法是，通过对隐变量$\theta$和$\phi$积分，得到边缘概率分布$p(w,z|\alpha ,\beta )$（也是联合分布），其中变量$w$是可观测的，变量$z$是不可观测的；对后验概率分布$p(z|w,\alpha ,\beta )$进行吉布斯抽样，得到分布$p(z|w,\alpha ,\beta )$的样本集合；再利用这个样本集合对参数$\theta$和$\phi$进行估计，最终得到LDA模型$p(z|w,\alpha ,\beta )$的所有参数估计。

  根据书中第20.3.2节的抽样分布的表达式：

$p(w,z|\alpha ,\beta )=p(w|z,\alpha ,\beta )p(z|\alpha ,\beta )=p(w|z,\beta )p(z|\alpha )$

1. 处理第一个因子$p(w|z,\beta )$

  该分布表示在给定话题$z$和话题-单词分布参数$\phi$下的文本的分布，是一个多项分布。

  根据书中第20.3.2节的公式(20.22)：

$$p(w|z,\phi )=\prod _{k=1}^K\prod _{v=1}^V{\phi }_{kv}^{n_{kv}}$$

其中${\phi }_{kv}$表示第$k$个话题生成单词集合第$v$个单词的概率，$n_{kv}$是数据中第$k$个话题生成第$v$个单词的次数。

  现引入该多项分布的共轭先验，有

$$p(w,\phi |z,\beta )\propto p(w|z,\phi )p(\phi |\beta )$$

  对上式两边求$\phi$的积分得：

$$\begin{array}{rl}p(w|z,\beta )& =\int p(w,\phi |z,\beta )\text{d}\phi \\ & =\int p(w|z,\phi )p(\phi |\beta )\text{d}\phi \\ & =\prod _{k=1}^K\frac{\text{B}(n_k+\beta )}{\text{B}(\beta )}\end{array}$$

  根据假设 $p(w|z,\beta ),p(\phi |\beta )$ 均是多项分布，其中第$k$个话题的单词分布参数${\phi }_k$的先验分布为

$$p({\phi }_k|\beta )=\text{Dir}({\phi }_k|\beta )$$

  根据书中第20.3.3节的公式（20.32）：

2. 参数$\phi =\{{\phi }_k\}$的估计
   后验概率满足

$$p({\phi }_k|w,z,\beta )=\text{Dir}({\phi }_k|n_k+\beta )$$

这里$n_k=\{n_{k1},n_{k2},\dots ,n_{kV}\}$是第$k$个话题的单词的计数。

2. 处理第二个因子 $p(z|\alpha )$

  该分布表示在给定文本的话题分布参数$\theta$下话题的分布，是一个多项分布。

  根据书中第20.3.2节的公式(20.24)：

$$p(z|\theta )=\prod _{m=1}^M\prod _{k=1}^K{\theta }_{mk}^{n_{mk}}$$

其中${\theta }_{mk}$表示第$m$个文本生成第$k$个话题的概率，$n_{mk}$是数据中第$m$个文本生成第$k$个话题的次数。

  现引入该多项分布的共轭先验，有

$$p(z,\theta |\alpha )\propto p(z|\theta )p(\theta |\alpha )$$

  对上式两边求$\theta$积分得：

$$\begin{array}{rl}p(z|\alpha )& =\int p(z,\theta |\alpha )\text{d}\theta \\ & =\int p(z|\theta )p(\theta |\alpha )\text{d}\theta \\ & =\prod _{m=1}^M\frac{\text{B}(n_m+\alpha )}{\text{B}(\alpha )}\end{array}$$

  根据假设 $p(z,\theta |\alpha ),p(\theta |\alpha )$均是多项分布，其中第$m$个文本的话题分布参数${\theta }_m$的先验分布为

$$p({\theta }_m|\alpha )=\text{Dir}({\theta }_m|\alpha )$$

  根据书中第20.3.3节的公式（20.30）：

1. 参数$\theta =\{{\theta }_m\}$的估计
   后验概率满足

$$\begin{array}{}\text{(20.30)}& p({\theta }_m|z_m,\alpha )=\text{Dir}({\theta }_m|n_m+\alpha )\end{array}$$

这里$n_m=\{n_{m1},n_{m2},\cdots ,n_{mK}\}$是第$m$个文本的话题的计数。

第2步：LDA的变分EM算法

  根据书中第20.4.3节的变分EM算法：

  将变分EM算法应用到LDA模型的学习上，首先定义具体的变分分布，推导证据下界的表达式，接着推导变分分布的参数和LDA模型的参数估计，最后给出LDA模型的变分EM算法。 为简单起见，一次只考虑一个文本，记作$w$。文本的单词序列$w=(w_1,\cdots ,w_n,\cdots ,w_N)$，对应的话题序列$z=(z_1,\cdots ,z_n,\cdots ,z_N)$，以及话题分布$\theta$，随机变量$w$，$z$和$\theta$的联合分布是

$$\begin{array}{}\text{(20.42)}& p(\theta ,z,w|\alpha ,\phi )=p(\theta |\alpha )\prod _{n=1}^{N}p(z_n|\theta )p(w_n|z_n,\phi )\end{array}$$

其中$w$是可观测变量，$\theta$和$z$是隐变量，$\alpha$和$\phi$是参数。

  根据第1步已证明$p(\theta |\alpha )$是狄利克雷分布，即

$$p(\theta |\alpha )=\text{Dir}(\theta |\alpha )$$

  另一部分${\displaystyle \prod _{n=1}^{N}p(z_n|\theta )p(w_n|z_n,\phi )}$表示给定话题分布$\theta$产生第$n$个文本的话题序列 $p(z_n|\theta )$的概率以及给定话题单词分布参数$\phi$和第$n$个文本的话题序列$z_n$下产生第$n$个文本$w_n$的概率，根据第1步已证明均为多项分布。

  根据书中第20.4.3节的公式(20.43)：

$$\begin{array}{}\text{(20.43)}& q(\theta ,\mathbf{\text{z}}|\gamma ,\eta )=q(\theta |\gamma )\prod _{n=1}^{N}q(z_n|{\eta }_n)\end{array}$$

其中$\gamma$是狄利克雷分布参数，$\eta =({\eta }_1,{\eta }_2,\cdots ,{\eta }_n)$是多项分布参数，变量$\theta$和$\mathbf{\text{z}}$的各个分量都是条件独立的。

  由于 $w$是可观测变量，变分EM算法就是利用上式的变分分布$q(\theta ,z|\gamma ,\eta )$来近似后验分布$p(\theta ,z|w,\alpha ,\phi )$。

  根据书中第20.4.3节的公式(20.44)：

  由此得到一个文本的证据下界

$$L(\gamma ,\eta ,\alpha ,\phi )=E_q[\log p(\theta ,z,w|\alpha ,\phi )]-E_q[\log q(\theta ,z|\gamma ,\eta )]$$

展开证据下界式

$$\begin{array}{rl}L(\gamma ,\eta ,\alpha ,\phi )& =E_q[\log p(\theta |\alpha )]+E_q[\log p(z|\theta )]\\ & +E_q[\log p(w|z,\phi )]-E_q[\log q(\theta |\gamma )]-E_q[\log q(z|\eta )]\end{array}$$

  每一项都是对应的概率分布关于变分分布$q(\theta ,z|\gamma ,\eta )$的期望，逐一分析：

1. 第一项包含$p(\theta |\alpha )$，含有狄利克雷分布；
2. 第二项$p(z|\theta )$为多项分布，其中参数$\theta$在变分分布中已假设为狄利克雷分布$q(\theta |\gamma )=\mathrm{Dir}(\theta |\gamma )$，含有狄利克雷分布；
3. 第三项$p(w|z,\phi )$为多项分布，其中参数$\theta$的分布可以根据归一化性质消掉，剩余部分均为多项分布，因此不含狄利克雷分布；
4. 第四项$q(\theta |\gamma )$ 已假设为狄利克雷分布；
5. 第五项$q(z|\eta )$，其中参数$\theta$ 的分布根据归一化性质消掉，剩余均为多项分布，故不含狄利克雷分布。

第3步：LDA中使用狄利克雷分布的重要性

  由于LDA话题模型天然具备多项分布的性质，因此在进行贝叶斯学习时，则不可避免地需要引入狄利克雷分布作为多项分布的共轭先验，没有狄利克雷分布就无法进行话题模型的贝叶斯估计。

习题20.4

  给出LDA的吉布斯抽样算法和变分EM算法的算法复杂度。

解答：

解答思路：

1. 计算LDA的吉布斯抽样算法的算法复杂度

2. 计算LDA的变分EM算法的算法复杂度：根据书上409页算法20.4和411页算法20.5，同样计算循环次数以及每次循环的代价，注意这里需要分别考虑 E 步和 M 步的代价，最后加起来即为 LDA 变分EM算法的复杂度。

解答步骤：

第1步：LDA的吉布斯抽样算法的算法复杂度

  根据书中第20章的算法20.2（LDA吉布斯抽样算法），假设迭代次数为$T$

1. 第1、2步均为初始化相关，对于所有文本 ${\mathbf{\text{w}}}_m$，需要进行$M$次循环，即文本数量；对于第$m$个文本中的所有单词 $w_{mn}$，需要进行$N_m$次循环，即第$m$个文本包含的单词数，因此一共进行$T\cdot M\cdot N_m$次循环，故复杂度为$T\cdot M\cdot N_m$

2. 第3步为迭代，每次循环过程中，第3.a、3.c和3.d步均为固定更新，算法复杂度均为1

3. 第3.b步涉及到随机采样，该公式第一个因子表示话题生成该位置单词的概率，计算复杂度为$V$，即词典大小；第二个因子表示该位置的文本生成话题的概率，计算复杂度为$K$，即话题数量；由于词典大小远大于话题数量，这里可以近似认为计算复杂度为$V$

  因此，总算法复杂度为$T\cdot M\cdot N_m\cdot V$。

第2步：LDA的变分EM算法

  根据书中第20章的算法20.5（LDA的变分EM算法），假设迭代次数为$T$

1. E步：估计变分参数$\gamma ,\eta$，假设重复 $N_E$ 次后收敛，根据书中第409页算法20.4（LDA的变分参数估计算法）中的第4、5步的循环次数为 $N_E\cdot N\cdot K$次；第6步在计算${\displaystyle \mathrm{\Psi }(\sum _{l=1}^K{\gamma }_l^{(t)})}$时的算法复杂度为$K$，其余计算的复杂度均为1
   则每次总算法复杂度为$N_E\cdot N\cdot K\cdot K$

2. M步：估计模型参数 $\alpha ,\phi$，根据书中公式(20.63)

$$L[{\phi }_w]=\sum _{m=1}^M\sum _{n=1}^{N_m}\sum _{k=1}^K\sum _{v=1}^V{\eta }_{mnk}{w}_{mn}^v\log {\phi }_{kv}+\sum _{k=1}^K{\lambda }_k(\sum _{v=1}^V{\phi }_{kv}-1)$$

对${\phi }_{kv}$求偏导并令其为零，归一化求解，得到参数${\phi }_{kv}$的估计值

$${\phi }_{kv}\propto \sum _{m=1}^M\sum _{n=1}^{N_m}{\eta }_{mnk}{w}_{mm}^v$$

其算法复杂度为$M\cdot N_m$；
再通过证据下界的最大化估计参数 $\alpha$，根据公式(20.65)

$$L_{|\alpha |}=\sum _{m=1}^M\{\log \mathrm{\Gamma }\left(\sum _{l=1}^K{\alpha }_l\right)-\sum _{k=1}^K\log \mathrm{\Gamma }({\alpha }_k)+\sum _{k=1}^K({\alpha }_k-1)[\mathrm{\Psi }({\gamma }_{mk})-\mathrm{\Psi }\left(\sum _{l=1}^K{\gamma }_{ml}\right)]\}$$

计算其关于$\alpha$的Hessian矩阵，然后应用牛顿法求该函数的最大化，假设牛顿法迭代次数为$N_M$，$\alpha$包含$K$个话题的狄利克雷分布参数，牛顿法的计算复杂度为$N_M\cdot K\cdot K$
则每次总算法复杂度为$M\cdot N_m+N_M\cdot K\cdot K$

  最后一共进行了 $T$ 次 EM 算法的迭代，因此，总算法复杂度为$T\cdot (N_E\cdot N\cdot K\cdot K+M\cdot N_m+N_M\cdot K\cdot K)$。

习题20.5

  证明变分EM算法收敛。

解答：

解答思路：

1. 给出变分EM算法
2. 给出EM算法的收敛性
3. 证明在变分EM算法中，证据下界随着迭代进行单调递增

解答步骤：

第1步：变分EM算法

  根据书中第20章的算法20.3：

  假设模型式联合概率分布$p(x,z|\theta )$，其中$x$是观测变量，$z$是隐变量，$\theta$是参数。目标是通过观测数据的概率（证据）$\log p(x|\theta )$的最大化，估计模型的参数$\theta$。使用变分推理，导入平均场$q(z)=\prod _{i=1}^{n}q(z_i)$，定义证据下界

$$\begin{array}{}\text{(20.39)}& L(q,\theta )=E_q[\log p(x,z|\theta )]-E_q[\log q(z)]\end{array}$$

通过迭代，分别以$q$和$\theta$为变量时对证据下界进行最大化，就得到变分EM算法。

算法20.3（变分EM算法）
循环执行以下E步和M步，直到收敛。
（1）E步：固定$\theta$，求$L(q,\theta )$对$q$的最大化。
（2）M步：固定$q$，求$L(q,\theta )$对$\theta$的最大化。
给出模型参数$\theta$的估计值。

第2步：EM算法的收敛性

  根据书中第9.2节的定理9.1：

  定理9.1 设$P(Y|\theta )$为观测数据的似然函数，${\theta }^{(i)}(i=1,2,\cdots )$为EM算法得到的参数估计序列，$P(Y|{\theta }^{(i)})(i=1,2,\cdots )$为对应的似然函数序列，则$P(Y|{\theta }^i)$是单调递增的，即

$$P(Y|{\theta }^{(i+1)})⩾P(Y|{\theta }^{(i)})$$

  根据书中第9.2节的定理9.2：

  定理9.2 设$L(\theta )=\log P(Y|\theta )$为观测数据的对数似然函数，${\theta }^{(i)}(i=1,2,\cdots )$为EM算法得到的参数估计序列，$L({\theta }^{(i)})(i=1,2,\cdots )$为对应的对数似然函数序列。
  （1）如果$P(Y|\theta )$有上界，则$L({\theta }^{(i)})=\log P(Y|{\theta }^{(i)})$收敛到某一值$L^{\ast }$；
  （2）在函数$Q(\theta ,{\theta }^{\prime })$与$L(\theta )$满足一定条件下，由EM算法得到的参数估计序列${\theta }^{(i)}$的收敛值${\theta }^{\ast }$是$L(\theta )$的稳定点。

第3步：证明在变分EM算法中，证据下界随着迭代进行单调递增

  对于证据下界

$$\begin{array}{}\text{(1)}& L(q,\theta )=E_q[\log p(x,z|\theta )]-E_q[\log q(z)]\end{array}$$

  假设迭代次数从$t$到$t-1$

1. E步：固定$\theta$，求$L(q,{\theta }^{(t-1)})$对$q$的最大化，由于 $\log p(x|\theta )$并不依赖于$q(z)$，所以 $L(q^{(t)},{\theta }^{(t-1)})$的最大值出现在KL散度等于0的时候，即$q(z)$与后验概率分布$q(z|x,{\theta }^{(t-1)})$相等的时候，此时有证据下界等于对数似然函数，即

$$\begin{array}{}\text{(2)}& \log p(x|{\theta }^{(t-1)})=L(q^{(t)},{\theta }^{(t-1)})\end{array}$$

  如下图所示

2. M步：固定$q(z)$求 $L(q,\theta )$对$\theta$的最大化，更新模型参数$\theta$，此时下界$L(q,\theta )$ 增⼤，有

$$\begin{array}{}\text{(3)}& L(q^{(t)},{\theta }^{(t-1)})⩽L(q^{(t)},{\theta }^{(t)})\end{array}$$

  由于 $q$ 由模型参数 ${\theta }^{(t-1)}$ 确定，并且在 M 步中固定，所以它不会等于更新新的后验概率分布 $p(z|x,{\theta }^{(t)})$， 因此 KL 散度⾮零，对数似然函数的增量⼤于证据下界的增量，即

$$\begin{array}{}\text{(4)}& L(q^{(t)},{\theta }^{(t)})⩽\log p(x|{\theta }^{(t)})\end{array}$$

  如下图所示

  结合公式(2)、(3)、(4)，可得

$$\log p(x|{\theta }^{(t-1)})⩽\log p(x|{\theta }^{(t)})$$

  即证得在变分EM算法中，$\log p(x|\theta )$每次迭代都保证单调递增，根据第2步的定理9.1和定理9.2，因此最后一定收敛。