基于烟花粒子滤波的电液舵机轴承剩余寿命预测

郭润夏，张国良，王 雨

（中国民航大学电子信息与自动化学院，天津 300300）

轴承的剩余寿命（RUL，remaining useful life）主要是通过提取振动信号的时域和频域特征来进行分析和判断的[1]，主要有两种方法：物理模型方法和数据驱动方法[2]。物理模型方法是指根据物理原则建立轴承状态数学模型[3]，再根据以往的观测值预测轴承未来的状态，如用自回归的方法来预测回转轴承的RUL[4]。数据驱动方法是通过大量的趋势数据学习轴承的状态，并预测其未来状态，如引入深度神经网络跟踪轴承状态并预测RUL[5-6]。然而，物理模型方法和数据驱动方法都有其自身的局限性，一般而言，物理模型难以建立，而数据驱动方法需要大量的数据，浪费计算资源。

粒子滤波算法适用于基于物理模型的非线性、非高斯系统的状态估计，其通过一组加权粒子来近似后验概率密度，主要包括3 个步骤[7]：重要性采样、权重计算、重采样。到目前为止，粒子滤波算法已成功应用于许多领域的RUL 预测[8-9]。但粒子滤波算法存在一个严重的问题，即粒子退化问题[10]，为了将粒子滤波算法有效地应用于实际工程中，首先要解决粒子退化问题。

为解决上述问题，提出了一种改进的粒子滤波算法——基于改进烟花算法的粒子滤波（FPF，fireworks algorithm-based particle filter）算法，并通过实验验证了该算法的有效性。FPF 采用了新的重采样过程，即对烟花算法进行改进，然后利用改进的烟花算法优化粒子权值，最后通过高斯近似方法产生新粒子。改进后的粒子滤波算法能有效抑制粒子退化。在此基础上，根据Paris-Erdogan 模型预测了轴承的RUL。

1 烟花算法

1.1 标准烟花算法

烟花算法[11]是一种新型的群智能寻优算法，每个初始烟花将被视为可行解空间的一个解, 通过爆炸操作使每个烟花爆炸产生爆炸火花，然后选择一部分烟花随机变异，产生变异火花，最后,从烟花、爆炸火花和变异火花的集合中选择下一代烟花。

烟花算法主要由4 个基本部分组成：爆炸操作、变异操作、映射操作和选择策略。不失一般性，假设待解决的优化问题目标函数为

式中：X 是函数的独立变量；Ω 是解空间的可行域。

1）爆炸操作

根据适应度函数计算每个烟花Xi的适应度值f（Xi）。适应度越好的烟花产生的火花越多，适应度越差的烟花产生的火花越少，从而达到局部搜索和全局搜索的平衡。烟花Xi的爆炸半径Ai和火花数Si的计算式分别为

式中：A 和M 分别是用来调整爆炸幅度、爆炸数目的两个常数；уmin和уmax是烟花的最小适应度值和最大适应度值；ε 是一个机器数。

同时，为了避免产生过多或过少的火花，采用以下规则来限制火花的数量，即

式中：a 和b 是两个常数，a ＜b ＜1；round()是取整函数。

根据式（2）～式（4）计算每个烟花Xi在第k 维度爆炸产生的火花为

2）变异操作

随机挑选部分烟花变异以增强种群多样性，烟花Xi在第k 维变异产生的火花计算如下

式中gauss()为高斯函数。

3）映射操作

当烟花Xi在第k 维进行爆炸、变异操作时，如果超出了烟花种群在该维度的界限，则进行如下操作

4）选择策略

从原始烟花、火花中选择出下一代烟花，在候选集合K 中，适应度值最小的原始烟花或火花必被选择成为下一代烟花，其余N-1 个烟花以轮盘赌的方式选择，每个烟花或火花被选择的概率为

1.2 改进的烟花算法

重采样过程是粒子优化过程，而可以用烟花算法来解决粒子优化问题。为了使烟花算法更好地设计重采样过程，将有针对性地对烟花算法进行改进，具体改进操作如下。

1）选择性的爆炸操作

除了适应度值最小的烟花之外，根据式（1）～式（5）对剩余的烟花执行爆炸操作。

2）选择初始下一代烟花

执行完上一步操作后，根据适应度值从初始烟花、火花中选出N 个初始下一代烟花，被选出的烟花需要满足以下条件

3）选择需变异的烟花

用烟花集合C 表示选出来的初始下一代烟花，将集合中的元素按适应度值大小进行排序，并按照黄金分割比将N 个烟花划分为两个部分，选择小比例部分的烟花，分割示意图如图1 所示。

图1 分割示意图Fig.1 Splitting diagram

假设烟花集合中有8 个烟花，先将8 个烟花按适应度值从小到大排序，根据黄金分割比将其划分为红、绿两部分，3 个红色的烟花即是需要变异的烟花。

4）变异操作

对步骤3）中选择出的烟花执行变异操作，变异后的烟花与未变异的烟花最终形成下一代烟花，变异操作如下

当爆炸火花、变异火花超出界限时，根据式（7）执行映射操作。

2 粒子滤波

本节主要介绍基础粒子滤波（Basic PF）[2，12-13]及引入改进的烟花算法来设计粒子滤波算法的重采样过程。

2.1 基础粒子滤波算法

一般而言，系统的空间模型如下

式中：xk，zk分别表示k 时刻的系统状态和测量数据；vk-1，wk分别表示已知分布的k-1 时刻的过程噪声和k 时刻的测量噪声, 并假设过程噪声与测量噪声独立不相关。

已知初始状态分布为p（x0），根据式（11）和式（12）分别计算系统的状态转移概率密度函数、似然概率密度函数为

式中：p（xk|xk-1）表示在状态xk-1下xk的状态转移概率密度函数；p（zk| xk）是在状态xk下zk的似然概率密度函数；pv()为过程噪声的概率密度函数；pw()为测量噪声的概率密度函数。

粒子滤波算法[12-13]基于蒙特卡罗模拟思想，从重要性概率密度函数q（xk| xk-1，zk）中采集N 个样本粒子，并计算这些粒子的权值集，然后用该组粒子及其对应的权值来近似系统状态的后验概率密度，每个粒子的对应权值计算式如下

根据式（15）和式（16），可以得到权值的迭代公式如下

计算得到粒子的权值后，归一化权值

再计算系统的后验概率密度

最后，可以得到系统的状态估计为

由式（17）可看出，随着迭代次数的增加，大部分粒子的权值会越来越小，这就是粒子退化现象，因此，有必要重新设计重采样过程解决粒子的退化问题。

2.2 基于改进烟花算法的粒子滤波算法

改进烟花算法设计重采样过程并给出了FPF 的具体操作步骤。

2.2.1 重采样

重采样是为了抑制粒子的退化，通常用有效粒子数Neff来衡量粒子的退化程度，有效粒子数越少，退化越严重，有效粒子数计算如下

重采样的基本过程是从离散近似的后验概率密度中重采样产生一组新的粒子。重采样的原则是在权重大的粒子附近产生新粒子，使新粒子也具有较大的权值。当粒子出现退化，即有效粒子数小于设置的重采样阀值Nthr时，经过改进的烟花算法进行操作，在权重大的粒子附近产生较大权重的新粒子，从而增加了有效粒子数，满足重采样的原则要求，解决了基础粒子滤波算法的粒子退化问题，提高了算法的估计精度。

粒子的权值称为烟花权值，当新的观测数据到达时，得到带有烟花权值的粒子集合当出现粒子退化现象时，根据改进的烟花算法对烟花权值进行优化，获得下一代烟花权值，并通过高斯近似产生新的粒子。

在新的重采样过程中，适应度函数为

重采样流程图如下。

图2 重采样步骤Fig.2 Resampling steps

2.2.2 操作步骤

图3 FPF 的执行过程Fig.3 The execution process of FPF

步骤1初始化：

（1）从初始概率密度p（x0）采样；

步骤2序惯重要性采样：

步骤3当Neff＜Nthr，重采样及状态估计：

2.3 FPF 仿真验证

在本节中，采用经典的非线性模型

验证FPF 的有效性和可靠性[13]。式中wk服从gauss（0，1）分布。

初始粒子服从gauss（0.5，1）分布。粒子数N=200，重采样阀值Nthr=N/2。同时，以均方根误差（RMSE）来评价算法的性能，即

假设系统模型是精确的，根据状态方程及初始真值x0=0.5 计算下一时刻的状态值。为了展示FPF 能有效抑制粒子退化，提高状态估计精度，将使用Basic PF 与之进行比较。图4 为状态估计仿真图，图5 和图6 分别为有效粒子数和均方根误差仿真图。

图4 状态估计Fig.4 State estimation

图5 有效粒子数Fig.5 Number of effective particles

由图4 可看出，FPF 的状态估计更接近真实值。图5 中，Basic PF 的有效粒子数几乎变为1，但FPF 的有效粒子数仍在100～200 之间。图6 中，与Basic PF相比，FPF 的均方根误差较小。综上，FPF 可以在抑制粒子退化的同时增加有效粒子的数量，状态估计更准确。

图6 均方根误差Fig.6 Root mean square error

3 电液舵机轴承的RUL 预测

RUL 定义为从预测起点开始到预测损伤指标达到某一阀值的时间。当损伤指标超过阀值时，认为该部件不能再继续长期使用。一般来说，在不同的领域，RUL的计算方法不同，因此，RUL 的计算定义如下

式中：d()是距离函数；xi是预测起点的状态；xj是预测终点的状态；Δt 是采样间隔；xthreshold是预定义的阀值；ξ是一个可调的较小的常数。利用FPF 算法对轴承状态进行估计后，可根据式（25）得到预测的RUL。

为了有效地预测轴承的RUL，提取振动加速度信号的时域特征均方根（RMS）作为损伤指标[2，14]（该指标也是轴承在k 时刻的健康状态xk）；定义为

根据提取出的轴承特征状态，建立轴承的退化模型[14]如下

式中β 和m 是模型参数。

当加速度信号均方根值达到1.5 g 时，轴承振动信号波动剧烈，说明轴承磨损严重，因此，设置xthreshold=1.5 g。

实验平台由加速疲劳装置、电控可编程逻辑控制器（PLC，programmable logic controller）装置、电机、齿轮箱和传感器组成，如图7 所示。传感器主要包括加速度传感器、速度传感器和扭矩传感器。信号（如振动）由传感器测量，并通过采集卡传递到计算机。

图7 实验平台Fig.7 Experimental platform

实验平台给出不同时间段的FPF 的预测结果，同时给出Basic PF、相对熵粒子滤波（MREIS-PF）的预测结果对比。图8 是FPF 对3 个时间点损伤指标的预测；图9 为FPF、Basic PF 和MREIS-PF 在同一时间起点的状态预测；图10 为FPF、Basic PF 和MREIS-PF在整个过程中多个时间起点预测的RUL，并与真实RUL 进行对比。

图8 FPF 在3 个时间起点的状态预测Fig.8 State prediction by FPF at three time starting points

图9 不同算法在同一时间起点的状态预测Fig.9 State prediction by different algorithms at the same time starting point

图10 剩余寿命预测结果Fig.10 Prediction of RUL

表1 给出了基于不同算法的多组轴承的RUL 预测和误差率指标，误差率指标定义如下

表1 不同算法的多组轴承的预测结果Tab.1 Predicted results of multiple sets of bearings obtained from different algorithms

由图8～图10 和表1 可看出：由于粒子的退化，Basic PF 在状态预测期间波动较大且不够准确，导致预测的轴承RUL 与真实值相差较大；由于引入了设计的重采样过程，FPF 在状态跟踪过程中更加稳定，比Basic PF 和MREIS-PF 有更好的预测结果；在预测的早期阶段，由于轴承退化过程缓慢，FPF 预测的RUL时间要比实际RUL 时间长，在预测的中间阶段，轴承的磨损速率比早期阶段快，FPF 预测的结果更接近真实值，当轴承进入预测的后期阶段，轴承退化加剧,由于模型包含了缓慢退化过程的信息，不能完全适应轴承的快速退化，导致预测值大于真实值。综合分析各阶段预测结果，FPF 能较好地跟踪轴承的磨损状态，预测轴承的RUL 与真实值比较接近，在预测轴承RUL 方面显示出较强的优势。

4 结语

针对粒子滤波算法中的粒子退化问题，采用改进的烟花算法设计重采样过程，提出了一种改进的粒子滤波算法，可以减小粒子权值方差，抑制粒子退化。根据Paris-Erdogan 模型建立的物理模型能较好地适应轴承的退化过程。同时，基于振动加速度信号提取的轴承状态的时域特征均方根值能较好地反映轴承的磨损程度。在此基础上，利用所提出的粒子滤波算法对电液舵机轴承的RUL 进行预测，取得了良好的效果。所提出的基于物理模型的FPF 预测电液舵机轴承RUL 的方法具有较好的实时性。然而，物理模型的参数会影响预测RUL 的精度，如何更好地辨识模型参数，从而获得更精确的轴承RUL 值，需要进一步的研究。