CN108153267A - 一种基于误差主元分析模型的工业过程监测方法 - Google Patents

Abstract

本发明公开一种基于误差主元分析模型的工业过程监测方法，该方法旨在将将原始数据转换成服从或近似服从高斯分布的误差信息，然后对该误差实施监测从而避免因非高斯数据造成的种种不便。具体来讲，本发明方法首先逐一假设过程对象中每个变量的测量数据缺失；其次，根据主元分析(PCA)模型中处理缺失数据的技巧推测出相应缺失变量的估计值；最后，利用假设的缺失数据实测值与估计值之间的误差作为被监测对象，再次建立基于PCA的故障检测模型实施在线故障检测。与传统方法相比，本发明方法不拘泥于原始数据是否满足于高斯分布假设，能较大幅度改善传统PCA方法的故障监测性能，是一种更为优选的过程监测方法。

Description

一种基于误差主元分析模型的工业过程监测方法

技术领域

本发明涉及一种数据驱动的工业过程监测方法，尤其涉及一种基于误差主元分析模型的工业过程监测方法。

背景技术

随着市场竞争的日趋激烈化，保证生产过程的安全性与产品质量的稳定性是企业稳健发展的根本途径。近几十年来，以故障检测为核心要素的过程监测方法与技术一直都是工业界与学术界重点关注的对象。在故障检测研究的起始阶段，基于机理模型与设计观测器是主流的技术手段。通常来讲，误差生成是基于机理模型方法的核心。这类方法大多利用过程对象中某些参数实际值与估计值之间的误差来反应过程运行状态，以达到实时故障检测的目的。这些被监测的误差信息在正常运行状态下一般都是服从高斯分布的，只有在异常状态下才会出现较大波动。然而，随着现代化工业对象规模的不断复杂化与大型化，获取符合一定精度要求的机理模型变得越来越不切实际。取而代之的是，基于数据的过程建模与故障检测方法。这主要得益于先进仪表与计算机技术在工业过程中的广泛应用，海量的生产过程数据可以很容易地被测量与储存。与基于机理模型的故障检测方法生成误差不同，基于数据的故障检测方法旨在挖掘过程数据中的潜藏的有用信息，并对这些有用信息进行监测。

在诸多的基于数据的方法中，多变量统计过程监测(Multivariate StatisticalProcess Monitoring，MSPM)是最常见的故障检测方法。其中，又以主元分析(PrincipalComponent Analysis，PCA)与独立元分析(Independent Component Analysis，ICA)为最主流的数据分析与建模方法，已被广泛用于工业过程监测中，各种延伸拓展形式层出不穷。PCA模型过程监测时，一般需要假设过程数据服从或近似服从高斯分布以方便地确定出监测统计量的控制上限。若不满足高斯分布假设，PCA故障检测模型就无法较好地描述正常数据的波动范围，相应的故障检测能力将会大幅下降。相比较之下，ICA方法不存在这个假设，能更好地适应于非高斯过程对象的建模与监测问题。值得一提的是，无论是PCA还是ICA，对应的故障检测模型通常都是采用基于距离的监测统计量来定义正常数据范围。从几何空间角度看，PCA与ICA在将原始数据进行转换后，定义了一个超椭球来限定其波动范围。可想而知，若是这些经PCA或ICA转换过后成分信息不服从高斯分布，那么定义的超椭球内部就存在较多的“空洞”，处于这些空洞中的点按理说应该是故障数据，但是PCA或ICA却无法将其有效识别出来。

由于实际过程对象的复杂与大规模化特性，采样数据很少符合或近似符合高斯分布。因此，基于距离的统计监测指标及其对应的控制上限无法精准地定义训练数据的正常波动范围。一种可行的解决方法是采用多维核密度估计法或单分类的支持向量机，对数据的分布情况加以描述，这样一来就可以准确地定义的正常范围。然而，参数选择确是这类方法面临的主要问题。参数选择不当，很容易出现过拟合现象。另一种可行的解决思路可以是将原始不服从高斯分布的数据通过某种方式转换成误差。与基于机理模型的故障检测方法类似，通过监测实际值与估计值之间的误差来实现对生产过程状态的有效监测。一般而言，实际值与估计值之间的误差通常是符合或近似符合高斯分布的。若能将原始数据转换成这种的类型的误差，那么基于距离的统计监测指标以及对应的上限就不会存在“空洞”，相应的故障检测模型的可靠性与有效性又可以得到保证。

发明内容

本发明所要解决的主要技术问题是：如何从数据的角度将原始数据转换成服从或近似服从高斯分布的误差信息，然后对该误差实施监测从而避免因非高斯数据造成的种种不便。为此，本发明提供一种基于缺失数据的分布式工业过程监测方法。首先，本发明方法逐一假设过程对象中每个变量的测量数据缺失；然后，根据PCA模型中处理缺失数据的技巧推测出相应缺失变量的估计值；最后，利用假设的缺失数据实测值与估计值之间的误差作为被监测对象，再次建立基于PCA的故障检测模型实施在线故障检测。

本发明解决上述技术问题所采用的技术方案为：一种基于误差主元分析模型的工业过程监测方法，包括以下步骤：

(1)收集生产过程正常运行状态下的数据样本，组成训练数据集X∈R^n×m，并对每个变量进行标准化处理，得到均值为0，标准差为1的新数据矩阵其中，n为训练样本数，m为过程测量变量数，R为实数集，R^n×m表示n×m维的实数矩阵。

(2)利用PCA算法对矩阵进行分解其中T∈R^n×m与P∈R^m×m分别为主元得分矩阵与载荷矩阵，并初始化变量下标号i＝1。

(3)假设矩阵中第i列数据缺失，为不失一般性，可将矩阵表述成其中，为假设缺失的数据(实为矩阵中第i列)，由矩阵中剩余的列组成，对应的载荷矩阵P也可以类似的表述成P＝[P_i ^＃T，P_i ^*T]^T，其中，为矩阵P中的第i行，P_i ^*∈R^(m-1)×d由矩阵P中除去第i行以外的元素组成，上标号T表示矩阵或向量的转置。

(4)按照如下所示公式推测出中的第i列数据的估计值

(5)计算缺失数据的估计误差并判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(3)～(5)；若否，将所有得到的估计误差向量组成误差矩阵F＝[F₁，F₂，…，F_m]∈R^n×m后执行步骤(6)。

(6)对误差矩阵F建立基于PCA的故障检测模型，保留模型载荷矩阵H∈R^m×d、主元得分矩阵G∈R^n×d的协方差矩阵Λ∈R^d×d、以及监测统计量的控制上限T_c ²与Q_c，其中，d为主元模型保留的主元个数。

(7)收集新采样时刻的数据样本x∈R^m×1，对其实施与步骤(1)中相同的标准化处理得到新数据向量后，初始化i＝1。

(8)假设向量中第i个数据缺失，同理，可表示成其中，x_i ^＃为第i个缺失的数据，x_i ^*由向量中除缺失数据以外的元素组成，按照如下所示公式计算中缺失变量实测值与估计值之间的误差e_i：

e_i＝x_i ^#-P^#ΛP^*T(P^*ΛP^*T)^-1x_i ^* (2)

(9)判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(8)～(9)；若否，将所有得到的估计误差组成误差向量e＝[e₁，e₂，…，e_m]^T后执行下一步骤(10)。

(10)调用步骤(6)中保留的故障检测模型参数，并根据如下所示公式分别计算监测统计量指标T²与Q的具体数值：

T²＝e^THΛH^Te (3)

Q＝e^T(I-HH^T)e (4)

(11)将T²与Q的具体数值与对应控制上限T_c ²与Q_c进行比较，判断当前监测样本是否来自于故障工况。

与传统方法相比，本发明方法的优势在于：

本发明方法以PCA模型为基础通过假设缺失数据以生成缺失数据实测值与估计值之间的估计误差，巧妙的将原始数据转换成服从或近似服从高斯分布的误差。本发明方法不拘泥于原始数据是否满足于高斯分布假设，而通过对估计误差实施基于PCA的故障检测流程来避免直接面对这一问题。本发明方法虽然同样采用的是传统PCA故障检测模型中的距离型监测统计量指标，但是由于监测对象为满足于高斯分布的估计误差，本发明方法定义的正常状态的波动区域是比较准确的，不会出现“空洞”现象。可以说，本发明方法能较大幅度改善传统PCA方法的故障监测性能，是一种更为优选的过程监测方法。

附图说明

图1为本发明方法中基于缺失数据的误差生成示意图。

图2为本发明方法的实施流程图。

图3为本发明方法能将非高斯原始数据转变成高斯分布的误差信息示意图。

图4为本发明方法与其他方法在监测TE过程故障5时监测详情对比图。

具体实施方式

下面结合附图对本发明方法进行详细的说明。

如图1所示，本发明方法提供一种基于误差主元分析模型的工业过程监测方法，具体包括以下步骤：

首先是离线建模部分，如图2中左图所示，实施流程如下所示：

步骤(1)：收集生产过程正常运行状态下的数据样本，组成训练数据集X∈R^n×m，并对每个变量进行标准化处理，得到均值为0，标准差为1的新数据矩阵其中，n为训练样本数，m为过程测量变量数，R为实数集，R^n×m表示n×m维的实数矩阵。

步骤(2)：利用PCA算法对矩阵进行分解其中T∈R^n×m与P∈R^m×m分别为主元得分矩阵与载荷矩阵，并初始化变量下标号i＝1，具体实施过程如下所示：

①计算的相关性矩阵

②求解C所有特征值λ₁≥λ₂≥…≥λ_m所对应的特征向量p₁，p₂…，p_m；

③将m个特征向量组成载荷矩阵P＝[p₁，p₂…，p_m]；

④根据公式计算主元矩阵T∈R^n×m。

步骤(3)：假设矩阵中第i列数据缺失，为不失一般性，可将矩阵表述成其中，为假设缺失的数据(实为矩阵中第i列)，由矩阵中剩余的列组成，对应的载荷矩阵P也可以类似的表述成P＝[P_i ^＃T，P_i ^*T]^T，其中，P_i ^#∈R^{1 ×d}为矩阵P中的第i行，P_i ^*∈R^(m-1)×d由矩阵P中除去第i行以外的元素组成，上标号T表示矩阵或向量的转置。

步骤(4)：按照如下所示公式推测出中的第i列数据的估计值：

步骤(5)：计算缺失数据的估计误差并判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(3)～(5)；若否，将所有得到的估计误差向量组成误差矩阵F＝[F₁，F₂，…，F_m]∈R^n×m后执行步骤(6)。

步骤(6)：对误差矩阵F建立基于PCA的故障检测模型，保留模型载荷矩阵H∈R^m×d、主元得分矩阵G∈R^n×d的协方差矩阵Λ∈R^d×d、以及监测统计量的控制上限T_c ²与Q_c。为F建立基于PCA的故障检测模型包括以下步骤：

①计算F的协方差矩阵S＝F^TF/(n-1)；

②求解S所有特征值γ₁≥γ₂≥…≥γ_m所对应的特征向量h₁，h₂…，h_m；

③设置保留的主成分个数d为满足如下所示条件的最小值，并将对应的d个特征向量组成载荷矩阵H＝[h₁，h₂…，h_d]；

④根据公式G＝FH计算主元得分矩阵G∈R^n×d，并计算其相应的协方差矩阵Λ＝G^TG/(n-1)；

⑤根据如下所示公式分别确定监测统计量T²与Q对应的控制上限T_c ²与Q_c：

上两式中，置信水平α＝99％，F_α(d，n-d)表示自由度为d与n-d的F分布，表示权重为g＝v/2b，自由度为h＝2a²/b的χ²分布，a与b分别是F对应的Q统计量的估计均值和估计方差。

其次，实施在线故障检测的流程如图2中右图所示，详细的过程如下所示：

步骤(7)：收集新采样时刻的数据样本x∈R^m×1，对其实施与步骤(1)中相同的标准化处理得到新数据向量后，初始化i＝1。

步骤(8)：假设向量中第i个数据缺失，同理，可表示成其中，x_i ^＃为第i个缺失的数据，x_i ^*由向量中除缺失数据以外的元素组成，按照如下所示公式计算中缺失变量实测值与估计值之间的误差e_i：

e_i＝x_i ^#-P^#ΛP^*T(P^*ΛP^*T)^-1x_i ^* (9)

步骤(9)：判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(8)～(9)：若否，将所有得到的估计误差组成误差向量e＝[e₁，e₂，…，e_m]^T后执行下一步骤(10)。

步骤(10)：用步骤(6)中保留的故障检测模型参数，并根据如下所示公式分别计算监测统计量指标T²与Q的具体数值：

T²＝e^THΛH^Te (10)

Q＝e^T(I-HH^T)e (11)

步骤(11)：判断T²与Q的具体数值是否大于对应控制上限T_c ²与Q_c？若否，则当前样本为正常工况采样；若是，则当前采样数据有可能来自故障工况，理应继续监测接下来的3～6个新样本，若都超限，则说明当前工况已出现故障，若都没超限，则说明当前工况仍旧处于正常状态。

下面结合一个具体的工业过程的例子来说明本发明方法相对于现有方法的优越性与可靠性。该过程数据来自于美国田纳西-伊斯曼(TE)化工过程实验，原型是伊斯曼化工生产车间的一个实际工艺流程。目前，TE过程因其流程的复杂性，已作为一个标准实验平台被广泛用于故障检测研究。整个TE过程包括22个测量变量、12个操作变量、和19个成分测量变量。所采集的数据分为22组，其中包括1组正常工况下的数据集与21组故障数据。而在这些故障数据中，有16个是已知故障类型，如冷却水入口温度或进料成分的变化、阀门粘滞、反应动力学漂移等，还有5个故障类型是未知的。为了对该过程进行监测，选取如表1所示的33个过程变量，接下来结合该TE过程对本发明具体实施步骤进行详细的阐述。

表1：TE过程监测变量。

序号	变量描述	序号	变量描述	序号	变量描述
						1	物料A流量	12	分离器液位	23	D进料阀门位置
2	物料D流量	13	分离器压力	24	E进料阀门位置
						3	物料E流量	14	分离器塔底流量	25	A进料阀门位置
4	总进料流量	15	汽提塔等级	26	A和C进料阀门位置
						5	循环流量	16	汽提塔压力	27	压缩机循环阀门位置
6	反应器进料	17	汽提塔底部流量	28	排空阀门位置
						7	反应器压力	18	汽提塔温度	29	分离器液相阀门位置
8	反应器等级	19	汽提塔上部蒸汽	30	汽提塔液相阀门位置
						9	反应器温度	20	压缩机功率	31	汽提塔蒸汽阀门位置
10	排空速率	21	反应器冷却水出口温度	32	反应器冷凝水流量
						11	分离器温度	22	分离器冷却水出口温度	33	冷凝器冷却水流量

首先，利用TE过程正常工况下的采样数据建立故障检测模型，包括以下步骤：

(1)采集正常工况下的过程数据，以其中960个正常数据组成训练数据矩阵X∈R^{960 ×33}，对其进行标准化处理得到

(2)利用PCA算法对矩阵进行分解并初始化变量下角标i＝1。

(3)假设矩阵中第i列数据缺失，为不失一般性，可将矩阵表述成对应的载荷矩阵P也可以类似的表述成P＝[P_i ^＃T，P_i ^*T]^T。

(4)按照如下所示公式推测出中的第i列数据的估计值

(5)计算缺失数据的估计误差并判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(3)～(5)；若否，将所有得到的估计误差向量组成误差矩阵F＝[F₁，F₂，…，F_m]∈R^n×m后执行步骤(6)。

(6)对误差矩阵F建立基于PCA的故障检测模型，保留模型载荷矩阵H∈R^m×d、主元得分矩阵G∈R^n×d的协方差矩阵Λ∈R^d×d、以及监测统计量的控制上限T_c ²与Q_c。

为验证本发明方法可以将非高斯的原始数据转变成高斯分布的误差信息，特对TE过程第10、19、和31号测量变量的采样数据执行高斯分布检验，相应的检验图显示于图3中的上半部分。可以很明显的发现这三个测量变量不服从高斯分布。因为，高斯分布检验图为一条直线时，所测试的数据向量才服从高斯分布。在经本发明方法转变成误差后，以主元估计误差F₁₀，F₁₉，F₃₁为例，相应的高斯分布检验图如图3中下半部分所示，都呈现出近乎直线的特性。因此，本发明方法在将其转换成误差信号后，该误差是服从高斯分布的。

最后，实施在线故障检测，为对比验证本发明方法与其他方法(如传统PCA方法与Modified ICA方法)的优越性，采集TE过程中第五种故障类型工况下的960个采样数据，其中前160个数据是正常的，后800个数据才处于异常工况。

(7)收集新采样时刻的数据样本x∈R^1×33，对其实施与步骤(1)中相同的标准化处理得到新数据向量后，初始化i＝1。

(8)假设向量中第i个数据缺失，同理，可表示成计算中缺失变量实测值与估计值之间的误差e_i。

(9)判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(8)～(9)；若否，将所有得到的估计误差组成误差向量e＝[e₁，e₂，…，e₃₃]^T后执行下一步骤(10)。

(10)计算监测误差向量e的统计量指标T²与Q的具体数值。

(11)判断T²与Q的具体数值是否大于对应控制上限T_c ²与Q_c？若否，则当前样本为正常工况采样；若是，则当前采样数据有可能来自故障工况，理应继续监测接下来的3～6个新样本，若都超限，则说明当前工况已出现故障，若都没超限，则说明当前工况仍旧处于正常状态。对应的故障监测详图显示于图4中，可以很明显的发现，本发明方法取得的效果明显优越于其他两种方法。

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

Claims (3)

1.一种基于误差主元分析模型的工业过程监测方法，其特征在于，包括以下步骤：

(1)收集生产过程正常运行状态下的数据样本，组成训练数据集X∈R^n×m，并对每个变量进行标准化处理，得到均值为0，标准差为1的新数据矩阵其中，n为训练样本数，m为过程测量变量数，R为实数集，R^n×m表示n×m维的实数矩阵；

(2)利用主元分析(PCA)算法对矩阵进行分解其中T∈R^n×m与P∈R^m×m分别为主元得分矩阵与载荷矩阵，并初始化变量下标号i＝1；

(3)假设矩阵中第i列数据缺失，为不失一般性，可将矩阵表述成其中，为假设缺失的数据(实为矩阵中第i列)，由矩阵中剩余的列组成，对应的载荷矩阵P也可以类似的表述成P＝[P_i ^#T，P_i ^*T]^T，其中，P_i ^#∈R^1×d为矩阵P中的第i行，P_i ^*∈R^(m-1)×d由矩阵P中除去第i行以外的元素组成，上标号T表示矩阵或向量的转置；

(4)按照如下所示公式推测出中的第i列数据的估计值

(5)计算缺失数据的估计误差并判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(3)～(5)；若否，将所有得到的估计误差向量组成误差矩阵F＝[F₁，F₂，…，F_m]∈R^n×m后执行步骤(6)；

(6)对误差矩阵F建立基于PCA的故障检测模型，保留模型载荷矩阵H∈R^m×d、主元得分矩阵G∈R^n×d的协方差矩阵Λ∈R^d×d、以及监测统计量的控制上限与Q_c，其中，d为主元模型保留的主元个数；

(7)收集新采样时刻的数据样本x∈R^m×1，对其实施与步骤(1)中相同的标准化处理得到新数据向量后，初始化i＝1；

(8)假设向量中第i个数据缺失，同理，可表示成其中，x_i ^#为第i个缺失的数据，x_i ^*由向量中除缺失数据以外的元素组成，按照如下所示公式计算中缺失变量实测值与估计值之间的误差e_i：

(9)判断是否满足条件i＜m？若是，置i＝i+1后，重复执行步骤(8)～(9)；若否，将所有得到的估计误差组成误差向量e＝[e₁，e₂，…，e_m]^T后执行下一步骤(10)；

(10)调用步骤(6)中保留的故障检测模型参数，并根据如下所示公式分别计算监测统计量指标T²与Q的具体数值：

T²＝e^THΛH^Te (3)

Q＝e^T(I-HH^T)e (4)

(11)判断T²与Q的具体数值是否大于对应控制上限与Q_c？若否，则当前样本为正常工况采样；若是，则当前采样数据有可能来自故障工况，理应继续监测接下来的3～6个新样本，若都超限，则说明当前工况已出现故障，若都没超限，则说明当前工况仍旧处于正常状态。

2.根据权利要求1所述的一种基于缺失数据的分布式工业过程监测方法，其特征在于，所述步骤(2)的具体实施过程如下所示：

①计算的相关性矩阵

②求解C所有特征值λ₁≥λ₂≥…≥λ_m所对应的特征向量p₁，p₂…，p_m；

③将m个特征向量组成载荷矩阵P＝[p₁，p₂…，p_m]；

④根据公式计算主元矩阵T∈R^n×m。

3.根据权利要求1所述的一种基于缺失数据的分布式工业过程监测方法，其特征在于，所述步骤(6)中建立基于PCA的故障检测模型的具体实施过程如下所示：

①计算F的协方差矩阵S＝F^TF/(n-1)；

②求解S所有特征值γ₁≥γ₂≥…≥γ_m所对应的特征向量h₁，h₂…，h_m；

③设置保留的主成分个数d为满足如下所示条件的最小值，并将对应的d个特征向量组成载荷矩阵H＝[h₁，h₂…，h_d]；

④根据公式G＝FH计算主元得分矩阵G∈R^n×d，并计算其相应的协方差矩阵Λ＝G^TG/(n-1)；

⑤根据如下所示公式分别确定监测统计量T²与Q对应的控制上限T_c ²与Q_c：

上两式中，置信水平α＝99％，F_α(d，n-d)表示自由度为d与n-d的F分布，表示权重为g＝v/2b，自由度为h＝2a²/b的χ²分布，a与b分别是F对应的Q统计量的估计均值和估计方差。