任小甜， 褚小立， 田松柏

(中國石化 石油化工科學(xué)研究院，北京 100083)

減壓餾分油(VGO)是目前重油加工最主要的原料之一。一般而言，VGO中含有各種類型的硫化物，其中噻吩類硫化物主要有苯并噻吩、二苯并噻吩、萘苯并噻吩等。目前，VGO中硫化物的組成分布主要采用GC-MS方法分析，操作復(fù)雜、時間長，造成生產(chǎn)工藝參數(shù)調(diào)整的延遲。實現(xiàn)VGO中不同類型噻吩硫化物組成分布的及時預(yù)測分析，有助于VGO加工工藝中各項參數(shù)的實時優(yōu)化和調(diào)整。相對而言，VGO的基本物性數(shù)據(jù)易于進行實時分析，且在生產(chǎn)過程中也積累了大量的物性分析數(shù)據(jù)，因此，構(gòu)建通過VGO基本物性數(shù)據(jù)預(yù)測其噻吩類硫化物組成的模型，從而實現(xiàn)VGO中噻吩硫化物的組成分布的實時預(yù)測具有重要的意義。

近年來，已有一些基于石油餾分基本物性數(shù)據(jù)計算和預(yù)測其烴類組成的方法，基本思路是由餾分油基本物性數(shù)據(jù)計算和轉(zhuǎn)換得到可以區(qū)分不同烴族的特征參數(shù)，然后以特征參數(shù)構(gòu)建線性回歸模型并進行預(yù)測。Riazi等[1-3]利用石油餾分的相對分子質(zhì)量、折光指數(shù)、黏度以及密度等物性構(gòu)建了不同的烴特征參數(shù)，并用其平均值代表對應(yīng)的烴族，建立各物性的線性方程組，通過求解方程組來計算其中烷烴、環(huán)烷烴、芳烴以及硫化物的含量。邢波等[4]根據(jù)VGO的常規(guī)物性構(gòu)建了5個新的特征參數(shù)，并選取80組VGO樣本進行關(guān)聯(lián)計算，通過求解帶約束的最小二乘規(guī)劃問題獲得各項特征參數(shù)的系數(shù)，確定了VGO中鏈烷烴，單環(huán)和多環(huán)環(huán)烷烴，單環(huán)、雙環(huán)和多環(huán)芳香烴含量的計算關(guān)聯(lián)式。

以上方法需要先確定能區(qū)分不同烴族的特征參數(shù)，其本質(zhì)是消除原始物性變量間的多元共線性，因此只能用于族組成種類較少的情況，對于詳細的烴組成則找不到相應(yīng)的特征參數(shù)。因此，除構(gòu)造有物理意義的特征參數(shù)以外，還可將基本物性兩兩組合來構(gòu)造新的特征參數(shù)。劉四斌等[5]收集了27個VGO樣品的基本物性和族組成數(shù)據(jù)，并由其中5個基本物性(密度、折光率、相對分子質(zhì)量、硫含量和氮含量)兩兩組合構(gòu)建了15個回歸因子，通過逐步回歸法從中選擇顯著變量，分別建立預(yù)測13種烴類化合物組成的線性回歸模型，其中總噻吩模型的預(yù)測精度最高，其顯著變量為硫含量和折光率的乘積，表明VGO的總噻吩含量與其硫含量和折光率有較強的相關(guān)性。

除了線性回歸模型外，利用非線性回歸模型也可以直接擬合VGO族組成與其基本物性之間的關(guān)系。孫仁金等[6]利用動量-自適應(yīng)學(xué)習(xí)率的BP人工神經(jīng)網(wǎng)絡(luò)(ANN)構(gòu)建了由VGO的平均沸點、密度、相對分子質(zhì)量和折光率預(yù)測其飽和烴含量的非線性預(yù)測模型；并進行了深入研究，添加VGO的運動黏度作為輸入特征，建立了預(yù)測其飽和烴含量的GA-ANN模型；此外，利用遺傳算法(GA)確定隱含層神經(jīng)元數(shù)量和學(xué)習(xí)率，可以提高模型預(yù)測的準(zhǔn)確度[7-8]。Stratiev等[9-10]則建立了由VGO的平均沸點、密度和結(jié)構(gòu)特征參數(shù)(由經(jīng)驗公式計算得到)分別預(yù)測其飽和烴加單環(huán)芳烴的總含量、多環(huán)芳烴含量、總芳烴含量的多項式回歸模型。非線性回歸模型不需要構(gòu)造新的物性特征參數(shù)，適用范圍更廣，預(yù)測準(zhǔn)確度更高。

目前，研究主要集中在VGO烴族組成的預(yù)測，還未見預(yù)測VGO中不同類型噻吩類硫化物組成的報道。因此，筆者從VGO的基本物性出發(fā)，采用隨機森林回歸算法構(gòu)建模型，分別預(yù)測其中苯并噻吩、二苯并噻吩、萘苯并噻吩以及總噻吩的質(zhì)量分數(shù)。

1 模型的構(gòu)建

1.1 算法原理

用RFR建立的模型訓(xùn)練時間短，不需要進行特征數(shù)據(jù)的預(yù)處理，而且還可以計算出每個特征的重要程度用于特征選擇；同時，模型對離群的異常樣本不敏感，穩(wěn)健性好，不易過擬合，有較強的泛化能力，能保證較高的準(zhǔn)確度。

1.2 構(gòu)建數(shù)據(jù)集

收集160個直餾VGO樣本，切割自不同基屬的國內(nèi)外原油(包括塔河原油、勝利原油、科威特原油和俄羅斯原油等)。測定其各項基本物性數(shù)據(jù)，包括20 ℃密度(ρ20)，碳、氫、硫元素的質(zhì)量分數(shù)(w(C)、w(H)、w(S))，70 ℃折光率n70，平均相對分子質(zhì)量M，體積平均沸點TV，40 ℃、80 ℃和100 ℃下的運動黏度(v40、v80、v100)，黏度指數(shù)VI和相關(guān)指數(shù)BMCI。依照標(biāo)準(zhǔn)方法ASTM D3239測定VGO中不同類型噻吩類硫化物(包括苯并噻吩、二苯并噻吩和萘苯并噻吩硫物)的含量，并計算得到總噻吩的含量。匯總這160個樣本的基本物性和噻吩硫化物的組成數(shù)據(jù)，構(gòu)造數(shù)據(jù)集。對數(shù)據(jù)集進行隨機劃分，其中75%的樣本(120個)作為訓(xùn)練集，用于隨機森林回歸模型的訓(xùn)練和調(diào)參；其余25%的樣本(40個)作為測試集，用來評估模型的準(zhǔn)確性和泛化能力。

1.3 模型構(gòu)建

利用Python進行模型的編寫，其中調(diào)用開源的機器學(xué)習(xí)庫Scikit-Learn[12]實現(xiàn)隨機森林回歸算法。RFR算法屬于非參數(shù)模型，其對樣本數(shù)據(jù)的總體分布沒有先驗的假設(shè)，直接通過算法擬合待測性質(zhì)與變量之間的非線性關(guān)系。構(gòu)建模型時，先對訓(xùn)練集樣本進行有放回的抽樣得到若干個樣本子集，然后在每個樣本子集上構(gòu)造回歸樹；對樹中的每個非葉節(jié)點進行分裂時，隨機選擇若干個特征數(shù)進行計算。當(dāng)每個子集都構(gòu)造起相應(yīng)的回歸樹時，即可建立起不同特征下相應(yīng)的劃分規(guī)則和取值，整個隨機森林回歸模型也就構(gòu)建完成。

將VGO的12項基本物性作為模型輸入特征X， 3種噻吩硫化物和總噻吩的質(zhì)量分數(shù)依次作為輸出變量y，構(gòu)建4種隨機森林回歸預(yù)測模型，分別為：苯并噻吩質(zhì)量分數(shù)的預(yù)測模型M1、二苯并噻吩質(zhì)量分數(shù)的預(yù)測模型M2、萘苯并噻吩質(zhì)量分數(shù)的預(yù)測模型M3和總噻吩質(zhì)量分數(shù)的預(yù)測模型M4。然后，確定4種模型的超參數(shù)，進行模型的訓(xùn)練和調(diào)參，并評估其預(yù)測的準(zhǔn)確性和泛化能力。

1.4 超參數(shù)尋優(yōu)

在對隨機森林回歸模型進行訓(xùn)練之前先要確定模型的超參數(shù)。有2種超參數(shù)需要尋優(yōu)，第一種是模型的框架參數(shù)，即抽樣得到的樣本子集個數(shù)，也就是構(gòu)造回歸樹的數(shù)量NT。一般來說，NT太小則模型訓(xùn)練不足而容易欠擬合，預(yù)測準(zhǔn)確度較差；NT太大則模型預(yù)測的準(zhǔn)確度又會降低，計算量也會增加。因此，需要選擇一個適中的數(shù)值。第二種需要尋優(yōu)的參數(shù)是回歸樹參數(shù)，包括樹的最大深度、葉節(jié)點中的最少樣本數(shù)、葉子節(jié)點的最高數(shù)量、每個非葉節(jié)點劃分使用的最大特征數(shù)MF等。由于本研究的樣本數(shù)量和特征數(shù)量都較小，只需要對回歸樹的每個非葉節(jié)點劃分使用的最大特征數(shù)MF進行調(diào)參，其他回歸樹的參數(shù)均設(shè)為默認值即可。本研究以訓(xùn)練集樣本的袋外估計為基礎(chǔ)進行調(diào)參，省去了繁瑣的交叉驗證計算，先對回歸樹的數(shù)量NT進行尋優(yōu)，然后再調(diào)節(jié)最大特征數(shù)MF。選擇不同的超參數(shù)，設(shè)定其取值范圍，然后計算對應(yīng)模型中每個袋外樣本的預(yù)測值，以所有袋外樣本的預(yù)測標(biāo)準(zhǔn)偏差(RMSE_OOB)為評價指標(biāo)，當(dāng)其取值最小時對應(yīng)的超參數(shù)即為最優(yōu)。

2 結(jié)果與討論

2.1 模型的超參數(shù)

對M1、M2、M3和M4依次進行超參數(shù)的尋優(yōu)計算，分別得到最優(yōu)的回歸數(shù)的數(shù)量NT和最大特征數(shù)MF。苯并噻吩質(zhì)量分數(shù)預(yù)測模型M1的調(diào)參過程如圖1所示。由圖1可知：M1袋外樣本的預(yù)測標(biāo)準(zhǔn)偏差(RMSE_OOB)隨著回歸樹數(shù)量(NT)的增加先減小后增大，當(dāng)NT為160時其取得最小值(圖1(a))；且RMSE_OOB隨著最大特征數(shù)(MF)的增大而減小，當(dāng)MF為12時其取得最小值(圖1(b))。由此，可以確定該模型的最優(yōu)超參數(shù)NT和MF分別為160和12。類似地，依次對其他3種模型進行調(diào)參，匯總4種模型的最優(yōu)超參數(shù)得到：二苯并噻吩質(zhì)量分數(shù)預(yù)測模型M2的最優(yōu)超參數(shù)為：NT=60，MF=6；萘苯并噻吩質(zhì)量分數(shù)預(yù)測模型M3的最優(yōu)超參數(shù)為：NT=130，MF=10；總噻吩質(zhì)量分數(shù)預(yù)測模型M4的最優(yōu)超參數(shù)為：NT=285，MF=12。

2.2 模型的訓(xùn)練

確定了各模型的最優(yōu)超參數(shù)之后，分別在對應(yīng)的120個訓(xùn)練集樣本上進行模型的訓(xùn)練計算，利用訓(xùn)練集樣本的校正標(biāo)準(zhǔn)偏差(RMSEC)和擬合決定系數(shù)R2來評價模型的準(zhǔn)確性和擬合效果，由袋外得分可初步評價模型的泛化能力。各模型中訓(xùn)練集樣本的噻吩硫化物的組成分布范圍和模型評價指標(biāo)見表1。由表1可知，4種模型中噻吩硫化物的組成分布范圍較寬，模型具有較好的代表性；利用隨機森林回歸算法可以實現(xiàn)對直餾VGO中3種噻吩硫化合物和總噻吩質(zhì)量分數(shù)的準(zhǔn)確預(yù)測，對于噻吩硫化物質(zhì)量分數(shù)為0的樣本，模型也能給出較準(zhǔn)確的預(yù)測。這4種模型中，總噻吩質(zhì)量分數(shù)預(yù)測模型M4的擬合效果最好，其擬合決定系數(shù)R2可達到0.993，萘苯并噻吩模型M3的擬合效果較差，擬合決定系數(shù)R2為0.976；RMSEC值的大小和噻吩硫化物質(zhì)量分數(shù)的大小呈正相關(guān)關(guān)系，因而總噻吩質(zhì)量分數(shù)預(yù)測模型的RMSEC值最大，但各模型的RMSEC值均小于標(biāo)準(zhǔn)分析方法(ASTM D3239)中的重復(fù)性要求，說明模型的擬合效果都較好。

圖1 苯并噻吩質(zhì)量分數(shù)預(yù)測模型的超參數(shù)優(yōu)化Fig.1 Optimization of hyper-parameters of the prediction model of benzothiophenes mass fractions(a) Number of trees in the forest； (b) Max features used for the split

從袋外得分OOB_score來看：總噻吩質(zhì)量分數(shù)預(yù)測模型M4的得分最高，達到0.948，說明其具有較好的泛化能力，能對訓(xùn)練集之外的樣本作出較準(zhǔn)確的預(yù)測；而萘苯并噻吩質(zhì)量分數(shù)預(yù)測模型M3的得分最低，只有0.826，同樣其對應(yīng)的訓(xùn)練集樣本的擬合效果也最差(R2=0.976)，說明該模型的準(zhǔn)確性和泛化能力都較低。

VGO中總噻吩質(zhì)量分數(shù)的線性回歸預(yù)測模型[5]的預(yù)測校正標(biāo)準(zhǔn)偏差RMSEC為0.580%，擬合決定系數(shù)R2為0.984。與之相比，本研究構(gòu)建的總噻吩質(zhì)量分數(shù)RFR預(yù)測模型的準(zhǔn)確性更高，RMSEC為0.264%。而且，其他3種模型預(yù)測3種不同類型噻吩硫化物的質(zhì)量分數(shù)的結(jié)果也較準(zhǔn)確。

表1 各模型的噻吩硫化物組成分布范圍和評價指標(biāo)Table 1 The distribution range and evaluation index ofthiophene sulfides composition of each model

2.3 特征重要度分析

隨機森林回歸算法可以計算出模型中每個特征的重要程度。具體來說，在構(gòu)造1顆回歸樹時，非葉節(jié)點按照分裂后的各個樣本子集方差最小的準(zhǔn)則進行分裂。這樣可以計算出每個特征在某一節(jié)點處分裂前后的方差減少量，進而得到每個特征在所有回歸樹的相應(yīng)節(jié)點上的平均方差減少量。將所有特征的方差減少量進行歸一化計算，即得到每個特征的重要度。重要度表示每個特征在模型中的貢獻值大小，可用于模型的特征篩選，即可從所有的12項物性特征中選擇出若干個重要度較大的特征作為相應(yīng)模型的主要特征，以這些特征為新的輸入可以得到更簡單的預(yù)測模型。

圖2為各模型對應(yīng)的特征重要度分析結(jié)果。由圖2可知，這4種模型中硫質(zhì)量分數(shù)特征的重要度最高，碳質(zhì)量分數(shù)特征次之，而其他物性特征的重要度都很低，說明這些物性指標(biāo)和噻吩硫化物組成的相關(guān)性很小。所以，各模型的重要特征均為硫質(zhì)量分數(shù)(w(S))和碳質(zhì)量分數(shù)(w(C))。在實際應(yīng)用中，為了構(gòu)建更加簡單和輕量的VGO中噻吩硫化物質(zhì)量分數(shù)的預(yù)測模型，可以只選擇硫質(zhì)量分數(shù)和碳質(zhì)量分數(shù)作為輸入特征來構(gòu)建模型，也可以達到一定的預(yù)測準(zhǔn)確度。

此外，由于減壓餾分油中的硫基本上都分布在噻吩硫化物中，所以各模型中硫質(zhì)量分數(shù)特征的重要度都最高；同時，VGO的各種噻吩硫化物中，苯并噻吩的含量最高，其占總噻吩含量的比例也最高；因此，對于苯并噻吩和總噻吩的2種模型(M1和M4)，w(S)的特征重要度要遠高于w(C)；而二苯并噻吩及萘苯并噻吩的預(yù)測模型(M2和M3)中w(C)的特征重要度則相對較大，因為與苯并噻吩相比，這2種噻吩硫化物分子中分別增加6和12個C(分別為苯并噻吩分子中增加了1個和2個苯環(huán))；其結(jié)構(gòu)中C的比例明顯增大，使其硫化物質(zhì)量分數(shù)與VGO中碳質(zhì)量分數(shù)的相關(guān)性更大。

圖2 各模型的特征重要度分布Fig.2 The feature importance distribution of each model(a) M1; (b) M2; (c) M3; (d) M4

2.4 測試集樣本的預(yù)測分析

將測試集的40個樣本數(shù)據(jù)分別代入這4種模型進行預(yù)測，利用測試樣本的預(yù)測標(biāo)準(zhǔn)偏差RMSEP和擬合決定系數(shù)R2來評價各模型預(yù)測效果和泛化能力，并主要通過比較R2的大小來對比4種模型的預(yù)測準(zhǔn)確性。測試集樣本的3種噻吩化合物各自的質(zhì)量分數(shù)及總噻吩質(zhì)量分數(shù)的實測值和預(yù)測值對比如圖3所示。從圖3可知：苯并噻吩質(zhì)量分數(shù)、二苯并噻吩質(zhì)量分數(shù)和總噻吩質(zhì)量分數(shù)的預(yù)測值和實測值基本吻合，表明模型M1、M2和M4的預(yù)測效果很好；萘苯并噻吩質(zhì)量分數(shù)的預(yù)測值和實測值吻合度不好，有5個樣本的預(yù)測偏差較大，說明模型M3的預(yù)測效果較差。同時，萘苯并噻吩模型的擬合決定系數(shù)最小，為R2=0.925，說明其擬合效果較差，泛化能力較弱；而苯并噻吩模型、二苯并噻吩模型和總噻吩模型的R2分別為0.963、0.988和0.981，說明苯并噻吩質(zhì)量分數(shù)、二苯并噻吩質(zhì)量分數(shù)和總噻吩質(zhì)量分數(shù)預(yù)測模型的預(yù)測效果較好、泛化能力強。

分析萘苯并噻吩質(zhì)量分數(shù)預(yù)測模型M3泛化能力較差的原因：由于隨機森林回歸算法不易產(chǎn)生過擬合現(xiàn)象，同時模型M3的訓(xùn)練集樣本的擬合決定系數(shù)R2均小于其他3種模型，說明模型M3的自變量并不能很好地解釋因變量的變化，即目前所選的12項物性特征不能充分反映VGO中萘苯并噻吩的組成信息，導(dǎo)致算法從現(xiàn)有數(shù)據(jù)集中學(xué)習(xí)得到的萘苯并噻吩質(zhì)量分數(shù)和物性特征的映射關(guān)系不夠準(zhǔn)確，因此，對于測試集樣本，其預(yù)測準(zhǔn)確度較低，模型的泛化能力也較差。

選擇模型的最重要特征硫質(zhì)量分數(shù)(w(S))進行

圖3 測試集樣本中噻吩硫化物組成的實測值和預(yù)測值對比Fig.3 Comparison of the measured and predicted values of thiophene sulfides’composition of the test set samples(a) Benzothiophenes; (b) Dibenzothiophenes; (c) Naphthathiophenes; (d) Total thiophenes

相關(guān)性可視分析，訓(xùn)練集樣本中萘苯并噻吩質(zhì)量分數(shù)與硫質(zhì)量分數(shù)的關(guān)系如圖4所示。由圖4可知，VGO的萘苯并噻吩質(zhì)量分數(shù)和(w(S))之間沒有明顯的相關(guān)關(guān)系，且存在硫質(zhì)量分數(shù)相差較大的樣本中萘苯并噻吩質(zhì)量分數(shù)相同的情況。進一步計算VGO中萘苯并噻吩質(zhì)量分數(shù)與12項物性特征(ρ20、w(C)、w(H)、w(S)、n70、M、TV、v40、v80、v100、VI和BMCI)的Spearman相關(guān)系數(shù)，結(jié)果分別為：0.58、-0.30、-0.69、0.80、0.62、-0.13、0.11、0.36、0.41、0.38、-0.32 和0.59。其中，w(S)的相關(guān)系數(shù)最高，為0.80；其他特征相關(guān)系數(shù)的絕對值均低于w(S)的。這說明VGO中萘苯并噻吩質(zhì)量分數(shù)與上述12項物性特征的相關(guān)性都較差。因此，萘苯并噻吩質(zhì)量分數(shù)預(yù)測模型的優(yōu)化，需要篩選添加與其相關(guān)性好的VGO物性特征來優(yōu)化預(yù)測模型，提高預(yù)測模型的準(zhǔn)確度和泛化能力。

12個特征預(yù)測模型的預(yù)測評價結(jié)果列于表2。由表2可以看出，模型M1、M2、M3和M4的預(yù)測標(biāo)準(zhǔn)偏差(RMSEP)分別為：0.268%、0.131%、0.111%和0.385%。各模型的預(yù)測RMSEP值與其訓(xùn)練集樣本中噻吩硫化物質(zhì)量分數(shù)為正相關(guān)關(guān)系，即預(yù)測效果最差的萘苯并噻吩模型M3的RMSEP值最小，而預(yù)測最準(zhǔn)確的總噻吩模型M4的RMSEP值則最大。3種噻吩硫化物預(yù)測模型的RMSEP值均小于標(biāo)準(zhǔn)分析方法(ASTM D3239)中數(shù)據(jù)的重復(fù)性要求：苯并噻吩0.8%、二苯并噻吩0.3%、萘苯并噻吩0.3%。這表明4種模型均能滿足生產(chǎn)過程中快速分析和過程分析的要求。

圖4 訓(xùn)練集樣本的苯并噻吩質(zhì)量分數(shù)和硫質(zhì)量分數(shù)的散點圖Fig.4 Scatter plot between the mass fraction ofnapathabenzthiophenes and sulfur of the train set samples

2.5 簡化特征模型的構(gòu)建與預(yù)測結(jié)果分析

為了提高模型的適用性，以硫質(zhì)量分數(shù)(w(S))和碳質(zhì)量分數(shù)(w(C))為輸入特征重新構(gòu)建模型，并對比不同特征數(shù)量模型的預(yù)測效果，結(jié)果如表2所示。由表2可以看出，對于3種噻吩硫化物和總噻吩的質(zhì)量分數(shù)，只含2個特征預(yù)測模型的各項預(yù)測指標(biāo)(袋外得分OOB_score、預(yù)測決定系數(shù)R2和預(yù)測標(biāo)準(zhǔn)偏差RMSEP)均略差于包含12個特征的預(yù)測模型，但也可以達到較高的預(yù)測精度，并且其RMSEP值同樣均小于標(biāo)準(zhǔn)分析方法(ASTM D3239)中的數(shù)據(jù)重復(fù)性要求。

表2 不同數(shù)量特征的模型的評價指標(biāo)Table 2 The evaluation index of the models with different number of features

上述結(jié)果表明，本研究構(gòu)建的12特征預(yù)測模型和2特征預(yù)測模型都具有較強的泛化能力，能夠準(zhǔn)確預(yù)測訓(xùn)練集之外的樣本。其中，2特征預(yù)測模型更加簡便，僅由VGO的硫質(zhì)量分數(shù)和碳質(zhì)量分數(shù)就可以較準(zhǔn)確預(yù)測3種噻吩硫化物和總噻吩的質(zhì)量分數(shù)，實用性更強。

3 結(jié) 論

利用隨機森林回歸方法分別構(gòu)建VGO中3種噻吩硫化物和總噻吩的質(zhì)量分數(shù)預(yù)測模型，由直餾VGO的12項基本物性數(shù)據(jù)可以快速計算得到VGO中苯并噻吩、二苯并噻吩、萘苯并噻吩和總噻吩的質(zhì)量分數(shù)，計算結(jié)果表明4種模型預(yù)測的準(zhǔn)確性和預(yù)測結(jié)果的重復(fù)性較好，達到標(biāo)準(zhǔn)方法ASTM D3239的數(shù)據(jù)重復(fù)性要求，具有較強的泛化能力。

根據(jù)特征重要性的計算結(jié)果，選擇VGO中硫質(zhì)量分數(shù)和碳質(zhì)量分數(shù)為輸入特征構(gòu)建了簡化的預(yù)測模型，其預(yù)測結(jié)果也較準(zhǔn)確，滿足標(biāo)準(zhǔn)方法ASTM D3239的數(shù)據(jù)重復(fù)性要求，且更加簡便、實用性更強。