CN104635724A - 基于动态核独立分量分析的天然气净化过程异常检测方法 - Google Patents

Abstract

本发明公开了一种基于动态核独立分量分析的天然气净化过程异常检测方法，该诊断方法按如下步骤进行：确定所监测的高含硫天然气净化脱硫生产过程中的工艺参数；对高含硫天然气净化脱硫生产过程进行数据采集；对数据进行预处理；分析数据的自回归模型，确定其动态时滞后阶次，得到其动态扩展矩阵；在核主元空间对扩展矩阵进行白化处理，提取核主元分量，并采用独立分量分析估计独立元；计算独立元对应的SPE统计量和T²统计量，并分析统计量是否超限，若超限则发生异常工况，否则正常，再结合T²贡献图法，进行异常工况参数追溯。本方法能够及时检测故障发生，并追溯工艺操作参数导致故障发生原因，从而为系统故障排查和恢复提供决策参考依据。

Description

基于动态核独立分量分析的天然气净化过程异常检测方法

技术领域

本发明属于高含硫天然气脱硫生产过程故障检测与诊断技术，涉及一种基于动态核独立分量分析的天然气净化过程异常检测方法。

背景技术

高含硫天然气净化脱硫工业流程复杂，过程工艺参数众多，受温度、压力、流量、设备老化和原料气处理量等不确定因素影响，是典型的复杂非线性动态特性化工系统。高含硫天然气净化脱硫过程主要包括以下部分：主吸收塔MDEA溶液吸收酸性组分H₂S和CO₂，水解反应器脱除(COS)，再生塔MDEA溶液的循环再生以及换热过程，具体工艺流过程如图2所示。一旦系统发生异常，如何及时检测故障发生，并追溯工艺操作参数导致故障发生原因，从而为系统故障排查和恢复提供决策参考依据，对指导实际工业生产具有重要的现实意义。

原料气处理量负荷波动、脱硫单元吸收塔发泡、硫磺回收装置硫收率下降是高含硫天然气过程常见的三种异常工况。目前高含硫天然气异常工况检测与诊断主要依赖现场生产经验，异常工况发生机理描述的不够透彻，缺乏对高含硫天然气异常工况的检测与诊断机制。

基于多元统计的故障检测与诊断方法是工业大数据技术在故障检测与诊断领域应用的代表，已成为工业过程监控的热点。基于多元统计过程的高含硫天然气净化脱硫过程故障检测与诊断采用基于T²和SPE统计量质量控制图和PCA、ICA技术分析DCS采集的相关工艺参数，从而为系统运行状态提供决策参考依据，并结合贡献图法实现异常工况的参数追溯。

核独立元分析在故障诊断中是一种常见的特征分析方法，这种方法是利用非线性映射Φ和Ψ将原始随机向量空间x_i和y_i映射到高维空间F_x和F_y，使得数据在特征空间中成为线性数据，再在特征空间中利用独立成分分析法对变换后的数据提取独立元，从而实现特征变量互相独立(不相关)。但现有的核独立分量分析没有考虑到高含硫天然气净化脱硫过程是一个动态过程，采集的原始数据在时序上呈现自相关性。

基于T²和SPE统计量的质量控制图通过判断独立元对应的T²和SPE统计量是否超出设定的控制限，判断系统运行状态。若超过控制限，系统处于故障状态；否则系统处于正常状态。贡献图法通过分析故障时刻不同参数对故障贡献程度实现异常参数追溯。

发明内容

本发明的目的在于克服现有技术不能及时检测故障发生的不足，提供一种基于动态核独立分量分析的天然气净化过程异常检测方法，它能够及时检测故障发生，并追溯工艺操作参数导致故障发生原因，从而为系统故障排查和恢复提供决策参考依据。

本发明的目的是这样实现的：

一种基于动态核独立分量分析的天然气净化过程异常检测方法，该诊断方法按如下步骤进行：

步骤1：确定所监测的高含硫天然气净化脱硫生产过程中的m个工艺参数，其中，m＝10，各工艺参数分别为：x₁表示脱硫吸收塔胺液入口流量，x₂表示尾气吸收塔胺液入口流量，x₃表示原料气处理量，x₄表示半富胺液循环量，x₅表示一级吸收塔胺液入口温度，x₆表示二级吸收塔胺液入口温度，x₇表示闪蒸罐压力，x₈表示重沸器A口蒸汽消耗量，x₉表示重沸器B口蒸汽消耗量，x₁₀表示蒸汽预热器流量；

步骤2：对高含硫天然气净化脱硫生产过程进行数据采集，所得的样本数据为[X_m×N,Y₁,Y₂]，其中：X为工艺参数，m为工艺参数数量，N为采集数据样本数量，Y₁为H₂S含量，Y₂为CO₂含量；

采集高含硫天然气净化脱硫生产过程中的工艺参数，并采集对应工艺参数下生产的天然气当中的H₂S含量和CO₂含量，用于后续数据分析。

步骤3：对步骤2所得的样本数据进行预处理：

3.1剔除采集数据中缺失参数的样本，并保证样本满足企业净化气技术指标，得到的数据为[X_m×n,Y₁,Y₂]，n为处理后采集样本数量，n＜N；

通过剔除采集数据中缺失参数的样本能够得到最能反映出生产过程实际特性的有效数据。

3.2进行数据归一化处理，得到的工业数据为[X′_m×n,Y′₁,Y′₂]；

数据归一化的目的是为了把不同来源的数据统一到一个参考系下，方便排查。

步骤4：分析步骤3所得的工业数据[X'_m×n,Y₁',Y₂']的自回归模型，确定其动态时滞后阶次h，得到其动态扩展矩阵；

计算时滞后阶次h的方法为：首先处理h＝0的静态情况，数据矩阵中的静态关系数等于变量数与主元个数之差；然后令h＝1,计算新的动态关系数，新的动态关系数等于变量数减去主元个数和第一步计算出的静态关系数，h函数逐渐增加，新关系函数按以下公式递推计算：

$r_{\mathrm{new}}\left(h\right)=r\left(h\right)-\underset{i=0}{\overset{h-1}{\mathrm{\Σ}}}(h-i+1)r_{\mathrm{new}}\left(i\right),$

直至r_new(h)<0，即持续到没有新的静态和动态关系为止，观测向量x_i＝[x_i1,x_i2,…x_im]，考虑其h时滞后阶次信息，则 $x_i^h=[x_{i1},x_{i-11},...,x_{i-h1},x_{i2},x_{i-12},...x_{i-l2},...,x_{\mathrm{im}},x_{i-1m},...x_{i-\mathrm{hm}}]$ 最终得到的动态扩展矩阵为： $X^h={[x_{h+1}^h,x_{h+2}^h,...{,x}_n^h]}^T\&Element;R^{(n-h)\&times;m(h+1)};$

引入自回归模型，计算动态时滞后阶次：利用动态特性确定算法确定自回归模型阶次，并对数据进行动态扩展处理，提取动态阶次信息。

步骤5：在核主元空间对扩展矩阵X^h进行白化处理，提取核主元分量Z^h，并采用独立分量分析估计独立元s^h；

步骤6：计算独立元s^h对应的SPE统计量和T²统计量，并分析统计量是否超限，若超限则发生异常工况，否则正常，再结合T²贡献图法，进行异常工况参数追溯。

步骤3.1中的企业净化气技术指标为H₂S含量低于6mg/m³即4ppmv(Y₁<4)，CO₂含量低于3％(Y₂<3％)。

步骤3.2中的具体归一化处理方法如下：

$x_i^{\&prime;}=\frac{x_i-x_{\min }}{x_{\max }-x_{\min }},y_i^{\&prime;}=\frac{y_i-y_{\min }}{y_{\max }-y_{\min }},$

其中，x_i为归一化前的操作变量，x′_i为归一化后的操作变量，x_min为归一化前x_i的最小值，x_max为归一化前x_i的最大值，y_i为归一化前的观测变量，y′_i为归一化后观测变量，y_min为归一化前y_i最小值，y_max为归一化前y_i最大值。

由于采用了上述技术方案，本发明具有如下有益效果：

本发明克服了动态、非高斯过程故障检测与诊断精度不高缺陷，它能够及时检测故障发生，并追溯工艺操作参数导致故障发生原因，从而为含硫天然气净化系统故障排查和恢复提供决策参考依据。本发明还可以推广到其他复杂动态化工系统，为其故障检测与诊断提供参考依据。

附图说明

图1为动态核独立分量分析算法流程；

图2为某高含硫天然气脱硫模拟工艺流程；

图3为本发明的流程图；

图4为基于T²统计量动态独立分量分析异常工况检测；

图5为基于SPE统计量动态独立分量分析异常工况检测；

图6为T²统计量对应的贡献图法进行工艺参数追溯。

附图标记

图2中，1水解反应器进料分离器；2水解反应器预热器；3水解反应器；4水解反应器进/出料换热器；5水解反应器后冷器；6二级吸收塔；7一级吸收塔；8贫胺液泵；9中间胺液泵；10中间胺液冷却器；11贫胺液后冷器；12再生塔；13再生塔顶空冷器；14胺液再生塔顶回流罐；15酸水回流泵；16再生塔底贫胺液泵；17贫富液换热器；18贫胺液空冷器；19胺液闪蒸罐；Feed gas：原料天然气；Treated gas：净化气；Acid gas：酸气。

具体实施方式

下面结合附图和实施例对本发明作进一步说明：

实施例1：

参见图3，一种基于动态核独立分量分析的天然气净化过程异常检测方法，该诊断方法按如下步骤进行：

步骤1：确定所监测的高含硫天然气净化脱硫生产过程中的m个工艺参数，其中，m＝10，各工艺参数分别为：x₁表示脱硫吸收塔胺液入口流量，x₂表示尾气吸收塔胺液入口流量，x₃表示原料气处理量，x₄表示半富胺液循环量，x₅表示一级吸收塔胺液入口温度，x₆表示二级吸收塔胺液入口温度，x₇表示闪蒸罐压力，x₈表示重沸器A口蒸汽消耗量，x₉表示重沸器B口蒸汽消耗量，x₁₀表示蒸汽预热器流量；

步骤2：对高含硫天然气净化脱硫生产过程进行数据采集，所得的样本数据为[X_m×N,Y₁,Y₂]，其中：X为工艺参数，m为工艺参数数量，N为采集数据样本数量，Y₁为H₂S含量，Y₂为CO₂含量；

采集高含硫天然气净化脱硫生产过程中的工艺参数，并采集对应工艺参数下生产的天然气当中的H₂S含量和CO₂含量，用于后续数据分析。

步骤3：对步骤2所得的样本数据进行预处理：

3.1剔除采集数据中缺失参数的样本，并保证样本满足企业净化气技术指标，得到的数据为[X_m×n,Y₁,Y₂]，n为处理后采集样本数量，n＜N；

企业净化气技术指标为H₂S含量低于6mg/m³即4ppmv(Y₁<4)，CO₂含量低于3％(Y₂<3％)。

通过剔除采集数据中缺失参数的样本能够得到最能反映出生产过程实际特性的有效数据。

3.2进行数据归一化处理，得到的工业数据为[X′_m×n,Y′₁,Y′₂]；

数据归一化的目的是为了把不同来源的数据统一到一个参考系下，方便排查。

具体归一化处理方法如下：

$x_i^{\&prime;}=\frac{x_i-x_{\min }}{x_{\max }-x_{\min }},y_i^{\&prime;}=\frac{y_i-y_{\min }}{y_{\max }-y_{\min }},$

其中，x_i为归一化前的操作变量，x′_i为归一化后的操作变量，x_min为归一化前x_i的最小值，x_max为归一化前x_i的最大值，y_i为归一化前的观测变量，y′_i为归一化后观测变量，y_min为归一化前y_i最小值，y_max为归一化前y_i最大值。步骤3.2中的具体归一化处理方法如下：

$x_i^{\&prime;}=\frac{x_i-x_{\min }}{x_{\max }-x_{\min }},y_i^{\&prime;}=\frac{y_i-y_{\min }}{y_{\max }-y_{\min }},$

其中，x_i为归一化前的操作变量，x′_i为归一化后的操作变量，x_min为归一化前x_i的最小值，x_max为归一化前x_i的最大值，y_i为归一化前的观测变量，y′_i为归一化后观测变量，y_min为归一化前y_i最小值，y_max为归一化前y_i最大值。

步骤4：分析步骤3所得的工业数据[X'_m×n,Y₁',Y₂']的自回归模型，确定其动态时滞后阶次h，得到其动态扩展矩阵；

计算时滞后阶次h的方法为：首先处理h＝0的静态情况，数据矩阵中的静态关系数等于变量数与主元个数之差；然后令h＝1,计算新的动态关系数，新的动态关系数等于变量数减去主元个数和第一步计算出的静态关系数，h函数逐渐增加，新关系函数按以下公式递推计算：

$r_{\mathrm{new}}\left(h\right)=r\left(h\right)-\underset{i=0}{\overset{h-1}{\mathrm{\&Sigma;}}}(h-i+1)r_{\mathrm{new}}\left(i\right),$

直至r_new(h)<0，即持续到没有新的静态和动态关系为止，观测向量x_i＝[x_i1,x_i2,…x_im]，考虑其h时滞后阶次信息，则 $x_i^h=[x_{i1},x_{i-11},...,x_{i-h1},x_{i2},x_{i-12},...x_{i-l2},...,x_{\mathrm{im}},x_{i-1m},...x_{i-\mathrm{hm}}]$ 最终得到的动态扩展矩阵为： $X^h={[x_{h+1}^h,x_{h+2}^h,...{,x}_n^h]}^T\&Element;R^{(n-h)\&times;m(h+1)};$

引入自回归模型，计算动态时滞后阶次：利用动态特性确定算法确定自回归模型阶次，并对数据进行动态扩展处理，提取动态阶次信息。

步骤5：在核主元空间对扩展矩阵X^h进行白化处理，提取核主元分量Z^h，并采用独立分量分析估计独立元s^h；

步骤6：计算独立元s^h对应的SPE统计量和T²统计量，并分析统计量是否超限，若超限则发生异常工况，否则正常，再结合T²贡献图法，进行异常工况参数追溯，如图1所示。

以某高含硫净化厂天然气脱硫脱硫装置的生产数据进行分析，监控数据由DCS系统即时显示。选取某一个月内744组数据进行影响因素分析，包括有脱硫吸收塔胺液入口管流量x₁(t/h)、尾气吸收塔胺液入口管流量x₂(t/h)、原料气处理量x₃(kNm3/h)、半富胺液的循环量x₄(t/h)、一级吸收塔胺液入口温度x₅(℃)、二级吸收塔胺液入口温度x₆(℃)、闪蒸罐压力x₇(MPa)、重沸器A入口蒸气消耗量x₈(kg/h)、重沸器B入口蒸气消耗量x₉(kg/h)和蒸气预热器流量x₁₀(t/h)。采样数据如表1所示。原料气处理量负荷波动是高含硫天然气净化脱硫过程典型的异常工况。

表1 12个原始特征的744组样本:

剔除采集数据中缺失参数的样本，并保证样本满足企业净化气技术指标，即H₂S含量低于6mg/m³(Y₁<6)，，CO₂含量低于3％(Y₂<3)，并进行归一化至[-1,1]得到[X′_m×n,Y′₁,Y′₂]，如表2所示。

表2归一化处理后数据:

表3动态时滞后阶次确定方法:

分析数据[X′_m×n,Y′₁,Y′₂]的自回归模型，确定其动态时滞后阶次h，如表3所示。得到其动态扩展矩阵，如表4所示。

核函数选用高斯函数K(x,y)＝exp(-||x-y||²/σ)，σ＝3600。所选取的独立元对应的特征值满足λ>0.0001。由步骤5，6可得不同观测时刻独立元对应的T²和SPE统计量，如图4,5所示。

表4动态扩展矩阵数据

当某一时刻T²和SPE统计量超出控制限，则该时刻工况处于异常状态；否则表明该时刻工况正常。异常工况检测效果如表5所示所示。

表5基于T²统计量和SPE统计量异常工况检测

对于发生异常工况的时刻，采用T²统计量的贡献图法进行工艺参数追溯，挖掘异常发生的原因，如图6所示。

通过分析T²统计量对应的贡献图法进行工艺参数追溯，发现第三个变量原料气处理量是导致异常发生的原因，属于原料气处理量负荷波动异常工况。

最后说明的是，以上优选实施例仅用以说明本发明的技术方案而非限制，尽管通过上述优选实施例已经对本发明进行了详细的描述，但本领域技术人员应当理解，可以在形式上和细节上对其作出各种各样的改变，而不偏离本发明权利要求书所限定的范围。

Claims (3)

1.一种基于动态核独立分量分析的天然气净化过程异常检测方法，其特征在于该诊断方法按如下步骤进行：

步骤1：确定所监测的高含硫天然气净化脱硫生产过程中的m个工艺参数，其中，m＝10，各工艺参数分别为：x₁表示脱硫吸收塔胺液入口流量，x₂表示尾气吸收塔胺液入口流量，x₃表示原料气处理量，x₄表示半富胺液循环量，x₅表示一级吸收塔胺液入口温度，x₆表示二级吸收塔胺液入口温度，x₇表示闪蒸罐压力，x₈表示重沸器A口蒸汽消耗量，x₉表示重沸器B口蒸汽消耗量，x₁₀表示蒸汽预热器流量；

步骤2：对高含硫天然气净化脱硫生产过程进行数据采集，所得的样本数据为[X_m×N,Y₁,Y₂]，其中：X为工艺参数，m为工艺参数数量，N为采集数据样本数量，Y₁为H₂S含量，Y₂为CO₂含量；

步骤3：对步骤2所得的样本数据进行预处理：

3.1剔除采集数据中缺失参数的样本，并保证样本满足企业净化气技术指标，得到的数据为[X_m×n,Y₁,Y₂]，n为处理后采集样本数量，n＜N；

3.2进行数据归一化处理，得到的工业数据为[X′_m×n,Y′₁,Y′₂]；

步骤4：分析步骤3所得的工业数据[X'_m×n,Y₁',Y₂']的自回归模型，确定其动态时滞后阶次h，得到其动态扩展矩阵；

计算时滞后阶次h的方法为：首先处理h＝0的静态情况，数据矩阵中的静态关系数等于变量数与主元个数之差；然后令h＝1,计算新的动态关系数，新的动态关系数等于变量数减去主元个数和第一步计算出的静态关系数，h函数逐渐增加，新关系函数按以下公式递推计算：

$r_{\mathrm{new}}\left(h\right)=r\left(h\right)-\underset{i=0}{\overset{h-1}{\mathrm{\&Sigma;}}}(h-i+1)r_{\mathrm{new}}\left(i\right),$

直至r_new(h)<0，即持续到没有新的静态和动态关系为止，观测向量x_i＝[x_i1,x_i2,…x_im]，考虑其h时滞后阶次信息，则 $x_i^h=[x_{i1},x_{i-11},...x_{i-h1},x_{i-h1},x_{i2},x_{i-12},...x_{i-l2},...,x_{\mathrm{im}},x_{i-1m},...x_{i-\mathrm{hm}}]$ 最终得到的动态扩展矩阵为： $X^h={[x_{h+1}^h,x_{h+2}^2,...,x_n^h]}^T\&Element;R^{(n-h)\&times;m(h+1)};$

步骤5：在核主元空间对扩展矩阵X^h进行白化处理，提取核主元分量Z^h，并采用独立分量分析估计独立元s^h；

步骤6：计算独立元s^h对应的SPE统计量和T²统计量，并分析统计量是否超限，若超限则发生异常工况，否则正常，再结合T²贡献图法，进行异常工况参数追溯。

2.根据权利要求1所述的基于动态核独立分量分析的天然气净化过程异常检测方法，其特征在于：步骤3.1中的企业净化气技术指标为H₂S含量低于6mg/m³，CO₂含量低于3％。

3.根据权利要求1所述的基于动态核独立分量分析的天然气净化过程异常检测方法，其特征在于：步骤3.2中的具体归一化处理方法如下：

$x_i^{\&prime;}=\frac{x_i-x_{\min }}{x_{\max }-x_{\mathrm{mun}}},y_i^{\&prime;}=\frac{y_i-y_{\min }}{y_{\max }-y_{\min }},$

其中，x_i为归一化前的操作变量，x′_i为归一化后的操作变量，x_min为归一化前x_i的最小值，x_max为归一化前x_i的最大值，y_i为归一化前的观测变量，y′_i为归一化后观测变量，y_min为归一化前y_i最小值，y_max为归一化前y_i最大值。