CN106909756A - 一种滚动轴承剩余寿命预测方法 - Google Patents

Abstract

本发明公开了滚动轴承剩余寿命预测方法，涉及滚动轴承剩余寿命预测领域。本发明提供一种滚动轴承剩余寿命预测方法，包括：一种滚动轴承剩余寿命预测方法，其特征在于，包括：步骤一、采集并得到表征滚动轴承运行状态的原始振动信号；步骤二、对所述步骤一得到的所述原始振动信号进行消噪处理，得到能够具体表征滚动轴承运行状态的有用信号；步骤三、提取所述步骤二得到的所述有用信号的时域特征和频域特征；步骤四、对所述时域特征与所述频域特征进行特征融合，得到表征滚动轴承运行趋势的特征指标；步骤五、构建剩余寿命预测模型，对所述特征指标进行剩余寿命预测。本发明提供的滚动轴承剩余寿命预测方法预测精度高，预测效果好。

Description

一种滚动轴承剩余寿命预测方法

技术领域

本发明涉及滚动轴承剩余寿命预测领域，具体而言，涉及滚动轴承剩余寿命预测方法。

背景技术

机械设备应用于人类生活、工作和生产的方方面面，并在其中扮演了举足轻重的角色。目前，机械设备正朝着大型化、高速化、精密化、系统化、连续化和自动化方向发展，机械设备的运行环境越来越复杂多变，更是为设备的健康管理提出了新的挑战。随着设备的运行，零件老化、可靠性降低、剩余寿命减少等问题，与设备能否持续安全高效地工作、维护能否及时有效地执行，亟待人类解决。轴承是机械设备重要的组成元件，在机械系统中发挥着十分重要的作用。轴承的性能退化趋势及寿命预测方法一直是机械设备健康管理研究的重点。为了有效防止设备工作精度下降，最大限度地利用以轴承为代表的机械设备关键零部件的工作能力，节约物料开支，减少事故发生，越来越有必要跟踪零部件的运行过程。对机械设备关键零部件进行寿命研究已成为现代设备健康管理中极其重要的一环。

发明内容

本发明的目的在于提供一种滚动轴承剩余寿命预测方法，其预测 精度高，预测效果好。

本发明提供一种技术方案：

一种滚动轴承剩余寿命预测方法，包括：步骤一、采集并得到表征滚动轴承运行状态的原始振动信号；步骤二、对所述步骤一得到的所述原始振动信号进行消噪处理，得到能够具体表征滚动轴承运行状态的有用信号；步骤三、提取所述步骤二得到的所述有用信号的时域特征和频域特征；步骤四、对所述时域特征与所述频域特征进行特征融合，得到表征滚动轴承运行趋势的特征指标；步骤五、构建剩余寿命预测模型，对所述特征指标进行剩余寿命预测。

进一步地，所述步骤五包括拟合特征指标步骤和计算误差权值步骤，所述拟合特征指标步骤：根据所述滚动轴承的历史工作状态信息确定寿命研究的共有时间段。并对该时间段内的各零部件的所述特征指标进行曲线拟合，为保证拟合的有效性，利用基于粒子群算法的最小二乘支持向量机模型进行拟合，得到进行寿命研究多个所述剩余寿命预测模型；

所述计算误差权值步骤：确定目标零部件某个已知的时间点及其对应特征指标值，将该指标值分别输入到多个所述剩余寿命预测模型中，每个模型都会有一个预测结果输出。将每个模型的预测时间值和目标零部件的实际的时间值相对照，通过e_ik＝(t_i-t_k)²计算各个预测模型的误差。计算每个预测模型的权重如下式：

其中，n为预测模型个数，ω_i为预测模型各自的权重，t_i是由各预测模型预测的当前运行时间点。t_pre为计算所得的零部件性能退化的总体时间，t_k为目标零部件测试点处的已知时间。

进一步地，所述特征融合的方法为主成分分析法，所述时域特征与所述频域特征组成多维特征矩阵，所述多维特征矩阵通过所述主成分分析法得到所述特征指标。

进一步地，所述步骤二利用小波消噪方法对所述原始振动信号进行消噪，并得到所述有用信号。

本发明的有益效果在于：依靠在线监测振动传感器采集的滚动轴承的振动信号，获取有用信息并在时域、频域、时频域等域提取信息的状态特征，得到能表征轴承运行状态的多维特征集，并进行滚动轴承剩余寿命预测，是轴承基于视情维修、提高可靠性的重要基础。对于寿命预测的研究，首先确定寿命研究的时间段，通过基于粒子群算法优化的最小二乘支持向量回归机对同型号的多个零部件退化性能特征指标进行拟合，得到多个寿命预测的预测模型。利用待测零部件的在线状态特征输入到各预测模型进行训练，得到误差权值。输入待测零部件当前在线状态监测信息，预测当前寿命和计算剩余寿命。其 拥有较高的预测精度和预测效果。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍。应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对范围的限定。对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。

图1为本发明的实施例提供的滚动轴承剩余寿命预测方法的流程示意图；

图2为本发明的实施例提供的误差累计和方法的流程示意图。

具体实施方式

下面结合附图，对本发明的具体实施方式进行详细说明。

请参阅图1，本实施例提供了一种滚动轴承剩余寿命预测方法。下面对该方法的具体步骤进行说明。

步骤(1)：在线监测振动传感器采集滚动轴承振动信号(包括了反映轴承运行性能状态的有用信号及噪声信号)。

步骤(2)：采用MATLAB软件自带的小波消噪函数(主要依据Matlab算法，其中基函数、分解层数、阈值方法及阈值等选取是关键)消除振动信号中的大部分噪声干扰，并保留能够反映轴承运行性能状态的有用信号。

步骤(3)：提取时域特征集(包括量纲特征及无量纲特征，详细见表1和表2)、频域特征集(详细见表3)、时频域特征集(经验模态分 解方法(Empirical ModelDecomposition,EMD)能有效地提取轴承状态监测信号中的非线性、非平稳信号特征，且具备较高时频分辨率，此处的时频域特征集即为EMD提取的轴承状态监测信号中的非线性特征)及威布尔参数特征集(威布尔分布的形状参数和尺度参数可作为反映轴承运行状态的特征量。由于希尔伯特变换能较好地提取轴承早期故障信息，因此状态监测信号先经过希尔伯特变换得到新的数据，而后采用极大似然方法即可得到威布尔分布的形状参数和尺度参数，即威布尔参数特征集)，得到能够表征轴承运行性能状态的多维特征集，并构建多维特征集矩阵。

表1含量纲时域特征指标

表2无量纲时域特征指标

表3频域特征

机械设备关键零部件的状态监测信号中通常包含非线性、非平稳信息，经验模态分解方法(Empirical Mode Decomposition,EMD)是一种非线性、非平稳信号分析方法，它将信号分解成一系列本征模函数(Intrinsic Mode Function,IMF)之和。对于非线性的特征，应用较为广泛的方法还有小波方法。相较于小波方法，经验模态分解方法不但能够有效的提取信号中的微弱特征，而且不用选取基函数，因此，自适应性较强，非常适合处理非线性信号。通过经验模态分解方法消除了人为因素对结果的影响，分解结果依赖分析信号本身，能有效地提取非线性、非平稳信号的微弱特征，并可得到较高的时频分辨率，具有良好的时频聚集性。所以本实施例利用经验模态分析方法来提取零部件监测信号中的非线性特征。

设通过经验模态分解，信号获得多个模态分量f_i(t)和余项r_n(t)，把r_n(t)看作信号的第n+1个分量f_n+1(t)，则第i(i＝1,2,...,n+1)个分量f_i(t)的能量为：

式中，N为IMF分量f_i(t)的数据长度。

根据经验模态分解完备正交性，可以得出：

E[x(t)]＝E[f₁(t)]+E[f₂(t)]+...+E[f_n+1(t)]

威布尔分布是在可靠性分析中常用的统计分布模型之一，可分为两参数威布尔分布和三参数威布尔分布。威布尔分布的形状参数和尺度参数也可作为反映零部件运行状态的特征量。本实施例考虑两参数威布尔分布，其概率密度函数为：

式中β是形状参数，η是尺度参数。

通常将参数估计方法分为图解法和解析法两大类，图解法包括经验分布图法、威布尔概率图法和风险率统计图法等；解析法包括极大似然估计法和回归估计法等。本实施例通过极大似然估计法求解参数。极大似然估计函数为：

式中，N为数据长度。由于希尔伯特变换对零部件早期故障信息的提取具有良好效果，x_i为状态检测信号经过希尔伯特变换后的数据。

关于上式，分别对β和η求偏导，并令其为零，得：

上式为非线性方程组，利用牛顿迭代法求解得到参数β和η的估计值。由于通过极大似然估计法得到的估计值通常会存在较大误差，修正形状参数β可以使估计误差减小，其修正公式为：

式中，N为数据长度，β_U是修正值，β是估计值。

步骤(4)：利用主成分分析法PCA(Principal Components Analysis)对滚动轴承消噪后振动信号的多维特征进行融合，在不减少原有数据所包含信息内容的基础上，得到能全面表征滚动轴承运行状态性能退化趋势的特征指标。

利用主成分分析进行多维特征融合得到零部件运行过程中的性 能退化特征指标是将原来的多维特征集通过某种线性映射，形成比原有维数更少的新的多维特征集。新的特征集间两两相互正交，且向量按照特征间的差异性从大到小排序，排序第一的向量反映了特征间的最大差异。这样由多维特征间最大差异量组成的新向量不仅可以涵盖特征中的有用信息，而且去除了原来多维特征集的冗余部分，实现了用主成分的方式来描述高维数据特征。

设某次采集到的数据中共有n个数据，每个数据有p个特征变量，则构成一个n×p阶的多维特征矩阵：

式中x_i＝(x_1i,x_2i,...,x_ni)^T，x_ij是第i个数据的第j个特征变量。

通过主成分变换得到由x₁,x₂,...,x_p表示的线性组合：

如果系数l_ij满足：1)l_i1 ²+l_i2 ²+...+l_ip ²＝1，i＝1,2,..,p；2)系数l_ij使线性组合y_i与y_j(i≠j)相互无关；3)y_i按照x₁,x₂,...,x_p的所有线性组合方差大小排序，则原多维特征集的第一、第二、第p主成分为y₁,y₂,...,y_p。

上式所表示的主成分分析过程就是通过对多维特征的协方差矩 阵求特征值与特征向量，并将按大小顺序排列后的特征值所对应的特征向量作为线性组合的系数。利用主成分分析方法对多维特征X进行融合得到零部件性能退化特征指标的步骤是：

(1)计算零部件多维特征矩阵X的均值向量:

将多维特征中心化：

(2)计算多维特矩阵的协方差矩阵：

计算多维特征矩阵的特征值与特征向量，将最大特征值所对应的特征向量作为对多维特征矩阵线性变化的系数，得到零部件性能退化特征指标。

为进行机械设备关键零部件的寿命研究需要解决如何获得与寿命密切相关的特征指标和如何预测剩余寿命。本文在第三章中主要针对第一个问题进行研究，提出了的基于小波消噪、多域特征提取和主成分分析方法的解决方案。该方案可以将零部件运行状态的监测信号转化为表征这一运行过程的特征指标，该特征指标极大的保留了状态监测信号中的有用信息。针对第二个关键问题，首先考虑到特征指标值存在一定波动，为避免可能因此产生的误差，本文利用曲线拟合特 征指标数据。至于预测模型的选择，在本文的讨论中，基于粒子群算法优化最小二乘支持向量回归机参数的预测模型表现了良好的预测能力，在对机械设备关键零部件进行寿命预测时本文仍然选用该模型。但在进行寿命预测时，由于单一的预测模型可能会存在较大误差，多个模型相互配合能减小误差发生的几率。而不同模型在一次预测中的预测结果往往各异，对待测目标值的影响也不一样。预测误差权值计算方法通过训练计算每个预测模型对于研究目标预测结果所占的权值，从而预测当前特征指标值所对应的工作寿命，进而得到剩余寿命。预测误差权值计算方法是根据每个训练误差值的大小来确定该模型在预测结果中所占的权值。这种由多个预测模型误差权值计算得到的被测目标寿命可以减少由于模型选取不当而引起误差产生的概率。误差累积和方法在某些情况下会对训练模型进行重新训练，而发生这一情况的根本原因是由于特征指标的实际值发生了较大变化，也就是说零部件的运行状态发生了变化，即发生了故障。根据这一事实，可以利用误差累积和方法估计零部件的潜在故障点和功能故障点，即确定进行寿命预测的时间段(P-F间隔段)。

综上所述，本文进行机械设备关键零部件寿命预测的研究步骤如下：

(1)状态监测及信号处理。对同型号的多个零部件在相同条件下进行加速度试验，分别采集相应的振动信号。对机械设备关键零部件的运行状态进行在线监测，对获得的每个零部件的振动信号进行滤波处理(小波消噪方法)，以保证振动数据包含的特征信息不被噪声淹没；

(2)从时域、频域和时频域提取每个零部件的监测数据的特征信息，利用特征融合的方式(主成分分析方法)获得机械设备关键零部件运行过程的特征指标。

(3)拟合特征指标。根据零部件历史工作状态信息确定寿命研究的共有时间段。并对该时间段内的各零部件的特征指标进行曲线拟合，为保证拟合的有效性，利用基于粒子群算法的最小二乘支持向量机模型进行拟合，得到进行寿命研究多个预测模型；

(4)计算误差权值。确定目标零部件某个已知的时间点及其对应特征指标值，将该指标值分别输入到由第三步获得的多个预测模型中，每个模型都会有一个预测结果输出。将每个模型的预测时间值和目标零部件的实际的时间值相对照，通过e_ik＝(t_i-t_k)²计算各个预测模型的误差。计算每个预测模型的权重如下式：

其中n为预测模型个数，ω_i为预测模型各自的权重，t_i是由各预测模型预测的当前运行时间点。t_pre为计算所得的零部件性能退化的总体时间，即P-F间隔时间，t_k为目标零部件测试点处的已知时间。

则目标零部件在此时刻预测的剩余寿命可以表示为：

t_remain＝t_pre-t_k (0-3)

寿命预测实例研究

步骤(5)：根据步骤(4)得到的结论，设计预测模型及参数优化。以最小二乘支持向量机作为预测模型，对参数模型优化中使用交叉验证法，并综合使用粒子群优化算法，来提高每次训练的准确率和精确度。并利用误差累计和方法(Cumulative Sum,CUSUM)对当前模型的预测能力进行检验，进一步提高预测的准确性。

通过粒子群算法优化得到的预测模型，虽然参数已经达到某种意义上的最优，但是却缺乏对模型的预测能力进行进一步的检验，本实施例利用用误差累积和方法(Cumulative Sum,CUSUM)对当前模型的预测能力进行检验，进一步提高预测的准确性。

请参阅图2，以下给出误差累计和方法的具体步骤。

假设预测模型的误差Z＝y_i-y_i是服从均值为μ，方差为σ²的正态分布，对服从正态分布的数组Z进行标准化后得到：

累积和方法利用两个指标来检验当前模型的预测能力：

UB_i＝max[0,(d_i-m)+UB_i-1]

LB_i＝max[0,(-d_i-m)+LB_i-1]

在上式中，UB_i表示误差的正偏差，LB_i表示误差的负偏差。初始 状态UB₀和LB₀都为0，通常情况下，参数m取0.5，阈值h取3。

本实施例提出的基于最小二乘支持向量回归机、粒子群算法和误差累积和方法进行机械设备关键零部件趋势预测的研究步骤如下：

(1)对机械设备关键零部件的运行状态进行在线监测，对获得的振动信号进行滤波处理(小波消噪方法)，以保证振动数据包含的特征信息不被噪声淹没；

(2)提取监测数据的多维特征信息，利用特征融合的方式(主成分分析方法)获得机械设备关键零部件运行状态性能退化趋势的特征指标；

(3)建立基于最小二乘支持向量回归机、粒子群算法和累积和方法的趋势预测模型；

(4)将特征指标输入预测模型，实现对机械设备关键零部件运行状态的退化趋势预测，并输出预测结果。

步骤(6)：趋势预测，将步骤4中得到的特征指标作为趋势预测模型的输入数据，并确定训练集和测试集；对滚动轴承运行状态性能退化趋势进行预测。

在本实施例中，训练集长度为100，本发明所提方法和对比方法进行单步预测，其中交叉验证法均选择5-Fold；粒子群算法的主要参数设置分别为：最大进化数量为200、种群最大数量为20、全局和局部搜索能力分别设置为1.7和1.5；遗传算法的主要参数设置为：最大进化代数为200，种群最大数量为20。

Claims (4)

1.一种滚动轴承剩余寿命预测方法，其特征在于，包括：

步骤一、采集并得到表征滚动轴承运行状态的原始振动信号；

步骤二、对所述步骤一得到的所述原始振动信号进行消噪处理，得到能够具体表征滚动轴承运行状态的有用信号；

步骤三、提取所述步骤二得到的所述有用信号的时域特征和频域特征；

步骤四、对所述时域特征与所述频域特征进行特征融合，得到表征滚动轴承运行趋势的特征指标；

步骤五、构建剩余寿命预测模型，对所述特征指标进行剩余寿命预测。

2.根据权利要求1所述的滚动轴承剩余寿命预测方法，其特征在于，所述步骤五包括拟合特征指标步骤和计算误差权值步骤，所述拟合特征指标步骤：根据所述滚动轴承的历史工作状态信息确定寿命研究的共有时间段。并对该时间段内的各零部件的所述特征指标进行曲线拟合，为保证拟合的有效性，利用基于粒子群算法的最小二乘支持向量机模型进行拟合，得到进行寿命研究多个所述剩余寿命预测模型；

所述计算误差权值步骤：确定目标零部件某个已知的时间点及其对应特征指标值，将该指标值分别输入到多个所述剩余寿命预测模型中，每个模型都会有一个预测结果输出。将每个模型的预测时间值和目标零部件的实际的时间值相对照，通过e_ik＝(t_i-t_k)²计算各个预测模型的误差。计算每个预测模型的权重如下式：

${\mathrm{\ω}}_i=\frac{(\underset{i=1}{\overset{n}{\Σ}}{e}_{ik}/e_{ik})}{\underset{i=1}{\overset{n}{\Σ}}(\underset{i=1}{\overset{n}{\Σ}}{e}_{ik}/e_{ik})}$

$t_{pre}=\underset{i=1}{\overset{n}{\Σ}}{\mathrm{\ω}}_i{t}_i$

其中，n为预测模型个数，ω_i为预测模型各自的权重，t_i是由各预测模型预测的当前运行时间点。t_pre为计算所得的零部件性能退化的总体时间，t_k为目标零部件测试点处的已知时间。

3.根据权利要求1所述的滚动轴承剩余寿命预测方法，其特征在于，所述特征融合的方法为主成分分析法，所述时域特征与所述频域特征组成多维特征矩阵，所述多维特征矩阵通过所述主成分分析法得到所述特征指标。

4.根据权利要求1所述的滚动轴承剩余寿命预测方法，其特征在于，所述步骤二利用小波消噪方法对所述原始振动信号进行消噪，并得到所述有用信号。