基于两阶段自适应Wiener 过程的剩余寿命预测方法

2022-03-10郑建飞胡昌华牟含笑

自动化学报 2022年2期

董青郑建飞胡昌华李冰牟含笑

随着高新技术的迅猛发展,现代工业设备正朝着大型化、复杂化和智能化趋势快速发展.这类设备在运行过程中由于受到内部和外部因素的随机影响,性能和健康状态不可避免地呈现下降趋势乃至退化失效,导致无法完成正常任务和功能,进而引发严重事故,造成环境破坏和人员伤亡[1-3].如果能在设备性能退化初期对其进行剩余寿命预测(Remaining useful life,RUL),并基于预测结果确定维修决策的最佳时机,制定相应的备件订购或替换策略,将有效提高设备运行可靠性、降低运行成本.近年来,预测与健康管理技术(Prognostics and health management,PHM)得到广泛关注和应用,已经成为可靠性领域的热点研究方向,而PHM 技术的关键在于预测运行设备的剩余寿命.因此,如何得到精确且符合实际情况的剩余寿命,对切实保障系统的运行安全性、可靠性与经济性具有重要的意义[4-7].

经过几十年的发展,RUL 预测取得了丰硕的理论成果并得到广泛应用.袁烨等[8]将寿命预测研究主要分为模型方法和数据驱动方法.数据驱动方法主要包括统计方法和机器学习方法,基于设备大量退化数据推导退化模型,进而判断超过失效阈值的时间预测剩余寿命,如:最小二乘方法、支持向量机方法和深度学习方法等.模型方法分为物理退化模型(机理建模)和经验退化模型.相比于物理退化模型,经验退化模型能通过随机模型对监测数据建模,进而得到寿命或剩余寿命的概率分布,便于量化寿命或剩余寿命的不确定性,从而为健康管理奠定基础,更加适用于现代复杂工业设备.而在经验退化模型中,Wiener 过程和Gamma 过程是两种最常用随机过程退化建模的方法.Gamma 过程是一种增量非负的单调过程,主要适用于单调退化过程,如磨损过程、疲劳扩展过程.但在工程实际中,设备的退化过程大多为非单调,如锂电池容量退化、惯性平台陀螺仪的退化等.Wiener 过程作为非单调退化过程,凭借良好的数学特性,在RUL 预测和健康管理领域得到广泛应用.

近年来,大多数基于Wiener 过程退化建模方法普遍假设系统在退化过程中是一种遵循单一阶段的随机退化模型.但在工程实际中,由于受到内部因素(如:退化机理突变)或外部因素(如:动态环境、状态切换等)的影响,许多设备的退化特性呈现出两阶段乃至多阶段退化特征[9].例如锂电池[10]开始时经历一个平稳退化期,随着充放电的进行,固体电解质层在电极上的生长以及副反应导致的活性材料的损失,导致锂电池容量在后一阶段迅速衰落;液力耦合器[11]开始时经历一个快速退化期,当到达某一时刻(变点)退化速度明显下降,与之类似的还有半导体激光器[12]、等离子显示板[13-14]等.

对于这种存在变点、呈现两阶段退化特性的设备进行退化建模和RUL 预测,已有不少学者进行了研究和拓展.Ng[12]根据退化数据的两阶段特性,提出一种基于单个变点独立增量的两阶段随机退化模型,并采用期望最大化(Expectation maximization,EM) 算法对模型参数进行估计.Yan 等[15]基于两阶段Wiener 过程模型对液力耦合器进行可靠性校验,并根据赤池信息准则(Schwarz information criterion,SIC)对变点进行辨识.Chen 等[16]改进两阶段线性对数模型来描述滚球轴承的分阶段退化过程,并用贝叶斯方法更新模型参数进行寿命估计.Wang 等[17]提出了一种两阶段退化模型用于轴承退化数据的建模,在第一段假设处于健康状态,在第二段结合卡尔曼滤波和EM 算法进行RUL 估计.Peng 等[18]为了提高RUL 预测的鲁棒性和效率,开发了一种半解析预测模型,该模型可以避免RUL 预测的大幅度波动,自动跟踪不同的退化阶段,并自适应地更新超参数.Zhang 等[19]在两阶段Wiener 过程退化模型的框架下,推导出基于首达时间意义的寿命分布,该模型优势在于充分考虑并量化变点处退化量的不确定性同时能够推广至更具有一般性的多阶段退化模型中.

尽管两阶段以及多阶段退化模型已经取得了一些理论与实际应用成果,但仍存在一些问题有待解决.目前大多数两阶段退化模型(如:Zhang 等[19])都是基于Wang 等[20]所提出的一阶自回归模型进行建模,但该模型存在三点不足:1) 假设噪声项是独立且均匀分布,并且仅适用于均匀测量间隔.由于不是自动测量或根据某些设计方案进行测量等原因,在工程实际中设备退化过程的测量间隔往往是不均匀的;2) 当使用多组同类型退化设备的历史数据或先验信息估计模型未知参数时,必须要求监测数据的测量频率与历史数据中的测量频率相同.否则,历史数据将不再适用;3) 该模型退化建模存在一个潜在假设,即在后一时刻估计的随机参数与前一时刻的随机参数的后验估计完全相等,并且当该模型用于RUL 预测时,使用最新的监测值来更新漂移系数,该漂移系数从最后监测点开始保持不变,直到系统发生故障.这意味着该模型假设可以根据实时监测数据自适应更新漂移系数,但在未来的RUL预测中忽略这种自适应漂移可变性.

针对上述问题,本文提出了一种基于自适应Wiener 过程的两阶段退化模型,突破测量间隔固定和采样频率一致的要求限制,同时考虑对表征退化个体差异性的漂移系数实现自适应更新.在首达时间意义下,推导出两阶段自适应Wiener 过程模型的RUL 分布解析式,结合EM 算法和Kalman 滤波技术对模型参数进行估计和更新,并基于SIC 实现退化变点辨识,最后通过锂电池的实例研究验证了本文所提方法可有效实现两阶段退化设备的RUL预测.

1 两阶段随机退化过程建模

建立两阶段退化模型主要针对退化过程表现为两阶段特征的设备,且整个退化时间内存在一个变点,变点前后的退化率呈现显著差异性.

1.1 两阶段线性Wiener 过程退化建模

令X(t) 表示设备在运行中的退化过程,常用的Wiener 过程模型具有以下形式[21]

其中,λ和σ分别表示退化过程的漂移系数和扩散系数,B(t) 是一个标准Brownian 运动,x0为退化初值,通常假设x0=0 .

基于上述模型和假设条件,对于存在变点的随机退化系统可以建立两阶段Wiener 过程退化模型:

其中,x0表示退化过程初值,xτ表示第二阶段退化初值即变点处的退化量,τ表示变点发生时间;λ1和σ1分别表示第一阶段退化过程的漂移系数和扩散系数,λ2和σ2分别表示第二阶段退化过程的漂移系数和扩散系数.

1.2 两阶段自适应Wiener 过程退化建模

由于上述的Wiener 过程模型存在测量间隔不均匀、测量频率不一致以及在RUL 预测中没有利用实时监测数据自适应更新漂移系数三点不足.因此,本文考虑基于Zhai 等[22]所提如下自适应Wiener过程模型建模:

其中,λ(t) 是一个遵循Wiener 过程且随时间变化的漂移系数.λ0是初始漂移率,k是自适应漂移率的扩散系数,W(t) 是一个独立于B(t) 的标准Brownian 运动.

根据式(1)、(3),可以推导出两者各自的离散状态空间模型:

其中,ηi～N(0,Q),ΔBi=B(ti)-B(ti-1),ΔWi=W(ti)-W(ti-1) .

从上式中,可以看到自适应Wiener 过程离散模型的噪声项比式(4)的Wiener 过程离散模型,增加了一个自适应漂移项,当使用最新的监测值Xi更新漂移系数时,自适应Wiener 过程的漂移项从最后监测点开始仍可以动态变化,直到系统发生故障.此外,如式(4)所示的Wiener 过程模型假定两个监测点之间的漂移系数λi为一个随机游走模型,并且依赖于前一个时刻的漂移系数,当测量间隔发生变化时,可能导致模型参数估计不准确.鉴于此,自适应Wiener 过程增加了一个随时间变化的漂移项.当测量间隔发生变化时,模型漂移部分能动态相应变化,适用于不等间隔测量下的RUL 预测.

为了进一步说明,本文给出式(1)和式(3)的两步预测模型:

从式(6)中观察到,Wiener 过程两步预测模型的噪声项中第一个元素ηi+ηi-1的方差是一步预测模型的2 倍,然而两步预测模型的噪声项中第二个元素除了σ(ΔBi+ΔBi-1) 还有一个附加项ηi-1Δti,与一步预测模型的噪声项第二个元素相比,方差不能构成2 倍关系.因此,当式(1)所示的Wiener 过程应用于测量间隔不均匀的数据时,可能造成估计值不准确,RUL 预测准确度随之下降.相比之下,自适应Wiener 过程的两步预测模型增加了一个漂移项,导致附加项kΔWi-1Δti的影响可忽略不计内容.因此,自适应Wiener 过程两步预测与一步预测可以相互兼容,当测量间隔变化时能解决此类问题带来的影响.

基于上述结论,本文对存在变点的随机退化系统建立两阶段自适应Wiener 过程退化模型:

其中,x0表示退化初值,xτ表示第二阶段退化初值即变点处的退化量,τ表示变点发生时间;λ1(s)和σ1分别表示第一阶段退化过程的漂移系数和扩散系数,λ2(s)和σ2分别表示第二阶段退化过程的漂移系数和扩散系数.

2 两阶段自适应Wiener 过程剩余寿命预测

为描述同批次设备中某一个体的退化过程,体现个体差异性,将退化模型的漂移系数随机化[23],即和.若发生变点时的退化量已知,根据文献[22] 中自适应Wiener过程的寿命分布,可以推导出两阶段自适应Wiener 过程寿命分布的概率密度函数(Probability density function,PDF),如式(9)所示.其中,X0表示退化初值,D表示设备退化的失效阈值,τ表示变点发生时间.

实际中,在变点出现前,变点处退化量准确值Xτ是未知的,为了得到寿命估计值,首先要得到首达时间意义下Xτ的分布形式,即在Xτ＜D条件下经过时间τ,退化量从0 到Xτ的转移概率gτ(Xτ) .因此,要计算退化过程在 (Xτ,∞) 失效概率,需保证退化过程在 (0,Xτ) 上未超过失效阈值.如果gτ(Xτ) 的解析表达式可以得到,则基于全概率公式可推导出寿命分布的PDF.

引理 1[19,22].若退化过程为两阶段自适应Wiener 过程模型,且漂移系数λ1、λ2服从正态分布.如果变点时间τ给定,那么在首达时间意义下的寿命分布的PDF,如式(10)、(11)所示.

推论 1.若已知当前时刻tk的退化状态xk,用lk表示设备剩余寿命,fL(lk) 表示设备RUL 分布的PDF,在随机退化速率λ1和λ2的影响下,可获得首达意义下两阶段自适应Wiener 过程模型RUL的PDF,其形式与首达意义下得到的寿命分布PDF,即与式(10)、(11)类似,具体可分为以下两种情况:

情况 1.当前时刻tk位于变点前,即tk＜τ,此时随机设备退化失效又存在两种情况:1) 失效阈值位于变点前,即lk+tk≤τ;2) 失效阈值位于变点后,即lk+tk＞τ,此时RUL 的PDF,如式(12)所示.

情况 2.当前时刻tk位于变点后,即tk＞τ,此时退化设备RUL 的PDF 为

在推论1 中,变点发生时间为某一固定值,只适用于事先预设情况.在实际中,利用监测数据进行预测时,通常两阶段变点位置在不同情况或不同个体下存在差异性.因此,假设变点时间τ为随机变量来描述这种差异性.在这种情况下,随机退化设备的寿命和剩余寿命PDF 为

其中,p(τ) 为变点发生时间的PDF.由于寿命与剩余寿命分布形式比较复杂,上述的积分难以得到具体解析表达式,故本文考虑采用数值积分方法求解.

3 模型参数估计与更新

3.1 潜在退化状态估计

当前时间可定义为tk,而当前运行设备从时间t0～tk获取的退化数据为x0:k={x0,x1,···,xk},如果此时变点未出现,即tk≤τ,也就是说退化设备仍处于第一阶段且尚无第二阶段退化数据,那么仅需根据收集的退化数据来更新第一阶段模型参数;反之,若变点已经出现,即tk＞τ,那么仅需要更新第二阶段模型参数.

类似地,若tk＞τ,可利用当前运行设备退化数据更新参数λ2,由于第一阶段数据与第二阶段模型无关,因此仅需要数据xτ:k={xτ+1,xτ+2,···,xk}用于更新.

当两阶段自适应Wiener 模型用于RUL 预测时,模型参数a10,a20,P10,P20,均是未知的.对此,本文采用EM 算法对参数自适应估计,使得估计的寿命更好地反映设备当前健康状态.

3.2 基于EM 算法的自适应估计

假设对某个退化设备进行状态监测,监测点为m个,即x={x1,x2,···,xm},其各自对应的监测时间为{t1,t2,···,tm}.同时,本文假设变点发生时间已知,即τ∈{t1,t2,···,tm},那么{x1,x2,···,xτ}表示设备第一阶段的退化数据,{xτ+1,xτ+2,···,xm}表示设备第二阶段的退化数据.

式中,p(x0:k|θ1) 为监测数据x0:k的联合PDF.

然后,基于监测数据x0:k的似然函数,θ1的极大似然估计值可通过最大化似然函数得到,表示为

在本文中,由于漂移系数λ1被视作一个 “隐含状态”.无法使未知参数θ1最大化.而EM 算法可通过最大化联合似然函数p(λ1k,x0:k|θ1) 来估计逼近参数的极大似然估计,估计值可以通过迭代以下两步实现.

通过不断迭代E-步骤和M-步骤直到满足某一收敛条件截止,由此得到对应的参数估计值,一般来说随着迭代次数增加,得到的参数估计值会越来越好.

对于第一阶段的随机退化模型,EM 算法中的联合对数似然函数可以表示为

第二阶段参数估计方法同上,不再赘述.在退化实验中,监测到的设备性能退化数据一般为离散值,变点τ的值通常未知.SIC 是Akaike 信息准则的改进,对变点检测效果良好,下面通过SIC,确定变点τ的值[15].

3.3 变点检测

SIC 是由Schwarz.于1978 年提出,可以判断模型是否存在变点问题.其原理是如果待检测序列存在变点,其样本的熵要大于不存在变点的样本的熵[25].利用SIC 来估计变点的个数和位置是较为简单的,对变点的检测效果良好.其定义为

原假设H0:各参数值相等,表示模型中不存在变点.

为了描述不同设备变点时间的个体差异性,本文假设变点时间τ服从随机分布.相比于其他分布,Gamma 分布能包含其他常见分布,如指数分布等,并且形状参数α越大,Gamma 分布越逼近正态分布,计算也较为容易.因此,本文假设变点时间τ服从Gamma 分布,且形状参数为α、尺度参数为β.

通过SIC 方法,利用退化先验信息可离线估计设备的变点发生时间,通过Gamma 分布的统计形式,进而得到变点时间τ的分布参数及概率分布函数p(τ) .

4 实例验证

在本节,通过引入蒙特卡洛仿真,验证本文基于自适应Wiener 过程所提出的模型优越性.然后,将所提出的模型应用于锂电池容量退化数据中.

4.1 模拟仿真

为了验证本文所提模型能够解决现有两阶段模型不能刻画测量间隔不均匀、测量频率不一致的问题,在此考虑将所提出模型与Zhang 等[19]所提的两阶段模型作比较.

首先,在这里主要考虑漂移参数的随机效应,利用两种退化模型分别对仿真得到的数据进行退化建模,并利用第四节模型参数估计方法,可求得两个阶段模型参数估计值,最后得到两种方法各自预测的RUL 结果.在此,本文增加一个仿真验证的例子,设定仿真参数值为 (λ1,λ2,σ1,σ2,k1,k2)=(2,1,1,1,1,0.5),并且生成一组间隔为1、次数为200 的监测数据.为了证明本文所提模型比现有模型适用于测量间隔不均匀、测量频率不一致情况,考虑将数据变为以下情况:存在测量间隔为1、2、4 的测量间隔不均匀的混合数据.其中,变点发生时间设为120,失效阈值为155.

图1 为本文所提模型和Zhang 等所提的两阶段模型在不均匀退化数据下的RUL 结果.从图中,可以看到RUL 预测的PDF 随着时间逐渐变窄,表明预测的不确定性越来越小.对于测量间隔不均匀的数据,所提模型的PDF 曲线更加窄而尖锐,说明在寿命预测方面可以取得更好的预测结果.

图1 两种模型RUL 预测结果Fig.1 RUL prediction results of the two models

为了进一步量化测量间隔不均匀条件下两种模型的预测结果,本文采用可靠性领域常用的性能指标:绝对误差(Absolute error,AE)和相对误差(Relative error,RE)指标.

从图2和表1 中,可以看出随着退化数据的逐渐累积,两种模型各自预测的误差值也在随之下降.总体上看,本文所提模型的预测准确度要优于Zhang 等所提的模型.当退化数据呈不均匀分布时,相比较Zhang 等所提的两阶段模型,本文所提模型能较好解决这种情况带来的影响,更为准确地估计参数,进而提高预测准确度.

图2 两种模型RUL 预测绝对误差Fig.2 Absolute error of RUL prediction of the two models

表1 不同监测点相对误差的比较结果Table 1 Comparison results of relative error at different monitoring points

综上,本文基于自适应Wiener 过程所建立的模型相比于Zhang 等[19]的两阶段退化模型,可以克服测量间隔不均匀、测量频率不一致问题,较为准确预测设备RUL.

4.2 锂电池实例验证

本节中,基于马里兰大学Pecht[26]课题组的锂电池容量退化数据进行实例验证.该数据是在室温条件下通过充放电实验得到的,记录了电池状态信息(包括容量)随充放电次数的变化.由于复杂的化学反应,锂电池容量开始时经历一个平稳退化期,随着充放电的进行,由于固体电解质层在电极上的生长以及副反应导致的活性材料的损失,导致锂电池容量在后一阶段迅速降低[9].编号为CS2-35、CS2-36、CS2-37、CS2-38 四组电池容量退化数据如图3所示,图中退化过程呈现出明显的两阶段特性,这里采用CS2-37 锂电池数据进行RUL 预测验证,其余三组(即CS2-35、CS2-36和CS2-38)用作变点时间离线参数估计.

图3 锂电池容量退化Fig.3 Lithium battery capacity degradation

首先,对CS2-37 锂电池利用SIC 进行变点检测,确定变点所在时刻.根据SIC,分别计算相应的SIC 值,对应的SIC 值变化趋势如图4 所示.

图4 CS2-37 锂电池SIC 值Fig.4 SIC value of CS2-37 lithium battery

从图中可以观测到,第96 监测点数值最小,由于本文采用两阶段随机退化建模,从起始退化到第96 监测点,监测时间较短不适宜作为变点,因而,对于此情况忽略不计.此时,SIC(1 006)=-16 735＞SIC(736),且S IC(736)=-26 891 最小,所以CS2-37 锂电池退化中存在变点,且变点发生在第736 监测点.同理,其余三组电池数据辨识到的变点分别为623、681和753.由于变点时间τ服从Gamma分布,可利用辨识得到的数据进行拟合,可求得形状参数和尺度参数分别为α=106、β=7.98 .将变点引入到参数估计中,再结合退化数据,可得两个阶段漂移系数的均值和方差估计值分别为ua=8.056×10-4、σa=2.94×10-5、ub=0.00221、σb=7.41×10-5.利用上述估计值,结合CS2-37 锂电池数据和卡尔曼滤波技术实现漂移系数在线更新.图5、图6 展示了隐含状态即漂移系数的在线更新过程.结果表明,两个阶段容量的退化率相差较大,如果只用单一阶段的退化过程进行退化建模,误差将会比较大.由于变点发生时刻在第736 监测点,因此第一阶段模型参数在变点出现后不再更新,第二阶段模型参数在变点之前不进行更新.

图5 第一阶段模型参数更新Fig.5 The first stage of model parameter updating

图6 第二阶段模型参数更新Fig.6 The second stage of model parameter updating

为了验证所提模型的有效性,用建立模型对CS2-37 锂电池退化数据进行拟合.图7 为符合Wiener 过程的锂电池容量退化预测情况,从图中可以看出所建立的模型能较好地反映锂电池退化过程.

图7 CS2-37 模型拟合效果Fig.7 Fitting effect of CS2-37 model

由于变点在第736 监测点,而一般电池容量的失效阈值定义为电容量损失到初始容量的百分比.在本文中设定失效阈值为初始电容量的45 %.为了说明本文方法的有效性和合理性,将估计的模型参数代入到式(8)中,可得到RUL 的PDF和预测值,如图8 所示.

图8 CS2-37 的RUL 预测结果Fig.8 RUL prediction results of CS2-37

从图8 中可以看出,对于锂电池退化数据的RUL 预测,所提模型与单一阶段线性模型相比较,前者更为准确.与Zhang 等提出的两阶段模型相比,本文的方法能取得更好的结果,且随着监测数据的增加,RUL 预测结果不确定性越来越小,精度越来越高.为了更加直观说明本文方法有效性,给出几种方法的RUL 预测绝对误差和α-β性能指标对预测结果进行验证.

从图9、图10 中可以看出,与单一阶段线性模型相比,本文方法考虑变点前后呈现两阶段特征,即变点前后的退化速率存在明显差异进行建模,且考虑同批产品个体差异性的影响,其模型更符合实际退化情况.与Zhang 等所提的两阶段模型相比,本文方法考虑Wiener 过程模型存在测量间隔不均匀、测量频率不一致以及在RUL 预测中忽略了自适应漂移的可变性等三点不足,结果表明在监测前期退化数据较少时,所提模型能取得较好的预测结果.其原因是本文所提模型与Zhang 等所提模型相比,噪声项增加了一个自适应漂移项它是一个随时间变化的过程.在设备监测初期,剩余寿命lk值比较大,无法忽略自适应漂移可变性的影响,因此所提模型预测结果优于现有模型.随着电池充放电循环在寿命将尽时,自适应漂移项的影响降到最低,此时两种退化模型结构相似,进而提供近似的预测结果.