基于不动点方程的非负矩阵分解算法

刘伟强,林 鹭

(厦门大学数学科学学院,福建 厦门 361005)

非负矩阵分解(nonnegative matrix factorization,NMF)[1-2]是近年来提出的一种对大规模非负数据降维的方法,通过非负的限制,让分解得到的特征矩阵以及系数矩阵具有明显的可解释性,使得一些潜在的信息被有效地挖掘出来.NMF已经广泛地应用在人脸识别[3]、图像分析[4]、文本聚类[5]、盲源信号分类[6]、生物基因探测[7],以及生物数据挖掘等领域[8].

近年来,学者从不同角度出发,提出相应的NMF算法.Lee等[2]于2001年首次提出交替方向的乘性迭代算法(MU).Hoyer[9-10]从罚函数的角度出发,提出了带稀疏约束条件的NMF算法.Lin[11-12]则提出了基于梯度投影的NMF算法.但对NMF问题算法分析较多,理论分析较少.本文中从线性互补问题出发,得到基于不动点方程的NMF(fixed point equation NMF,FP-NMF)算法,并从理论上证明算法的收敛性.根据下降方向的不同,分别提出了不动点方程的最速下降(FP-steepest descent,FP-SD)算法与不动点方程的最小梯度(FP-minimal gradient,FP-MG)算法,证明了这两种算法的收敛性.

1 NMF问题与不动点方程的关系

s.t.W≥0,H≥0.

(1)

若将目标函数看成是W与H两个变量的函数,则式(1)变为一个非凸优化问题,但目标函数对于单变量W或者H而言均为凸函数.因此,常用交替迭代法求解.

固定W,此时的目标函数可改写为

其中Vj代表矩阵V的第j列,下同.式(1)等价于求解n个非负最小二乘问题[13],即

(2)

同理,对式(1)固定H,可转化为求解:

(3)

(4)

由定理1可知,利用交替方向迭代求解NMF问题时,其KKT条件正好为一个线性互补问题;又因问题(2)关于H是凸函数,故求解问题(2)等价于求解一个线性互补问题(4).

定理2[14]问题(4)与下列不动点方程同解:

Hj=(Hj-ηjgrad(Hj))+,ηj>0,

(5)

其中H+=max(H,0),下同.

类似地,对于问题(3),有以下结论:

(6)

定理4[14]式(6)与下列不动点方程同解:

(7)

2 FP-SD与FP-MG算法

用梯度下降法求解问题(2)与问题(3).固定W,按列计算新的H:

Hj∶=Hj-ηj(grad(H))j,

其中ηj是负梯度方向上的步长,j=1,2,…,n,grad(H)=WT(WH-V).

同理,固定H,按行计算新的W:

其中ξi是负梯度方向上的步长,i=1,2,…,m,grad(W)=(WH-V)HT.

分别记问题(2)与问题(3)的Hesse矩阵为G(H)与G(W),则

G(H)=WTW,G(W)=HHT.

现在考虑如何得到最优的步长ηj与ξi,使得算法具有更快的收敛速度.为此,用SD[15]优化以下式子:

j=1,2,…,n.

(8)

同理,固定H,定义:

(9)

算法1FP-SD

对k=0,1,…，

1) 计算(Wk)TWk的谱半径ρ((Wk)TWk);

2) 计算grad(Hk)=(Wk)T(WkHk-V),

G(Hk)=(Wk)TWk;

3) 对j=1,2,…,n,

4) 计算Hk(Hk)T的谱半径ρ(Hk(Hk)T);

5) 计算grad(Wk)=(WkHk-V)(Hk)T,

G(Wk)=Hk(Hk)T;

6) 对i=1,2,…,m,

注意到式(8)与(9)在选取步长时所用的方法为SD方法.当然,也可以考虑使得梯度值最小所对应的步长的MG方法[16],为此,考虑以下优化式子:

(10)

式(10)是关于步长ηj的一个非负最小二乘问题,解之可得:

同理,固定H,定义:

i=1,2,…,m.

(11)

式(11)是关于步长ξi的一个非负最小二乘问题,解之可得:

类似于Lin的方法[18-19]将式(10)与(11)应用于问题(1)中,得到如下基于投影梯度方法的MG算法.

算法2FP-MG

对k=0,1,…，

1) 计算(Wk)TWk的谱半径ρ((Wk)TWk);

2) 计算grad(Hk)=(Wk)T(WkHk-V),

G(Hk)=(Wk)TWk;

3) 对j=1,2,…,n,

5) 计算Hk(Hk)T的谱半径ρ(Hk(Hk)T);

6) 计算grad(Wk)=(WkHk-V)(Hk)T,

G(Wk)=Hk(Hk)T;

7) 对i=1,2,…,m,

3 FP-SD与FP-MG算法的收敛性

证明由已知可设WTW的特征值为0<λ1≤λ2≤…≤λr=ρ(WTW),则

又有

故有

再由引理1知:

表1 MU、FP-SD、FP-MG算法在不同秩上关于迭代次数、相对误差和运行时间的比较

Tab.1 Comparison of MU,FP-SD and FP-MG on the number of iterations,reletive error,and running time in different ranks

r迭代次数相对误差运行时间/s MUFP-SDFP-MGMUFP-SDFP-MGMUFP-SDFP-MG 53 0326979350.029 30.028 00.028 2671203302 104 2481 2462 6670.024 20.022 10.022 41 031375907 206 9891 1741 5790.020 50.018 00.018 71 905423613 308 5611 0651 4210.019 00.016 10.017 22 781385597 409181 5220.015 10.016 55751 023 605 2208411 4310.016 20.014 00.016 02 9925701 096 804 2767861 3550.015 20.013 40.015 72 643406840

4 数值实验

4.1 实验数据

选取剑桥大学的ORL人脸图库[20]做数值实验.ORL图库包含40个人脸的图片,每个人脸从不同角度拍摄10张相片,共400张图片,每张图片的像素为112×92,其元素值介于0到255之间.因此V为10 304×400矩阵.

4.2 KKT剩余量

W≥0,

grad(W)≥0,

W.*grad(W)=0,

和

H≥0,

grad(H)≥0,

H.*grad(H)=0.

上述条件可以等价简化为

min(W,grad(W))=0,

min(H,grad(H))=0，

(12)

其中式(12)左边为按元素取最小所得的矩阵.令

则式(12)成立的充分必要条件是δ=0.显然δ与W及H的阶数有关.为排除阶数的影响,实际计算中,取如下的标准KKT剩余量

作为度量,其中

δW=#(min(W,grad(W))≠0),

δH=#(min(H,grad(H))≠0).

4.3 实验结果与分析

表1是以KKT剩余量Δ≤0.01或者迭代步长上限10 000作为停机准则.对每一组r,所有的算法均取相同的初始矩阵对,迭代次数、相对误差及运行时间均为50组初始矩阵对计算后的平均值.

表中空白部分代表迭代次数超过上限10 000.观察表1可知:

1) 从迭代次数上看,FP-SD以及FP-MG算法均优于MU算法,且FP-SD算法表现最佳.当r取20,30,40,60时,MU算法未能在5 000次迭代内达到停机准则,FP-SD和FP-MG均在1 000步左右达到停机准则,并且FP-SD算法几乎全在1 000步内达到.当r为80时,FP-SD以及FP-MG算法的迭代次数均不足MU算法的1/3,其中FP-SD不足MU的1/5.除MU算法外,我们发现一个有趣的现象,即随着r的增长,迭代步数呈先增后减的倒U型趋势.

2) 从相对误差上看,随着r的增加,3种算法均呈现下降的趋势.且FP-SD以及FP-MG算法均较MU算法有较好的表现,其中相对误差平均下降了7%以上,有的下降了15%.

3) 从运行时间上看,随着r的增加,3种算法均呈现先增长后下降的趋势.且FP-SD以及FP-MG算法较MU算法都有明显的提高,运行效率大多提高了2倍以上,最多的提高了7倍.

下面固定r,对照MU,FP-SD及FP-MG算法的平均相对误差以及平均KKT剩余量.由于ORL图库包含40个人的人脸图片,所以取r=40,对比结果如图1.

图1 三种NMF算法的平均相对误差对比Fig.1 Comparison of the average relative error of the three NMF algorithms

图1对比了3种算法的平均相对误差.从图中可以看到,FP-SD与FP-MG算法的相对误差下降较快,且当迭代次数在20到40之间时,这2种算法较MU算法的优势更为显著.MU算法迭代近100次后,相对误差达到0.02,而FP-SD与FP-MG算法只需50次迭代,相对误差即可达到0.02,步长减少了一半.图1从侧面验证了所有算法的相对误差是单调下降的,这与前面的理论分析相吻合.平均KKT剩余量对比结果如图2.

图2 三种NMF算法的平均KKT剩余量对照Fig.2 Comparison of the average KKT residual of the three NMF algorithms

由图2可以看出,迭代15次后,FP-SD以及FP-MG算法的KKT剩余量均小于MU算法.图中还看到,要使KKT剩余量达到0.2,FP-SD以及FP-MG算法大约需要迭代40次,而MU迭代100次后仍未达到.进一步,FP-SD算法迭代70次时,KKT剩余量已低于0.1,而FP-MG需要迭代100次.总之,从KKT剩余量看,这3种算法中,FP-SD以及FP-MG算法效果较好.

5 结 论

本研究将运筹学中的线性互补问题与NMF问题联系起来,利用NMF问题与线性互补问题的关系,将NMF问题转化为求解不动点方程问题,从而设计出了基于不动点方程NMF算法,讨论了算法的收敛性,并从理论上给予证明.借助这个思路,利用最速下降法以及最小梯度法选取下降步长,分别得到基于不动点方程FP-SD与FP-MG算法.数值实验的结果表明,本文中提出的两种算法较MU算法在收敛速度上有明显的改进.