CN102779232A - 基于样条插值改进多变量灰色模型的故障预测方法 - Google Patents

Abstract

基于样条插值改进多变量灰色模型的故障预测方法，属于故障预测领域，本发明为解决现有多变量灰色模型在计算灰色系数和卷积积分过程中，采用代数精度较低的梯形公式或复化梯形公式，故障预测精度低的问题。本发明在不改变原有多变量灰色模型故障预测框架的前提下，当计算灰色系数和卷积积分时，分别用特定的样条函数代替原方法来逼近相应变量，推导一套全新的多变量灰色模型计算方法，有利于提高故障预测精度。

Description

基于样条插值改进多变量灰色模型的故障预测方法

技术领域

本发明涉及基于样条插值改进多变量灰色模型的故障预测方法，属于故障预测领域。

背景技术

故障预测技术能对设备状态进行估计并为决策提供服务与参考，是视情维修的关键技术之一。但实际设备往往由于其结构的复杂性和严重的不确定性，致使故障预测过程较为困难。灰色系统理论是1982年华中科技大学邓聚龙教授提出的研究少数据、不确定性的方法，能通过灰生成(包括累加生成，逆累加生成和均值生成)，降低原始数据序列的随机性，是故障预测的有效方法。

作为灰色系统理论的一个重要组成部分，灰色模型通过建立特殊的微分方程来揭示系统内部事物连续发展变化的规律，得到对原始数据较好的预测结果。事实上，灰色模型可以统称为GM(P，Q)，其中G代表Grey，M代表Model，P代表灰色模型是P阶方程，Q代表灰色模型有Q个变量。最常用的灰色模型是一维情形的GM(1，1)(P＝1，Q＝1)，而本发明要探讨的多变量灰色模型可以表示为GM(1，N)(P＝1，Q＝N)。与其他故障预测方法(如，回归分析预测，最小方差预测，马尔科夫模型预测，指数平滑预测，统计趋势预测等方法)相比，基于灰色模型的故障预测方法有以下优点：①允许少数据预测，最少只需用4个(对GM(1，1))或N+3个(对GM(1，N))数据即可进行预测；②计算过程简单，快速便捷；③无需知道太多的关联因素，数据比较容易得到；④尤其适合进行短期预测，也可以用于中长期预测。灰色模型的这些优良特性，促使其在很多领域，特别是数据贫乏和不确定性显著的故障预测情形中，得到了广泛应用。

然而，目前对灰色模型的应用大都局限于一维情形的GM(1，1)，而对二维乃至多维情形的GM(1，N)，应用却很少。究其原因，主要是由于GM(1，N)的理论还很不完善，模型模拟能力和故障预测能力不强。近年来，很多学者致力于GM(1，N)的研究，其中最为突出的是台湾国立虎尾科技大学的Tzu Li Tien教授研究的一种基于卷积积分的多变量模型(Convolution Integral based GM(1，N)，以下简称GMC(1，N))。在故障预测过程中，GMC(1，N)不再把原始GM(1，N)(Original GM(1，N)，以下简称OGM(1，N))中难以求积分的变量粗糙地设定为常量，而是把这个变量定义为卷积积分，试图通过对该卷积积分进行准确计算，取得较高的预测精度。此外，对灰色系数的估计也是非常重要的，它关系到模型模拟是否精确。遗憾的是，在OGM(1，N)中，计算灰色系数时采用了代数精度较低的梯形公式，而GMC(1，N)分别采用代数精度较低的梯形公式和复化梯形公式来计算灰色系数和卷积积分，不利于故障预测精度的提高。

发明内容

本发明目的是为了解决现有多变量灰色模型在计算灰色系数和卷积积分过程中，采用代数精度较低的梯形公式或复化梯形公式，故障预测精度低的问题，提供了一种基于样条插值改进多变量灰色模型的故障预测方法。

本发明所述基于样条插值改进多变量灰色模型的故障预测方法，该方法包括以下步骤：

步骤一、输入的待测设备的故障原始特征数据序列为

$X_1^{\left(0\right)}=\{x_1^{\left(0\right)}\left(1\right),x_1^{\left(0\right)}\left(2\right),...,x_1^{\left(0\right)}\left(k\right),...,x_1^{\left(0\right)}\left(n\right)\},$ k＝1，2，...，n    (1)

且 $x_1^{\left(0\right)}\left(k\right)\≥0;$

输入的待测设备的故障相关数据序列为

$X_i^{\left(0\right)}=\{x_i^{\left(0\right)}\left(1\right),x_i^{\left(0\right)}\left(2\right),...,x_i^{\left(0\right)}\left(t\right),...,x_i^{\left(0\right)}(n+m)\},$ i＝2，3，...，N，t＝1，2，...，n+m    (2)

且 $x_i^{\left(0\right)}\left(t\right)\≥0,$

其中，m为待预测数据个数，

由

i＝2，3，...，N决定，

通过公式 $x_1^{\left(1\right)}\left(k\right)=\underset{l=1}{\overset{k}{\mathrm{\Σ}}}{x}_1^{\left(0\right)}\left(l\right)---\left(3\right)$

和

公式 $x_i^{\left(1\right)}\left(t\right)=\underset{l=1}{\overset{t}{\mathrm{\Σ}}}{x}_i^{\left(0\right)}\left(l\right)---\left(4\right)$

计算i＝1，2，...，N的一阶累加生成数

步骤二、以

为初值，建立如下多变量灰色模型SF-GMC(1，N)

$\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}+{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)=\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N---\left(5\right)$

令灰色系数为κ＝[a b₁ b₂...b_N]^T，计算R-GMC(1，N)的灰色系数κ，其中

Y＝[u(1)u(2)…u(k)…u(n)]^T    (6)

$u\left(k\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(k\right)}{\mathrm{dk}}---\left(7\right)$

步骤二中κ的推导过程为：

将公式(5)中的

用三次样条或B样条表示，则可以近似地表示为多项式形式，设为 $x_1^{\left(1\right)}\left(t\right)\≈f{x}_1^{\left(1\right)}\left(t\right),$

直接对

求导，则有 $\frac{d{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}\≈\frac{\mathrm{df}{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$

令 $\frac{{\mathrm{dfx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}=u\left(t\right),$ 则 $u\left(t\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$ t＝1，2，，...，n，

将公式(7)代入公式(5)，则有

$-{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)+\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N=u\left(t\right),$ t＝1，2，....，n+m    (13)

公式(13)化简为Bκ＝Y，

其中κ＝[λ μ₁ μ₂...μ_N]^T，Y和B分别由公式(6)和公式(8)求得，

则由最小二乘法，可以得出κ＝(B^TB)^-1B^TY。

步骤三、利用步骤二求得的灰色系数κ并按公式(9)对

进行预测：

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_N}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\tilde{I}\left(t\right)---\left(9\right)$

其中，

为的预测值，为与步骤二选定的样条函数有关的卷积积分，步骤四、由步骤三的

可求得

${\hat{x}}_1^{\left(0\right)}\left(t\right)={\hat{x}}_1^{\left(1\right)}\left(t\right)-{\hat{x}}_1^{\left(1\right)}(t-1),$ t＝2，3，...，n+m    (10)

其中， ${\hat{x}}_1^{\left(0\right)}\left(1\right)=x_1^{\left(0\right)}\left(1\right),$ 为

的预测值；

步骤五、由步骤一至四得到了

t＝1，2，...，n+m，就可以通过下式计算SF-GMC(1，N)的预测精度

$\mathrm{RMSPEPR}=\sqrt{\frac{1}{n}\underset{t=1}{\overset{n}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(11\right)$

$\mathrm{RMSPEPO}=\sqrt{\frac{1}{m}\underset{t=n+1}{\overset{n+m}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(12\right)$

其中，RMSPEPR为由前部分样本，即t＝1，2，...，n组成的根平均平方误差百分比，而RMSPEPO为由后部分样本，即t＝n+1，n+2，...，n+m组成的根平均平方误差百分比；

步骤六、判断条件

是否满足，其中ε₁为故障阈值上界，ε₂为故障阈值下界，

若满足，则说明待测设备无故障；否则，判定待测设备有故障。

本发明的优点：

1)在计算灰色系数的过程中，本发明方法采用特定的样条函数来代替梯形公式，对微分方程进行有效的离散化，最终得到准确的灰色系数；

2)在计算卷积积分的过程中，本发明方法采用与1)相同的样条函数来代替复化梯形公式，最终得到精确的故障预测值；

3)本发明不改变原有多变量灰色模型故障预测的框架，计算过程简单，易于实现，有利于推广应用。

附图说明

图1是本发明所述一种基于样条插值改进多变量灰色模型的故障预测流程图；

图2是本发明方法与GMC(1，N)求得的

t＝1，2，...，10的故障预测值及原始值对比图；

图3是本发明方法与GMC(1，N)求得的

t＝1，2，...，10的故障预测值及原始值对比图；

图4是本发明方法与GMC(1，N)的故障预测精度(包括RMSPEPR和RMSPEPO)对比图。

具体实施方式

具体实施方式一：下面结合图1说明本实施方式，本实施方式所述基于样条插值改进多变量灰色模型的故障预测方法，该方法包括以下步骤：

步骤一、输入的待测设备的故障原始特征数据序列为

$X_1^{\left(0\right)}=\{x_1^{\left(0\right)}\left(1\right),x_1^{\left(0\right)}\left(2\right),...,x_1^{\left(0\right)}\left(k\right),...,x_1^{\left(0\right)}\left(n\right)\},$ k＝1，2，...，n    (1)

且 $x_1^{\left(0\right)}\left(k\right)\≥0;$

输入的待测设备的故障相关数据序列为

$X_i^{\left(0\right)}=\{x_i^{\left(0\right)}\left(1\right),x_i^{\left(0\right)}\left(2\right),...,x_i^{\left(0\right)}\left(t\right),...,x_i^{\left(0\right)}(n+m)\},$ i＝2，3，...，N，t＝1，2，...，n+m    (2)

且 $x_i^{\left(0\right)}\left(t\right)\≥0,$

其中，m为待预测数据个数，

由

i＝2，3，...，N决定，

通过公式 $x_1^{\left(1\right)}\left(k\right)=\underset{l=1}{\overset{k}{\mathrm{\Σ}}}{x}_1^{\left(0\right)}\left(l\right)---\left(3\right)$

和

公式 $x_i^{\left(1\right)}\left(t\right)=\underset{l=1}{\overset{t}{\mathrm{\Σ}}}{x}_i^{\left(0\right)}\left(l\right)---\left(4\right)$

计算i＝1，2，...，N的一阶累加生成数

步骤二、以为初值，建立如下多变量灰色模型SF-GMC(1，N)

$\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}+{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)=\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N---\left(5\right)$

令灰色系数为κ＝[a b₁ b₂...b_N]^T，计算R-GMC(1，N)的灰色系数κ，其中

Y＝[u(1)u(2)…u(k)…u(n)]^T    (6)

$u\left(k\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(k\right)}{\mathrm{dk}}---\left(7\right)$

步骤三、利用步骤二求得的灰色系数κ并按公式(9)对进行预测：

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_N}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\tilde{I}\left(t\right)---\left(9\right)$

其中，

为

的预测值，

为与步骤二选定的样条函数有关的卷积积分，步骤四、由步骤三的

可求得

${\hat{x}}_1^{\left(0\right)}\left(t\right)={\hat{x}}_1^{\left(1\right)}\left(t\right)-{\hat{x}}_1^{\left(1\right)}(t-1),$ t＝2，3，...，n+m    (10)

其中， ${\hat{x}}_1^{\left(0\right)}\left(1\right)=x_1^{\left(0\right)}\left(1\right),$

为

的预测值；

步骤五、由步骤一至四得到了

t＝1，2，...，n+m，就可以通过下式计算SF-GMC(1，N)的预测精度

$\mathrm{RMSPEPR}=\sqrt{\frac{1}{n}\underset{t=1}{\overset{n}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(11\right)$

$\mathrm{RMSPEPO}=\sqrt{\frac{1}{m}\underset{t=n+1}{\overset{n+m}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(12\right)$

其中，RMSPEPR为由前部分样本，即t＝1，2，...，n组成的根平均平方误差百分比，而RMSPEPO为由后部分样本，即t＝n+1，n+2，...，n+m组成的根平均平方误差百分比；

步骤六、判断条件

是否满足，其中ε₁为故障阈值上界，ε₂为故障阈值下界，

若满足，则说明待测设备无故障；否则，判定待测设备有故障。

具体实施方式二：本实施方式对实施方式一作进一步说明，步骤二中κ的推导过程为：

将公式(5)中的用三次样条或B样条表示，则

可以近似地表示为多项式形式，设为 $x_1^{\left(1\right)}\left(t\right)\≈f{x}_1^{\left(1\right)}\left(t\right),$

直接对

求导，则有 $\frac{d{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}\≈\frac{\mathrm{df}{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$

令 $\frac{{\mathrm{dfx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}=u\left(t\right),$ 则 $u\left(t\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$ t＝1，2，...，n，

将公式(7)代入公式(5)，则有

$-{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)+\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N=u\left(t\right),$ t＝1，2，...，n+m    (13)

公式(13)化简为Bκ＝Y，

其中κ＝[λ μ₁ μ₂...μ_N]^T，Y和B分别由公式(6)和公式(8)求得，

则由最小二乘法，可以得出κ＝(B^TB)^-1B^TY。

具体实施方式三：本实施方式对实施方式一作进一步说明，步骤三中公式(9)的推导过程为：

以

为初值，求解公式(5)所表示的SF-GMC(1，N)，得

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\left(\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right){\mathrm{\μ}}_N\right)\mathrm{d\τ}$

$=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_N}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ}---\left(14\right)$

其中，

为

的预测值。

与选定的样条函数有关的卷积积分 $\tilde{I}\left(t\right)=\underset{1}{\overset{t}{\∫}}\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}=\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ},$ 其中 $\tilde{g}\left(\mathrm{\τ}\right)=e^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right),$ 则公式(14)等价于公式(9)，

用选定的样条函数来逼近

设为

那么

可表示为 $\tilde{I}\left(t\right)=\underset{1}{\overset{t}{\∫}}\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}\≈\underset{1}{\overset{t}{\∫}}f\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}.$

具体实施方式四：下面结合图1至图4给出一个具体的实施例，建立基于样条插值的多变量灰色模型(Spline Function based GMC(1，N)，以下简称SF-GMC(1，N)，本实施例中取S-GMC(1，2)，N＝2)的实施例，来阐述本发明的具体实施方式。

执行步骤一、假设输入的故障原始特征数据序列和相关数据序列分别为 $X_1^{\left(0\right)}=\{x_1^{\left(0\right)}\left(1\right),x_1^{\left(0\right)}\left(2\right),...,x_1^{\left(0\right)}\left(5\right)\}$ 及 $X_2^{\left(0\right)}=\{x_2^{\left(0\right)}\left(1\right),x_2^{\left(0\right)}\left(2\right),...,x_2^{\left(0\right)}\left(5\right),...,x_2^{\left(0\right)}\left(10\right)\},$ 具体数值如表1所示，其中

可由

反映，待预测数据个数为5，

k＝1，2，...，5且t＝1，2，...，10， $X_1^{\left(0\right)}=\{x_1^{\left(0\right)}\left(6\right),x_1^{\left(0\right)}\left(7\right),...,x_1^{\left(0\right)}\left(10\right)\}$ 为

的原始值，用于计算预测精度。通过

k＝1，2，...，5和

k＝1，2，...，10，计算i＝1，2的一阶累加生成数(the First-order Accumulated Generating Operation，以下简称1-AGO)。

执行步骤二、以为初值，建立如下SF-GMC(1，2)

$\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}+{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)=\underset{i=1}{\overset{2-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_2,$ t＝1，2，...，10    (15)

由κ＝[λ μ₁ μ₂]^T＝(B^TB)^-1B^TY，计算SF-GMC(1，2)的灰色系数，其中

Y＝[u(1)u(2)…u(5)]^T    (16)

$u\left(t\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},t=\mathrm{1,2},...,5---\left(17\right)$

$B=\left[\begin{array}{ccc}-x_1^{\left(1\right)}\left(1\right)& x_2^{\left(1\right)}\left(1\right)& 1\\ -x_1^{\left(1\right)}\left(2\right)& x_2^{\left(1\right)}\left(2\right)& 1\\ \·& \·& \·\\ \·& \·& \·\\ \·& \·& \·\\ -x_1^{\left(1\right)}\left(5\right)& x_2^{\left(1\right)}\left(5\right)& 1\end{array}\right]---\left(18\right)$

事实上，以上步骤二中κ的推导过程如下：

选定某种特定的样条(如，三次样条，B样条等，本实施例以三次样条为例)，然后把(11)中的用该样条表示，则

可以近似地表示为如下多项式形式，设为

则

$x_1^{\left(1\right)}\left(t\right)\≈{\mathrm{fx}}_1^{\left(1\right)}\left(t\right)=\underset{j=1}{\overset{5}{\mathrm{\Σ}}}[x_1^{\left(1\right)}\left(j\right){\mathrm{\α}}_j\left(t\right)+{\mathrm{fx}}_1^{\′\left(1\right)}\left(j\right){\mathrm{\β}}_j\left(t\right)]---\left(19\right)$

其中，

为

j＝1，2，...，5的一阶导数，

设

则

和

可表示为

${\mathrm{fx}}_1^{\left(1\right)}\left(t\right)=M\left(j\right)\frac{{[(j+1)-t]}^3}{6}+M(j+1)\frac{{(t-j)}^3}{6}+(x_1^{\left(1\right)}\left(j\right)-\frac{M\left(j\right)}{6})[(j+1)-t]$ (22)

$+(x_1^{\left(1\right)}(j+1)-\frac{M(j+1)}{6})(t-j),$ 如果t∈[j，j+1]

$\frac{{\mathrm{dfx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}=-3M\left(j\right)\frac{{[(j+1)-t]}^2}{6}+3M(j+1)\frac{{(t-j)}^2}{6}-(x_1^{\left(1\right)}\left(j\right)-\frac{M\left(j\right)}{6})$ (23)

$+(x_1^{\left(1\right)}(j+1)-\frac{M(j+1)}{6}),$ 如果t∈[j，j+1]

依据三次样条理论，公式(23)中的M(t)可以由带三对角矩阵的方程求得。特别地，对于本实施例给定的数据，利用上述三次样条函数及最小二乘法，有κ＝[-0.5982 -0.7978 4.0622]。

执行步骤三、利用步骤二求得的灰色系数κ＝[λ μ₁ μ₂]^T＝[-0.5982 -0.7978 4.0622]，就可以由下式对t＝1，2，...，10进行预测

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_2}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\tilde{I}\left(t\right),$ t＝1，2，...，10    (24)

其中，t＝1，2，...，10为

t＝1，2，...，10的预测值，t＝1，2，...，10为与步骤二选定的样条函数(本实施例中为三次样条)有关的卷积积分。

事实上，步骤三中公式(24)的推导过程如下：

以

为初值，求解公式(15)所表示的SF-GMC(1，2)，得

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\left(\underset{i=1}{\overset{2-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right){\mathrm{\μ}}_2\right)\mathrm{d\τ}$

$=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_2}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\underset{i=1}{\overset{2-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ}---\left(25\right)$

其中，

t＝1，2，...，10为

t＝1，2，...，10的预测值。

为了求公式(25)中的

t＝2，3，...，10，最重要的是要计算卷积积分 $\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{2-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ},$ 定义 $\tilde{I}\left(t\right)=\underset{1}{\overset{t}{\∫}}\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}=\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{2-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ},$ 其中

则公式(25)等价于(24)。更进一步地，可以用步骤二选定的三次样条来逼近

设为

那么

可表示为 $\tilde{I}\left(t\right)=\underset{1}{\overset{t}{\∫}}\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}\≈\underset{1}{\overset{t}{\∫}}f\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}.$ 对于本实施例，有

其中，M(t)与步骤二的相同。将公式(26)代入公式(24)，则能得到

t＝1，2，...，10的预测值

t＝1，2，...，10。图2所示为本发明方法与基于卷积积分的多变量模型(Convolution Integral based GM(1，N)，以下简称GMC(1，N)，本实施例中取GMC(1，2)，N＝2)求得的

t＝1，2，...，10的预测值及原始值对比图。可以看出，本发明得到的t＝1，2，...，10与原始值较为贴近，预测结果比较理想。

执行步骤四、由步骤三的

t＝1，2，...，10，可求得

${\hat{x}}_1^{\left(t\right)}\left(t\right)={\hat{x}}_1^{\left(1\right)}\left(t\right)-{\hat{x}}_1^{\left(1\right)}(t-1),$ t＝2，3，...，10    (27)

其中，

t＝1，2，...，10为

t＝1，2，...，10的预测值。图3显示了本发明方法SF-GMC(1，2)与GMC(1，2)求得的

t＝1，2，...，10的预测值及原始值对比图。与图2类似，可以发现SF-GMC(1，2)得到的

t＝1，2，...，10与原始值较为贴近，预测结果较为理想。

执行步骤五、由步骤一至四得到了

t＝1，2，...，10，就可以通过下式计算SF-GMC(1，2)的预测精度

$\mathrm{RMSPEPR}=\sqrt{\frac{1}{5}\underset{t=1}{\overset{5}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(28\right)$

$\mathrm{RMSPEPO}=\sqrt{\frac{1}{5}\underset{t=6}{\overset{10}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(29\right)$

其中，RMSPEPR为由前部分样本(即t＝1，2，...，5)组成的根平均平方误差百分比(RootMean Squared Percentage Error of PRiori-sample period)，而RMSPEPO为由后部分样本(即t＝6，7，...，10)组成的根平均平方误差百分比(Root Mean Squared Percentage Error ofPOst-sample period)。图4所示为SF-GMC(1，2)与GMC(1，2)的预测精度(包括RMSPEPR和RMSPEPO)对比图。从图4可以看出，本发明方法得到的RMSPEPR和RMSPEPO比GMC(1，2)分别降低了3.63％和18.72％，预测精度大大提高。

步骤六、判断条件

t＝6，7，...，10是否满足，其中ε₁＝0，ε₂＝45分别为故障上界和下界。从图3知，GMC(1，2)在序列6之后开始严重偏离原始值，并在序列9至序列10判断出故障，产生了虚警，而S-GMC(1，2)与原始值吻合较好，能准确地对设备故障进行预测。

总之，图2至图4充分说明了本发明SF-GMC(1，N)方法能有效提高故障预测精度，准确反映原始数据变化趋势，是对原有方法的巧妙改进。

表1实施例用到的数据序列

表1中给出的本实施方式提供的实施例涉及的相关数据。

Claims (3)

1.基于样条插值改进多变量灰色模型的故障预测方法，其特征在于，该方法包括以下步骤：

步骤一、输入的待测设备的故障原始特征数据序列为

$X_1^{\left(0\right)}=\{x_1^{\left(0\right)}\left(1\right),x_1^{\left(0\right)}\left(2\right),...,x_1^{\left(0\right)}\left(k\right),...,x_1^{\left(0\right)}\left(n\right)\},$ k＝1，2，...，n    (1)

且 $x_1^{\left(0\right)}\left(k\right)\≥0;$

输入的待测设备的故障相关数据序列为

$X_i^{\left(0\right)}=\{x_i^{\left(0\right)}\left(1\right),x_i^{\left(0\right)}\left(2\right),...,x_i^{\left(0\right)}\left(t\right),...,x_i^{\left(0\right)}(n+m)\},$ i＝2，3，...，N，t＝1，2，...，n+m    (2)

且 $x_i^{\left(0\right)}\left(t\right)\≥0,$

其中，m为待预测数据个数，

由

i＝2，3，...，N决定，

通过公式 $x_1^{\left(1\right)}\left(k\right)=\underset{l=1}{\overset{k}{\mathrm{\Σ}}}{x}_1^{\left(0\right)}\left(l\right)---\left(3\right)$

和

公式 $x_i^{\left(1\right)}\left(t\right)=\underset{l=1}{\overset{t}{\mathrm{\Σ}}}{x}_i^{\left(0\right)}\left(l\right)---\left(4\right)$

计算

i＝1，2，...，N的一阶累加生成数

步骤二、以

为初值，建立如下多变量灰色模型SF-GMC(1，N)

$\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}+{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)=\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N---\left(5\right)$

令灰色系数为κ＝[a b₁ b₂...b_N]^T，计算R-GMC(1，N)的灰色系数κ，其中

Y＝[u(1)u(2)…u(k)…u(n)]^T    (6)

$u\left(k\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(k\right)}{\mathrm{dk}}---\left(7\right)$

步骤三、利用步骤二求得的灰色系数κ并按公式(9)对

进行预测：

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_N}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\tilde{I}\left(t\right)---\left(9\right)$

其中，

为

的预测值，

为与步骤二选定的样条函数有关的卷积积分，步骤四、由步骤三的

可求得

${\hat{x}}_1^{\left(0\right)}\left(t\right)={\hat{x}}_1^{\left(1\right)}\left(t\right)-{\hat{x}}_1^{\left(1\right)}(t-1),$ t＝2，3，...，n+m    (10)

其中， ${\hat{x}}_1^{\left(0\right)}\left(1\right)=x_1^{\left(0\right)}\left(1\right),$

为

的预测值；

步骤五、由步骤一至四得到了

t＝1，2，...，n+m，就可以通过下式计算SF-GMC(1，N)的预测精度

$\mathrm{RMSPEPR}=\sqrt{\frac{1}{n}\underset{t=1}{\overset{n}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(11\right)$

$\mathrm{RMSPEPO}=\sqrt{\frac{1}{m}\underset{t=n+1}{\overset{n+m}{\mathrm{\Σ}}}{[{\hat{x}}_1^{\left(0\right)}\left(t\right)-x_1^{\left(0\right)}\left(t\right)]}^2/{\left[x_1^{\left(0\right)}\left(t\right)\right]}^2}\×100\%---\left(12\right)$

其中，RMSPEPR为由前部分样本，即t＝1，2，...，n组成的根平均平方误差百分比，而RMSPEPO为由后部分样本，即t＝n+1，n+2，...，n+m组成的根平均平方误差百分比；

步骤六、判断条件

是否满足，其中ε₁为故障阈值上界，ε₂为故障阈值下界，

若满足，则说明待测设备无故障；否则，判定待测设备有故障。

2.根据权利要求1所述基于样条插值改进多变量灰色模型的故障预测方法，其特征在于，步骤二中κ的推导过程为：

将公式(5)中的

用三次样条或B样条表示，则

可以近似地表示为多项式形式，设为 $x_1^{\left(1\right)}\left(t\right)\≈f{x}_1^{\left(1\right)}\left(t\right),$

直接对

求导，则有 $\frac{d{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}\≈\frac{\mathrm{df}{x}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$

令 $\frac{{\mathrm{dfx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}}=u\left(t\right),$ 则 $u\left(t\right)\≈\frac{{\mathrm{dx}}_1^{\left(1\right)}\left(t\right)}{\mathrm{dt}},$ t＝1，2，...，n，

将公式(7)代入公式(5)，则有

$-{\mathrm{\λx}}_1^{\left(1\right)}\left(t\right)+\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\η}}_i{x}_{i+1}^{\left(1\right)}\left(t\right)+{\mathrm{\μ}}_N=u\left(t\right),$ t＝1，2，...，n+m    (13)

公式(13)化简为Bκ＝Y，

其中κ＝[λ μ₁ μ₂...μ_N]^T，Y和B分别由公式(6)和公式(8)求得，

则由最小二乘法，可以得出κ＝(B^TB)^-1B^TY。

3.根据权利要求1所述基于样条插值改进多变量灰色模型的故障预测方法，其特征在于，步骤三中公式(9)的推导过程为：

以

为初值，求解公式(5)所表示的SF-GMC(1，N)，得

${\hat{x}}_1^{\left(1\right)}\left(t\right)=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\left(\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right){\mathrm{\μ}}_N\right)\mathrm{d\τ}$

$=x_1^{\left(1\right)}\left(1\right)e^{-\mathrm{\λ}(t-1)}+\frac{{\mathrm{\μ}}_N}{\mathrm{\λ}}(1-e^{-\mathrm{\λ}(t-1)})+\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\λ}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\μ}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ}---\left(14\right)$

其中，

为

的预测值。

与选定的样条函数有关的卷积积分 $\tilde{I}\left(t\right)=\underset{1}{\overset{t}{\∫}}\tilde{g}\left(\mathrm{\τ}\right)\mathrm{d\τ}=\underset{1}{\overset{t}{\∫}}{e}^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right)\mathrm{d\τ},$ 其中 $\tilde{g}\left(\mathrm{\τ}\right)=e^{-\mathrm{\α}(t-\mathrm{\τ})}\underset{i=1}{\overset{N-1}{\mathrm{\Σ}}}{\mathrm{\β}}_i{x}_{i+1}^{\left(1\right)}\left(\mathrm{\τ}\right),$ 则公式(14)等价于公式(9)。

Images