本发明涉及一种故障检测方法，具体是涉及一种基于多采样概率核主成分模型的化工过程故障检测方法。
背景技术：
在现代流程工业中，随着集散控制系统(DCS)的应用和计算机技术的进步，工业现场采集和存储了大量的在线和离线测量数据，基于多元统计分析的过程监测(MSPM)技术得到了快速发展，它具有基于数据、降维、易于可视化与易于实际应用等优点，在化工、制药、半导体制造等多个工业领域得到了广泛应用。其中，主成分分析(PCA)和偏最小二乘估计(PLS)及其扩展方法是MSPM技术的代表性模型。基于PCA的过程监测技术通过检测过程变量的波动以有效监测生产过程是否发送了异常工况。然而，该类方法无法判别生产过程的波动是否影响了最终的产品质量。传统的PLS及PLVR模型均假设过程与质量变量一一对应，以提取过程与质量变量的相关关系，并对生产过程进行相应的故障检测。
然而对于一些化工过程，特别是合成某些化工原料(比如合成氨生产过程)的过程中，有些质量数据的采样频率很低，诸如组分、浓度、分子量等质量变量难以测量，需要通过实验室化验，一般以小时或天为测量单位。在此期间，诸如温度、压力、流量等过程变量则以分钟甚至是秒为测量单位被DCS大量采集并记录。因此，需要提出一种针对化工生产过程多采样率样本的数据建模与故障检测方法。同时，由于化工生产过程具有强烈的耦合性和非线性，因此设计的故障检测方法需要考虑到此类复杂数据特性。
技术实现要素：
本发明的目的在于针对现有技术的不足，提供一种基于多采样概率核主成分分析模型的化工生产过程故障检测方法。
本发明同时提供一种基于多采样概率核主成分分析模型的合成氨生产过程故障检测方法，以合成氨生产过程不同采样率的过程变量、质量变量为建模样本，提取了这些样本间的非线性相关关系，并在此模型的基础上建立了故障检测方法，以实现合成氨生产过程的过程监控。
本发明的目的是通过以下技术方案来实现的：
一种基于多采样概率核主成分模型的化工过程故障检测方法，包括如下步骤：
(1)收集待检测化工过程在正常运行时的至少三种不同采样率的过程变量，作为样本数据组成建模用的训练样本集；
(2)对训练样本集进行预处理，将预处理后的数据集投影到高维特征空间，得到训练样本集对应的高维映射值数据集；
(3)求取高维映射值数据集的核函数，结合当前模型参数的初始值，利用期望最大化(EM)算法，估计得到模型参数的更新值；
(4)当所有模型参数的更新值满足收敛要求时，得到模型参数的最终值，进而得到多采样概率核主成分分析模型，进入步骤(5)；否则，将模型参数的更新值作为模型参数的初始值，返回步骤(3)；
(5)根据建立的多采样概率核主成分分析模型，得到对应的T2和SPE统计量检测控制限；
(6)在线收集待检测化工过程实际运行过程中的工艺参数，得到测试样本集，对测试样本进行步骤(2)的预处理，然后将预处理后的数据集投影到高维特征空间，得到测试样本集对应的高维映射值数据集，根据得到的多采样概率核主成分分析模型，计算测试样本的和SPEtest统计量，然后与步骤(5)得到检测控制限对比，得出化工过程的在线检测结果。
本发明中，所述的化工过程可以是各种化工过程，特别是需要多参数控制的化工过程，比如可是造纸废水处理过程、钢铁冶炼过程或者合成氨生产过程等等。所述的过程变量可以是温度、压强、浓度、重量、pH值、流量等等，其采样率可以相同，也可以不同。步骤(1)中，相同时间内，采集三种不同采样率的过程变量。对于某一具体的化工生产过程，每一种特定采样率下的过程变量可能是一个或者多个，所以，三个不同采样率的过程变量可能是三个或者三个以上。根据采样率的不同，一个过程变量一般由多个过程参数组成。
在合成氨生产过程中，有些过程参数的采样频率很低，比如组分、浓度、分子量等质量变量难以测量，需要通过实验室化验，一般以小时或天为测量单位。而温度、压力、流量等过程变量则以分钟甚至是秒为测量单位被DCS大量采集并记录，为提高生产过程中的检测精度，作为优选，可将本发明的检测方法应用于合成氨生产过程中，即：所述化工过程为合成氨生产过程。
作为优选，步骤(1)中，利用集散控制系统收集至少三种不同采样率的过程变量。
步骤(1)中，在化工过程正常运行时，采集训练样本集，保证了检测的准确度。
步骤(1)中，收集三种不同采样率的过程变量，得到三个样本集，样本数量分别为K,N,J,K，N，J均为大于等于1的正整数，且J≤N≤K；三个样本集分别为X，Y，Z：
X∈RM×K,X＝{x1,x2,...,xK}
Y∈RV×N,Y＝{y1,y2,...,yN}
Z∈RU×J,Z＝{z1,z2,...,zJ}
其中：R表示实数集；M为第一种采样率过程变量数，V为第二种采样率过程变量数，U为第三种采样率过程变量数。
步骤(2)中，所述预处理主要包括归一化处理，当然可以根据需要增加去噪声等预处理过程。通过步骤(2)的预处理，使得每个过程变量的均值为零，方差为1。通过步骤(2)的预处理步骤，使得不同数量级的过程参数进一步简化，转换为统一数量级的数据集，进一步方便了后续的建模。
步骤(2)中，对于此步骤得到的三部分的二维数据矩阵集X∈RM×K,Y∈RV×N,Z∈RR×J，每一部分的采样率不同，首先将归一化后的数据投影到高维特征空间中，三个样本X，Y，Z的高维映射值分别为Φ(x),Φ(y),Φ(z)，分别组成高维映射值数据集Φ(X),Φ(Y),Φ(Z)，假设Φ(x),Φ(y),Φ(z)与潜隐变量间存在线性相关关系，则：
Φ(x)＝Wt+e
Φ(y)＝Pt+f (1)
Φ(z)＝Qt+g
其中F1,F2,F3分别为Φ(x),Φ(y),Φ(z)的维数；t∈RD×1为该模型的潜隐变量，D是潜隐变量的维度；为模型负载矩阵；e,f,g分别为过程的测量噪声，它们分别服从各向同性的高斯分布为分别为高斯分布的方差。
步骤(3)中，利用期望最大化(EM)算法来估计模型参数，在进行参数估计时，虽然无法直接获得模型的高维映射值数据集Φ(X),Φ(Y),Φ(Z),但是可以通过计算其对应的核函数Φ(X)ΦT(X),Φ(Y)ΦT(Y),Φ(Z)ΦT(Z)的值来进行估计；
首先，由于该模型的模型参数W,P,Q无法直接估计，定义模型的辅助参数为ΨX,ΨY,ΨZ,CX,CY,CZ，其被定义为：
ΨZ＝Ψ3∈RJ×D
CX＝C3+C2+C1
CY＝C3+C2 (3)
CZ＝C3
其中Ψ3,Ψ2,Ψ1为潜隐变量的后验概率的期望值组成的矩阵，其被定义为：
C3,C2,C1分别为潜隐变量的后验概率的二阶矩的和，其被定义为：
其中：表示第i个样本对应的高维映射值Φ(xi),Φ(yi),Φ(zi)的潜隐变量ti的后验概率的期望值，i＝1～K；
表示第i个样本对应的高维映射值Φ(xi),Φ(yi),Φ(zi)的潜隐变量的后验概率的二阶矩，i＝1～K。
对模型参数随机进行初始化；在模型参数估计的E步，根据当前的模型参数，获取模型潜隐变量估计的更新值，具体公式为：
其中分别为Ψ3,Ψ2,Ψ1,C3,C2,C1的估计值；
其中辅助变量定义如下：
其中：HX，HY，HZ分别为Φ(X),Φ(Y),Φ(Z)的核函数，即HX，HY，HZ分别为Φ(X)ΦT(X),Φ(Y)ΦT(Y),Φ(Z)ΦT(Z)；
HX1、HX2、HX3分别为HX的子矩阵；HY2、HY3分别为HY的子矩阵；
在M步，根据E步的更新结果，获取模型参数{σ32,σ22,σ12}的更新值如下：
其中：trace()表示矩阵的迹。
步骤(4)中，利用模型参数的更新值Θnew与其原模型参数Θold的差异，如果对于所有模型参数，满足||Θnew-Θold||2＜ε，则进入到步骤(5)，否则，返回步骤(3)，其中，ε为模型收敛的阈值，对于多个模型参数，其收敛的阈值可以相同也可以不同。
步骤(5)中，根据建立的多采样概率核主成分分析模型估计训练样本的潜隐变量的期望值t，潜隐变量的方差var-1(t|x)和模型预测误差{eX,eY,eZ}，构建相应的T2和SPE统计量监测统计限；具体为：根据步骤(3)模型收敛后计算得到的{Ψ3,Ψ2,Ψ1}获得训练样本潜隐变量的期望值tnormal、潜隐变量的方差var-1(tnormal|xnormal)，利用潜隐变量的期望值tnormal和方差var-1(tnormal|xnormal)，可构建T2统计量如下：
T2＝tnormalTvar-1(tnormal|xnormal)tnormal (7)；
其中，{Ψ3,Ψ2,Ψ1}与tnormal的对应关系如下：
T2统计量的控制限由χ2分布估计如下：
其中，D为潜隐变量的维度个数。
同时，基于模型的预测误差，还可构建SPE统计量以反应模型残差空间的变化：
其中eX，eY，eZ分别为Φ(x),Φ(y),Φ(z)的预测误差，分别为不同采样率数据的模型预测值；
进一步推倒可得：
其中：H(x,x)为Φ(x)与Φ(x)的核矩阵，即H(x,x)＝<Φ(x),Φ(x)>＝Φ(x)Φ(x)T，其中<>为计算内积；H(z,z)为Φ(z)与Φ(z)的核矩阵，H(x,X)为Φ(x)与Φ(X)的核矩阵，H(y,Y)为Φ(y)与Φ(Y)的核矩阵，H(z,Z)为Φ(z)与Φ(Z)的核矩阵，H(X,X)为Φ(X)与Φ(X)的核矩阵，H(Y,Y)为Φ(Y)与Φ(Y)的核矩阵，H(Z,Z)为Φ(Z)与Φ(Z)的核矩阵；其中H(X,X)＝HX，H(Y,Y)＝HY，H(Z,Z)＝HZ；
SPE统计量的控制限估计方法为：即SPE服从卡方分布，其中，
其中mean()为求均值，var()为求方差，代表卡方分布，g和h分别为卡方分布的系数和自由度；通过式(11)可以求得g和h，进而求得SPE统计量的控制限。
作为优选，所述步骤(7)中，
所述测试样本集中必定包括采用率为K的测试样本集Xtest，有可能包括采用率为N或J的测试样本集Ytest或Ztest，Xtest由一个或多个测试样本xtest组成，当测试样本xtest被采集后，其潜隐变量的期望值ttest为：
其中H(X,xtest)为Φ(X)与Φ(xtest)的核函数；
利用潜隐变量的期望值ttest，构建测试样本的统计量如下：
其中var-1(ttest|x1test)为潜隐变量的方差；然后基于模型的预测误差，构建测试样本的SPEtest统计量；
其中：H(xtest,X)表示Φ(xtest)与Φ(X)的核函数；H(ytest,Y)表示Φ(ytest)与Φ(Y)的核函数；H(ztest,Z)表示Φ(ztest)与Φ(Z)的核函数；
最后，判断其是否超过统计限，得出合成氨生产过程的在线检测结果。
在检测过程中，一般以采用率最大的采样率作为检测过程的检测采用率，当然也可以设置更大的采样率；这样一来，测试样本集中，必定包含采用率为K的测试样本xtest，在某些采样数据中，也有可能包括ytest或者ztest。我们直接采用采样最多的测试样本xtest，计算其潜隐变量的期望值ttest。计算SPEtest统计量时，可根据实际是否包括ytest或者ztest，选择式(14)中对应的公式进行计算，进而进行最终的故障判断。
本发明中，H()表示两个矩阵或者向量的核函数，以H(A,B)为例,H(A,B)＝<Φ(A),Φ(B)>＝Φ(A)Φ(B)T。出现的上标T均表示对应矩阵或者向量的转置。上标＾表示对应参数的估计值。
本发明的有益效果是：
本发明利用基于多采样概率核主成分分析模型的核学习方法，建立了一个有效的非线性故障检测模型，并克服了化工生产过程采样速率不同所造成的问题，提高了该过程的在线检测效率和性能，从而使得化工生产过程更加可靠，产品质量监测更加稳定。
通过将合成氨生产过程中采集到的不同频率的数据矩阵沿着时间点方向展开为新的二维数据矩阵，并对其顺序进行重新排列，并建立了多采样率的模型，克服了过程中采样率不同的问题。基于核方法，对于合成氨生产过程的非线性数据特性进行特征提取。相比目前的其它合成氨过程的检测方法，本发明不仅可以大大提高合成氨过程非线性多采样率建模以及故障检测的效果，减小故障误报率和漏报率，而且很大程度上改善了模型的预测能力，基于此的故障检测方法也更有效。
附图说明
图1为本发明的基于多采样概率核主成分分析模型的合成氨生产过程故障检测方法流程图。
具体实施方式
以合成氨生产过程为例，对本发明做进一步说明：
一种基于多采样概率核主成分分析模型的合成氨生产过程故障检测方法，该方法针对合成氨过程的故障检测问题，首先利用集散控制系统收集正常工作状态下不同采样率的数据，建立多采样概率核主成分分析模型。该模型结构由期望最大化算法估计得到。在此基础上，利用模型的潜隐变量和预测误差构建了两个检测统计量T2和SPE及其对应的统计限和SPElim。对在线合成氨生产过程进行检测，获得测试样本，然后可利用已有的模型结构估计测试样本的潜隐变量及预测误差，并计算出相应的统计量，并得到最终的故障检测结果。
参见图1，本发明一种基于多采样概率核主成分分析模型的合成氨生产过程故障检测方法，包括以下步骤：
第一步：利用集散控制系统收集的合成氨生产过程正常运行的不同采样率的数据组成数据建模用的训练样本集,假设采集到了三种采样率的正常样本，其样本数量分别为K,N,J,且J≤N≤K：
X∈RM×K,X＝{x1,x2,...,xK}
Y∈RV×N,Y＝{y1,y2,...,yN}
Z∈RU×J,Z＝{z1,z2,...,zJ}
R表示实数集；M为第一种采样率过程变量数，V为第二种采样率过程变量数，U为第三种采样率过程变量数，将这些数据存入历史数据库；
第二步：对数据集X,Y,Z进行预处理和归一化，使得各个过程变量的均值为零，方差为1，将归一化后的数据投影到高维空间中，其映射结果为Φ(X),Φ(Y),Φ(Z)，在高维空间中的投影值Φ(x),Φ(y),Φ(z)被假设与隐变量间存在线性相关关系，则
Φ(x)＝Wt+e
Φ(y)＝Pt+f
Φ(z)＝Qt+g(1)
其中F1,F2,F3分别为Φ(x),Φ(y),Φ(z)的维数；t∈RD×1为该模型的潜隐变量，D是潜隐变量的维度；为模型负载矩阵；e,f,g分别为过程的测量噪声，它们分别服从各向同性的高斯分布为
第三步：利用期望最大化(EM)算法来更新模型参数，在E步利用当前模型参数估计潜隐变量的后验概率；之后，在M步中分别针对各个参数求似然函数的一阶偏导数，以获取模型参数的更新值。最后，反复迭代E步和M步直至达到模型收敛条件。
首先，定义模型的辅助参数为ΨX,ΨY,ΨZ,CX,CY,CZ，其被定义为：
ΨZ＝Ψ3∈RJ×D
CX＝C3+C2+C1
CY＝C3+C2 (3)
CZ＝C3
其中Ψ3,Ψ2,Ψ1为潜隐变量的后验概率的期望值组成的矩阵，其被定义为：
C3,C2,C1分别为潜隐变量的后验概率的二阶矩的和，其被定义为：
其中：表示第i个样本对应的高维映射值Φ(xi),Φ(yi),Φ(zi)的潜隐变量ti的后验概率的期望值，i＝1～K；
表示第i个样本对应的高维映射值Φ(xi),Φ(yi),Φ(zi)的潜隐变量的后验概率的二阶矩，i＝1～K；
对模型参数随机进行初始化；在模型参数估计的E步，根据当前的模型参数的初始值，获取模型潜隐变量估计的更新值，具体公式为：
其中分别为Ψ3,Ψ2,Ψ1,C3,C2,C1的估计值；其中辅助变量定义如下：
HX＝[HX3 HX2 HX1]∈RK×K，HX3∈RK×J,HX2∈RK×(N-J),HX1∈RK×(K-N)
HY＝[HY3 HY2]∈RN×N，HY3∈RN×J,HY2∈RN×(N-J)
其中：HX，HY，HZ分别为Φ(X),Φ(Y),Φ(Z)的核函数，即HX，HY，HZ分别为Φ(X)ΦT(X),Φ(Y)ΦT(Y),Φ(Z)ΦT(Z)；
其中HX,HY根据公式需要，将得到的不同采样率的原始数据，进行分类后，求出其核函数，因此，HX包含了三种采样率的数据，HY包含了两种采样率的数据。
在M步，根据E步的更新结果，获取模型参数{σ32,σ22,σ12}的更新值如下
第四步：利用模型参数的更新值Θnew与其原模型参数Θold的差异，如果||Θnew-Θold||2＜ε，则进入到第五步，否则，返回第三步，其中，ε为模型收敛的阈值；
第五步：根据建立的多采样概率核主成分分析模型估计训练样本的潜隐变量的期望值t，潜隐变量的方差var-1(t|x)和模型预测误差{eX,eY,eZ}，构建相应的T2和SPE统计量监测统计限；
模型收敛后计算得到的{Ψ3,Ψ2,Ψ1}获得训练样本潜隐变量的期望值tnormal，其对应关系为：
利用潜隐变量的期望值tnormal，可构建T2统计量如下：
T2＝tnormalTvar-1(tnormal|xnormal)tnormal (7)
T2统计量的控制限由χ2分布估计如下：
其中，D为潜隐变量的个数。
同时，基于模型的预测误差，还可构建SPE统计量以反应模型残差空间的变化：
进一步推倒可得：
其中：H(x,x)为Φ(x)与Φ(x)的核矩阵；H(z,z)为Φ(z)与Φ(z)的核矩阵，H(x,X)为Φ(x)与Φ(X)的核矩阵，H(y,Y)为Φ(y)与Φ(Y)的核矩阵，H(z,Z)为Φ(z)与Φ(Z)的核矩阵，H(X,X)为Φ(X)与Φ(X)的核矩阵，H(X,X)＝HX，H(Y,Y)为Φ(Y)与Φ(Y)的核矩阵，H(Y,Y)＝HY，H(Z,Z)为Φ(Z)与Φ(Z)的核矩阵，H(Z,Z)＝HZ；
SPE统计量的控制限估计方法为：SPE:g·其中，
gh＝mean(SPE)
2g2h＝var(SPE) (11)
其中mean()为求均值，var()为求方差，代表卡方分布，g和h分别为卡方分布的系数和自由度；通过式(11)可以求得g和h，进而求得SPE统计量的控制限。
第六步：在线收集新的合成氨生产过程中，与训练样本集过程变量对应的不同采样率样本数据，采样时间间隔由采用率最高的过程变量的采用率确定(可以等于或者大于最高的采用率，当然也可以小于，本实施例中，选择等于)，得到测试样本集Xtest(也可能同时包括Ytest、Ztest)，并进行归一化；然后将预处理后的数据集投影到高维特征空间，得到测试样本集对应的高维映射值数据集Φ(Xtest)(Φ(Ytest)、Φ(Ztest))；
Xtest∈RM×K',Xtest＝{xtest1,xtest2,...,xtestK'}
Ytest∈RV×N',Ytest＝{ytest1,ytest2,...,ytestN'}
Ztest∈RU×J',Ztest＝{ztest1,ztest2,...,ztestJ'}
其中K',N',J'分别采样时间间隔内，采集的样本数量；
第七步：采用多采样概率核主成分分析模型对其进行检测，计算测试样本的和SPEtest统计量，当测试样本xtest被采集后，其潜隐变量的期望值为
利用潜隐变量的期望值ttest，可构建T2统计量如下：
其中var-1(ttest|x1test)为潜隐变量的方差。基于模型的预测误差，构建SPEtest统计量
其中：H(xtest,X)表示Φ(xtest)与Φ(X)的核函数；H(ytest,Y)表示Φ(ytest)与Φ(Y)的核函数；H(ztest,Z)表示Φ(ztest)与Φ(Z)的核函数；
判断其是否超过步骤(5)得到的控制限，得出合成氨生产过程的在线检测结果：如果超过所述的控制限，则判断为发生故障；否则则判断生产过程正常。