一种非线性过程故障诊断方法

专利名称：一种非线性过程故障诊断方法
技术领域：
本发明属于故障检测与诊断技术领域，特别涉及一种基于改进的核独立元分析和支持向量机相结合的非线性过程诊断方法。

背景技术：
随着计算机和电子测量技术的飞速发展，现代工业过程大都具有完备的传感测量装置，可以在线获得大量的过程数据。对这些数据的统计分析可以帮助操作人员及时发现过程故障，避免重大事故的发生，这就促进了人们研究基于数据分析的过程监测方法。目前该类方法中使用最多的是主元分析(PCA)、偏最小二乘(PLS)及独立元分析(ICA)等。然而这些方法都是基于线性变换的统计方法，即假定过程变量之间满足线性相关关系。但是对于复杂的工业过程，像化工业的田纳西-伊斯曼过程、废水处理过程以及矿业中的连续采煤机采煤过程等等，变量之间往往呈现出非线性的相关关系特性，利用传统的统计方法将给过程监测带来较大的误差，增大过程故障误报、漏报的概率。
针对这种情况，Kramer在《AIChE Journal》期刊中提出了一种基于自联想的五层(输入，映射，瓶颈，解映射，及输出层)神经网络的非线性PCA。Dong和McAvoy在《Computers and ChemicalEngineering》期刊中提出了一种基于主曲线和神经网络的非线性PCA，并把它应用到非线性过程监测。基于遗传编程(genetic programming)以及输入训练神经网络的非线性PCA方法也已经分别在《Computers and Chemical Engineering》和《AIChE Journal》期刊中被提出。但是大部分现有的非线性PCA方法都是基于神经网络，必须离线训练，要计算主元的话也必须解决非线性最优化问题，并且在训练神经网络之前必须事先明确主元的个数。而且这些方法使用的非线性变换函数一般难以获得，神经网络的训练也较为困难。目前，核理论在非线性工业过程中的应用逐渐增多。对于非线性过程监测和故障诊断，引入了核主元分析(KPCA)和核独立元分析(KICA)。KICA首次是由Lee和Qin等人提出的。它是一种KPCA和ICA结合的新的非线性扩展。KPCA在一个与输入空间非线性相关的高维特征空间中计算主元。ICA提取高阶统计量以及分解观测数据。因此，独立元(ICs)比主元(PCs)从观测数据中展露了更多动态信息。同其它的非线性方法相比，KICA结合了KPCA和ICA的优点，成为一种在线监测故障的非线性动态方法。目前，把KICA方法用于工业过程存在的问题是形成的核矩阵的维数是采样次数的平方，这样对于过程非常不利。ICA方法是用于动态过程中提取独立元，但是由于输入空间和特征空间不能像PCA那样自如地转换，使得故障诊断变得困难。


发明内容
针对现有技术存在的问题，本发明提供一种非线性过程故障诊断方法。对核独立元分析(KICA)方法进行了改进，并同最小二乘支持向量机(LS-SVM)结合，提出了一种基于改进的KICA和LS-SVM的非线性过程故障诊断方法。
针对故障检测与诊断，本发明在KICA的基础上进行了如下改进 1)能够表达所有采样数据信息的维数适当的新的核矩阵的形成技术，提出了一个组群中数据的相似性概念。经过相似性分析之后，提取了特征空间中非线性信息来表达所有观测数据，这样当训练样本个数很大时，就降低了KPCA的计算量。
2)基于LS-SVM提出了一种故障分类策略。
核理论已经越来越多的应用于非线性工业过程，表明了核函数在工业过程对观测数据非线性化的优越性。建模之前必须进行数据相似性分析，原因如下 1)映射数据进入特征空间变得线性冗余； 2)假如数据是线性的，那么当使用核技巧非线性化时将会产生更大的误差； 3)在KPCA训练过程中，核矩阵的大小是样本个数的平方。当样本数变大，特征值和特征向量的计算将会十分耗费时间。
本发明通过分别探索输入空间和特征空间样本的线性相关性来分析观测数据。并引入了相似性概念来定义新数据和节点集之间的相似性因子。
改进的KICA和LS-SVM均利用mercer核，两者很容易进行结合。本发明用LS-SVM进行故障诊断，一旦检测到故障，来自改进的KICA的核变换得分将直接被引入，作为LS-SVM的输入对故障进行分类。
本发明的故障诊断方法包括如下步骤 步骤一、采集数据 采集过程中相关变量的数据，对于每个故障，产生两组数据，即训练数据和实时工况数据。训练数据用于建立模型，实时工况数据用于在线监测。并用均值和标准偏差规范化采集的数据。
步骤二、进行相似性分析 相似性分析包括输入空间的相似性分析和特征空间的相似性分析。本发明根据相似性因子对数据进行分析，去除掉相似性较强的数据。既解决了特征空间中数据的线性冗余，又降低了计算负载，还减小了使用核技巧带来的误差。
1、输入空间的相似性分析 在输入空间中，假如三个数据点x，y，z满足 那么x，y，z之间是线性关系。
学习的开始，只有两个数据点，即第i个节点集被表示为i＝1，...，t1，且满足ni＜t1，t1为训练数据的个数。根据方程(1)定义相似性因子为 式中 xi+1-训练数据中的第i+1个的数据点；

-节点集Ni中任意两个数据点； S1满足0≤S1≤1。当S1＝1时相似性最强。因此，将xi+1同节点集Ni中所有数据都两两进行相似性计算，假如所有其中γ0是满足的预定义小变量，新数据xi+1将被引入建立新扩展节点集Ni+1＝{Ni，xi+1}；否则，这个数据点将被拒绝，即Ni+1＝Ni。在输入空间相似性分析后，获得第一个数据子集Ω1。
2、特征空间中的相似性分析 得到子集Ω1后，将Ω1中的数据通过非线性函数Ф映射到特征空间，进行特征空间的相似性分析。令x，y及z是Ω1中输入空间的三个数据点。映射三个数据点到特征空间，他们变为Ф(x)，Ф(y)及Ф(z)。
假如 那么 θ＝0 也就是，Ф(x)，Ф(y)及Ф(z)在特征空间都是线性的。
学习的开始，只有两个数据，即特征空间第j个节点集被表示为j＝1，...，t2，且满足nj＜t2，t2为子集Ω1中数据的个数。根据方程(3)定义特征空间的相似性因子为 式中Ф(xj+1)-数据子集Ω1中第j+1个数据点xj+1映射后的数据；

-节点集Lj中的任意两个数据； 因而，S2满足0≤S2≤1。当S2＝1时相似性是最强的。
S2的可选择相似性因子能按下面过程来决定，方程(3)等价于方程(4) 其中 f(x，y，z)＝(Ф(x)，Ф(y))-(Ф(y)，Ф(z))-(Ф(x)，Ф(z))+(Ф(z)，Ф(z)) 方程(5)利用核技巧 k(a，b)＝(Ф(a)，Ф(b))(6) 其中a，b为输入空间的两个数据点，将方程(4)变为 简化方程(6)，本文采用径向基核函数有 方程(7)等价于方程(8) 令 h＝1-k(y，z)-k(x，z) 因此假如数据点x，y，z满足 那么θ＝0，也就是，数据点x，y，z在特征空间是线性的。
且假如 k(x，z)＝k(x，y)(12) k(y，z)＝4k(x，y)-3 那么满足方程(10)。
当核函数选取径向基核函数，那么当数据点x，y，z满足方程(12)时x，y，z在特征空间是线性的。因此，特征空间中相似性因子可以转化为在输入空间进行计算，定义如下 且 其中S21，S22满足0≤S21≤ 1，0≤S22 ≤ 1，当S21＝1且S22＝1时相似性是最强的。因此不需要对子集Ω1中的数据进行非线性映射就可在特征空间对数据进行相似性分析，重新定义第j个节点集假如且其中γ1是一个满足的预定义小变量，xj+1将被引入建立新扩展节点集Lj+1＝{Lj，xj+1}。否则这个数据点将被拒绝，即Lj+1＝Lj。在特征空间对Ω1中的数据进行相似性分析后，保留的数据点个数表示为N；且获得了第二个数据集Ω2。
步骤三、利用核主元分析(KPCA)对数据做白化处理，求解白化后的观测变量z 通过非线性映射将输入空间映射到一个特征空间接着在此特征空间对观测数据进行白化处理，得到白化后的观测变量z。通过白化处理后，减少了ICA需要估计的参数，从而简化计算过程。
子集Ω2中包含m个变量的观测数据xk∈Rm，k＝1，...，N，其中N是观测的个数。利用非线性映射ФRm→F，原始空间中的观测数据就扩展到高维特征空间F，Ф(xk)∈F。此特征空间的协方差结构是一个单位矩阵。因而，特征空间中协方差矩阵将为 其中Ф(xj)，j＝1，...，N，假定为零均值和单位方差。令Θ＝[Ф(x1)，...，Ф(xN)]，因而CF能表示为定义一个N×N维的Gram核矩阵K [K]ij＝Kij＝<Ф(xi)，Ф(xj)>＝k(xi，xj)(16) 有K＝ΘTΘ。核函数k(xi，xj)的应用可以在不执行非线性映射的情况下在F中计算内积。即可以避免执行非线性映射，并且在特征空间通过引入一个核函数k(x，y)＝<Ф(x)，Ф(y)>计算内积。常用的核函数有径向基核函数多项式核函数k(x，y)＝<x，y>r，Sigmoid(S形)核函数k(x，y)＝tanh(β0<x，y>+β1)。核函数的选择潜在的决定了映射Ф以及特征空间F。从核矩阵K可知，高维空间中Ф(xk)的中心化能按以下方法执行，即将Ф(xk)的中心化转化为对K的中心化处理。中心化核矩阵

可由下式获得 其中


特征值分解 式中λ为

的特征根，α为λ对应的特征向量。由式(18)能够获得同

的d个最大正特征值λ1≥λ2≥...≥λd相应的标准正交特征向量α1，α2，...，αd。因而CF的d个最大正特征值是相应的特征向量v1，v2，...，vd能够表示为 j＝1，...，d(19) 特征向量矩阵V＝[v1，v2，...，vd]可以简单的由下式表达 V＝ΘHΛ-1/2(20) 其中Λ＝diag(λ1，λ2，...，λd)和H＝[α1，α2，...，αd]分别为

的d个最大特征值的对角阵及对应的 特征向量矩阵。特征空间中的白化矩阵P和映射数据的白化变换如下 z＝PTФ(x)(21) 具体为， 式中

k＝[k(x1，x)，...，k(xN，x)]T，1t＝(1/I)[1，...，1]∈R1×N； x1，...，xN—Ω2中的数据； x-需要白化处理的数据； 步骤四、利用修正独立元分析(ICA)提取独立元 在KPCA变换空间提取独立元s。修正ICA进一步从白化后的观测变量提取一组独立元。并使独立元各变量之间相互统计独立。找出p(≤d)个独立元，即s＝{s1，s2，...，sp}，满足E(ssT)＝D＝diag{λ1，...，λp}以使s的元素变得尽可能彼此独立，利用 s＝CTz(23) 其中C∈Rd×p为得分转换矩阵和CTC＝D。
定义归一化的独立元sn为 sn＝D-1/2s＝D-1/2CTz＝CnTz(24) Cn为标准得分转换矩阵，很明显D-1/2CT＝CnT，CnTCn＝I，以及E(snsnT)＝I。因此作为正交的结果，把原始ICA寻找一个任意分离矩阵的问题简化为寻找更少参数估计的矩阵Cn的简单问题。因而，所提出算法的目标改为从z∈Rd寻找sn∈Rp和Cn，使sn的元素变得尽可能彼此独立且满足E(snsnT)＝I。设定z的第一个数据中p个元为sn的初始元。令CnT＝[Ip0]，其中Ip是p维的单位矩阵且0是p×(d-p)的零矩阵。根据修正ICA算法计算标准得分转换矩阵Cn (1)选择p，估计独立元的个数，设定计数器i←1； (2)获取初始向量ci； (3)令ci←E{zg(ciTz)}-E{g′(ciTz)}ci，其中g′是g的一阶导数，Hyvrinen给出了三种g的表达函数为g1(u)＝tanh(a1u)，g2(u)＝u exp(-a2u2/2)，g3(u)＝u3； (4)执行正交化这种正交化排除了现有方案所包含的信息； (5)归一化 (6)假如ci还没有收敛，返回步骤(3)； (7)假如ci收敛，输出向量ci。若i≤p，则设定i←i+1且返回步骤(2)；一旦找到Cn，从下式可得到独立元 s＝D1/2CnTz(25) 步骤五、利用T2和SPE统计量以及LS-SVM进行故障检测与诊断 采用T2和SPE统计量进行在线故障检测，当观测数据的统计量没有超出统计量规定的控制限，则属于正常数据，反之属于异常数据，表明出现故障。一旦故障发生，将实时工况数据的独立元作为LS-SVM的输入进行故障诊断，根据训练数据建立的决策函数，即可判断故障的类型。
从Lee和Qin(2006)那，知道T2统计量和SPE统计量定义如下 T2＝sTD-1s(26) z＝PTФ(x)(27) 其中

能够由下式得来 因为s不服从高斯分布，T2的控制限通过F分布来决定。
SPE的控制限由下面权重χ2分布来计算(Qin，2003) SPE～μχh2，μ＝b/2a，h＝2a2/b(29) 其中a和b分别是标准操作数中SPE的估计均值和方差。
要对故障进行分类，首先要建立各个故障的决策函数，具体如下 LS-SVM分类器在F中构造一个最大边界(margin)超平面。在输入空间，这与个非线性决策边界相对应。
其中w表示权重向量且b为偏置(bias)。两种情况都能合成为单一框架如下 yi[<w，Ф(si>+b]≥+1-ξi 其中ξi≥0，i＝1，2，...，N是正松弛变量，用于允许不等式集的错误分类。因为两个超平面之间的边界宽是2/‖w‖，目标函数被定义为‖w‖2和惩罚项的最小值 其中C∈R+是一整定超参数，它决定边界最大化和训练误差最小化之间的权衡；松弛变量ξi≥0，i＝1，2，...，N是用于放宽硬边界LS-SVM分类器的约束以便允许存在一些被错误分类的数据；在LV-SVM中，当引入核技巧后，上式(31)经过一定的运算后，上式最优化问题能重新描述为 式中If-IfT＝[1，...，1]∈R1×N；
Ω-Ωij＝yiyjФ(si)TФ(sj)＝yiyjk(si，sj)； α-拉格朗日因子； 根据方程(32)可以求出偏置b和拉格朗日因子α，根据这两个参数可以确定相应的决策函数即分类器如下 式中s-需要判断故障类型数据的独立元； si-训练数据的独立元，i＝1，...，N； 决策函数建立后，当检测到故障，输入实时工况数据的独立元，根据对应的f(s)的值判断故障的类型。选取其中一个故障的决策函数，当f(s)＝+1时，发生本决策函数对应的故障，反之f(s)＝-1时，没有发生本决策函数对应的故障，继续根据其他决策函数进行故障诊断。
本发明的具体操作过程如图11所示的流程图，其中定义训练数据的子集Ω2中的观测变量为xk∈Rm，白化后的观测变量为z，独立元为s。实时工况数据为xt∈Rm，白化后的观测变量为zt，独立元为st。上面步骤中的计算公式(22)，(25)，(26)以及(33)等属于通用公式，均可用于训练数据和实时工况数据的计算。
本发明——基于改进的核独立元分析和支持向量机相结合的非线性过程故障诊断方法的优越性表现为 提出了非线性动态过程故障诊断技术，结合Kernel、ICA和LS-SVM三者的优点，即，发挥Kernel对非线性的表达能力，同时发挥ICA对动态特性的把握能力和LS-SVM的分类能力。



图1三个数据点在特征空间的角度； 图2田纳西-伊斯曼过程示意图； 图中1，2，3，4，5，6，7，8，9，10，11，12，13分别表示为流1，流2，流3，流4，流5，流6，流7，流8，流9，流10，流11，流12，流13；FC流量控制；FI流量指示器；PI压力指示器；PHL压力控制；SC同步回旋加速器；TC温度控制；TI温度指示器；LI液位指示器；LC液位控制；XC成分控制；XA，XB，XD，XE分别是成分A分析，成分B分析，成分D分析，成分E分析； 图3 a)，b)分别是田纳西-伊斯曼过程中故障8基于KPCA的监测结果和两个主元(nPC)；图4特征空间中的成分分析； 图5 a)，b)分别是田纳西-伊斯曼过程中故障8基于KICA的监测结果和两个独立元(IC)； 图6 a)，b)分别是田纳西-伊斯曼过程中故障8基于改进的KICA的监测结果和两个IC； 图7废水处理过程WWTP示意图； 图中Qin流入量；Qa内部循环流量；Qe流出量；Qf进料层流量；Qr返回污泥的流量；Qu沉淀池底部流流量；Qw处理后的水的输入流量； 图8标准操作条件数据的变量模式； 图9 a)，b)分别是WWTP在两场暴雨情况下KICA的监测结果和KICA的两个IC； 图10 a)，b)分别是WWTP在两场暴雨情况下基于γ0＝0.05的改进的KICA的监测结果和两个IC； 图11本发明方法实现的流程图。
以美国ROCKWELL公司的可编程控制器(PLC)实现基础控制，监控程序用RSView32提供的VBA应用软件编制。监控软件在单独的计算机上运行，该计算机上装有RSLinx通讯软件，负责与PLC和上位机进行数据通讯，RSLinx与监控程序之间通过DDE方式进行双向通讯。把监控结果输出到计算机的系统管理画面，同时把监控结果保存到实时数据库中，为操作者或相关技术工人进行监督操作提供参考指导作用。

具体实施例方式 实例1、田纳西-伊斯曼(Tennessee Eastman)过程 本发明提出的方法应用到了田纳西-伊斯曼过程仿真数据中，并同KPCA及原KICA的监测结果进行对比。田纳西-伊斯曼过程是一个复杂非线性过程，它是由Eastman化学品公司所创建，其目的是为评价过程控制和监测方法提供一个真实的工业过程。控制结构如图2所示。过程包括五个主要单元反应器，冷凝器，压缩机，分离器，汽提塔；而且，它包含了八种成分A，B，C，D，E，F，G和H。四种反应物A，C，D以及E和惰性B一起被加进反应器里，形成产品G和H，还有副产品F。田纳西-伊斯曼过程包括21个预设定的故障，如表1所示。包括22个连续过程测量，12个控制变量，和19个成分测量。如表2所示。除了反应器的搅拌器的搅动速度(因为没有对它进行控制)，共52个观测变量用于本研究的监测。
表1.田纳西-伊斯曼过程的过程故障描述 表2.田纳西-伊斯曼过程中的监测变量 步骤一、采集数据 对于训练数据和实时工况数据采用了三分钟的采样间隔采集数据。每个故障的训练数据由480个观测数构成，实时工况数据由960个观测数构成。开始时都是没有故障的，训练数据中，故障都是在第20次采样时引入，实时工况数据中，故障都是在第160次采样引入的。训练数据和实时工况数据中的数据都包含了52个观测变量。本实例主要针对故障8进行建模分析，随机选取了故障8的训练数据和实时工况数据中控制变量的5组数据分别如表3和表4所示 表3.故障8训练数据中控制变量的五组数据

表4.故障8实时工况数据中控制变量的五组数据
步骤二、进行相似性分析 对训练数据的数据进行输入空间的相似性分析，取γ0＝0.95。学习的开始只有两个数据点，即节点集N1＝{x1，x2}。当引入数据点x3时，根据公式(1)进行计算，当时， x3加入节点集N2中，得到N2＝{x1，x2，x3}，否则N2＝N1。当引入数据点x4时，也根据公式(1)进行计算，当时，x4加入节点集N3中，得到N3＝{N2，x4}，否则N3＝N2。
依次类推，每引入一个新的数据点，都根据公式(1)进行计算，直到引入x480为止，得到第一个数据子集Ω1。
再对数据子集Ω1中的数据进行特征空间的相似性分析，取γ1＝0.95，核函数选径向基核函数。由于Ω1中数据的编号不再连续所以要对Ω1中的数据重新编号。学习的开始也只有两个数据点，即节点集L1＝{X1，x2}。当引入数据点x3时，根据公式(13)、(14)进行计算，当且时，x3加入节点集L2，得到L2＝{x1，x2，x3}，否则L2＝L1。当引入数据点x4时，根据公式(13)、(14)进行计算，当且时，x4加入节点集L3，得到L3＝{L2，x4)，否则L3＝L2。依次类推，直到Ω1中最后一个数据点进行相似性分析，得到第二个数据子集Ω2。
步骤三、利用KPCA对数据做白化处理 选取径向基核函数作为核函数，其中σ＝0.5。对第二个数据子集Ω2中的数据根据公式(16)进行核计算，得到Gram核矩阵K。接着根据公式(17)对Gram核矩阵K进行中心化处理，得到中心化的核矩阵

。通过

计算出了它的特征值，并选取了12个最大的正特征值为λ1≥λ2≥...≥λ12。根据公式(18)得到12个特征值对应的特征向量矩阵H。最后在特征空间进行白化，根据公式(22)依次对Ω2中的数据进行白化处理得到白化后的观测变量z。
根据公式(22)对实时工况数据进行白化处理得到白化后的观测变量z，。
步骤四、利用修正ICA提取独立元 取p＝8，根据ICA算法计算训练数据的标准得分转换矩阵Cn，其中g(u)＝tanh(1.5u)。接着利用上一步得到白化后的观测变量z，根据公式s＝D1/2CnTz计算得到训练数据的独立元s。再利用上一步得到的白化后的观测变量zt，根据公式st＝D1/2CnTzt计算得到实时工况数据的独立元st。
步骤五、故障检测与诊断 将训练数据的独立元

作为LS-SVM的输入，确定相应的决策函数。建立训练样本S{si，y}，i＝1，2，…，n，其中yi∈{+1，-1)，当si为故障8的数据时，yi＝+1，否则yi＝-1。核函数仍然选取径向基核函数。由此就建立了故障8相应的决策函数为。决策函数建立后，根据公式(26)、(27)计算实时工况数据的T2和SPE统计量。从统计量图中看出检测到故障，此时将实时工况数据的独立元st代入决策函数进行分类，根据f(st)的值判断实时工况数据中发生故障的故障类型。具体监测结果见对比实例。
对比实例 对所有的数据都在KPCA、原始KICA和改进的KICA应用之前进行规范化。对于统计量，每幅图中的点线是控制限，控制限上方的统计量表示检测到一个故障(统计量用实线表示)。对于故障8，基于KPCA的监测结果如图3所示。非线性主元分析如图4所示。基于原始KICA的监测结果如图5所示。基于改进的KICA的监测结果如图6所示。在故障8的情况下，流4的A，B，C原料成分是随机改变的。在训练的过程中LS-SVM需要正常训练样本和故障训练样本。对于本次故障，KPCA能够在大约第185个样本检测到故障(图3)且延时时间是25s。然而，尽管发生故障，但是有些样本(样本210-220，370-390，630-650，750-760，790-820以及895-905)仍低于95％的控制限，这就给过程操作员带来一幅错误的过程状态图。原始KICA监测图显示统计量成功地从第180个样本检测到故障(图5)延时时间是20s。且比起低于95％的控制限的KPCA需要更少的样本。从改进的KICA监测图中可以看出，它和原始的KICA监测结果基本相同，但通过相似性分析，改进的KICA降低了计算负载。由表5可知，改进的KICA的整个监测率明显高于KPCA的。其中故障3，9和15的检测率比所有的故障检测率都要低很多。对于故障3，9，15，任何一个统计量都将导致低的检测率，换句话说，故障3，9以及15从实时工况数据中是不可观测的。
表5.田纳西-伊斯曼过程中每种方法的故障诊断率 分类器的LS-SVM输入变量的选择对它的分类性能有着很重要的影响。输入变量是根据过程知识和操作特性进行选定的。(Chiang，2001)尽管考虑了有着不同类型和位置的三种故障，这三种故障是故障4(反应器冷水入口温度)，故障6(A进给丢失(流1))和故障11(反应器冷水入口温度)。然而，只有变量1对故障6是重要的，只有变量9和变量51对故障4和故障11是重要的。其余的49个变量就分类而言不提供重要信息，也就是，对于故障4，6和11的训练数据和实时工况数据观测不到这49个变量有什么不同。故障4的数据集同故障11重叠。从表6和7可以看出，通过使用主元nPCs和独立元ICs，从原始信息中提取出了有用的特征。且使用ICs的分类率高于nPCs的分类率。原因是改进的KICA中的负熵能比KPCA更好地解释高阶原始输入信息。从表8中，发现改进的KICA和LS-SVM结合的训练时间比原始KICA与LS-SVM结合的更短。
表6.田纳西-伊斯曼过程中LS-SVM利用nPCs后的故障非类率 表7.田纳西-伊斯曼过程中LS-SVM利用ICs后的故障非类率 表8.田纳西-伊斯曼过程的训练时间 实例2、废水处理过程(WWTP) 把基于改进的KICA的监测方法用于废水处理过程。它包括活性污泥模型No.1(ASMl)和一个十层沉淀池模型(settler model)，活性污泥模型No.1和一个十层沉淀池模型分别用于仿真生物反应和沉淀过程。WWTP系统的工艺布局如图7所示。生物反应器的第一二个分室不充气，其它的充气。所有分室都被认为是理想混合的，而二级沉淀池用一系列一维十层进行建模。对于过程监测，选定了8个变量如表9所列，因为他们都很重要且在实际WWTP系统监测中比较典型。这种过程的变量往往一个周期有很大的波动，他们的均值和方差不能够一直保持为常量。根据这种情况，常用的多变量统计过程监测(MSPM)，隐含地假定了一个稳定的潜在过程，会导致大量的错误警报及故障丢失。所以提出了本发明方法对此过程进行监测。
表9.WWTP基准的监测变量 步骤一、数据采样 采集如表9所示变量的数据。包括一周的标准数据(训练数据)和两周的实时工况数据，其中采样周期为15分钟，由此获得了672组标准数据和1344组实时工况数据。在实时工况数据中，此过程经历了长期干燥的天气后突然下了两场暴雨，第一场暴雨样本850-865，第二场暴雨样本1050-1110。因为WWTP受到每日流量大的波动及进给流的成分的影响，数据有一定的周期特性如图8所示。分别随机选取了标准数据和实时工况数据中的十组数据如表10和表11所示 表10.WWTP标准数据中的十组数据
表11.WWTP实时工况数据中的十组数据
步骤二、进行相似性分析 对标准数据进行输入空间的相似性分析，取γ0＝0.95 。学习的开始只有两个数据点，即节点集N1＝{x1，x2}。当引入数据点x3时，根据公式(1)进行计算，当时，x3加入节点集N2中，得到N2＝{x1，x2，x3}，否则N2＝N1。依次类推，每引入一个新的数据点，都根据公式(1)进行计算，直到引入x672为止，得到第一个数据子集Ω1。
再对数据子集Ω1中的数据进行特征空间的相似性分析，取γ1＝0.95，核函数选径向基核函数。学习的开始也只有两个数据点，即节点集L1＝{x1，x2}。当引入数据点x3时，根据公式(13)、(14)进行计算，当且时，x3加入节点集L2，得到L2＝{x1，x2，x3}，否则L2＝L1。依次类推，直到Ω1中最后一个数据点进行相似性分析，得到第二个数据子集Ω2。
步骤三、利用KPCA对数据做白化处理 选取径向基核函数作为核函数，其中σ＝0.5。对第二个数据子集Ω2中的数据根据公式(16)进行核计算，得到Gram核矩阵K。接着根据公式(17)对Gram核矩阵K进行中心化处理，得到中心化的核矩阵

通过

计算出了它的特征值，并选取了4个最大的正特征值为λ1≥λ2≥...≥λ4。根据公式(18)得到4个特征值对应的特征向量矩阵H。通过标准数据得到公式(22)所需参数后，就可以在特征空间利用公式(22)对实时工况数据进行白化处理，得到白化后的观测变量zt。
步骤四、利用修正ICA提取独立元 取p＝3，根据ICA算法计算训练数据的标准得分转换矩阵Cn，其中g(u)＝tanh(1.5u)。接着利用上一步得到的白化后的观测变量zt，根据公式st＝D1/2CnTzt计算得到实时工况数据的独立元st。
步骤五、故障检测 利用上一步得到的实时工况数据的独立元st根据公式(26)和(27)计算实时工况数据的T2统计量和SPE统计量。得到改进的KICA监测结果如图10所示，其中点线为它们的控制限，实线为统计量的分布曲线。通过实时工况数据的统计量分布，可以检测到故障发生。对于原始的KICA，得到如图9所示的监测结果。从图9和图10可以看出，对于这两场暴雨，原始的KICA能够成功的对其进行监测，改进的KICA也成功地对两场暴雨进行了监测。但是改进的KICA降低了计算负载，而且达到的监测结果与原始的KICA基本一致。由此表明改进的KICA比原始的KICA更有优势。此实例由于故障就是两场暴风雨的引入，所以不再对故障的种类进行诊断。
本发明——基于改进的核独立元分析和支持向量机相结合的非线性过程故障诊断方法还能对其他的非线性过程像连续采煤机采煤过程、轧钢过程等进行故障诊断。通过对他们的监测，及时检测出故障，避免了系统崩溃及由此造成的物质损失和人员伤亡，提高了工业过程的生产效益。
权利要求
1.一种非线性过程故障诊断方法，其特征在于该方法包括以下步骤
步骤一、采集数据
采集过程中相关变量的数据，对于每个故障，产生两组数据，即训练数据和实时工况数据；训练数据用于建立模型，实时工况数据用于在线监测；并用均值和标准偏差规范化采集的数据；
步骤二、进行相似性分析
相似性分析包括输入空间的相似性分析和特征空间的相似性分析，根据相似性因子对数据进行分析，去除掉相似性较强的数据；
步骤三、利用核主元分析对数据进行白化处理，求解白化后的观测变量z
通过非线性映射将输入空间映射到一个特征空间接着在此特征空间对观测数据进行白化处理，得到白化后的观测变量z；
步骤四、利用修正独立元分析ICA提取独立元
在核主元分析KPCA变换空间提取独立元s；利用修正ICA方法从白化后的观测变量z中提取一组独立元；并使独立元各变量之间相互统计独立；
步骤五、利用T2和SPE统计量和LS-SVM进行故障检测与诊断
采用T2和SPE统计量进行在线故障检测，当观测数据的统计量没有超出统计量规定的控制限，则属于正常数据，反之属于异常数据，表明出现故障；一旦故障发生，将实时工况数据的独立元作为LS-SVM的输入进行故障诊断，根据训练数据建立的决策函数，即可判断故障的类型。
2.根据权利要求1所述的一种非线性过程故障诊断方法，其特征在于所述的步骤二中的相似性分析的具体过程如下
1)输入空间的相似性分析
对采集的数据进行相似性分析；
学习的开始，只有两个数据点，即假设第i个节点集被表示为i＝1，...，t1，且满足ni＜t1，t1为训练数据的个数；定义相似性因子为
式中xi+1-训练数据中的第i+1个数据点；
-节点集Ni中任意两个数据点；
S1满足0≤S1≤1；当S1＝1时相似性最强；将xi+1同节点集Ni中所有数据都两两进行相似性计算，假如所有其中γ0是满足的预定义小变量，新数据xi+1将被引入建立新扩展节点集Ni+1＝{Ni，xi+1}；否则，这个数据点将被拒绝，即Ni+1＝Ni；在输入空间对所有训练数据进行相似性分析后，获得第一个数据子集Ω1；
2)特征空间的相似性分析
得到子集Ω1后，将Ω1中的数据通过非线性函数Ф映射到特征空间，进行特征空间的相似性分析；学习的开始，只有两个数据，即假设特征空间第j个节点集被表示为j＝1，...，t2，且满足nj＜t2，t2为子
集Ω1中数据的个数；定义特征空间的相似性因子为
式中Ф(xj+1)-数据子集Ω1中第j+1个数据点xj+1映射后的数据；
-节点集Lj中的任意两个数据；
由于非线性函数Ф难以直接确定，利用核技巧k(a，b)＝<Ф(a)，Ф(b)>，其中a，b为输入空间的数据，对此问题进行了解决，核函数采用径向基核函数
经过一定转换，特征空间中相似性因子转化为在输入空间进行计算，定义如下
且
其中S21，S22满足0≤S21≤1，0≤S22≤1，当S21＝1且S22＝1时相似性最强；因此不需要对子集Ω1中的数据进行非线性映射就可在特征空间对数据进行相似性分析，重新定义第j个节点集假如且其中γ1是一个满足的预定义小变量，xj+1将被引入建立新扩展节点集Lj+1＝{Lj，xj+1}；否则这个数据点将被拒绝，
即Lj+1＝Lj；在特征空间对Ω1中的数据进行相似性分析后，保留的数据点个数表示为N；且获得了第二个数据集Ω2。
3.根据权利要求1所述的一种非线性过程故障诊断方法，其特征在于所述的步骤三中求解白化后的观测变量z的具体过程如下
子集Ω2中包含m个变量的观测数据xk∈Rm，k＝1，...，N，其中N是观测的个数；利用非线性映射ФRm→F，原始空间中的观测数据就扩展到高维特征空间F，Ф(xk)∈F；
此特征空间的协方差结构是一个单位矩阵；因而，特征空间中协方差矩阵将为
其中Ф(xj)，j＝1，...，N，假定为零均值和单位方差；令Θ＝[Ф(x1)，...，Ф(xN)]，因而CF能表示为定义一个N×N维的Gram核矩阵K
[K]ij＝Kij＝<Ф(xi)，Ф(xj)>＝k(xi，xj)(16)
有K＝ΘTΘ；核函数k(xi，xj)的应用可以在不执行非线性映射的情况下在F中计算内积；即可以避免执行非线性映射，并且在特征空间通过引入一个核函数k(x，y)＝<Ф(x)，Ф(y)>计算内积；从核矩阵K可知，高维空间中Ф(xk)的规范化能按以下方法执行，即将Ф(xk)的中心化转化为对K的中心化处理；中心化核矩阵
可有下式获得
其中
把
特征值分解
式中λ为
的特征根，α为λ对应的特征向量；由式(18)能够获得同
的d个最大正特征值λ1≥λ2≥...≥λd相应的标准正交特征向量α1，α2，...，αd；因而CF的d个最大正特征值是相应的特征向量v1，v2，...，vd能够表示为
j＝1，...，d(19)
特征向量矩阵V＝[v1，v2，...，vd]可以简单的由下式表达
V＝ΘHΛ-1/2(20)
其中Λ＝diag(λ1，λ2，...，λd)和H＝[α1，α2，...，αd]分别为
的d个最大特征值的对角阵及对应的特征向量；特征空间中的白化矩阵P和映射数据的白化变换如下
z＝PTФ(x)(21)
具体为，
式中
k＝[k(x1，x)，...，k(xN，x)]T，1t＝(1/I)[1，...，1]R1×N；
x1，...，xN，—Ω2中的数据；
x-需要白化处理的数据。
4.根据权利要求1所述的一种非线性过程故障诊断方法，其特征在于所述的步骤四求解独立元s的具体过程如下
找出p(≤d)个独立元，即s＝{s1，s2，...，sp}，满足E(ssT)＝D＝diag{λ1，...，λp}以使s的元素变得尽可能彼此独立，利用
s＝CTz(23)
其中C∈Rd×p为得分转换矩阵和CTC＝D；定义归一化的独立元为
sn＝D-1/2s＝D-1/2CTz＝CnTz(24)
Cn为标准得分转换矩阵，很明显D-1/2CT＝CnT，CnTCn＝I，以及E(snsnT)＝I；因此算法的目标改为从z∈Rd寻找sn∈Rp和Cn，使sn的元素变得尽可能彼此独立且满足E(snsnT)＝I；设定z的第一个数据中p个元为sn的初始元；令CnT＝[Ip0]，其中Ip是p维的单位矩阵且0是p×(d-p)的零矩阵；根据修正ICA算法计算标准得分转换矩阵Cn
(1)选择p，估计独立元的个数，设定计数器i←1；
(2)获取初始向量ci；
(3)令ci←E{zg(ciTz)}-E{g′(ciTz)}ci，其中g′是g的一阶导数，Hyvrinen给出了三种g的表达函数为g1(u)＝tanh(a1u)，g2(u)＝uexp(-a2u2/2)，g3(u)＝u3；
(4)执行正交化这种正交化排除了现有方案所包含的信息；
(5)归一化
(6)假如ci还没有收敛，返回步骤(3)；
(7)假如ci收敛，输出向量ci，若i≤p，则设定i←i+1且返回步骤(2)；
一旦找到Cn，根据下式得到独立元
s＝D1/2CnTz(25)。
5.根据权利要求1所述的一种非线性过程故障诊断方法，其特征在于所述的步骤五中检测统计量的计算以及决策函数的建立过程如下
T2统计量和SPE统计量定义如下
T2＝sTD-1s(26)
z＝PTФ(x)(27)
其中
能够由下式得来
因为s不服从高斯分布，T2的控制限通过F分布来决定；SPE的控制限由下面权重χ2分布来计算(Qin，2003)
SPE～μχh2，μ＝b/2a，h＝2a2/b(29)
其中a和b分别是标准操作数中SPE的估计均值和方差；
要对故障进行分类，首先要建立各个故障的决策函数，具体如下
LS-SVM分类器在F中构造一个最大边界(margin)超平面；在输入空间，这与一个非线性决策边界相对应；
其中w表示权重向量且b为偏置(bias)；两种情况都能合成为单一框架如下
yi[<w，Ф(si)>+b]≥+1-ξi
其中ξi≥0，i＝1，2，...，N是正松弛变量，用于允许不等式集的错误分类；因为两个超平面之间的边界宽是2/‖w‖，目标函数被定义为‖w‖2和惩罚项的最小值
其中C∈R+是一整定超参数，它决定边界最大化和训练误差最小化之间的权衡；松弛变量ξi≥0，i＝1，2，...，N是用于放宽硬边界LS-SVM分类器的约束以便允许存在一些被错误分类的数据；在LV-SVM中，当引入核技巧后，方程(31)经过一定的运算后，能重新描述为
式中If-IfT＝[1，...，1]∈R1×N；
Ω-Ωij＝yiyjФ(si)TФ(sj)＝yiyjk(si，sj)；
α-拉格朗日因子；
根据方程(32)求出偏置b和拉格朗日因子α，根据这两个参数确定相应的决策函数即分类器如下
式中s-需要判断故障类型数据的独立元；
si-训练数据的独立元，i＝1，...，N；
决策函数建立后，当检测到故障，输入实时工况数据的独立元，根据f(s)的值判断故障的类型，选取其中一个故障的决策函数，当f(s)＝+1时，发生本决策函数对应的故障，反之f(s)＝-1时，没有发生本决策函数对应的故障，继续根据其他决策函数进行故障诊断。
全文摘要
一种非线性过程故障诊断方法，包括采集数据、进行相似性分析、利用核主元分析对数据进行白化处理，求解白化后的观测变量z、利用修正独立元分析ICA提取独立元、利用T2和SPE统计量和LS-SVM进行故障检测与诊断步骤。本发明提出了非线性动态过程故障诊断技术，结合Kernel、ICA和LS-SVM三者的优点，即，发挥Kernel对非线性的表达能力，同时发挥ICA对动态特性的把握能力和LS-SVM的分类能力。
文档编号G05B23/02GK101158873SQ20071001295
公开日2008年4月9日 申请日期2007年9月26日 优先权日2007年9月26日
发明者张颖伟, 秦泗钊, 滢 王 申请人:东北大学