一种基于粒子滤波的性能状态百分位值估计方法

一种基于粒子滤波的性能状态百分位值估计方法
【技术领域】
1.本发明属于可靠性领域，具体涉及一种基于粒子滤波的性能状态百分位值估计方法。


背景技术：

2.在工业和工程领域中，随着产品和设备运行时间的增加，其零部件和结构不可避免会发生退化或损耗，最终导致设备的失效。而在一些关键的任务中，设备的失效将引发重大的事故，带来巨大的经济甚至人员损失，因此为保障产品和设备的正常运行，对它们的健康状态和可靠性进行准确评估尤为重要。
3.传统的可靠性分析方法往往基于失效时间数据来分析设备的可靠性，通过选取适量样本进行寿命试验，得到大量失效时间数据，然后根据失效时间的分布建立相应的模型，从而实现设备可靠性的评估。而近些年来，随着结构的精细化和设计的智能化，以及制造技术和材料工艺的发展，设备越来越呈现出高成本，长寿命的特点。对于这类设备而言，受时间、费用等因素限制，难以获得大量失效时间数据，导致传统的可靠性分析方法不再适用。随着传感器应用的不断增加，基于状态监测的可靠性实时评估技术越来越受到人们的关注。该方法基于设备性能状态的在线监测信息，通过判断性能状态是否到达失效阈值来评估设备的可靠性，能够实现单机设备的可靠性实时评估。但是在实际中，在线监测信息往往受到量测误差的影响，并且这个量测误差还往往是随着时间变化的，给性能状态的估计精度带来了重大影响。
4.因此，本发明提出了一种基于粒子滤波(particle filter,pf)的性能状态百分位值估计方法，通过滤波融合的方法克服了量测误差带来的影响，并据此给出了性能状态百分位值的估计方法。


技术实现要素：

5.本发明的目的是提供一种基于粒子滤波的性能状态百分位值估计方法。通过pf实现量测误差和性能状态的同时估计，从而提高性能状态的估计精度，并且基于估计得到的量测误差信息，推导了性能状态百分位值，为当前设备可靠性的估计提供了重要依据。
6.本发明为一种基于粒子滤波的性能状态百分位值估计方法，它包含以下六个步骤：
7.步骤一：针对量测噪声标准差随时间增加的系统，建立其状态方程和量测方程，并基于状态扩维改写系统状态方程和量测方程。
8.对于一个设备系统，基于其物理原理而建立的反映其性能状态变化的状态方程如式(1)所示：
9.xk＝f(x
k-1
)+w
k-1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
10.式中，k为时间步数；xk表示系统的性能状态，该参数是反映系统健康状态的指标。f(
·
)是一个确定性的函数，w
k-1
表示系统噪声，其均值为0，协方差为q
k-1
。
11.系统的量测方程表示为：
12.yk＝h(xk)+vkꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
13.式中，yk表示系统性能状态的量测值，h(
·
)是一个确定性的函数，vk表示量测噪声，其均值为0，在实际中，量测噪声的标准差往往随着时间增加而增大，可以表示为：
14.σk＝a
·kꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
15.式中，σk表示量测噪声的标准差，a为未知参数。
16.基于式(3)，对状态进行扩维，将系统方程(1)-(2)改写为：
[0017][0018]
式中，ak表示在k时刻对参数a的估计，w
s,k-1
是为了保证算法多样性而设置的虚拟噪声，一般根据实际情况调整，vs为均值为0、方差为1的高斯白噪声。
[0019]
步骤二：滤波初始化
[0020]
根据性能状态和参数的初始概率密度p(x0)和p(a0)随机产生n个初始粒子，记为
[0021]
对于k＝1,2,
…
，执行步骤三到步骤六。
[0022]
步骤三：先验估计，基于状态方程，对所有粒子进行一步时间更新，得到当前时刻性能状态和参数粒子的先验估计
[0023]
性能状态和参数的先验估计如下：
[0024][0025]
式中，每一个噪声w
k-1,i
和w
s,k-1,i
都是分别基于噪声w
k-1
和w
s,k-1
的概率密度函数随机产生的。
[0026]
步骤四：计算似然概率密度，根据当前时刻最新量测值计算每个粒子的似然概率密度并归一化
[0027]
以系统性能状态量测值yk为条件计算每个粒子的似然概率密度qi，这个可以通过以下公式进行计算：
[0028][0029]
随后将qi进行归一化为：
[0030]
[0031]
步骤五：基于粒子的似然概率密度进行重采样，得到性能状态和参数的后验估计
[0032]
s51、在[0,1]的均匀分布上产生一个随机数r，累加归一化后的似然概率密度bi，直到累加的总和大于r，即但此时新的后验粒子被设置为等于先验粒子
[0033]
s52、对于i＝1,
…
,n，重复进行上述步骤s51的操作，直至生成n个后验粒子，完成重采样。此时，性能状态和参数的后验估计可以表示为：
[0034][0035]
步骤六：性能状态百分位值估计
[0036]
对于服从均值为u、方差为σ2的正态分布的系统性能状态，其概率为p的百分位值可以表示为：
[0037]
x
p,k
＝u+u
p
σ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
[0038]
基于步骤五的估计，假设k时刻，系统的性能状态xk服从均值为方差为的正态分布。则性能状态概率为p的百分位值的估计为：
[0039][0040]
式中，u
p
为标准正态偏量。
[0041]
通过上述六个步骤，可以实时地估计系统的性能状态的百分位值进而与给定的失效阈值x
th
进行比较评估系统的健康状态。若不超过阈值x
th
，则说明当前系统可靠的概率不低于p，否则说明当前系统可靠的概率已经低于p。
[0042]
本发明的有益效果是：通过基于粒子滤波的方法，提高了性能状态的估计精度，并据此提出了一种实时的性能状态百分位值估计方法，从而能够实时评估系统的健康状态，避免因为设备失效引发的灾难性事故。
【附图说明】
[0043]
图1为本发明方法流程示意图。
【具体实施方式】
[0044]
下面结合附图和实施例对本发明作详细说明。
[0045]
本发明提出了一种基于粒子滤波的性能状态百分位值估计方法，其流程图如图1所示，它包括以下六个步骤：
[0046]
步骤一：建立系统方程
[0047]
对于可在线监测的设备或系统来说，其健康状态往往可以通过某个性能状态来表征，并且，通过基于系统的物理原理或者经验模型及其失效特性可以建立反映其性能状态
变化轨迹的状态方程。比如对于一个变速箱的轴承来说，冲击和磨损导致的裂纹扩展是其主要失效模式，因此裂纹长度是该系统的主要性能状态，而裂纹长度逐渐扩展这一过程可以用paris模型来表示。状态方程一般如式(11)所示。
[0048]
xk＝f(x
k-1
)+w
k-1
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0049]
式中，k为时间步数，xk表示系统的性能状态，该参数是反映系统健康状态的指标。f(
·
)是一个确定性的函数，w
k-1
表示系统噪声，其均值为0，协方差为q
k-1
。
[0050]
传感器对系统性能状态的监测方式可以分为直接监测和间接监测，根据监测数据与性能状态之间的关系，从而可以得到系统的量测方程，一般表示为：
[0051]
yk＝h(xk)+vkꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(12)
[0052]
式中，yk表示系统性能状态的量测值，h(
·
)是一个确定性的函数，vk表示量测噪声，其均值为0，在实际中，随着系统健康状态的恶化，性能状态和量测的方差往往逐渐增加，这一特征可以用量测噪声标准差随时间增加而增大来表征，如下式所示：
[0053]
σk＝a
·kꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(13)
[0054]
式中，σk表示量测噪声的标准差，a为未知参数。
[0055]
为了同时估计系统的性能状态和量测的噪声特性，对式(11)和(12)中的状态进行扩维，从而将系统方程改写为：
[0056][0057]
式中，ak表示在k时刻对参数a的估计，w
s,k-1
是为了保证算法多样性而设置的虚拟噪声，一般根据实际条件和计算收敛情况调整其大小，vs为均值为0、方差为1的高斯白噪声。
[0058]
步骤二：滤波初始化
[0059]
在初始时刻，根据性能状态和参数的初始概率密度p(x0)和p(a0)随机产生n个初始粒子，记为粒子数目n选取得越多，结果越精确，但计算量越大，一般根据精度和时间需求的平衡进行选取。
[0060]
对于时间步数k从1开始，直到最终时刻，重复执行步骤三到步骤六。
[0061]
步骤三：先验估计
[0062]
在前一时刻的后验估计结果的基础上，对所有粒子进行更新，得到当前时刻性能状态和参数的先验估计，如式(15)所示：
[0063][0064]
式中，每一个噪声w
k-1,i
和w
s,k-1,i
都是分别基于噪声w
k-1
和w
s,k-1
的概率密度函数随机产生的。
[0065]
步骤四：计算似然概率密度
[0066]
在得到当前时刻量测值yk的基础上，通过式(16)计算每个粒子的似然
概率密度：
[0067][0068]
随后通过式(17)将似然概率密度qi进行归一化：
[0069][0070]
步骤五：重采样
[0071]
s51、在[0,1]的均匀分布上产生一个随机数r，累加归一化后的似然概率密度bi，直到累加的总和大于r，即但记录此时先验粒子的编号j。此时新的后验粒子被设置为等于第j个先验粒子
[0072]
s52、对于i＝1,
…
,n，重复进行上述步骤s51的操作，直至生成所有n个后验粒子，完成重采样。基于重采样后的后验粒子，性能状态和参数的后验估计为：
[0073][0074]
步骤六：性能状态百分位值估计
[0075]
基于步骤五的估计，可以假设在当前时刻，系统的性能状态xk服从均值为方差为的正态分布。则性能状态概率为p的百分位值的估计为：
[0076][0077]
式中，u
p
为标准正态偏量。
[0078]
重复运行步骤三到步骤六，直至到达最终时刻。