一种基于局部特征异常因子的初馏塔采样数据粗差判别方法

本发明涉及一种采样数据粗差判别方法，特别涉及一种基于局部特征异常因子的初馏塔采样数据粗差判别方法。

背景技术：

初馏塔是炼油厂常减压装置的龙头，其控制的好坏直接影响原油收率的提高和后续工序的操作稳定。初顶石脑油干点是初馏塔重要的控制指标，但该指标无合适的在线分析仪进行测定，只能离线分析且时间长。因此初顶石脑油干点软测量模型的建立可以为生产操作条件及时调整提供指导。而影响初顶石脑油干点的各个因素数据的准确性又直接关系到软测量模型的性能。因此确保这些数据的正确无误，去除采样数据中的粗差数据是非常重要的。

因测量仪表和计算机技术的广泛应用，初馏塔以及其他化工生产单元会每隔几分钟会采样并存储过程数据。在数据采集的过程中，因内在或外在因素导致采样数据存在较大的偏差，但是通过单个测量变量来看又在其允许的变化内。这些特殊的采样数据被称为粗差数据，采样数据中存在粗差数据会严重影响后续的基于数据的分析与建模，实际使用采样数据时是绝对不能采用粗差数据的。

在初馏塔生产运行过程中，能够影响初顶馏分干点的几个主要的测量变量包括：初馏塔处理量、塔顶温度、塔顶压力、顶回流带出能量、回流比、初顶石脑油流量、初顶循带出能量、初中段带出能量、及进料温度。由于生产过程积累了海量的、高维样本数据对所有数据逐个人为判断是否为粗差。然而，建立初顶馏分干点软测量模型时，又要求这些采样数据不能存在粗差数据。此外，初馏塔的采样数据中是否存在粗差数据，以及到底存在多少个粗差数据是无法预知和确认的，这给初馏塔采样数据粗差判别任务带来了严峻的挑战。

虽然现有科研文献与专利技术材料中，有不同种类的粗差判别方法技术。这些粗差判别方法要么依赖于平方马氏距离，要么依赖于计算采样数据的协方差矩阵。因此，若是采样数据之间线性相关性较强，其相应的协方差矩阵会存在不可逆的情况。而且，初馏塔采样数据运行的动态性，会导致采样数据分布情况复杂，不能简单的依靠马氏距离所定义的椭圆来一概而论。

技术实现要素：

本发明所要解决的主要技术问题是：如何通过局部特征异常因子来判别出初馏塔采样数据中的粗差数据。具体来讲，本发明方法通过设计一种局部特征异常因子来量化各个样本数据的异常程度，从而可以将局部特征异常因子较大的样本数据判别为粗差。本发明方法涉及的局部特征异常因子，需要先通过为各个样本数据找到多个近邻样本数据，然后以近邻样本数据为参考数据集，通过优化得到能区分该样本数据与其近邻之间差异的判别型特征，从而计算得到局部特征异常因子。

本发明方法解决上述问题所采用的技术方案为：一种基于局部特征异常因子的初馏塔采样数据粗差判别方法，包括以下所示步骤。

步骤(1)：确定初馏塔中影响初顶馏分干点的9个变量，具体包括：初馏塔进料流量，进料温度，塔顶温度，塔顶压力，塔顶回流带出能量，回流比、初顶石脑油流量，初顶循带出能量，初中断带出能量；采集初馏塔这9个变量的n组样本数据x1，x2，…，xn，并将其组成为一个n×9维的数据矩阵x∈r^n×9；其中，r^n×9表示n×9维的实数矩阵，r表示实数集，xi∈r^1×9表示第i个样本数据，r^1×9表示1×9维的实数向量，i∈{1，2，…，n}。

步骤(2)：根据公式分别为数据矩阵x中的列向量z1，z2，…，z9实施抗差标准化处理，得到抗差标准化处理后的数据矩阵其中，zk和分别表示x和中第k列的列向量，k∈{1，2，…，9}表示变量的序号，μ(zk)表示计算列向量zk中所有元素的中位数，表示计算列向量zk中所有元素的中位绝对差，列向量表示计算列向量中所有元素的中位数，|zk-μ(zk)|表示计算zk-μ(zk)的绝对值。

步骤(3)：按照如下所示步骤(3.1)至步骤(3.5)依次计算当i分别等于1，2，…，n时的局部特征异常因子l1，l2，…，ln。

步骤(3.1)：初始i＝1，并将中的第i行向量记为

步骤(3.2)：根据公式计算行向量与数据矩阵中第j行向量之间的平方距离，并将中与之间平方距离最小的c个行向量组成参考矩阵其中，j∈{1，2，…，n}且j≠i，上标号t为矩阵或向量的转置符号。

步骤(3.3)：求解广义特征值问题中最大特征值λi对应的特征向量pi后，再计算变换向量

步骤(3.4)：计算局部特征后，再根据公式li＝(si)²计算第i个局部特征异常因子li。

步骤(3.5)：判断是否满足i＜n；若是，则设置i＝i+1后返回步骤(3.2)；若否，则得到n个局部特征异常因子l1，l2，…，ln。

步骤(4)：计算n个局部特征异常因子l1，l2，…，ln的中位数μ(l)以及中位绝对差δ(l)后，再确定局部特征异常因子的变化上限llim＝μ(l)+3×δ(l)；其中，l＝[l1，l2，…，ln]。

步骤(5)：根据如下所示步骤(5.1)至步骤(5.3)逐个判别初馏塔采样数据中的粗差数据。

步骤(5.1)：初始化i＝1。

步骤(5.2)：判断是否满足条件li＞llim；若是，则将第i个样本数据xi判别为粗差数据；若否，则第i个样本数据不是粗差数据。

步骤(5.3)：判断是否满足i＜n；若是，则设置i＝i+1后返回步骤(5.2)；若否，则结束初馏塔采样数据的粗差判别。

值得说明的是，上述步骤(3.3)求解广义特征值问题的实施过程实际上是为寻找到能最大化程度区分与其近邻组成的参考矩阵之间差异的变换向量wi。为了最大化程度的实现区分，需要通过寻找变换向量wi，将相应的行向量xi及其参考矩阵经过wi变换后的局部特征从距离上最大化的分开来，即可得到如下所示的最大问题：

上式①通过将参考矩阵对应的特征向量约束成单位长度的向量，并同时最大化经wi变换后的局部特征与原点之间的距离。上式①的求解可使用经典的朗格朗日乘子法，即通过拉格朗日乘子λi构造如下所示的朗格朗日函数j。

其中，求解j相对于wi的偏微分，并设置其等于0即可得到如下所示的等式关系

进而转换为广义特征值问题：若在该广义特征值问题的等式两边分别同左乘则可得到因此，λi等于上式①中的最大化目标，即广义特征值问题需要计算最大的特征值。

在本发明方法的实施步骤(3.3)中，先通过求解广义特征值问题的特征向量pi后，再计算变换向量wi＝pi/||pi||。不仅完成了在线特征分析技术所涉及的广义特征问题的求解，而且还进一步将变换向量wi归一化成单位长度。

通过以上所述实施步骤，本发明方法的优势介绍如下。

本发明方法未曾涉及计算马氏距离，而是寻找各个样本数据的近邻样本数据，并通过变换向量最大化两者之间的差异，从而通过局部特征与原点之间的距离来判别各个样本是否为粗差。从这个角度上看，本发明方法避免了协方差矩阵不可逆的问题，并同时给出样本数据是否为粗差数据及其采样数据集中有多少个粗差数据的判别。

附图说明

图1为本发明方法的实施流程示意图。

图2为粗差识别可视化图像。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

如图1所示，本发明公开了一种基于局部特征异常因子的初馏塔采样数据粗差判别方法。下面结合一个具体应用实例来说明本发明方法的具体实施方式。

采集某炼油厂常减压装置初馏塔中能影响初顶馏分干点的9个变量的200组样本数据，即n＝200。这200组样本数据中是否存在粗差数据，以及存在几个粗差数据是未知的。通过本发明方法实施粗差判别。

步骤(1)：确定初馏塔中影响初顶馏分干点的9个变量，具体包括：初馏塔进料流量，进料温度，塔顶温度，塔顶压力，塔顶回流带出能量，回流比、初顶石脑油流量，初顶循带出能量，初中断带出能量；采集初馏塔这9个变量的n组样本数据x1，x2，…，xn，并将其组成为一个n×9维的数据矩阵x∈r^n×9。

步骤(2)：根据公式分别为数据矩阵x中的列向量z1，z2，…，z9实施抗差标准化处理，得到抗差标准化处理后的数据矩阵

步骤(3)：按照如下所示步骤(3.1)至步骤(3.5)依次计算当i分别等于1，2，…，n时的局部特征异常因子l1，l2，…，ln。

步骤(3.1)：初始i＝1，并将中的第i行向量记为

步骤(3.2)：根据公式计算行向量与数据矩阵中其它行的行向量之间的平方距离，并将中与之间平方距离最小的c个行向量组成参考矩阵其中，表示数据矩阵中第j行的行向量，j∈{1，2，…，n}且j≠i，上标号t为矩阵或向量的转置符号。

步骤(3.3)：求解广义特征值问题中最大特征值λi对应的特征向量pi后，再计算变换向量

步骤(3.4)：计算局部特征后，再根据公式li＝(si)²计算第i个局部特征异常因子li。

步骤(3.5)：判断是否满足i＜n；若是，则设置i＝i+1后返回步骤(3.2)；若否，则得到n个局部特征异常因子l1，l2，…，ln。

步骤(4)：计算n个局部特征异常因子l1，l2，…，ln的中位数μ(l)以及中位绝对差δ(l)后，再确定局部特征异常因子的变化上限llim＝μ(l)+3×δ(l)；其中，l＝[l1，l2，…，ln]。

步骤(5)：根据如下所示步骤(5.1)至步骤(5.3)逐个判别初馏塔采样数据中的粗差数据。

步骤(5.1)：初始化i＝1。

步骤(5.2)：判断是否满足条件li＞llim；若是，则将第i个样本数据xi判别为粗差数据；若否，则第i个样本数据不是粗差数据。

步骤(5.3)：判断是否满足i＜n；若是，则设置i＝i+1后返回步骤(5.2)；若否，则结束初馏塔采样数据的粗差判别。

将得到的200个局部特征异常因子l1，l2，…，l200以及变化上限llim绘制于图2中。在图2中，横线表示局部特征异常因子的变化上限llim，根据步骤(5)可知超出变化上限的样本数据即为粗差数据。因此，从图2中可以发现，哪几个样本数据通过本发明方法判别为粗差。