王玺，胡昌华,*，任子强，熊薇

1. 火箭军工程大学 导弹工程学院，西安 710025 2. 国防大学 联合勤务学院，北京 100089

航空发动机作为飞机的“心脏”，不仅为飞机提供动力，同时也为航空事业的发展注入强大动力。然而，航空发动机由于长期运行在复杂多变的恶劣环境下，其性能会随之产生衰减，并在累积到某种量级时，导致发动机失效[1]。目前，航空发动机在航空领域发挥着重要的作用，一旦发生由失效引起的事故，所造成的后果将难以估计。因此，通过监测的航空发动机性能衰减数据，建立相应的性能衰减模型，进而实现航空发动机的剩余寿命(Remanining Useful Life,RUL)预测是近年来热门的一项关键技术[2]。然而，航空发动机的寿命有多种定义，如首翻寿命、大修寿命、总寿命等，同时又受多种因素约束，如性能衰减、关键件寿命等。鉴于此，本文主要研究的是在首翻寿命的定义下，基于性能衰减建模的航空发动机剩余寿命预测。

目前，基于性能衰减建模的剩余寿命预测方法已经广泛地应用于航空发动机[3-5]、导弹[6-8]、电机[9]等。例如，黄亮等[3]基于多阶段Wiener过程建立了航空发动机的性能衰减模型，通过融合发动机的历史性能衰减数据和实时监测数据，实现对发动机的剩余寿命预测。王浩伟等[6]针对导弹多个部件间的竞争退化，综合考虑了退化失效和突发失效2种竞争模式，提高了对导弹剩余寿命预测的精度。袁庆洋等[9]针对电机的退化特性，建立了多段Wiener过程模型，实现了电机的剩余寿命预测。

然而，对于航空发动机这类复杂的机械系统，在性能衰减过程中普遍存在着三源不确定性，即时变不确定性[10]、个体差异不确定性[11]、测量不确定性[12]。时变不确定性是指性能衰减过程随着时间变化的固有不确定性；个体差异不确定性是指同类型设备之间的性能衰减差异，即同类型设备的不同个体在性能衰减过程中虽然具有相似的衰减路径，但实际的衰减率是不相同的；测量不确定性是指测量的性能衰减数据与实际的性能衰减状态之间的差异。在实际中，由于受噪声、干扰、仪器精度等影响，想要对性能衰减状态的精确测量几乎是不可能的，因此测量中难免会引入测量误差，造成测量的不确定性。

目前，有关同时考虑三源不确定性的研究较少。为此，文献[13]基于线性的Wiener过程建立了同时考虑三源不确定性的性能衰减模型，并使得预测的剩余寿命分布可以随着性能衰减数据在线更新。然而，对于航空发动机而言，在性能衰减过程中不仅具有三源不确定性，还具有强非线性。鉴于此，文献[14]在非线性下考虑了三源不确定性，并且证明了同时考虑三源不确定性可以有效地提高剩余寿命预测的精度。文献[15]在此基础上，利用粒子滤波实现了对锂电池的剩余寿命预测。但是，以上有关考虑三源不确定性的工作中存在一个共性问题，即假定有多组同类型设备的历史性能衰减数据来离线地辨识模型初始参数，并且这些参数一旦确定，就不再随着新的性能衰减数据进行更新。因此，剩余寿命预测的精度受限于数据量，并且预测的结果在很大程度上依赖于选择的模型初始参数。然而，针对新研发的航空发动机，往往是缺乏此类相关的历史性能衰减数据和先验信息，因此，以上方法具有一定的局限性，且目前相关的研究较少。

对于目前航空发动机剩余寿命预测存在的问题，本文建立了同时考虑三源不确定性和非线性的性能衰减模型，并在首达时间下推导出剩余寿命的分布。然后，针对新研发的航空发动机，提出了一种基于Kalman滤波和条件期望最大化(Expectation Conditional Maximization,ECM)算法的参数估计方法，以克服缺乏历史性能衰减数据和先验信息的问题，并当新的监测数据可用时，能够自适应地估计模型参数，进而更新剩余寿命的分布，实现了对新研发航空发动机剩余寿命的在线预测。最后，通过实例验证了本文方法可以有效地提高预测的精度并降低预测的不确定性，且更加适用于新研发的航空发动机。

1 模型描述

针对航空发动机在性能衰减过程中存在的非线性和不确定性，本文采用一类一般性的非线性Wiener过程来建模航空发动机的性能衰减过程。具体地，性能衰减过程{X(t),t≥0}由标准Brown运动B(t)驱动，可以表示为

(1)

Y(t)=X(t)+ε

(2)

在以上框架下，本文的主要目标是针对缺乏历史数据和先验信息的新研发航空发动机，利用实时监测的性能衰减数据来自适应地预测航空发动机的剩余寿命，并且在获得新的性能衰减数据后实现对剩余寿命的在线更新。具体地，首先假定获得的性能衰减数据是在离散时间点0=t0

(3)

2 剩余寿命自适应预测

考虑到表征设备健康水平的潜在性能衰减状态{X(t),t≥0}一旦首次达到设定的失效阈值ω(一般由工业标准设定)时，则认为设备失效。因此，基于首达时间的概念[1-5]，设备在tn时刻的剩余寿命Ln定义为

Ln=inf{ln>0:X(ln+tn)≥ω}

(4)

相应地，剩余寿命Ln的概率密度函数(PDF)为fLn|Y0:n(ln|Y0:n)。因此，预测航空发动机剩余寿命的关键是基于监测的性能衰减数据Y0:n={y0,y1,…,yn}，求解fLn|Y0:n(ln|Y0:n)，并在获得新的性能衰减数据后，实现剩余寿命分布的在线更新。

为了实现剩余寿命的自适应预测，首先将式(3)进一步改写为

(5)

(6)

式中：δn～N(0,φn)，且有

定义基于监测数据Y0:n估计的zn的期望和协方差分别为

(7)

同样地，定义一步前向的期望和协方差为

(8)

状态估计

(9)

协方差更新

Pn|n=Pn|n-1-KnCPn|n-1

(10)

滤波初值

(11)

根据建立的状态空间模型式(6)，基于监测的性能衰减数据Y0:n，计算剩余寿命的PDF为

fLn|Y0:n(ln|Y0:n)≅

(12)

3 基于ECM算法的模型参数估计

(13)

根据极大似然估计(MLE)理论，性能衰减数据Y0:n关于模型参数向量Θ的对数似然函数为

(14)

(15)

然而，由于将潜在性能衰减状态xn和随机参数λn视作“隐含”状态，因此无法直接最大化式(15)。为解决这一问题，自然会采用期望最大化(EM)算法进行参数估计，通过最大化联合似然函数p(zn,Y0:n|Θ)去逼近模型参数的MLE。EM算法的估计过程可以分为以下2步。

E步：

(16)

M步：

(17)

首先，为了计算对数联合似然函数，可以通过条件概率公式得到

2lnp(zn,Y0:n|Θ)=

2[lnp(Y0:n|zn,Θ)+lnp(zn|Θ)]=

(18)

然后，通过E步对式(18)求条件期望，有

(19)

式中：

(20)

显然，为了推导出E步结果，需要先计算式(20)。这里，本文采用Rauch-Tung-Striebel(RTS)平滑滤波算法[19]计算Ω、Ξ、Φ：

平滑滤波

(21)

协方差初值

Mn|n=(I2×2-KnC)AnPn-1|n-1

(22)

后向迭代

(23)

式中：Mi|n=cov(zi,zi-1|Y0:n)。基于RTS算法，条件期望Ω、Ξ、Φ可以计算为

(24)

(25)

(26)

(27)

(28)

CM2：为了基于性能衰减数据Y0:n估计参数θ，首先将式(25)和式(28)一起代入到式(19)中，便可以通过最大化来估计参数即

(29)

由此可见，与传统的EM算法相比，本文的ECM算法通过CM1和CM2两步，可以将原本需要优化的3个参数简化为一个参数θ，降低了计算复杂度，提高了在线性能。此外，本文提出的参数估计方法可以在监测到一个新的性能衰减数据时，实现对模型参数的自适应估计和在线更新，克服了新研发航空发动机历史数据和先验信息不足的限制，这对于工程实际具有重要意义。

4 实例验证

根据相关文献可知，航空发动机的性能衰减和故障通常表现在排气温度裕度(EGTM)的衰减[22-23]。因此，EGTM是表征航空发动机健康状态的一个重要指标。EGTM的定义为航空发动机出厂前测定的排气温度(EGT)限制值与全推力(或全功率)起飞时实际的EGT之间的差值。随着航空发动机工作时间的增加，EGT会逐渐升高，反映到EGTM就会逐渐衰退，因此准确预测EGTM对提高航空发动机的可靠性至关重要。

实际中，由于外部环境和内部因素的随机性，航空发动机在某一时刻的性能衰减量也是随机的，而且同类型航空发动机EGTM的衰减轨迹虽然具有相似性，但实际上不同个体之间仍然存在着差异性。因此，监测的EGTM数据不仅会受到时变不确定性的影响，还会受到个体差异不确定性的影响。此外，EGT数据来源于航空发动机上安装的传感器，经过空地数据链、报文解析程序等步骤进入到性能监测系统，在这整个过程中不可避免地会引入测量误差，造成测量的不确定性。根据EGTM的定义可知，EGTM是基于EGT计算的，因此EGTM数据中必然存在测量误差。通过上述分析发现，监测的EGTM数据中存在着三源不确定性共存的现象，因此同时考虑三源不确定性是必要的，也是更加符合实际需求的。下面以某型号航空发动机在6年内监测的EGTM数据来验证本文方法的有效性，收集的性能衰减数据如图1所示。根据工程经验，设定EGTM的失效阈值ω=0 ℃，并且考虑以下幂函数形式的性能衰减模型：

图1 航空发动机的性能衰减路径

X(t)=x0+λ·tθ+σBB(t)

(30)

这里，根据EGTM数据的特点，选择式(1)中非线性函数的形式为ρ(t;θ)=θ·tθ-1。到目前为止，此类模型已经在航空发动机、轴承、电池等复杂系统的性能衰减建模中得到了广泛的应用[3,12,17]。

如图1所示，监测间隔为100次循环，在监测到第4 000次循环时，EGTM首次达到失效阈值ω，期间共监测了41次。显然，随着循环次数的增加，航空发动机的性能衰减路径整体呈下降趋势。然而，本文所提出的剩余寿命预测方法是基于性能衰减路径呈上升趋势的情况下推导出来的，因此需要对原始性能衰减数据进行转换。这里，本文利用EGTM的初始值减去其他所有的性能衰减数据，相应地，失效阈值转换为初始值减去ω。下面，基于转换后的性能衰减数据来验证本文方法的有效性和优越性。随着性能衰减数据的积累，采用所提出的参数估计方法可以对模型参数Θ在每一个监测点进行自适应地估计和更新，结果如图2所示。

图2 模型参数的自适应估计过程

图2的结果表明，模型参数可以随着性能衰减数据的积累很快地收敛，并且在获得一个新的性能衰减数据时，模型参数和剩余寿命的PDF都可以进行自适应地估计和更新。因此，克服了历史数据和先验信息不足的限制，能够适用于新研发航空发动机的剩余寿命预测。

为了进一步证明本文方法的优越性，选择以下2种方法进行对比研究：

图3给出了3种方法在最后6个监测点预测的剩余寿命的PDF曲线，从图中可以看出，3种方法得到的PDF曲线均包含了实际的剩余寿命，且随着循环次数的增加逐渐变窄，这意味着预测的剩余寿命的不确定性越来越低。此外，由于本文方法在性能衰减建模中同时考虑了三源不确定性，并且在每一个监测点可以对模型参数进行自适应地估计和在线更新，因此本文方法得到的PDF曲线更加窄而尖锐，说明本文方法在剩余寿命预测方面具有更高的精度和更好的性能。

图3 三种方法的比较结果

为了进一步量化比较结果，本文采用可靠性领域中常用的均方误差(MSE)指标，该指标同时考虑了剩余寿命预测的精度和不确定性，在tn时刻的定义式为

(31)

图4 所有监测点的MSE

图5 4个不同监测点的箱形图(点划线表示实际剩余寿命)

图4给出了3种方法在所有监测点的MSE比较结果。由于本文方法和M2方法均不依赖于历史数据和先验信息，且模型的初始参数是随机给定的，因此在监测初期性能衰减数据较少时，剩余寿命预测的MSE高于M1，但随着性能衰减数据的积累，本文方法的MSE明显低于M1和M2。为进一步验证，图5给出了4个不同监测点的箱形图。其中，箱子的大小表示预测的剩余寿命的不确定性，若箱子越大，则预测的不确定性越大，精度越低；点划线表示当前监测点的实际剩余寿命；箱子中的实线表示预测的剩余寿命的期望，实线越接近点划线，说明预测的精度越高，结果越准确。从图5(a)中可以看出，此时由于性能衰减数据较少，M1优于本文方法和M2，但是随着性能衰减数据的积累，本文方法的箱子更小且更加接近实际的剩余寿命，如图5(b)～图5(d)所示。此外，引入相对误差(RE)指标来进一步量化剩余寿命预测的精度，在寿命的20%、45%、70%和95%分位点上分别给出了3种方法RE的比较结果，如表1所示。

通过表1的结果可以看出，本文方法可以有效减少剩余寿命预测的RE，进而提高了剩余寿命预测的精度，尤其在寿命的95%分位点处，剩余寿命预测的RE仅为1.68%，明显小于其他2种方法。以上的实验结果表明，本文方法具有更好的预测能力，剩余寿命预测的精度更高，且不依赖于历史数据和先验信息，更加适用于新研发的航空发动机。

表1 不同寿命分位点相对误差的比较结果

5 结 论

1) 本文在三源不确定性下建立了非线性的性能衰减模型，不仅克服了现有方法的潜在假设，而且可以较好地建模航空发动机的性能衰减轨迹，并能够将三源不确定性的影响纳入到预测的剩余寿命分布中。

2) 利用Kalman滤波技术可以对潜在的性能衰减状态进行实时估计，并基于首达时间的概念，得到了同时考虑三源不确定性的剩余寿命分布。

3) 在获得一个新的性能衰减数据时，基于RTS平滑和ECM算法对模型参数进行自适应地估计和在线更新，克服了新研发航空发动机缺乏历史数据和先验信息的问题，降低了参数估计的不确定性，显著提高了剩余寿命预测的精度，尤其在寿命的95%分位点处，剩余寿命预测的RE仅为1.68%。