基于BP神經(jīng)網(wǎng)絡(luò)和多元線性回歸的辛烷值預(yù)測

摘要為降低硫、烯烴含量及辛烷值損失，保證汽油清潔化生產(chǎn)，基于S Zorb裝置運行積累的數(shù)據(jù)，首先利用Lasso算法初步篩選建模變量，并基于BP神經(jīng)網(wǎng)絡(luò)計算指標因子貢獻度，進一步篩選出15個主要變量用于建立辛烷值損失預(yù)測模型；其次對比分析4種模型，得出BP神經(jīng)網(wǎng)絡(luò)預(yù)測精度更優(yōu)，更適合作為辛烷值損失預(yù)測模型，并經(jīng)過10折交叉驗證得到均方誤差（MSE）均值為0.027 193，R2均值為0.904 87，驗證了該模型的可靠性；最后在控制油品硫質(zhì)量分數(shù)不大于5 μg/g的前提下，結(jié)合多元線性回歸對主要變量進行優(yōu)化調(diào)控.結(jié)果表明，需同時改變多個變量才能使辛烷值損失降幅大于30%，多元線性回歸模型預(yù)測精度較好，能按照一定比例對主要變量進行正反向調(diào)控.本文還可視化展示了優(yōu)化過程中辛烷值和硫含量的變化軌跡.

關(guān)鍵詞BP神經(jīng)網(wǎng)絡(luò)；多元線性回歸；Lasso算法；辛烷值損失預(yù)測；優(yōu)化調(diào)控

中圖分類號TP183；TP273 文獻標志碼A

0 引言

汽車的普及使用增加了人們交通出行的便捷性，但也隨之加重了環(huán)境受污染的程度.現(xiàn)如今日益嚴峻的健康和環(huán)境問題促使各國都在重新考慮汽油中各種化合物的質(zhì)量標準，清潔措施的重中之重就是降低汽油中的硫、烯烴含量.在催化裂化為核心的重油輕質(zhì)化工藝進行脫硫和降烯烴過程中，普遍會降低汽油辛烷值.辛烷值（RON）是反映交通工具所使用燃料（汽油）燃燒性能的最重要指標.一般來說，提高汽油中辛烷值的占比，將會有效提高其抵抗震爆的性能，而抗爆性能的高低是體現(xiàn)汽油燃燒性能的主要指標.而煉油生產(chǎn)過程中辛烷值的損失將導(dǎo)致油品經(jīng)濟效益的驟然下跌，每增加1個單位的精制汽油辛烷值損失，每噸汽油銷售價格將降低150元.例如在一個每年生產(chǎn)100萬t精制汽油的S Zorb裝置中，如果能夠使辛烷值損失降低0.3個單位，銷售經(jīng)濟效益將提升4 500萬元.因此，降低催化裂化汽油精制脫硫裝置中辛烷值的損失具有重要的理論價值和現(xiàn)實意義.

在現(xiàn)有的研究文獻中，針對汽油精制過程中辛烷值的研究可總結(jié)為三種：一是討論使用汽油物理特性測試數(shù)據(jù)通過線性建模分析來快速測定汽油辛烷值.丁怡曼等[1]利用紅外光譜法結(jié)合偏最小二乘法構(gòu)建PLS模型來對113個樣品進行汽油辛烷值的預(yù)測；Kardamakis等[2]利用近紅外光譜法收集了249個汽油數(shù)據(jù)樣本，構(gòu)建辛烷值和苯含量分析的定量預(yù)測模型，同時對比分析了C_PLS法和D_PLS法的模型求解效果.二是針對汽油各種組成成分或者原子團對辛烷值貢獻進行定性和定量分析，并構(gòu)建其關(guān)系之間的擬合模型.黃水望等[3]利用氣相色譜法和偏最小二乘法構(gòu)建汽油詳細組分和辛烷值之間關(guān)系的數(shù)學模型，模型測定結(jié)果與實際偏差范圍在0～1.1個單位之間，預(yù)測性能和精度較好.Ghosh等[4]利用氣相色譜技術(shù)測定不同型號汽油的辛烷值，構(gòu)建烴類貢獻度大小不同的多元線性回歸模型，明確烴類組成對汽油生產(chǎn)調(diào)和辛烷值的影響.但不管是根據(jù)汽油物理特性數(shù)據(jù)還是各組成成分或原子團進行預(yù)測分析，都不可避免地要使用到相關(guān)實驗測試設(shè)備，而且儀器的運轉(zhuǎn)和維護需要耗費較高的費用，實驗測試效率較低.三是利用主成分分析、相關(guān)分析、偏最小二乘法、逐步線性回歸等方法研究汽油的其他理化質(zhì)量指標和抗爆性能指標辛烷值之間的關(guān)系，得出相關(guān)關(guān)系式來測定不同汽油的辛烷值.熊春華等[5]以乙醇汽油和車用汽油為研究對象，分析其抗爆性能指標與其他理化指標的關(guān)聯(lián)，采用逐步線性回歸和相關(guān)分析來建立方程關(guān)系式，研究表明各理化指標產(chǎn)生的影響均不相同，需遵循相關(guān)規(guī)律來配合實際生產(chǎn).

截至目前，利用各種理論測量方法來進行化工過程建模預(yù)測汽油辛烷值的相關(guān)研究工作取得了一定的進展，但這些傳統(tǒng)研究方式一般都是簡單的數(shù)據(jù)關(guān)聯(lián)或機理建模，所構(gòu)建模型中涉及的操作變量數(shù)量相對較少，相關(guān)變量分析存在明顯缺陷，使得辛烷值的預(yù)測分析結(jié)果與實際相比存在較大的誤差.而現(xiàn)階段催化裂化精制汽油所需的生產(chǎn)設(shè)備多種多樣，加工技術(shù)及工藝過程非常復(fù)雜，各操作變量之間具有嚴重非線性和影響因素相互強耦聯(lián)的特點，若繼續(xù)采用以往的研究方式將出現(xiàn)過程優(yōu)化響應(yīng)時效性較差、變量優(yōu)化成效一般的情況.

因此，本文考慮到現(xiàn)有化工過程監(jiān)測和控制硬件設(shè)備技術(shù)的發(fā)展，將利用采集到的中國石化公司多年來催化裂化生產(chǎn)精制汽油保留下來的大量歷史數(shù)據(jù)，結(jié)合數(shù)據(jù)挖掘技術(shù)發(fā)現(xiàn)隱藏在深層的重要信息.同時，因為神經(jīng)網(wǎng)絡(luò)技術(shù)強大的函數(shù)映射能力和高度非線性描述能力等優(yōu)點，它已經(jīng)被廣泛應(yīng)用于化工過程非線性系統(tǒng)建模領(lǐng)域.本文在利用Lasso算法初步篩選變量后，選擇利用BP神經(jīng)網(wǎng)絡(luò)計算指標因子貢獻度，得出主要變量用于構(gòu)建汽油辛烷值損失預(yù)測模型.經(jīng)過模型效果分析對比決定使用BP神經(jīng)網(wǎng)絡(luò)模型來預(yù)測辛烷值，并把BP神經(jīng)網(wǎng)絡(luò)和多元線性回歸模型相結(jié)合來對主要變量進行優(yōu)化調(diào)控，對優(yōu)化調(diào)整過程中辛烷值和硫含量變化進行可視化展示，以提升汽油品質(zhì)為企業(yè)的煉油生產(chǎn)過程提供可靠的理論操作借鑒.

1 數(shù)據(jù)收集和清洗

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

本文采集到的實驗數(shù)據(jù)來自于中國石化上海高橋石油化工有限公司催化裂化汽油精制脫硫裝置運行多年保留下來的歷史數(shù)據(jù) [5-7].關(guān)于汽油精制生產(chǎn)工作，該公司建立了2個規(guī)模可觀的數(shù)據(jù)庫，即PHD和LIMS實驗數(shù)據(jù)庫.有關(guān)實驗建模的原料、產(chǎn)品和催化劑等性質(zhì)數(shù)據(jù)均可以每周2次的頻率從這2個實時數(shù)據(jù)庫中采集得到.為了確保實驗和分析等方面的有效性與準確率，本文將數(shù)據(jù)采集時間跨度范圍擴充到3年，即從LIMS數(shù)據(jù)庫中獲取2017-04—2019-09和2019-10—2020-05這2個時間段的生產(chǎn)信息.而操作變量數(shù)據(jù)可從PHD數(shù)據(jù)庫中獲取，2次數(shù)據(jù)采集時間的頻次不同，分別是每3 min 1次和每6 min 1次.采集到的原始實驗數(shù)據(jù)樣本一共有325個，每個數(shù)據(jù)樣本中都有367個特征變量，包括7個原料性質(zhì)、2個待生吸附劑性質(zhì)、2個再生吸附劑性質(zhì)、2個產(chǎn)品性質(zhì)等不可操作變量以及另外354個操作變量.為了使數(shù)據(jù)處理和分析更加系統(tǒng)化，可以按時間戳將其進行降序排列，并結(jié)合數(shù)據(jù)挖掘技術(shù)應(yīng)用到化工過程建模中，得到隱藏在這些數(shù)據(jù)之后的更多更重要的信息.

1.2 數(shù)據(jù)清洗

1）刪除全部為空值的位點.遍歷全部數(shù)據(jù)樣本，篩選出共有19個操作變量數(shù)據(jù)全部為空值，故刪除上述位點.

2）對樣本數(shù)據(jù)節(jié)點進行過濾，刪除2個由于殘缺數(shù)據(jù)過半而無法補充完整的變量位點.

3）對于樣本中其余值為零的數(shù)據(jù)，使用控制前后2 h的數(shù)據(jù)平均值替代.

4）根據(jù)汽油精制工藝要求及操作經(jīng)驗可以得出操作變量的最大最小取值范圍，對于7個超出此范圍的變量數(shù)據(jù)進行剔除.

5）運用拉依達準則（3σ準則）判斷數(shù)據(jù)是否存在粗大誤差，剔除此類誤差數(shù)據(jù).利用MATLAB計算誤差并比較發(fā)現(xiàn)樣本原始數(shù)據(jù)中不存在異常值，無需剔除.

6）由于研究需以辛烷值作為目標數(shù)據(jù)，而辛烷值的測量比較麻煩，一周僅2次，樣本較少，無法與操作數(shù)據(jù)樣本相匹配，且測量結(jié)果存在滯后，故最終樣本數(shù)據(jù)取為辛烷值測量時間前2 h的數(shù)據(jù)平均值.計算完成數(shù)據(jù)預(yù)處理后，替換原始收集到的樣本數(shù)據(jù).

2 二步法篩選建模的主要變量

2.1 基于Lasso的變量初步篩選模型

為了更有效地對工程技術(shù)應(yīng)用效果進行分析，需要根據(jù)實際催化裂化汽油精制過程所得到的325個樣本數(shù)據(jù)，先對367個變量進行降維，剔除一些次要的對辛烷值損失影響不大的標量，篩選出一定數(shù)量的主要變量，以便后續(xù)能精準有效地建立辛烷值的損失預(yù)測模型.Lasso線性模型在變量選擇方面的精度比逐步回歸法和嶺回歸法等能更精確、更全面地篩選出主要影響變量，其最大優(yōu)點在于可以直接將不重要變量的系數(shù)直接壓縮為0，而不保留所有變量[6-10].因此，采用Lasso對367個影響變量進行初步篩選，通過比較解釋變量與被解釋變量之間相關(guān)性的大小，刪除不重要變量并保留主要變量，降低影響變量之間的多重共線性，使其相互獨立，從而提高解的空間穩(wěn)定性并進一步使得模型的泛化能力增強.

2.1.1 Lasso線性回歸模型的建立

運用Lasso方法初步篩選變量的步驟如下：

1） 特征標準化.觀察367個變量發(fā)現(xiàn)其量綱并不一致，為了避免量綱對研究結(jié)果的影響，需要利用極差標準化法對數(shù)據(jù)進行標準化處理.

2） 建立Lasso線性回歸模型.由于研究樣本中有367個自變量和1個因變量，可根據(jù)線性回歸模型建立以下關(guān)系式：

式中，α0為常數(shù)項，α1，α2，…，α367為回歸系數(shù)，ε為隨機擾動項.

該模型中未知參數(shù)α的Lasso估計的定義為

其中，t≥0為調(diào)和參數(shù).

3） 選擇最佳的調(diào)整參數(shù).對于控制回歸系數(shù)壓縮量的問題，在估計時可以通過調(diào)和參數(shù)t來實現(xiàn)，經(jīng)過若干步驟之后可得不同t值下的所有Lasso估計值.在研究中，可選擇10折交叉驗證法對參數(shù)進行調(diào)整從而選擇出最佳的調(diào)整參數(shù).

4） 篩選出重要變量.對367個變量進行Lasso回歸之后，得出一定數(shù)量回歸系數(shù)不為0的變量，即為第1步所要篩選出的重要變量，其余回歸系數(shù)為0的變量就會在Lasso線性回歸模型中刪除掉，被刪除掉的可能是導(dǎo)致變量之間存在多重共線性的不重要變量.

2.1.2 Lasso模型初步篩選變量的結(jié)果分析

運用Stata軟件對辛烷值損失的367個影響因素進行初步篩選，通過lassopack命令實現(xiàn)Lasso回歸，并使用10折交叉驗證方法對模型的參數(shù)進行調(diào)整以達到最佳的狀態(tài).根據(jù)運行結(jié)果可以發(fā)現(xiàn)當Lamda = 1.373 349 8時，模型的均方預(yù)測誤差（MSPE）值達到最小，如圖1所示.

根據(jù)Lasso所估計出的變量系數(shù)是否非零來篩選變量，其中系數(shù)非零的影響變量被保留下來作為主要變量，且它們之間的多重共線性已經(jīng)得到一定的削減.經(jīng)降維篩選過后的主要變量如表1所示，按照從左到右對辛烷值損失影響程度大小進行排列.

2.2 基于BP神經(jīng)網(wǎng)絡(luò)-指標因子貢獻度排名的變量篩選優(yōu)化

2.2.1 BP神經(jīng)網(wǎng)絡(luò)-指標因子貢獻度模型的建立

采用Lasso線性回歸模型對367個變量進行初步篩選后得出n個主要影響因素，而這些因素對辛烷值損失影響程度的大小可以通過權(quán)重大小來決定.因此，需要進一步優(yōu)化辛烷值損失的有效變量操作方案.鑒于BP神經(jīng)網(wǎng)絡(luò)模型可以通過數(shù)據(jù)逼近任意線性連續(xù)的函數(shù)，這一特點與原料性質(zhì)、待生吸附劑性質(zhì)、再生吸附劑性質(zhì)、產(chǎn)品性質(zhì)、操作變量對辛烷值損失影響方式的特點相吻合.所以選擇BP神經(jīng)網(wǎng)絡(luò)模型進行n個主要變量的計算，并將其指標數(shù)值作為神經(jīng)網(wǎng)絡(luò)的輸入層，將損失的辛烷值作為輸出層.此外，輸入層主要變量對輸出層辛烷值損失的影響，是由輸入層對隱含層的影響和隱含層對輸出層的影響這2個部分組成的，如圖2所示.因此，輸入層n個主要影響變量指標對輸出層辛烷值損失的影響權(quán)重需綜合上述2個部分來進一步計算得出，過程如下：

1） 假設(shè)各輸入變量對各隱含層變量都有一定程度的影響，這些影響的程度可通過輸入層作用到隱含層的權(quán)重反映.權(quán)重計算公式為

式中：a0表示每個輸入層主要變量對隱含層的影響權(quán)重之和；aij表示輸入層n個主要變量指標對隱含層中各個變量的權(quán)重，其中i和j分別表示輸入層和隱含層的節(jié)點；∑njwij表示單個輸入層節(jié)點對隱含層中所有變量的影響大小之和，i取整數(shù).

2） 計算隱含層對輸出層辛烷值損失的影響程度.由于輸出層僅有辛烷值損失一個指標，所以隱含層對輸出層的影響權(quán)重等價于隱含層中各個變量與輸出層節(jié)點的權(quán)重之比，其計算公式為

式中：bkl為隱含層對輸出層影響權(quán)重的比例，k和l分別表示隱含層和輸出層的節(jié)點；bk為隱含層中單個變量對輸出層的權(quán)重；∑mkbk表示隱含層全部變量對輸出層權(quán)重之和，k取整數(shù).

3） 計算輸入層中n個具體主要影響變量指標對輸出層的影響權(quán)重，計算公式為

式中，Si為輸入層中單個的具體主要影響變量對輸出層的影響權(quán)重.

4） 計算輸入層中每個指標對輸出層影響大小之比例，計算公式為

式中，P0表示所有輸入層各個指標權(quán)重之和，Pi為輸入層各個具體主要影響變量指標對輸出層辛烷值損失影響的占比.

經(jīng)過對BP神經(jīng)網(wǎng)絡(luò)模型的不斷調(diào)整，最終設(shè)置n個輸入層、1個隱含層，將隱含層節(jié)點數(shù)設(shè)為10、輸出層為1，此設(shè)置方式可以使得誤差達到最低，如圖2所示.

2.2.2 BP神經(jīng)網(wǎng)絡(luò)-指標因子貢獻度排名結(jié)果分析

基于BP神經(jīng)網(wǎng)絡(luò)的指標因子貢獻度計算模型，得出基于Lasso方法初步篩選的61個有用變量對辛烷值損失的影響貢獻度.根據(jù)貢獻度的大小確定最終留下的主要變量，并用于分析辛烷值損失.而選擇的個數(shù)可以參考確定獨立篩選法SIS中的方法[11]，選取n/log n個，其中n為61.結(jié)合實際情況計算選取15個最主要影響變量，它們對辛烷值損失貢獻度排序如圖3所示，圖3中序號對應(yīng)的影響變量名稱參見表1.圖3中CD表示影響變量對辛烷值損失貢獻度大小.其中，原料辛烷值（a2）和產(chǎn)品辛烷值（a9）2個變量對辛烷值損失的影響貢獻度最大，遙遙領(lǐng)先于其他變量，而辛烷值損失值就是由這2個值計算得出的[12].由此可見，基于BP神經(jīng)網(wǎng)絡(luò)的指標因子貢獻度計算模型來反映變量對辛烷值損失的影響貢獻程度，符合工程的實際應(yīng)用，具有一定的科學性.

3 基于BP神經(jīng)網(wǎng)絡(luò)的辛烷值損失預(yù)測模型

3.1 構(gòu)建辛烷值損失預(yù)測模型

利用篩選得出的15個主要變量進行建模預(yù)測辛烷值損失，考慮使用多元線性回歸、灰色預(yù)測、隨機森林回歸、BP神經(jīng)網(wǎng)絡(luò)4種常見建模方法，計算不同方法的均方誤差（MSE）與決定系數(shù)（R2），對比不同方法的模型精度，進而選擇誤差最小的模型作為辛烷值損失預(yù)測模型.

1） 多元線性回歸的擬合過程通常利用最小二乘法來逼近，一般假設(shè)變量間存在線性關(guān)系，如式（10）所示.回歸分析時，需要計算回歸系數(shù)bk，使得計算的因變量與原始數(shù)據(jù)間的誤差最小.

2） GM（1，N）模型與GM（1，1）模型類似，區(qū)別在于輸入變量個數(shù)為N個.基于篩選得到的15個主要影響變量，則N為15，輸出變量1個，因此適用于GM（1，N）預(yù)測模型.首先對特征數(shù)據(jù)序列X（0）1和相關(guān)因素序列X（0）2，…，X（0）N進行計算，生成一次累加序列X（1）i，i=1，2，…，N，再對該序列中兩兩鄰近數(shù)取平均生成值序列Z（1）1，然后根據(jù)上述序列建立灰色微分方程：

3） 隨機森林回歸是一個用CART算法構(gòu)建的沒有剪枝的分類決策樹的集合，輸出采用單棵樹輸出結(jié)果的均值.根據(jù)每棵決策樹的權(quán)重ωi（x，θt）（t=1，2，…，k）取每棵決策樹觀測值的均值作為最終的結(jié)果.

4） BP神經(jīng)網(wǎng)絡(luò)具有較好的非線性映射能力，RON損失預(yù)測中影響因素眾多，盡管通過Lasso回歸及貢獻度排序篩選了一部分變量，但剩余15個主要變量間的關(guān)系依舊不明確，對RON損失的影響也不清晰，利用該方法的非線性映射可以較好地擬合絕大多數(shù)變量與辛烷值損失的關(guān)系.此外，神經(jīng)網(wǎng)絡(luò)的魯棒性較佳、容錯能力強，變量數(shù)據(jù)的部分缺失或者異常并不會導(dǎo)致模型產(chǎn)生誤差[13].

利用BP神經(jīng)網(wǎng)絡(luò)映射15個主要變量與辛烷值損失的關(guān)系，通過對325個樣本訓(xùn)練的擬合預(yù)測判斷模型精度，進而決定該算法是否可以用于辛烷值損失預(yù)測.利用Matlab自帶工具箱nftool對數(shù)據(jù)進行處理，把降維得出的15個主要變量作為輸入層變量，輸出層即為辛烷值損失.設(shè)置隱含層個數(shù)為10，訓(xùn)練算法選用最小二乘優(yōu)化算法（Levenberg-Marquardt）.BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)如圖4所示.

3.2 模型求解與比較分析

利用多元線性回歸、灰色預(yù)測、隨機森林回歸、BP神經(jīng)網(wǎng)絡(luò)4種方法建立辛烷值預(yù)測模型，使用Matlab軟件進行求解從而找出預(yù)測精度最優(yōu)的模型.通過多次測試設(shè)置、調(diào)試好各個模型參數(shù)后，對樣本進行模型訓(xùn)練，計算均方誤差（MSE）以及決定系數(shù)（R2）：

通過上述模型預(yù)測的RON損失檢驗結(jié)果如表2所示，BP神經(jīng)網(wǎng)絡(luò)預(yù)測模型的MSE最小，決定系數(shù)R2最接近于1，模型精度要優(yōu)于其他3個模型，故采用BP神經(jīng)網(wǎng)絡(luò)模型預(yù)測辛烷值損失.

3.3 BP神經(jīng)網(wǎng)絡(luò)辛烷值預(yù)測模型驗證分析

通過對比4種預(yù)測模型后，選用BP神經(jīng)網(wǎng)絡(luò)模型作為辛烷值損失預(yù)測模型，并將對此使用10折交叉驗證的方法來判斷結(jié)果的可靠性.即隨機抽取10份樣本，把其中1份樣本作為測試集，其余9份作為訓(xùn)練集，如此循環(huán)10次，使得每份樣本都能作為一次測試集.對隨機建立的10個訓(xùn)練集來計算10個模型，對測試集分別得到10個均方誤差（MSE）、決定系數(shù)（R2），再求出10個模型的MSE均值以及R2均值.

從表3的10折交叉驗證結(jié)果可知，MSE均值為0.027 193（接近0），R2均值為0.904 87（接近1），故認為BP神經(jīng)網(wǎng)絡(luò)的預(yù)測精度較高，符合預(yù)期.但上述10份樣本的模型訓(xùn)練精度低于用全部325個樣本進行訓(xùn)練所得模型，因樣本數(shù)較少，模型訓(xùn)練時不可控因素較多，故10個訓(xùn)練模型僅用于驗證BP神經(jīng)網(wǎng)絡(luò)預(yù)測方法的適用性，而不用于最終的辛烷值損失預(yù)測模型.

本文的辛烷值損失預(yù)測模型將原始樣本中的70%數(shù)據(jù)序列（227個樣本）用于模型訓(xùn)練，15%的數(shù)據(jù)序列（49個樣本）用于模型驗證，剩余15%的數(shù)據(jù)序列（49個樣本）作為模型測試集.基于BP神經(jīng)網(wǎng)絡(luò)的辛烷值損失預(yù)測模型誤差結(jié)果如圖5所示.

根據(jù)BP神經(jīng)網(wǎng)絡(luò)建模過程可知，該模型在第8次迭代時達到最優(yōu)精度，此時驗證集的MSE為0.001 182，訓(xùn)練集和測試集的誤差都較小.圖6是對每個樣本實際值與模型預(yù)測輸出值誤差建立的直方圖，圖中誤差集中在0附近，有76.3%的樣本誤差位于-0.01和0.01之間.除142號樣本外，其余樣本誤差絕對值均小于0.1，所有樣本誤差平均值為-0.001 009 44.模型訓(xùn)練結(jié)果如表4所示，訓(xùn)練集、驗證集、測試集的MSE都足夠小，接近于0，且R2都大于0.99，接近于1，故認為該預(yù)測模型可用于汽油辛烷值損失預(yù)測工程中.

圖7是各數(shù)據(jù)集的輸出值與辛烷值損失實際數(shù)據(jù)的比較，可以發(fā)現(xiàn)模型輸出值能較好地擬合原始數(shù)據(jù)，同樣可驗證利用BP神經(jīng)網(wǎng)絡(luò)模型在15個主要輸出變量基礎(chǔ)上預(yù)測辛烷值損失的精度較高.

4 主要變量操作方案的優(yōu)化

根據(jù)已建立的辛烷值損失預(yù)測模型，并分析了S Zorb裝置產(chǎn)品的歷史脫硫數(shù)據(jù)，計算得知所有油品的硫質(zhì)量分數(shù)平均值是4.75 μg/g.為保證汽油產(chǎn)品脫硫效果，歐VI和國VI標準中汽油產(chǎn)品硫質(zhì)量分數(shù)不得大于10 μg/g，但為了給企業(yè)裝置操作留有空間[14]，要求在實際生產(chǎn)中產(chǎn)品硫質(zhì)量分數(shù)不大于5 μg/g.在此前提下，利用已知數(shù)據(jù)樣本對主要變量進行優(yōu)化，從而使辛烷值損失降幅達到理想狀態(tài).

首先探索單一變量的改變對基于BP神經(jīng)網(wǎng)絡(luò)建立的辛烷值預(yù)測模型的影響.使用控制變量法使其他變量保持不變，只改變單一變量取值，進而預(yù)測辛烷值的損失.在此基礎(chǔ)上，為使辛烷值損失降幅達到30%，利用上文所構(gòu)建的多元線性回歸模型對需要調(diào)控的操作變量進行處理.在優(yōu)化處理過程中原料性質(zhì)、待生吸附劑性質(zhì)、再生吸附劑性質(zhì)保持不變，因此主要調(diào)節(jié)15個主要影響變量中的11個操作變量，通過添加調(diào)節(jié)系數(shù)進行變量調(diào)控.最后提出具體的優(yōu)化操作條件，使辛烷值損失降到最低.

4.1 基于BP神經(jīng)網(wǎng)絡(luò)的主要變量優(yōu)化調(diào)控模型

使用控制變量法探究單一變量對BP神經(jīng)網(wǎng)絡(luò)建立辛烷值預(yù)測模型的影響.首先篩選出硫質(zhì)量分數(shù)不大于5 μg/g的樣本共268個，再隨機抽取10個樣本進行優(yōu)化.預(yù)測模型的主要變量共15個，其中操作變量11個，其取值范圍以及調(diào)整幅度值如表5所示，使用Matlab建立循環(huán)語句，保持其他變量不變，將單一變量由取值范圍的最低值以調(diào)整幅度值為間隔變化至取值范圍的最大值.利用前文建立的BP神經(jīng)網(wǎng)絡(luò)模型預(yù)測辛烷值損失，具體過程如圖8所示.

圖9是對10份樣本分別控制9個操作變量的取值：再生器頂?shù)撞顗骸-101輻射室出口壓力、穩(wěn)定塔底出口溫度、D-201含硫污水排量、加氫裂化輕石腦油進裝置流量、干氣出裝置流量、S_ZORB AT-0010、D-109吸附劑料位、R-102床層吸附劑料位密度等.基于前文模型預(yù)測出的辛烷值損失變化曲線，由于D-109壓力、非凈化風進裝置壓力2個變量取值范圍較小，調(diào)整幅度值相對較大，導(dǎo)致循環(huán)次數(shù)過少，因此未在圖中表示.分析圖9發(fā)現(xiàn)，單一控制穩(wěn)定塔底出口溫度升高，辛烷值（RON）損失預(yù)測值明顯下降，且溫度低于120 ℃左右時，下降幅度較大，此后溫度繼續(xù)升高對RON損失的影響較小.除了該變量外，單一控制其他變量并不能對RON損失預(yù)測值產(chǎn)生明顯影響，或影響不穩(wěn)定，如改變D-201含硫污水排量時，4個樣本的損失值無明顯變化，3個樣本損失值明顯上升，一個樣本損失值小幅度波動，還有2個樣本損失值呈下降趨勢.

基于以上分析，研究認為僅僅改變單一變量取值并不能影響RON損失值，更不能達到辛烷值損失降幅大于30%的生產(chǎn)要求，所以需要考慮同時改變多個變量，優(yōu)化操作條件，從而達到盡量降低辛烷值損失的目的.由于神經(jīng)網(wǎng)絡(luò)算法屬于“黑箱方法”，在研究時只能得出輸入輸出變量，無法了解內(nèi)部結(jié)構(gòu)，若利用該算法同時對多個變量進行優(yōu)化較為復(fù)雜.而多元線性回歸模型的預(yù)測精度僅次于BP神經(jīng)網(wǎng)絡(luò)，且算法簡單、易操作，更適合同時對多個變量進行優(yōu)化，故將利用該模型預(yù)測辛烷值損失并優(yōu)化主要操作變量.

4.2 基于多元線性回歸的主要變量優(yōu)化調(diào)控模型

多元線性回歸模型體現(xiàn)的是多個解釋變量和被解釋變量之間的關(guān)系問題，利用篩選得出的15個最主要影響變量和1個因變量，建立以下關(guān)系式：

式中：α0為常數(shù)項；α2，α9，…，α11為回歸系數(shù).

此外選擇擬合優(yōu)度檢驗法對該模型進行計算驗證.擬合優(yōu)度的含義是樣本具體觀察數(shù)值在回歸線附近聚集的緊密程度.一般選用R2判斷多元線性回歸擬合優(yōu)度.而擬合優(yōu)度是在分解總離差平方的基礎(chǔ)上測算出來的.其中總離差平方和計算公式為

其中，SSE為殘差平方和，SSR為回歸平方和，SST為總離差平方.計算公式如下：

其中，是樣本觀察值均值，是估計值.決定系數(shù)R2便是通過回歸平方和占總離差平方和的比例，計算公式為

決定系數(shù)R2反映線性回歸方程的擬合程度，表示解釋變量和被解釋變量之間存在的回歸關(guān)系可以用來解釋所有偏差中的百分比.R2取值范圍為[0，1]，趨近于1，擬合效果越好，越趨近于0效果越差[15].

為使辛烷值損失降幅達到30%，對需要調(diào)控的操作變量進行處理.在優(yōu)化處理過程中原料性質(zhì)、待生吸附劑性質(zhì)、再生吸附劑性質(zhì)保持不變，即在11個操作變量前添加調(diào)節(jié)系數(shù)來進行調(diào)控，基于多元線性回歸模型的變量調(diào)控模型表達式為

在最小二乘法準則的指導(dǎo)下，對各個調(diào)控系數(shù)進行參數(shù)估計.將產(chǎn)品硫質(zhì)量分數(shù)不大于5 μg/g的樣本數(shù)據(jù)代入式（20），得到方程組：

通過多次測試設(shè)置并調(diào)試好各個模型參數(shù)后，對樣本進行模型訓(xùn)練，然后利用所得模型預(yù)測辛烷值損失，運用上文所述相關(guān)誤差檢驗公式，計算得出MSE為0.000 26，R2為0.995 1，模型精度較好，故可采用多元線性回歸模型預(yù)測辛烷值損失和優(yōu)化操作變量.

根據(jù)最小二乘法估算出來的調(diào)節(jié)系數(shù)數(shù)值大小和正負性來對變量進行調(diào)控.對再生器頂?shù)撞顗汉蚏-102床層吸附劑料位密度按照一定的比例進行反向調(diào)控，對D-109壓力、F-101輻射室出口壓力、穩(wěn)定塔底出口溫度、D-201含硫污水排量、非凈化風進裝置壓力、S_ZORB AT-0010和D-109吸附劑料位按照一定的比例進行正向調(diào)控，那么辛烷值損失降幅就會達到甚至超過30%.

5 優(yōu)化調(diào)整過程的可視化展示

工業(yè)生產(chǎn)在對汽油進行精制過程中，需要對催化裂化汽油進行脫硫和降烯烴，進而提高汽油的燃燒性能.然而由于生產(chǎn)條件、機器設(shè)備各方面的限制，無法將各操作變量一次性提高至最優(yōu)值，故需要逐步平穩(wěn)調(diào)整.由于所收集到的133號樣本與其他企業(yè)的大部分S Zorb裝置生產(chǎn)數(shù)據(jù)相似，所以主要針對133號樣本數(shù)據(jù)進行分析[16].為充分了解并分析操作變量優(yōu)化過程中汽油辛烷值與硫含量的變化趨勢，首先需要建立硫含量預(yù)測模型，由上文中預(yù)測模型誤差對比發(fā)現(xiàn)，BP神經(jīng)網(wǎng)絡(luò)算法精度最高，故可利用該方法構(gòu)建硫含量預(yù)測模型.操作條件的優(yōu)化通過操作變量在取值范圍內(nèi)按照一定幅度的逐漸變化實現(xiàn)，并在主要操作變量基礎(chǔ)上預(yù)測辛烷值損失及硫含量，繪制出調(diào)整過程中產(chǎn)品性質(zhì)的變化軌跡[17].

5.1 構(gòu)建硫含量預(yù)測模型

與基于BP神經(jīng)網(wǎng)絡(luò)辛烷值損失預(yù)測模型相類似，采用Matlab軟件工具箱nftool處理變量數(shù)據(jù)，輸入層為15個主要影響變量，輸出層為硫含量，隱含層節(jié)點數(shù)為10，變量數(shù)據(jù)為原有數(shù)據(jù)中硫質(zhì)量分數(shù)小于等于5 μg/g的樣本，一共268個.其中70%樣本數(shù)據(jù)作為訓(xùn)練集，15%樣本數(shù)據(jù)作為驗證集，剩余15%樣本數(shù)據(jù)為測試集，從而建立基于BP神經(jīng)網(wǎng)絡(luò)的硫含量預(yù)測模型（圖10）.模型的MSE為0.059 3，與辛烷值損失預(yù)測模型相比誤差較大，但仍在允許誤差范圍內(nèi)，因此該模型可用于硫含量預(yù)測[18].

5.2 主要變量優(yōu)化的可視化分析

關(guān)于操作變量的優(yōu)化條件，已在本文第4章中詳細闡述，可以得出再生器頂?shù)撞顗海ㄓ米兞?代替）、F-101輻射室出口壓力（用變量2代替）、非凈化風進裝置壓力的大幅度變化，D-109壓力（用變量3代替）、穩(wěn)定塔底出口溫度（用變量4代替）、S_ZORB AT-0010（用變量5代替）的小幅度優(yōu)化會對降低產(chǎn)品的辛烷值損失產(chǎn)生影響[19].

圖11為改變變量1和變量2所預(yù)測得出的產(chǎn)品硫質(zhì)量分數(shù)變化軌跡.由圖11可以發(fā)現(xiàn)，樣本的硫質(zhì)量分數(shù)預(yù)測結(jié)果波動較大.在變量1位于[30，35]，變量2位于[-0.3，-0.2]范圍內(nèi)時，硫質(zhì)量分數(shù)取極小值；變量1接近于上限，變量2接近于下限時，硫質(zhì)量分數(shù)取極大值；同時變量1極小，變量2極大時，硫質(zhì)量分數(shù)也是小范圍的極大值.

圖12為改變變量1和變量2所預(yù)測得出的產(chǎn)品辛烷值損失變化軌跡.三維曲面圖較為平滑，符合上文的分析結(jié)果，同時可以發(fā)現(xiàn)變量1為極小值，變量2為極大值時，辛烷值的損失預(yù)測結(jié)果最小.

圖13表示其他變量不變時，133號樣本的D-109壓力（變量3）分別為取值范圍的上、下限時，改變變量4和變量5取值預(yù)測出的硫含量變化軌跡.當變量3取下限時，產(chǎn)品硫含量較低，但兩者變化趨勢大體相同，都在變量4取值最大時達到極小值.

與圖13相類似，圖14是不同D-109取值下辛烷值損失預(yù)測值的變化軌跡，可知變量3的不同取值不會影響辛烷值損失的變化范圍，但改變了辛烷值損失的變化軌跡.變量3取下限時，辛烷值損失在變量4較低時達到極小；而變量3取上限時，恰恰相反辛烷值損失在變量4較低時達到極大[20].

6 結(jié)論

通過對S Zorb裝置運行積累的數(shù)據(jù)進行建模分析，得出以下幾個結(jié)論：

1） 對于催化裂化汽油精制脫硫裝置運行積累的大量數(shù)據(jù)，挑選建模變量時利用了二步篩選法.首先基于Lasso線性模型對多個影響變量進行初步篩選，基于此使用BP神經(jīng)網(wǎng)絡(luò)的指標因子貢獻度計算模型確定參與化工建模的主要變量.綜合利用了2個篩選模型的優(yōu)勢，第1步去除變量的多重共線性，第2步進一步保證留下的影響變量對辛烷值損失有足夠的貢獻度，兩者結(jié)合使得篩選得出的變量不僅更精簡而且與因變量關(guān)系更加緊密.

2） 在構(gòu)建辛烷值預(yù)測模型時把神經(jīng)網(wǎng)絡(luò)和多元線性回歸模型相結(jié)合，首先將基于最小二乘優(yōu)化算法（Levenberg-Marquardt）的BP神經(jīng)網(wǎng)絡(luò)模型用于辛烷值損失預(yù)測，其精度和可操作性較高，能夠更加準確地預(yù)測辛烷值損失.而在優(yōu)化調(diào)控階段，將神經(jīng)網(wǎng)絡(luò)模型和多元線性回歸模型同時用于影響變量的調(diào)控，兩者優(yōu)勢互補并探索出一個有效的變量優(yōu)化模式.

3） 為充分了解并分析操作變量優(yōu)化過程中汽油辛烷值與硫含量的變化趨勢，分析誤差對比得知神經(jīng)網(wǎng)絡(luò)算法精度最高，即建立基于神經(jīng)網(wǎng)絡(luò)的硫含量預(yù)測模型.操作條件的優(yōu)化通過操作變量在取值范圍內(nèi)按照一定幅度的逐漸變化實現(xiàn)，并在主要操作變量基礎(chǔ)上預(yù)測辛烷值損失及硫含量，可視化展示了優(yōu)化調(diào)整過程中產(chǎn)品性質(zhì)的變化軌跡.

4） 研究所構(gòu)建的降低汽油精制過程中的辛烷值損失模型能明顯降低石油企業(yè)煉油生產(chǎn)過程辛烷值的損失，并保證產(chǎn)品的脫硫效果，有效地解決了普通化工建模操作變量少、嚴重非線性及相互強耦聯(lián)的問題，可以為石油企業(yè)的S Zorb裝置生產(chǎn)實際操作提供合理的指導(dǎo).而后續(xù)研究中可根據(jù)實際生產(chǎn)的實時性，在預(yù)測汽油辛烷值時利用自適應(yīng)在線極限學習機模型，從而更有效地保證油品質(zhì)量，提高生產(chǎn)效率，增加經(jīng)濟效益[21].

參考文獻

References

［1］ 丁怡曼，薛曉康，范賓，等.基于PLS-紅外光譜的汽油辛烷值測定方法研究[J].化學研究與應(yīng)用，2021，33（5）：863-867

DING Yiman，XUE Xiaokang，F(xiàn)AN Bin，et al.Determination of gasoline octane number based on PLS-infrared spectroscopy[J].Chemical Research and Application，2021，33（5）：863-867

[2] Kardamakis A A，Pasadakis N.Autoregressive modeling of near-IR spectra and MLR to predict RON values of gasolines[J].Fuel，2010，89（1）：158-161

[3] 黃水望，趙曉鋒，郭振，等.氣相色譜法計算汽油的研究法辛烷值[J].廣州化工，2018，46（1）：145-146，186

HUANG Shuiwang，ZHAO Xiaofeng，GUO Zhen，et al.Determination of research octane number of gasoline by gas chromatography[J].Guangzhou Chemical Industry，2018，46（1）：145-146，186

[4] Ghosh P，Hickey K J，Jaffe S B.Development of a detailed gasoline composition-based octane model[J].Industrial amp; Engineering Chemistry Research，2006，45（1）：337-345

[5] 熊春華，田高友.汽油抗爆性能指標與其他理化指標關(guān)聯(lián)研究[C]//中國汽車工程學會燃料與潤滑油分會第十四屆年會論文集.沈陽：遼海出版社，2010：126-132

XIONG Chunhua，TIAN Gaoyou.The correlations study of antiknocking properties with other physical-chemical properties of gasoline[C]//Proceedings of the 14th Annual Meeting of Fuel and Lubricant Branch of Chinese Society of Automotive Engineering.Shenyang：Liaohai Publishing House，2010：126-132

[6] 楊帆，周敏，金繼民，等.智能優(yōu)化算法及人工神經(jīng)網(wǎng)絡(luò)在催化裂化模型分析中的應(yīng)用進展[J].石油學報（石油加工），2020，36（4）：878-888

YANG Fan，ZHOU Min，JIN Jimin，et al.Research progress on application of intelligent optimization algorithm and artificial neural network in FCC model analysis[J].Acta Petrolei Sinica （Petroleum Processing Section），2020，36（4）：878-888

[7] 歐陽福生，游俊峰，方偉剛.BP神經(jīng)網(wǎng)絡(luò)結(jié)合遺傳算法優(yōu)化MIP工藝的產(chǎn)品分布[J].石油煉制與化工，2018，49（8）：98-104

OUYANG Fusheng，YOU Junfeng，F(xiàn)ANG Weigang.Optimizing product distribution of MIP process by BP neural network combined with genetic algorithm[J].Petroleum Processing and Petrochemicals，2018，49（8）：98-104

[8] 楊帆，周敏，戴超男，等.基于人工智能算法的催化裂化裝置汽油收率預(yù)測模型的構(gòu)建與分析[J].石油學報（石油加工），2019，35（4）：807-817

YANG Fan，ZHOU Min，DAI Chaonan，et al.Construction and analysis of gasoline yield prediction model for FCC unit based on artificial intelligence algorithm[J].Acta Petrolei Sinica （Petroleum Processing Section），2019，35（4）：807-817

[9] 張玉瑞，陳微微，周曉龍，等.一種改進的調(diào)合辛烷值模型預(yù)測汽油研究法辛烷值[J].石油煉制與化工，2016，47（6）：42-46

ZHANG Yurui，CHEN Weiwei，ZHOU Xiaolong，et al.RON prediction by improved model for blended gasoline[J].Petroleum Processing and Petrochemicals，2016，47（6）：42-46

[10] 王天宇，劉忠保，黃明富，等.采用人工神經(jīng)網(wǎng)絡(luò)方法建立加氫裂化反應(yīng)體系模型[J].石油煉制與化工，2015，46（8）：90-95

WANG Tianyu，LIU Zhongbao，HUANG Mingfu，et al.Modeling VGO hydrocracking process by BP-ANN technology[J].Petroleum Processing and Petrochemicals，2015，46（8）：90-95

[11] 朱燚丹，陳興榮，李秋萍.基于信息增益率的超高維變量選擇[J].統(tǒng)計與決策，2021，37（22）：18-21

ZHU Yidan，CHEN Xingrong，LI Qiuping.Ultra-high dimensional variable selection based on information gain rate[J].Statistics and Decision，2021，37（22）：18-21

[12] 楊軼男，任曄，毛安國，等.影響催化裂化裝置汽油辛烷值變化的技術(shù)因素分析[J].煉油技術(shù)與工程，2019，49（6）：32-35

YANG Yinan，REN Ye，MAO Anguo，et al.Analysis of technical factors affecting gasoline octane number in catalytic cracking unit[J].Petroleum Refinery Engineering，2019，49（6）：32-35

[13] 劉禹含，曹萃文.基于LightGBM的催化重整裝置產(chǎn)品預(yù)測及操作優(yōu)化相關(guān)性分析[J].石油學報（石油加工），2020，36（4）：756-766

LIU Yuhan，CAO Cuiwen.Product prediction technology and optimal operation correlation analysis for catalytic reforming unit based on LightGBM[J].Acta Petrolei Sinica （Petroleum Processing Section），2020，36（4）：756-766

[14] 胡碧霞，張紅光，盧建剛，等.汽油辛烷值近紅外光譜檢測的改進極限學習機建模方法[J].南京理工大學學報，2017，41（5）：660-665

HU Bixia，ZHANG Hongguang，LU Jiangang，et al.Novel modeling method based on improved extreme learning machine algorithm for gasoline octane number detection by near infrared spectroscopy[J].Journal of Nanjing University of Science and Technology，2017，41（5）：660-665

[15] 鄭斌，孫洪霞，王維民.基于隨機森林回歸的汽油研究法辛烷值預(yù)測[J].石油煉制與化工，2020，51（12）：69-75

ZHENG Bin，SUN Hongxia，WANG Weimin.Prediction of gasoline research octane number based on random forest regression[J].Petroleum Processing and Petrochemicals，2020，51（12）：69-75

[16] Kubic W L，Jenkins R W，Moore C M，et al.Artificial neural network based group contribution method for estimating cetane and octane numbers of hydrocarbons and oxygenated organic compounds[J].Industrial amp; Engineering Chemistry Research，2017，56（42）：12236-12245

[17] Abdul J A G，van Oudenhoven V，Emwas A H，et al.Predicting octane number using nuclear magnetic resonance spectroscopy and artificial neural networks[J].Energy amp; Fuels，2018，32（5）：6309-6329

[18] Razavi-Far R，Hallaji E，F(xiàn)arajzadeh-Zanjani M，et al.Information fusion and semi-supervised deep learning scheme for diagnosing gear faults in induction machine systems[J].IEEE Transactions on Industrial Electronics，2019，66（8）：6331-6342

[19] Yang X B，Zheng X L，Gao H J.SGD-based adaptive NN control design for uncertain nonlinear systems[J].IEEE Transactions on Neural Networks and Learning Systems，2018，29（10）：5071-5083

[20] Yuan X F，Huang B，Wang Y L，et al.Deep learning-based feature representation and its application for soft sensor modeling with variable-wise weighted SAE[J].IEEE Transactions on Industrial Informatics，2018，14（7）：3235-3243

[21] 韓萬龍，范崢，薛崗，等.利用BP人工神經(jīng)網(wǎng)絡(luò)預(yù)測天然氣中重組分對凈化裝置的影響[J].石油與天然氣化工，2018，47（6）：1-6

HAN Wanlong，F(xiàn)AN Zheng，XUE Gang，et al.Effects prediction of heavy components in natural gas on purification unit by BP artificial neural network[J].Chemical Engineering of Oil amp; Gas，2018，47（6）：1-6

Octane number prediction based on BP neural network and multiple linear regression

XU Meixian1 ZHENG Yan1 ZHOU Ruolan1 ZHANG Ruyi1

1College of Automobile and Traffic Engineering，Nanjing Forestry University，Nanjing 210037

Abstract In order to reduce the sulfur and olefin and the loss of octane number so as to promote the clean production of gasoline，an octane number loss prediction model is established based on data accumulated by the S Zorb device.First，the Lasso is used to screen out the modeling variables，then the index factor contributions are calculated by the BP neural network，based on which 15 main variables are screened out to build the model.Second，four modeling approaches are compared and analyzed，which shows that the BP neural network has better prediction accuracy thus is more suitable to model the octane number loss.The ten-fold cross-validation produces the average MSE value of 0.027 193 and the average R2 value of 0.904 87，verifying the reliability of the model.Furthermore，the main variables are optimized and adjusted by multiple linear regression under the premise that the sulfur content is not greater than 5 μg/g.The results show that multiple variables need to be adjusted simultaneously to reduce the octane number loss by more than 30%.The multiple linear regression model has good prediction accuracy and can adjust main variables positively or negatively according to a certain proportion.The trajectories of octane number and sulfur content are also visualized in the paper.

Key words BP neural network;multiple linear regression;Lasso algorithm;octane number loss prediction;optimized regulation

收稿日期2022-04-26

資助項目國家自然科學基金（71701099，71501090）；江蘇省高等學校自然科學研究項目（17KJB580008）

作者簡介許美賢，女，碩士生，主要從事復(fù)雜工業(yè)過程建模、數(shù)據(jù)挖掘、人工智能算法的研究.xumeixian3210@163.com

鄭琰（通信作者），女，博士，副教授，主要從事復(fù)雜工業(yè)過程建模、機器學習的研究.yzheng_x@163.com