12. 隐变量线性模型

12.1 因子分析(Factor analysis)

混合模型的局限性在于它只能使用一个隐变量来生成观察值。特别地,每个观察值只能属于$K$个原型之一。我们可以认为混合模型中使用的隐变量是某个聚簇身份的one-hot编码,其中包含了$K$个二元隐变量。但因为这些变量之间是互斥的,所以模型的表达能力仍然是有限的。

一种替代方案是使用一个实数隐变量$\mathbf{z}_i\in \mathbb{R}^L$,并使用最简单的先验分布:高斯分布(我们会在下文考虑其他的先验分布),即: $$ p(\mathbf{z}_i)=\mathcal{N}(\mathbf{z}_i|\mathbf{\mu}_0,\mathbf{\Sigma}_0) \tag{12.1} $$ 如果观察值是连续值,即$\mathbf{x}_i \in \mathbb{R}^D$,我们或许可以使用一个高斯分布作为似然函数。就像在线性回归中的那样,假设期望是(隐藏)输入的线性函数,从而得到: $$ p(\mathbf{x}_i|\mathbf{z}_i,\mathbf{\theta})=\mathcal{N}(\mathbf{W}\mathbf{z}_i+\mathbf{\mu},\mathbf{\Psi}) \tag{12.2} $$ 其中$\mathbf{W}$是一个$D\times L$的矩阵,又被称为因子载入矩阵(factor loading matrix),$\mathbf{\Psi}$是一个$D \times D$的协方差矩阵。我们令$\mathbf{\Psi}$是一个对角矩阵,因为该模型的一个主要特点是“强迫”$\mathbf{z}_i$来解释变量之间的相关性,而不是将“这件事”交给观察值的协方差矩阵。整个模型被称为因子分析(factor analysis,FA)。一种特殊的情况是$\mathbf{\Psi}=\sigma^2\mathbf{I}$,此时模型被称为概率性主元分析(probabilistic principal components,PPCA),这个名字的来源我们将在后文介绍。

图12.1展示了当$L=1,D=2,\mathbf{\Sigma}_0$为球状协方差矩阵,$\mathbf{\Psi}$为对角矩阵时模型的示意图。我们使用一个各向同性的高斯分布作为$\mathbf{z}_i$的先验,然后将它沿着$1$维直线$\mathbf{w}z_i+\mathbf{\mu}$滑动(译者注:注意这里是滑动,不是采样得到单一的$z_i$值)。最终导致一个在$2$维空间中被拉长的(也就是相关的)高斯分布,该分布包含了沿对角线分布的噪声。 $$ figure 12.1 $$

12.1.1 FA是MVN的一个低秩版的参数化

FA可以认为是一种使用较少参数来描述变量$\mathbf{x}$的联合分布的方式。为了说明这一点,根据式4.126,推导出来的边缘分布$p(\mathbf{x}_i|\mathbf{\theta})$是一个高斯分布: $$ \begin{align} p(\mathbf{x}_i|\mathbf{\theta}) & =\int\mathcal{N}(\mathbf{x}_i|\mathbf{W}\mathbf{z}_i+\mathbf{\mu},\mathbf{\Psi})\mathcal{N}(\mathbf{z}_i|\mathbf{\mu}_0,\mathbf{\Sigma}_0)d\mathbf{z}_i \tag{12.3} \ & = \mathcal{N}(\mathbf{x}_i|\mathbf{W}\mathbf{\mu}_0+\mathbf{\mu},\mathbf{\Psi}+\mathbf{W}\mathbf{\Sigma}_0\mathbf{W}^T) \tag{12.4} \end{align} $$ 因此$\mathbb{E}[\mathbf{x}]=\mathbf{W}\mathbf{\mu}_0+\mathbf{\mu},\rm{cov}[\mathbf{x}]=\mathbf{W}\mathbb{E}[\mathbf{z}\mathbf{z}^T]\mathbf{W}^T+\mathbf{\Psi}=\mathbf{W}\mathbf{\Sigma}_0\mathbf{W}^T+\mathbf{\Psi}$。不失一般性,我们可以设置$\mathbf{\mu}_0=\mathbf{0}$,因为我们总可以使用$\mu$来包含$\mathbf{W}\mathbf{\mu}_0$的影响。类似地,不失一般性,我们可以令$\mathbf{\Sigma}_0=\mathbf{I}$,因为通过定义一个新的权重矩阵$\tilde{\mathbf{W}}=\mathbf{W}\mathbf{\Sigma}_0^{-\frac{1}{2}}$,我们总可以模拟出一个可以模拟相关性的高斯分布,因为: $$ \rm{cov}[\mathbf{x}]=\mathbf{W}\Sigma_0\mathbf{W}^T+\mathbf{\Psi}=\tilde{\mathbf{W}}\tilde{\mathbf{W}}^T+\mathbf{\Psi} \tag{12.5} $$ 因此我们发现FA使用一个低秩的分解,近似地模拟了观测向量的协方差矩阵: $$ \mathbf{C}\triangleq\rm{cov}[\mathbf{x}]=\mathbf{W}\mathbf{W}^T+\mathbf{\Psi} \tag{12.6} $$ 上式需要的参数数量的数量级为$O(LD)$,从而实现了在包含$O(D^2)$个参数的全协方差矩阵与包含$O(D)$个参数的对角协方差矩阵之间的灵活压缩。值得注意的是,如果我们不强制$\mathbf{\Psi}$为对角矩阵,反而将其设置为全协方差矩阵,这样的话就可以令$\mathbf{W}=\mathbf{0}$,此时,模型也就不再需要隐变量$\mathbf{z}$。

每个可观测变量的边缘方差定义为${\rm{var}}[xj]=\sum{k=1}^Lw_{jk}^2+\psi_j$,其中第一项是共享因子导致的方差,第二项$\psi_j$被称为唯一性(uniqueness),它特指那个维度的方差项。

12.1.2 隐藏因子的推理

尽管FA只是被认为是一种定义关于变量$\mathbf{x}$概率分布的方式,但我们经常使用它,是因为我们希望隐藏因子$\mathbf{z}$可以揭露一些关于数据的有趣的规律。为了实现这一点,我们需要计算关于隐藏因子的后验分布。基于贝叶斯法则,可以得出: $$ \begin{alignat}{3} p(\mathbf{z}_i|\mathbf{x}_i,\mathbf{\theta}) & = \mathcal{N}(\mathbf{z}_i|\mathbf{m}_i,\mathbf{\Sigma}_i) \tag{12.7} \ \mathbf{\Sigma}_i & \triangleq (\mathbf{\Sigma}_0^{-1}+\mathbf{W}^T\mathbf{\Psi}^{-1}\mathbf{W})^{-1} \tag{12.8} \ \mathbf{m}_i & \triangleq \mathbf{\Sigma}_i(\mathbf{W}^T\mathbf{\Psi}^{-1}(\mathbf{x}_i-\mathbf{\mu})+\mathbf{\Sigma}_0^{-1}\mathbf{\mu}_0) \tag{12.9} \end{alignat} $$ 值得注意的是,在FA模型中,$\mathbf{\Sigma}_i$事实上与$i$无关,所以我们可以将其表示为$\mathbf{\Sigma}$。计算这个矩阵所花费的时间复杂度为$O(L^3+L^2D)$,计算每个$\mathbf{m}_i=\mathbb{E}[\mathbf{z}_i|\mathbf{x}_i,\mathbf{\theta}]$的时间复杂度为$O(L^2+LD)$。$\mathbf{m}_i$有时又被称为隐藏分值(scores)或者隐藏因子(factors)

让我们给出一个简单的例子。考虑一个数据集,其中变量的数量$D=11$,样本数量$N=387$,这些变量描述了不同车辆的各种方面,比如发动机大小,气缸数量,每加仑公里路(miles per gallon,MPG),价格等等。我们首先训练一个$L=2$的模型。我们可以在二维空间中绘制$\mathbf{m}_i$来实现数据的可视化,结果如图12.2所示。

为了进一步加深对隐藏因子的理解,我们可以将每个特征维度对应的单位向量$\mathbf{e}_1=(1,0,…,0),\mathbf{e}_2=(0,1,0,…,0)$等投影到一个低维度的空间。这些结果如图中蓝色直线所示,这被称为双标图(biplot)。我们发现水平轴代表价格,对应于特征”经销商”和”零售”,且贵的车在右侧。垂直轴对应于燃油效率(以MPG进行衡量)与尺寸的比例:重型车辆效率更低且分布在上面,然而轻型车辆效率更高且分布在下面。我们可以通过点击某些点来验证这种解释,找到在训练集中距离最近的样例,然后打印出它们的名字。然而,通常情况下,对隐藏因子进行解释是非常困难的,这一点我们会在12.1.3节进行解释。

12.1.3 无法辨别性

FA模型的参数具有不可辨别性。为了说明这一点,假设$\mathbf{R}$是一个任意的正交旋转矩阵,满足$\mathbf{R}\mathbf{R}^T=\mathbf{I}$。定义$\tilde{\mathbf{W}}=\mathbf{WR}$,基于这个变体矩阵的似然函数与没变化之前的结果是一样的,因为 $$ \begin{align} {\rm{cov}}[\mathbf{x}] & = \tilde{\mathbf{W}}\mathbb{E}[\mathbf{z}\mathbf{z}^T]\tilde{\mathbf{W}}^T+\mathbb{E}[\mathbf{\epsilon\epsilon}^T] \tag{12.10} \ & = \mathbf{WR}\mathbf{R}^T\mathbf{W}^T+\mathbf{\Psi}=\mathbf{W}\mathbf{W}^T+\mathbf{\Psi} \tag{12.11} \end{align} $$ 几何上,使用一个正交矩阵与$\mathbf{W}$相乘等价于在生成$\mathbf{x}$之前旋转$\mathbf{z}$;但因为$\mathbf{z}$是从各向同性的高斯分布中采样得到的,所以导致最终的似然函数没有任何变化。因此,我们没有办法对$\mathbf{W}$进行唯一确定,也因此我们没有办法唯一确定隐藏因子。

为了确定一个唯一的解,考虑到正交矩阵$\mathbf{R}$的大小为$L×L$,所以我们需要移除$L(L-1)/2$个自由度。整体上,$FA$模型有$D+LD-L(L-1)/2$个自由参数(不包括期望),第一项为$\mathbf{\Psi}$的参数量。显然我们需要模型参数量小于或者等于$D(D+1)/2$,该值就是一个非受限(但是对称)的协方差矩阵的参数量。从而我们得到一个关于$L$的上限值: $$ L_{max}=\left \lfloor D+0.5(1-\sqrt{1+8D}) \right \rfloor \tag{12.12} $$ 举例来说,如果$D=6$,则$L\le3$。但通常情况下我们并不会选择这个上界,因为这会导致过拟合(12.3节将讨论如何选择$L$)。

不幸的是,即使我们令$L<L_{max}$,我们依然不能唯一的确定参数,因为旋转的不确定性依然存在。不可辨识性并不影响模型的预测性能。然而,它对载入矩阵存在影响,从而导致对潜在因子的解释存在不确定性。因为因子分析通常是用来揭示数据的结构,所以这个问题需要被解决。下面是一些常见的解决方案:

  • 强制要求$\mathbf{W}$为正交矩阵 可能解决不可辨识性问题的最直接的方案是强制$\mathbf{W}$为正交矩阵,然后按照对应的潜在因子的方差大小来对矩阵的列进行排序。这是PCA中所使用的方法,这一点我们将在12.2节讨论。其结果不一定更具备可解释性,但至少结果是唯一确定的。

  • 强制要求$\mathbf{W}$为下三角矩阵 一种在贝叶斯领域十分著名的方法是强制要求第一个观测变量只由第一个潜在因子产生。第二个观测变量仅由前两个潜在因子产生,以此类推。举例来说,如果$L=3,D=4$,对应的因子载入矩阵为: $$ \mathbf{W}= \begin{pmatrix} w{11} & 0 & 0 \ w{21} & w{22} & 0 \ w{31} & w{32} & w{33} \ w{41} & w{42} & w{43} \end{pmatrix} \tag{12.13} $$ 同时我们要求$w{jj}>0,j=1:L$。在这个受约束的矩阵中整体的参数量为$D+DL-L(L-1)/2$,与唯一可辨识的参数量相等。该方法的缺点在于,前$L$个被称为founder 变量的可观测变量,会影响潜在因子的解释,所以必须小心地进行选择。

  • 为权重增加稀疏控制先验 相较于提前指定矩阵$\mathbf{W}$中的哪些元素应该为0,我们可以使用$\mathcal{l}_1$正则,ARD或者spike-and-slab先验来鼓励矩阵中的某些元素为0。这被称为稀疏因子分析。这种方法不一定能确保一个唯一的MAP估计,但它鼓励具备可解释性的解。13.8节将介绍有关稀疏编码的细节。

  • 选择一个信息旋转矩阵 存在很多启发式的方式,来尝试去找到一个旋转矩阵$\mathbf{R}$,来修正矩阵$\mathbf{W}$以及潜在因子,从而增加可解释性,典型情况下是鼓励他们变得稀疏。一种著名的方式是$\mathbf{varimax}$(方差最大化旋转)。

  • 对潜在因子使用非高斯先验 如果我们将潜在因子的先验分布$p(\mathbf{z}_i)$替换成一个非高斯分布,有的时候我们可以唯一确定矩阵$\mathbf{W}$以及潜在因子。这个技术被称为独立因素分析或者ICA,这一点将在12.6节进行解释。

12.1.4 混合因子分析

FA模型假设数据服从一个低维的线性流型。在现实场景中,大部分数据更适合使用一些低维“曲线”流型来建模。我们可以使用多个局部线性流型来近似一个曲线流型。这将导致如下的模型:令第$k$个线性子空间的维度为$L_k$,由$\mathbf{W}_k$表示。其中$k=1:K$。假设我们有一个潜在因子指示变量$q_i \in{1,…,K}$,用来指定使用哪一个子空间来生成数据。我们从一个高斯先验分布中采样$\mathbf{z}_i$,并将它传递给矩阵$\mathbf{W}_k(k=q_i)$,并且添加噪声。更加精确的模型表示为: $$ \begin{align} p(\mathbf{x}_i|\mathbf{z}_i,q_i=k,\mathbf{\theta}) & = \mathcal{N}(\mathbf{x}_i|\mathbf{\mu}_k+\mathbf{W}_k\mathbf{z}_i,\mathbf{\Psi}) \tag{12.14} \ p(\mathbf{z}_i|\mathbf{\theta}) & = \mathcal{N}(\mathbf{z}_i|\mathbf{0},\mathbf{I}) \tag{12.15} \ p(q_i|\mathbf{\theta}) & = {\rm{Cat}}(q_i|\mathbf{\pi}) \tag{12.16} \end{align} $$ 这被称为混合因子分析(mixture of factor analysers,MFA)。图12.3表示了其中的条件独立假设。

我们也可以将这个模型看做是混合高斯模型的低秩版本,特别地,该模型的空间复杂度为$O(KLD)$,而非$O(KD^2)$,后者对应全协方差矩阵。这将有利于减少过拟合。事实上,对于高维的实数变量数据而言,MFA是一个很好的通用密度模型。

12.1.5 因子分析模型的EM算法

使用第4章的结果,我们可以直接推导出训练一个FA模型所使用的EM算法。只需要一些额外的工作,我们就可以训练一个混合FA模型。接下来,我们将直接展示结果,而不给出证明。

为了获得单个因子分析器的结果,只需要在下面的公式中设置$r_{ic}=1$和$c=1$。在12.2.5节,我们将看到这些公式的进一步的简化版,这个版本就是训练PPCA模型的算法,其结果将具备一个特别简单而且优雅的解释。

回到最一般的情况:在步骤$\rm{E}$中,我们计算数据点$i$属于聚簇$c$的后验概率: $$ r{ic}\triangleq p(q_i=c|\mathbf{x}_i,\mathbf{\mathbf{\theta}}) = \pi_c\mathcal{N}(\mathbf{x}_i|\mathbf{\mu}_c,\mathbf{W}_c\mathbf{W}_c^T+\mathbf{\Psi}) \tag{12.17} $$ 关于$\mathbf{z}_i$的条件后验分布为: $$ \begin{align} p(\mathbf{z}_i|\mathbf{x}_i,q_i=c,\mathbf{\theta}) & = \mathcal{N}(\mathbf{z}_i|\mathbf{m}{ic},\mathbf{\Sigma}{ic}) \tag{12.18} \ \mathbf{\Sigma}{ic} &=(\mathbf{I}L+\mathbf{W}_c^T\mathbf{\Psi}_c^{-1}\mathbf{W}_c)^{-1} \tag{12.19} \ \mathbf{m}{ic} &= \mathbf{\Sigma}_{ic}(\mathbf{W}_c^T\mathbf{\Psi}_c^{-1}(\mathbf{x}_i-\mathbf{\mu}_c)) \tag{12.20} \end{align} $$ 待定

12.1.6 含缺失值的FA模型训练

12.2 主元分析(Principal components analysis,PCA)

考虑在FA模型中添加约束:$\mathbf{\Psi}=\sigma^2\mathbf{I}$,$\mathbf{W}$为正交矩阵。结果表明,当$\sigma^2\rightarrow 0$时,该模型将退化为典型的(非概率性的)主元分析(PCA),又被称为Karhunen Loeve变换。如果$\sigma^2 \gt 0$,则该模型被称为概率性PCA(probabilistic PCA)或者敏感性PCA(sensible PCA)

为了使上述结果更加直观,我们首先需要了解经典的PCA。接着将PCA与SVD相联系。最后我们回到对PPCA的讨论。

12.2.1 经典PCA:定理阐述

经典PCA的综合观点可以总结为如下理论。

定理12.2.1. 假设我们希望找到一个由$L$个向量$\mathbf{w}_j\in\mathbb{R}^D$组成的正交基,以及每个样本所对应的分值$\mathbf{z}_i\in \mathbb{R}^L$,从而使得平均重构误差(reconstruction error)达到最小 $$ J(\mathbf{W},\mathbf{Z})=\frac{1}{N}||\mathbf{x}_i-\hat{\mathbf{x}}_i||^2 \tag{12.26} $$ 其中$\hat{\mathbf{x}}_i=\mathbf{W}\mathbf{z}_i$,$\mathbf{W}$为正交矩阵。等价地,我们可以将上述目标写成如下形式: $$ J(\mathbf{W},\mathbf{Z})=||\mathbf{X}-\mathbf{WZ}^T||^2_F \tag{12.27} $$ 其中$\mathbf{Z}$为一个$N\times L$的矩阵,$\mathbf{z}_i$为它的每一行。$||\mathbf{A}||_F$为矩阵$\mathbf{A}$的Frobenius范数。