CN106769052B - 一种基于聚类分析的机械系统滚动轴承智能故障诊断方法 - Google Patents

Abstract

本发明为一种基于聚类分析的机械系统滚动轴承智能故障诊断方法，先训练诊断模型：包括采集外圈、内圈、滚动体和保持架故障及正常五种轴承状态标准振动信号样本；信号分解，原始振动信号和分解后分量的时域与频域特征提取，得原始特征集；自权重算法和AP聚类去除冗余得Z个优选特征；再次用AP聚类，完成样本状态分类，得到训练好的诊断模型。故障诊断时采集轴承实时振动信息，信号分解，提取模型确定的优选特征，导入AP按诊断模型参数聚类，与模型中已知的Z个特征比较，得到当前未知信号类别，完成故障诊断。本发明同时使用EEMD与WPT分解振动信号，获取更多精细的轴承状态信息，自权重法和AP聚类增加诊断的智能化，保证精确诊断。

Description

一种基于聚类分析的机械系统滚动轴承智能故障诊断方法

技术领域

本发明属于装备系统故障监测、诊断领域，具体涉及一种基于聚类分析的机械系统滚动轴承智能故障诊断方法。

背景技术

随着科技的发展和社会的进步，各类型机械装备在工程中得到了广泛应用。滚动轴承是机械装备中的关键部件也是使用较多的部件。在运转过程中多种原因可能会造成滚动轴承损坏，如装配不当、润滑不良、水分和异物侵入、腐蚀和过载等都可能会导致滚动轴承过早损坏。即使在安装、润滑和使用维护都正常的情况下，经过一段时间运转，滚动轴承也会出现疲劳剥落、磨损、点蚀等故障，导致轴承不能正常工作。一般的轴承故障分为外圈故障、内圈故障、滚动体故障、保持架故障等类别。滚动轴承是机械装备中的重要零部件，一旦发生故障，有可能造成财产损失，严重的甚至会造成人员伤亡。所以对轴承故障的及时发现和诊断是非常重要的。

传统的滚动轴承诊断方式，在诊断对象复杂时通常需要专家人工干预，不仅增加了人力成本，而且容易出现漏诊和误诊。因此，需要发展以人工智能为核心的智能故障诊断方法，辅助一般的操作人员，即使他们不掌握系统原理，也可借助人工智能进行可靠的决策。

现有的智能故障诊断技术大多在诊断时只用单一的指标或几个指标组合来判断故障类型。但在对复杂对象的故障进行识别时，少数几个指标有时并不能完全表现故障特点，因此故障诊断的准确率低。

虽然随着信号处理和特征提取技术的快速发展，越来越多的特征可以参与到对故障的诊断之中，但现有的诊断方法难以对信号分类提取故障信息特征，增加的特征使识别运算效率低，无法完成实时的故障诊断。

因此需要一种改进的机械系统滚动轴承智能故障诊断方法，能够采用更多的信号特征，且运算效率高，故障诊断的准确率高。

发明内容

本发明目的在于提供一种基于聚类分析的机械系统滚动轴承智能故障诊断方法，该方法通过采集滚动轴承的振动信号；采用EEMD和WPT算法分解振动信号得到更精细轴承状态信息；提取时域，频域特征，组合成为原始特征矩阵；将原始特征矩阵输入自适应特征选择算法剔除非敏感和冗余的特征，从而得到优选特征矩阵；最后优选特征矩阵输入AP算法，得到故障聚类结果，完成识别不同故障的诊断模型的训练，实际使用时，当前采集的未知振动信号导入完成训练的诊断模型，用AP算法与其匹配，诊断当前机械系统的轴承状态。

本发明设计的一种基于聚类分析的机械系统滚动轴承智能故障诊断方法包括两大部分，第一部分是训练诊断模型，第二部分是采用训练好的诊断模型对机械系统中的滚动轴承进行故障诊断。

第一部分训练诊断模型的主要步骤如下：

第一步：振动信号采集

机械系统中的关键部件待检测的滚动轴承的轴承座上安装振动传感器，采集机械装备的振动信号x(t)。振动传感器的信号线连接上位机，各轴承实时振动信息储存于上位机中为原始振动信号。所述振动传感器为压电加速度传感器、位移传感器和电涡流位移传感器中的任一种。

本发明轴承故障分为外圈故障、内圈故障、滚动体故障、保持架故障，还有正常状态，即无故障状态，C＝5是轴承状态类别数，包括4种故障状态以及一种正常状态。在本训练诊断模型部分，采集外圈故障、内圈故障、滚动体故障、保持架故障和正常状态的五种轴承状态标准振动信号，上位机采集各振动信号后，将每种振动信号分成信号片段，每一段信号作为一个样本，Num是每种轴承状态的样本个数；每种轴承状态都有Num个样本，样本总数T＝Num*C。

第二步：信号分解

上位机对采集到的信号x(t)进行集合经验模式分解(EEMD)得到各个频带的本征模态函数(IMF)，即得到多个本征模式分量(IMF分量)。集合经验模式分解EEMD中设置的参数包括添加高斯白噪声标准差(Nstd)和添加噪声的次数(Ensemble number,NE)。本发明Nstd值设置为信号标准差的0.1～0.2倍；NE值设置为信号点数的0.08～0.15。由于集合经验模式分解是自适应分解，得到的本征模式分量的具体个数由当前采集信号决定。

上位机还对采集到的信号x(t)进行小波包分解(WPT)得到分解后的子频带。本发明的小波包分解采用DB(Daubechies)系列小波基函数。DB系列小波基函数，具有N阶消失矩，又具有紧支撑性的正交小波基，并且是近似对称的，满足本发明小波包分解所需的紧支撑性和至少1阶以上的消失矩，以及连续可微性和对称性，可有效地提取振动信号的故障特征。

本发明的小波包分解采用5～10阶消失矩的DB系列小波基函数。

小波分解层数L为2～10层。根据采样频率大小确定分析频率范围，从而确定小波分解层数。

第三步：特征提取

上位机分别对第一步采集的原始振动信号和第二步集合经验模式分解后的本征模式分量进行时域与频域特征提取。

分别提取如下14个特征：

31、均值(Mean)F_m:

式中x(n)是信号的时域序列，n＝1,2,……N；N为振动样本点数；实际采集到的是连续信号x(t)的离散点，用x(n)表示采集的离散信号的时域序列，x(n)相当于x(t)，其中n＝1,2，……N。

32、标准差(Standard deviation)F_std:

33、均方根值(Root mean square)F_rms:

34、峰值指标(Peak)F_p:

F_p＝max|x(n)|

35、波形指标(Shape Factor)F_sf:

36、峰值因子(Crest Factor)P_cf:

P_cf＝F_p/F_rms

37、脉冲指标(Impulse Factor)F_if:

38、裕度指标(Clearance Factor)F_clf:

39、峭度(Kurtosis)F_k:

310、偏斜度(Skewness)F_s:

311、频率均值(Mean frequency)F_mf:

式中s(k)是信号x(n)的频谱，k＝1,2,……K；K是谱线数；

312、中心频率(Center frequency)F_cf:

式中f_k是第k条谱线的频率值。

313、频率均方根(Root mean square frequency)F_rmfs:

314、频率标准偏差(Standard deviation frequency)F_stdf:

对原始振动信号提取上述时域与频域指标，得到14个特征，设为特征集1；

第二步原始振动信号集合经验模式分解后得到a个本征模式分量，对每个本征模式分量分别提取上述时域与频域指标，得到14×a个特征，设为特征集2；

第二步原始振动信号小波包分解后得到B个子频带，分解L层中第l层的第j个子频带能量用公式(1)表示：

是第j个子频带的成分，提取每个子频带能量作为特征得到B个特征，设为特征集3。

以上三个特征集综合为原始特征集X，其中特征总数为14+14×a+B。

第四步：优选特征

在第三步中所得原始信号和其集合经验模式分解后的本征模式分量信号进行特征提取，以及原始信号小波包分解后的子频带能量特征构成了原始特征集X，特征数量巨大，其中包含了一定量的冗余特征与非敏感特征，因上述所列出的14种时频域特征中，存在着对故障不敏感的特征。如果人工一一检查每个特征的敏感性，显然费时费力。本发明提出了一种自适应特征选择方法，对每个特征进行评价，用于去除原始特征集X中的冗余特征与非敏感特征，挑选出最优特征。本自适应特征选择方法不需要先验知识参照，而是根据特征本身的特性来决定该特征是否被去除。该方法包含两个主要步骤，一是敏感特征选择，二是冗余特征去除。

41、敏感特征选择

假设包含了C种轴承状态类别的特征矩阵为

{f_g，vg＝1，2，...，Num；i＝1，2，...，I}_C，

本发明样本总数T＝Num*C，为上述原始特征矩阵中的行数；I是特征种数，是上述原始特征矩阵中的列数，由第三步可知I＝14+14×a+B；从各样本得到的特征总个数为Num×C×I＝T×I，为上述原始特征矩阵的元素总数；f_g，i是第n个样本的第i个特征值。

所有的特征集合定义为{v_t，i}t＝1，2，...，T；i＝1，2，...，I自权重算法具体过程如下：

411、v_t，i进行正规化，得到正规化后的特征集{x_t，i}

412、计算自相似因子(self-similarity factor)SF_mt，i

SF_mt，i＝||x_m，i-x_t，i||²，m，t＝1，2...，T；i＝1，2，...，I (3)

413、对每个特征建立权重矩阵

414、计算每个特征的自权重Sw_i

根据步骤413得到的权重，计算每个特征的自权重Sw_i

Sw_i＝mean(W_i) (5)

显然，权重矩阵W_i和自权重值Sw_i是根据每个特征X_t，i的特点自适应得到的，计算完成后，每个特征都有一个与之对应的自权重值Sw_i。不同轴承状态样本之间的差异性越大，则得到的特征的自权重值Sw_i越大，故用特征的自权重值可以评价该特征的故障样本差异性，故障样本差异大的特征的自权重值Sw_i大。将所有特征的自权重值Sw_i从大到小排列，选出前R个自权重相对较大的特征作为去冗余的基础，此R个特征为敏感特征。原始特征矩阵{x_t，i}中选出了R个特征后，得到敏感特征矩阵{x_t，r}。R值由用户自己视具体振动信号的情况根据经验决定，R推荐取值为特征总数的20％～30％。

42、冗余特征去除

步骤41的自权重算法只选取了R个敏感特征，并没有考虑冗余特征的问题，所以得到的敏感特征矩阵{x_t，r}中仍然包含冗余特征。为了进一步降低特征矩阵维数，提高算法效率，本发明使用AP聚类算法去除冗余特征。一般情况下，AP聚类算法是按照样本数划分类别的，本发明将特征矩阵转置，把特征作为样本输入AP聚类算法，得到特征之间的聚类结果，用聚类所得的每一类别的聚类中心作为此类特征的代表，此类中的其他特征作为冗余特征。聚类后所得的各类中最有代表性特征组成优选特征集合。利用AP去冗余步骤如下：

421、敏感特征矩阵的转置

步骤41自权重选出的敏感特征矩阵{x_t，r}＝X_t，r

化为转置阵

T是样本总个数，R是在步骤41中被选出的敏感特征个数。

422、AP算法的偏向参数

用欧式距离建立X_r，t的相似矩阵Sim₄(Similarity)，由公式7设置AP算法中的偏向参数P₄。

Sim₄＝-||X_r，t-X_r+i，t||² (6)

P₄＝2·median(Sim₄) (7)

423、优选特征矩阵X_t，z

将X_r，t导入AP聚类算法，自动为每个特征打上类标签，聚类得到Z个类别，则选取每一类的聚类中心作为最终的优选特征，得到优选特征矩阵X_t，z。

第五步：故障识别

第四步所得的优选特征矩阵中，每一个特征都是和轴承状态息息相关的，大量冗余的特征已经被去除。本步骤将第四步所得优选特征矩阵再次输入AP聚类算法，自动识别不同轴承状态。

AP聚类算法不同于其他一般的聚类算法，AP聚类算法同时把所有的数据点当作潜在的聚类中心。数据点之间的信息实值在整个网络边缘传输，直到一个良好的聚类中心和其对应集群出现为止。AP算法将每个样本间的相似度之和最大化，使得聚类中心逐渐出现，见公式(8)：

h是数据点p的潜在聚类中心，Sim(p,h)是数据点p和h之间的相似度值。AP算法中传递两种类型的消息，即吸引度(Responsibility)和归属度(Availability)。吸引度re(p,h)表示从点p发送到候选聚类中心h的数值消息，反映h点是否适合作为p点的聚类中心；归属度av(p,h)表示从候选聚类中心h发送到p的数值消息，反映p点是否选择h作为其聚类中心。re(p,h)与a av(p,h)越大，则h点作为聚类中心的可能性就越大，并且p点隶属于以h点为聚类中心的聚类的可能性也越大。

将第四步获得的优选特征矩阵X_t，z导入AP聚类算法中，识别故障的主要步骤如下：

51、建立相似矩阵Sim₅；相似度Sim₅度量采用欧氏距离，偏向参数P₅设置为Sim₅均值的2倍，如下公式所示。

建立T个数据点的相似矩阵Sim₅，设置偏向参数P₅的值；

Sim₅＝-||X_t，z-X_t+1，z||²，t＝1，2，...，T；z＝1，2，...，Z (9)

P₅＝2·media(Sim₅) (10)

52、初始化归属度矩阵av＝0；

53、更新吸引度的值re(p,h)

(11)式中h′表示另外一个数据点。h′∈{1，2，...，T}，但h′≠h。按(11)式循环，直至得到{av(p，h′)+Sim₅(p，h′)}的最大值。

为了避免在迭代过程中出现数值震荡的情况，在每次re更新时加入一个阻尼系数λ，0≤λ≤1，re按下式更新

re_p＝(1-λ)re_p+λre_p-1 (12)

54、更新归属度的值av(p，h)

(13)式中p′表示另外一个数据点。p′∈{1，2，...，T}，但p′≠h同时p′≠p；

re(h，h)＝P(h)-max{av(p，h)+sim₅(p，h)}

为了避免在迭代过程中出现数值震荡的情况，在每次av更新时加入一个相同的阻尼系数λ，av按下式更新，

av_p＝(1-λ)av_p+λav_p-1 (14)

55、确定潜在聚类中心

同时考虑吸引度和归属度的值，当满足av(h，h)+re(h，h)＞0时，h点被选择成为一个潜在聚类中心；

56、如果在数次迭代后，聚类中心点仍然不变，那么进入步骤57；否则，返回步骤53；

57、使用最邻近规则为剩余的数据点安排归属类，同一类中每个数据点的相似度均为最大值。

第四步用AP算法识别出聚类中心有多少个，再将剩余的样本点归属于与之最近的聚类中心，选择得到Z个优选特征。本步再次使用AP算法识别故障，即完成分类。

消息传递过程的终止是达到了用户设置的固定迭代次数。

本发明设置迭代次数为800～1200。

经上述五个步骤完成本发明诊断模型训练，确定各步骤的集合经验模式(EEMD)分解中的添加高斯白噪声标准差(Nstd)值，添加噪声的次数(NE)，本征模态函数(IMF)个数a；小波包分解所采用的小波基函数，子频带个数B，分解层数L；自权重法选择的敏感特征个数R，优选特征个数Z以及AP算法中的相似度，偏向参数，阻尼系数和迭代次数，训练好的故障样本和正常样本均保留在模型中，作为用于对某机械系统的滚动轴承故障进行识别的训练好的诊断模型存储于上位机。

第二部分、故障诊断

第1步：振动信号采集

采用第一步的振动传感器安装于机械系统中的每个待检测的滚动轴承的轴承座上，采集机械装备的振动信号x(t)。振动传感器的信号线连接上位机，各轴承实时振动信号储存于上位机中为未知振动信号。

高速重载轴承、多受冲击的轴承或暴露在高温、高压、易腐蚀环境中的轴承为恶劣工况下易发生故障的轴承，对此类轴承需要对同一部位定期进行振动检测，并按时间先后进行比较，根据振动状态的变化，粗略判断该轴承是否需要进一步精细检测，即采用本方法检测。

第2步：信号分解

上位机取完成训练的诊断模型中确定的参数对采集到的未知信号x(t)进行集合经验模式分解(EEMD)得到各个频带的本征模态函数(IMF)，并进行小波包分解(WPT)得到分解后的B个子频带。

第3步：特征提取

在未知振动信号、第2步分解所得的本征模态函数以及B个子频带中提取按完成训练的诊断模型确定的Z个优选特征；

第4步：故障诊断

第3步提取的Z个优选特征，导入AP算法，按完成训练的诊断模型中AP算法中的相似度，偏向参数，阻尼系数和迭代次数进行聚类，未知信号的优选特征与诊断模型中已知的Z个特征比较，根据未知信号的优选特征归属于完成训练的诊断模型中的哪个聚类中心，即判断得到当前未知信号的类别，即完成故障识别诊断。

与现有技术相比，本发明一种基于聚类分析的机械系统滚动轴承智能故障诊断方法的优点为：1、同时使用EEMD与WPT技术分解原始振动信号，再对其分解后分量进行大量的特征提取，获取了更多精细的轴承状态信息，参与到对故障的诊断之中；2、利用自权重算法选择对故障分类有用的敏感特征，降低特征矩阵维度，减小模式识别算法的负担，保证精确诊断；3、使用AP聚类算法去除冗余特征，大量不相关特征自动去除，增加了诊断模型的智能化程度，不仅减小后续智能算法的负担，而且增加了轴承状态识别精确度；4、识别核心采用AP聚类(仿射传播聚类AP cluster)，属于无监督聚类算法，快速而准确地识别不同种类的故障；5、适用于大型复杂机械系统中滚动轴承的在线或离线智能故障识别，减少人力成本，避免突发性事故发生，减小经济损失；6、本方法一般操作人员均可掌握，无需理解系统原理即可完成准确可靠的故障诊断，提高故障诊断的自动化及智能化程度。

附图说明

图1为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例流程图；

图2为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例第一步不同类别的原始振动信号及其频谱图；

图3为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例第二步中集合经验模式分解所得的本征模式分量图；

图4为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例第二步中小波包分解所得的子频带图；

图5为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例第四步的优选特征结果示意图；

图6为本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例第五步中智能诊断结果示意图。

具体实施方式

下面结合附图和具体实施例对本发明的内容作进一步详细说明。

本实施例采用本发明基于聚类分析的机械系统滚动轴承智能故障诊断方法，根据单一轴承的振动信号，对外圈故障、内圈故障、滚动体故障、保持架故障的振动信号和正常轴承的振动信号进行分析，最后给出诊断结果。

本基于聚类分析的机械系统滚动轴承智能故障诊断方法实施例如图1所示，包括两大部分，第一部分是训练诊断模型，第二部分是采用训练好的诊断模型对机械系统中的滚动轴承进行故障诊断。

第一部分训练诊断模型主要步骤如下：

第一步：振动信号采集

压电式加速度传感器安装于待测轴承的轴承座上，采集滚动轴承振动信号。首先采集包括外圈故障、内圈故障、滚动体故障、保持架故障的标准振动信号以及正常状态的标准振动信号。本实施例采用的是LW149536型压电式加速度传感器，其特性参数如表1所示。各类振动信号及其频谱如图2所示。

表1 LW149536型压电式加速度传感器特性参数表

上位机与压电式加速度传感器相连接，采集振动信号。

上位机先将采集的原始振动信号分成信号片段，本例每个信号片段点数为4096。

本例轴承故障分为外圈故障、内圈故障、滚动体故障、保持架故障以及正常状态，C是轴承状态类别数，C＝5。本例每种轴承状态类别取50段振动信号片段作为样本，每种轴承状态类别的样本个数Num＝50。本例5类轴承状态的样本总数为T＝Num×C＝250。

第二步：信号分解

本例上位机对采集到的信号x(t)进行集合经验模式分解(EEMD)得到各个频带的本征模态函数(IMF)，即得到a＝8个本征模式分量。

上位机还对采集到的信号x(t)进行小波包分解(WPT)，本例采用10阶消失矩的DB(Daubechies)系列小波基函数。小波包分解每一段信号后得到各个子频带。

本实施例集合经验模式分解、小波包分解以及自权重算法的参数设置如下表2所示。

表2 EEMD、WPT分解以及自权重算法输入参数表

集合经验模式分解所得的本征模式分量如图3所示，小波包分解所得的子频带图如图4所示。

第三步：特征提取

上位机分别对第一步采集的原始振动信号按表3提取14种时域与频域特征，作为特征集1。

对集合经验模式分解后的8个本征模式分量按表3各提取14种时域与频域特征，共得14×8个特征，作为特征集2。

并提取第二步原始振动信号小波包分解后得到8个子频带能量作为特征，分解3层的第j个子频带能量按下式计算

是第j个子频带的成分，共得到B个特征，设为特征集3。

本例三个特征集合成的原始特征集X的特征总数I＝14+14×8+8＝134。

表3时域与频域特征指标

第四步：优选特征

在第三步提取了描述信号状态的大量特征，但其敏感性大有不同。

41、敏感特征选择

本例包含了5种轴承状态类别的特征矩阵为

{f_n，i，n＝1，2，...，Num；i＝1，2，...，I}_C，

本例样本总数T＝Num×C＝250，为上述原始特征矩阵中的行数；特征总数I＝134，是上述原始特征矩阵中的列数；特征总个数即为上述原始特征矩阵的元素总数，为Num×C×I＝T×I＝250×134；f_n，i是第n个样本的第i个特征值。

将所有的特征集合，也就是250×134的特征矩阵，定义为{v_t，i}t＝1，2，...，T；i＝1，2，...，I。自权重算法具体过程如下：

411、v_t，i进行正规化，得到正规化后的特征集{x_t，i}

412、计算自相似因子(self-similarity factor)SF_mt，i

SF_mt，i＝||x_m，i-x_t，i||²，m，t＝1，2...，T；i＝1，2，...，I (3)

413、对每个特征建立权重矩阵

414、计算每个特征的自权重Sw_i

根据步骤413得到的权重，计算每个特征的自权重Sw_i

Sw_i＝mean(W_i) (5)

本例所有特征的自权重值Sw_i从大到小排列，选出前R＝34个自权重较大的特征作为敏感特征，得到敏感特征矩阵{X_t，r}。

42、冗余特征去除

421、敏感特征矩阵的转置

步骤41自权重选出的敏感特征矩阵{x_t，r}＝X_t，r

化为转置阵

T是特征样本总个数，R＝34，是在步骤41中被选出的敏感特征个数。

422、AP算法的偏向参数

用欧式距离建立X_r，t的相似矩阵Sim₄(Similarity)，由公式7设置AP算法中的偏向参数P₄。

Sim₄＝-||X_r，t-X_r+i，t||² (6)

P₄＝2·median(Sim₄) (7)

423、优选特征矩阵X_t，z

将X_r，t导入AP聚类算法，自动为每个特征打上类标签，本例聚类得到Z＝11个类别，则选取每一类的聚类中心作为最终的优选特征，得到优选特征矩阵X_t，z。

经过层层筛选，无用、冗余的特征已经被去除，剩下的特征均为对5种轴承状态具有良好区分度且无冗余的优选特征值。所得优选特征结果如图5所示，图5中的○即为自适应特征选择算法选出的11个优选特征。

第五步：智能故障识别

将第四步得到的优选特征矩阵X_t，z，导入AP聚类算法，识别不同种类的故障。主要步骤如下：

51、建立相似矩阵Sim₅

建立T个数据点的相似矩阵Sim₅，设置偏向参数P₅的值；

Sim₅＝-||X_t，z-X_t+1，z||²，t＝1，2，...，T；z＝1，2，...，Z (9)

P₅＝2·media(Sim₅) (10)

52、初始化归属度矩阵av＝0；

53、更新吸引度的值re(p，h)

(11)式中h′表示另外一个数据点。h′∈{1，2，...，T}，但h′≠h。按(11)式循环，直至得到{av(p，h′)+Sim₅(p，h′)}的最大值。

在每次re更新时加入一个阻尼系数λ，本例λ＝0.8，re按下式更新

re_p＝(1-λ)re_p+λre_p-1 (12)

54、更新归属度的值av(p，h)

(13)式中p′表示另外一个数据点。p′∈{1，2，...，T}，但p′≠h，同时p′≠p；

re(h，h)＝P(h)-max{αv(p，h)+sim₅(p，h)}

在每次av更新时本例也加入一个阻尼系数λ＝0.8，av按下式更新，

av_p＝(1-λ)av_p+λav_p-1 (14)

55、确定潜在聚类中心

当满足av(h，h)+re(h，h)＞0时，h点被选择成为一个潜在聚类中心；

56、如果在数次迭代后，聚类中心点仍然不变，那么进入步骤57；否则，返回步骤53；

57、使用最邻近规则为剩余的数据点安排归属类，同一类中每个数据点的相似度均为最大值。

本例设置迭代次数为1000。

本例模型训练结果如图6所示，用主成分分析法把聚类结果可视化，把Z个优选特征转化为2个主成分，即第一主成分与第二主成分，其中每个主成分都能够反映原始变量的大部分信息，且所含信息互不重复。图6的x轴与y轴分别为第一主分量与第二主分量，展示最后聚类效果，不同种类的信号被识别，包括内圈故障、外圈故障、滚动体故障、保持架故障以及正常信号5种类别，分类正确率达到了99.3％。

本例上述五个步骤完成模型训练，确定各步骤的集合经验模式(EEMD)分解中的添加高斯白噪声标准差值Nstd＝02，添加噪声的次数NE＝200，本征模态函数(IMF)个数a＝8；小波包分解采用的10阶消失矩的DB系列小波基函数，子频带个数B＝8，分解层数L＝3；自权重法选择的敏感特征个数R＝32，优选特征个数Z＝11以及AP算法中的相似度，偏向参数P为相似矩阵Sim均值的二倍，阻尼系数λ＝0.8和迭代次数为1000，训练好的故障样本和正常样本均保留在模型中，作为用于对某机械系统的滚动轴承故障进行识别的训练好的诊断模型存储于上位机。

第二部分、故障诊断

第1步：振动信号采集

采用与第一步相同的压电式加速度传感器安装于机械系统中的每个待检测的滚动轴承的轴承座上，采集机械装备的振动信号x(t)。振动传感器的信号线连接上位机，各轴承实时振动信息储存于上位机中为未知振动信号。

第2步：信号分解

上位机取完成训练的诊断模型中确定的参数对采集到的未知信号x(t)进行集合经验模式分解(EEMD)得到各个频带的本征模态函数(IMF)8个，并进行小波包分解(WPT)得到分解后的B＝8个子频带。

第3步：特征提取

在未知振动信号、第2步分解所得的本征模态函数以及B个子频带中提取按完成训练的诊断模型确定的Z＝11个优选特征；

第4步：故障诊断

第3步提取的Z＝8个优选特征，导入AP算法，按完成训练的诊断模型中AP算法中的相似度，偏向参数，阻尼系数和迭代次数进行聚类，未知信号的优选特征与诊断模型中已知的Z个特征比较，根据未知信号的优选特征归属于完成训练的诊断模型中的哪个聚类中心，即判断得到当前未知信号的类别，即完成故障识别诊断。

上述实施例，仅为对本发明的目的、技术方案和有益效果进一步详细说明的具体个例，本发明并非限定于此。凡在本发明的公开的范围之内所做的任何修改、等同替换、改进等，均包含在本发明的保护范围之内。

Claims (6)

1.一种基于聚类分析的机械系统滚动轴承智能故障诊断方法，包括两大部分，第一部分是训练诊断模型，第二部分是采用训练好的诊断模型对机械系统中的滚动轴承进行故障诊断；

第一部分训练诊断模型的主要步骤如下：

第一步：振动信号采集

机械系统中的关键部件待检测的滚动轴承的轴承座上安装振动传感器，采集机械装备的振动信号x(t)；振动传感器的信号线连接上位机，各轴承实时振动信息储存于上位机中为原始振动信号；

本发明轴承故障分为外圈故障、内圈故障、滚动体故障、保持架故障，还有正常状态；C＝5是轴承状态类别数，包括4种故障状态以及一种正常状态；在本训练诊断模型部分，采集外圈故障、内圈故障、滚动体故障、保持架故障和正常状态的五种轴承状态标准振动信号，上位机采集各振动信号后，将每种振动信号分成信号片段，每一段信号作为一个样本；每种轴承状态都有Num个样本，样本总数T＝Num*C；

第二步：信号分解

上位机对采集到的信号x(t)进行集合经验模式分解得到各个频带的本征模态函数，即得到多个本征模式分量；集合经验模式分解中设置的参数包括添加高斯白噪声标准差和添加噪声的次数；

上位机还对采集到的信号x(t)进行小波包分解得到分解后的子频带；小波分解层数L为2～10层；

第三步：特征提取

上位机分别对第一步采集的原始振动信号和第二步集合经验模式分解后的本征模式分量进行时域与频域特征提取；

分别提取如下14个特征：

31、均值F_m:

式中x(n)是信号的时域序列，n＝1,2,……N；N为振动样本点数；实际采集到的是连续信号x(t)的离散点，用x(n)表示采集的离散信号的时域序列，x(n)相当于x(t)，其中n＝1,2，……N；

32、标准差F_std:

33、均方根值F_rms:

34、峰值指标F_p:

F_p＝max|x(n)|

35、波形指标F_sf:

36、峰值因子P_cf:

P_cf＝F_p/F_rms

37、脉冲指标F_if:

38、裕度指标F_clf:

39、峭度F_k:

310、偏斜度F_s:

311、频率均值F_mf:

式中s(k)是信号x(n)的频谱，k＝1,2,……K；K是谱线数；

312、中心频率F_cf:

式中f_k是第k条谱线的频率值；

313、频率均方根F_rmfs:

314、频率标准偏差F_stdf:

对原始振动信号提取上述时域与频域指标，得到14个特征，设为特征集1；

第二步原始振动信号集合经验模式分解后得到a个本征模式分量，对每个本征模式分量分别提取上述时域与频域指标，得到14×a个特征，设为特征集2；

第二步原始振动信号小波包分解后得到B个子频带，分解L层中第l层的第j个子频带能量用公式(1)表示：

是第j个子频带的成分，提取每个子频带能量作为特征得到B个特征，设为特征集3；

以上三个特征集综合为原始特征集X，其中特征总数为14+14×a+B；

第四步：优选特征

41、敏感特征选择

包含了C种轴承状态类别的特征矩阵为

{f_g，i，g＝1，2，...，Num；i＝1，2，...，I}_C，

本发明样本总数T＝Num*C，为上述原始特征矩阵中的行数；I是特征种数，是上述原始特征矩阵中的列数，由第三步可知I＝14+14×a+B；从各样本得到的特征总个数为Num×C×I＝T×I，为上述原始特征矩阵的元素总数；f_g，i是第g个样本的第i个特征值；

所有的特征集合定义为{v_t，i}t＝1,2,…,T；i＝1,2,…,I；自权重算法具体过程如下：

411、v_t，i进行正规化，得到正规化后的特征集{x_t，i}

412、计算自相似因子(self-similarity factor)SF_mt，i

SF_mt，i＝||x_m，i-x_t，i||²，m，t＝1，2...，T；i＝1，2，...，I (3)

413、对每个特征建立权重矩阵

414、计算每个特征的自权重Sw_i

根据步骤413得到的权重，计算每个特征的自权重Sw_i

Sw_i＝mean(W_i) (5)

将所有特征的自权重值Sw_i从大到小排列，选出前R个自权重相对较大的特征作为去冗余的基础，此R个特征为敏感特征，得到敏感特征矩阵{x_t，r}；

42、冗余特征去除

利用AP算法去冗余步骤如下：

421、敏感特征矩阵的转置

步骤41选出的敏感特征矩阵{x_t，r}＝X_t，r

化为转置阵

422、AP算法的偏向参数

用欧式距离建立X_r，t的相似矩阵Sim₄，由公式7设置AP算法中的偏向参数P₄；

Sim₄＝-||X_r，t-X_r+i，t||² (6)

P₄＝2·median(Sim₄) (7)

423、优选特征矩阵X_t，z

将X_r，t导入AP聚类算法，自动为每个特征打上类标签，聚类得到Z个类别，选取每一类的聚类中心作为最终的优选特征，得到优选特征矩阵X_t，z；

第五步：故障识别

AP算法将每个样本间的相似度之和最大化，见公式(8)：

式(8)中的h是数据点p的潜在聚类中心，Sim(p，h)是数据点p和h之间的相似度值；AP算法中传递吸引度和归属度两种类型的消息；吸引度re(p,h)表示从点p发送到候选聚类中心h的数值消息，反映h点是否适合作为p点的聚类中心；归属度av(p,h)表示从候选聚类中心h发送到p的数值消息，反映p点是否选择作为其聚类中心；re(p,h)与a av(p,h)越大，则h点作为聚类中心的可能性就越大，并且p点隶属于以h点为聚类中心的聚类的可能性也越大；

将第四步获得的优选特征矩阵X_t，z导入AP聚类算法中识别故障，主要步骤如下：

51、建立相似矩阵Sim₅；相似度Sim₅度量采用欧氏距离，偏向参数P₅设置为Sim₅均值的2倍，如下公式所示；

建立T个数据点的相似矩阵Sim₅，设置偏向参数P₅的值；

Sim₅＝-||X_t，z-X_t+1，z||²，t＝1，2，...，T；z＝1，2，...，Z (9)

P₅＝2·medi(Sim₅) (10)

52、初始化归属度矩阵av＝0；

53、更新吸引度的值re(p,h)

(11)式中h′表示另外一个数据点；h′∈{1，2，...，T}，但h′≠h；按(11)式循环，直至得到{av(p，h′)+Sim5(p，h′)}的最大值；

在每次re更新时加入一个阻尼系数λ，0≤λ≤1，re按下式更新

re_p＝(1-λ)re_p+λre_p-1 (12)

54、更新归属度的值av(p,h)

(13)式中p′表示另外一个数据点；p′∈{1，2，...，T}，但p′≠h同时p′≠p

re(h，h)＝P(h)-max{av(p，h)+sim₅(p，h)}

在每次av更新时加入一个相同的阻尼系数λ，av按下式更新

av_p＝(1-λ)av_p+λav_p-1 (14)

55、确定潜在聚类中心

当满足av(h，h)+re(h，h)＞0时，h点被选择成为一个潜在聚类中心；

56、如果在数次迭代后，聚类中心点不变，那么进入步骤57；否则，返回步骤53；

57、使用最邻近规则为剩余的数据点安排归属类，同一类中每个数据点的相似度均为最大值；

吸引度和归属度两类消息传递过程的终止是达到了用户设置的固定迭代次数；

经上述五个步骤完成诊断模型训练，确定各步骤的集合经验模式分解中的添加高斯白噪声标准差值，添加噪声的次数，本征模态函数个数a；小波包分解所采用的小波基函数，子频带个数B，分解层数L；自权重法选择的敏感特征个数R，优选特征个数Z以及AP算法中的相似度，偏向参数，阻尼系数和迭代次数，训练好的故障样本和正常样本均保留在模型中，作为用于对某机械系统的滚动轴承故障进行识别的训练好的诊断模型存储于上位机；

第二部分、故障诊断

第1步：振动信号采集

采用第一步的振动传感器安装于机械系统中的每个待检测的滚动轴承的轴承座上，采集机械装备的振动信号x(t)；振动传感器的信号线连接上位机，各轴承实时振动信息储存于上位机中为未知振动信号；

第2步：信号分解

上位机取完成训练的诊断模型中确定的参数对采集到的未知信号x(t)进行集合经验模式分解得到各个频带的本征模态函数，并进行小波包分解得到分解后的B个子频带；

第3步：特征提取

在未知振动信号、第2步分解所得的本征模态函数以及B个子频带中提取按完成训练的诊断模型确定的Z个优选特征；

第4步：故障诊断

第3步提取的Z个优选特征，导入AP算法聚类，采用完成训练的诊断模型中AP算法中的相似度，偏向参数，阻尼系数和迭代次数，未知信号的优选特征与诊断模型中已知的Z个特征比较，根据未知信号的优选特征归属于完成训练的诊断模型中的哪个聚类中心，判断得到当前未知信号的类别，即完成故障识别诊断。

2.根据权利要求1所述的基于聚类分析的机械系统滚动轴承智能故障诊断方法，其特征在于：

所述第一步振动信号采集中的振动传感器为压电加速度传感器、位移传感器和电涡流位移传感器中的任一种。

3.根据权利要求1所述的基于聚类分析的机械系统滚动轴承智能故障诊断方法，其特征在于：

所述第二步的信号分解中添加高斯白噪声标准差值设置为信号标准差的0.1～0.2倍；添加噪声值设置为信号点数的0.08～0.15。

4.根据权利要求1所述的基于聚类分析的机械系统滚动轴承智能故障诊断方法，其特征在于：

所述第二步的信号分解中的小波包分解采用5～10阶消失矩的DB系列小波基函数。

5.根据权利要求1所述的基于聚类分析的机械系统滚动轴承智能故障诊断方法，其特征在于：

所述第四步优选特征中的R值由用户视具体振动信号的情况根据经验决定，R取值为特征总数的20％～30％。

6.根据权利要求1所述的基于聚类分析的机械系统滚动轴承智能故障诊断方法，其特征在于：

所述第五步故障识别时设置迭代次数为800～1200。