羅維平, 曹長昕

(1.武漢紡織大學(xué) 機(jī)械工程與自動化學(xué)院， 湖北 武漢 430200；2.湖北省數(shù)字化紡織裝備重點實驗室 湖北 武漢 430200)

0 引言

在工業(yè)成品汽油生產(chǎn)過程中，反應(yīng)汽油燃燒性能的最重要指標(biāo)就是辛烷值，辛烷值越高，原油利用率越高。目前國內(nèi)大部分汽油都是對原油進(jìn)行催化裂化生產(chǎn)得到，而在工業(yè)進(jìn)行催化裂化原油過程中普遍降低了辛烷值，在一定程度上造成了原油的資源浪費。在工業(yè)提純過程中，工廠各個工序只是簡單地進(jìn)行線性搭配，而機(jī)械與用料之間并不只是簡簡單單的線性關(guān)系，這種粗略的線性關(guān)系導(dǎo)致了原油辛烷值降低，增加了企業(yè)成本。部分企業(yè)通過實驗來測定調(diào)試汽油樣品辛烷值，不僅需要花費大量的時間與人力，還需要價格不菲的實驗儀器與樣本試劑，因此，需要一種準(zhǔn)確的辛烷值預(yù)測模型在為工業(yè)提純機(jī)械與用料之間的搭配建立基礎(chǔ)。

對于汽油辛烷值的預(yù)測目前主要有3種方法：第一種方法是傳統(tǒng)方法，主要對汽油從分子層面進(jìn)行研究[1]。第二種方法是利用儀器對汽油組成成分進(jìn)行分析。王拓等[2]使用偏最小二乘法提取汽油拉曼光譜特征譜段，結(jié)合研究法對辛烷值進(jìn)行定量分析，但模型整體精度偏低；孫中奇等[3]基于拉曼分析技術(shù)設(shè)計了汽油調(diào)和過程中辛烷值的測定系統(tǒng)，能有效測量汽油辛烷值，對辛烷值的預(yù)測提供一定的幫助；丁怡曼等[4]將紅外光譜法結(jié)合偏最小二乘法建立了辛烷值快速預(yù)測方法，該方法操作簡單，速度快，但預(yù)測準(zhǔn)確率不高。第三種方法是利用汽油辛烷值數(shù)據(jù)進(jìn)行模型預(yù)測。蔣偉等[5]通過構(gòu)建隨機(jī)森林模型實現(xiàn)對汽油辛烷值的預(yù)測，但只采用了單一模型進(jìn)行預(yù)測，沒有進(jìn)行多種模型的比較選優(yōu)過程；高萍等[6]通過構(gòu)建BP神經(jīng)網(wǎng)絡(luò)能夠較好地預(yù)測提純工業(yè)步驟中的辛烷值損失，但不能實現(xiàn)對最后辛烷值的預(yù)測功能。

針對上述問題，考慮到各機(jī)械與用料之間的非線性關(guān)系，以及根據(jù)人工和通過傳統(tǒng)的統(tǒng)計方法所選取的特征對模型精度有較大的影響，造成現(xiàn)有的辛烷值預(yù)測模型精度較低的現(xiàn)狀，本文提出一種基于特征工程[7]作為模型特征選擇，以及通過對Xgboost[8]和支持向量機(jī)[9]進(jìn)行比較去劣擇優(yōu)的方法，以得到精度較高的辛烷值預(yù)測模型。

1 數(shù)據(jù)來源

本文的數(shù)據(jù)來源于某一大型石化企業(yè)。該企業(yè)通過對工業(yè)生產(chǎn)中的催化裂化汽油精制脫硫裝置進(jìn)行實時跟蹤測量，每隔3 min完成1次記錄，對2 h的數(shù)據(jù)進(jìn)行取平均值的操作，獲得325個樣本的原始數(shù)據(jù)。

原始數(shù)據(jù)中的每一個樣本包含了原油的辛烷值、硫含量等7個主要組成成分，催化裂化汽油精制脫硫裝置的354個操作變量主要信息，經(jīng)過工業(yè)脫硫和降烯烴過程后的成品汽油的辛烷值與硫含量以及待生吸附劑性質(zhì)和再生吸附劑性質(zhì)共計367個變量，是目前工業(yè)中對原油進(jìn)行處理且包含裝置信息的較為全面、頗具權(quán)威的數(shù)據(jù)集。

2 數(shù)據(jù)預(yù)處理

通過對該企業(yè)所提供的數(shù)據(jù)集進(jìn)行整體觀察，發(fā)現(xiàn)數(shù)據(jù)集中存在因人工記錄失誤或者計算錯誤等而造成缺失值與異常值等問題，需要對數(shù)據(jù)進(jìn)行清理篩選，整理出一份較為完美的數(shù)據(jù)集，有助于后續(xù)的建模分析等流程。數(shù)據(jù)清洗流程如圖1所示。

2.1 數(shù)值范圍檢測

數(shù)據(jù)清洗的第一步是對數(shù)據(jù)集的缺失值進(jìn)行處理。觀察數(shù)據(jù)發(fā)現(xiàn)有不少變量的記錄值為0，但不排除變量本身的值為0。在缺失值填補之前，需要根據(jù)該石化企業(yè)通過實驗獲得的變量取值范圍，對每一個變量的數(shù)值進(jìn)行對比篩選，部分變量取值范圍見表1。對不在其變量取值范圍內(nèi)的值，采用均值填補的方法。

圖1 數(shù)據(jù)清洗流程Fig.1 Data cleaning process

表1 變量取值范圍Tab.1 Variable value range

2.2 數(shù)據(jù)篩選

范圍檢測完成之后，再次觀察數(shù)據(jù)，進(jìn)行缺失值填補。當(dāng)某一變量的缺失值過多時，無論采用任意一種填補方式，其數(shù)值都不能較好地表達(dá)變量的變化趨勢，將其視為無用的特征變量進(jìn)行刪除；當(dāng)缺失值在合理范圍之內(nèi)時，需要對其進(jìn)行填補。目前缺失值填補主要分為單變量缺失值插補和多變量缺失值填補兩大類，主流的單變量缺失值插補方法主要有眾數(shù)填補、平均值填補、中值填補以及上、下值填補，但這些方法或多或少會改變數(shù)據(jù)的原始分布，造成一定量的模型誤差。本文采用多變量缺失值填補方法對辛烷值數(shù)據(jù)的缺失值進(jìn)行填補。

隨機(jī)森林是目前較為主流的機(jī)器學(xué)習(xí)集成算法，分為分類和回歸2種用法，其主要思想為內(nèi)部使用多個弱決策樹，結(jié)合Bagging(裝袋算法)方法以及特征子空間，通過隨機(jī)抽樣的方法在建模過程中抽取樣本，通過投票的方式?jīng)Q定模型輸出結(jié)果。由于模型的處理機(jī)制特殊，隨機(jī)森林回歸算法對缺失值以及噪聲點具有較好的容錯率，不容易過擬合，因此采用隨機(jī)森林填補缺失值。

對于擁有n個特征的辛烷值數(shù)據(jù)集，其中多列存在缺失值，遍歷所有特征，從缺失值最少的開始填補，將其他特征的缺失值暫時用0代替，每完成一次模型預(yù)測，就將預(yù)測值放入到特征矩陣，再進(jìn)行下一次的缺失值預(yù)測填補，逐漸完成所有缺失值的預(yù)測填補，得到完整的數(shù)據(jù)集。

2.3 異常值檢測

完成對缺失值的填補后，通過對變量繪制核密度圖觀察數(shù)據(jù)分布情況，部分變量核密度圖如圖2所示，分別以汽油中的硫含量和飽和烴為例。

(a)硫含量 (b)飽和烴

由圖中可以看出，變量中依舊存在離群值，需要對其進(jìn)行處理。本文采用孤立森林的方式進(jìn)行離散值處理。孤立森林算法隨機(jī)遞歸分割數(shù)據(jù)集，直到所有的數(shù)據(jù)均為離群值則停止運行，而離群值相較于整體數(shù)據(jù)來說偏離數(shù)據(jù)分離點，更容易被分割，所以使用較少次數(shù)分割出來的點即為離群值。采用孤立森林算法后，硫含量分布情況如圖3所示。圖中方形點即為硫含量中的離群值，對于孤立森林捕捉到的離群值，算法中所自帶的接口可以對其索引進(jìn)行緩存幫助修改。

通過數(shù)據(jù)檢測到異常值檢測的操作，對該數(shù)據(jù)有一個較為完整的清洗流程，得到一個共計346變量的、相對于之前較好的辛烷值數(shù)據(jù)集，有利于后續(xù)的數(shù)據(jù)降維與建模分析。

圖3 硫含量分布情況Fig.3 Sulfur content distribution

3 降維

主成分分析(principal component analysis)是一種用來探索高維數(shù)據(jù)結(jié)構(gòu)的技術(shù)，一般被用來實現(xiàn)高維的數(shù)據(jù)集的探索與可視化以及用來進(jìn)行數(shù)據(jù)壓縮。

本文中的辛烷值數(shù)據(jù)集，除去標(biāo)簽共計345個特征變量，在高維數(shù)據(jù)的情況下會出現(xiàn)數(shù)據(jù)樣本稀疏，距離計算困難等問題，更容易導(dǎo)致模型的過擬合問題，被統(tǒng)稱為維度災(zāi)難。PCA的本質(zhì)就是找一些投影方向，且這些投影方向是相互正交的，使得數(shù)據(jù)在這些投影方向上的方差是最大的。方差越大，說明其在對用正交基上包含更多的信息量，證明原始數(shù)據(jù)協(xié)方差矩陣的特征值越大，對應(yīng)所包含的信息越多。

根據(jù)PCA降維的原理，可以分為以下幾個計算步驟：

①對數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理，轉(zhuǎn)化為無量綱的純數(shù)值，去除單位的限制，便于不同量級的指標(biāo)之間能夠相互比較。

(1)

②計算相關(guān)系數(shù)矩陣。

(2)

③對相關(guān)系數(shù)矩陣進(jìn)行求解，得出其特征向量和特征值。

④進(jìn)行主要成分提取。根據(jù)累積貢獻(xiàn)率進(jìn)行特征提取，主成分zi的方差貢獻(xiàn)率為

(3)

圖4 主要成分和累積解釋方差Fig.4 Main components and cumulative explanatory variance

利用PCA降維對原始辛烷值數(shù)據(jù)集進(jìn)行降維，其方差貢獻(xiàn)率如圖4所示。從圖中可以明顯看出，當(dāng)提取36個主要成分之后，剩余提取特征對原本數(shù)據(jù)的解釋能力幾乎為0，單解釋方差曲線趨向于平滑，所以決定提取36個主要成分，由原來的345個特征變量通過PCA降維，成功降到36個特征變量。

4 模型建立

本文要解決的問題為回歸預(yù)測問題，所以根據(jù)問題選擇目前較為主流的支持向量機(jī)回歸以及Xgboost2種回歸算法進(jìn)行預(yù)測對比，使用k折交叉驗證(k-fold cross validation)對模型的穩(wěn)定性及其對數(shù)據(jù)集外的數(shù)據(jù)的泛化性進(jìn)行驗證，再通過模型評估指標(biāo)MSE(平均均方誤差)以及R2(決定系數(shù))對模型進(jìn)行評估，選取最優(yōu)模型。模型構(gòu)建流程如圖5所示。

圖5 模型構(gòu)建流程Fig.5 Model building process

4.1 支持向量機(jī)

不同于上文介紹過的隨機(jī)森林回歸，支持向量回歸(SVR)是由支持向量機(jī)(SVM)從分類問題推廣到回歸問題，SVM分類是找到超平面，讓2個分類集合的支持向量或者所有數(shù)據(jù)離分類平面最遠(yuǎn)；SVR回歸不同于SVM的是使得集合內(nèi)所有數(shù)據(jù)到平面的距離最近，具體為回歸預(yù)測值f(x)與真實值y之間得偏離程度不大時，即可以認(rèn)為預(yù)測正確，損失不用計算。

關(guān)于SVR的問題可以簡化為

(4)

式中：C是正則化參數(shù)；ε為ε-的一個不敏感損失函數(shù)，

(5)

SVR能較好地解決局部最優(yōu)問題，并且有較優(yōu)的抗過擬合能力。使用該模型時，核函數(shù)的選取也是一個重要的影響因素，選擇合適的核函數(shù)能更好地擬合數(shù)據(jù)，得到更好的回歸預(yù)測值，有較低的模型損失。

4.2 Xgboost模型

Xgboost模型(eXtreme gradient boosting)是一種tree boosting的可拓展機(jī)器學(xué)習(xí)系統(tǒng)，旨在通過結(jié)合k個精度較低的回歸樹組成一個高精度的模型，使得樹群的預(yù)測值盡量接近真實值，具有較高的準(zhǔn)確率，并且有很強(qiáng)的泛化能力。其目標(biāo)函數(shù)為傳統(tǒng)的損失函數(shù)加上模型復(fù)雜度：

(6)

(7)

4.3 評價標(biāo)準(zhǔn)

對于辛烷值的預(yù)測，本文更關(guān)注辛烷值的預(yù)測值與原本的真實值之間的誤差，所以采用均方誤差(mean squared error，MSE)作為模型評價指標(biāo)，

(8)

R2也可以作為另外一種參考的指標(biāo)，

(9)

5 實驗結(jié)果與分析

5.1 實驗結(jié)果

將通過有效性處理后的數(shù)據(jù)輸入本文所選取的2種算法，再分別對其進(jìn)行10折交叉驗證，利用模型評價指標(biāo)對各模型效果進(jìn)行評估，結(jié)果如圖6所示。

由圖中數(shù)據(jù)可得，在未進(jìn)行參數(shù)調(diào)優(yōu)之前，在R2指標(biāo)上，SVR的值明顯大于Xgboost的；而在均方誤差上，SVR的值小于Xgboost的，模型的損失較小，SVR模型性能整體高于Xgboost，因此本文基于SVR構(gòu)建辛烷值預(yù)測模型。

5.2 模型改進(jìn)

SVR主要核函數(shù)有4種，分別為linear、poly、rbf、sigmoid，使用不同核函數(shù)得到的模型精度如圖7所示。由圖中數(shù)據(jù)可知，rbf核函數(shù)的R2值最大，平均誤差最小，所以選取rbf核函數(shù)作為本文模型的核函數(shù)。

SVR模型的內(nèi)置參數(shù)有2個，分別為C和γ。C，即為懲罰系數(shù)，表示對誤差的寬容度，C值越大，表示對誤差的容忍度越低，越不能容忍出現(xiàn)誤差，過擬合的風(fēng)險越大；反之，則表示容忍度越高，容易出現(xiàn)欠擬合,懲罰系數(shù)過大或者過小，模型的泛化能力都會變差。γ是選擇rbf核函數(shù)之后函數(shù)自帶的參數(shù)，隱含的表示了數(shù)據(jù)映射到新的特征空間的分布，γ越小，支持向量越多；反之，則支持向量越少。

圖6 2種模型評估指標(biāo)對比效果Fig.6 Comparison of evaluation indexes of two models

圖7 不同核函數(shù)模型精度對比Fig.7 Comparison of accuracy of different kernel function models

使用學(xué)習(xí)曲線選取最優(yōu)參數(shù)，模型學(xué)習(xí)曲線如圖8所示。

由圖8數(shù)據(jù)可知，當(dāng)參數(shù)γ、C分別取0.08、17時，模型效果達(dá)到峰值。通過調(diào)整參數(shù)后，模型精度達(dá)到84.36%，平均誤差降到0.169。結(jié)合表2可知，當(dāng)采用隨機(jī)森林結(jié)合孤立森林進(jìn)行數(shù)據(jù)篩選后，通過改進(jìn)模型得到的模型精度相對于普通數(shù)據(jù)處理方式使用未改進(jìn)模型后的精度有較大的提升。

表2 模型精度對比Tab.2 Comparison of model accuracy

圖9、10分別為調(diào)參前后模型預(yù)測效果圖，三角號為真實值，星號為預(yù)測值，兩者的覆蓋率表示模型預(yù)測的準(zhǔn)確度。相比于圖9，圖10中的預(yù)測值與真實值兩者相交范圍明顯擴(kuò)大，表明準(zhǔn)確率有明顯提升。相應(yīng)的模型調(diào)優(yōu)前，存在部分偏離正常范圍的數(shù)值點，形成孤立現(xiàn)象；調(diào)優(yōu)后，部分偏離的真實值與預(yù)測值之間產(chǎn)生交集或者兩者之間距離縮短，表明誤差減小，可以很好地預(yù)測數(shù)據(jù)。

圖9 調(diào)參前模型預(yù)測效果Fig.9 Rrediction effect of model before parameter adjustment

圖10 調(diào)參后模型預(yù)測效果Fig.10 Rrediction effect of model after parameter adjustment

5.3 模型比較

以目前國內(nèi)普遍的辛烷值預(yù)測方法紅外光譜測量法以及拉曼光譜結(jié)合偏最小二乘法作為比較對象，紅外光譜測量法即通過紅外光譜分析技術(shù)，分析已知辛烷值的汽油產(chǎn)品的特征基團(tuán)和表征結(jié)構(gòu)的近紅外光譜參數(shù)，使用多元線性回歸的方法結(jié)合馬達(dá)法、研究法以及主成分分析來預(yù)測辛烷值。比較結(jié)果見表3。

表3 模型預(yù)測誤差比較Tab.3 Comparison of model prediction error

結(jié)合表3發(fā)現(xiàn)，在最小誤差上，本文所建立的模型預(yù)測結(jié)果與馬達(dá)法以及研究法相差不大，但在最大誤差以及平均絕對誤差上均小于其他2種方法；而在模型R2的比較上拉曼光譜結(jié)合偏最小二乘法略高于主成分分析法以及本文所建立模型，但本文所建立的辛烷值預(yù)測模型的平均絕對誤差遠(yuǎn)小于另外兩者方法，可見本文所建立的模型表現(xiàn)要優(yōu)于其他4種方法。

6 結(jié)語

本文將特征工程應(yīng)用到辛烷值預(yù)測模型特征處理與選取上，不同于一般的數(shù)據(jù)處理方式，采用隨機(jī)森林與孤立森林結(jié)合的方法，對缺失值與異常值進(jìn)行有效處理，使得樣本數(shù)據(jù)更加精確合理。使用PCA降維的方式，分析特征對原始數(shù)據(jù)的解釋能力，通過特征累積解釋方差進(jìn)行特征提取，使得所選取的數(shù)據(jù)特征契合不同算法，較好地提升了模型對數(shù)據(jù)的預(yù)測能力與泛化效果，對模型的預(yù)測精度有較大的提升。通過對Xgboost與SVR的模型評估指標(biāo)對比與參數(shù)調(diào)優(yōu)，進(jìn)行棄劣留優(yōu)操作選出效果最優(yōu)模型，使得模型整體精度達(dá)到84.36%，平均誤差降到0.169。通過本文所構(gòu)建的預(yù)測模型可以對汽油辛烷值進(jìn)行更好的預(yù)測，降低工業(yè)辛烷值提純的損失，為后續(xù)工業(yè)辛烷值提純設(shè)備需求數(shù)據(jù)的完備提供可能。

由于辛烷值樣本數(shù)據(jù)量不夠完備，因此預(yù)測模型在辛烷值預(yù)測性能中未能達(dá)到理想的效果。在以后將進(jìn)一步對辛烷值預(yù)測模型進(jìn)行深入研究，提升其預(yù)測精度。