EEMD-四线性平行因子的机械多故障欠定盲分离方法

EEMD-四线性平行因子的机械多故障欠定盲分离方法

朱亚静，李泽东，李志农,谷士鹏，马亚平（1.南昌航空大学无损检测教育部重点实验室，南昌330063；2.中国飞行试验研究院，西安710089）平行因子（ParallelFactor,PARAFAC）在宽松条件下具有分解唯一性的优势，在机械故障诊断中获得了成功应用。李晓明等[1]将传统奇异谱分解结合标准PARAFAC张量分解方法进行扩展，对机械故障特征进行提取。苗育茁等[2]利用连续小波变换构造3阶张量，通过PARAFAC对其进行分解，提出了一种多尺度平行因子分析算法，可以有效地从机械非线性多故障状态中提取特征，进行故障诊断。杨诚等[3]结合Volterra模型与PARAFAC构建了故障预测模型，主要解决了PARAFAC在分解机械非线性系统信号时估计参数较多的问题。李志农等[4]结合变分模态分解(VariationalModeDecomposition，VMD)和PARAFAC各自优点，提出了基于VMDPARAFAC欠定盲分离方法，并将其应用到滚动轴承的复合故障特征提取中。孙宜权等[5]采用自适应PARAFAC和BSS算法对柴油机喷油多故障进行诊断。然而，上述方法中的PARAFAC分解都是通过构造3阶张量完成的，采用三线性PARAFAC进行机械故障诊断时，只能对振动信号的部分信号建模，存在信息不完整的问题。在平行因子的发展过程中，四线性PARAFAC的提出具有重要的意义。四线性PARAFAC具有三线性平行因子所不具有的优势，且采用四线性平行因子盲分离算法可有效地对机械多故障进行源数估计，更好地分解信号，提取故障特征[6]。然而，在研究机械多故障欠定盲分离时，通常情况下对四线性PARAFAC模型采用QALS迭代拟合，而在欠定条件下，收敛过程中会出现陷入局部收敛的问题，得到的估计源信号往往不是很理想。集成经验模态分解(EnsembleEmpiricalModeDecomposition，EEMD)是针对经验模态分解的不足而提出的1种噪声辅助数据分析方法，该方法充分利用白噪声均值为零的特点，使得由于干扰导致的任何尺度的信号所出现间断，都能得到白噪声的补充，使间断信号在相应尺度形成连续信号，该方法很好抑制了由于EMD（EmpiricalModeDecomposition，EMD）分离所导致各分量之间出现混频现象[7]。文献[8]结合EEMD和具有自适应增长特性的输出隐藏反馈Elman网络对滚动轴承进行诊断。文献[9]将EEMD应用于机械故障中的复杂信号处理，提取故障特征。文献[10]提出1种改进的EEMD方法，以消除在EEMD方法中只添加噪声而导致的局限性，结果表明所提方法可以抑制机械信号模态混淆并且在故障诊断方面效果明显。文献[11]将EEMD和散布熵结合，构建故障信号的高维特征，在结合LPP-KNN进行故障诊断。文献[12]将集合经验模态分解(EEMD)和低相干K-SVD相结合进行齿轮故障特征提取。上述研究均表明EEMD在进行机械振动信号分解时具有独特的优势。在此，结合四线性平行因子盲分离算法和集合经验模态分解算法(EEMD)的的独特优势，提出了一种基于EEMD分解的四线性PARAFAC欠定盲源分离方法。该方法首先利用EEMD分解传感器所采集的观测信号，接着从分解得到的IMF分量中选出恰当相关IMF分量信号。根据选出的相关IMF分量与原观测信号重新构造信号，使得经四线性分解处理的重组信号数大于或等于源信号数，这样盲分离中的欠定就被转换为正定或超定。最后利用四线性交替最小二乘法迭代(QuadrilinearAlternatingLeastSquares，QALS)进行拟合迭代，得到混合矩阵估计，用最短路径法得到源信号估计。然后，通过仿真和实验研究以验证该方法的有效性。1EEMD-四线性平行因子盲分离方法信号线性混叠模式可以表示为如下形式：式中：S(t)=[s1(t),…,sM(t)]T为源信号矩阵，每个源信号包含L个数据点，A(t)是N×M维的混合矩阵(N<p>在观测信号X(t)中加入白噪声信号，利用EEMD方法将含噪的观测信号分解成多个IMF信号分量。具体步骤如下：第一步：首先在处理观测信号X(t)=[x1(t),…,xN(t)]T前设定平均处理次数M，其中i=1,2…,M；第二步：接着在观测信号X(t)添加白噪声ni(t)，将加入白噪声的信号作为一个整体，然后与原观测信号组成新的观测信号Xi(t)=X(t)+ni(t)；第三步：用EMD分解Xi(t)得到IMF，得到子分量ci,n(t)；其中：ri,n(t)是得到的残余分量；第四步：重复前面3个步骤，通过分解得到新的一系列的IMF；第五步：对得到的IMF分量进行M次平均获得cn(t)；其中：i=1,2…,M；n=1,2…M。第六步：利用相似系数公式求通过分解得到的IMF分量相关系数，然后选取相关性较大的有效IMF。选取相关性较强的IMF信号分量，并将其与原来的传感器数据组合，构建新的观测信号X(t)。这样，实现了将欠定条件下的盲源分离转化为正定或超定条件下盲源分离。然后将经中心化处理后的观测信号进行分段处理，其被均匀划分成J个不会重叠的数据段，经计算可知每段的数据点数为L/J个，记其为K，即K=L/J，具体段数可以任意选择，若最后一段的点数不够，可以补零。每一个数据段对应的采集时间为h，则其共被分为t/h个时间段，记其为I，此时：信号块的时滞协方差矩阵为：其中：k=1,2,…,K，序号tk表示t内第k个数据点。传感器在h段采集到的各源信号数据时滞协方差Rn可以通过计算得到，有：式(6)中：tk表示t内第k个源信号。将块时滞协方差矩阵转换成为对角阵，并且使该矩阵等于矩阵B，将源时滞协方差矩阵转换成为对角矩阵并且使该矩阵等于矩阵C，则式(5)、式(6)可以分别改写为式(7)和式(8)。将源时滞协方差矩阵以及块时滞协方差矩阵叠加，叠加后所得的3阶矩阵计为RJ,K其中：n=1,2,…,N,k=1,2,…,K。相邻对角阵相乘可以换位，则式(9)可表示为：将式中N×K个切片累积成N×J×K×I维的R（四维数据集）。其中：n=1,2,…,N,k=1,2…,K，j=1,2,…,J，i=1,2,…I，式(11)即为所构造的EEMD-四线性盲分离模型。如上实现了集合经验模态四线性平行因子盲源分离算法的构建，其分解过程中采用QALS得到混合矩阵估计A′。分解算法的主要思想为让上一次估计更新的矩阵参与下一次的矩阵求解，同时使根据最终估计得到的4个载荷矩阵与被分解矩阵的差值尽量最小，迭代过程步骤分为4步。第二步：将A重新代入式(12)，接着最小化式(12)得到估计的，令其等于B，更新矩阵B。第三步：将求得的矩阵B重新代入式(12)，接着最小化式(12)得到矩阵估计，令其等于C，更新矩阵C。第四步：将按照以上循环步骤得到的载荷矩阵代入式(17)，若满足式(17)，则退出，反之继续。一般认为ε为1×10-6，ςς3是残余误差的平方和，表示含噪的切片，‖∙‖F为Frobenius范数，i=1,2,…,I，j=1,2,…,J，上标+为Moore-Penrose伪逆，和表示A、B和C的上一次估计，m为迭代次数。当式（17）所示的停止准则得到满足即表示矩阵更新结束算法收敛，得到最终的混合矩阵估计A′，令A′=A。令ai为A第i个列向量，si(t)表示源信号的第i个源向量。线性混叠模型X(t)=AS(t)可以表示为：其中：X(t)是观测信号，它可以表示成a1s1(t)、a2s2(t)、…、ansn(t)的累加。ai表示矢量方向，si(t)表示向量长度。对于欠定混合的情况，采用最短路径法[9]恢复源信号。推断源的问题可以用式（19）表示：由式(18)以及式(19)即可得到源信号估计。其主要思想是求每个时刻最优解，然后得到恢复的源信号。2仿真研究构造的两个仿真信号s1(t)、s2(t)如式(20)所示，式中：f1=20Hz，f2=30Hz，f3=50Hz，f4=60Hz。选取采样频率fs=5000Hz，采样点数N=20000。绘制两个源信号s1(t)、s2(t)的时域波形如图1所示。相应的频谱如图2所示。图1源信号时域波形图图2源信号频谱图用1×2维的随机混合矩阵A，将两个源信号s1(t)、s2(t)进行随机混合，混合后得到1路的混合信号如图3所示。已知源信号的个数为2，混合后观测信号的个数为1，因此，该盲源分离问题为欠定求解问题。图3混叠信号时域波形图利用四线性平行因子盲分离算法将该混合信号分离，得到恢复源信号如图4所示。相应的频谱如图5所示。图4恢复信号时域波形图图5恢复信号的频谱对比图4和图1可知，采用四平行因子盲分离算法得到的恢复信号与源信号相似度不是特别高，分离效果不理想。由图5可知，分离得到的恢复信号的特征频率分别为20Hz、30Hz、50Hz、60Hz，与源信号特征频率不一致，分离得到的频谱图效果不是很理想，因此，在欠定条件下直接采用四线性平行因子盲源分离算法分解效果不理想，需对其进一步改进。为了改进四线性平行因子盲源分离算法的这种不足，用集合经验模态平行因子盲源分离算法分解混合信号。该算法首先用EEMD算法分解混合信号，通过分离得到14个子分量，然后选择相关子分量。选择相关分量的具体做法是利用相关系数求各子分量与观测信号的相关度，一般认为信号与观测信号越相似，则相关度的值越接近为1。根据相关度的值定义其在0到0.3之间为微相关，0.3到0.5之间为实相关，0.5到0.8之间为显著相关，0.8到1之间为高度相关。根据EEMD算法分解得到的14个子分量的相关系数如表1所示。由表1可知，子分量IMF1、IMF6、IMF7、IMF8的相关系数值相对较高，其他分量均为不相关，所以选择IMF1、IMF6、IMF7、IMF8为相关分量并将其与原观测信号组合得到新观测信号，这样就将欠定转化为超定。绘制相关分量IMF1、IMF6、IMF7、IMF8的时域波形图和频谱图，得到的图形分别如图6以及图7所示。图6根据EEMD分解得到相关子分量时域波形图图7根据EEMD分解得到相关子分量的幅值谱图表1各子分量相关系数将通过分解得到相关子分量和混叠信号组合，得到新观测信号，用核一致算法估计振源数，如图9所示。图9组件数由图9可知，振源数为2。将得到的观测信号用四线性平行因子盲分离算法处理，用QALS算法拟合迭代，当达到迭代终止条件时得到4个载荷矩阵A、B、C、D。由图8可知，4个载荷矩阵中，载荷矩阵B、D分解形式相同，表明建立的模型符合式(11)所示的四线性模型，并未偏离四线性模型。图8载荷矩阵用集合经验模态四线性平行因子盲源分离算法处理通过组合得到的观测信号，依据模型求得混合矩阵估计，以及根据式(19)求得估计源信号，绘制恢复信号的时域波形图如10图所示。相应的恢复信号频谱如图11所示。由图11可知，图11(a)中的特征频率分别为50Hz、60Hz，图11(b)中的特征频率分别为20Hz、30Hz，除排列顺序有所不同外，其分别与两源信号的特征频率相一致，这种排列顺序改变并不影响分离结果的正确性。对比图10与源信号时域波形，恢复信号的时域波形图和源信号时域波形图非常相似。由此可知，提出的方法能够有效分离混合信号，分离效果非常令人满意。对比提出的方法和传统的四线性平行因子盲源分离算法，发现欠定条件下采用四线性平行因子盲源分离算法并不能得到很好的分离效果，而采用集合经验模态四线性平行因子盲源分离算法可以很好分离混合信号并且分离效果比四线性平行因子盲源分离算法好很多，证明了该算法的有效性。图10恢复信号时域波形图11恢复信号的频谱图3实验研究为了验证所提算法的有效性，将集合经验模态四线性平行因子盲分离算法应用于滚动轴承复合故障盲源分离中。实验所需要的数据由滚动轴承诊断试验台采集，试验台由旋转电机、滚子轴承和加速度传感器组成。轴承的代号是N205EM。电机转速约为1300r/min，转动频率约为21.7Hz，采样频率fs=100kHz。数据中包含一组内滚道和滚子的复合故障信号。使用线切割机技术人为制造缺陷以产生故障源，内圈故障大小为0.3mm×0.15mm，滚子故障大小为0.3mm×0.15mm，利用轴承故障频率机理公式可计算得到滚子的故障特征频率fb=102.3Hz，内圈故障特征频率fi=145.7Hz。所用信号的采样点数为1001000，绘制采集的混合观测信号时域波形图，如图12所示。图12观测信号时域波形图将观测信号用EEMD算法分解，得到16个IMF分量。利用相关系数法选择相关性较大的子分量，各子分量的相关系数如表2所示。由表2可知，子分量IMF1、IMF2、IMF3、IMF4、IMF5的相关系数均大于0.3，为有效子分量，其他子分量信号的相似系数值均在0.3以下，所以选择IMF1、IMF2、IMF3、IMF4、IMF5为相关分量。将选择的相关子分量与原观测信号组合得到新的观测信号。绘制IMF1、IMF2、IMF3、IMF4、MF5的时域图，得到的图形如图13所示。表2各子分量相关系数图13采用EEMD分解观测信号得到的有效子分量时域波形图将分解得到的这5个相关的子分量信号和混叠信号重新构造成新的观测信号，计算观测信号段的时滞协方差矩阵，并将所有数据段的时滞协方差叠加成4阶张量形式，得到四线性模型，用QALS拟合，当满足式(17)时迭代结束，求得估计矩阵。利用所求的估计矩阵和最短路径法得到恢复信号，绘制恢复信号的时域图和频谱如图13、图14所示。图14恢复信号时域波形