CN105912843A - 基于两视角和粒子滤波的非监督故障预报方法 - Google Patents

Abstract

本发明针对新设备或新系统没有过往的设备或部件的失效数据（历史故障数据）情况，提出一种两视角方法结合粒子滤波进行非监督的故障预报方法。一个视角是分析观测数据的变化范围来给出健康指标，该指标评估系统的当前数据与已知的健康数据之间变化范围的一致性；另一个视角关注系统的预测数据对健康数据的偏离，并设计基于这种偏离程度的潜在故障指标。预测数据的获取是通过状态空间模型描述系统的故障演化过程，然后利用粒子滤波算法求出未来的系统运行数据。在当前监测时刻，若潜在故障指标连续大于健康指标，则给出故障即将来临的警报（即故障预报）。本发明可实现故障的早期预报，是一种有效的故障预报方法。

Description

基于两视角和粒子滤波的非监督故障预报方法

技术领域

本发明涉及一种基于两视角和粒子滤波的非监督故障预报方法。

背景技术

故障预报对于避免事故发生和保障设备或系统的正常运行具有重要的现实意义。对新设备或新系统而言，通常没有过往的设备或部件的失效数据(历史故障数据)，此种情况下对即将发生故障进行预报是一个非监督的任务。具有潜在故障的动态系统可以用状态空间模型描述，而模型中系统的当前和未来状态则采用粒子滤波算法求解，基于当前观测和未来预测数据的故障预报方法给出故障预报。由于粒子滤波方法具有系统是非线性或非高斯噪声情况都适用的优点，目前国内外研究人员已将粒子滤波方法用于故障检测和预报之中，例如采用随机摄动粒子的粒子滤波方法、基于高斯混合模型的粒子滤波方法和使用转移概率核的粒子滤波方法已用于三容水箱系统的故障预报，粒子滤波结合状态空间方程的方法已用于飞行器齿轮盘裂纹故障的检测，粒子滤波结合衰减状态模型的方法已用于轴承故障的检测，并行计算的粒子滤波方法已用于移动机器人的故障检测，以及粒子滤波结合触点温升模型的方法已用于电力设备载流故障预测等。对于从两个视角设计方法再合作进行故障预报，即视角一分析系统运行的当前数据与已知的健康数据之间变化范围的一致性，视角二考虑系统运行的预测数据对健康数据的偏离，这样的两视角合作的非监督故障预报，是上述故障预报/检测方法和现有的文献都没有论及的。

发明内容

有鉴于此，本发明的目的在于提供一种基于两视角和粒子滤波的非监督故障预报方法，可成功地实现故障的早期预报。

为实现上述目的，本发明采用如下技术方案：一种基于两视角和粒子滤波的非监督故障预报方法，其特征在于包括以下步骤：

步骤S1：系统初始化：采用状态空间模型描述由状态变量x和观测变量y构成的系统，每个时刻t测量并获得观测数据D₀＝{y_t}，其中t＝1,2,3,…,m；m为正整数；移除D₀＝{y_t}中最大的5％数据及最小的5％数据后计算其均值μ₀和方差σ₀，在所述状态空间模型的基础上采用粒子滤波算法估计所述系统的状态并获得和其中t＝1,2,3,…,m；设定步长为q，并令k＝m+1；

步骤S2：健康指标的计算：视健康系统的观测数据近似于高斯分布p(y)＝N(μ₀,σ₀ ²)，从时间点k-q+1到k，将观测数据代入分布函数p(y)并求p(y_t)的均值：

$E_k({\mathrm{\μ}}_0,{\mathrm{\σ}}_0)=\frac{1}{q\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_0}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{{\mathrm{\σ}}_0}^2}\}$

在时刻k，计算所述系统的健康指标：

$H_k=\frac{E_k({\mathrm{\μ}}_0,2.55{\mathrm{\σ}}_0)}{p\left({\mathrm{\μ}}_0\right)}=\frac{1}{q}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{\left(2.55{\mathrm{\σ}}_0\right)}^2}\}$

步骤S3：潜在故障指标的计算：在时刻k，在所述状态空间模型的基础上采用粒子滤波算法计算并计算预测数据和其中t＝k+1,k+2,…,k+q；计算预测数据对健康范围的偏离程度，偏离程度的计算公式为：

设Z表示第一次发生故障时的偏离，β为指数函数的均值参数，则在某个偏离z下的无故障概率为：

P(Z＞z)＝exp{-z/β}

在时刻k，当用已获得的预测数据代入偏离公式时，计算所述系统的潜在故障指标：

$F_k=1-\frac{1}{q}\underset{j=k+1}{\overset{k+q}{\Σ}}\exp \{-z\left({\hat{y}}_j\right)/\mathrm{\β}\}$

步骤S4：故障预报：在时刻k，其中k>m+2，如果F_k>H_k且F_k-1>H_k-1且F_k-2>H_k-2，给出即将来临故障的预报，否则不预报；k＝k+1，转至步骤S2。

进一步的，所述步骤S3中的健康范围为M₀＝[μ₀-3σ₀,μ₀+3σ₀]。

本发明与现有技术相比具有以下有益效果：

1.本发明从两个角度分析系统运行的异常征兆，即健康指标评估当前观测数据属于系统运行的健康范围的程度，潜在故障指标评估系统运行的预测数据对健康运行数据的偏离程度；

2.本发明基于预测数据的潜在故障指标与基于观测数据的健康指标这两个角度合作进行故障预报。

附图说明

图1是本发明算法流程图。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

请参照图1，本发明提供一种基于两视角和粒子滤波的非监督故障预报方法，其特征在于包括以下步骤：

步骤S1：系统初始化：采用状态空间模型描述由状态变量x和观测变量y构成的系统，在每个时刻t测量并获得观测数据D₀＝{y_t}，其中t＝1,2,3,…,m；m为正整数；移除D₀＝{y_t}中最大的5％数据及最小的5％数据后计算其均值μ₀和方差σ₀，在所述状态空间模型的基础上采用粒子滤波算法估计所述系统的状态并获得和其中t＝1,2,3,…,m；设定步长为q，并令k＝m+1；

步骤S2：健康指标的计算：从第一个视角，我们聚焦于观测数据的变化范围来分析系统的健康状态。假设健康系统是一个近似的平稳随机过程，即观测时间序列具有近似固定的均值和方差，则健康状态下的观测数据将在一个近似固定的范围内变动(称为健康范围)，设计一个健康指标来测量当前观测数据属于健康范围的程度。

设系统的初始阶段是健康阶段，其观测数据为D₀＝{y_i}(i＝1,2,3,…,m)；当{y_i}围绕其均值μ₀变化时，D₀可视为近似来自于高斯分布p(y)＝N(μ₀,σ₀ ²)，y的概率密度函数为：

$p\left(y\right)=\frac{1}{\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_0}\exp \{-\frac{{(y-{\mathrm{\μ}}_0)}^2}{2{{\mathrm{\σ}}_0}^2}\}$

从离散时间点k-q+1到k(共q个时间点)，将所述观测数据{y_t}代入分布函数p(y)并求p(y_t)的均值：

$E_k({\mathrm{\μ}}_0,{\mathrm{\σ}}_0)=\frac{1}{q\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_0}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{{\mathrm{\σ}}_0}^2}\}$

我们知道满足分布N(μ₀,σ₀ ²)的数据D₀中的99.7％落入范围M₀＝[μ₀-3σ₀,μ₀+3σ₀]，这里将M₀视为健康范围。这样，系统健康运行时的数据D₀被设置为基准，则我们可以通过与基准数据的对比来评估当前系统的健康状态。具体地，当前观测数据落入M₀时，系统运行在健康状态，其健康程度用数值0.5至1来描述(观测数据越靠近μ₀，值越接近1)；与此对照，当一些或所有观测数据跑出M₀时，其健康程度用数值0至0.5(小于0.5)来描述；特别地，当观测数据都远离M₀时，其健康程度值应当靠近0。

这种健康程度的评估适合用概率密度函数p(y)来实现，我们可以发现：离M₀的中心越远，p(y)的值越小；p(y)/p(μ₀)具有值域(0,1]；为在M₀的边界μ₀±3σ₀处满足健康值0.5，p(y)/p(μ₀)中的σ₀被扩大为2.55σ₀，于是p(μ₀±3σ₀)/p(μ₀)＝0.5006。按照上述思路，我们利用概率密度函数p(y)和当前观测数据{y_i}(i＝k-q+1,k-q+2,…,k)设计一个健康指标来评估一个系统的健康状况。

当E_k中的σ₀被扩大为2.55σ₀且E_k用p(μ₀)规范化之后，在时刻k，所述系统的健康指标H_k∈(0,1]被定义为：

$H_k=\frac{E_k({\mathrm{\μ}}_0,2.55{\mathrm{\σ}}_0)}{p\left({\mathrm{\μ}}_0\right)}=\frac{1}{q}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{\left(2.55{\mathrm{\σ}}_0\right)}^2}\}$

步骤S3：潜在故障指标的计算：从另一个视角，我们关注预测数据对健康数据的偏离，并设计基于这种偏离程度的潜在故障指标。为更早地发现这种偏离(故障征兆)，一种可行的方法是对系统的未来状态进行预测并对预测数据进行分析。为此，第一项任务是在设计潜在故障指标之前获得预测数据。因此在时刻k，在所述状态空间模型的基础上采用粒子滤波算法计算并计算预测数据和其中t＝k+1,k+2,…,k+q；具体的：通过粒子滤波算法估计系统状态求得在当前时刻k对系统状态方程使用粒子滤波算法，可进行q-步预测给出未来状态依据系统状态方程可获得预测数据

当前观测数据偏离健康范围M₀，这是动态系统中异常行为的征兆。假设这异常行为来自于系统部件的磨损，则偏离越多，预示部件的磨损越严重，故障即将来临的概率越大。

在可靠性理论中，设备或部件的失效常常符合一个泊松过程，直到第一个事件发生(例如故障发生)时的时间长度(T)超过时间t的概率可以用指数函数P(T＞t)＝exp{-t/β}来描述，这是在时刻t不发生事件的概率。在系统连续出现异常行为之后，随着时间的增长，部件的磨损程度具有增长的趋势。这样，在异常行为连续出现之后，时间增长近似对应于偏离的增长，或者说时间与数据偏离之间有正相关关系。因此，在某个偏离值下无事件(例如无故障)的概率可以近似采用相同的指数函数来描述。

当观测数据D₀＝{y_t}偏离健康范围M₀时，这种偏离被设计为观测值y与μ₀±3σ₀之间的差值z，并用σ₀进行规范化后的偏离z为：

上述偏离公式也适用于计算预测数据对健康范围的偏离程度，即将预测数据取代y代入偏离公式计算

设Z表示第一次发生故障时的偏离，β为指数函数的均值参数，则在某个偏离z下的无故障概率为：

P(Z＞z)＝exp{-z/β}

考虑到3σ₀为M₀边界到中心μ₀的距离，特做如下设计：指数分布从M₀边界出发，伸展到3σ₀作为指数分布的均值，即β＝3。在当前时刻k使用预测数据和我们设计一个基于指数函数和偏离数据的潜在故障指标。

在当前时刻k，当用已获得的预测数据代入偏离公式z(y)时，潜在故障指标F_k∈[0,1)被定义为：

$F_k=1-\frac{1}{q}\underset{j=k+1}{\overset{k+q}{\Σ}}\exp \{-z\left({\hat{y}}_j\right)/\mathrm{\β}\}$

步骤S4：故障预报：为了能更早给出即将到来故障的警告，在当前时刻k，潜在故障指标F_k使用了预测数据以使该指标具有超前性。健康指标H_k和潜在故障指标F_k一起合作进行故障预测，即在时刻k，其中k>m+2，如果F_k>H_k且F_k-1>H_k-1且F_k-2>H_k-2，给出即将来临故障的预报，否则不预报；k＝k+1，转至步骤S2。

于本实施例中，所述步骤S1和步骤S3中在所述状态空间模型的基础上采用粒子滤波算法的具体计算方法如下：

所述系统的故障演化过程采用如下的状态空间模型描述：

x_k＝f(x_k-1,O_k-1,θ_k)+u_k-1 (1)

y_k＝h(x_k,θ_k)+v_k (2)

公式(1)为系统的状态方程，公式(2)为系统的观测方程；式中，x_k和y_k分别为系统的状态变量和观测变量，O_k为系统的输入变量，θ_k为所述状态空间模型的参数，u_k和v_k分别为系统噪声和观测噪声；

在所述状态空间模型中设计状态变量x与系统的故障演化过程相关，当可获得系统从起始时刻t＝0至当前时刻t＝k的观测序列{y₀,y₁,...,y_k}，即所述观测数据D₀＝{y_t}时，由粒子滤波算法实现对模型中未知的状态变量x的求解，即求出跟踪序列以及预测值和其中j>0；

令公式(1)简记为x_k＝f(x_k-1)，设u_k和v_k为高斯噪声，用于求解涉及故障演化的状态空间模型的粒子滤波算法的简要步骤如下：

1)在初始时刻k＝0，抽样生成服从先验分布p(x₀)的N_s个粒子(样本)其中i＝1,2,…,N_s，每个粒子对应的权值设为1/N_s，抽样生成对应的高斯系统噪声令k＝k+1；

2)在时刻k，抽样生成服从参考分布的N_s个粒子(样本)其中i＝1,2,…,N_s；例如，选定参考分布则由产生N个粒子；抽样生成对应的高斯系统噪声

3)给每个粒子分配对应的权值当选定参考分布为时，有

其中

4)对每个权值进行归一化：

5)估计时刻k的状态：

6)q步前向预测：

7)令k＝k+1，转步骤2)。

以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本发明的涵盖范围。

Claims (2)

1.一种基于两视角和粒子滤波的非监督故障预报方法，其特征在于包括以下步骤：

步骤S1：系统初始化：采用状态空间模型描述由状态变量x和观测变量y构成的系统，每个时刻t测量并获得观测数据D₀＝{y_t}，其中t＝1,2,3,…,m；m为正整数；移除D₀＝{y_t}中最大的5％数据及最小的5％数据后计算其均值μ₀和方差σ₀，在所述状态空间模型的基础上采用粒子滤波算法估计所述系统的状态并获得和其中t＝1,2,3,…,m；设定步长为q，并令k＝m+1；

步骤S2：健康指标的计算：视健康系统的观测数据近似于高斯分布p(y)＝N(μ₀,σ₀ ²)，从时间点k-q+1到k，将观测数据代入分布函数p(y)并求p(y_t)的均值：

$E_k({\mathrm{\μ}}_0,{\mathrm{\σ}}_0)=\frac{1}{q\sqrt{2\mathrm{\π}}{\mathrm{\σ}}_0}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{{\mathrm{\σ}}_0}^2}\}$

在时刻k，计算所述系统的健康指标：

$H_k=\frac{E_k({\mathrm{\μ}}_0,2.55{\mathrm{\σ}}_0)}{p\left({\mathrm{\μ}}_0\right)}=\frac{1}{q}\underset{i=k-q+1}{\overset{k}{\Σ}}\exp \{-\frac{{(y_i-{\mathrm{\μ}}_0)}^2}{2{\left(2.55{\mathrm{\σ}}_0\right)}^2}\}$

步骤S3：潜在故障指标的计算：在时刻k，在所述状态空间模型的基础上采用粒子滤波算法计算并计算预测数据和其中t＝k+1,k+2,…,k+q；计算预测数据对健康范围的偏离程度，偏离程度的计算公式为：

设Z表示第一次发生故障时的偏离，β为指数函数的均值参数，则在某个偏离z下的无故障概率为：

P(Z＞z)＝exp{-z/β}

在时刻k，当用已获得的预测数据代入偏离公式时，计算所述系统的潜在故障指标：

$F_k=1-\frac{1}{q}\underset{j=k+1}{\overset{k+q}{\Σ}}\exp \{-z\left({\hat{y}}_j\right)/\mathrm{\β}\}$

步骤S4：故障预报：在时刻k，其中k>m+2，如果F_k>H_k且F_k-1>H_k-1且F_k-2>H_k-2，给出即将来临故障的预报，否则不预报；k＝k+1，转至步骤S2。

2.根据权利要求1所述的基于两视角和粒子滤波的非监督故障预报方法，其特征在于：所述步骤S3中的健康范围为M₀＝[μ₀-3σ₀,μ₀+3σ₀]。