多模态化工过程动态多点故障监测方法

胡瑾秋， 罗 静， 郭 放

(中国石油大学(北京) 油气资源与工程国家重点实验室 机械与储运工程学院， 北京 102249)

实际生产中，多模态化工过程除了可能存在非高斯特性以外，还会受到环境和操作干扰，产生噪声。当生产过程处于不同生产模态时，正常的过程数据的均值、方差和相关关系等特征变量会出现较大的波动。目前，在多模态化工过程监测领域应用最多的是多元统计过程故障监测技术[1]。多元统计方法不需要根据过程机理建立复杂模型，并且可以降低数据维度。近年来，国内外学者也开展了相应的研究，主要分为以下两类：基于主元分析方法的技术和基于独立成分分析方法的技术。大多针对多模态化工过程的特点采取主元分析与其他方法相结合的手段对过程进行监控。卢春红等[2]采用主元分析结合核函数的方法处理高斯数据的非线性问题，从概率角度刻画数据集的多个局部分量特征，在提取的核主元分量内获得测试样本的后验概率，结合贝叶斯推理进行故障检测。Zhang等[3]采用递归核主元分析法应用于连续退火工艺的动态过程监测，能够有效地监测出过程故障。部分学者用主元分析结合聚类算法和奇异值识别算法来解决过程的多模态问题[4-8]。基于独立成分分析方法的技术，主要处理过程的非高斯问题，也与其他方法结合处理多模态过程中的其他问题。赵小强等[9]采用独立成分分析(Independent component analysis, ICA)结合核函数解决非高斯数据的非线性问题，进行化工过程的故障检测。还有一些将ICA[10]结合聚类算法、奇异值算法和其他方法解决过程的多模态问题[11-15]。

然而，虽然已有的多元统计监测方法解决了多模态化工过程的非高斯性，但是这些方法均采用静态控制限对统计量进行监测，鲁棒性差，会因噪声产生误报。因此，针对上述问题，笔者提出多模态化工过程动态多点故障监测方法。基于粒子群优化的ICA算法和自回归(Autoregressive, AR)模型构造不同平稳过程的非高斯监测模型，计算平稳过程的单点监测统计量和多点异常统计量。基于粒子群优化的ICA算法构造过渡过程的非高斯监测模型。平稳过程和过渡过程均采用动态控制限实施监控。案例分析中，分别将传统的统计量控制限求解方法3σ法与动态多点监测方法应用到丙烯计量罐装置进行对比。

1 动态多点故障监测方法基本理论

1.1 AR模型

为了解决多模态化工过程监测数据的非高斯性和高维度，引入粒子群优化的ICA算法来计算不同过程模态的非高斯统计量，基于非高斯统计量建立平稳过程AR监测模型。AR模型要求时间序列{Xt}趋势平稳、均值为0。AR模型是p阶自回归模型，记为AR(p)，其中p是模型的阶数，N是时间窗宽，如公式(1)所示。

Xt=φ1Xt-1+φ2Xt-2+…+φpXt-p+εt

(1)

式中：φ1,φ2,…,φp是AR(p)的自回归系数，εt是均值为0、方差为σ2的独立同分布高斯随机白噪声。AR模型的阶数p与窗口N满足约束条件：

0≤p≤0.1N

(2)

(3)

(4)

(5)

1.2 单点监测统计量

参考网络流量异常定义准则，用观测值与AR模型预测值的残差来定义化工过程异常，而过程的最终监测统计量则采用单点监测统计量和多点异常统计量。

(1)残差e

定义零均值化后的观测值序列为{…,x(t+1),x(t+2),x(t+3),…}，由AR模型拟合所得的预测值序列为{…,y(t+1),y(t+2),y(t+3),…}，那么，残差序列{…,e(t+1),e(t+2),e(t+3),…}按照下式计算：

e(t+i)=x(t+i)-y(t+i)

(6)

(2)单点监测统计量W

(7)

其中，ξ2=(e2(t+1)+e2(t+2)+…+e2(t+N+1))/(N+1)。通过单点监测统计量W(t+N+1)判断预测值y(t+N+1)是否正常，当W(t+N+1)>U(t+N+1)时，y(t+N+1)是异常的，否则，y(t+N+1)是正常的。U(t+N+1)代表单点监测统计量W(t+N+1)在t+N+1时刻的控制限，数值由下式计算：

U(t+N+1)=μ+k×σ

(8)

其中，μ和σ分别是正常历史观测值对应的残差正序列的均值和标准差，k取值为2或3，初值为2。若前一时刻的预测值y(t+N)被判断为正常状态，则令当前时刻k=2。当y(t+N)指示异常时，为了减少过程干扰产生的误判，此时令k=3，即适当降低检测点的监测标准。

1.3 多点异常统计量

多模态化工过程由于其自身的复杂特性和员工的频繁操作，即使处于正常状态，过程检测点数值也有可能具有短暂的较大波动。化工过程异常时，不会只有某个孤立的检测点异常，势必会引发链式效应。因此，为了确保有效报警率，有必要设置多点异常统计量λ，根据多个连续检测点的异常情况来判定是否需要报警。定义λ为当前检测点距离前一个异常检测点的时间间隔a和一定时间内异常发生次数的函数。在一定时间内，过程发生异常的次数越多，λ值越大，过程异常程度越大。记λt为检测点在单点时刻t的异常统计量，初值为0。在初次检测到异常发生时，根据单点监测统计量超出控制限的部分计算λt和λ，如式(9)和式(10)所示：

λt=W(t+N+1)-U(t+N+1),

W(t+N+1)>U(t+N+1)

(9)

(10)

当多点异常统计量λ超过其控制限Uλ，说明在当前时间窗口中有多个连续的检测点发生异常，报告异常发生。λ的控制限Uλ为发生异常的点数n与k=3时的单点监测统计量控制限的积，计算如下式：

Uλ=n×(μ+3σ)

(11)

为避免误报警，需进一步计算多点异常统计量λ来判断异常。增大k值，减小λt，放缓多点异常统计量λ的增长。如果有2个以上连续相邻的点的单点异常统计量λt超过阈值，会造成多点异常统计量λ以2倍指数的形式呈现爆炸式增长，并迅速越过报警限。如果没有连续异常点出现，单点监测统计量W(t+N+1)必须远远超出阈值，λ才会越过控制限。若单点异常统计量λt不大，且相邻的点没有检测到异常，则多点异常统计量λ会随着时间的推移逐渐递减为0。

2 动态多点故障监测方法

针对目前多元统计监测方法中所采用的统计量静态控制限、鲁棒性差、易因噪声产生误报的问题，提出多模态化工过程动态多点故障监测方法。基于自回归(Autoregressive, AR)模型和粒子群优化的ICA算法，构造平稳模态的单点监测统计量和多点异常统计量，建立起平稳模态的非高斯监测模型。基于粒子群优化的ICA算法构造过渡模态的非高斯监测模型，平稳模态监测模型和过渡模态监测模型均采用动态监控策略，实现在线故障监测。以Hotelling’s T2统计量为例展示动态多点故障监测方法的具体步骤，平方预测误差(Squared Prediction Error, SPE)统计量同样适用于该方法。

2.1 离线建模

步骤1：挑选模型训练数据。对多模态化工过程变量的正常历史数据做Lilliefors test正态性检验，由检验标准挑出n个非高斯变量，如表1所示。取这n个非高斯变量的M个时刻的正常历史数据作为训练样本X∈Rn×M；

表1 正态分布检验标准Table 1 Normal inspection standards

步骤2：模态划分。用模糊C均值(Fuzzy C-means, FCM)算法将训练样本x离线划分为c个不同的平稳模态和c-1个过渡过程；

(12)

(13)

2.2 在线监测

Us,i(k)=μei+d1×σei

(14)

步骤5：判断单点监测统计量Wi(k)超限与否。若Wi(k)>Us,i(k)，d1=3，则继续步骤8；否则，d1=2，λi(k)=0，k=k+1，返回步骤1。

步骤6：根据式(9)和(10)计算当前时刻的多点异常统计量λi(k)，并剔除当前观测值。

步骤7：根据式(11)计算当前时刻多点异常统计量控制限Uλ,i(k)。

步骤8：判断多点异常统计量λi(k)超限与否。若λi(k)>Uλ,i(k)，则报告异常发生，ni=ni+1；否则，k=k+1，返回步骤1。

步骤9：根据式(15)计算当前时刻过渡过程统计量控制限Ut(k)

Ut(k)=μt(k)+d2×σt(k)

(15)

2.3 方法流程图

图1为动态多点故障监测方法步骤图。

3 案例分析

丙烯投料计量控制系统的作用是准确计量投入聚合釜中的丙烯量，为聚合反应控制提供前提，如图2所示。平时丙烯经气动三通球阀直接回丙烯中间罐循环，由可编程控制器控制气动三通球阀进行投放料。丙烯投料计量控制系统的监控变量如表2所示。

每5 s采集1组数据，每种状态下各采集537组数据。过程运行状态描述：(1)正常状态下，前1020 s，进料压力PI6113控制在1.7 MPa以内，进料流量FT6101控制在0.01 t/h，其他变量保持稳定，过程处于平稳状态。在第1025 s到1760 s，调小进料流量FT6101，使其他变量发生变化，过程处于过渡状态。在1765 s以后，进料流量FT6101控制在(10.65±0.35) t/h，过程再次处于平稳状态。(2)故障状态下，前1270 s与正常状态同样的操作条件，第1275 s至第1430 s，进料阀门开度过大，进料流量FT6101大幅升高，经调节进料阀门开度，在第1495 s过程恢复正常状态。

图1 动态多点故障监测方法步骤Fig.1 Dynamic multi-point fault monitoring method steps

图2 丙烯投料计量控制点流程图Fig.2 Propylene feed metering control point flow chart

3.1 动态多点故障监测方法

挑选正常状态下所有监测变量的537组数据作为监测模型的训练样本X∈R5×537，537组故障数据作为测试样本Test∈R5×537。用FCM算法将训练数据划分为2个平稳模态和1个过渡过程(如图3所示)，并针对不同模态的训练数据分别建立各个模态下的动态多点故障监测非高斯模型。

图3 模态划分结果Fig.3 Result of mode division

采用事先建立的故障监测模型来监测测试样本数据，并绘制监测曲线，如图4、图5和图6所示。图4(a)为平稳模态1的Hotelling’s T2单点监测统计量监测曲线。可以看出，曲线在第765 s至第880 s间有较大波动，在个别点控制限有所提高，但在第840 s和第860 s仍有数据超过控制限，这两点是连续异常点，初步判断此时系统异常。图4(b)为平稳模态1的Hotelling’s T2多点异常统计量监测曲线，只有在第840 s和第860 s，多点异常统计量有数值，但均未超过控制限。故平稳模态1的Hotelling’s T2监测结果为全过程处于正常状态。图4(c)是平稳模态1的SPE单点监测统计量监测曲线，可以看出，曲线随着时刻增长随机波动，整体上无明显趋势。控制限数值几乎固定在2.7，只有在第740 s、第870 s和第970 s控制限数值稍有提高，但无数据点超过控制限，系统全过程处于正常状态。

图4 平稳模态1的监测结果Fig.4 Monitoring results of the stationary mode 1(a) Hotelling’s T2 single point monitoring statistic; (b) Hotelling’s T2 multiple points anomaly statistic; (c) SPE single point monitoring statistic

图5(a)为过渡过程的Hotelling’s T2统计量监测曲线。可以看出，前1270 s曲线随着时刻增长有下滑趋势，在第1275 s发生阶跃式增长，并越过动态控制限，系统出现异常，在第1465 s发生突降，回落至动态控制限以下，此后的过渡阶段曲线趋于平稳。图5(b)为过渡过程的SPE统计量监测曲线。可以看出，前1270 s曲线随着时刻增长有小幅随机波动，在第1275 s越过动态控制限，系统出现异常，在第1465 s发生突降，回落至动态控制限以下，此后的过渡阶段曲线在动态控制限以下维持小幅波动。

图5 过渡模态监测结果Fig.5 Transition mode monitoring results(a) Hotelling’s T2 statistic; (b) The SPE statistic

图6(a)为平稳模态2的Hotelling’s T2单点监测统计量监测曲线。可以看出，曲线在第1865 s至第2530 s间虽有波动，但整体趋势较平稳。在第2535 s，统计量数值略微增长，此后维持小幅随机波动。在第2435 s和第2655 s，曲线存在突然增长点，控制限也相应提高。除了这两点，整个平稳模态2，控制限几乎保持平稳，而Hotelling’s T2单点监测统计量没有超过控制限，系统处于正常状态。图6(b)为平稳模态2的SPE单点监测统计量监测曲线。可以看出，曲线随着时刻增长上下大幅波动，没有明显趋势。在第2565 s至第2595 s间有个别点波动较大，控制限也相应提高。过程的所有数据点均未超过控制限，说明系统维持在正常状态。另外，平稳模态2无连续异常点出现，故多点异常统计量始终为0。

图6 平稳模态2的监测结果Fig.6 Monitoring results of the stationary mode 2(a) Hotelling’s T2 single point monitoring statistic; (b) SPE single points monitoring statistic

3.2 基于非高斯模型的3σ法

为了验证本文所提方法的有效性，与3σ法的监测结果对比。根据3σ法，统计量控制限为统计量均值与其3倍标准差的和。应用3σ法监测Hotelling’s T2统计量和SPE统计量，监测结果如图7所示。

3.3 结果分析

将动态多点故障监测方法与传统的基于非高斯模型的3σ法的误报率和漏报率统计于表3。

图7 3σ法监测结果Fig.7 Monitoring results of 3σ(a) Hotelling’s T2 index of the whole process; (b) SPE index of the whole process; (c) Hotelling’s T2 index of transition process; (d) SPE index of transition process

MethodT2SPEFalse alarm rate/%False negative rate/%False alarm rate/%False negative rate/%Dynamic multi-point fault monitoring method0.190.191.120.743σmethod0.930.197.450.74

结果表明，2种方法的漏报率相同，但相比较基于非高斯模型的3σ法，动态多点故障监测方法的T2统计量和SPE统计量的误报率分别降低了0.74百分点和6.33百分点。

4 结 论

(1) 针对传统方法所采用的静态控制限因不能排除噪声的干扰而产生误报警的问题，提出多模态化工过程动态多点故障监测方法。用粒子群优化的ICA算法计算不同过程模态的非高斯统计量，平稳过程基于自回归模型构造非高斯统计量的单点监测统计量和多点异常统计量，采用动态控制限监测，过渡过程直接采用动态控制限对非高斯统计量进行监测。

(2) 案例分析中，将动态多点监测方法应用到丙烯计量罐装置，离线划分过程模态，计算不同过程模态的非高斯统计量，构造平稳过程的单点监测统计量和多点异常统计量，分别采用动态控制限进行平稳模态和过渡模态的监测。

(3) 结果表明，两种方法的漏报率相同，但相比较基于非高斯模型的3σ法，动态多点故障监测方法的T2统计量和SPE统计量的误报率分别降低了0.74百分点和6.33百分点。因此，新方法能够提高多模态化工过程故障监测的准确率。

符号说明：

a———时间间隔,s；

c——模态种类；

d1、d2——系数；

e——残差；

f——样本批次；

M——时刻，s；

m——过渡过程异常数据个数；

N——时间窗宽度，s；

n——平稳模态异常数据个数；

p——模型阶数；

t——时间，s；

U——监测统计量控制限；

W——单点监测统计量；

X——时间序列矩阵；

xt——训练样本；

Y——模型预测序列；

y——预测样本；

εt——白噪声；

φi——模型系数；

σ——标准差；

λ——多点异常统计量；

μ——均值。