劉 可,尹 成,趙 虎,張運(yùn)龍

(1.西南石油大學(xué)地球科學(xué)與技術(shù)學(xué)院,四川成都610500;2.中國石油新疆油田分公司勘探開發(fā)研究院,新疆克拉瑪依834000)

神經(jīng)網(wǎng)絡(luò)作為一種由數(shù)據(jù)驅(qū)動的非線性算法,因其“端對端”的特點(diǎn)以及強(qiáng)大的學(xué)習(xí)能力,廣泛應(yīng)用于圖像處理、語音識別、自然語言處理等領(lǐng)域[1]。根據(jù)神經(jīng)網(wǎng)絡(luò)所承擔(dān)任務(wù)的數(shù)學(xué)特點(diǎn)可分為:聚類任務(wù)、分類任務(wù)、回歸任務(wù)、重構(gòu)任務(wù)等。

在地球物理領(lǐng)域,聚類任務(wù)主要包括地震相劃分[2],常使用自組織映射(SOM)神經(jīng)網(wǎng)絡(luò);分類任務(wù)包括儲層識別[3]、初至波識別[4]、斷層識別[5]等,常使用概率神經(jīng)網(wǎng)絡(luò)(PNN)和具有二維空間適應(yīng)能力的卷積神經(jīng)網(wǎng)絡(luò)(CNN)等;回歸任務(wù)包括巖石物理參數(shù)預(yù)測[6]、儲層參數(shù)預(yù)測[7]、測井曲線預(yù)測[8]等,往往使用前饋神經(jīng)網(wǎng)絡(luò)(BP)、徑向基神經(jīng)網(wǎng)絡(luò)(RBF)、循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)等;重構(gòu)任務(wù)包括疊前、疊后地震數(shù)據(jù)插值重構(gòu)[9]、地震數(shù)據(jù)去噪[10]等,主要使用卷積神經(jīng)網(wǎng)絡(luò)(CNN)、深度殘差網(wǎng)絡(luò)(RESnet)以及U型神經(jīng)網(wǎng)絡(luò)(U-net)等。

在回歸任務(wù)中的儲層參數(shù)預(yù)測研究方面,除神經(jīng)網(wǎng)絡(luò)方法外,樹模型類方法[11]、支持向量機(jī)類方法[12]、多元回歸類方法[13]等方法也被廣泛應(yīng)用。在傳統(tǒng)方法基礎(chǔ)上,近年提出了以深度神經(jīng)決策森林(DNDF)[14]、自適應(yīng)神經(jīng)決策樹(ANT)[15]為代表的多方法融合深度學(xué)習(xí)模型。雖然這些模型多適用于分類任務(wù),但只須將模型進(jìn)行少許變形即可適用于回歸任務(wù),如由深度神經(jīng)決策森林(DNDF)發(fā)展而來的深度回歸森林(DRF)[16]。本質(zhì)而言,這些方法或模型均可視為一個復(fù)雜的、非線性的函數(shù),即每個輸入都有確定唯一的輸出與之對應(yīng),因此這些方法或模型僅適用于復(fù)雜唯一性問題。

地震數(shù)據(jù)反演通常具有非唯一解,并且受噪聲等因素影響較大。由地震數(shù)據(jù)估計儲層參數(shù)具有明顯的不確定性,因此直接應(yīng)用上述方法或模型估計儲層參數(shù)會存在一定程度的誤差。利用概率分布解決不確定性問題效果較好。結(jié)合神經(jīng)網(wǎng)絡(luò)和概率分布兩者優(yōu)勢的混合密度網(wǎng)絡(luò)是目前解決不確定性問題的較好選擇。

混合密度網(wǎng)絡(luò)(MDN)最早由BISHOP等[17]和CORNFORD等[18]提出,該網(wǎng)絡(luò)不僅擁有神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)能力,還擁有概率分布描述不確定性問題能力,已經(jīng)成功應(yīng)用于具有不確定性的海面風(fēng)場反演。陳坤堂等[19]和馮倩等[20]聯(lián)合BP神經(jīng)網(wǎng)絡(luò)和混合密度網(wǎng)絡(luò)成功完成海面風(fēng)場反演,證明了該網(wǎng)絡(luò)的有效性。SHAHRAEENI等[21]將該網(wǎng)絡(luò)引入巖石物理反演。殷習(xí)容等[22]通過該網(wǎng)絡(luò),利用橫縱波速度反演含水飽和度、泥質(zhì)含量和孔隙度等巖石物理參數(shù),同時測試該網(wǎng)絡(luò)抗噪聲能力,證明該網(wǎng)絡(luò)對巖石物理參數(shù)反演的有效性。但混合密度網(wǎng)絡(luò)在儲層參數(shù)預(yù)測方面的研究文獻(xiàn)較少。

由于地震數(shù)據(jù)與儲層參數(shù)的不確定性關(guān)系,通過混合密度網(wǎng)絡(luò)得到的儲層參數(shù)“觀測”概率分布往往呈現(xiàn)“多峰”,為提升儲層參數(shù)預(yù)測的精確度,需要對該概率分布進(jìn)行約束。由于地下介質(zhì)具有連續(xù)分布的、空間橫向漸變性的特征,因而可以假設(shè)相鄰地震道儲層參數(shù)之間具有一定的相似性[23]。以該假設(shè)為基礎(chǔ),可以從已知井點(diǎn)的準(zhǔn)確值出發(fā),建立相鄰道的“估計”概率分布來約束“觀測”概率分布。相鄰道“估計”概率分布的約束,使得“觀測”概率分布的“多峰”性得到一定程度壓制。以往研究表明受“過擬合”等影響,單獨(dú)使用神經(jīng)網(wǎng)絡(luò)得到的結(jié)果往往分布雜亂,空間橫向連續(xù)性差[24]。因此反演中需考慮相鄰道儲層參數(shù),使得相鄰道之間儲層參數(shù)過渡更為平滑。

本文在混合密度網(wǎng)絡(luò)的基礎(chǔ)上,以相鄰道橫向漸變?yōu)榧s束條件,提出了基于空間近似概率約束的混合密度網(wǎng)絡(luò)儲層參數(shù)預(yù)測方法。以砂體厚度預(yù)測為例,分別利用本文方法對理論模型和實際資料進(jìn)行了試驗,同時與BP神經(jīng)網(wǎng)絡(luò)和混合密度網(wǎng)絡(luò)的預(yù)測結(jié)果進(jìn)行了對比,并結(jié)合預(yù)測結(jié)果分析了3種方法在地震儲層參數(shù)預(yù)測方面的效果。

1 算法原理

1.1 混合密度網(wǎng)絡(luò)

1.1.1 混合模型

將多個不同頻率的正弦或余弦信號進(jìn)行組合可以得到一個復(fù)雜信號。同理,將多個簡單概率分布進(jìn)行組合也可以得到一個復(fù)雜概率分布,即混合概率分布模型。假設(shè)簡單概率分布為高斯概率分布(N),在地震數(shù)據(jù)反演中地震屬性(s)已知,則儲層參數(shù)(r)的概率分布可表述如下:

(1)

1.1.2 混合密度網(wǎng)絡(luò)結(jié)構(gòu)

混合密度網(wǎng)絡(luò)數(shù)學(xué)模型來自于混合概率分布模型。與BP神經(jīng)網(wǎng)絡(luò)相似,混合密度網(wǎng)絡(luò)一般由4部分組成:輸入層,隱藏層,預(yù)處理層和輸出層(圖1),但其輸出層輸出的值并不是儲層參數(shù),而是一系列與混合概率分布模型相關(guān)的參數(shù)[25]。

圖1 混合密度網(wǎng)絡(luò)的一般網(wǎng)絡(luò)結(jié)構(gòu)

對于上述混合密度網(wǎng)絡(luò)結(jié)構(gòu),其前向過程可分為:輸入層-隱藏層、隱藏層-預(yù)處理層以及預(yù)處理層-輸出層3個部分。

輸入層-隱藏層與BP神經(jīng)網(wǎng)絡(luò)一致,采取全連接方式計算:

(2)

隱藏層-預(yù)處理層也采取全連接方式計算,但不使用激活函數(shù):

(3)

根據(jù)高斯混合復(fù)雜概率分布的限制條件,對預(yù)處理層-輸出層不同類別節(jié)點(diǎn)采取不同處理方式:各高斯分布權(quán)值為正數(shù),權(quán)值之和為1以及方差為正數(shù)。則由混合密度網(wǎng)絡(luò)得到的“觀測”概率分布(Pobs)可表述如下:

(4)

根據(jù)最大似然估計,混合密度網(wǎng)絡(luò)的誤差函數(shù)L一般為:

(5)

式中:Xd和Yd分別為訓(xùn)練集中第d個訓(xùn)練樣本的輸入和輸出;D為參與訓(xùn)練樣本數(shù)目。

根據(jù)誤差函數(shù),使用BP算法進(jìn)行訓(xùn)練?；旌厦芏染W(wǎng)絡(luò)可視為BP神經(jīng)網(wǎng)絡(luò)的一種概率上的推廣[25],在單獨(dú)使用混合密度網(wǎng)絡(luò)時,往往將預(yù)測區(qū)間內(nèi)概率最大處作為預(yù)測輸出。

1.2 相鄰道估計模型

以往神經(jīng)網(wǎng)絡(luò)僅僅單獨(dú)應(yīng)用于每一地震道,未考慮地層介質(zhì)在橫向上具有的漸變性,由此可假設(shè)相鄰地震道儲層參數(shù)之間具有一定相似性[23],在此基礎(chǔ)上建立相鄰道“估計”數(shù)學(xué)模型。對于二維地震數(shù)據(jù),每一地震道均具有唯一的一維坐標(biāo)n,其中n為該地震道的inline號。若第n道儲層參數(shù)rn未知,而第n-1道儲層參數(shù)rn-1已知,則rn的“估計”概率分布(Pest)可表述如下:

(6)

(7)

式中:sn和sn-1為第n道和第n-1道的地震屬性;C為整個工區(qū)地震屬性組成矩陣的協(xié)方差矩陣;const為常數(shù),表示整個工區(qū)的基礎(chǔ)差異。

當(dāng)?shù)趎-1道儲層參數(shù)rn-1不確定,僅知其概率分布為bel(rn-1)時,則rn的“估計”概率分布(Pest)可表述如下:

(8)

實際地震數(shù)據(jù)多為三維數(shù)據(jù),對此每一地震道均具有唯一的二維坐標(biāo)(n,m),其中n和m分別代表其inline號和crossline號。仿照卡爾曼濾波相鄰估計方法,將其相鄰估計方法應(yīng)用于二維空間(圖2),其中藍(lán)色表示第一輪相鄰估計,紅色表示第二輪相鄰估計。卡爾曼濾波相鄰估計多用于一維時間域,是單向、單點(diǎn)且不交叉的(圖2a),而二維空間域相鄰估計是多向、集體且存在交叉的(圖2b)。

1.3 儲層參數(shù)預(yù)測

對于儲層參數(shù)未知的(n,m)道,使用混合密度網(wǎng)絡(luò)可以得到該道儲層參數(shù)rn,m的“觀測”概率分布(Pobs),同時根據(jù)相鄰道“估計”關(guān)系((6)式和(8)式)得到(n,m)道儲層參數(shù)rn,m的“估計”概率分布(Pest)。由于二維空間域相鄰估計具有交叉的特點(diǎn)(圖2b),為簡便表達(dá),設(shè)(n,m)道的相鄰道隸屬于集合w,將兩種概率分布進(jìn)行融合,得:

圖2 兩種相鄰估計a 卡爾曼濾波相鄰估計; b 本文二維相鄰估計

(9)

式中;η為歸一化參數(shù);K′為融合后概率分布的高斯分布數(shù)目。

顯然,Pmerge擁有的高斯分布數(shù)目遠(yuǎn)大于Pest和Pobs擁有的高斯分布數(shù)目。為運(yùn)算簡便,采取近似方法,即在Pmerge中選取具有最大的K″個權(quán)值的高斯分布,重新歸一化后組合,得到(n,m)道儲層參數(shù)rn,m的概率分布bel(rn,m):

(10)

以井位處作為起始位置,根據(jù)二維空間外推方案(圖2b),以及公式(6)、公式(8)、公式(9)及公式(10)求取相鄰道儲層參數(shù)概率分布,再以這些相鄰道為起始位置重復(fù)以上過程,直到完成全三維工區(qū)地震道的計算,即完成了全工區(qū)儲層參數(shù)概率預(yù)測。單井儲層參數(shù)概率預(yù)測流程如圖3所示。

圖3 單井儲層參數(shù)概率預(yù)測流程

實際工區(qū)有多口井,往往會產(chǎn)生多個不同初始點(diǎn)的全工區(qū)儲層參數(shù)概率預(yù)測結(jié)果。概率理論下的相鄰道估計過程中會不斷積累誤差,為了在減少誤差的同時融合不同初始點(diǎn)的反演結(jié)果,根據(jù)預(yù)測道與各初始位置的距離進(jìn)行加權(quán)平均,得到加權(quán)后的預(yù)測道儲層參數(shù)概率分布belavg(rn,m):

(11)

(12)

2 理論模型試驗

設(shè)計的理論模型如圖4a所示,包括砂泥巖薄互層以及上、下兩套圍巖,其中目的層厚度50m。砂體分為類菱形、類平行四邊形以及凹凸形3類,地震波速度2350m/s,密度2.1g/cm3;泥巖地震波速度2420m/s,密度2.2g/cm3;圍巖地震波速度2500m/s,密度2.0g/cm3。圖4b為采用主頻為50Hz的零相位雷克子波進(jìn)行正演模擬的地震數(shù)據(jù)。

對于該地震數(shù)據(jù)(圖4b),首先根據(jù)頂、底界面提取了振幅類、瞬時類、頻譜類等八大類60種地震屬性[26],同時采用粗糙集(RS)[27]、三參數(shù)敏感屬性分析(SDC)[28]等方法,聯(lián)合挑選出7種敏感地震屬性。然后隨機(jī)選取其中30%,20%,10%以及5%道數(shù)據(jù)作為訓(xùn)練數(shù)據(jù),最終得到如圖5至圖8所示的不同訓(xùn)練樣本數(shù)下砂體厚度(d)的預(yù)測結(jié)果及絕對誤差曲線。

圖4 理論模型(a)及正演地震數(shù)據(jù)(b)

圖5 3種方法的砂體厚度預(yù)測(a)及絕對誤差曲線(b)(30%訓(xùn)練數(shù)據(jù))

圖6 3種方法的砂體厚度預(yù)測(a)及絕對誤差曲線(b)(20%訓(xùn)練數(shù)據(jù))

對比本文方法(紅色曲線)、BP神經(jīng)網(wǎng)絡(luò)(藍(lán)色曲線)以及混合密度網(wǎng)絡(luò)(綠色曲線)這3種方法的預(yù)測結(jié)果可以看出,隨著訓(xùn)練樣本的減少,3種方法預(yù)測曲線越來越偏離實際厚度曲線。但是在相同條件下,本文方法與實際厚度曲線更為接近,橫向過渡更為連續(xù)。

考慮到單次隨機(jī)選取訓(xùn)練數(shù)據(jù)不具有代表性,進(jìn)行50次訓(xùn)練數(shù)據(jù)的Monte Carlo隨機(jī)挑選,進(jìn)行3種方法對比預(yù)測試驗。選擇預(yù)測砂體厚度曲線與實際砂體厚度曲線的平均誤差和相關(guān)系數(shù)作為預(yù)測效果的評價指標(biāo)。采用3種方法分別計算50次預(yù)測試驗的平均誤差和平均相關(guān)系數(shù),結(jié)果如圖9所示。

由圖9a可見,隨著訓(xùn)練樣本的減少,BP神經(jīng)網(wǎng)絡(luò)預(yù)測的平均誤差由0.6m上升到3.1m;混合密度網(wǎng)絡(luò)平均誤差由1.3m上升到3.9m;本文方法平均誤差由0.2m上升到2.1m。由圖9b可見,隨著訓(xùn)練樣本的減少,BP神經(jīng)網(wǎng)絡(luò)平均相關(guān)系數(shù)由0.98下降至0.83;混合密度網(wǎng)絡(luò)平均相關(guān)系數(shù)由0.89下降至0.67;本文方法平均相關(guān)系數(shù)由0.99下降至0.92。綜合圖5至圖9可見,在相同訓(xùn)練樣本下,單純混合密度網(wǎng)絡(luò)預(yù)測效果最差,BP神經(jīng)網(wǎng)絡(luò)較好,本文方法最好。

圖7 3種方法的砂體厚度預(yù)測(a)及絕對誤差曲線(b)(15%訓(xùn)練數(shù)據(jù))

圖8 3種方法的砂體厚度預(yù)測(a)及絕對誤差曲線(b)(5%訓(xùn)練數(shù)據(jù))

圖9 50次Monte Carlo試驗結(jié)果統(tǒng)計a 50次試驗平均誤差; b 50次試驗平均相關(guān)系數(shù)

3 實際數(shù)據(jù)應(yīng)用

基于上述實驗結(jié)果,將本文方法應(yīng)用于X工區(qū)某實際地震數(shù)據(jù)。從工區(qū)提取的多種屬性中,同樣通過粗糙集(RS)、三參數(shù)敏感屬性分析(SDC)等方法聯(lián)合優(yōu)選出8種敏感屬性(圖10至圖13)。將本工區(qū)的56口井分為41口訓(xùn)練井(圖10至圖13中含Tr前綴井)和15口測試井(圖10至圖13中含Te前綴井)。

圖10 基于敏感性分析得到的包絡(luò)差(a)和累計能量比(b)

圖11 基于敏感性分析得到的工區(qū)能量半時間(a)和變異系數(shù)(b)

圖12 基于敏感性分析得到的波形對稱度(a)和有效帶寬(b)

圖13 基于敏感性分析得到的平均頻率(a)和平均振幅(b)

在訓(xùn)練井位置處抽取地震道敏感屬性與訓(xùn)練井組合成訓(xùn)練集,訓(xùn)練BP神經(jīng)網(wǎng)絡(luò)和混合密度網(wǎng)絡(luò),將訓(xùn)練完成的網(wǎng)絡(luò)應(yīng)用于整個工區(qū),得到工區(qū)的砂體厚度分布預(yù)測結(jié)果(圖14a和圖14b)。再采用本文方法對整個工區(qū)進(jìn)行預(yù)測,得到工區(qū)的砂體厚度分布(圖14c)。從圖14可以看出,BP神經(jīng)網(wǎng)絡(luò)和混合密度網(wǎng)絡(luò)預(yù)測結(jié)果砂體厚度分布較為雜亂,橫向不連續(xù),而本文方法預(yù)測的砂體厚度橫向更為連續(xù),分布集中。

圖14 3種方法預(yù)測的砂體厚度平面分布a BP神經(jīng)網(wǎng)絡(luò)預(yù)測結(jié)果; b 混合密度網(wǎng)絡(luò)預(yù)測結(jié)果; c 本文方法預(yù)測結(jié)果

選取該工區(qū)測試井位處砂體預(yù)測厚度與已知的測試井位處的砂體厚度值進(jìn)行對比,其結(jié)果如圖15 所示。統(tǒng)計各個測試井?dāng)?shù)據(jù),BP神經(jīng)網(wǎng)絡(luò)預(yù)測的平均誤差為1.08m,誤差大于1m的井有7口,最大誤差為3.6m,位于井Te2處;混合密度網(wǎng)絡(luò)平均誤差為1.59m,誤差大于1m的井有7口,最大誤差為5.8m,位于井Te13處;本文方法平均誤差為0.7m,誤差大于1m的井有3口,最大誤差為2.5m,位于Te3處。由此可見本文方法預(yù)測的誤差最低。

圖15 測試井點(diǎn)三種方法的預(yù)測結(jié)果a 3種方法預(yù)測的砂體厚度與井上實際厚度; b 3種方法的絕對誤差

綜上所述,綜合判斷本文方法預(yù)測結(jié)果能夠反映工區(qū)砂體厚度的分布情況。

4 結(jié)論

傳統(tǒng)神經(jīng)網(wǎng)絡(luò)可視為一個復(fù)雜的、非線性的函數(shù),即每個輸入都有確定的輸出與之對應(yīng)。但是利用地震數(shù)據(jù)估計儲層參數(shù)具有明顯的不確定性,因而傳統(tǒng)神經(jīng)網(wǎng)絡(luò)不能很好適應(yīng)。對地震預(yù)測的不確定性問題,最好使用概率分布來間接解決。結(jié)合了神經(jīng)網(wǎng)絡(luò)和概率分布優(yōu)勢的混合密度網(wǎng)絡(luò),雖然從理論上能很好地解決多解性,但實際上的表現(xiàn)卻不如BP神經(jīng)網(wǎng)絡(luò),這是由于混合密度網(wǎng)絡(luò)需要訓(xùn)練的參數(shù)過多,以及概率分布本身所具有的不確定性造成的。

本文提出的基于空間近似概率約束的混合密度網(wǎng)絡(luò)砂體厚度預(yù)測方法,結(jié)合了混合密度網(wǎng)絡(luò)和橫向漸變假設(shè),借鑒地質(zhì)統(tǒng)計學(xué)從井出發(fā)的觀點(diǎn),實現(xiàn)了概率分布下的儲層參數(shù)最優(yōu)估計。通過理論數(shù)據(jù)和實際資料的應(yīng)用測試,對比了BP神經(jīng)網(wǎng)絡(luò)、混合密度網(wǎng)絡(luò)以及本文方法的預(yù)測效果。結(jié)果表明,本文方法具有更高的精確性和較大的應(yīng)用潛力。

本文方法中,地震道儲層參數(shù)概率分布需要由各已知井外推得到,而各已知井外推是相互獨(dú)立的運(yùn)算,所以需要較大次數(shù)的迭代運(yùn)算才能實現(xiàn)所有已知井的外推,如何實現(xiàn)多井同步外推是下一步研究的重點(diǎn)。