基于聚类优化的非负矩阵分解方法及其应用

栗茂林　梁　霖　陈元明　徐光华　何康康

1.西安交通大学工程坊，西安，7100492.西安交通大学机械工程学院，西安，710049

0　引言

随着信息获取技术的不断进步，描述系统状态的信息量越来越大，导致特征空间维数不断增加，从而引发维数灾难等问题。在传统的特征约简方法中，主分量分析、独立成分分析以及矢量量化等方法要求信号平稳、满足高斯条件，限制了应用范围。

作为一种新兴的多元数据处理方法，非负矩阵分解(non-negative matrix factorization，NMF)可在高维空间中获得原始数据的局部特征，其纯加性的表达方式符合“局部构成整体”的认知规律，成为信号处理、模式识别等领域的热门工具。张培林等[1]采用NMF技术对发动机故障信号进行特征参数提取，获得了更高的分类精度。李兵等[2]采用二维NMF技术来提取时频分布矩阵特征参数，得到了较好的诊断效果。除了特征提取外，NMF在复杂工业过程和机电设备的监测诊断中也有应用[3-4]。

根据NMF的问题模型，可将其求解归结为一个优化问题，即通过目标函数来刻画它的逼近程度。实际应用中，一般采用交替迭代方法获得局部最优解。常用的乘性迭代算法、梯度下降算法和交替最小二乘算法是基于不同思路提出的，而对于设备诊断的特征提取，需要选择出合适的分解模型。

研究表明，NMF进行数据约简的目的是估计出原数据中的结构，而其分解的基向量与K均值聚类中的“类”有相通之处，即具有软聚类特性。因此，本文从分类性能和迭代效率角度出发，基于NMF的聚类特性，将交替非负最小二乘算法用于故障诊断，并通过基矩阵的聚类性优化出约束参数与嵌入维数，最后，通过测试数据和特征选择实例应用验证了其有效性。

1　非负矩阵分解的迭代准则

1.1　NMF原理

NMF是一种非负性约束下的矩阵分解方法，具体可描述如下：给定一个非负矩阵V∈Rm×n(m为样本数，n为特征个数)，将矩阵V分解成非负的基矩阵W∈Rm×r和系数矩阵H∈Rr×n(通常情况下，要求r

非负矩阵分解可以通过优化迭代来解决，即通过目标函数来刻画逼近程度，并在非负约束条件下进行迭代求解。距离目标函数中，基于欧氏距离目标函数的优化问题应用较广泛[5]，其表达式如下：

(1)

1.2　迭代准则

在目标函数(式(1))中，当W和H同时为变量时，求解最小化问题是非凸的，因此采用交替迭代更新W和H以获得局部解。常用的迭代算法有乘性迭代(multiplicative iterative，MI)、Lin’s投影梯度(Lin’s projected gradient，LPG)及交替最小二乘法(alternating least squares，ALS)。其中，MI算法是针对欧氏距离目标函数的优化问题提出的；LPG算法用投影梯度法进行更新迭代，采用Armijo规则来搜索每次迭代步长；ALS算法根据库恩塔克一阶最优性条件，采用直接估计“驻点”的固定点方法获得局部解[6]。在这三种典型迭代算法中，ALS算法形式更简单，理论上分解结果优于MI算法，在处理高维数据时的效率高于LPG算法，缺点是对噪声敏感，但可以通过增加约束条件来提高优化效果。标准ALS算法的更新规则为

(2)

式中，[·]+表示非负矩阵。

2　面向聚类的非负矩阵分解算法

2.1　NMF的聚类特性

NMF数据约简的目的就是估计出原始数据的本质结构，即获得基矩阵W。事实上，基矩阵W也可理解为模式识别理论中的“模式”，基向量空间就形成了K均值聚类中的“类”[7]。因此，NMF就是要挖掘出原始特征中的本质结构，发现其中的“类”，将相似的局部特征聚集为一类，然后用提取出的“类”代表原始特征集。

NMF分解过程中没有正交性约束的要求，即具有软聚类性，有助于克服一些硬聚类所遇到的问题。同时，NMF不受样本数据的限制，无需预先确定聚类数，有效提高了其应用能力。

2.2　施加约束的改进型交替最小二乘法

为了获得更稳定、优质的分解结果，需要对W和H施加其他约束条件，可同时对W和H或只对其中一个施加约束。从故障诊断的角度来说，希望基矩阵W的聚类效果尽可能好，相关性尽可能小，也希望权矩阵H具有一定的稀疏度。因此，对基矩阵W施加相关性约束，对权矩阵H施加稀疏性约束，最优化问题转化为[8]

(3)

式中，αW为相关系数，αW≥0；JW(W)为W的迹，JW(W)=tr(WBWT)；矩阵B为k阶全1矩阵；αH为稀疏化系数；JH(H)为H的迹，JH(H)=tr(HBHT)。

约束条件下，改进型ALS的迭代更新准则为

W←[VHT(HHT+αWB)-1]+

(4)

H←[(WTW)-1(WTV-αHB)]+

(5)

诊断应用中，约束条件可根据情况进行选择，可同时施加两个约束，也可只施加其一。约束参数的选取可通过基矩阵W的聚类效果进行确定。

2.3　初始化方法的选择

迭代更新方式对W、H的初始值很敏感，常用的初始化方法[9-11]中，随机初始化法、随机Acol初始化法以及模糊C均值初始化法的初始值不唯一，分解结果不稳定，而基于奇异值分解的初始化方法可以获得唯一的初始值，使得分解结果稳定，具有可复现性。

2.4　嵌入维数的选择

3　实验分析

3.1　测试数据集

选取4个常用的UCI数据集(电离层数据集、垃圾邮件分类数据集、钢板缺陷数据集以及肺癌数据集)对迭代算法的分解结果进行比较分析，详情见表1。

表1　测试数据集的描述Tab.1　Description of test datasets

3.2　迭代算法的效果对比

为了验证标准ALS的效果，将其与MI和LPG算法进行对比，其中，嵌入维数r的选取见表1。由于随机初始化的多次处理有可能包含最优解，因此，3种算法均采用随机初值，统计10次的迭代结果。图1所示为不同迭代算法在4个数据集上应用的分类正确率。

由图1可知，MI算法的分类正确率总体上落后于其他2种算法，仅仅是针对肺癌数据集的分类率优于LPG算法；MI算法不稳定且易收敛至局部点，如图1a中第10次和图1c中第3次实验结果，分类正确率都处于最低点。LPG算法在三者中的稳定性最好，波动较小。标准ALS算法可以获得比MI算法和LPG算法更优的解，如图1a中第2次和图1d中的所有结果，但波动性比较大，主要是对噪声比较敏感。此外，即使是同一种迭代方法，每次实验结果都有较大的差距，说明随机初始值对NMF影响较大。标准ALS算法在分类性能方面表现良好，而通过施加约束条件还有提高分解效果的可能。另外，从数据集的分解效果来看，垃圾邮件分类数据的3种分解效果基本相同，这与垃圾邮件分类较明确有关。

3.3　初始化方法的影响

采用电离层和垃圾邮件分类数据集对初始化方法进行测试，目标函数为欧氏距离，迭代规则为标准ALS，嵌入维数r分别为5和4，评判标准不变。不同的初始化方法得到的基矩阵W的分类结果如图2所示。

(a)电离层数据集的分类正确率

(b)垃圾邮件分类数据集的分类正确率

(c)钢板缺陷数据集的分类正确率

(d)肺癌数据集的分类正确率图1　三种迭代算法的分类正确率Fig.1　Accuracy of three iterative algorithms

对于电离层和垃圾邮件分类两种数据集来说，最优和最差的分解结果皆由随机法和随机Acol法获得，这两种方法的波动性也最大；模糊C均值法与奇异值分解法虽然最终的分类率不如上述两种方法，但其波动性很小，且与最优值差距不大。这四种方法当中，奇异值分解法的结果最稳定，平均分类率也较高，且具有可重复性，是最为理想的初始化方法。

3.4　约束条件的影响

(a)电离层数据集

(b)垃圾邮件分类数据集图2　两种数据集的初始化方法对比结果Fig.2　Comparison with the initial methods for two datasets

图3　相关性约束对电离层数据集的影响Fig.3　Effects of the Ionosphere dataset with the correlation constraints

相关性约束可以减少基矩阵W中的冗余结构，有利于控制基向量的冗余性。图3所示为相关性约束对电离层数据集的影响结果，其中，初始化方法和评判标准不变。显然，施加约束下的分类正确率略有提高。当分类结果较差(图3中第10次实验)时，施加相关性约束后的效果较为明显。当分类结果已较为理想(图3中第7次实验)时，约束效果反而变差，究其原因，主要是当已较为准确地提取出数据的本质结构时，施加约束反而无法迭代出更好的结果，为此，在后续应用中不施加相关性约束条件。

特征约简的重点是获取稀疏的权矩阵，即希望用极少维度的向量表征数据的突出特征。因此，对钢板缺陷数据集进行稀疏性约束，其实验结果如图4所示。由图4可知，在10次随机初始优化搜索中，在施加稀疏性约束的情况下，NMF分解得到基矩阵的聚类效果基本都得到了改善，即分类正确率得到了提高。这表明了稀疏约束对提高特征约简效果的有效性。

图4　稀疏性约束对钢板缺陷数据集的影响Fig.4　Effects of the Steel plates faults dataset with the sparsity constraints

图5　不同稀疏度参数对电离层数据集的影响Fig.5　Effects of the Ionosphere dataset with the different sparsitys

稀疏度约束能够提高分解质量，但稀疏度系数的不同对分类效果也有一定的影响。图5所示为电离层数据集在不同稀疏度系数下的分类正确率变化情况，其中，嵌入维数r=5，稀疏度系数在0～0.5内变化，KNN分类器的邻域为5。曲线变化情况表明，当αH=0.4时，分解得到的基矩阵W的聚类效果最好，对应着较大的分类正确率。参数选择不佳，则导致约束下的分类效果变差，因此在应用中需要根据分类情况来优选稀疏度系数。

4　在特征选择中的应用

4.1　实例对象数据

采用由伊斯曼化学公司创建的田纳西-伊斯曼过程(Tennessee-Eastman process，TEP)进行实例分析。TEP过程包含了1种正常状态IDV(0)和21种故障状态(IDV(1)～IDV(21))，涉及22个过程测量变量(XMEAS(1)～XMEAS(22))、19个成分测量值(XMEAS(23)～XMEAS(41))和12个控制变量(XMV(1)～XMV(12))。22个过程测量变量包括进料、反应器压力以及汽提器的流量等信息；19个成分测量值为反应过程中各种气体成分测量值，都包含高斯噪声；12个控制变量描述了进料量、分离器罐液流量以及冷却水流量等信息。其中，搅拌速度变量恒定不变，在应用中不包含该项，即一共52个观测变量。

从TEP中选择IDV(2)、IDV(3)、IDV(4)、IDV(5)四类故障数据(每类故障480个样本)的52个特征进行特征子集的选择。因此，待分解数据为1920×52的高维矩阵V。

4.2　约束参数和嵌入维数的选择

特征选择的目的是从原始52维特征中选出对故障分类敏感的特征。因此，首先在原始特征中，通过特征间的相关系数剔除掉矩阵V中相关性较高的冗余特征，将矩阵特征维数从52维降到16维，则待分解矩阵X为1 920×16维矩阵。其次，评估不同参数下分解矩阵X得到的基矩阵W的分类性能，将其中分类性能好的基矩阵作为特征选择的基础。在评估中，低维嵌入维数r在3～5内变化，稀疏度约束αH以0.05的步长在0～0.5内变化，通过分解基矩阵的聚类情况确定出优化的维数和稀疏度。

图6所示为不同参数下基矩阵的分类率曲线，其中，KNN分类器的近邻数为5。由图6可知，在r=5，αH=0.5时分解得到的基矩阵W的分类性能最好。根据各基向量在分类性能上的互补性和组合性，剔除掉基矩阵W中的第2和第4冗余基向量，保留第1、3、5基向量，使得彼此间具有结构上的互补性，组合在一起则具有良好的分类性，其投影空间中的样本分布如图7所示。

图6　不同参数下基矩阵的分类正确率Fig.6　Accuracy curve of basis matrix with different parameters

4.3　特征选择的效果

根据优化参数确定的权矩阵H分布能够进行特征选择，即在最优的基矩阵W对应的系数矩阵H中找出对应行幅值最大的元素，该元素所在的列即为有效特征。

在优化参数确定的系数矩阵H中，找出对应行幅值最大的元素，该元素所在的列依次为特征16、15和2，对应原始集合特征子集为{52，51，10}。其中，特征52、51和10分别为控制变量XMV(11)、XMV(10)和过程测量变量XMEAS(10)。该子集的KNN分类率达到89.74%。

(a)第1基向量空间中的样本分布

(b)第3基向量空间中的样本分布

(c)第5基向量空间中的样本分布图7　三个基向量空间中的样本分布Fig.7　Distribution of three projection vectors

图8　子集{52，51，10}的三维样本分布Fig.8　Distribution of the feature subset{52，51，10}

图8所示为{52，51，10}所张成的三维空间投影，由样本分布可知，除IDV(3)故障样本(“*”)和IDV(5)故障样本(“□”)之间有一定的重合外，其他故障样本间都区分得较为明显，显然是一个较理想的特征子集。

若不考虑稀疏度的约束条件，在r=5的情况下，稀疏度和相关系数均取零，选出的特征依次为12、14和15，对应原始集合中的分类特征子集为{41，51，52}，该子集的KNN分类率只能达到81.87%。由此可见，不考虑约束条件的话，识别效果会有一定程度的下降。

基于主成分分析的特征选择结果为{10，28，47}，样本分布的正确率仅为65.62%，显然提取的特征子集更不理想。

5　结语

结合NMF的聚类性能，在NMF的迭代求解优化问题基础上，结合相关性约束和稀疏性约束，提出了一种面向NMF聚类性的迭代优化算法，通过基矩阵空间中的分类性能够选择出具有最优聚类性的嵌入维数及相关参数。通过测试数据验证了最小二乘迭代算法的有效性以及约束条件的效果，通过复杂机电系统的应用实例验证了特征选择的有效性。

参考文献：

[1]张培林，王怀光，张磊，等.非负矩阵分解在发动机故障特征提取中的应用[J].振动工程学报，2013，26(6)：944-950.

ZHANG Peilin，WANG Huaiguang，ZHANG Lei，et al. Feature Extraction for Engine Fault Diagnosis by Utilizing Adaptive Multi-scale Morphological Gradient and Non-negative Matrix Factorization[J]. Journal of Vibration Engineering，2013，26(6)：944-950.

[2]李兵，米双山，刘鹏远，等.二维非负矩阵分解在齿轮故障诊断中的应用[J].振动测试与诊断，2012，32(5)：836-840.

LI Bing，MI Shuangshan，LIU Pengyuan，et al. Application of Two-dimensional Non-negative Factorization for Gear Fault Diagnosis[J]. Journal of Viaration, Measurement & Diagnosis，2012，32(5)：836-840.

[3]陈霖，邹金慧.基于NMF-SVM的复杂化工过程故障诊断[J].计算机与应用化学，2014，31(8)：1015-1018.

CHEN Lin，ZOU Jinhui.Fault Diagnosis of Complex Chemical Process Based on NMF-SVM [J].Computers and Applied Chemistry，2014，31(8)：1015-1018.

[4]唐曦凌，梁霖，高慧中，等.结合连续小波变换和多约束非负矩阵分解的故障特征提取方法[J].振动与冲击，2013，32(19)：7-11.

TANG Xiling，LIANG Lin，GAO Huizhong，et al.Fault Feature Extraction Method Combining Continuous Wavelet Transformation with Multi-consitraint Nonnegative Matrix Factorization[J]. Journal of Vibration and Shock，2013，32(19)：7-11.

[5]WANG Yuxiong，ZHANG Yujin.Nonnegative Matrix Factorization: a Comprehensive Review[J]. IEEE Transactions on Knowledge and Data Engineering，2013，25(6)：1336-1353.

[6]BERRY M W，BROWNE M，LANGVILLE A N，et al. Algorithms and Applications for Approximate Nonnegative Matrix Factorization[J]. Computational Statistics & Data Analysis，2007，52(1)：155-173.

[7]LI T，DING C.The Relationships among Various Nonnegative Matrix Factorization Methods for Clustering[C]//Sixth International Conference on Data Mining. Hong Kong，2006：362-371.

[8]CICHOCKI A，ZDUNEK R，PHAN A H，et al. Nonnegative Matrix and Tensor Factorizations: Applications to Exploratory Multi-way Data Analysis and Blind Source Separation [M]. Hoboken：Wiley，2009：203-213.

[9]YU Shaohui，ZHANG Yujun，LIU Wenqing，et al. A Novel Initialization Method for Nonnegative Matrix Factorization and Its Application in Component Recognition with Three-dimensional Fluorescence Spectra[J]. Spectrochimica Acta Part A: Molecular and Biomolecular Spectroscopy，2012，86：315-319.

[10]XUE Yun，TONG C S，CHEN Ying，et al. Clustering-based Initialization for Non-negative Matrix Factorization[J]. Applied Mathematics and Computation，2008，205 (2)：525-536.

[11]JANECEK A，TAN Y. Using Population Based Algorithms for Initializing Nonnegative Matrix Factorization[J]. Lecture Notes in Computer Science, 2011，6729：307-316.

(编辑张洋)