1 前言

安全生产与产品质量在几乎所有的工业过程中都有着举足轻重的地位，对工业生产过程进行实时监测是保证生产安全与维持产品质量稳定性的最直接途径。无论是工业界还是学术界，对过程监测方法的研究一直伴随着整个工业发展的进程。从最开始的基于机理模型的故障检测方法，到目前逐步发展成数据驱动的工业过程监测技术，几十年间的研究与应用成果层出不穷^[1-6]。基于机理模型的故障检测方案主要依赖于生成某些过程参数实际值与估计值之间的误差，并通过误差的变化情况来反映过程对象是否出现异常状态^[1-2]。因此，基于机理模型的故障检测方法需要被监测对象精确的物理模型，限制了该类方法技术在现代大规模工业过程中的应用范围^[7]。值得指出的是，基于机理模型的故障检测方法核心在于生成的误差^[8]，并以误差为被监测对象是值得借鉴学习的。数据驱动的过程监测方法体系旨在通过采样数据建立描述正常数据变化特征的数学模型，任何有悖于该数学模型的采样数据都被定义为非正常工况下的数据样本^[9]。现代工业过程广泛采用计算机与先进的仪表技术，海量的采样数据为数据驱动的工业过程研究提供了坚实的数据基础。

在数据驱动的过程监测领域，研究与应用最多的是多变量统计分析算法，其中以主成分分析(principal component analysis, PCA)与独立成分分析(independent component analysis, ICA)两种算法为主^{[5, 9]}。PCA与ICA用于故障检测的相同之处在于将高维原始数据投影变换成低维的潜藏有用成分信息，然后再对该潜藏信息与模型拟合残差建立基于距离的统计量监测指标实施在线故障检测。两者的不同之处在于提取潜藏信息的方式与出发点不一样，PCA在于挖掘不相关的主成分信息并最大化地保留方差，而ICA在高阶统计量的指引下提取相互独立的独立成分信息，挖掘的是潜藏的非高斯信息，所以可以应对非高斯的过程监测问题^[10]。然而，无论是PCA还是ICA算法，其基本思想都是对数据实施投影变换，将原始数据转换成少数几个特征成分。从某种程度上讲，这些投影变换方向向量上的系数各不相同，会对不同变量赋予了不同的权重值。换句话说，PCA或ICA模型中所有变量因投影变换向量的转换并非得到了同等程度的对待^[11]。

由于基于PCA或是ICA的过程监测方法只需利用正常工况下的采样数据进行离线建模，所有测量变量在故障检测中的重要性程度是相同的。因此，若是不等同对待某些变量，相应的故障检测性能就无法得到充分保证。借鉴分散式建模的思路以及误差生成的策略，本文提出一种基于缺失变量估计误差的工业过程监测方法。一般来讲，缺失变量有可能会出现在实际采样过程中，导致统计模型面临如何进行模型预测的问题^[12]。现已有很多种可以处理包含缺失变量数据的模型建立与模型预测问题^[13]。以PCA算法为例，FOLCH- FORTUNY等人详细比较了多种处理缺失变量的方法。本文借鉴的缺失变量处理方法正是由FOLCH-FORTUNY等^[12]提出的，其优势在于未假设训练数据服从高斯分布，直接通过迭代算法(iterative algorithm, IA)获取缺失变量的估计值。值得强调的是，实际的变量缺失问题并不是本文研究的关注点，只是单纯地假设测量变量数据缺失，从而提出一种新颖的建模与故障检测方法。具体来讲，所提出方法是在已有PCA模型的基础上，通过逐一假设采样数据各个变量的测量值缺失，再利用IA计算出相应缺失变量的估计值。根据缺失变量的估计值，原始数据就可以一一对应地转变成实测值与估计值之间的误差了。本文逐一假设各个测量变量缺失的主要动机包括：首先，每个测量变量在故障检测上的同等重要性，逐一假设可保证公允性；其次，只在一个变量缺失的假设前提下，缺失变量的估计值更接近于实际值，估计精度更高。可以说，这种逐一假设各变量缺失的方法，不仅等同地对待了所有的过程测量变量，而且通过生成的误差就能反映出在缺失该变量的条件下对特征成分的影响程度。

2 经典过程监测方法 2.1 基于PCA的过程监测方法

PCA算法旨在将经标准化后的训练数据X∈R^n×m(n为样本数，m为变量数)进行分解：

$ \mathit{\boldsymbol{X}} = \mathit{\boldsymbol{T}}{\mathit{\boldsymbol{P}}^{\rm{T}}} + \mathit{\boldsymbol{E}} $

(1)

从而得到得分矩阵T=XP∈R^n×d、载荷矩阵P∈R^m×d，以及模型残差E∈R^n×m。其中，参数d为PCA模型中保留的主成分个数，一般可采用累计方差贡献率(cumulative percentage variance, CPV)或交叉验证(cross validation, CV)原则确定。可以发现，PCA旨在将原始数据通过载荷矩阵P投影变换后，分解成两个子空间：主元子空间与残差子空间。载荷矩阵P中各列向量即是投影变换的各个方向，这些向量各元素的差异体现了各测量变量的权重大小。因此，经PCA投影变换后，原始数据各测量变量被不同程度地加权了。

实施在线监测时，通常采用的T²与Q统计量分别监测主元子空间与残差子空间，即：

$ {T^\mathit{\boldsymbol{2}}} = \mathit{\boldsymbol{xP \boldsymbol{\varLambda} }}{\mathit{\boldsymbol{P}}^{\rm{T}}}{\mathit{\boldsymbol{x}}^{\rm{T}}} \le \frac{{d\left( {n - 1} \right)}}{{n - d}}{F_\alpha }\left( {d, n - d} \right) $

(2)

$ Q = \mathit{\boldsymbol{x}}\left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{P}}{\mathit{\boldsymbol{P}}^{\rm{T}}}} \right){\mathit{\boldsymbol{x}}^T} \le g\chi _{h, \alpha }^2 $

(3)

分别用来监测新样本数据x∈R^1×m经PCA投影变换后的主成分与残差的变化情况。式(2)和(3)中，α为置信限，Λ = T^TT/(n-1)，F_α(d, n-d)表示自由度为d与n-d的F分布，gχ_h,α²表示权重为g=v/2b，自由度为h=2a²/b的χ²分布，a与b分别是训练数据集对应的Q统计量的估计均值和估计方差^[14]。

2.2 基于ICA的过程监测方法

ICA算法起源于盲源信号分离，旨在从数据集X中分理处出k < m个相互独立的成分，即：

$ \mathit{\boldsymbol{S}} = \mathit{\boldsymbol{W}}{\mathit{\boldsymbol{X}}^{\rm{T}}} $

(4)

式中，S∈R^k×n为独立成分矩阵；W∈R^k×m为分离矩阵，通常使用Fast ICA算法迭代求取。然而，传统Fast ICA算法针对不同初始值最后的收敛结果不一致。为此，LEE等^[15]提出一种新的Modified ICA算法并将之应用于过程监测。该算法首先利用PCA对原始数据进行白化处理，并假设经白化处理后的数据可作为一个“较好”独立成分初值。然后，经过迭代过程求取ICA模型的混合矩阵A∈R^m×k与分离矩阵W∈R^k×m，那么，相应的ICA模型可表述成：

$ {\mathit{\boldsymbol{X}}^{\rm{T}}} = \mathit{\boldsymbol{AS}} + \mathit{\boldsymbol{F}} $

(5)

其中，F为模型残差；S包含k个独立成分且都是非高斯分布的，这是与PCA挖掘特征最为显著的不同之处。此外，监测统计量控制上限不能按照高斯分假设来确定，只能依据单变量的核密度估计法确定^[15]。

2.3 分散式过程监测方法

分散式的过程监测旨在将过程对象分解成多个子块并建立多个模型，这不仅能降低多过程分析的复杂程度，而且多模型的效果通常优越于单个模型^[16]。因此，分散式的过程监测方法于近年来得到了较多的青睐^{[17, 18]}。Ge等^[17]曾利用PCA模型中各载荷向量对过程变量实施分解，累计得到d+1个子变量块，并对各个子块分别建立基于PCA的过程监测模型。这种Distributed PCA方法不需要过程的机理知识，而完全是从数据统计角度出发实施过程分解的。通过贝叶斯概率融合的方式，可将d+1个PCA模型的T²与Q监测指标分别融合成概率型指标BIC_T2与BIC_Q以方便最终故障决策。

类似地，在实施变量分解后同样可以建立多个ICA过程监测模型，从而实施分布式的非高斯过程监测。例如，LI等^[18]提出将测量变量分解成相关与弱相关的两个子块，然后分别建立基于ICA算法的过程监测模型。相比于单个基于ICA的非高斯过程监测模型，这种基于two-block ICA的分散式过程监测模型能够展现出更优越的故障检测性能。

3 基于误差的工业过程监测方法

如图 1所示，在已知PCA模型的基础上提出一种基于IA缺失变量(missing variable, MV)估计误差的过程监测方法，简称为IAMV方法。从图 1中可以发现，本文所提出的IAMV方法将分散建模发挥得更加彻底，即为每个测量变量建立其相应的估计误差生成模型。值得强调的是，图 1中所示建模过程并未考虑真实意义上的变量缺失问题，只是在逐个假设各测量变量数据缺失的前提下使用IA方法得出缺失变量的估计值。

3.1 基于IA的误差生成方式

为不失一般性，采样数据向量x可描述成：

$ \mathit{\boldsymbol{x = }}\left[ {\mathit{\boldsymbol{x}}_i^\# , \mathit{\boldsymbol{x}}_i^*} \right] $

(6)

其中，i = 1, 2, …, m; x_i^#∈R^1×1为向量x中假设的第i个缺失变量，x_i^*为向量x中剩余的可测量数据。载荷矩阵P可表述成对应的形式，即：

$ \mathit{\boldsymbol{P}} = \left[ \begin{array}{l} \mathit{\boldsymbol{P}}_i^\# \\ \mathit{\boldsymbol{P}}_i^* \end{array} \right] $

(7)

式中，P_i^#∈R^1×d与P_i^*∈R^(m-1)×d分别对应x中缺失变量与可观测变量。利用IA生成缺失变量x_i^#的估计值$\hat x_i^\# $的流程如图 2所示，具体为：

(1) 初始化$\hat x_i^\# $=0；

(2) 由$ \mathit{\boldsymbol{\hat x = }}\left[ {\mathit{\boldsymbol{\hat x}}_i^\# , \mathit{\boldsymbol{x}}_i^*} \right]$计算出相应主成分向量$\mathit{\boldsymbol{\hat t}} = \mathit{\boldsymbol{\hat xP}} = \mathit{\boldsymbol{\hat x}}_i^\# \mathit{\boldsymbol{P}}_i^\# + \mathit{\boldsymbol{x}}_i^*\mathit{\boldsymbol{P}}_i^* $；

(3) 更新缺失变量的估计值$\hat x_i^\# = \mathit{\boldsymbol{\hat tP}}_i^{\# {\rm{T}}} $；

(4) 判断$ \hat x_i^\# $是否收敛？若否，则返回(2)继续迭代操作；若是，则得到缺失变量估计值$ \hat x_i^\# $；

(5) 计算估计误差${e_i} = x_i^\# - \hat x_i^\# $。

从上述迭代步骤中可以看出，IA方法不假设原始数据服从高斯分布，这是与文献中经常使用的期望最大化(expectation maximization, EM)法估计缺失数据的不同之处。由于本文对采样数据不做任何分布假设，IA方法更适合用来估计缺失变量。

3.2 基于IAMV的过程监测方法

通过逐一假设各个变量数据缺失，可以将训练数据X中每一列都对应转换成误差Θ∈R^n×m中的每一列。由于这些误差都是通过同一个已知的PCA模型转换来的，存在较多的冗余信息。因此，可对误差Θ实施奇异值分解，即：

$ \Theta = \mathit{\boldsymbol{UD}}{\mathit{\boldsymbol{V}}^{\rm{T}}} $

(8)

其中，对角矩阵D中包含了l个非零奇异值。那么，从误差Θ到U的变换过程为：

$ \mathit{\boldsymbol{U}} = \Theta \mathit{\boldsymbol{V}}{\mathit{\boldsymbol{D}}^{ - 1}} = \Theta \mathit{\boldsymbol{H}} $

(9)

式中，变换矩阵H=VD^-1。

至此，整个基于IAMV的过程监测方法的实施流程可以归纳为离线建模和在线故障检测流程。

离线建模流程：

1) 对正常工况下的训练数据X∈R^n×m实施标准化处理，以消除变量间量纲差异的影响；

2) 利用PCA算法对标准化后的数据X建立PCA模型X=TP^T+E；

3) 依次假设X中第i个变量的测量数据缺失，并根据3.1节介绍的基于IA的误差生成方法，将X中第i列数据变换成估计误差向量Θ_i；

4) 对生成的误差矩阵Θ=[Θ₁, Θ₂, ..., Θ_m]实施奇异值分解，从而得到变换矩阵H；

在线故障检测流程：

5) 对新采样数据实施与步骤1)中一样的标准化处理；

6) 依次假设x中第i个数据缺失，并利用IA推测出相应的估计误差e_i；

7) 在将x中全部元素变换成估计误差e=[e₁, e₂, ..., e_m]后，计算监测统计量Q，即：

$ Q\mathit{\boldsymbol{ = eH}}{\mathit{\boldsymbol{H}}^{\rm{T}}}{\mathit{\boldsymbol{e}}^{\rm{T}}} $

(10)

8) 监测统计量Q的控制上限为自由度是m的χ²分布在置信水平α下的取值(即χ_α,1²)，决策当前采样数据是否来自异常工况。

4 仿真实验研究

TE仿真模型因其结构的复杂性，已成为不同控制方法和过程监测策略的标准实验平台^[19-20]。TE过程主要有连续搅拌反应器，产品冷凝器、气液分离塔、汽提塔和离心式压缩机等5个生产单元组成，可连续测量22个过程变量和12个操作变量，还可以仿真模拟如表 1所列18种不同的故障。选取如参考文献[15]中所列的33个变量作为监测变量。离线建模阶段，利用正常工况下960个样本分别建立基于PCA、Distributed PCA^[17]、Modified ICA、Two-block ICA^[18]、以及IAMV的故障检测模型。其中，置信限α统一取值99%，PCA模型中参数d采用CPV > 85%的准则确定，Distributed PCA通过贝叶斯融合得到两个概率型监测指标BIC_T²与BIC_Q。非高斯过程监测模型(Modified ICA与Two-block ICA)的模型参数分别与参考文献[15]与[18]一致。

对比4种不同故障检测方法(PCA、Distributed PCA、Modified ICA和IAMV)用于监测TE过程不同故障工况时的故障检测效果，将故障漏报率详情列于表 1中。其中故障3、9和15由于对过程数据的影响甚微，已被很多文献证实它们很难被可靠地检测出^[21]，在本文的研究中未加考虑。针对TE过程剩余18种异常工况的监测中，取得最小漏报率的数值已用黑体标出。可以发现，IAMV方法在绝大多数故障的监测上能取得最优的监测效果。相比于PCA、Distributed PCA以及非高斯的Modified ICA与Two-block ICA方法，IAMV方法能显著改善针对TE过程故障5、10、16、19、和20的故障检测效果，相应的故障检测漏报率得到了较大幅度的降低。另外，Distributed PCA方法与Two-block ICA方法由于占据分散式多模型的优势，过程监测效果分别明显优越于传统的PCA方法与Modified ICA方法。虽然，在故障2、8、17和18上的检测上，IAMV方法不是最佳的，但是故障漏报率差别仅在2%左右。最后，将故障20的过程监测详情展示于图 3中。从图 3的对比中，可以发现IAMV方法能持续不断地触发故障警报，而PCA与Modified ICA都存在较明显的漏报情况。从以上的对比实验中，可以总结出本文所提出的IAMV方法由于等同对待各测量变量，且利用了分散式建模的优势，充分地保证了过程监测性能。

5 结论

在假设变量缺失的前提下，提出了一种基于缺失变量估计误差的故障检测方法。该方法借助可处理缺失变量的IA方法，通过逐一假设采样数据中各个变量的测量值缺失，通过PCA模型计算出每个缺失变量对应的估计值，从而将原始数据转变成估计误差。一定程度上讲，利用了分散式建模策略等同地对待每个测量变量。最后，通过对比实验验证了IAMV方法相比于其他方法的优越性与可靠性，是一种更为优选的过程监测实施方案。然而，当前的研究还未涉及非线性模型，也未能据此进一步提出相应的故障诊断方法。在未来的研究中，如何将这种缺失变量生成误差的思路应用在非线性模型上，以及如何实施故障分类诊断，都是值得深入研究的。