王 杰，程 順，劉 松，歐陽福生，趙明洋

(1.華東理工大學(xué)化工學(xué)院石油加工研究所，上海 200237；2.中國(guó)石化上海高橋分公司)

催化裂化(FCC)汽油約占中國(guó)汽油總量的70%[1]。成品汽油中的硫主要來自FCC汽油，因此，如何降低FCC汽油硫含量從而滿足環(huán)保要求成為了一個(gè)亟需解決的問題。S Zorb工藝[2-4]以其良好的脫硫能力而成為FCC汽油脫硫的主要技術(shù)。該技術(shù)能在較低的辛烷值損失下生產(chǎn)硫質(zhì)量分?jǐn)?shù)低于10 μg/g的汽油，但在工業(yè)裝置運(yùn)行過程中運(yùn)行結(jié)果存在一定的波動(dòng)，例如某S Zorb裝置4年運(yùn)行的汽油辛烷值(RON)損失平均值為1.45個(gè)單位，大大高于同類裝置的平均指標(biāo)(RON損失0.6個(gè)單位)[5]，從而影響了企業(yè)經(jīng)濟(jì)效益。由于來自FCC裝置汽油性質(zhì)的波動(dòng)，以及回?zé)挼牟缓细衩摿蚱彤a(chǎn)品比例和性質(zhì)的不確定性，使得S Zorb裝置的汽油原料性質(zhì)經(jīng)常會(huì)發(fā)生變化。在實(shí)際操作過程中，該裝置對(duì)原料油性質(zhì)的變化考慮較少，未能及時(shí)針對(duì)不同的汽油原料調(diào)整操作條件，因此較難降低汽油辛烷值損失。

聚類分析[6]是機(jī)器學(xué)習(xí)、數(shù)據(jù)挖掘、模式識(shí)別等領(lǐng)域的重要組成內(nèi)容，它在無標(biāo)記樣本的條件下根據(jù)不同的準(zhǔn)則對(duì)數(shù)據(jù)進(jìn)行分類，找到這些數(shù)據(jù)信息的內(nèi)部結(jié)構(gòu)和規(guī)律。由于同一類別的數(shù)據(jù)具有相似性，在數(shù)據(jù)分析過程中，可以將一個(gè)類別中的數(shù)據(jù)對(duì)象作為一個(gè)整體來處理，從而達(dá)到簡(jiǎn)化數(shù)據(jù)、精確分析目標(biāo)的目的。傳統(tǒng)的聚類算法主要包括基于模型、劃分、密度和層次的聚類[7]。此外，針對(duì)更復(fù)雜多樣的數(shù)據(jù)來源，還有模糊聚類、綜合聚類等新的聚類算法[8]。目前未見針對(duì)S Zorb工藝原料油聚類分析研究的文獻(xiàn)報(bào)道。

本研究以某S Zorb裝置運(yùn)行4年的FCC汽油原料性質(zhì)數(shù)據(jù)為基礎(chǔ)，以產(chǎn)品辛烷值損失最小為目標(biāo)，使用聚類算法建立原料油的聚類評(píng)價(jià)模型，將原料油性質(zhì)最為相近的聚為一類，并對(duì)每一類原料油的特征進(jìn)行描述，以期為針對(duì)不同種類的S Zorb裝置原料油分別建立產(chǎn)品辛烷值預(yù)測(cè)模型，從而為尋找最佳操作條件來降低產(chǎn)品辛烷值損失提供指導(dǎo)。

1 數(shù)據(jù)收集

S Zorb工藝的目的是脫除FCC汽油中的微量硫化物，但在脫硫反應(yīng)過程中常常伴隨著烴類的化學(xué)反應(yīng)，其中以烯烴的臨氫飽和反應(yīng)為主，烴類的裂化、異構(gòu)化和芳構(gòu)化等反應(yīng)為輔[9]，因此產(chǎn)品汽油中烴類組成會(huì)發(fā)生變化，從而影響到產(chǎn)品汽油的辛烷值?？梢?，原料油中的飽和烴、芳烴和烯烴含量的變化均會(huì)影響到汽油的辛烷值；在研究如何降低辛烷值損失的同時(shí)，必須以保證S Zorb工藝的脫硫效果為前提；密度反映原料油的輕重，而原料油中的重組分會(huì)發(fā)生裂化反應(yīng)，從而影響產(chǎn)品辛烷值；溴值是有機(jī)化合物不飽和度的重要指標(biāo)，在臨氫條件下，不飽和烴會(huì)發(fā)生加成、異構(gòu)化反應(yīng)，從而影響產(chǎn)品汽油的辛烷值。因此，本研究采用飽和烴含量、烯烴含量、芳烴含量、硫含量、密度和溴值6個(gè)變量作為原料油聚類分析的變量。

本研究的聚類分析所需的原料油性質(zhì)數(shù)據(jù)采集自某S Zorb裝置的實(shí)驗(yàn)室信息管理系統(tǒng)(LIMS)，采集時(shí)間為2015年6月1日至2019年4月16日，通過整理共得到117組樣本，原料油性質(zhì)分布見表1。

表1 原料油性質(zhì)分布

2 原料油聚類模型的建立

K-means算法理論可靠、算法簡(jiǎn)單、收斂速率快并且對(duì)大數(shù)據(jù)集有較高的效率[10]。模糊聚類算法可以得到每個(gè)樣本屬于各個(gè)類別的不確定程度，應(yīng)用范圍廣泛[11]。本研究利用MATLAB編程平臺(tái)，分別采用K-means和模糊C均值聚類算法(FCM)來建立原料油性質(zhì)的聚類模型。

2.1 K-means算法聚類分析

K-means是一種被廣泛使用的基于劃分的聚類算法，它通過迭代不斷更新聚類中心，最終獲得最小化的誤差平方和(SSE)，其計(jì)算式為：

(1)

式中：Yi為第i類的樣本點(diǎn)；x為Yi中的一個(gè)數(shù)據(jù)點(diǎn)；mi為Yi的聚類中心。該算法的具體步驟是：

(1)確定最佳聚類數(shù)目kopt，即把數(shù)據(jù)分為k類。

(2)隨機(jī)選取k個(gè)數(shù)據(jù)點(diǎn)，作為k類數(shù)據(jù)的初始聚類中心。

(3)利用歐幾里得函數(shù)式(2)分別計(jì)算每個(gè)數(shù)據(jù)到各初始聚類中心的距離，根據(jù)就近原則進(jìn)行數(shù)據(jù)分類，即數(shù)據(jù)離哪個(gè)聚類中心最近就被分為哪一類。

(2)

式中，xik和xjk分別表示第i個(gè)和第j個(gè)數(shù)據(jù)對(duì)象在屬性k上的取值。

(4)所有數(shù)據(jù)分好類后，重新計(jì)算每一類新的聚類中心。

(5)重復(fù)上述過程，直到滿足終止條件：沒有聚類中心再發(fā)生變化。

本研究先確定最佳聚類數(shù)目的上下限，令kmin=2，根據(jù)經(jīng)驗(yàn)式(3)求得kmax=11(N為樣本數(shù)目)，然后通過綜合手肘法和輪廓系數(shù)法的結(jié)果來確定最佳聚類數(shù)目kopt。

(3)

手肘法的基本思路是：做出SSE隨k變化的圖形，找到圖中斜率突變的節(jié)點(diǎn)，即為kopt。由于節(jié)點(diǎn)處的圖形類似于手的肘部，所以稱之為手肘法。

輪廓系數(shù)法的基本思路是：使組內(nèi)距離最小，組間距離最大。求出不同k值所對(duì)應(yīng)的平均輪廓系數(shù)(S)后，找出最大S對(duì)應(yīng)的k值即為最佳聚類數(shù)目。根據(jù)Kaufman等提出的輪廓系數(shù)的相關(guān)概念[12]，平均輪廓系數(shù)的定義式為：

(4)

式中：Si為第i個(gè)樣本的輪廓系數(shù)；n為樣本總數(shù)。Si的定義式為：

(5)

式中：a(i)為樣本i和與其同屬一類的其他所有樣本之間的平均距離，b(i)為樣本i和與其不屬于一類的其他各個(gè)類中所有樣本的平均距離的最小值。a(i)和b(i)的定義式為：

(6)

(7)

式中：假設(shè)樣本i屬于第t類；nt表示第t類的樣本總數(shù)；p表示非t類；Yt、Yp分別表示第t類和非t類的樣本點(diǎn)；d(i，j)表示歐式空間中樣本i和j的歐式距離。

圖1為SSE隨k的變化。圖2為S隨k的變化。由圖1可見，當(dāng)k=3和k=5時(shí)，曲線的斜率變化較為明顯，但是均無顯著的肘部位置，需要結(jié)合平均輪廓系數(shù)法的結(jié)果來確定kopt。由圖2可見，k=9時(shí)所對(duì)應(yīng)的S最大。但當(dāng)k=9時(shí)，圖1的曲線斜率變化趨于平緩，說明k值已經(jīng)遠(yuǎn)遠(yuǎn)偏離kopt。綜合手肘法和輪廓系數(shù)法的分析結(jié)果，k=3時(shí)，圖1曲線的斜率變化較為明顯，圖2中對(duì)應(yīng)的S也較大，因此最終取kopt=3。

圖1 SSE隨k的變化

圖2 S隨k的變化

表2為kopt=3時(shí)每一類原料油聚類中心的原料性質(zhì)。由表2可以看出，聚類中心相互之間區(qū)分度較大的原料性質(zhì)是硫含量、飽和烴含量、烯烴含量和溴值，而溴值是表示有機(jī)化合物中不飽和程度的一種指標(biāo)，與烯烴含量成正相關(guān)關(guān)系，因此本研究選取硫含量、飽和烴含量和烯烴含量的高低來度量不同類別原料油性質(zhì)的差異。第1類原料油的烯烴含量最低，硫和飽和烴含量較高，因此可命名為低烯烴高硫高飽和烴原料油；第2類原料油的硫和飽和烴含量最低，烯烴含量最高，因此命名為低硫低飽和烴高烯烴原料油；第3類原料油的硫和飽和烴含量最高，烯烴含量較低，因此命名為高硫高飽和烴低烯烴原料油。

表2 kopt=3時(shí)每一類聚類中心(K-means)的原料性質(zhì)

圖3為溴值與芳烴含量之間的二維關(guān)系。由圖3可見：3類原料油數(shù)據(jù)點(diǎn)分布集中，分布區(qū)域相對(duì)獨(dú)立且具有較明顯界限，同時(shí)聚類中心分布也較分散；溴值在3類原料油中的分布存在較為明顯的梯度，可見溴值對(duì)最終聚類結(jié)果影響較大；芳烴含量從低到高在3類原料油中均有數(shù)據(jù)點(diǎn)分布，可見芳烴含量對(duì)最終的聚類結(jié)果影響不大。

圖3 K-means算法溴值與芳烴含量之間的二維關(guān)系●—第1類原料油； ▲—第2類原料油； ■—第3類原料油。較大標(biāo)記為聚類中心。圖4～圖8同

圖4為密度與烯烴含量之間的二維關(guān)系。由圖4可見：3類原料油數(shù)據(jù)點(diǎn)分布集中，分布區(qū)域相對(duì)獨(dú)立且具有較明顯界限，同時(shí)聚類中心分布也較分散；烯烴在3類原料油中的分布存在明顯的梯度，說明烯烴含量對(duì)最終的聚類結(jié)果影響較大；密度從低到高在3類原料油都有數(shù)據(jù)點(diǎn)，可見密度對(duì)最終的聚類結(jié)果影響不大。

圖4 K-means算法密度與烯烴含量之間的二維關(guān)系

圖5為烯烴含量、飽和烴含量與硫含量之間的三維關(guān)系。由圖5可見：與二維圖相比，三維圖中各類原料油數(shù)據(jù)點(diǎn)之間的界限更為明顯，數(shù)據(jù)點(diǎn)分布更為集中，聚類中心更加分散；烯烴分布結(jié)果與圖4一致；飽和烴含量與硫含量在3類原料油中的分布存在明顯的梯度，說明飽和烴含量與硫含量對(duì)最終的聚類結(jié)果影響均較大。

圖5 K-means算法烯烴含量、飽和烴含量與硫含量之間的三維關(guān)系

2.2 FCM算法聚類分析

傳統(tǒng)的聚類方法(如K-means算法)是一種非此即彼的、清楚的分類算法，但是自然界中有許多概念是模糊的，例如大小、老少、冷暖等，它們不能用簡(jiǎn)單的是或否去表示。模糊聚類就是用以表示這種具有模糊屬性的數(shù)據(jù)集合。它將是或否，也就是集合{0，1}擴(kuò)充為[0，1]，即將二值0和1推廣到無窮多值[13]。最常見的模糊聚類方法是FCM算法，該方法設(shè)計(jì)簡(jiǎn)單、解決問題范圍廣，并且最終可以轉(zhuǎn)化為優(yōu)化問題而借助經(jīng)典數(shù)學(xué)的分線性規(guī)劃理論求解，且易于實(shí)現(xiàn)，因此，F(xiàn)CM算法得到了廣泛的應(yīng)用[8]。FCM算法的具體步驟為：

(1)初始化參數(shù)，包括模糊參數(shù)、聚類數(shù)目、隸屬度矩陣和終止條件。

(2)計(jì)算聚類中心ci。

(8)

式中：n為數(shù)據(jù)元素的個(gè)數(shù)；ci為第i類的聚類中心；μ為某數(shù)據(jù)元素對(duì)某一組的隸屬度；m為模糊加權(quán)指數(shù)；xj為第j個(gè)數(shù)據(jù)點(diǎn)。

(3)按式 (9) 計(jì)算目標(biāo)函數(shù)的值，如果目標(biāo)函數(shù)滿足終止條件則終止計(jì)算，否則繼續(xù)執(zhí)行步驟(4)。

(9)

式中：J為目標(biāo)函數(shù)；dik為第i個(gè)數(shù)據(jù)元素到第k個(gè)聚類中心的歐幾里得距離。

(4)按式 (10) 重新計(jì)算隸屬度矩陣，然后返回到步驟(2)。

(10)

最小聚類數(shù)cmin=2，最大聚類數(shù)也采用式(3)進(jìn)行計(jì)算，得cmax=11。FCM算法的聚類參數(shù)設(shè)置值見表3。

表3 FCM算法主要參數(shù)

采用吳成茂等[14]提出的歸一劃分系數(shù)式(11)來確定最佳聚類數(shù)目copt。該準(zhǔn)則是在Bezdek準(zhǔn)則[式(12)]基礎(chǔ)上得到的。

(11)

(12)

式中：U為隸屬度矩陣；F(U；c)為劃分系數(shù)；GF(U；c)為歸一劃分系數(shù)。

按c從2到11運(yùn)行MATLAB程序，得到相應(yīng)的GF(U；c)的值以及收斂迭代次數(shù)，結(jié)果見表4。由表4可以發(fā)現(xiàn)：隨c的增加，程序收斂次數(shù)也不斷增加；當(dāng)GF(U；c)取最大值時(shí)，其所對(duì)應(yīng)的c值即為最佳聚類數(shù)目，copt=3。

表4 FCM算法聚類結(jié)果

分別計(jì)算每個(gè)數(shù)據(jù)樣本對(duì)每一類原料油的隸屬度，結(jié)果見表5。以樣本1為例，它對(duì)第1，2，3類原料油的隸屬度分別為0.472 2，0.161 6，0.366 2，隸屬度之和為1，其中，對(duì)第1類原料油的隸屬度最大，因此樣本1被歸為第1類原料油。

表5 原料油樣本對(duì)于不同類別的隸屬度

表6為copt=3時(shí)每一類原料油聚類中心的原料性質(zhì)，與表2 中K-means聚類算法得到的聚類中心相比雖然略有差別，但是總體類似。同樣地，第1類原料油可命名為低烯烴高硫高飽和烴原料油；第2類原料油可命名為低硫低飽和烴高烯烴原料油；第3類原料油可命名為高硫高飽和烴低烯烴原料油。

表6 copt=3時(shí)每一類聚類中心(FCM)的原料性質(zhì)

圖6～圖8分別為FCM算法的溴值與芳烴含量，密度與烯烴含量，以及烯烴含量、飽和烴含量與硫含量的分布。由圖6～圖8可見，3類原料油數(shù)據(jù)點(diǎn)分布集中，分布區(qū)域相對(duì)獨(dú)立且具有較明顯界限，同時(shí)聚類中心分布也比較分散。溴值、烯烴含量、飽和烴含量和硫含量在3類原料油中的分布均存在較為明顯的梯度，可見溴值、烯烴含量、飽和烴含量和硫含量對(duì)最終聚類結(jié)果影響較大；芳烴含量和密度從低值到高值在3類原料油中均有數(shù)據(jù)點(diǎn)分布，可見芳烴含量和密度對(duì)最終的聚類結(jié)果影響不大。該結(jié)果與K-means算法的聚類結(jié)果一致。

圖6 FCM算法溴值與芳烴含量之間的二維關(guān)系

圖7 FCM算法密度與烯烴含量之間的二維關(guān)系

圖8 FCM算法烯烴含量、飽和烴含量與硫含量之間的三維關(guān)系

2.3 兩種聚類算法的結(jié)果比較

對(duì)比表2和表6列出的K-means和FCM聚類算法在最佳聚類數(shù)目時(shí)各類聚類中心的原料性質(zhì)可以看出，不同算法所對(duì)應(yīng)的每一類聚類中心重合度均較高。K-means和FCM聚類算法均將117組樣本分為3類，K-means算法的第1至第3類原料油的樣本數(shù)分別為38，38，41；FCM算法的第1至第3類原料油的樣本數(shù)分別為43，38，36。兩種聚類算法僅在5組樣本的分類結(jié)果上有差異，其余112組樣本的分類結(jié)果完全一致，這進(jìn)一步說明了兩種聚類算法結(jié)果的一致性。

3 結(jié) 論

以某S Zorb裝置的原料油性質(zhì)數(shù)據(jù)為基礎(chǔ)，通過對(duì)工藝和反應(yīng)機(jī)理的分析，選取原料油性質(zhì)中的硫含量、飽和烴含量、烯烴含量、芳烴含量、密度和溴值6個(gè)聚類變量，采用MATLAB編程平臺(tái)，分別使用K-means和FCM聚類算法進(jìn)行原料油聚類研究。結(jié)果表明，兩種聚類法均將原料油分為邊界比較明顯的3類，分類結(jié)果高度重合，聚類分析效果均較好，為針對(duì)每一類原料油建立相應(yīng)的汽油辛烷值預(yù)測(cè)模型，從而為尋找到使汽油辛烷值損失最小的操作條件打下良好基礎(chǔ)。