非负矩阵分解

原论文 Algorithms for Non-negative Matrix Factorization


摘要

非负矩阵(NMF)对于多元数据的分解来说很有用,这里分析两个不同的NMF算法,其中一个最小化二乘误差,另一个最小化广义Kullback-Leibler分歧。两种算法的单调收敛可以使用类似于用于证明期望最大化算法的收敛的辅助功能来证明。算法也可以解释为对角重新缩放的梯度下降,其中最佳地选择重定标因子以确保收敛。

非负矩阵的分解

我们考虑以下形式化问题的求解: NMF 给定一个非负矩阵\(V\),找到非负矩阵\(W\)\(H\)使得有以下式子: \[V \approx WH \tag{1}\]

NMF可以被应用到多元数据的统计分析当中,给定一个n维的数据向量集,由一个 \(n*m\) 的矩阵\(V\)来表示,其中\(m\)表示数据集的例子个数。这个矩阵可以被分解成\(n*r\)的矩阵\(W\)以及\(r*m\)的矩阵\(H\),通常\(r\)都是要比\(n\)\(m\)要小的,所以\(W\)\(H\)也比原始矩阵要小一些。也就是说结果可以看做是原始数据矩阵的一种压缩。

那么公式(1)有什么意义呢?他可以携程列的形式\(v \approx Wh\)\(v\)\(h\)都是对应矩阵\(V\)\(H\)的某一列。也就是说,每一列的数据向量\(v\)都是对于矩阵\(W\)的一个线性变换,权重取值为\(h\)。因此\(W\)可以被视为包含针对V中的数据的线性近似而优化的基础。因为很多基本列都可以变换得到其他的数据向量,因此一个好的变换就是可以通过这些基本列来实现。

现在所关注的不是NMF的如何应用,而是关注与找到非负矩阵如何分解。当然,其他类型的矩阵分解已经在数值线性代数中被广泛研究了,但是由于非负的限制,使得之前的工作与现在的例子不完全兼容。

现在我们考虑两种基于不断迭代更新\(W\)\(H\)的NMF算法,因为这些算法容易实现且可以收敛,在实际情况中也可以应用。其他算法可能也有效但是比较难实现。

在我们算法的每一次迭代中,\(W\)\(H\)的新值可以通过公式(1)将当前值乘以质量因子而得。我们证明得出这种近似估计的质量随着应用的多次迭代,单调递增。也就是说这种重复迭代通过更新规则保证了矩阵分解的局部优化的收敛。

成本函数

为了找到近似分解\(V \approx WH\),我们首先需要确定成本函数来量化估计质量。成本函数可以通过构造两个非负矩阵\(A\)\(B\)的距离来确定。一个有效的测量就是简单计算两个矩阵的欧氏距离 \[ \left \| A-B \right \|^2=\sum_{ij}(A_{ij}-B_{ij})^2 \tag{2} \] 根据上式可以明显看出当\(A=B\)时,取下界0。另外一个有效的测量方法是 \[ D(A||B)=\sum_{ij}(A_{ij}log\frac{A_{ij}}{B_{ij}})-A_{ij}+B_{ij} \tag{3} \] 虽然上式的下界也是0,但是却不能称之为距离,因为并没有根据A和B对称,所以我们将其称之为A与B的区分度,它简化了Kullback-Leibler分歧,或者说是相对熵。当\(\sum_{ij}A_{ij}=\sum_{ij}B_{ij}=1\),A和B可以被看作是归一化的概率分布。 现在我们考虑两个可能的形式化NMF优化问题: Problem 1 最小化关于W和H的式子\(\left \|V-WH \right \|^2\),其中\(W,H\geq0\) Problem 2最小化关于W和H的式子\(D(V \|WH)\),其中\(W,H\geq0\)

尽管上式对于\(WH\)都分别是凸函数,但是二者一起来看就不是凸函数了,所以在全局对于\(WH\)找到一个最小值是不现实的,然而有很多技术可以进行数值优化找到局部最小值。 梯度下降法可能是最简单的一种方法,但是收敛速度可能比较慢,另外的如共轭梯度可以快速收敛找到一个局部最小值,但是实现复杂。除此之外,梯度下降法对于步长还比较敏感,对于大规模的应用来说适应性比较低。

乘法更新规则

我们发现"乘法更新规则"对于解决Problem 1&2来说有较高的效率提升且易于实现。

Theorem 1 欧氏距离\(\left \| V -WH \right \|\)在以下的更新规则中是非增长的当W和H达到平稳点的时候,欧氏距离不再发生变化。 \[ H_{a \mu} \leftarrow H_{a \mu} \\frac{(W^TV)\_{a\mu}}{(W^TWH)_{a\mu}} \]

\[ W_{ia} \leftarrow W_{ia} \frac{(VH^T)\_{ia}}{(WHH^T)_{ia}} \tag{4} \]

Theorem 2 \(D(V \|WH)\)的分歧子啊以下更新规则中不再增加。 \[ H_{a \mu} \leftarrow H_{a \mu} \frac{\sum_iW_{ia}V_{i\mu}/(WH)_{i\mu}}{\sum_kW_{ka}} \\\\ W_{ia} \leftarrow W_{ia} \frac{\sum_{\mu}H_{a\mu}V_{i\mu}/(WH)_{i\mu}}{\sum_vH_{av}} \tag{5} \] 当W和H达到平稳点时,区分度会保持不变

乘法与加法更新规则

拿乘法更新规则与梯度下降法相比较是很有用的。特别的有,一个对于H的简单加法更新可以减小平方距离,可以被写作: \[ H_{a\mu} \leftarrow H_{a\mu}+ \eta_{a\mu}[(W^TV)\_{a\mu}-(W^TWH)_{a\mu}] \tag{6} \]

如果\(\eta_{a\mu}\)都等于小的正整数,那么就是常规的梯度下降法,只要这个数足够小那么该式就可以规约成\(\left \| V-WH \right \|\) 如果我们规定变量

\[ \eta_{a\mu} = \frac{H_{a\mu}}{(W^TWH)_{a\mu}} \tag{7} \]

那么我们就得到了Theorem 1的更新规则。注意到重缩放的结果是由乘法因子与分母中的梯度的正分量以及分子中的负分量的绝对值来决定的。

至于区分度,斜向重缩放的梯度下降表示形式为 \[ H_{a\mu} \leftarrow H_{a\mu}+\eta_{a\mu}[\sum_iW_{ia}\frac{V_{i\mu}}{(WH)_{i\mu}}-\sum_iW_{ia}] \tag{8} \] 同样的,如果\(\eta_{a\mu}\)是小的正数,这个更新应该会减少\(D(V||WH)\),如果我们设定 \[ \eta_{a\mu}=\frac{H_{a\mu}}{\sum_iW_{ia}} \tag{9} \] 我们将会得到对于H的Theorem 2更新规则。这个重缩放能够被解释成一个在分母中有正因子梯度下降的乘法规则和分子中有负因子的乘法因子。 由于我们选择的\(\eta_{a\mu}\)并不非常小,所以它并不保证这样的梯度下降会造成成本函数的收敛。但是下一节将会对此证明。

收敛证明

为了证明\(Theorems 1,2\),我们会重复使用一个辅助函数来近似期望最大化的算法。

定理 1 \(G(h,h')\)是一个对于\(F(h)\)的辅助函数,当且仅当以下条件满足时 \[ G(h,h') \geq F(h),G(h,h)=F(h) \tag{10} \]

由于Lemma 1,辅助函数是一个非常有帮助的函数,在Fig.1中有所阐释。

如果G是一个辅助函数,那么F在更新过程中是一个非增长的。 \[ h^{t+1}=arg min_h\,G(h,h^t) \tag{11} \]

证明: \(F(h^{t+1})\leq G(h^{t+1},h^t)\leq G(h^t,h^t) = F(h^t)\)

注意到当\(h^t\)是函数\(G(h,h^t)\)的局部最小点的时候,有\(F(h^{t+1})=F(h^t)\)。如果\(F\)存在分歧且在\(h^t\)的很小一个邻域内连续,就是说\(\bigtriangledown F(h^t)=0\),因此,由(11)的更新规则,我们能得到一个序列,在局部对目标函数的最小值连续的收敛估计\(h_{min}=arg\,min_hF(h)\) \[ F(h_{min}) \leq ...F(h^{t+1}) \leq F(h^t) ... \leq F(h_2) \leq F(h_1) \leq F(h_0) \tag{12} \]

定理 2 如果\(K(h^t)\)是一个对角矩阵 \[ K_{ab}(h^t) = \delta_{ab}(W^TWh^t)_a/h_a^t/h_a^t \tag{13} \] 然后 \[ G(h,h^t) = F(h^t)+(h-h^t)^T\bigtriangledown F(h^t) + \frac{1}{2}(h-h^t)^TK(h^t)(h-h^t) \tag{14} \] 作为以下函数的辅助函数 \[ F(h)=\frac{1}{2}\sum_i(v_i-\sum_aW_{ia}h_a)^2 \tag{15} \] 证明有: 因为G(h,h) = F(h),我们需要证明\(G(h,h^t)\geq F(h)\) 我们比较下式与(14) \[ F(h) = F(h^t)+(h-h^t)^T\bigtriangledown F(h^t) + \frac{1}{2}(h-h^t)^T(W^TW)(h-h^t) \tag{16} \] 发现要满足\(G(h,h^t)\geq F(h)\)等价于 \[ 0 \leq (h-h^t)^T\[K(h^t)-W^TW\]\(h-h^t) \tag{17} \] 为了证明半定,考虑以下矩阵:

\[ M_{ab}(h^t) = h_a^t(K(h^t)-W^TW)_{ab}h_b^t \tag{18} \]

这样,当且仅当\(M\)是一个半正定矩阵时,\(K-W^TW\)也是一个半正定的矩阵。有:

\[ \begin{aligned} \upsilon^TM\upsilon = \frac{1}{2}\sum_{ab}(W^TW)_{ab}h_a^th_b^t(\upsilon_a - \upsilon_b)^2 \end{aligned} \] 上式大于0,我们可以证明Theorem1的收敛 Theorem 1证明 将公式(14)的结果代入公式(11)中得到更新规则: \[ h^{t+1} = h^t - K(h^t)^{-1} \bigtriangledown F(h^t) \tag{21} \] 因为公式(14)是一个辅助函数,F在此更新规则下是非增长的,根据Lemma 1,可以明确得到下式: \[ h_a^{t+1} = h_a^t\frac{W^T\upsilon_a}{W^TWh^t}_a \tag{22} \] 翻转W与H的角色,F也可以在W的更新规则下保持不增长。 现在我们考虑对于分歧开销函数的估计函数。 Lemma 3 定义 \[ G(h,h^t) = \sum_i(\upsilon_i log\upsilon_i - \upsilon_i)+\sum_{ia}W_{ia}h_a \\ -\sum_{ia}\upsilon_i \frac{W_{ia}h_a^t}{\sum_bW_{ib}h_b^t}(logW_{ia}h_a-log\frac{W_{ia}h_a^t}{\sum_b}W_{ib}h_b^t) \tag{23} \] 这是对于下式的辅助函数 \[ F(h) = \sum_i \upsilon_ilog(\frac{\upsilon_i}{\sum_aW_{ia}h_a}) - \upsilon_i + \sum_aW_{ia}h_a \tag{24} \] 证明 可以很直接证明\(G(h,h)=F(h)\),为了证明\(G(h,h) \geq F(h)\),有 \[ -log \sum_aW_{ia}h_a \leq -\sum_a\alpha_a log\frac{W_{ia}h_a}{\alpha_a} \tag{25} \] 总的来说这里的\(\alpha_a\)是非负的。设定有 \[ \alpha_a = \frac{W_{ia}h_a^t}{\sum_bW_{ib}h_b^t} \tag{26} \] 我们可以由此得到 \[ -log\sum_aW_{ia}h_a \leq - \sum_a\frac{W_{ia}h_a^t}{\sum_bW_{ib}h_b^t}(log W_{ia}h_a - log\frac{W_{ia}}{\sum_bW_{ib}h_b^t}) \tag{27} \] 从上述不等式可以看出\(F(h)\leq G(h,h^t)\)

Theorem 2证明 关于h的\(G(h,h^t)\)的最小值是由梯度下降为0时而决定的: \[ \frac{dF(h,h^t)}{dh_a}=-\sum_i \upsilon_i\frac{W_{ia}h_a^t}{\sum_bW_{ib}h_b^t}W_{ia} \tag{28} \] 因为G是一个辅助函数,在公式(24)里面,F在更新下是一个非增长的函数。重写该矩阵可以得到公式(5)。翻转H和W的角色,关于W的更新规则也可以看做是非增长的。