楊占偉，姜振學(xué)，梁志凱，吳 偉，王軍霞，宮厚健，李維邦，蘇展飛，郝綿柱

（1.中國石油大學(xué)（北京）油氣資源與探測國家重點實驗室，北京 102249；2.中國石油大學(xué)（北京）非常規(guī)油氣科學(xué)技術(shù)研究院，北京 102249；3.中國石油西南油氣田公司頁巖氣研究院，成都 610051；4.中國石油大學(xué)（北京）理學(xué)院，北京 102249）

0 引言

近年來，隨著非常規(guī)油氣領(lǐng)域探索的不斷深入，頁巖油氣勘探開發(fā)已然成為全球油氣資源勘探的熱點方向。有機質(zhì)碳含量（TOC）作為評價烴源巖生烴潛力及預(yù)測儲層油氣“甜點區(qū)”的關(guān)鍵性指標，對于頁巖氣勘探開發(fā)以及資源評價具有指導(dǎo)意義，但其值準確性受到地質(zhì)條件復(fù)雜程度及測定方法差異等一系列因素的影響。目前，盡管地球化學(xué)方法被認為是測量TOC 含量最受歡迎的方法，但有些缺陷可能會影響結(jié)果的可靠性，例如，不能覆蓋大范圍的深度，成本高且耗時。此外，由于樣品長時間暴露在空氣中，測量往往不準確，這增加了游離有機物被氧化和逸出的機會。利用測井資料預(yù)測TOC 含量也是常用手段，其常規(guī)方法有多元回歸法和ΔlogR法，但這2 種方法都存在一定缺陷［1-3］。多元回歸法將TOC 含量與一系列測井參數(shù)之間建立多元回歸關(guān)系，但因TOC含量往往受多種地質(zhì)因素的影響，且與各因素間存在著復(fù)雜的非線性關(guān)系，常規(guī)的回歸方法難以表達其內(nèi)在聯(lián)系，預(yù)測效果難以達到要求［4-7］。ΔlogR法通過將對數(shù)坐標的電阻率曲線與算術(shù)坐標的聲波時差曲線疊加，以非烴源巖段調(diào)準基線計算幅度差Δ logR［8-9］，該方法需要設(shè)定有機碳背景值（泥巖普遍含有一定量的有機碳）、成熟度參數(shù)及人為對準多個非烴源巖基線，操作復(fù)雜、誤差較大，對于異常點處理并沒有系統(tǒng)的標準，人為主觀性較強［10-11］。如今隨著人工智能領(lǐng)域與各個學(xué)科相互交叉融合，機器學(xué)習(xí)方法在地質(zhì)領(lǐng)域得到了一系列的應(yīng)用，通過引用不同的人工智能算法實現(xiàn)對烴源巖TOC 含量的準確預(yù)測［12-14］。常用的機器學(xué)習(xí)方法主要包括神經(jīng)網(wǎng)絡(luò)、支持向量機、極限學(xué)習(xí)等［15-17］。為了進一步提高預(yù)測精度，國內(nèi)外學(xué)者通過結(jié)合具體研究區(qū)特點，不斷優(yōu)化算法，降低預(yù)測誤差［18］。王貴文等［19］使用BP 神經(jīng)網(wǎng)絡(luò)模型對塔里木盆地臺盆區(qū)寒武系—奧陶系烴源巖進行了TOC 含量預(yù)測，取得了較好的效果；蔣德鑫等［20］對珠江口盆地陸豐凹陷文昌組烴源巖與TOC 含量建立多元回歸模型、人工神經(jīng)網(wǎng)絡(luò)模型和曲線疊合模型，并分析了3 種模型的預(yù)測效果及適用條件。Johnson 等［21］使用改進的神經(jīng)網(wǎng)絡(luò)模型對澳大利亞坎寧盆地?zé)N源巖TOC 含量進行了預(yù)測并達到了基本精度要求。上述方法雖取得了一定效果，但預(yù)測精度還有待提高，且在不同地區(qū)適用性不同，針對四川盆地尚未提出有效預(yù)測方法。

以川南五峰組—龍馬溪組為例，首先采集該地區(qū)多口井的測井曲線及實測TOC 含量數(shù)據(jù)進行建模，然后利用主成分分析法對收集的測井資料進行預(yù)處理，建立并訓(xùn)練基于BP 神經(jīng)網(wǎng)絡(luò)模型和基于梯度提升決策樹（Gradient Boosting Decision Tree，GBDT）的預(yù)測模型，并將這2種模型應(yīng)用于生產(chǎn)井的TOC 含量預(yù)測中驗證其準確性，以期為該區(qū)TOC含量預(yù)測提供新思路，為后續(xù)勘探開發(fā)提供依據(jù)。

1 地質(zhì)概況

四川盆地位于上揚子臺地西北部，在龍門山—大巴山臺緣坳陷和滇黔川鄂臺褶帶中間，是經(jīng)歷多期復(fù)雜構(gòu)造運動后形成的疊合深層盆地。在奧陶紀晚期，受到南部黔中隆起、西部川中隆起以及東部雪峰古隆起的相互作用，四川盆地逐漸由早中奧陶世廣闊海域逐漸轉(zhuǎn)變?yōu)楸粬|南西三面隆起包圍的半局限性海域，沉積主體也演變?yōu)榘刖窒逌\海相和深水—半深水陸棚亞相［22-23］。

根據(jù)川南長寧、瀘州等地區(qū)頁巖氣勘探開發(fā)資料顯示，該地區(qū)下古生界海相黑色富有機質(zhì)頁巖廣泛發(fā)育，4 500 m 以淺的五峰龍馬溪組成為黑色泥頁巖具有分布范圍廣、有機質(zhì)（TOC）含量高、成熟度高等一系列優(yōu)點，已經(jīng)成為了我國南方海相頁巖氣主力開發(fā)層系，2020 年產(chǎn)量已達200 億m3，該層位頁巖氣資源量達到3.7 萬億m3，可采資源量達200 億m3，具有極高的勘探開發(fā)價值。本次研究所涉及的區(qū)域主要為長寧和瀘州地區(qū)，所涉及的地層主要為五峰組—龍馬溪組頁巖層系。

2 TOC 含量預(yù)測模型

本次研究中井徑和地層電阻率等一系列測井?dāng)?shù)據(jù)、實測TOC 含量數(shù)據(jù)均來源于中國石油西南油氣田公司頁巖氣研究院，測試儀器為碳硫分析儀KLT-005，測試溫度為16～18 ℃。為改進常規(guī)地球化學(xué)方法測定TOC 含量時樣品長時間暴露在空氣中導(dǎo)致測量不準確的弊端，本次樣品密閉取心，封閉保存，并通過調(diào)研前人對于研究區(qū)TOC 含量的測定結(jié)果，與實測值進行比對，確保了測試的準確性。為了減弱個別參數(shù)不準確對模型精度的影響，采用主成分分析法對收集到的各項測井?dāng)?shù)據(jù)進行優(yōu)選，得到新的輸入變量，并分別代入BP 和GBDT預(yù)測模型中進行訓(xùn)練，調(diào)整各項參數(shù)使模型均方根誤差達到最小，利用建立的預(yù)測模型對不同頁巖氣井進行TOC 含量的預(yù)測。

2.1 主成分分析

主成分分析法的主要原理是利用降維的方法，以丟失較少的信息為前提，把多個相互關(guān)聯(lián)的參數(shù)轉(zhuǎn)化為幾個綜合性參數(shù)的一種多元統(tǒng)計方法。利用此方法轉(zhuǎn)化、生成的綜合性參數(shù)稱為主成分，其中每個主成分都是原始變量的線性組合，且各個主成分之間互不相關(guān)，利用此方法得到的主成分相比初始變量數(shù)據(jù)具有更為簡潔準確的性能［24-27］。

進行主成分分析，首先需要將原始數(shù)據(jù)進行標準化處理，以此來消除量綱所帶來的影響。假設(shè)進行主成分分析的指標有m個：x1，x2，x3，…，xm，評價對象共有n個，其中第1 個評價對象的第j個指標的取值為x1j，通過此方法可將原始數(shù)據(jù)排列為矩陣

將各個指標值xij轉(zhuǎn)換為標準指標：

其中

其次建立變量之間的相關(guān)系數(shù)矩陣R，其中

式中：rij=1，rij=rji，rij是第i 個指標與第j個指標的相關(guān)系數(shù)。

計算相關(guān)系數(shù)矩陣R的特征值λ1≥λ2…≥λm≥0，及對應(yīng)的特征向量u1，u2，um，其中

由特征向量組成m個新的指標變量。式中y1是第1 主成分，y2是第2 主成分，…ym是第m主成分。

最后確定主成分個數(shù)：

式中：αk為各個主成分的貢獻率；αP為主成分y1，y2，…，yp的累計貢獻率。

當(dāng)累計貢獻率αP大于或等于80%～85% 時，選擇其相對應(yīng)的前p個指標變量作為主成分分析得到的主成分，用這p個主成分替代m個初始變量指標，更加簡便準確。

收集川南長寧、瀘州等地區(qū)多口井的測井曲線及龍馬溪組17 口井627 組實測TOC 含量數(shù)據(jù)，分別將一系列不同測井參數(shù)和其對應(yīng)的TOC 含量進行相關(guān)性分析，優(yōu)選出相關(guān)性較高的7 個測井參數(shù)：聲波時差（AC）、井徑（CAL）、補償中子（CNL）、補償密度（DEN）、自然伽馬（GR）、無鈾伽馬（KTH）和地層深電阻率（RLLD），作為評價分析TOC 含量的影響因子［28］。

在創(chuàng)建預(yù)測模型時，僅使用任何單一測井參數(shù)均無法準確預(yù)測TOC 含量，為了達到準確得到TOC含量與各個測井曲線間的非線性映射關(guān)系的目的，對以上7 組測井參數(shù)進行主成分分析，計算相關(guān)系數(shù)矩陣，其中相關(guān)性系數(shù)越靠近1，則表明測井參數(shù)間的相關(guān)性越好，反之則相關(guān)性越差，結(jié)果如圖1所示。

圖1 輸入特征變量關(guān)聯(lián)度分析表Fig.1 Correlation analysis of input characteristic variables

分析計算得到各成分對于7 個測井參數(shù)信息的貢獻率和累計貢獻率（表1）。

表1 主成分特征值及貢獻率Table 1 Eigenvalues and contribution rate of the principal components

將貢獻率從高至低排列，前4 個主成分的累計貢獻率達到了86.505%，依據(jù)累計貢獻率大于85%的原則，選取前4 個主成分作為新的輸入變量，應(yīng)用于建立TOC 含量計算模型中。表2 為主成分的因子載荷矩陣，反映各指標對主成分載荷的相對大小和作用方向。

表2 主成分的因子載荷矩陣Table 2 Factor loading matrix of the principal components

利用表2 和4 個主成分初始特征值，計算得到新的輸入變量y1，y2，y3，y4公式：

經(jīng)過主成分分析得到了新的4 個綜合變量y1，y2，y3，y4，各個綜合變量之間相互獨立，每個綜合變量都包含前述7 個測井參數(shù)信息，但各有側(cè)重。根據(jù)各個系數(shù)大小可以看出，y1主要反映了AC和CAL的信息；y2主要反映了DEN和GR的信息；y3主要反映了KTH的信息；y4主要反映了CNL和RLLD的信息。

2.2 BP 神經(jīng)網(wǎng)絡(luò)模型

BP 神經(jīng)網(wǎng)絡(luò)是一種按照反向誤差進行傳播的多層前饋網(wǎng)絡(luò)，利用學(xué)習(xí)信號正向傳播和誤差逆向傳播雙重作用機制來對數(shù)據(jù)進行訓(xùn)練。BP 算法是以網(wǎng)絡(luò)誤差平方為目標函數(shù)、采用梯度下降法來計算目標函數(shù)的最小值［29-30］。通過正向傳播，隱含層可通過內(nèi)部計算獲取輸入向量的輸出值，輸出層則會在隱含層的基礎(chǔ)上進行計算，進而得到輸出值。當(dāng)進行反向傳播時，首先需要計算出輸出層所得到的輸出值，若獲得的輸出值達不到預(yù)期效果，則會重新將數(shù)據(jù)導(dǎo)入輸入層，再通過不斷修改輸入層與隱含層，隱含層與輸出層的連接強度和閾值，直到使誤差控制住預(yù)期范圍內(nèi)，訓(xùn)練則會停止［31］。經(jīng)過此方法訓(xùn)練好的神經(jīng)網(wǎng)絡(luò)則具有很強的非線性映射能力和柔性的網(wǎng)絡(luò)結(jié)構(gòu)。

本文建立的BP 神經(jīng)網(wǎng)絡(luò)主要包含輸入層、隱含層和輸出層三部分結(jié)構(gòu)。模型借助Matlab軟件編寫模型程序，并使用工具箱中newff 函數(shù)創(chuàng)建神經(jīng)網(wǎng)絡(luò)，再利用經(jīng)過主成分分析法降維后得到的4 個新的綜合性參數(shù)y1，y2，y3，y4作為模型的輸入變量，將實測得到的TOC 含量作為模型的輸出變量，隱含層的傳遞函數(shù)選擇tansig 函數(shù)，輸出層傳遞函數(shù)則選擇pureline 線性函數(shù)，隱函數(shù)節(jié)點數(shù)計算公式如下：

式中：k為樣本數(shù)，個；M為隱含層神經(jīng)元的數(shù)量，個；n為輸入層神經(jīng)元的個數(shù)；i為取值0～n的正整數(shù)，當(dāng)i＞M=0。

樣本數(shù)為627，選取隱含層神經(jīng)元個數(shù)M=15，最終確立的網(wǎng)絡(luò)結(jié)構(gòu)為4×15×1 的3 層BP 神經(jīng)網(wǎng)絡(luò)模型，為使算法達到較快的收斂速度和較高的訓(xùn)練精度，網(wǎng)絡(luò)訓(xùn)練采用rainlm 函數(shù)，神經(jīng)網(wǎng)絡(luò)模型具體相關(guān)參數(shù)如表3 所列。利用學(xué)習(xí)后的神經(jīng)網(wǎng)絡(luò)模型對測試集進行識別，通過計算誤差率，分析訓(xùn)練模型準確率，具體拓撲圖如圖2 所示。

表3 BP 神經(jīng)網(wǎng)絡(luò)模型相關(guān)參數(shù)Table 3 Relevant parameters of BP neural network model

圖2 BP 神經(jīng)網(wǎng)絡(luò)模型拓撲圖Fig.2 Topological diagram of BP neural network model

2.3 梯度提升決策樹（GBDT 模型）

GBDT 是一種具有出色預(yù)測能力的機器學(xué)習(xí)模型。通過回歸樹將計算值與目標值之間的殘差進行快速歸類分析，并利用逐步提升算法不斷減小殘差，使得計算值逐漸逼近目標值［32-34］。這種方法不僅能靈活處理各種類型的數(shù)據(jù)，還能在極短的調(diào)參時間下，達到較高的預(yù)測精度。由于回歸樹會對不同的殘差值進行不同處理，即使樣本中有錯誤樣點，訓(xùn)練結(jié)果也不會受到太大影響。此外，模型還可通過調(diào)整參數(shù)和使用更多回歸樹來挖掘出數(shù)據(jù)之間的對應(yīng)關(guān)系，使其具備較好的泛化能力。GBDT算法的核心是在多次基礎(chǔ)模型中，利用損失函數(shù)的負梯度值作為該輪基礎(chǔ)模型損失值的近似值，并通過這一近似值來構(gòu)建下一輪的基礎(chǔ)模型，能使目標函數(shù)的求解更為簡便。GBDT 算法的實現(xiàn)步驟如圖3 所示。

此模型利用Python 語言編寫模型程序，并借助機器學(xué)習(xí)庫內(nèi)Sklearn 模塊來進行建模分析。模型主要調(diào)節(jié)參數(shù)為Boosting 框架參數(shù)和弱學(xué)習(xí)器參數(shù)［35］。其中Boosting 框架的重要參數(shù)主要包括：最大迭代次數(shù)、權(quán)重縮減系數(shù)、損失函數(shù)等；弱學(xué)習(xí)器的主要參數(shù)包括：最大特征數(shù)、決策樹最大深度、葉子節(jié)點最少樣本數(shù)等。為了使得預(yù)測模型取得較好的效果，在建立模型前，需對模型參數(shù)進行調(diào)整。選擇參數(shù)時，迭代次數(shù)過小，容易出現(xiàn)欠擬合的情況；學(xué)習(xí)速率過小，則需要更為復(fù)雜的迭代過程和更大的計算量；過大的葉子節(jié)點深度會出現(xiàn)模型過擬合現(xiàn)象。因此，采用交叉驗證的方法來確定合理的參數(shù)體系，參數(shù)選取原則主要是利用準確率判斷模型擬合的好壞，通過不斷調(diào)整合適的參數(shù)以達到較高的準確率，模型具體參數(shù)如表4 所列。

表4 GBDT 模型相關(guān)參數(shù)Table 4 Relevant parameters of GBDT model

具體建模操作流程如下：

（1）輸入。

訓(xùn)練數(shù)據(jù)集：

損失函數(shù)：

式中：Y為實測TOC 質(zhì)量分數(shù)，%。

（2）初始化弱分類器。

對損失函數(shù)進行求導(dǎo)，令其導(dǎo)數(shù)為零，得到使損失函數(shù)達到極小值的常數(shù)c。因此初始化時，c值可取訓(xùn)練樣本標簽值的均值，即

計算損失函數(shù)的負梯度值，即：

式中：i為樣本數(shù)，個，i=1，2，…，627；m為迭代次數(shù)，即生成的基礎(chǔ)模型的個數(shù)，m=1，2，…，M。

利用上式得到的數(shù)據(jù)集（yi，rmi）來擬合下一輪的基礎(chǔ)模型，訓(xùn)練得到模型的回歸樹ft（y），計算每個葉子節(jié)點Rmj，j表示葉子節(jié)點數(shù)，j=1，2，…，J，并計算每個葉子節(jié)點Rmj的最佳擬合值rmj，使得其損失函數(shù)達到最小。

則

結(jié)合前m-1 輪基礎(chǔ)模型，可得到最終模型：

3 預(yù)測效果分析

對川南五峰組—龍馬溪組17 口井627 組實測數(shù)據(jù)分別按照訓(xùn)練樣本70%、驗證樣本15%、測試樣本15%的比例進行劃分，即439 組數(shù)據(jù)作為訓(xùn)練樣本，94 組數(shù)據(jù)作為驗證驗本、94 組數(shù)據(jù)作為測試樣本。BP 神經(jīng)網(wǎng)絡(luò)模型與GBDT 模型對于TOC含量的預(yù)測結(jié)果如圖4 所示，2 種模型對于TOC 含量的預(yù)測效果均較好，GBDT 模型的綜合效果要好于BP 神經(jīng)網(wǎng)絡(luò)，其預(yù)測值整體偏離TOC 含量實測值程度較小，相關(guān)性更高，達0.90 以上。

圖4 2 種模型訓(xùn)練及預(yù)測結(jié)果Fig.4 Training and prediction results of two models

分別計算兩者對于627組數(shù)據(jù)預(yù)測與真實TOC含量之間的誤差（圖5），GBDT模型整體誤差更小，絕大部分誤差值均在0.5 范圍內(nèi)，由此認為GBDT模型在預(yù)測TOC含量方面更具優(yōu)越性。

圖5 2 種模型預(yù)測誤差結(jié)果Fig.5 Prediction errors of two models

4 實際應(yīng)用

根據(jù)上述的2 種TOC 含量預(yù)測模型，利用川南長寧CNX202 井五峰組—龍馬溪組131組測井?dāng)?shù)據(jù)對頁巖氣TOC 含量進行預(yù)測。將經(jīng)過主成分分析的數(shù)據(jù)輸入訓(xùn)練好的模型中，如圖6 所示，BP 模型測試樣本的預(yù)測值與實測值相關(guān)性為0.761，GBDT模型測試樣本的預(yù)測值與實測值相關(guān)性為0.970。

圖6 2 種模型預(yù)測預(yù)測結(jié)果Fig.6 Prediction results of two models

以長寧CNX202 井為例，綜合對比上述2 種預(yù)測模型與傳統(tǒng)ΔlogR法對于TOC 含量的預(yù)測效果（圖7）可知，2 種預(yù)測模型的預(yù)測精度均高于傳統(tǒng)ΔlogR方法。

圖7 川南長寧地區(qū)CNX202 井五峰組—龍馬溪組TOC含量綜合對比圖Fig.7 Comparison between predicted and measured TOC content of Wufeng-Longmaxi Formation of well CNX202 in Changning area，southern Sichuan Basin

為了更加準確地評價3 種預(yù)測方法的優(yōu)劣，采用平均誤差、平均誤差率、均方根誤差等3 種誤差指標進行對比（表5）。根據(jù)對比結(jié)果顯示，BP 神經(jīng)網(wǎng)絡(luò)和GBDT 預(yù)測模型各誤差指標參數(shù)均小于傳統(tǒng)ΔlogR法，GBDT 預(yù)測模型表現(xiàn)出更高的準確度。

表5 3 種模型誤差指標對比Table 5 Comparison of error indexes of three models

5 結(jié)論

（1）優(yōu)選出相關(guān)性較高的聲波時差（AC）、井徑（CAL）、補償中子（CNL）、補償密度（DEN）、自然伽馬（GR）、無鈾伽馬（KTH）和地層深電阻率（RLLD）等7個測井參數(shù)作為評價TOC 含量的重要影響因子，利用主成分分析方法，將這7 個測井參數(shù)轉(zhuǎn)化為4個綜合指標，增加了計算精度和計算效率，減少計算冗余。

（2）在川南長寧地區(qū)五峰組—龍馬溪組采用BP 神經(jīng)網(wǎng)絡(luò)模型和GBDT 模型預(yù)測TOC 含量，訓(xùn)練樣本與測試樣本的預(yù)測值與實測值相關(guān)性均超過0.80，擬合效果均良好，GBDT 模型預(yù)測精度更高，其預(yù)測值與真實TOC 含量值的相關(guān)性達0.90以上，預(yù)測誤差更小，627 組樣本數(shù)據(jù)中絕大多數(shù)的誤差值在0.5 以內(nèi)，表現(xiàn)出更穩(wěn)定的優(yōu)勢。

（3）BP 神經(jīng)網(wǎng)絡(luò)模型和GBDT 模型均能較好的預(yù)測川南長寧地區(qū)五峰組—龍馬溪組TOC 含量，預(yù)測效果均好于常規(guī)ΔlogR方法，其中GBDT模型預(yù)測結(jié)果與實際測試值偏差更小，該方法在瀘州、長寧等地區(qū)有廣闊的應(yīng)用前景。