陳志禮 ，寧正福 ，杜華明 ，黃亮 ，葉洪濤 ，張文通

（1.中國(guó)石油大學(xué)（北京）油氣資源與探測(cè)國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京102249；2.中國(guó)石油大學(xué)（北京）石油工程學(xué)院，北京102249）

0 引言

頁(yè)巖吸附氣資源豐富，占頁(yè)巖氣總量的20%～85%，頁(yè)巖吸附氣量的預(yù)測(cè)是儲(chǔ)量評(píng)價(jià)的基礎(chǔ)［1］。目前，學(xué)者們主要依靠物理吸附實(shí)驗(yàn)以及理論吸附模型來(lái)研究頁(yè)巖吸附能力和吸附機(jī)理，常用的等溫吸附模型包括 Langmuir，F(xiàn)reundlich，TothL，BET，B-BET 和 TBET，Langmuir-Freundlich 等［2-5］。然而，頁(yè)巖吸附量受TOC（總有機(jī)碳質(zhì)量分?jǐn)?shù)）、Ro（鏡質(zhì)體反射率）、T（溫度）、p（壓力）等多重因素的共同作用［3］，常規(guī)的等溫吸附模型僅考慮恒定溫度下，壓力變化的情況，而無(wú)法對(duì)Ro，TOC，T，p的同時(shí)變化情況進(jìn)行分析。因此，建立多因素作用下的頁(yè)巖吸附量預(yù)測(cè)模型尤為必要。

人工神經(jīng)網(wǎng)絡(luò)具有很強(qiáng)的非線性映射能力，近年來(lái)在油氣勘探與開發(fā)領(lǐng)域得到了廣泛應(yīng)用［6-7］。目前，實(shí)際應(yīng)用最為成功的神經(jīng)網(wǎng)絡(luò)當(dāng)屬BP神經(jīng)網(wǎng)絡(luò)。BP神經(jīng)網(wǎng)絡(luò)主要用于模式識(shí)別、函數(shù)逼近及預(yù)測(cè)等［8］，但BP神經(jīng)網(wǎng)絡(luò)存在收斂速度慢、搜索能力弱、易陷入局部最小的不足，且預(yù)測(cè)結(jié)果隨訓(xùn)練而波動(dòng)［9］。遺傳算法（GA）［10］和核密度估計(jì)法［11］是改進(jìn) BP 網(wǎng)絡(luò)的可行方法。GA優(yōu)化BP網(wǎng)絡(luò)權(quán)值和閾值，加速收斂，跳出局部最優(yōu)。核密度估計(jì)法通過(guò)循環(huán)預(yù)測(cè)和核密度估計(jì)尋求BP網(wǎng)絡(luò)最優(yōu)預(yù)測(cè)值，能在一定程度上穩(wěn)定網(wǎng)絡(luò)預(yù)測(cè)。

筆者在總結(jié)分析頁(yè)巖吸附能力影響因素的基礎(chǔ)上，構(gòu)建BP神經(jīng)網(wǎng)絡(luò)，利用試錯(cuò)法確定隱含層最優(yōu)節(jié)點(diǎn)數(shù)，運(yùn)用GA和核密度估計(jì)法優(yōu)化權(quán)值、閾值，提高預(yù)測(cè)穩(wěn)定性，最終建立考慮多因素共同作用的頁(yè)巖吸附量預(yù)測(cè)模型，并對(duì)影響因素的敏感性進(jìn)行了分析。

1 頁(yè)巖吸附能力的主要影響因素

總有機(jī)碳質(zhì)量分?jǐn)?shù)是控制頁(yè)巖氣聚集的關(guān)鍵因素。國(guó)內(nèi)外學(xué)者［12-14］研究表明：頁(yè)巖的TOC越高，吸附量也相應(yīng)越大，兩者存在良好的正相關(guān)關(guān)系。Ro與頁(yè)巖的生氣特征與吸附能力密切相關(guān)。在生油窗階段與濕氣、凝析油生成階段，頁(yè)巖吸附能力與Ro呈正相關(guān)；進(jìn)入干氣階段后，Ro對(duì)頁(yè)巖吸附能力具有抑制作用，且Ro越大，抑制作用越顯著［15］。溫度和壓力是影響頁(yè)巖吸附能力的2個(gè)重要外部因素。頁(yè)巖吸附甲烷屬于放熱過(guò)程，隨溫度升高，頁(yè)巖吸附能力顯著降低。Chalmers等［16］通過(guò)吸附實(shí)驗(yàn)發(fā)現(xiàn)，溫度與氣體吸附能力呈負(fù)相關(guān)關(guān)系，且其影響遠(yuǎn)大于TOC。壓力對(duì)頁(yè)巖吸附的影響可用Langmuir模型［17］解釋，當(dāng)實(shí)際壓力低于飽和吸附壓力時(shí)，壓力增加，促進(jìn)甲烷吸附；反之，甲烷吸附達(dá)到平衡，壓力對(duì)甲烷吸附量貢獻(xiàn)不大。

2 BP神經(jīng)網(wǎng)絡(luò)

2.1 BP神經(jīng)網(wǎng)絡(luò)原理

BP神經(jīng)網(wǎng)絡(luò)是一種信號(hào)向前傳遞、誤差反向傳播的多層前饋網(wǎng)絡(luò)，包括輸入層、隱含層和輸出層，具體的算法流程包括4部分：數(shù)據(jù)的選擇與歸一化處理、BP神經(jīng)網(wǎng)絡(luò)構(gòu)建、網(wǎng)絡(luò)訓(xùn)練、預(yù)測(cè)與反歸一化處理。

為消除數(shù)量級(jí)引起的誤差，采用最大最小法需對(duì)原始數(shù)據(jù)進(jìn)行歸一化處理：

式中：Xi為原始訓(xùn)練、測(cè)試自變量數(shù)據(jù)；Xmax，Xmin分別為原始訓(xùn)練、測(cè)試自變量數(shù)據(jù)中最大值、最小值；Xin為歸一化后自變量數(shù)據(jù)。

BP神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)輸出時(shí)，需要對(duì)結(jié)果數(shù)據(jù)進(jìn)行反歸一化處理：

式中：Ymax，Ymin分別為原始因變量數(shù)據(jù)中最大值、最小值；Yon為BP神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)的歸一化數(shù)據(jù)；Y0為反歸一化處理后的最終預(yù)測(cè)結(jié)果值。

2.2 BP神經(jīng)網(wǎng)絡(luò)的改進(jìn)

雖然BP神經(jīng)網(wǎng)絡(luò)具有誤差反向傳遞的特征，但仍存在3個(gè)顯著不足：1）學(xué)習(xí)率太小，導(dǎo)致收斂速度過(guò)慢。2）存在局部最小值，可能導(dǎo)致網(wǎng)絡(luò)陷入局部最小點(diǎn)。3）結(jié)果不唯一，每次運(yùn)行對(duì)應(yīng)的權(quán)值和閾值發(fā)生波動(dòng)，結(jié)果可能存在差異。針對(duì)上述問(wèn)題，本文采用遺傳算法［10］、核密度估計(jì)法［11］對(duì) BP 神經(jīng)網(wǎng)絡(luò)進(jìn)行改進(jìn)，改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)算法流程見圖1。

圖1 改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)算法流程

本文應(yīng)用GA優(yōu)化BP神經(jīng)網(wǎng)絡(luò)的初始化權(quán)值和閾值，以彌補(bǔ)其收斂速度慢、搜索能力弱且易陷入局部最小的不足。GA優(yōu)化的一般步驟：

1）種群初始化（編碼）。二進(jìn)制編碼和實(shí)數(shù)編碼是常用的個(gè)體編碼方法，當(dāng)編碼串較短時(shí)，二進(jìn)制編碼存在映射誤差，編碼串較長(zhǎng)時(shí)對(duì)應(yīng)的搜索空間膨脹式增加［18］。因此，對(duì)輸入層和隱含層的連接權(quán)值、隱含層和輸出層的連接權(quán)值、隱含層和輸出層的閾值進(jìn)行實(shí)數(shù)編碼。

2）適應(yīng)度函數(shù)?；趥€(gè)體解碼，對(duì)BP網(wǎng)絡(luò)進(jìn)行訓(xùn)練和預(yù)測(cè)，預(yù)測(cè)結(jié)果誤差絕對(duì)值的總和作為個(gè)體適應(yīng)度F，計(jì)算式為

式中：m為輸出層節(jié)點(diǎn)數(shù)；Yi，yi分別為第i個(gè)節(jié)點(diǎn)的預(yù)測(cè)值、期望值。

3）選擇操作。采用基于適應(yīng)度比例的輪盤賭法對(duì)樣本個(gè)體進(jìn)行選擇，個(gè)體i對(duì)應(yīng)的選擇概率Pi的計(jì)算公式為

式中：Fi為第i個(gè)個(gè)體適應(yīng)度值；w為種群規(guī)模。

4）交叉操作。采用實(shí)數(shù)交叉法進(jìn)行交叉操作，染色體Ak與Al在j處的交叉，產(chǎn)生的新個(gè)體Akj和Alj計(jì)算公式為

5）變異操作。利用隨機(jī)數(shù)控制基因變異，第i個(gè)個(gè)體的第j個(gè)基因Aij的變異操作，計(jì)算公式為

式中：Amax，Amin分別為 Aij基因上、下界限；r為 0～1 的隨機(jī)數(shù)；d，dmax分別為當(dāng)前、最大迭代次數(shù)。

在BP神經(jīng)網(wǎng)絡(luò)訓(xùn)練和預(yù)測(cè)環(huán)節(jié)設(shè)定循環(huán)，利用非參數(shù)核密度估計(jì)法得到核密度最大值對(duì)應(yīng)的橫軸坐標(biāo)，即為最優(yōu)預(yù)測(cè)值。

核密度估計(jì)函數(shù)為

式中：xi為一元連續(xù)總體的樣本；n為樣本總數(shù)；K為核函數(shù)（本文選用Gaussian函數(shù)）；h為窗寬。

h是影響核密度估計(jì)函數(shù)的重要因素，可利用綜合均方誤差函數(shù)MISE求取最佳窗寬：

式中：M為綜合均方誤差函數(shù)MISE。

求解公式（11）可得最佳窗寬：

式中：f（x）為總體的真實(shí)分布密度；σ可視為樣本的標(biāo)準(zhǔn)差。

3 頁(yè)巖吸附量預(yù)測(cè)模型

3.1 樣本數(shù)據(jù)的選擇與處理

基于前人實(shí)驗(yàn)數(shù)據(jù)［3］，累計(jì)4個(gè)巖樣，80組不同T，p，TOC和Ro下的吸附量實(shí)驗(yàn)數(shù)據(jù)，其中75組數(shù)據(jù)用于訓(xùn)練學(xué)習(xí)，剩余5組數(shù)據(jù)用于檢驗(yàn)，利用最大最小法進(jìn)行處理。

3.2 BP神經(jīng)網(wǎng)絡(luò)構(gòu)建與初始化

頁(yè)巖吸附能力主要與TOC，Ro，T和p有關(guān)，因此輸入、輸出層的神經(jīng)元個(gè)數(shù)分別設(shè)置為4和1?？紤]到訓(xùn)練速度不宜過(guò)慢，本文采用單隱含層，通過(guò)調(diào)整隱含層節(jié)點(diǎn)數(shù)來(lái)提高預(yù)測(cè)精度。筆者參考經(jīng)驗(yàn)公式確定隱含層最優(yōu)節(jié)點(diǎn)數(shù)為

式中：m，n分別為輸入、輸出層節(jié)點(diǎn)數(shù)；a為0～10的常數(shù)。

通過(guò)試錯(cuò)法［19］求取隱含層最優(yōu)節(jié)點(diǎn)數(shù)大致范圍為2～12，即利用訓(xùn)練數(shù)據(jù)依次對(duì)2～12節(jié)點(diǎn)情況進(jìn)行預(yù)測(cè)檢驗(yàn)，對(duì)比分析誤差以得到隱含層最優(yōu)節(jié)點(diǎn)數(shù)，各節(jié)點(diǎn)對(duì)應(yīng)的誤差數(shù)據(jù)如表1所示。由表1可知，隨著節(jié)點(diǎn)數(shù)的增加，平均相對(duì)誤差、均方誤差均呈現(xiàn)下降趨勢(shì)。當(dāng)節(jié)點(diǎn)數(shù)等于12時(shí)，平均相對(duì)誤差為0.019 2，均方誤差為2.480%，預(yù)測(cè)精度最高。為保證預(yù)測(cè)精度，本文隱含層節(jié)點(diǎn)數(shù)選用12，構(gòu)建的BP神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)為4—12—1（見圖 2）。

表1 誤差數(shù)據(jù)

圖2 構(gòu)建的BP神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)

遺傳算法、BP神經(jīng)網(wǎng)絡(luò)基本參數(shù)的設(shè)定是神經(jīng)元訓(xùn)練和學(xué)習(xí)的基礎(chǔ)，直接影響預(yù)測(cè)效果，筆者結(jié)合實(shí)際使用情況進(jìn)行參數(shù)設(shè)置（見表2）。

表2 遺傳算法與BP網(wǎng)絡(luò)參數(shù)設(shè)置

3.3 改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)訓(xùn)練

根據(jù)改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)基本原理和上述初始化參數(shù)，采用MATLAB編程?；谶z傳算法得到最優(yōu)個(gè)體適應(yīng)度變化曲線（見圖3）和最優(yōu)權(quán)值、閾值，利用BP神經(jīng)網(wǎng)絡(luò)進(jìn)行訓(xùn)練，相對(duì)誤差基本處于±2%，平均相對(duì)誤差為2.480%，基本滿足訓(xùn)練精度要求。

圖3 遺傳算法優(yōu)化過(guò)程

4 預(yù)測(cè)與結(jié)果分析

利用訓(xùn)練完畢的改進(jìn)BP神經(jīng)網(wǎng)絡(luò)，對(duì)5組檢驗(yàn)數(shù)據(jù)進(jìn)行100次循環(huán)預(yù)測(cè)，采用核密度估計(jì)法確定最優(yōu)預(yù)測(cè)結(jié)果。以檢驗(yàn)1數(shù)據(jù)為例。通過(guò)核密度分布圖確定最大核密度估計(jì)值坐標(biāo)為（0.796，3.176），預(yù)測(cè)值最大分布區(qū)間為0.750～0.800，占整體的73%，0.796即為檢驗(yàn)1組數(shù)據(jù)對(duì)應(yīng)的最優(yōu)預(yù)測(cè)值（見圖4）。

圖4 檢驗(yàn)1組核密度分布

同理，可得到其他4組的最優(yōu)預(yù)測(cè)結(jié)果，見表3。利用表3數(shù)據(jù)作出預(yù)測(cè)值與實(shí)測(cè)值間的關(guān)系曲線 （見圖5）。由表3可知，5組檢驗(yàn)數(shù)據(jù)的相對(duì)誤差分布區(qū)間為-10%～10%，平均相對(duì)誤差為4.999%。由圖5可知，預(yù)測(cè)值與實(shí)際值的相關(guān)系數(shù)為0.978 9，兩者關(guān)系曲線與標(biāo)準(zhǔn)線基本重合，表明該方法構(gòu)建的BP神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)精度較高，可用于頁(yè)巖吸附能力的預(yù)測(cè)。

表3 檢驗(yàn)數(shù)據(jù)對(duì)應(yīng)的預(yù)測(cè)結(jié)果

圖5 預(yù)測(cè)值與實(shí)測(cè)值的關(guān)系曲線

5 多因素敏感性分析

為探究不同因素對(duì)頁(yè)巖吸附量的影響，筆者采用正交試驗(yàn)方法［20］開展了多因素敏感性分析。首先，設(shè)定p，T，TOC和Ro的因素水平，而后利用正交設(shè)計(jì)軟件建立試驗(yàn)計(jì)劃表，按本文的預(yù)測(cè)模型確定對(duì)應(yīng)的吸附氣量，進(jìn)而得到直觀結(jié)果（見表4）。

表4 直觀結(jié)果

通過(guò)正交試驗(yàn)，同時(shí)求取預(yù)測(cè)模型的計(jì)算方差，進(jìn)行顯著性排序。由表5可知，4個(gè)因素的顯著性排序?yàn)閜＞T＞TOC＞Ro，即頁(yè)巖吸附量主要受壓力影響，其次受溫度和總有機(jī)質(zhì)碳質(zhì)量分?jǐn)?shù)的影響，而有機(jī)質(zhì)成熟度對(duì)頁(yè)巖吸附量的影響最弱。

表5 方差分析

6 結(jié)論

1）利用遺傳算法和核密度估計(jì)改進(jìn)了BP神經(jīng)網(wǎng)絡(luò)，建立了頁(yè)巖吸附量預(yù)測(cè)模型。基于80組實(shí)驗(yàn)數(shù)據(jù)對(duì)BP網(wǎng)絡(luò)進(jìn)行訓(xùn)練和檢驗(yàn)，75組訓(xùn)練數(shù)據(jù)預(yù)測(cè)結(jié)果的平均相對(duì)誤差為2.480%，5個(gè)檢驗(yàn)組的平均相對(duì)誤差為4.999%。預(yù)測(cè)值與實(shí)際值的相關(guān)系數(shù)為0.978 9，表明改進(jìn)的BP神經(jīng)網(wǎng)絡(luò)模型可用于預(yù)測(cè)多因素共同作用的頁(yè)巖吸附量。

2）基于構(gòu)建的多因素預(yù)測(cè)模型，通過(guò)正交試驗(yàn)分析了p，T，TOC和Ro對(duì)頁(yè)巖吸附量的敏感性，發(fā)現(xiàn)吸附量主要受壓力影響，其次受溫度和總有機(jī)質(zhì)碳質(zhì)量分?jǐn)?shù)的影響，而有機(jī)質(zhì)成熟度對(duì)其的影響最弱。

［1］張金川，金之鈞，袁明生.頁(yè)巖氣成藏機(jī)理和分布［J］.天然氣工業(yè)，2004，24（7）：15-18.

［2］熊健，劉向君，梁利喜.頁(yè)巖中超臨界甲烷等溫吸附模型研究［J］.石油鉆探技術(shù)，2015，43：（3）：96-102.

［3］李武廣，楊勝來(lái)，徐晶，等.考慮地層溫度和壓力的頁(yè)巖吸附氣含量計(jì)算新模型［J］.天然氣地球科學(xué)，2012，23（4）：189-194.

［4］岳長(zhǎng)濤，李術(shù)元，李林玥，等.荊門地區(qū)頁(yè)巖儲(chǔ)層吸附特征及吸附勢(shì)理論應(yīng)用［J］.石油實(shí)驗(yàn)地質(zhì)，2016，38（3）：346-353.

［5］GASPARIK M，GHANIZADEH A，BERTIER P，et al.High pressure methane sorption isotherms of black shales from the netherlands［J］.Energy&Fuels，2012，26（8）：4995-5004.

［6］樊靈，趙孟孟，殷川，等.基于BP神經(jīng)網(wǎng)絡(luò)的油田生產(chǎn)動(dòng)態(tài)分析方法［J］.斷塊油氣田，2013，20（2）：204-206.

［7］趙廷峰，周啟成，趙春艷，等.基于PCA-BP的鉆頭優(yōu)選方法與應(yīng)用［J］.斷塊油氣田，2016，23（2）：243-247.

［8］王小川.MATLAB神經(jīng)網(wǎng)絡(luò)43個(gè)案例分析［M］.北京：北京航空航天大學(xué)出版社，2013：20-33.

［9］童飛.基于BP神經(jīng)網(wǎng)絡(luò)的水上交通事故預(yù)測(cè)及MATLAB實(shí)現(xiàn)［D］.武漢：武漢理工大學(xué)，2005.

［10］王銀年.遺傳算法的研究與應(yīng)用［D］.無(wú)錫：江南大學(xué)，2009.

［11］謝中華.MATLAB統(tǒng)計(jì)分析與應(yīng)用［M］.北京：北京航空航天大學(xué)出版社，2010：208-216.

［12］肖鋼，唐穎.頁(yè)巖氣及其勘探開發(fā)［M］.北京：高等教育出版社，2012：25-37.

［13］張烈輝，唐洪明，陳果，等.川南下志留統(tǒng)龍馬溪組頁(yè)巖吸附特征及控制因素［J］.天然氣工業(yè)，2014，34（12）：63-69.

［14］BOWKER K A.Barnett shale gas production，fort worth basin：issues and discussion［J］.AAPG Bulletin，2007，91（4）：523-533.

［15］MONTGOMERY S L，JARVIEDM，BOWKERKA，etal.Mississippian barnett shale，fort worth basin，north-central Texas：gas-shale play with multitrillion cubic foot potential［J］.AAPG Bulletin，2005，90（6）：963-966.

［16］CHALMERS G R L，BUSTIN R M.Lower cretaceous gas shales in northeastern British Columbia，partⅡ：evaluation of regional potential gas resources ［J］.Bulletin of Canadian Petroleum Geology，2008，56（1）：22-61.

［17］武瑾，王紅巖，拜文華，等.渝東南龍馬溪組頁(yè)巖儲(chǔ)層特征及吸附影響因素分析［J］.斷塊油氣田，2013，20（6）：713-718.

［18］張超群，鄭建國(guó)，錢潔.遺傳算法編碼方案比較［J］.計(jì)算機(jī)應(yīng)用研究，2011，28（3）：819-822.

［19］劉春艷，凌建春，寇林元，等.GA-BP神經(jīng)網(wǎng)絡(luò)與BP神經(jīng)網(wǎng)絡(luò)性能比較［J］.中國(guó)衛(wèi)生統(tǒng)計(jì)，2013，30（2）：173-176.

［20］吳曉東，張玉豐，劉彥輝.蒸汽吞吐井注汽工藝參數(shù)正交優(yōu)化設(shè)計(jì)［J］.石油鉆探技術(shù)，2007，35（3）：1-4.