CN104655425A - 基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法 - Google Patents

Abstract

本发明提供了一种基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，该方法克服了现有单通道机械复合故障诊断方法中信号分解不完备、重构信号难以较好保持观测信号特征等不足，通过完备总体经验模态分解法实现信号从一维到高维的转化，保证了分解的完备性，抑制了模态混叠现象；同时，在盲源信号的特征提取处理过程中，引入了基于稀疏表示的维数约简方法，通过稀疏表示对数据进行稀疏重构，从全局数据中提取出数据特征信息，使得重构信号能够更好的保持了观测信号数据特征；并且，在对待测轴承机型故障类型分类处理过程中，引入了大间隔分布学习分类方法，借助大间隔分布学习分类方法的泛化能力，能够提高轴承故障诊断的准确性和有效性。

Description

基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法

技术领域

本发明属于机械故障诊断和计算机人工智能技术领域，尤其涉及一种基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法。

背景技术

盲信号处理技术是在20世纪80年代后期发展起来的一项新技术，具有优异的盲分离、盲辨识、特征提取能力。它的主要思想是：在源信号先验信息匮乏、传输通道参数均未知的情况下，仅通过观测信号就可对源信号和传输通道参数进行估计。近些年，盲源分离的应用逐渐延伸到机械领域，为故障诊断的研究提供了新的方法。但是，目前应用于机械故障振动源分离的盲源分离方法仍存在两点不足：其一，常规盲源分离算法为获得较好的分离效果，假设信号的观测通道数不能少于源信号数，但工程中由于各种因素的限制，有时存在一个观测通道的特殊情况，此时常规盲源分离方法已不再适用；其次，实际的机械振动信号中常伴有大量的非平稳信号，但传统盲源分离方法通常将采集到的振动信号视作平稳信号进行处理，从而造成分离结果的精度误差较大。

针对盲源分离中的极端欠定(单通道)这一难点问题，单通道盲源分离方法被提出。单通道盲源分离可以有效解决常规盲源分离方法在机械诊断过程中要求观测通道数不得少于源信号数的限制，保证分离结果的精确性。单通道盲源分离方法为解决单通道下传统盲源分离算法失效问题首先要对观测信号进行升维处理，以满足盲源分离假设。当前使用的升维方法主要包括：Hoonbin等在《Separation of fault features from a single-channel mechanical signal mixture using wavelet decomposition》中提出的利用小波分解单通道机械故障信号的盲分离方法；周晓峰等在《相关机械振源的盲源分离方法》中提出的基于小波包分解的机械振源盲源分离方法；毋文峰等在《基于经验模式分解的单通道机械信号盲分离》中提出的经验模式分解法以及孟宗等在《基于EEMD子带提取相关机械振动信号单通道盲源分离》中提出的总体经验模态分解法。上述方法中小波分解缺乏自适应性且带有能量的泄漏。小波包分解方法将观测信号用小波包分解成不同频率子带信号，将单通道的一维信号转变为多维信号，对于相关振源分离具有一定优势，但该方法存在小波包分解层数确定的问题，若分解的层数过小重构的局部子带信号中仍然是相关信号的复合，分离效果欠佳；若分解的层数过大会导致算法性能退化。经验模式分解法实现了单观测通道下机械振动信号盲分离，亦能克服源信号稀疏性限制。但该方法存在模态混叠较严重的缺陷。

总体经验模态分解方法为改进的经验模式分解法，其本质是把非平稳信号进行平稳化处理，对信号中不同尺度的波动由高频到低频的顺序依次进行分解，进而得到一系列包含有不同尺度信号的本征模态函数。通过分解，单通道观测信号实现从一维到高维的转化。该方法能有效抑制混叠现象，但在使用中仍然存在添加白噪声残留较大和分解不完备等不足。因此，现有方法均难以实现单通道机械信号的精确分离。

信号升维之后是信号重构过程，在当前的信号重构选取过程中，大多通过求取互信息的方法来实现。由于单一子带分量统计独立性越高该子带越能清晰体现源信号独立子带成分，越适合组成子带观测信号进行盲分离，而互信息是衡量一组随机变量间独立性的重要尺度，故将其作为子带选取准则。通过计算各个子带内信号的互信息，选择互信息较小的几个子带信号重构成新的观测信号。此外，黄书华等在《一种自适应单入多出盲源分离方法》中通过主成分分析(PCA)法，针对多子带信息进行降维，得到主要特征量，以此来实现信号选取。然而上述方法重构的信号难以很好保持观测信号特征，从而难以取得较好的分离效果。而以稀疏理论为基础的稀疏表示降维算法，能够从全局数据中提取出反映某类数据特征的信息，在数据特性保持方面更具优势。通过升维和信号重构单通道问题转化为正定的盲源分离问题，从而可以通过独立分量分析方法估计观测信号的分离矩阵；最后用该矩阵分离原始观测信号，实现相关机械振源信号的分离。

单通道盲源信号实现有效分离后，在对分离信号的故障诊断环节，谱分析和支持向量机智能分类是最常用的故障诊断方法。谱分析法对分离的信号求取FFT谱或包络谱，通过对计算所得的故障频率和求取得频谱图像进行分析，判断故障类型。但是在大多数情况下，谱分析法需要人工判断。所得频谱有时十分复杂，难以直观识别，需要分析者具有丰富的分析经验。支持向量机(Support Vector Machine，SVM)是在统计学习理论的基础上提出的一种新的机器学习方法。SVM能很好的解决小样本学习问题，与神经网络相比有更好的泛化能力，且训练收敛速度更快，可以实现故障的智能分类。然而SVM是基于最大化间隔的机器学习方法，仅仅考虑了单点间隔，并未考虑间隔分布。

轴承是旋转机械中最为重要的机械零件之一，广泛应用于化工、冶金、电力、航空等各个重要部门，同时它也是最易损坏的元件之一。轴承性能与工况的好坏直接影响到与之相关联的轴以及安装在转轴上的齿轮乃至整台机器设备的性能，其缺陷会导致设备产生异常振动和噪声，甚至造成设备损坏。因此，对滚动轴承故障进行诊断，尤其是对于早初期故障的分析，避免事故的发生，在生产实际中尤为重要。

然而，由于现有的机械故障诊断手段都存在一定的缺陷，因此采用现有的机械故障诊断手段都难以对轴承故障进行较为准确的识别和诊断。

发明内容

针对现有技术中存在的上述问题，为了克服现有的机械故障诊断手段中存在的信号分解不完备、难以较好保持观测信号特征、导致故障诊断准确性不高等问题，本发明提供了一种基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，该方法通过完备总体经验模态分解法实现信号从一维到高维的转化，保证了分解的完备性，抑制了模态混叠现象，以提高轴承故障诊断的准确性和有效性。

为实现上述目的，本发明采用了如下技术手段：

基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，包括如下步骤：

1)在待测轴承转动工作时，通过单个加速度传感器采集待测轴承的振动加速度信号，作为待测轴承的单通道的盲源信号；

2)将获取到的盲源信号使用完备总体经验模态分解法进行分解，将盲源信号分解成一系列的本征模函数，并由分解后得到的全部本征模函数分量构成盲源信号的观测矩阵X；

3)根据盲源信号的观测矩阵X，采用贝叶斯信息准则估算出盲源信号所含的信号源数目；

4)以盲源信号所含的信号源数目作为目标维数，利用稀疏保持投影方法将盲源信号的观测矩阵X进行维数约简，约简至目标维数，构成投影特征矩阵；

5)将投影特征矩阵投入特征矩阵近似联合对角化分离算法，计算出分离矩阵，通过分离矩阵与观测信号矩阵求得盲源信号的分离信号；

6)分别计算历史故障信号和待测轴承的盲源信号的分离信号的时频域特征，构成历史故障信号和待测轴承盲源信号各自的特征矩阵；

7)以历史故障信号的特征矩阵作为分类基准数据，将待测轴承的盲源信号的特征矩阵作为分类识别特征，利用大间隔分布学习分类方法对盲源信号进行分类识别，确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步骤2)具体为：

21)分别对获取到的待测轴承的单通道盲源信号添加幅值为E₀的不同高斯白噪声，得到多个添加不同噪声的盲源噪声信号，分别对各个盲源噪声信号采用完备总体经验模态分解法进行分解；

22)对各个盲源噪声信号采用完备总体经验模态分解法分解后得到的所有分量进行总体平均，并将得到的总体平均分量视为盲源信号的第一阶的本征模函数分量I₁(t)，从盲源 信号中将第一阶的本征模函数分量I₁(t)分离出来后得到的剩余信号作为第一阶剩余信号r₁(t)，t表示时间；

23)将幅值为E₁的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W₁(t)；分别将各个噪声的幅值E₁与其第一阶本征模函数分量W₁(t)的积E₁W₁(t)添加到第一阶剩余信号r₁(t)，得到多个不同的第一阶重构信号S^* ₁(t)，并分别对各个第一阶重构信号S^* ₁(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第二阶本征模函数分量I₂(t)，再从第一阶剩余信号r₁(t)中将第二阶本征模函数分量I₂(t)分离出来后得到的剩余信号作为第二阶剩余信号r₂(t)；

24)对于得到的第b阶本征模函数分量和I_b(t)和第b阶剩余信号r_b(t)，b>1，将幅值为E_q的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W_b(t)；分别将各个噪声的幅值E_b与其第一阶本征模函数分量W_b(t)的积E_bW_b(t)添加到第b阶剩余信号r_b(t)，得到多个不同的第b阶重构信号S^* _b(t)，并分别对各个第b阶重构信号S^* _b(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第b+1阶本征模函数分量I_b+1(t)，再从第b阶剩余信号r_b(t)中将第b+1阶本征模函数分量I_b+1(t)分离出来后得到的剩余信号作为第b+1阶剩余信号r_b+1(t)；

25)重复执行步骤24)，直到剩余信号的极值点个数不超过两个时终止分解，由分解后得到的全部本征模函数分量，构成盲源信号的观测矩阵X。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步骤3)具体为：

31)计算观测矩阵X的相关矩阵R_x：

R_x＝XX^H；

其中H表示复共轭转置；

32)对相关矩阵R_x进行奇异值分解后，得到非零特征值λ₁,…,λ_j,…,λ_L，j∈{1,2,…,L}，L表示非零特征值的个数；对于序号k从1到L的取值，根据贝叶斯信息准则，分别计算各个取值的序号k的BIC代价函数值：

$\mathrm{BIC}\left(k\right)={\left(\underset{j=1}{\overset{k}{\mathrm{\Π}}}{\mathrm{\λ}}_j\right)}^{-N/2}{\mathrm{\σ}}_k^{-N(L-k)/2}{N}^{-(u+k)/2};$

其中，N为盲源信号的采样点数；u＝Lk-k(k+1)/2；

34)根据各个取值的序号k的BIC代价函数值，寻找使BIC代价函数值最大的序号k的取值，则盲源信号所含的信号源数目d的取值即为使BIC代价函数值最大的序号k的取值。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步骤4)具体为：

41)计算稀疏系数：根据稀疏理论构建目标函数：

${\min }_{S_i}{\left|\right|S_i\left|\right|}_0,$ s_·t_·x_i＝Xs_i

其中s_i表示盲源信号的观测矩阵X中第i个本征模函数分量x_i由其它本征模函数分量进行稀疏表示时对应的稀疏系数；

42)计算投影矩阵A：根据计算得到的观测矩阵X中各个本征模函数分量对应的稀疏系数，分别对观测矩阵X中的各个本征模函数分量进行稀疏重构，并使其重构误差最小，目标函数如下：

$\min \mathrm{\ϵ}\left(A\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left|\right|A^T{x}_i-A^{T}X{s}_i\left|\right|}^2=\mathrm{tr}\left\{A^{T}X(I-S-S^T+S^T)X^{T}A\right\};$

其中，S＝[s_1,s_2,,...,s_N]为稀疏系数矩阵，N为盲源信号的采样点数；增加约束A^TXX^TA＝I，I为单位矩阵，并求解投影矩阵A的最大值：

maxε(A)＝A^TXS_βX^TA/A^TXX^TA；

其中S_β＝S+S^T-S^TS，将上式转化为XS_βX^TA＝λXX^TA，求解广义特征值和特征向量，取出其中最大的d个特征值λ对应的特征向量a_1,a_2,,...,a_d，得到投影矩阵A＝[a_1,a_2,,...,a_d]，d为盲源信号所含的信号源数目；

43)将观测矩阵X在投影矩阵A上投影，得到一个d×N维的投影特征矩阵X₀。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步骤5)具体为：

51)计算投影特征矩阵X₀的自相关函数R_XR：

R_XR＝E(X₀X₀ ^*)；

其中上标*表示复共轭；对自相关函数R_XR进行广义特征值分解，得到的d个最大的特 征根和相应的特征向量分别用λ₁,…,λ_d和h₁,…,h_d表示；

52)对投影特征矩阵X₀进行预白化处理：

白化信号z(t)为：

z(t)＝[z₁(t),z₂(t),…,z_i(t),…,z_d(t)]^T；

其中1≤i≤d，d为盲源信号所含的信号源数目；

白化矩阵W为：

$W={[{\left({\mathrm{\λ}}_1\right)}^{-\frac{1}{2}}{h}_1,...,{\left({\mathrm{\λ}}_d\right)}^{-\frac{1}{2}}{h}_d]}^H;$

其中，上标H表示复共轭转置；

由盲源分离模型X₀＝BS_s可得：

z＝WX₀＝WBS_s＝US_s；

其中，B为混合矩阵，S_s为分离信号矩阵，U为待求解的酉矩阵；

53)联合对角化白化后的观测信号：

对于d阶的一个任意矩阵M，z(t)的四阶累积向量矩阵的第o行第p列元素[Q_z(M)]_op为：

${\left[Q_z\left(M\right)\right]}_{\mathrm{op}}\stackrel{\mathrm{def}}{=}\underset{q=1}{\overset{d}{\mathrm{\Σ}}}\underset{r=1}{\overset{d}{\mathrm{\Σ}}}\mathrm{cum}(z_o,{z_p}^{*},z_q,{z_r}^{*})m_{\mathrm{qr}};$

式中，um(z_o,z_p ^*,z_q,z_r ^*)表示白化信号z(t)中第o个分量z_o、第p个分量z_p、第q个分量z_q、第r个分量z_r的四阶累积量矩阵，1≤o≤d，1≤p≤d，1≤q≤d，1≤r≤d；m_qr为矩阵M的第q行第r列元素；

54)用U^TQ_z(M_i)U中对角元素的平方和作为对角化过程参照函数，即：

$F\left(U\right)=\underset{M_i\∈M}{\mathrm{\Σ}}{\left|\right|\mathrm{diag}\left(U^T{Q}_z(M_i\right)U\left|\right|}^2;$

通过优化算法使得参照函数F(U)最小化获得酉矩阵U的解，即获得盲源信号的分离信号S_s：

S_s＝U^TWX₀。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步骤6)具体为：

针对d×N维的盲源信号的分离信号S_s，d为盲源信号所含的信号源数目，N为盲源 信号的采样点数，通过计算盲源信号每周期采样点数R将分离信号S_s的每个分离信号分割为g个周期数据y_i，i＝1,2,...,g，在时域分别按周期提取分离信号S_s的峰值、平均幅值、标准差、均方根值、峭度指标、峰值指标、波形指标、裕度指标8个时域特征来反映信号变化；其中：

峰值为P＝max|y_i|；      平均幅值为 $\stackrel{\‾}{y_i}=\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}\left|y_i\right|;$

标准差为 $T=\sqrt{\frac{1}{R-1}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}|y_i-\stackrel{\‾}{y_i}|};$       均方根值为 $y_{\mathrm{rms}}=\sqrt{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{y_i}^2};$

峭度指标为 $K=\frac{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{\left(y_i\right)}^4}{{y_{\mathrm{rms}}}^4};$       峰值指标为 $C_f=\frac{P}{{y_{\mathrm{rms}}}^4};$

波形指标为 $S_f=\frac{y_{\mathrm{rms}}}{\left|\stackrel{\‾}{y_i}\right|};$       裕度指标为 ${\mathrm{CL}}_f=\frac{P}{{\left(\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}\sqrt{\left|y_i\right|}\right)}^2};$

然后，采用db4小波包函数对分离信号S_s进行3层正交小波包分解用以进行频域特征提取，得到分离信号S_s对应的由全频带均匀划分8个子频带的滤波信号，通过重构每个子频带的节点小波包系数，确保重构信号和原分离信号S_s的长度一样，然后计算每个子频带滤波信号的各个采样点的幅值平方和作为其能量，将各个子频带滤波信号的能量与全频带重构信号总能量的比值作为其频域统计特征，即(E₀/E),(E₁/E),...,(E₇/E)，其中E₀,E₁,...,E₇分别为8个子频带滤波信号的能量，E为全频带重构信号的总能量；由此，针对分离信号S_s，共构造了8个时域特征和8个频域统计特征，由该16个特征构成盲源信号的分离信号S_s的特征矩阵；

针对待测轴承盲源信号的分离信号，构造上述的8个时域特征和8个频域统计特征，构成待测轴承盲源信号的特征矩阵；

所述历史轴承故障数据中的每个故障信号数据，是对已知故障工况类型的轴承进行振动加速度信号采样后得到的振动加速度信号，并标记有其所属故障工况的类别标签；根据历史轴承故障数据包含的各种不同故障工况，从中分别在每种故障工况下选取多个故障信号数据，分别作为历史轴承故障信号样本，并针对每个历史轴承故障信号样本求取其分离信号，并同样构造上述的8个时域特征和8个频域统计特征，构成每个历史轴承故障信号样本的特征矩阵。

上述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法中，具体而言，所述步 骤7)具体为：

71)将待测轴承的盲源信号的特征矩阵作为测试样本，将从历史轴承故障数据中选取的全部历史轴承故障信号样本的特征矩阵作为训练样本，训练样本集表示为：

V＝{(v₁,c₁),(v₂,c₂),…,(v_i,c_i),…,(v_m,c_m)}；

其中，m为训练样本集中训练样本的个数，v_i表示第i个训练样本的特征矩阵，c_i表示第i个训练样本的类别标签；

72)利用一阶和二阶统计量，即均值和方差来描述间隔分布：

$\stackrel{\‾}{\mathrm{\γ}}=\frac{1}{m}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{c}_i{w}^T\mathrm{\φ}\left(c_i\right)=\frac{1}{m}{\left(V_{\mathrm{\φ}}c\right)}^{T}w;$

$\widehat{\mathrm{\γ}}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{m}{\mathrm{\Σ}}}{(c_i{w}^T\mathrm{\φ}\left(v_i\right)-c_j{w}^T\mathrm{\φ}\left(v_j\right))}^2=\frac{2}{m^2}\left(m{w}^T{V}_{\mathrm{\φ}}c{c}^T{V}_{\mathrm{\φ}}^{T}w\right);$

式中，V_φ＝[φ(v₁),φ(v₂),…,φ(v_i),…,φ(v_m)]，列向量c＝[c₁,c₂,…,c_i,…,c_m]^T；大间隔分布学习模型的线性系数向量σ＝[σ₁，σ₂，…，σ_m]为线性系数；φ(v_i)和φ(b_j)分别为第i个训练样本的特征矩阵v_i和第j个训练样本的特征矩阵v_j在核函数K上的特征映射，即有：

K(v_i,v_j)＝φ(v_i)^Tφ(v_j)；

其中，i＝1,2,...,m，j＝1,2,...,m，且i≠j；

在训练样本可以实现零误差的精确区分情况下，同时最大化间隔均值和最小化间隔方差引出：

$\begin{array}{cc}\underset{w}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}\end{array};$

s.t.c_iw^Tφ(v_i)≥1，i＝1,2,...,m；

式中：β₁和β₂分别为均值、方差的权值参数；

在训练样本不能实现零误差精确区分情况下，引出：

$\begin{array}{cc}\underset{w,\mathrm{\ϵ}}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}+C_c\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_i;\end{array}$

s.t.c_iw^Tφ(v_i)≥1-ε_i，ε_i≥0，i＝1,2,...，m；

C_c为惩罚参数，ε＝[ε₁,ε₂,…,ε_m]^T为松弛变量；

把带入上式，根据拉格朗日乘法确定系数σ＝[σ₁,σ₂,…,σ_i,…,σ_m]；

73)训练样本集中包含的各种不同故障工况的训练样本，针对其中第z类故障工况的训练样本看作正类，z∈{1,2,…,Z}，Z表示训练样本集中包含的故障工况的类别总数，将其它Z-1类故障工况的训练样本看作负类；将正类训练样本的类别标签取值为1，正类训练样本的类别标签取值为-1，然后按下式计算输入量α在第z类故障工况下对应的故障类型预测函数：

$\mathrm{sgn}\left(w^T\mathrm{\φ}\left(\mathrm{\α}\right)\right)=\mathrm{sgn}\left(\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_{i}K(v_i,\mathrm{\α})\right);$

其中，sgn(·)表示符号函数；

74)将测试样本分别作为各种工况故障工况下对应的故障类型预测函数的输入量，计算出测试样本作为输入量在各个故障工况类别下对应的故障类型预测函数值，将故障类型预测函数值为1时对应的故障工况类别判定为测试样本的故障工况类别，从而确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

相比于现有技术，本发明具有如下有益效果：

1、本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，克服了现有单通道机械复合故障诊断方法中信号分解不完备、重构信号难以较好保持观测信号特征等不足，通过完备总体经验模态分解法实现信号从一维到高维的转化，保证了分解的完备性，抑制了模态混叠现象。

2、本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，在盲源信号的特征提取处理过程中，引入了基于稀疏表示的维数约简方法，即稀疏保持投影，通过稀疏表示对数据进行稀疏重构，从全局数据中提取出数据特征信息，使得重构信号能够更好的保持了观测信号数据特征。

3、本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，在对待测轴承机型故障类型分类处理过程中，引入了大间隔分布学习分类方法，借助大间隔分布学习分类方法的泛化能力，能够提高轴承故障诊断的准确性和有效性。

附图说明

图1为本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法的流程框图。

图2为本发明实验示例中待测轴承的单通道盲源信号示例图。

图3为本发明实验示例中图2所示的盲源信号经过总体经验模态分解得到的9个本征模函数分量示例图。

图4为本发明实验示例中图2所示的盲源信号的分离信号示例图。

具体实施方式

本发明提出了一种基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，该方法通过完备总体经验模态分解法(Complete Ensemble Empirical Mode Decomposition，CEEMD)实现信号从一维到高维的转化，保证了分解的完备性，抑制了模态混叠现象，得到了更好的分解效果。针对分解后的多维子带信息，引入一种基于稀疏表示的维数约简方法完成高维到目标维数的转化，使的重构信号能够更好的保持观测信号数据特征。利用特征矩阵联合对角化盲源分离算法分离出源信号。对分离出源信号以及历史故障信号进行时频特征提取，根据所提取特征利用大间隔分布学习(Large Margin Distribution Machine，缩写为LDM)分类方法实现轴承故障的分类诊断，以期获得较好的分类效果。

本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法的处理流程如图1所示，其包括如下步骤：

1)在待测轴承转动工作时，通过单个加速度传感器采集待测轴承的振动加速度信号，作为待测轴承的单通道的盲源信号。

轴承在正常、轴承内圈故障、轴承滚动体故障、轴承外圈故障等不同故障工况下，其转动工作的振动加速度信号相互之间存在一定的差异，因此可以基于滚动轴承在不同工况下的振动加速度信号数据作为轴承的单通道盲源信号，用以进行故障情况的识别。

2)将获取到的盲源信号使用完备总体经验模态分解法进行分解，将盲源信号分解成一系列的本征模函数，并由分解后得到的全部本征模函数分量构成盲源信号的观测矩阵X。

为消除模态混叠现象，该步骤通过添加不同的高斯白噪声，对盲源信号继续拧总体经验模态分解(Empirical Mode Decomposition，缩写为EMD)处理，其具体实现方式如下：

21)分别对获取到的待测轴承的单通道盲源信号添加幅值为E₀的不同高斯白噪声，得到多个添加不同噪声的盲源噪声信号，分别对各个盲源噪声信号采用完备总体经验模态分解法进行分解；

22)对各个盲源噪声信号采用完备总体经验模态分解法分解后得到的所有分量进行总体平均，并将得到的总体平均分量视为盲源信号的第一阶的本征模函数分量I₁(t)，从盲源信号中将第一阶的本征模函数分量I₁(t)分离出来后得到的剩余信号作为第一阶剩余信号r₁(t)，t表示时间；

23)将幅值为E₁的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W₁(t)；分别将各个噪声的幅值E₁与其第一阶本征 模函数分量W₁(t)的积E₁W₁(t)添加到第一阶剩余信号r₁(t)，得到多个不同的第一阶重构信号S^* ₁(t)，并分别对各个第一阶重构信号S^* ₁(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第二阶本征模函数分量I₂(t)，再从第一阶剩余信号r₁(t)中将第二阶本征模函数分量I₂(t)分离出来后得到的剩余信号作为第二阶剩余信号r₂(t)；

24)对于得到的第b阶本征模函数分量和I_b(t)和第b阶剩余信号r_b(t)，b>1，将幅值为E_q的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W_b(t)；分别将各个噪声的幅值E_b与其第一阶本征模函数分量W_b(t)的积E_bW_b(t)添加到第b阶剩余信号r_b(t)，得到多个不同的第b阶重构信号S^* _b(t)，并分别对各个第b阶重构信号S^* _b(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第b+1阶本征模函数分量I_b+1(t)，再从第b阶剩余信号r_b(t)中将第b+1阶本征模函数分量I_b+1(t)分离出来后得到的剩余信号作为第b+1阶剩余信号r_b+1(t)；

25)重复执行步骤24)，直到剩余信号的极值点个数不超过两个时终止分解，由分解后得到的全部本征模函数分量，构成盲源信号的观测矩阵X。所得到的观测矩阵X是一个M×N的，其中M表示分解后所得到的所有分量个数，N表示盲源信号的采样点数。

3)根据盲源信号的观测矩阵X，采用贝叶斯信息准则估算出盲源信号所含的信号源数目。

该步骤的具体实现方式如下：

31)计算观测矩阵X的相关矩阵R_x：

R_x＝XX^H；

其中H表示复共轭转置；

32)对相关矩阵R_x进行奇异值分解后，得到非零特征值λ₁,…,λ_j,…,λ_L，j∈{1,2,…,L}，L表示非零特征值的个数；对于序号k从1到L的取值，根据贝叶斯信息准则，分别计算各个取值的序号k的BIC代价函数值：

$\mathrm{BIC}\left(k\right)={\left(\underset{j=1}{\overset{k}{\mathrm{\Π}}}{\mathrm{\λ}}_j\right)}^{-N/2}{\mathrm{\σ}}_k^{-N(L-k)/2}{N}^{-(u+k)/2};$

其中，N为盲源信号的采样点数；u＝Lk-k(k+1)/2；

34)根据各个取值的序号k的BIC代价函数值，寻找使BIC代价函数值最大的序号k的取值，则盲源信号所含的信号源数目d的取值即为使BIC代价函数值最大的序号k的取值。

4)以盲源信号所含的信号源数目作为目标维数，利用稀疏保持投影方法将盲源信号的观测矩阵X进行维数约简，约简至目标维数，构成投影特征矩阵。

该步骤的具体实现方式如下：

41)计算稀疏系数：根据稀疏理论构建目标函数：

${\min }_{S_i}{\left|\right|S_i\left|\right|}_0,$ s_·t_·x_i＝Xs_i

其中s_i表示盲源信号的观测矩阵X中第i个本征模函数分量x_i由其它本征模函数分量进行稀疏表示时对应的稀疏系数；

42)计算投影矩阵A：根据计算得到的观测矩阵X中各个本征模函数分量对应的稀疏系数，分别对观测矩阵X中的各个本征模函数分量进行稀疏重构，并使其重构误差最小，目标函数如下：

$\min \mathrm{\ϵ}\left(A\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left|\right|A^T{x}_i-A^{T}X{s}_i\left|\right|}^2=\mathrm{tr}\left\{A^{T}X(I-S-S^T+S^T)X^{T}A\right\};$

其中，S＝[s_1,s_2,,...,s_N]为稀疏系数矩阵，N为盲源信号的采样点数；增加约束A^TXX^TA＝I，I为单位矩阵，并求解投影矩阵A的最大值：

maxε(A)＝A^TXS_βX^TA/A^TXX^TA；

其中S_β＝S+S^T-S^TS，将上式转化为XS_βX^TA＝λXX^TA，求解广义特征值和特征向量，取出其中最大的d个特征值λ对应的特征向量a_1,a_2,,...,a_d，得到投影矩阵A＝[a_1,a_2,,...,a_d]，d为盲源信号所含的信号源数目；

43)将观测矩阵X在投影矩阵A上投影，得到一个d×N维的投影特征矩阵X₀。

5)将投影特征矩阵投入特征矩阵近似联合对角化分离算法，计算出分离矩阵，通过分离矩阵与观测信号矩阵求得盲源信号的分离信号。

该步骤的具体实现方式如下：

51)计算投影特征矩阵X₀的自相关函数R_XR：

R_XR＝E(X₀X₀ ^*)；

其中上标*表示复共轭；对自相关函数R_XR进行广义特征值分解，得到的d个最大的特 征根和相应的特征向量分别用λ₁,…,λ_d和h₁,…,h_d表示；

52)对投影特征矩阵X₀进行预白化处理：

白化信号z(t)为：

z(t)＝[z₁(t),z₂(t),…,z_i(t),…,z_d(t)]^T；

其中1≤i≤d，d为盲源信号所含的信号源数目；

白化矩阵W为：

$W={[{\left({\mathrm{\λ}}_1\right)}^{-\frac{1}{2}}{h}_1,...,{\left({\mathrm{\λ}}_d\right)}^{-\frac{1}{2}}{h}_d]}^H;$

其中，上标H表示复共轭转置；

由盲源分离模型X₀＝BS_s可得：

z＝WX₀＝WBS_s＝US_s；

其中，B为混合矩阵，S_s为分离信号矩阵，U为待求解的酉矩阵；这样，经过白化过程，一个混合矩阵B的求解问题转化为一个酉矩阵U的求解问题；

53)联合对角化白化后的观测信号：

对于d阶的一个任意矩阵M，z(t)的四阶累积向量矩阵的第o行第p列元素[Q_z(M)]_op为：

${\left[Q_z\left(M\right)\right]}_{\mathrm{op}}\stackrel{\mathrm{def}}{=}\underset{q=1}{\overset{d}{\mathrm{\Σ}}}\underset{r=1}{\overset{d}{\mathrm{\Σ}}}\mathrm{cum}(z_o,{z_p}^{*},z_q,{z_r}^{*})m_{\mathrm{qr}};$

式中，um(z_o,z_p ^*,z_q,z_r ^*)表示白化信号z(t)中第o个分量z_o、第p个分量z_p、第q个分量z_q、第r个分量z_r的四阶累积量矩阵，1≤o≤d，1≤p≤d，1≤q≤d，1≤r≤d；m_kl为矩阵M的第k行第l列元素；

54)用U^TQ_z(M_i)U中对角元素的平方和作为对角化过程参照函数，即：

$F\left(U\right)=\underset{M_i\∈M}{\mathrm{\Σ}}{\left|\right|\mathrm{diag}\left(U^T{Q}_z(M_i\right)U\left|\right|}^2;$

通过优化算法使得参照函数F(U)最小化获得酉矩阵U的解，即获得盲源信号的分离信号S_s：

S_s＝U^TWX₀。

6)分别计算历史故障信号和待测轴承的盲源信号的分离信号的时频域特征，构成历史故障信号和待测轴承盲源信号各自的特征矩阵。

该步骤的具体实现方式如下：

针对d×N维的盲源信号的分离信号S_s，d为盲源信号所含的信号源数目，N为盲源信号的采样点数，通过计算盲源信号每周期采样点数R将分离信号S_s的每个分离信号分割为g个周期数据y_i，i＝1,2,...,g，在时域分别按周期提取分离信号S_s的峰值、平均幅值、标准差、均方根值、峭度指标、峰值指标、波形指标、裕度指标8个时域特征来反映信号变化；其中：

峰值为P＝max|y_i|；平均幅值为 $\stackrel{\‾}{y_i}=\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}\left|y_i\right|;$

标准差为 $T=\sqrt{\frac{1}{R-1}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}|y_i-\stackrel{\‾}{y_i}|};$       均方根值为 $y_{\mathrm{rms}}=\sqrt{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{y_i}^2};$

峭度指标为 $K=\frac{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{\left(y_i\right)}^4}{{y_{\mathrm{rms}}}^4};$       峰值指标为 $C_f=\frac{P}{{y_{\mathrm{rms}}}^4};$

波形指标为 $S_f=\frac{y_{\mathrm{rms}}}{\left|\stackrel{\‾}{y_i}\right|};$       裕度指标为 ${\mathrm{CL}}_f=\frac{P}{{\left(\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}\sqrt{\left|y_i\right|}\right)}^2};$

然后，采用db4小波包函数对分离信号S_s进行3层正交小波包分解用以进行频域特征提取，得到分离信号S_s对应的由全频带均匀划分8个子频带的滤波信号，通过重构每个子频带的节点小波包系数，确保重构信号和原分离信号S_s的长度一样，然后计算每个子频带滤波信号的各个采样点的幅值平方和作为其能量，将各个子频带滤波信号的能量与全频带重构信号总能量的比值作为其频域统计特征，即(E₀/E),(E₁/E),...,(E₇/E)，其中E₀,E₁,...,E₇分别为8个子频带滤波信号的能量，E为全频带重构信号的总能量；由此，针对分离信号S_s，共构造了8个时域特征和8个频域统计特征，由该16个特征构成盲源信号的分离信号S_s的特征矩阵；

针对待测轴承盲源信号的分离信号，构造上述的8个时域特征和8个频域统计特征，构成待测轴承盲源信号的特征矩阵；

所述历史轴承故障数据中的每个故障信号数据，是对已知故障工况类型的轴承进行振动加速度信号采样后得到的振动加速度信号，并标记有其所属故障工况的类别标签；根据历史轴承故障数据包含的各种不同故障工况，从中分别在每种故障工况下选取多个故障信号数据，分别作为历史轴承故障信号样本，并针对每个历史轴承故障信号样本求取其分离信号，并同样构造上述的8个时域特征和8个频域统计特征，构成每个历史轴承故障信号 样本的特征矩阵。

7)以历史故障信号的特征矩阵作为分类基准数据，将待测轴承的盲源信号的特征矩阵作为分类识别特征，利用大间隔分布学习分类方法对盲源信号进行分类识别，确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

最新的理论研究证明间隔分布相对于最大间隔具有更好的泛化能力。TengZhang等在《Large Margin Distribution Machine》中提出了大间隔分布学习分类方法。本发明就利用大间隔分布学习分类方法对待测轴承盲源信号所属的故障工况类别进行分类诊断，具体实现方式如下：

71)将待测轴承的盲源信号的特征矩阵作为测试样本，将从历史轴承故障数据中选取的全部历史轴承故障信号样本的特征矩阵作为训练样本，训练样本集表示为：

V＝{(v₁,c₁),(v₂,c₂),…,(v_i,c_i),…,(v_m,c_m)}；

其中，m为训练样本集中训练样本的个数，v_i表示第i个训练样本的特征矩阵，c_i表示第i个训练样本的类别标签；

72)利用一阶和二阶统计量，即均值和方差来描述间隔分布：

$\stackrel{\‾}{\mathrm{\γ}}=\frac{1}{m}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{c}_i{w}^T\mathrm{\φ}\left(c_i\right)=\frac{1}{m}{\left(V_{\mathrm{\φ}}c\right)}^{T}w;$

$\widehat{\mathrm{\γ}}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{m}{\mathrm{\Σ}}}{(c_i{w}^T\mathrm{\φ}\left(v_i\right)-c_j{w}^T\mathrm{\φ}\left(v_j\right))}^2=\frac{2}{m^2}\left(m{w}^T{V}_{\mathrm{\φ}}c{c}^T{V}_{\mathrm{\φ}}^{T}w\right);$

式中，V_φ＝[φ(v₁),φ(v₂),…,φ(v_i),…,φ(v_m)]，列向量c＝[c₁,c₂,…,c_i,…,c_m]^T；大间隔分布学习模型的线性系数向量σ＝[σ₁,σ₂,…,σ_m]为线性系数；φ(v_i)和φ(b_j)分别为第i个训练样本的特征矩阵v_i和第j个训练样本的特征矩阵v_j在核函数K上的特征映射，即有：

K(v_i,v_j)＝φ(v_i)^Tφ(v_j)；

其中，i＝1,2,...,m，j＝1,2,...,m，且i≠j；

在训练样本可以实现零误差的精确区分情况下，同时最大化间隔均值和最小化间隔方差引出：

$\begin{array}{cc}\underset{w}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}\end{array};$

s.t.c_iw^Tφ(v_i)≥1，i＝1,2,...,m；

式中：β₁和β₂分别为均值、方差的权值参数；

在训练样本不能实现零误差精确区分情况下，引出：

$\begin{array}{cc}\underset{w,\mathrm{\ϵ}}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}+C_c\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_i;\end{array}$

s.t.c_iw^Tφ(v_i)≥1-ε_i，ε_i≥0，i＝1,2,…,m；

C_c为惩罚参数，ε＝[ε₁,ε₂,…,ε_m]^T为松弛变量；

把带入上式，根据拉格朗日乘法确定系数σ＝[σ₁，σ₂，…，σ_i，…，σ_m]；

73)训练样本集中包含的各种不同故障工况的训练样本，针对其中第z类故障工况的训练样本看作正类，z∈{1,2,…,Z}，Z表示训练样本集中包含的故障工况的类别总数，将其它Z-1类故障工况的训练样本看作负类；将正类训练样本的类别标签取值为1，正类训练样本的类别标签取值为-1，然后按下式计算输入量α在第z类故障工况下对应的故障类型预测函数：

$\mathrm{sgn}\left(w^T\mathrm{\φ}\left(\mathrm{\α}\right)\right)=\mathrm{sgn}\left(\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_{i}K(v_i,\mathrm{\α})\right);$

其中，sgn(·)表示符号函数；

74)将测试样本分别作为各种工况故障工况下对应的故障类型预测函数的输入量，计算出测试样本作为输入量在各个故障工况类别下对应的故障类型预测函数值，将故障类型预测函数值为1时对应的故障工况类别判定为测试样本的故障工况类别，从而确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

为了更好地说明本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法的技术方案和效果，下面通过实验案例对本发明进行进一步的说明。

实验案例： 

本实验利用旋转机械故障模拟实验台采集的轴承数据进行复合故障诊断，实验台通过加速度传感器来采集轴承的振动加速度信号，作为轴承的盲源信号，用以进行轴承故障分类诊断。本实验使用轴承型号为SKF6205，实验台的电机转速为1730r/m，加速度传感器放置于电机驱动端，信号的采样频率为12000次每秒。实验中，对正常、轴承内圈故障、轴承滚动体故障、轴承外圈故障四种不同故障工况下已知故障工况类型的轴承分别进行振动加速度信号采样后得到的振动加速度信号，作为故障信号数据，每种故障工况类型下采样25组故障信号数据，每组故障信号数据采集至少10500个采样点，由这些故障信号数据构成历史轴承故障数据。然后，对于故障工况类型未知的待测轴承，采用本发明的轴承 故障分类诊断方法进行故障分类诊断处理，其处理流程如下：

1)在待测轴承转动工作时，通过单个加速度传感器采集待测轴承的振动加速度信号，作为待测轴承的单通道的盲源信号。

在本实验中，采集到的一个待测轴承的单通道盲源信号如图2所示。

2)将获取到的盲源信号使用完备总体经验模态分解法进行分解，将盲源信号分解成一系列的本征模函数，并由分解后得到的全部本征模函数分量构成盲源信号的观测矩阵X。

本实验中，对如图2所示的盲源信号进行总体经验模态分解，得到9个本征模函数分量，构成观测矩阵X；9个本征模函数分量如图3所示，分别表示为IMF0～IMF8。

3)根据盲源信号的观测矩阵X，采用贝叶斯信息准则估算出盲源信号所含的信号源数目。

本实验中，通过贝叶斯信息准则计算，确定序号k＝2时BIC代价函数取得最大值，因此确定盲源信号所含的信号源数目d＝2。

4)以盲源信号所含的信号源数目作为目标维数，利用稀疏保持投影方法将盲源信号的观测矩阵X进行维数约简，约简至目标维数，构成投影特征矩阵。

由于盲源信号所含的信号源数目d＝2，因此在维数简约求解广义特征值特征向量的过程中，取出其中最大的2个对应的特征向量a₁,a₂构成投影矩阵A＝[a₁,a₂]，再将观测矩阵X在投影矩阵A上投影，得到一个2×N维的投影特征矩阵X₀。

5)将投影特征矩阵投入特征矩阵近似联合对角化分离算法，计算出分离矩阵，通过分离矩阵与观测信号矩阵求得盲源信号的分离信号。

本实验中，求得盲源信号的分离信号如图4所示。为定量评价分离效果，用源信号与分离信号相似系数作为分离的性能指标。相似系数为：

${\mathrm{\ρ}}_{\mathrm{ij}}=\frac{\mathrm{cov}(X_i,S_{s,j})}{\sqrt{\mathrm{cov}(X_i,X_i)\mathrm{cov}(S_{s,j},S_{s,j})}};$

式中，X_i为盲源信号的第i个分量，S_s,j为经盲源分离的与X_i对应的分离量，cov(·)为协方差运算符；0≤|ρ_ij|≤1，ρ_ij越接近于1，则说明盲分离算法得到的分离信号与之对应的盲源信号越相似。本实验中，通过计算得盲源信号相对于分离信号的相似系数ρ₁₂＝0.78，分离信号相对于盲源信号的相似系数ρ₂₁＝0.94，表明取得了较好的分离效果。

6)分别计算历史故障信号和待测轴承的盲源信号的分离信号的时频域特征，构成历史故障信号和待测轴承盲源信号各自的特征矩阵。

针对2×N维的盲源信号的分离信号，首先将其分割为周期数据。由于实验台的电机转速为1730转/分钟，每秒转动频率电机转动周期T＝1/f_r，信号的采样频率fs为12000次每秒，采样间隔t＝1/fs；因此，求得每周期采样点数R＝T/t＝416.2，取R＝420。分割为周期数据后，在时域分别按周期提取分离信号的峰值、平均幅值、标准差、均方根值、峭度指标、峰值指标、波形指标、裕度指标8个时域特征来反映信号变化，然后，采用db4小波包函数对分离信号S_s进行3层正交小波包分解用以进行频域特征提取，得到分离信号S_s对应的由全频带均匀划分8个子频带的滤波信号，通过重构每个子频带的节点小波包系数，确保重构信号和原分离信号的长度一样，然后计算每个子频带滤波信号的各个采样点的幅值平方和作为其能量，将各个子频带滤波信号的能量与全频带重构信号总能量的比值作为其频域统计特征；由此，针对分离信号，共构造了8个时域特征和8个频域统计特征，由该16个特征构成盲源信号的分离信号的特征矩阵。

针对待测轴承盲源信号的分离信号，构造上述的8个时域特征和8个频域统计特征，构成待测轴承盲源信号的特征矩阵。

将历史轴承故障数据中正常、轴承内圈故障、轴承滚动体故障、轴承外圈故障四种不同故障工况下的全部故障信号数据分别作为历史轴承故障信号样本，并针对每个历史轴承故障信号样本求取其分离信号，求取历史轴承故障信号样本分离信号的方法与求取盲源信号的分离信号的方法相一致，并同样构造上述的8个时域特征和8个频域统计特征，构成每个历史轴承故障信号样本的特征矩阵。

7)以历史故障信号的特征矩阵作为分类基准数据，将待测轴承的盲源信号的特征矩阵作为分类识别特征，利用大间隔分布学习分类方法对盲源信号进行分类识别，确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

上述流程介绍了针对一个待测轴承确定其故障类型的分类诊断流程。在本实验中，一共针对50个故障工况未知的待测轴承，采用本发明方法进行了轴承故障分类诊断，然后对该50个待测轴承通过人工检验的方式确定其各自的实际故障情况，并将本发明方法对待测轴承的轴承故障分类诊断结果与其实际故障情况进行比对，其中43组的诊断结果与实际故障情况一致，测试结果的分类精度达到43/50＝86％，表明本发明方法已具备非常好的轴承故障诊断准确性。

综上所述，本发明基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，克服了现有单通道机械复合故障诊断方法中信号分解不完备、重构信号难以较好保持观测信号特征等不足，通过完备总体经验模态分解法实现信号从一维到高维的转化，保证了分解的完 备性，抑制了模态混叠现象；同时，在盲源信号的特征提取处理过程中，引入了基于稀疏表示的维数约简方法，即稀疏保持投影，通过稀疏表示对数据进行稀疏重构，从全局数据中提取出数据特征信息，使得重构信号能够更好的保持了观测信号数据特征；并且，在对待测轴承机型故障类型分类处理过程中，引入了大间隔分布学习分类方法，借助大间隔分布学习分类方法的泛化能力，能够提高轴承故障诊断的准确性和有效性。

最后说明的是，以上实施例仅用以说明本发明的技术方案而非限制，尽管参照较佳实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的宗旨和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (7)

1.基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，包括如下步骤：

1)在待测轴承转动工作时，通过单个加速度传感器采集待测轴承的振动加速度信号，作为待测轴承的单通道的盲源信号；

2)将获取到的盲源信号使用完备总体经验模态分解法进行分解，将盲源信号分解成一系列的本征模函数，并由分解后得到的全部本征模函数分量构成盲源信号的观测矩阵X；

3)根据盲源信号的观测矩阵X，采用贝叶斯信息准则估算出盲源信号所含的信号源数目；

4)以盲源信号所含的信号源数目作为目标维数，利用稀疏保持投影方法将盲源信号的观测矩阵X进行维数约简，约简至目标维数，构成投影特征矩阵；

5)将投影特征矩阵投入特征矩阵近似联合对角化分离算法，计算出分离矩阵，通过分离矩阵与观测信号矩阵求得盲源信号的分离信号；

6)分别计算历史故障信号和待测轴承的盲源信号的分离信号的时频域特征，构成历史故障信号和待测轴承盲源信号各自的特征矩阵；

7)以历史故障信号的特征矩阵作为分类基准数据，将待测轴承的盲源信号的特征矩阵作为分类识别特征，利用大间隔分布学习分类方法对盲源信号进行分类识别，确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。

2.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤2)具体为：

21)分别对获取到的待测轴承的单通道盲源信号添加幅值为E₀的不同高斯白噪声，得到多个添加不同噪声的盲源噪声信号，分别对各个盲源噪声信号采用完备总体经验模态分解法进行分解；

22)对各个盲源噪声信号采用完备总体经验模态分解法分解后得到的所有分量进行总体平均，并将得到的总体平均分量视为盲源信号的第一阶的本征模函数分量I₁(t)，从盲源信号中将第一阶的本征模函数分量I₁(t)分离出来后得到的剩余信号作为第一阶剩余信号r₁(t)，t表示时间；

23)将幅值为E₁的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W₁(t)；分别将各个噪声的幅值E₁与其第一阶本征模函数分量W₁(t)的积E₁W₁(t)添加到第一阶剩余信号r₁(t)，得到多个不同的第一阶重构信号S^* ₁(t)，并分别对各个第一阶重构信号S^* ₁(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第二阶本征模函数分量I₂(t)，再从第一阶剩余信号r₁(t)中将第二阶本征模函数分量I₂(t)分离出来后得到的剩余信号作为第二阶剩余信号r₂(t)；

24)对于得到的第b阶本征模函数分量和I_b(t)和第b阶剩余信号r_b(t)，b>1，将幅值为E_q的多个不同高斯白噪声分别采用总体经验模态分解法进行分解，得到各个不同噪声的第一阶本征模函数分量W_b(t)；分别将各个噪声的幅值E_b与其第一阶本征模函数分量W_b(t)的积E_bW_b(t)添加到第b阶剩余信号r_b(t)，得到多个不同的第b阶重构信号S^* _b(t)，并分别对各个第b阶重构信号S^* _b(t)用总体经验模态分解法进行分解，分解后再将全部分解得到的分量进行总体平均，并将得到的总体平均分量视为盲源信号的第b+1阶本征模函数分量I_b+1(t)，再从第b阶剩余信号r_b(t)中将第b+1阶本征模函数分量I_b+1(t)分离出来后得到的剩余信号作为第b+1阶剩余信号r_b+1(t)；

25)重复执行步骤24)，直到剩余信号的极值点个数不超过两个时终止分解，由分解后得到的全部本征模函数分量，构成盲源信号的观测矩阵X。

3.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤3)具体为：

31)计算观测矩阵X的相关矩阵R_x：

R_x＝XX^H；

其中H表示复共轭转置；

32)对相关矩阵R_x进行奇异值分解后，得到非零特征值λ₁,…,λ_j,…,λ_L，j∈{1,2,…,L}，L表示非零特征值的个数；对于序号k从1到L的取值，根据贝叶斯信息准则，分别计算各个取值的序号k的BIC代价函数值：

$\mathrm{BIC}\left(k\right)={\left(\underset{j=1}{\overset{k}{\mathrm{\Π}}}{\mathrm{\λ}}_j\right)}^{-N/2}{\mathrm{\σ}}_k^{-N(L-k)/2}{N}^{-(u+k)/2};$

其中，N为盲源信号的采样点数；u＝Lk-k(k+1)/2；

34)根据各个取值的序号k的BIC代价函数值，寻找使BIC代价函数值最大的序号k的取值，则盲源信号所含的信号源数目d的取值即为使BIC代价函数值最大的序号k的取值。

4.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤4)具体为：

41)计算稀疏系数：根据稀疏理论构建目标函数：

${\min }_{s_i}{\left|\right|s_i\left|\right|}_0,s.t.x_i={\mathrm{Xs}}_i$

其中s_i表示盲源信号的观测矩阵X中第i个本征模函数分量x_i由其它本征模函数分量进行稀疏表示时对应的稀疏系数；

42)计算投影矩阵A：根据计算得到的观测矩阵X中各个本征模函数分量对应的稀疏系数，分别对观测矩阵X中的各个本征模函数分量进行稀疏重构，并使其重构误差最小，目标函数如下：

$\min \mathrm{\ϵ}\left(A\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\left|\right|A^T{x}_i-A^T{\mathrm{Xs}}_i\left|\right|}^2=\mathrm{tr}\left\{A^{T}X(I-S-S^T+S^T)X^{T}A\right\};$

其中，S＝[s_1,s_2,,…,s_N]为稀疏系数矩阵，N为盲源信号的采样点数；增加约束A^TXX^TA＝I，I为单位矩阵，并求解投影矩阵A的最大值：

maxε(A)＝A^TXS_βX^TA/A^TXX^TA；

其中S_β＝S+S^T-S^TS，将上式转化为XS_βX^TA＝λXX^TA，求解广义特征值和特征向量，取出其中最大的d个特征值λ对应的特征向量a_1,a_2,,...,a_d，得到投影矩阵A＝[a_1,a_2,,...,a_d]，d为盲源信号所含的信号源数目；

43)将观测矩阵X在投影矩阵A上投影，得到一个d×N维的投影特征矩阵X₀。

5.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤5)具体为：

51)计算投影特征矩阵X₀的自相关函数R_XR：

R_XR＝E(X₀X₀ ^*)；

其中上标*表示复共轭；对自相关函数R_XR进行广义特征值分解，得到的d个最大的特征根和相应的特征向量分别用λ₁,…,λ_d和h₁,…,h_d表示；

52)对投影特征矩阵X₀进行预白化处理：

白化信号z(t)为：

z(t)＝[z₁(t),z₂(t),…,z_i(t),…,z_d(t)]^T；

其中1≤i≤d，d为盲源信号所含的信号源数目；

白化矩阵W为：

$W=[{\left({\mathrm{\λ}}_1\right)}^{-\frac{1}{2}}{h}_1,...,{\left({\mathrm{\λ}}_d\right)}^{-\frac{1}{2}}{h}_d{]}^H;$

其中，上标H表示复共轭转置；

由盲源分离模型X₀＝BS_s可得：

z＝WX₀＝WBS_s＝US_s；

其中，B为混合矩阵，S_s为分离信号矩阵，U为待求解的酉矩阵；

53)联合对角化白化后的观测信号：

对于d阶的一个任意矩阵M，z(t)的四阶累积向量矩阵的第o行第p列元素[Q_z(M)]_op为：

${\left[Q_z\left(M\right)\right]}_{\mathrm{op}}\stackrel{\mathrm{def}}{=}\underset{q=1}{\overset{d}{\mathrm{\Σ}}}\underset{r=1}{\overset{d}{\mathrm{\Σ}}}\mathrm{cum}(z_o,{z_p}^{*},z_q,{z_r}^{*})m_{\mathrm{qr}};$

式中，um(z_o,z_p ^*,z_q,z_r ^*)表示白化信号z(t)中第o个分量z_o、第p个分量z_p、第q个分量z_q、第r个分量z_r的四阶累积量矩阵，1≤o≤d，1≤p≤d，1≤q≤d，1≤r≤d；m_qr为矩阵M的第q行第r列元素；

54)用U^TQ_z(M_i)U中对角元素的平方和作为对角化过程参照函数，即：

$F\left(U\right)=\underset{M_i\∈M}{\mathrm{\Σ}}{\left|\right|\mathrm{diag}(U^T{Q}_z\left(M_i\right)U\left|\right|}^2;$

通过优化算法使得参照函数F(U)最小化获得酉矩阵U的解，即获得盲源信号的分离信号S_s：

S_s＝U^TWX₀。

6.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤6)具体为：

针对d×N维的盲源信号的分离信号S_s，d为盲源信号所含的信号源数目，N为盲源信号的采样点数，通过计算盲源信号每周期采样点数R将分离信号S_s的每个分离信号分割为g个周期数据y_i，i＝1,2,…,g，在时域分别按周期提取分离信号S_s的峰值、平均幅值、标准差、均方根值、峭度指标、峰值指标、波形指标、裕度指标8个时域特征来反映信号变化；其中：

峰值为P＝max|y_i|；     平均幅值为

标准差为 $T=\sqrt{\frac{1}{R-1}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}|y_i-\stackrel{\‾}{y_i}|};$    均方根值为 $y_{\mathrm{rms}}=\sqrt{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{y_i}^2};$

峭度指标为 $K=\frac{\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{\left(y_i\right)}^4}{{y_{\mathrm{rms}}}^4};$     峰值指标为 $C_f=\frac{P}{{y_{\mathrm{rms}}}^4};$

波形指标为 $S_f=\frac{y_{\mathrm{rms}}}{\left|\stackrel{\‾}{y_i}\right|};$        裕度指标为 ${\mathrm{CL}}_f=\frac{P}{{\left(\frac{1}{R}\underset{i=1}{\overset{R}{\mathrm{\Σ}}}\sqrt{\left|y_i\right|}\right)}^2};$

然后，采用db4小波包函数对分离信号S_s进行3层正交小波包分解用以进行频域特征提取，得到分离信号S_s对应的由全频带均匀划分8个子频带的滤波信号，通过重构每个子频带的节点小波包系数，确保重构信号和原分离信号S_s的长度一样，然后计算每个子频带滤波信号的各个采样点的幅值平方和作为其能量，将各个子频带滤波信号的能量与全频带重构信号总能量的比值作为其频域统计特征，即(E₀/E),(E₁/E),...,(E₇/E)，其中E₀,E₁,...,E₇分别为8个子频带滤波信号的能量，E为全频带重构信号的总能量；由此，针对分离信号S_s，共构造了8个时域特征和8个频域统计特征，由该16个特征构成盲源信号的分离信号S_s的特征矩阵；

针对待测轴承盲源信号的分离信号，构造上述的8个时域特征和8个频域统计特征，构成待测轴承盲源信号的特征矩阵；

所述历史轴承故障数据中的每个故障信号数据，是对已知故障工况类型的轴承进行振动加速度信号采样后得到的振动加速度信号，并标记有其所属故障工况的类别标签；根据历史轴承故障数据包含的各种不同故障工况，从中分别在每种故障工况下选取多个故障信号数据，分别作为历史轴承故障信号样本，并针对每个历史轴承故障信号样本求取其分离信号，并同样构造上述的8个时域特征和8个频域统计特征，构成每个历史轴承故障信号样本的特征矩阵。

7.根据权利要求1所述基于稀疏表示和大间隔分布学习的轴承故障分类诊断方法，其特征在于，所述步骤7)具体为：

71)将待测轴承的盲源信号的特征矩阵作为测试样本，将从历史轴承故障数据中选取的全部历史轴承故障信号样本的特征矩阵作为训练样本，训练样本集表示为：

V＝{(v₁,c₁),(v₂,c₂),…,(v_i,c_i),…,(v_m,c_m)}；

其中，m为训练样本集中训练样本的个数，v_i表示第i个训练样本的特征矩阵，c_i表示第i个训练样本的类别标签；

72)利用一阶和二阶统计量，即均值和方差来描述间隔分布：

$\stackrel{\‾}{\mathrm{\γ}}=\frac{1}{m}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{c}_i{w}^T\mathrm{\φ}\left(c_i\right)=\frac{1}{m}{\left(V_{\mathrm{\φ}}c\right)}^{T}w;$

$\widehat{\mathrm{\γ}}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\underset{j=1}{\overset{m}{\mathrm{\Σ}}}{(c_i{w}^T\mathrm{\φ}\left(v_i\right)-c_j{w}^T\mathrm{\φ}\left(v_j\right))}^2=\frac{2}{m^2}\left({\mathrm{mw}}^T{V}_{\mathrm{\φ}}{\mathrm{cc}}^T{V_{\mathrm{\φ}}}^{T}w\right);$

式中，V_φ＝[φ(v₁),φ(v₂),…,φ(v_i),…,φ(v_m)]，列向量c＝[c₁,c₂,…,c_i,…,c_m]^T；大间隔分布学习模型的线性系数向量σ＝[σ₁,σ₂,…,σ_m]为线性系数；φ(v_i)和φ(b_j)分别为第i个训练样本的特征矩阵v_i和第j个训练样本的特征矩阵v_j在核函数K上的特征映射，即有：

K(v_i,v_j)＝φ(v_i)^Tφ(v_j)；

其中，i＝1,2,…,m，j＝1,2,…,m，且i≠j；

在训练样本可以实现零误差的精确区分情况下，同时最大化间隔均值和最小化间隔方差引出：

$\begin{array}{cc}\underset{w}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}\end{array};$

s.t.c_iw^Tφ(v_i)≥1，i＝1,2,…,m；

式中：β₁和β₂分别为均值、方差的权值参数；

在训练样本不能实现零误差精确区分情况下，引出：

$\begin{array}{cc}\underset{w,\mathrm{\ϵ}}{\min }& \frac{1}{2}{w}^{T}w+{\mathrm{\β}}_1\widehat{\mathrm{\γ}}-{\mathrm{\β}}_2\stackrel{\‾}{\mathrm{\γ}}\end{array}+C_c\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\ϵ}}_i;$

s.t.c_iw^Tφ(v_i)≥1-ε_i，ε_i≥0，i＝1,2,…,m；

C_c为惩罚参数，ε＝[ε₁,ε₂,…,ε_m]^T为松弛变量；

把带入上式，根据拉格朗日乘法确定系数σ＝[σ₁,σ₂,…,σ_i,…,σ_m]；

73)训练样本集中包含的各种不同故障工况的训练样本，针对其中第z类故障工况的训练样本看作正类，z∈{1,2,…,Z}，Z表示训练样本集中包含的故障工况的类别总数，将其它Z-1类故障工况的训练样本看作负类；将正类训练样本的类别标签取值为1，正类训练样本的类别标签取值为-1，然后按下式计算输入量α在第z类故障工况下对应的故障类型预测函数：

$\mathrm{sgn}\left(w^T\mathrm{\φ}\left(\mathrm{\α}\right)\right)=\mathrm{sgn}\left(\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_{i}K(v_i,\mathrm{\α})\right);$

其中，sgn(·)表示符号函数；

74)将测试样本分别作为各种工况故障工况下对应的故障类型预测函数的输入量，计算出测试样本作为输入量在各个故障工况类别下对应的故障类型预测函数值，将故障类型预测函数值为1时对应的故障工况类别判定为测试样本的故障工况类别，从而确定待测轴承的盲源信号所属的故障工况类别，作为待测轴承故障类型的分类诊断结果。