一种新型多模态化工过程异常状态检测方法

本发明涉及一种化工过程运行状态监测方法，特别涉及一种新型多模态化工过程异常状态检测方法。

现代化工过程在运行中，通常会不间断地实时地通过测量仪器仪表采集生产过程中的温度，压力，流量等数据信息。这些海量存储的采样数据以及实时更新的在线采样数据为当前化工“大数据”建设与智能化工应用奠定了坚实的数据基础。近十几年来，利用采样数据来检测化工过程运行中出现的异常状态在安全化工生产领域受到了越来越多的重视，其思想在于通过采样数据来检测化工过程运行是否进入异常状态，从而保证安全生产与维持产品质量的稳定性。无论是学术界还是工业界，都投入了大量的人力与财力研究数据驱动的异常状态检测方法与技术，应用对象除涵盖化工领域外，还涉及机械制造业与工业控制系统安全领域。数据驱动的化工过程监测经过多年的研究发展，出现了许多以主成分分析(Principal Component Analysis，缩写：PCA)和偏最小二乘(Partial Least Square，缩写：PLS)为基础的异常状态检测方法。这些主流的过程监测方法实施的核心主要关注于数据潜在特征的挖掘。换句话说，所建立的数据驱动模型都是旨在提取采样数据数据中潜藏的特征。

一个典型的化工过程对象的核心生产单元都包含有：反应设备、分离设备、冷却设备、和回流装置。这些设备之间相互链接，运行机理与控制系统回路较为复杂。此外，由于化工厂的订单要求波动、能耗调度问题等原因，化工过程对象并非一直运行在某个固定的生产状态下。即使不考虑各种异常状态情况，正常的运行状态就具备多种生产模态。例如，因化工企业的高能耗原因，通常会根据工业用电早晚价格的差异，适当调整白天与黑夜的生产产量。这种运行在多个正常工况下的化工过程称之为多模态化工过程，相应的异常状态检测就称为多模态的化工过程异常状态检测。在现有的科研文件与专利技术材料中，针对多模态化工过程异常状态检测的研究主要是为每个生产模态都各自建立相应的检测模型。例如，最经典的方法就是：首先区分出各个生产模态下的采样数据，然后利用PCA算法为各个模态分别建立相应的异常状态检测模型，最后分别利用多模态PCA模型对多模态化工过程实施在线的异常状态检测。

值得指出的是，即使化工生产过程对象有多个正常运行的生产模态，每个生产模态都是在同一套生产流程上完成的。因此，各个生产模态下的采样数据之间是存在潜在关系的。虽然多模态PCA模型很好的区分了各个模态，但是模态之间的关联特征却忽略了。正如前面所提到的，数据驱动的异常状态检测方法的核心在与如何充分有效地挖掘采样数据的潜在特征。从这个角度讲，经典的多模态PCA方法未曾考虑各个模态之间的关系，因此多模态PCA方法从实施路线上来看是无法挖掘多模态采样数据之间的关联特征。因此，发明出一种可以挖掘多模态之间相关特征的多模态化工过程建模算法并将之应用于异常状态检测是一个丞待解决的问题。

本发明所要解决的主要技术问题是：发明出一种可以挖掘多模态之间相关特征的多模态化工过程建模算法并将之应用于多模态化工过程异常状态检测。具体来讲，本发明方法首先根据多模态采样数据的两个特征：其一，多模态采样数据之间是存在相关性的；其二，多模态采样数据具备各自的独特性，来构造新型多模态建模算法的目标函数；其次，本发明方法在朗格朗日乘子法的基础上结合迭代算法来优化求解出能挖掘多模态之间相关特征的多个投影向量；最后，实时在线异常状态检测时，本发明方法对多模态相关特征与独特特征分别进行检测，从而实现对多模态化工过程对象的异常状态检测。

本发明方法解决上述问题所采用的技术方案为：一种新型多模态化工过程异常状态检测方法，包括以下步骤：

步骤(1)：利用化工过程对象中的测量仪器仪表采集该化工过程正常运行在C个生产模态下的样本数据，并按照生产模态的归属将样本数据区分成C个数据矩阵X1，X2，…，XC，其中，为第k个生产模态下的数据矩阵、下标号k＝1，2，…，C、表示m×Nk维的实数矩阵、Nk表示第k个生产模态下采集的样本数据总数、m为该化工过程对象的测量变量总数、数据矩阵Xk中的列向量为样本数据。

步骤(2)：确定出N1，N2，…，NC中的最小值，并将该最小值记为n后，分别从数据矩阵X1，X2，…，XC中随机选择n个列向量对应组成矩阵其中表示从第k个生产模态下的数据矩阵Xk中随机选择n个列向量组成的矩阵。

步骤(3)：将C个矩阵合并成一个矩阵后，对矩阵X实施标准化处理得到新矩阵并将表示成其中，N＝C×n，在新矩阵中的位置分别与在矩阵X中的位置相对应，即分别是经标准化处理后的矩阵。

步骤(4)：利用新型多模态建模算法得到相关特征提取矩阵W1∈Rm×D与独立特征提取矩阵W2∈Rm×(m-D)，其中D表示相关特征成分的个数。

如前所述，新型多模态建模算法需要提取多模态采样数据之间的相关性，其目标在于：通过优化搜寻一个投影向量wd∈Rm×1，使得特征成分向量g1，g2，…，gC两两之间的相关性最大化，其中表示经投影向量wd转换后得到的特征成分向量。

因此，新型多模态建模算法的优化目标函数即可量化为如下所示形式：

约束条件：

上式中，标量Hij的取值如下所示：

若是令wd＝(XTX)-1/2ud，则可将上式①等价转换成如下所示形式：

约束条件：

上式中，Φ＝(XTX)-1/2，上式③中的约束优化问题可通过拉格朗日乘子法进行求解。

首先，利用拉格朗日乘子λ构造拉格朗日函数然后，计算拉格朗日函数L相对于ud的偏导数：

当上式④中的偏导数等于0时，上式②中的目标函数取得最大值，即得到如下所示的等式关系：

从上式⑤中可以发现，这是一个非常典型的特征值求解问题；若是在上式④等号两边同时左乘再结合上式②中的约束条件，即可得到：

从上式⑥中可以看出，拉格朗日乘子λ等于目标函数值，因此上式④中的特征值问题需求解最大的特征值。

然而，在求解式⑤中的特征值问题前，需要确定Hij的具体数值，而Hij的具体数值在i≠j的情况下是待求解特征向量ud的函数。

为此，本发明方法设计出一种迭代循环往复求解特征向量的思路：初始化ud为任意的实数向量后，计算Hij的具体数值，从而求解上式⑤中最大特征值所对应的特征向量ud；然后根据ud再次计算Hij的具体数值，并求解最大特征值问题；如此往复循环即可最终得到ud。

根据上述推理过程，即可归纳总结出利用新型多模态建模算法得到多模态之间的相关特征成分矩阵S∈RD×N和特征提取矩阵W∈Rm×D的具体实施步骤。

步骤(4.1)：初始化下标号d＝1与初始化ud为任意m×1维的实数向量。

步骤(4.2)：根据如下所示公式计算得到对称矩阵Θ：

上式⑦中，下标号i＝1，2，…，C与下标号j＝1，2，…，C。

步骤(4.3)：若d＝1，则求解特征值问题：Θud＝λdud中最大特征值λd所对应的特征向量ud；若d＞1，则求解如下所示特征值问题中最大特征值λd所对应的特征向量ud：

上式⑧中，I表示m×m维的单位矩阵。

步骤(4.4)：判断ud是否收敛，收敛的标准为ud中元素基本不再发生变化；若是，则执行步骤(4.5)；若否，则返回步骤(4.2)。

步骤(4.5)：判断是否满足条件：d＜m；若是，则设置d＝d+1与初始化ud为任意m×1维的实数向量后返回步骤(4.2)；若否，则得到m个最大特征值λ1，λ2，…，λm及其对应的m个特征向量u1，u2，…，um。

步骤(4.6)：对这m个最大特征值λ1，λ2，…，λm进行降序排列后，将最大的D个特征值所对应的D个特征向量组成矩阵U1∈Rm×D，并将剩余的m-D个特征向量组成矩阵U2∈Rm×(m-D)。

步骤(4.7)：分别根据公式W1＝ΦU1与W2＝ΦU2计算相关特征提取矩阵W1∈Rm×D与独立特征提取矩阵W2∈Rm×(m-D)后，再根据公式S＝W1TX计算相关特征成分矩阵S∈RD×N。

步骤(5)：根据公式Q＝diag{STS}计算异常状态检测指标向量Q∈RN×1后，利用核密度估计法确定出在置信限α条件下的控制限Qlim，其中diag{}表示将大括号内的矩阵对角线元素转变成列向量的操作，置信限α的取值范围是95％至99％之间。

步骤(6)：根据公式计算第k个生产模态下的独立特征成分矩阵Vk，再计算Vk∈R(m-D)×n中所有列向量的均值向量βk以及Vk的协方差矩阵Λk。

步骤(7)：根据公式ψk＝diag{(Vk-Bk)TΛk-1(Vk-Bk)}计算异常状态检测指标向量ψk∈RN×1后，利用核密度估计法确定出在置信限α条件下的控制限ψk，lim，其中Bk∈R(m-D)×n由n个均值向量βk组成。

离线建模过程至此已全部结束，接下来就是利用在线采样数据不间断地实施多模态化工过程异常状态的检测。

步骤(8)：在线采集多模态化工过程对象最新采样时刻的样本数据xt∈Rm×1，并确定出样本数据xt所归属的生产模态b，其中下标号t表示最新采样时刻，b∈{1，2，…，C}。

值得说明的是，样本数据xt所归属的生产模态b是按照步骤(1)中的C个生产模态进行划分，因此生产模态b的取值来源于集合{1，2，…，C}中的元素，即{1，2，…，C}。

步骤(9)：对xt实施与步骤(3)中相同的标准化处理得到向量后，根据公式st＝W1Txt与分别计算出相关特征成分向量st∈RD×1与独立特征成分向量vt。

步骤(10)：根据步骤(8)中确定的生产模态归属b，利用如下所示公式计算异常状态检测指标Qt与ψt：

Qt＝stTst ⑨

ψt＝(vt-βb)TΛk-1(vt-βb) ⑩

步骤(11)：判断是否满足条件：Qt≤Qlim且ψt≤ψb，lim；若是，则当前采样时刻该化工该过程运行未出现异常状态，返回步骤(8)继续实施对下一采样时刻的异常状态检测；若否，则执行步骤(12)从而决策是否出现异常状态。

步骤(12)：返回步骤(8)继续实施对下一采样时刻的异常状态检测，若连续3个采样时刻的异常状态检测指标都不满足步骤(11)中的判断条件，则触发异常状态警报；否则，不触发异常状态警报，并返回步骤(8)继续实施对下一采样时刻的异常状态检测。

与经典的多模态化工过程异常状态检测方法(例如：多模态PCA)相比，本发明方法的优势在于：

首先，本发明方法所涉及的新型多模态建模算法能够挖掘出多模态采样数据的两类特征：其一，多模态之间的相关特征；其二，多模态各自的独立特征；相比之下，多模态PCA方法为各个模态采样数据单独建立模型，只考虑到了多模态的独立特征。

其次，本发明方法分别利用多模态之间的相关特征和多模态各自的独立特征实施在线异常状态检测，较全面地利用了挖掘出的潜在特征成分；最后，通过具体的实施案例，对比验证了本发明方法在多模态的连续搅拌反应釜运异常运行状态检测上的优越性。

图1为本发明方法的实施流程示意图。

图2为连续搅拌反应釜的流程图及其相应的测量仪表。

图3为本发明方法与多模态PCA方法在检测异常状态上的实时效果对比图。

下面结合附图和具体实施方式对本发明进行详细说明。

如图1所示，本发明公开了一种新型多模态化工过程异常状态检测方法，下面结合一个具体应用实例来说明本发明方法的具体实施方式。

如图2所示，连续搅拌反应釜(Continuous Stirred Tank Reactor，缩写：CSTR)的流程图及其相应的测量仪表。CSTR生产单元是化工厂最为常见的生产设备，本实施案例中的应用对象是一个涉及放热反应过程的CSTR设备。因此，该CSTR设备需配备一个冷凝器对反应物的出口温度进行降温处理。从图2中可以发现，该连续搅拌反应釜所涉及的测量变量有m＝7个，分别是：进料流量、反应器压力、反应器液位、反应器温度、反应器进料阀门开度、反应器冷凝水流量，以及冷凝器冷却水流量。

由于产量调度原因，该CSTR设备的进料流量会发生变化，进料流量会按照：低产量、中产量和高产量进行相应的调整，因此该化工过程对象可以正常运行在C＝3个生产模态下，是一个典型的多模态化工过程。

步骤(1)：利用CSTR设备的测量仪器仪表采集该化工过程正常运行在C＝3个生产模态下的样本数据，并按照生产模态的归属将样本数据区分成C＝3个数据矩阵X1，X2，X3。

步骤(2)：确定出各生产模态下采集的样本数据总数N1，N2，N3，并确定出最小值记为n后，分别从数据矩阵X1，X2，X3中随机选择n个列向量对应组成矩阵

步骤(3)：将C＝3个矩阵合并成一个矩阵后，对矩阵X实施标准化处理得到新矩阵并将表示成其中N＝3×n。

步骤(4)：利用新型多模态建模算法得到相关特征提取矩阵W1∈Rm×D与独立特征提取矩阵W2∈Rm×(m-D)，再根据公式S＝W1TX计算相关特征成分矩阵S∈RD×N。

步骤(4)中还需涉及如何确定相关特征成分个数D的实施过程，在本实施案例中，具体按照如下所示步骤确定相关特征成分个数D。

首先，对m＝7个最大特征值λ1，λ2，…，λ7进行降序排列后，得到η1≥η2≥…≥η7。

其次，根据公式Δi＝(ηi-ηi+1)/ηi计算差差值变化率Δ1，Δ2，…，Δ6，其中，下标号i＝1，2，…，6。

最后，确定出Δ1，Δ2，…，Δ6中的最大值Δj，那么相关特征成分个数为D＝3。

步骤(5)：根据公式Q＝diag{STS}计算异常状态检测指标向量Q∈RN×1后，利用核密度估计法确定出在置信限α＝99％条件下的控制限Qlim。

步骤(6)：根据公式计算第k个生产模态下的独立特征成分矩阵Vk，再计算Vk∈R(m-D)×n中所有列向量的均值向量βk以及Vk的协方差矩阵Λk。

步骤(7)：根据公式ψk＝diag{(Vk-Bk)TΛk-1(Vk-Bk)}计算异常状态检测指标向量ψk∈RN×1后，利用核密度估计法确定出在置信限α＝99％条件下的控制限ψk，lim，其中Bk∈R(m-D)×n由n个均值向量βk组成。

离线建模阶段至此完成，接下来进入在线异常状态检测阶段，需要实时利用CSTR化工过程对象的在线采样数据。

在本实施案例中，进料流量的调整可直接反映出生产模态的切换，在线采样数据的模态归属能够直接通过进料流量来确定出来。

刚开始实施在线异常状态检测时，该多模态CSTR过程对象运行在中产量生产模态下(即C＝2)，经过一段时间的正常运行后，CSTR设备冷凝器冷却水的温度发生了异常，由20摄氏度跳变至40摄氏度。

步骤(8)：在线采集多模态化工过程对象最新采样时刻的样本数据xt∈Rm×1，并确定出样本数据xt所归属的生产模态b，其中下标号t表示最新采样时刻，b∈{1，2，3}。

步骤(9)：对xt实施与步骤(3)中相同的标准化处理得到向量后，根据公式st＝W1Txt与分别计算出相关特征成分向量st∈RD×1与独立特征成分向量vt；

步骤(10)：根据步骤(8)中确定的生产模态b，利用上述公式⑨和⑩计算异常状态检测指标Qt与ψt：

步骤(11)：判断是否满足条件：Qt≤Qlim且ψt≤ψb，lim；若是，则当前采样时刻该化工该过程运行未出现异常状态，返回步骤(8)继续实施对下一采样时刻的异常状态检测；若否，则执行步骤(12)从而决策是否出现异常状态；

步骤(12)：返回步骤(8)继续实施对下一采样时刻的异常状态检测，若连续3个采样时刻的异常状态检测指标都不满足步骤(11)中的判断条件，则触发异常状态警报；否则，不触发异常状态警报，并返回步骤(8)继续实施对下一采样时刻的异常状态检测。

经典的多模态PCA方法使用两个异常状态检测指标T2与Q，并将T2与Q分别于对应的控制上限进行对比；如图3所示，本发明方法与多模态PCA方法应用于检测该多模态化工过程异常状态上的实时图像。

从图3中对比可以很明显地发现，本发明方法在异常状态检测灵敏度上明显优越于经典多模态PCA方法。因此，可以说本发明方法具有更可靠的性能。

上述实施案例只用来解释说明本发明的具体实施，而不是对本发明进行限制。在本发明的精神和权利要求的保护范围内，对本发明做出的任何修改，都落入本发明的保护范围。