王 杰，程 顺，刘 松，欧阳福生，赵明洋

(1.华东理工大学化工学院石油加工研究所，上海 200237；2.中国石化上海高桥分公司)

催化裂化(FCC)汽油约占中国汽油总量的70%[1]。成品汽油中的硫主要来自FCC汽油，因此，如何降低FCC汽油硫含量从而满足环保要求成为了一个亟需解决的问题。S Zorb工艺[2-4]以其良好的脱硫能力而成为FCC汽油脱硫的主要技术。该技术能在较低的辛烷值损失下生产硫质量分数低于10 μg/g的汽油，但在工业装置运行过程中运行结果存在一定的波动，例如某S Zorb装置4年运行的汽油辛烷值(RON)损失平均值为1.45个单位，大大高于同类装置的平均指标(RON损失0.6个单位)[5]，从而影响了企业经济效益。由于来自FCC装置汽油性质的波动，以及回炼的不合格脱硫汽油产品比例和性质的不确定性，使得S Zorb装置的汽油原料性质经常会发生变化。在实际操作过程中，该装置对原料油性质的变化考虑较少，未能及时针对不同的汽油原料调整操作条件，因此较难降低汽油辛烷值损失。

聚类分析[6]是机器学习、数据挖掘、模式识别等领域的重要组成内容，它在无标记样本的条件下根据不同的准则对数据进行分类，找到这些数据信息的内部结构和规律。由于同一类别的数据具有相似性，在数据分析过程中，可以将一个类别中的数据对象作为一个整体来处理，从而达到简化数据、精确分析目标的目的。传统的聚类算法主要包括基于模型、划分、密度和层次的聚类[7]。此外，针对更复杂多样的数据来源，还有模糊聚类、综合聚类等新的聚类算法[8]。目前未见针对S Zorb工艺原料油聚类分析研究的文献报道。

本研究以某S Zorb装置运行4年的FCC汽油原料性质数据为基础，以产品辛烷值损失最小为目标，使用聚类算法建立原料油的聚类评价模型，将原料油性质最为相近的聚为一类，并对每一类原料油的特征进行描述，以期为针对不同种类的S Zorb装置原料油分别建立产品辛烷值预测模型，从而为寻找最佳操作条件来降低产品辛烷值损失提供指导。

1 数据收集

S Zorb工艺的目的是脱除FCC汽油中的微量硫化物，但在脱硫反应过程中常常伴随着烃类的化学反应，其中以烯烃的临氢饱和反应为主，烃类的裂化、异构化和芳构化等反应为辅[9]，因此产品汽油中烃类组成会发生变化，从而影响到产品汽油的辛烷值。可见，原料油中的饱和烃、芳烃和烯烃含量的变化均会影响到汽油的辛烷值；在研究如何降低辛烷值损失的同时，必须以保证S Zorb工艺的脱硫效果为前提；密度反映原料油的轻重，而原料油中的重组分会发生裂化反应，从而影响产品辛烷值；溴值是有机化合物不饱和度的重要指标，在临氢条件下，不饱和烃会发生加成、异构化反应，从而影响产品汽油的辛烷值。因此，本研究采用饱和烃含量、烯烃含量、芳烃含量、硫含量、密度和溴值6个变量作为原料油聚类分析的变量。

本研究的聚类分析所需的原料油性质数据采集自某S Zorb装置的实验室信息管理系统(LIMS)，采集时间为2015年6月1日至2019年4月16日，通过整理共得到117组样本，原料油性质分布见表1。

表1 原料油性质分布

2 原料油聚类模型的建立

K-means算法理论可靠、算法简单、收敛速率快并且对大数据集有较高的效率[10]。模糊聚类算法可以得到每个样本属于各个类别的不确定程度，应用范围广泛[11]。本研究利用MATLAB编程平台，分别采用K-means和模糊C均值聚类算法(FCM)来建立原料油性质的聚类模型。

2.1 K-means算法聚类分析

K-means是一种被广泛使用的基于划分的聚类算法，它通过迭代不断更新聚类中心，最终获得最小化的误差平方和(SSE)，其计算式为：

(1)

式中：Yi为第i类的样本点；x为Yi中的一个数据点；mi为Yi的聚类中心。该算法的具体步骤是：

(1)确定最佳聚类数目kopt，即把数据分为k类。

(2)随机选取k个数据点，作为k类数据的初始聚类中心。

(3)利用欧几里得函数式(2)分别计算每个数据到各初始聚类中心的距离，根据就近原则进行数据分类，即数据离哪个聚类中心最近就被分为哪一类。

(2)

式中，xik和xjk分别表示第i个和第j个数据对象在属性k上的取值。

(4)所有数据分好类后，重新计算每一类新的聚类中心。

(5)重复上述过程，直到满足终止条件：没有聚类中心再发生变化。

本研究先确定最佳聚类数目的上下限，令kmin=2，根据经验式(3)求得kmax=11(N为样本数目)，然后通过综合手肘法和轮廓系数法的结果来确定最佳聚类数目kopt。

(3)

手肘法的基本思路是：做出SSE随k变化的图形，找到图中斜率突变的节点，即为kopt。由于节点处的图形类似于手的肘部，所以称之为手肘法。

轮廓系数法的基本思路是：使组内距离最小，组间距离最大。求出不同k值所对应的平均轮廓系数(S)后，找出最大S对应的k值即为最佳聚类数目。根据Kaufman等提出的轮廓系数的相关概念[12]，平均轮廓系数的定义式为：

(4)

式中：Si为第i个样本的轮廓系数；n为样本总数。Si的定义式为：

(5)

式中：a(i)为样本i和与其同属一类的其他所有样本之间的平均距离，b(i)为样本i和与其不属于一类的其他各个类中所有样本的平均距离的最小值。a(i)和b(i)的定义式为：

(6)

(7)

式中：假设样本i属于第t类；nt表示第t类的样本总数；p表示非t类；Yt、Yp分别表示第t类和非t类的样本点；d(i，j)表示欧式空间中样本i和j的欧式距离。

图1为SSE随k的变化。图2为S随k的变化。由图1可见，当k=3和k=5时，曲线的斜率变化较为明显，但是均无显著的肘部位置，需要结合平均轮廓系数法的结果来确定kopt。由图2可见，k=9时所对应的S最大。但当k=9时，图1的曲线斜率变化趋于平缓，说明k值已经远远偏离kopt。综合手肘法和轮廓系数法的分析结果，k=3时，图1曲线的斜率变化较为明显，图2中对应的S也较大，因此最终取kopt=3。

图1 SSE随k的变化

图2 S随k的变化

表2为kopt=3时每一类原料油聚类中心的原料性质。由表2可以看出，聚类中心相互之间区分度较大的原料性质是硫含量、饱和烃含量、烯烃含量和溴值，而溴值是表示有机化合物中不饱和程度的一种指标，与烯烃含量成正相关关系，因此本研究选取硫含量、饱和烃含量和烯烃含量的高低来度量不同类别原料油性质的差异。第1类原料油的烯烃含量最低，硫和饱和烃含量较高，因此可命名为低烯烃高硫高饱和烃原料油；第2类原料油的硫和饱和烃含量最低，烯烃含量最高，因此命名为低硫低饱和烃高烯烃原料油；第3类原料油的硫和饱和烃含量最高，烯烃含量较低，因此命名为高硫高饱和烃低烯烃原料油。

表2 kopt=3时每一类聚类中心(K-means)的原料性质

图3为溴值与芳烃含量之间的二维关系。由图3可见：3类原料油数据点分布集中，分布区域相对独立且具有较明显界限，同时聚类中心分布也较分散；溴值在3类原料油中的分布存在较为明显的梯度，可见溴值对最终聚类结果影响较大；芳烃含量从低到高在3类原料油中均有数据点分布，可见芳烃含量对最终的聚类结果影响不大。

图3 K-means算法溴值与芳烃含量之间的二维关系●—第1类原料油； ▲—第2类原料油； ■—第3类原料油。较大标记为聚类中心。图4～图8同

图4为密度与烯烃含量之间的二维关系。由图4可见：3类原料油数据点分布集中，分布区域相对独立且具有较明显界限，同时聚类中心分布也较分散；烯烃在3类原料油中的分布存在明显的梯度，说明烯烃含量对最终的聚类结果影响较大；密度从低到高在3类原料油都有数据点，可见密度对最终的聚类结果影响不大。

图4 K-means算法密度与烯烃含量之间的二维关系

图5为烯烃含量、饱和烃含量与硫含量之间的三维关系。由图5可见：与二维图相比，三维图中各类原料油数据点之间的界限更为明显，数据点分布更为集中，聚类中心更加分散；烯烃分布结果与图4一致；饱和烃含量与硫含量在3类原料油中的分布存在明显的梯度，说明饱和烃含量与硫含量对最终的聚类结果影响均较大。

图5 K-means算法烯烃含量、饱和烃含量与硫含量之间的三维关系

2.2 FCM算法聚类分析

传统的聚类方法(如K-means算法)是一种非此即彼的、清楚的分类算法，但是自然界中有许多概念是模糊的，例如大小、老少、冷暖等，它们不能用简单的是或否去表示。模糊聚类就是用以表示这种具有模糊属性的数据集合。它将是或否，也就是集合{0，1}扩充为[0，1]，即将二值0和1推广到无穷多值[13]。最常见的模糊聚类方法是FCM算法，该方法设计简单、解决问题范围广，并且最终可以转化为优化问题而借助经典数学的分线性规划理论求解，且易于实现，因此，FCM算法得到了广泛的应用[8]。FCM算法的具体步骤为：

(1)初始化参数，包括模糊参数、聚类数目、隶属度矩阵和终止条件。

(2)计算聚类中心ci。

(8)

式中：n为数据元素的个数；ci为第i类的聚类中心；μ为某数据元素对某一组的隶属度；m为模糊加权指数；xj为第j个数据点。

(3)按式 (9) 计算目标函数的值，如果目标函数满足终止条件则终止计算，否则继续执行步骤(4)。

(9)

式中：J为目标函数；dik为第i个数据元素到第k个聚类中心的欧几里得距离。

(4)按式 (10) 重新计算隶属度矩阵，然后返回到步骤(2)。

(10)

最小聚类数cmin=2，最大聚类数也采用式(3)进行计算，得cmax=11。FCM算法的聚类参数设置值见表3。

表3 FCM算法主要参数

采用吴成茂等[14]提出的归一划分系数式(11)来确定最佳聚类数目copt。该准则是在Bezdek准则[式(12)]基础上得到的。

(11)

(12)

式中：U为隶属度矩阵；F(U；c)为划分系数；GF(U；c)为归一划分系数。

按c从2到11运行MATLAB程序，得到相应的GF(U；c)的值以及收敛迭代次数，结果见表4。由表4可以发现：随c的增加，程序收敛次数也不断增加；当GF(U；c)取最大值时，其所对应的c值即为最佳聚类数目，copt=3。

表4 FCM算法聚类结果

分别计算每个数据样本对每一类原料油的隶属度，结果见表5。以样本1为例，它对第1，2，3类原料油的隶属度分别为0.472 2，0.161 6，0.366 2，隶属度之和为1，其中，对第1类原料油的隶属度最大，因此样本1被归为第1类原料油。

表5 原料油样本对于不同类别的隶属度

表6为copt=3时每一类原料油聚类中心的原料性质，与表2 中K-means聚类算法得到的聚类中心相比虽然略有差别，但是总体类似。同样地，第1类原料油可命名为低烯烃高硫高饱和烃原料油；第2类原料油可命名为低硫低饱和烃高烯烃原料油；第3类原料油可命名为高硫高饱和烃低烯烃原料油。

表6 copt=3时每一类聚类中心(FCM)的原料性质

图6～图8分别为FCM算法的溴值与芳烃含量，密度与烯烃含量，以及烯烃含量、饱和烃含量与硫含量的分布。由图6～图8可见，3类原料油数据点分布集中，分布区域相对独立且具有较明显界限，同时聚类中心分布也比较分散。溴值、烯烃含量、饱和烃含量和硫含量在3类原料油中的分布均存在较为明显的梯度，可见溴值、烯烃含量、饱和烃含量和硫含量对最终聚类结果影响较大；芳烃含量和密度从低值到高值在3类原料油中均有数据点分布，可见芳烃含量和密度对最终的聚类结果影响不大。该结果与K-means算法的聚类结果一致。

图6 FCM算法溴值与芳烃含量之间的二维关系

图7 FCM算法密度与烯烃含量之间的二维关系

图8 FCM算法烯烃含量、饱和烃含量与硫含量之间的三维关系

2.3 两种聚类算法的结果比较

对比表2和表6列出的K-means和FCM聚类算法在最佳聚类数目时各类聚类中心的原料性质可以看出，不同算法所对应的每一类聚类中心重合度均较高。K-means和FCM聚类算法均将117组样本分为3类，K-means算法的第1至第3类原料油的样本数分别为38，38，41；FCM算法的第1至第3类原料油的样本数分别为43，38，36。两种聚类算法仅在5组样本的分类结果上有差异，其余112组样本的分类结果完全一致，这进一步说明了两种聚类算法结果的一致性。

3 结 论

以某S Zorb装置的原料油性质数据为基础，通过对工艺和反应机理的分析，选取原料油性质中的硫含量、饱和烃含量、烯烃含量、芳烃含量、密度和溴值6个聚类变量，采用MATLAB编程平台，分别使用K-means和FCM聚类算法进行原料油聚类研究。结果表明，两种聚类法均将原料油分为边界比较明显的3类，分类结果高度重合，聚类分析效果均较好，为针对每一类原料油建立相应的汽油辛烷值预测模型，从而为寻找到使汽油辛烷值损失最小的操作条件打下良好基础。