基于一維卷積神經網絡的地層智能識別方法

曹茂俊,崔欣鋒

(東北石油大學 計算機與信息技術學院,黑龍江 大慶 163318)

0 引 言

如今,隨著社會經濟的快速發展,石油的消費不斷增加。為了滿足社會發展需要,需要進一步勘探開發油氣資源,而地層智能劃分可以為油氣資源的持續開發節省人力和物力。

測井曲線包含了隨井深變化的儲層性質,其直接反映了地下的地質特征。因此,測井曲線的地層劃分對于確定地層位置具有重要意義[1]。長期以來,地層劃分的人工解釋方法以測井解釋人員的主觀判斷為主。人工解釋方法需要人工觀察測井曲線圖形形狀的變化以尋找對應地層的半幅點位置,這些方法比較依賴專業領域知識和專家經驗[2],而且由于不同專家的經驗和熟練程度不同,也會導致地層劃分結果的巨大差異,不僅需要大量的人力和物力,而且分層誤差也相對較高。隨著計算機技術的發展和劃分任務的日益繁重,地層劃分逐漸由定性識別轉向定量識別[3],數理統計、活度函數、小波變換、模糊識別、人工神經網絡等方法[4]已經被應用到利用測井曲線進行地層劃分的任務中。然而,這些方法大都只考慮個別曲線進行劃分,而且大都只考慮點與點之間的關系,而忽略了測井曲線隨深度變化上下連續的特點,即局部特征,而深度學習方法為該問題提供了非常具有前景的解決方案。地層劃分問題是根據測井曲線的形態識別出不同深度所屬的地層,因此可轉化為多元序列的分類問題[5]。卷積神經網絡(Convolutional Neural Network,CNN)在圖像識別領域已經是相對成熟的技術了,由于其強大的圖像分類能力也被應用于測井曲線地層的劃分[6],但考慮到測井曲線一般為具有一定規律但不規則的形態,對于地層劃分來說,輸出數據為某一個分層類別,無法實現測井曲線的地層劃分。測井曲線序列是一種二維的多元序列,不能直接采用傳統的CNN進行地層劃分,一維卷積神經網絡(1DCNN)在處理序列數據方面效果比較突出[7]。因此,該文選擇1DCNN作為地層劃分的基礎模型,通過對測井曲線進行IMK特征工程處理,即對原始測井曲線分別采用INPEFA技術[8]、中值濾波處理以及K-Means聚類[9]處理,更好地提取了原始曲線的特征,為1DCNN網絡提供了更加容易識別的特征,從而增加了模型的地層劃分效率。該文結合已有的網絡模型及測井曲線的特點,搭建了基于IMK-1DCNN的地層劃分網絡模型,從而達到對地層劃分的智能化方法研究的目的。

1 IMK-1DCNN模型

如圖1所示,基于IMK-1DCNN的測井曲線地層識別模型流程大體分為4個部分:

圖1 IMK-1DCNN地層識別流程

(1)特征構造。利用INPEFA技術、中值濾波、K-means算法分別構造INPEFA曲線特征、中值濾波特征、聚類特征,簡稱IMK特征構造,與原有的測井曲線共同構成新的特征集合。

(2)數據集劃分。將新的特征集合劃分為訓練集、驗證集和測試集。

(3)模型訓練。搭建一維卷積神經網絡模型,使用訓練集訓練IMK-1DCNN模型,使用驗證集進行誤差計算和權值更新。

(4)模型驗證。最終在測試集上驗證IMK-1DCNN的分層效率。

1.1 IMK特征構造

圖2展示了IMK特征構造的流程,為了處理測井曲線在采集過程中產生的噪聲,以及更好地提取原始測井曲線的地層特征,在模型訓練之前利用特征程對原始曲線進行了處理,采用IMK特征構造方法構造了新的特征,與原始測井曲線特征共同構成了新的特征矩陣。

圖2 IMK特征構造流程

上述特征構造過程中的INPEFA技術是在最大熵頻譜分析法(Maximum Entropy Spectral Analysis,MESA)基礎上形成的預測誤差濾波分析(Prediction Error Filter Analysis,PEFA)通過進一步發展形成合成預測誤差濾波分析(Integrated Prediction Error Filter Analysis,INPEFA)技術[10]。其原理是通過最大熵頻譜變換,可以由已知點的值在熵最大的原則下推算出下一個點的數值,這個值稱為最大熵頻譜分析估計值(MESA)。預測誤差濾波分析(PEFA)是通過計算每一個深度點的MESA預測值和對應的測井曲線實際值而得到數據差值,如下式所示:

PEFA=RV-MESA

(1)

其中,RV代表實際值,MESA代表最大熵譜分析估計值,得到的較大正值或負值為不連續點,可能為地層不連續點或旋回界面。對PEFA曲線做一個特定積分處理,就得到了合成預測誤差濾波分析曲線(INPEFA),其本質是反映泥巖含量的多少且基于測井曲線分析,因此它能夠顯示通常在原始測井曲線中顯示不出來的趨勢和模式以及反映沉積的旋回性變化,通常被用來輔助測井解釋人員識別隱藏的地層界面[11]。該文運用INPEFA技術利用數據集中的GR曲線構造INPEFA曲線,并將其作為模型的特征之一來構造樣本集。

1.2 數據集劃分

經過IMK特征構造之后得到了新的特征集合,數據集的合理劃分能夠很好地加速神經網絡的集成,這里將其按照8∶1∶1的比例劃分為訓練集、驗證集和測試集,其中訓練集用來訓練模型,驗證集用來調整并尋找最佳參數,最后在測試集上進行測試以驗證模型的泛化性能及準確率。

1.3 模型訓練及驗證

在準備好數據集后,需要搭建模型進行訓練,該文選用了善于處理序列數據的一維卷積神經網絡作為模型的基礎架構,處理數據集使之滿足模型的輸入并進行訓練,訓練時需要設置網絡的批大小、訓練輪次、早停策略、Dropout比例等參數以加速網絡的訓練。在訓練完成之后,要在測試集上進行測試來驗證模型的識別效果,并通過繪圖來對比模型的識別結果與地質專家的標記結果。

2 基于IMK-1DCNN的地層識別方法

2.1 1DCNN網絡基本原理

(2)

式中,hi表示輸入序列中(i)th個感受野的輸出[13]。可以利用1DCNN對測井曲線序列切片做一維卷積運算,提取該深度區間的層位特征信息,從而實現劃分地層類型的目的。

2.2 模型總體框架

設計的1DCNN網絡架構如圖 3所示。網絡中各層的含義解釋如下。

(1)輸入層:將輸入的與地層密切相關的測井曲線通過IMK特征構造新的特征矩陣,從而輸入到網絡中進行特征提取并學習。

(2)卷積層:卷積層使用的是適合處理序列數據的一維卷積(Conv1d),一維卷積常被用以在一維輸入信號上進行鄰域濾波。為了減少過擬合現象,還在每個卷積層后面引入了BN(Batch Normalization)層[14],它可以加快網絡的訓練和收斂的速度,控制梯度爆炸且防止梯度消失。

(3)池化層:池化層的本質是一個下采樣,因為數據經過卷積之后維度越來越高,而且特征圖沒有多大改變,在連續多個卷積之后,會產生很大的參數量,不僅會大大地增加網絡訓練的難度,還容易造成過擬合的現象,所以在卷積層之后放一個池化層對數據進行壓縮,降低維度,減少參數量。該文使用了兩種池化層,一種是最大池化層(Max Pooling),可以減少卷積層參數誤差造成估計均值的偏移,另一種是全局平均池化層(Global Average Pooling),可以減小特征提取時鄰域大小受限導致估計值方差增大造成的誤差[15]。

(4)Dropout層:為了減少參數量,加快模型的訓練速度,緩解過擬合現象,該文引入了Dropout層,通過隨機選擇其中的一些神經元并將其臨時隱藏,然后再進行本次訓練和優化。在下一次迭代中,繼續隨機隱藏一些神經元,如此直至訓練結束。由于是隨機丟棄,故而每一個mini-batch都在訓練不同的網絡。設置Dropout層的隱藏比例為0.3(隨機將30%的神經元臨時隱藏)。

(5)Dense層:即全連接層,在模型中起到“分類器”的作用,它會將學到的“分布式特征表示”映射到樣本標記空間。通過將經過網絡訓練的輸出向量輸入到包含8個神經元的Dense層,輸出最終的地層分類。

2.3 數據歸一化

由于不同測井曲線的量綱不同,為了消除量綱的影響,需要對數據進行歸一化。如果對未歸一化的數據直接進行建模,可能會導致模型對數值大的變量學習過多,而對數值小的變量訓練不夠充分,往往模型效果不佳。常用的數據歸一化的方法有最大最小值歸一化、小數定標法、定量特征二值化及均值方差標準化等。

選用的歸一化方法是最小值-最大值(Min-Max)方法,顧名思義,該方法是利用數據列中的最大值和最小值進行歸一化處理,處理完之后可將數據范圍縮放到[0,1]范圍內。計算方式為數據列中的數據與該列的最小值作差,再除以其極差。設一條測井曲線為X=(x1,x2,…,xn),則歸一化后的數據x'可由式(3)計算得出:

(3)

其中,x'代表數據值x歸一化的結果,max(X)和min(X)分別代表測井曲線數據列中的最大值和最小值。

2.4 滑動窗口采樣

深度學習需要大量樣本進行訓練,而單純地將測井曲線構造成獨立的測井圖像不能滿足深度學習網絡需求。在傳統的深度學習網絡中,往往會采用數據增強的方式來擴充數據集,使得網絡能夠學習到足夠的特征,比如旋轉、平移、縮放、隨機遮擋、水平翻轉和噪聲擾動等[16]。但是由于測井曲線具有時序性,其垂向的形態變化往往代表著某種沉積環境,從而導致對測井曲線矩陣利用傳統的數據增強方法會破壞其本身所代表的地層和沉積環境的含義。

為了增加樣本量,該文提出一種使用滑動窗口對測井曲線進行采樣的方法。設滑動窗口的長度為L,步長為S,數據生成的方式如圖4所示,將特征矩陣劃分為多個長度為L,寬度為特征數量的切片,每個切片對應的地層標簽yi由每個切片中的出現次數最多的標簽決定,同時為了避免切片中的地層標簽差異過大,步長S不宜設置過大,然后利用樣條插值[17]對每個切片進行處理以保證生成的輸入矩陣長度大于特征數目。這里使用樣條插值是為了在擴充切片長度的同時保證訓練曲線的平滑,最大程度地保留原始曲線的地層趨勢。

圖4 利用滑動窗口生成測井曲線切片示意圖

3 實際應用與結果分析

3.1 數據來源

實驗使用的數據來自四川盆地川中地區,共包含20口井的測井數據(數據分布見圖5)。四川盆地位于四川省龍門山斷裂處以東,四面環山,是中新生代以來發展起來的大型沉積盆地。該文研究的主要地層為茅口組和棲霞組,其中茅口組頂部與上二疊統龍潭組呈平行不整合接觸,上覆于棲霞組之上。茅口組厚度為119 m～208 m,平均厚237 m,巖性特征為灰色—深灰色亮晶生物碎屑灰巖、泥晶生物碎屑灰巖,泥質灰巖夾硅質結核[18]。棲霞組厚度為23 m～313 m,平均125 m,巖性特征為:下部為深灰色、灰黑色泥質瀝青質灰巖,底部常見眼球狀微晶灰巖夾炭質瀝青質頁巖;上部為淺灰色、灰褐色微—亮晶生屑灰巖。

3.2 樣本集構造

所用的樣本集共20口井,共92 284條數據。將其中18口井作為訓練集,剩余2口井作為測試集,驗證模型對地層的識別性能。在構造樣本集之前要對數據進行預處理。在實際測井中由于測井數據受多因素影響,可能導致數值缺失或出現異常值。因此,對于原始數據,首先要進行缺失值和異常值的處理,該文選擇了用眾數填充空值及異常值,實驗證明該方法簡單有效。

結合地質領域知識,以及1.1節的IMK特征構造方法,構造新的特征。最終選取自然伽馬(GR)、自然電位(SP)、聲波時差(AC)和INPEFA曲線、中值濾波曲線及K-Means聚類曲線共6種測井曲線作為訓練特征。其中,GR曲線主要測量地層的放射性,SP曲線主要測量在地層電化學作用下產生的電位,AC曲線反映了地層聲波速度,常被用來進行地層劃分。另外,在構造中值濾波曲線特征時為了更好地提取測井曲線的邊緣信息,使用了20以內的奇數濾波窗口構造了9條新的中值濾波曲線特征。另外以地質專家人工標注地層作為標簽,確定了8種地層,分別是石炭系、龍潭組、梁山組、茅三段、茅二段、茅一段、棲二段、棲一段。

圖6展示了某單井的樣本數據測井解釋圖,測井數據的采樣間隔為0.125 m,包含部分測井曲線和地層標簽。為了避免不同量綱的影響,利用2.1節中的最大-最小值方法對數據進行了歸一化,將所有測井曲線的數值范圍縮放到[0,1]內。

圖6 well52井的測井曲線及其測井解釋

數據預處理完成后,將歸一化后的數據排列到一起形成測井曲線矩陣,然后利用2.4節提到的滑動窗口采樣方法將其轉化為測井曲線切片,通過對長度L及步長S取3～50內的奇數進行對比實驗,發現L為5,S=5時識別效果最佳,故設置長度L為5,即每個曲線切片的實際長度為0.625 m,步長S為5。利用該方法最終得到16 635個長度為5,寬為14的二維灰度測井曲線矩陣切片,經統計,樣本地地層平均厚度為160個厚度點,因此對每個切片進行樣條插值時,選擇160作為目標樣本長度,最終將其維度變為160×14以保證生成的輸入矩陣長度大于特征數目,其中有14 971張屬于訓練井,1 664張屬于測試井,最終處理完成的部分訓練樣本如表1所示。其中INPEFA、GR_Med_3、GR_Kmeans分別表示INPEFA曲線、濾波窗口為3的中值濾波曲線以及K-Means聚類曲線,Label代表將原來的8種地層標簽映射為0～7的數值。

表1 well52井處理后的訓練樣本數據

3.3 實驗環境和網絡參數配置

經過數據集的構造,最終由20口井得到了83 018個實際測井曲線切片及其對應的地層標簽,每個切片大小為160×14的二維矩陣,然后將數據集按照9∶1的比例劃分為訓練集和測試集,得到74 716個訓練樣本和8 301個測試樣本。然后對模型進行訓練和驗證,并尋找最優參數。

對測井曲線進行地層劃分的1DCNN網絡模型環境配置為Window10操作系統,使用Python3.6作為編程語言,深度學習框架為Tensorflow,GPU使用GTX1650,顯存大小為4 GB。

詳細的網絡配置和參數設置如表2所示。

表2 網絡配置和參數設置

實驗設置批次大小(Batch Size)為512,訓練輪次(Epoch)為100并運用早停減少模型訓練時間,優化器使用AdamOptimizer,學習率設為0.001。

3.4 評價指標

引入多個評價指標,利用這些指標對預測結果做全方位評價,綜合反映模型分類效果。包括每一類的混淆矩陣(如表3所示)、整體的準確率(A)、精確率(P)、召回率(R)和F1值。

表3 混淆矩陣

(1)混淆矩陣。

(2)準確率。

準確率反映了對于給定樣本,模型預測正確的樣本數占總樣本的比例:

(4)

(3)精確率。

精確率反映了模型預測為正例的結果中真正例所占的比重:

(5)

(4)召回率。

召回率又稱查全率,反映了模型預測中樣本的正例有多少被正確預測正確:

(6)

(5)F1值。

F1值為精確率和召回率的調和平均值,由于精確率(P)和召回率(R)是一對相互矛盾的量,所以為了更好地評價分類器的性能,一般使用F1-Score作為評價標準來衡量分類器的綜合性能:

(7)

公式(4)～公式(7)中:Tp、FN、TN和FP分別為真正例、假負例、真負例和假正例樣本的數量,多分類問題可以轉換為二分類,對于正例對應的某一類別,其余的類別均為負例。

3.5 模型訓練結果

圖7是文中模型訓練過程中對應的訓練集和驗證集的準確率上升曲線和損失值下降曲線。在準確率方面,訓練集最高為0.98,驗證集最高為0.92;損失值方面,訓練集最低為0.17,驗證集最低為0.28。表4顯示了某測試井的分類結果對應的混淆矩陣,可以直觀地看出模型在不同地層上的分類效果,其中最高的精確率為茅三段,達到了1,表示全部識別正確,最低的精確率為梁山組,為0.57。表5顯示了各類評價指標。由結果可以看出,文中方法對茅口組地層的識別效果較佳,可能與其灰巖巖性有關,在曲線數值上表現為GR數值較低、SP數值適中,比較易于識別。對于其它地層,各類評價指標略低,可能是由于曲線形態具有相似性,但總體的識別準確率都在0.85以上。總之,模型在驗證集上的整體預測能力表現較為優異,對于所有測試井的樣本而言,整體平均準確率達到了92.82%。

表4 混淆矩陣

表5 各類評價指標

圖7 訓練過程中模型準確率-損失值曲線

同時,為了驗證文中方法的優越性,將其與雙向長短期記憶方法[19]和隨機森林方法[20]進行對比實驗,實驗結果如表6所示。結果表明,文中方法的分類準確率為92.82%,高于雙向長短期記憶網絡模型的88.23%和隨機森林方法的84.43%。

表6 各類地層識別方法結果對比

3.6 盲井的地層識別及劃分結果

為了驗證模型在未經過模型訓練的盲井上的地層識別及劃分結果,選取某單井數據進行處理并輸入模型驗證,該盲井的測試樣本構造過程同訓練集相同,在經過IMK特征構造、滑動窗口采樣(長度L=5,步長S=5)及樣條插值后得到1 920個160×14的測試樣本作為模型輸入。模型的輸出為1 920×1的一維向量,其中向量中的每個值代表對應測試樣本的預測值。圖8給出了模型對某測試井的地層劃分可視化結果。前4列為訓練特征,其中包含GR自然伽馬曲線,INPEFA曲線、AC(聲波時差)曲線、SP(自然電位)曲線。后兩列的柱狀圖分別代表人工標定的地層標簽以及文中模型對該井的地層劃分結果。對比發現,文中方法的地層劃分結果,與地質專家人工標定的結果比較相似,具有良好的泛化性能,能夠很好地識別地層進而劃分出不同的地層。

圖 8 某單井的測井解釋及模型地層識別結果

4 結束語

提出的基于IMK-1DCNN的地層智能識別方法,構建了3種特征,即INPEFA曲線特征、中值濾波特征和K-Means聚類特征。INPEFA曲線能很好地反映不同地層之間的變化趨勢;中值濾波特征在保留有效邊緣信息的同時很好地去除了測井曲線中的高頻噪聲;K-Means聚類特征則能很好地反映測井曲線的時空相關性,區分不同地層之間的邊界。將文中方法與雙向長短期記憶方法和隨機森林方法進行對比,文中方法的分類準確率為92.82%,高于雙向長短期記憶網絡模型的88.23%和隨機森林方法的84.43%。文中方法可以實現測井曲線的端到端學習。實驗結果表明,文中方法有助于提取測井曲線序列的隱藏特征,準確地劃分不同地層。該方法能有效識別地層,具有較強的泛化能力和魯棒性,為測井曲線地層劃分提供了一種新的途徑。