CN109085816A - 一种基于动态正交的全局局部间歇过程故障检测方法 - Google Patents

Abstract

本发明提供了一种基于动态正交的全局局部间歇过程故障检测方法，包括以下步骤：(1)收集间歇过程正常运行下的各关键变量数据，形成正常运行状态下的训练样本X∈R^I×J×K；(2)将训练样本X首先沿着批次方向展开为二维数据X∈R^I×KJ，然后将展开的二维数据进行标准化，将标准化后的二维数据重新排列成X∈R^KI×J；(3)在二维数据X∈R^KI×J基础上，建立一个时滞矩阵X_D来消除过程变量的时序自相关；(4)对建立的时滞矩阵X_D构造动态正交全局局部模型；(5)分别在特征空间和残差空间建立T²和SPE统计模型，求取控制限；(6)收集在线的过程数据，进行标准化处理；(7)利用已经建立的动态正交全局局部模型将在线数据进行投影，通过T²和SPE统计模型，判断故障的发生。

Description

一种基于动态正交的全局局部间歇过程故障检测方法

技术领域

本发明属于工业过程监控技术领域，涉及一种基于动态正交的全局局部间歇过程故障检测方法。

背景技术

随着现代工业过程的规模不断扩大以及复杂度的增加，间歇过程作为一种重要的生产方式广泛的运用于高质量和小容量的过程生产中，如特殊化工、食品生产、医药生产和半导体生产等。间歇过程生产是由许多生产环节组成，每一环节对最终的产品质量和安全生产都至关重要，并且前一环节的完成情况直接影响下一环节的实施。由于多个环节的嵌套、生产规模不断扩大以及复杂度的增加，这就使得间歇过程生产更容易发生故障。如果故障不能及时的检测和排除，对整个生产过程的影响会越来越大，以至于最终影响整个生产。因此，这就迫切需要对间歇过程产品质量和生产安全建立监控系统。

由于分布式控制系统的广泛运用和数据分析技术的快速发展，多元统计的方法逐渐成为过程监测的主流。多元统计的方法不受对象模型的限制，只需要对生产过程的数据进行分析。对于日渐复杂的间歇过程来说，对其很难建立数学模型，因此基于多元统计的方法已经被深入研究并广泛用于间歇过程的故障检测。其中，多向主元分析和多向偏最小二乘已经成功的用于间歇过程监控中，这两种方法将间歇过程数据降维到一个低维空间，然后在低维空间中建立统计模型进行分析。一些学者对以上算法进行了改进，在间歇过程监控中考虑了过程的动态性，如动态多向主元分析等。传统的方法在特征提取时，只考虑了过程数据的全局结构，忽略了数据的局部结构，因此在降维过程中不能充分的保持过程数据的全局和局部结构，并且在数据重构时没有考虑广义正交问题，使得在过程监控中造成部分有用信息的丢失，导致监控效果不佳。

发明内容

本发明的目的在于针对现有技术存在的问题，提供一种基于动态正交的全局局部间歇过程故障检测方法，该方法可显著提高间歇过程故障检测的效果，能在故障发生时及时检测到故障，有利于对产品质量和生产安全进行严格把控。

一种基于动态正交的全局局部间歇过程故障检测方法，包括以下步骤：

(1)利用集散控制系统收集间歇过程正常运行下的各关键变量数据，形成正常运行状态下的训练样本X∈R^I×J×K，其中R为实数集并满足(I×J×K)的三维分布，I表示批次，J表示变量，K为采样点；

(2)将训练样本X首先沿着批次方向展开为二维数据X∈R^I×KJ，然后将展开的二维数据进行标准化，使得各个变量的均值为0，方差为1，然后将标准化后的二维数据重新排列成X∈R^KI×J；

(3)在混合展开的二维数据X∈R^KI×J基础上，建立一个时滞矩阵X_D来消除过程变量的时序自相关；

(4)对建立的时滞矩阵X_D构造动态正交全局局部模型，且其过程是首先建立全局局部特征提取的目标函数J_G-L，然后引入广义正交向量以求取目标函数J_G-L的投影矩阵；

(5)分别在特征空间和残差空间建立T²和SPE统计模型，求取控制限；

(6)收集在线的过程数据，进行标准化处理；

(7)利用已经建立的动态正交全局局部模型将在线数据进行投影，通过T²和SPE统计模型，判断故障的发生。

进一步地，步骤(3)中：时滞矩阵X_D如下所示：

上式中，d为窗口的宽度，T为总的采样数，X(t)＝[x_1,t,x_2,t,…,x_J,t]^T为在t采样时刻的测量向量。

进一步地，步骤(4)中：建立全局局部特征提取的目标函数J_G-L的具体过程如下：

首先建立全局特征提取的目标函数J_global：

上式中，Y＝[y₁,y₂,.....,y_n]＝A^TX，A＝[a₁,a₂,…,a_d]为J_global的投影矩阵，X＝[x₁,x₂,…,x_n]为采样矩阵，为局部均值矩阵，为全局结构矩阵；

局部特征提取的目标函数J_local通过邻域保持嵌入算法得到，具体过程如下：

首先，利用k近邻算法，选择近邻点，如果x_i为x_j的近邻点，则邻域图中节点i到节点j通过一条边连接，反之，如果不是近邻点，则没有连接，通过这种方法构成邻域图，邻域图的权重系数W＝[w₁,w₂,...,w_n]通过下式求得：

∑_jw_ij＝1,j＝1,2,…n

则局部特征提取的目标函数J_local可以表示为：

上式中，I为单位向量，M＝(I-W)^T(I-W)，限定条件为a^TXX^Ta＝1，则全局局部特征提取的目标函数J_G-L如下所示：

则上式可以转化为求解广义特征值问题，即：

Ga＝λMa

设{a₁，a₂...a_d}为正交基向量，第一个向量a₁为(M)^-Ga对应的最小的特征值所对应的特征向量，当d＞2时，如下式：

上式中限定条件为：引入拉格朗日乘子，则可转化为下式：

对C^(d)求取偏导则为：

对上式乘以可得：

通过连续的乘以可以获得如下d‐1个式子：

上式可以表示为：

式中，μ^(d-1)＝[μ₁,…,μ_d-1]，A^(d-1)＝[a₁,…,a_d-1]，B^(d-1)＝[A^(d-1)]^T(XX^T)^-1A^(d-1)；

通过乘以(M)^-1，则上式可以转化为：

2(M)^-1Ga_d-2λa_d-(M)^-1A^(d-1)μ^(d-1)＝0

{I-(M)^-1A^(d-1)[B^(d-1)]^-1[A^(d-1)]^T}·(M)^-1Ga_d＝λa_d

上式中，I为单位向量，λ为最小化判别因子，a^d为最小的特征值所对应的特征向量；A＝{a₁,a₂,…a_d}为(M)^-1Ga_d通过迭代计算的最小特征值所对应的特征向量，即正交全局局部特征提取的目标函数J_G-L的投影矩阵，也即J_G-L的正交解。

进一步地，步骤(5)为：设Y(y₁,…,y_n)∈R^n×d为原始变量空间X(x₁,…,x_m)∈R^n×m的低维表示，y_i＝A^Tx_i，则特征空间和残差空间可以描述为：

Y＝A^TX＝(B^TB)^-1B^TX

E＝X-BY

上式中，A^T＝(B^TB)^-1B^T为转移矩阵，表示主元空间，表示残差空间，E＝[e₁，e₂,.....e_n]为残差矩阵，T²和SPE统计模型为：

T²＝y^TΛ^-1y

SPE＝e^Te

在求得正常数据的SPE和T²统计后，通过核密度估计算法求取控制限SPE_α和求取过程为：

对一个给定的数据集j，则j的密度分布为：

上式中，j_i为观测值，σ为窗口宽度，n为观测值的个数，K为核函数，因此通过选择高斯核函数设检验水平为0.95下的SPE和T²的控制限为SPE_α和则有:

可以求取控制限为SPE_α和

步骤(6)中，对在线采集的测试数据X_test进行标准化，使得各个变量的均值为0，方差为1。

步骤(7)中，通过已经求取的全局局部特征提取的目标函数J_G-L的投影矩阵A，将采样的测试数据X_test进行投影，则Y_test＝A^TX_test

Y_test＝A^TX_test＝(B^TB)^-1B^TX_test

E_test＝(I-BA^T)X_test＝(I-B(B^TB)^-1B^T)X_test

因此在线数据的T²和SPE统计可以描述为：

如果SPE_test＞SPE_α或者则说明发生故障，否则没有故障发生。

本发明的有益效果在于：本发明通过对间歇过程数据进行建模统计分析，在提取过程数据时，首先对过程数据建立时滞窗来消除动态时序相关性，消除过程数据的动态性；之后考虑数据的全局和局部结构特性，并在求取全局局部结构投影时引入正交矩阵来充分的保持数据的几何结构；最后基于该模型建立相应的统计监测量，判断故障的发生。相比于目前的其他方法，本发明可显著提高间歇过程故障检测率并减少故障误报率，能在故障发生时及时检测到故障，有利于对产品质量和生产安全进行严格把控。

附图说明

图1是MPCA法对青霉素发酵过程搅拌功率添加0.03斜坡故障的故障检测图；

图2是MNPE法对青霉素发酵过程搅拌功率添加0.03斜坡故障的故障检测图；

图3是MGNPE法对青霉素发酵过程搅拌功率添加0.03斜坡故障的故障检测图；

图4是本发明方法对青霉素发酵过程搅拌功率添加0.03斜坡故障的故障检测图。

具体实施方式

下面结合具体实施例对本发明方法进行进一步的说明。

青霉素生产过程是一个典型的动态、非线性、时变、多阶段的间歇过程。本发明通过Pensim2.0青霉素发酵过程的标准仿真平台产生出间歇过程数据，Pensim2.0是美国Illinois州立理工学院为了更加方便地研究典型间歇过程而开发的，它可产生出不同初始条件和不同工况下青霉素发酵过程中各变量每个时刻的数据用以分析研究。在青霉素发酵模型里，充分的考虑了发酵过程的温度变化、pH值、空气流量的变化、底物流加速率、搅拌率等对菌体合成的影响，可以接近真实地模拟青霉素发酵的实际过程。接下来结合具体过程对本发明的方法步骤作一阐述：

(1)本发明将每批次的反应时间设定为400h，每隔1h采样一次，在不同初始条件且没有故障的情况下共产生30个批次正常工况下数据。从产生的18个变量数据中选择其中10个过程变量作为监控变量(如表1所示)，构成三维数据矩阵X(30×10×400)并作为训练样本。

表1受监控的过程变量

变量序号	变量名称
		1	通风速率
2	搅拌速率
		3	底物流加速率
4	补料温度
		5	基质浓度
6	溶解氧浓度
		7	反应器体积
8	排气二氧化碳浓度
		9	pH值
10	发酵罐温度

(2)对得到的三维标准数据X(30×10×400)沿着批次方向展开为X_b(30×4000)，然后进行标准化，最后对标准化后的X_b(30×4000)重新排列成沿变量方向的二维数据X_b-v(12000×10)。

(3)对混合展开的二维数据X_b-v(12000×10)，建立一个时滞矩阵来消除过程变量的时序自相关，时滞矩阵如下所示：

上式中，d为窗口的宽度，T为总的采样数，X(t)＝[x_1,t,x_2,t,…,x_J,t]^T为在t采样时刻的测量向量。

(4)对建立的时滞矩阵构造动态正交全局局部模型，求取全局局部特征提取的目标函数J_G-L的投影矩阵A，具体地：

首先建立如下的全局特征提取的目标函数J_global：

上式中，Y＝[y₁，y₂,.....,y_n]＝A^TX，A＝[a₁,a₂,…,a_d]为J_global的投影矩阵，X＝[x₁,x₂,…,x_n]为采样矩阵，为局部均值矩阵，为全局结构矩阵；

局部特征提取的目标函数J_local通过邻域保持嵌入算法得到，具体过程如下：

首先，利用k近邻算法，选择近邻点，如果x_i为x_j的近邻点，则邻域图中节点i到节点j通过一条边连接，反之，如果不是近邻点，则没有连接，通过这种方法构成邻域图，邻域图的权重系数W＝[w₁,w₂,...,w_n]通过下式求得：

∑_jw_ij＝1,j＝1,2,…n

则局部特征提取的目标函数J_local可以表示为：

上式中，I为单位向量，M＝(I-W)^T(I-W)，限定条件为a^TXX^Ta＝1，则全局局部特征提取的目标函数J_G-L如下所示：

通过构造正交基向量来保持数据结构，{a₁,a₂,…a_d}为正交基向量，第一个向量a₁为(M)^-1Ga对应的最小的特征值所对应的特征向量，当d＞2时，如下式：

上式中限定条件为：引入拉格朗日乘子，则可转化为下式：

对C^(d)求取偏导则为：

对上式乘以可得：

通过连续的乘以可以获得如下d-1个式子：

上式可以表示为：

式中，μ^(d-1)＝[μ₁,…,μ_d-1]，A^(d-1)＝[a₁,…,a_d-1]，B^(d-1)＝[A^(d-1)]^T(XX^T)^-1A^(d-1)；

通过乘以(M)^-1，则上式可以转化为：

2(M)^-1Ga_d-2λa_d-(M)^-1A^(d-1)μ^(d-1)＝0

{I-(M)^-1A^(d-1)[B^(d-1)]^-1[A^(d-1)]^T}·(M)^-1Ga_d＝λa_d

上式中，I为单位向量，λ为引入的最小化判别因子，a_d为最小的特征值所对应的特征向量；A＝{a₁,a₂,…a_d}为(M)^-1Ga_d通过迭代计算的最小特征值所对应的特征向量，即正交全局局部特征提取的目标函数J_G-L的投影矩阵。

(5)通过求取的投影矩阵A，求取T²和SPE的控制限

在特征空间和残差空间分别建立T²和SPE统计模型，设Y(y₁,…,y_n)∈R^n×d为建立的时滞矩阵的低维表示，则特征空间和残差空间可以描述为：

Y＝A^TX＝(B^TB)^-1B^TX

E＝X-BY

上式中，A^T＝(B^TB)^-1B^T为转移矩阵，表示主元空间，表示残差空间，E为残差矩阵，T²和SPE统计模型为：

T²＝y^TΛ^-1y

SPE＝e^Te

在求得正常数据的SPE和T²统计后，通过核密度估计算法求取控制限SPE_α和求取过程为：

对一个给定的数据集j，则j的密度分布为：

上式中，j_i为观测值，σ为窗口宽度，n为观测值的个数，K为核函数。因此通过选择高斯核函数设检验水平为0.95下的SPE和T²的控制限为SPE_α和则有:

可以求取控制限为SPE_α和

(6)Pensim2.0仿真平台不仅可以产生正常工况下数据，还提供了三种故障类型，分别为：通风率故障，搅拌功率故障，底物流加速率故障。为了验证算法的效果，本发明中引入故障类型2，即变量2搅拌功率(Agitator power)故障，在采样时间200～400h时加入+0.03的斜坡信号作为故障信号，产生测试数据X_test(400×10)，所产生的数据进行标准化后为X′_test(400×10)，将其作为故障样本以供在线检测。

分别用多向主成分分析算法(MPCA)、多向邻域保持嵌入算法(MNPE)、多向全局邻域保持嵌入算法(MGNPE)和本发明提出的算法来对故障样本进行监控，得到图1-4的监控图。由于本发明能够消除过程数据的动态时序相关性，并在降维过程中保持了过程数据的正交全局局部结构，因此基于此方法的故障检测效果得到明显的改善。

(7)通过已经求取的正交全局局部特征提取的目标函数J_G-L的投影矩阵A，将新的测试数据X_test(400×10)进行投影，Y_test＝A^TX_test。

Y_test＝A^TX_test＝(B^TB)^-1B^TX_test

E_test＝(I-BA^T)X_test＝(I-B(B^TB)^-1B^T)X_test

因此在线数据的T²和SPE统计可以描述为：

如果SPE_test＞SPE_α或者则说明发生故障，否则没有故障发生。

Claims (6)

1.一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，包括以下步骤：

(1)利用集散控制系统收集间歇过程正常运行下的各关键变量数据，形成正常运行状态下的训练样本X∈R^I×J×K，其中R为实数集并满足(I×J×K)的三维分布，I表示批次，J表示变量，K为采样点；

(2)将训练样本X首先沿着批次方向展开为二维数据X∈R^I×KJ，然后将展开的二维数据进行标准化，使得各个变量的均值为0，方差为1，然后将标准化后的二维数据重新排列成X∈R^KI×J；

(3)在混合展开的二维数据X∈R^KI×J基础上，建立一个时滞矩阵X_D来消除过程变量的时序自相关；

(4)对建立的时滞矩阵X_D构造动态正交全局局部模型，且其过程是首先建立全局局部特征提取的目标函数J_G-L，然后引入广义正交向量以求取目标函数J_G-L的投影矩阵；

(5)分别在特征空间和残差空间建立T²和SPE统计模型，求取控制限；

(6)收集在线的过程数据，进行标准化处理；

(7)利用已经建立的动态正交全局局部模型将在线数据进行投影，通过T²和SPE统计模型，判断故障的发生。

2.根据权利要求1所述的一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，步骤(3)中：时滞矩阵X_D如下所示：

上式中，d为窗口的宽度，T为总的采样数，X(t)＝[x_1,t,x_2,t,…,x_J,t]^T为在t采样时刻的测量向量。

3.根据权利要求1所述的一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，步骤(4)中：建立全局局部特征提取的目标函数J_G-L的具体过程如下：

首先建立全局特征提取的目标函数J_global：

上式中，Y＝[y₁,y₂,.....,y_n]＝A^TX，A＝[a₁,a₂,…,a_d]为J_global的投影矩阵，X＝[x₁,x₂,…,x_n]为采样矩阵，为局部均值矩阵，为全局结构矩阵；

局部特征提取的目标函数J_local通过邻域保持嵌入算法得到，具体过程如下：

首先，利用k近邻算法，选择近邻点，如果x_i为x_j的近邻点，则邻域图中节点i到节点j通过一条边连接，反之，如果不是近邻点，则没有连接，通过这种方法构成邻域图，邻域图的权重系数W＝[w₁,w₂,...,w_n]通过下式求得：

∑_jw_ij＝1,j＝1,2,…n

则局部特征提取的目标函数J_local表示为：

上式中，I为单位向量，M＝(I-W)^T(I-W)，限定条件为a^TXX^Ta＝1，则全局局部特征提取的目标函数J_G-L如下所示：

则上式转化为求解广义特征值问题，即：

Ga＝λMa

设{a₁，a₂...a_d}为正交基向量，第一个向量a₁为(M)^-Ga对应的最小的特征值所对应的特征向量，当d＞2时，如下式：

上式中限定条件为：引入拉格朗日乘子，则转化为下式：

对C^(d)求取偏导则为：

对上式乘以可得：

通过连续的乘以可以获得如下d‐1个式子：

上式可表示为：

式中，μ^(d-1)＝[μ₁,…,μ_d-1]，A^(d-1)＝[a₁,…,a_d-1]，B^(d-1)＝[A^(d-1)]^T(XX^T)^-1A^(d-1)；

通过乘以(M)^-1，则上式转化为：

2(M)^-1Ga_d-2λa_d-(M)^-1A^(d-1)μ^(d-1)＝0

{I-(M)^-1A^(d-1)[B^(d-1)]^-1[A^(d-1)]^T}·(M)^-1Ga_d＝λa_d

上式中，I为单位向量，λ为最小化判别因子，a_d为最小的特征值所对应的特征向量；A＝{a₁,a₂,…a_d}为(M)^-1Ga_d通过迭代计算的最小特征值所对应的特征向量，即正交全局局部特征提取的目标函数J_G-L的投影矩阵，也即J_G-L的正交解。

4.根据权利要求1所述的一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，步骤(5)为：

设Y(y₁,…,y_n)∈R^n×d为原始变量空间X(x₁,…,x_m)∈R^n×m的低维表示，y_i＝A^Tx_i，则特征空间和残差空间可以描述为：

Y＝A^TX＝(B^TB)^-1B^TX

E＝X-BY

上式中，A^T＝(B^TB)^-1B^T为转移矩阵，表示主元空间，表示残差空间，E＝[e₁，e₂,.....e_n]为残差矩阵，T²和SPE统计模型为：

T²＝y^TΛ^-1y

SPE＝e^Te

在求得正常数据的SPE和T²统计后，通过核密度估计算法求取控制限SPE_α和求取过程为：

对一个给定的数据集j，则j的密度分布为：

上式中，j_i为观测值，σ为窗口宽度，n为观测值的个数，K为核函数，因此通过选择高斯核函数设检验水平为0.95下的SPE和T²的控制限为SPE_α和则有:

可以求取控制限为SPE_α和

5.根据权利要求1所述的一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，步骤(6)中，对在线采集的测试数据X_test进行标准化，使得各个变量的均值为0，方差为1。

6.根据权利要求1所述的一种基于动态正交的全局局部间歇过程故障检测方法，其特征在于，步骤(7)中，通过已经求取的全局局部特征提取的目标函数J_G-L的投影矩阵A，将采样的测试数据X_test进行投影，则Y_test＝A^TX_test

Y_test＝A^TX_test＝(B^TB)^-1B^TX_test

E_test＝(I-BA^T)X_test＝(I-B(B^TB)^-1B^T)X_test

因此在线数据的T²和SPE统计可以描述为：

如果SPE_test＞SPE_α或者则说明发生故障，否则没有故障发生。