隋京坤，陳勝，鄭曉東*，胡天躍

(1.北京大學(xué)地球與空間科學(xué)學(xué)院，北京 100871； 2.中國(guó)石油勘探開(kāi)發(fā)研究院，北京 100083)

0 引言

由于不同時(shí)期沉積地層的巖石性質(zhì)變化一般較大，形成強(qiáng)反射界面；而同期地層內(nèi)不同小層的巖石性質(zhì)變化通常較小，形成弱反射界面。當(dāng)弱反射界面距強(qiáng)反射界面較近時(shí)，地震數(shù)據(jù)中的弱反射信號(hào)往往會(huì)被其附近的強(qiáng)反射信號(hào)所掩蓋，使地震數(shù)據(jù)的分辨能力降低，導(dǎo)致小層劃分困難。

為了識(shí)別弱反射界面的地震反射信號(hào)，傳統(tǒng)方法最常見(jiàn)思路為通過(guò)對(duì)地震記錄進(jìn)行分解后將強(qiáng)反射界面對(duì)應(yīng)的分量剔除，以突出小層之間界面對(duì)應(yīng)的弱反射信號(hào)。匹配追蹤算法是研究最多的技術(shù)之一，它對(duì)信號(hào)進(jìn)行自適應(yīng)的分解，實(shí)現(xiàn)強(qiáng)反射背景的匹配追蹤和識(shí)別[1-5]；子波分解和重構(gòu)技術(shù)[6-12]通過(guò)將地震數(shù)據(jù)分解為不同主頻的地震子波集合，將強(qiáng)反射信號(hào)對(duì)應(yīng)分量去除后，再進(jìn)行信號(hào)重構(gòu)，實(shí)現(xiàn)去除強(qiáng)反射的目的；基于奇異值分解的波形分解技術(shù)[13-15]也是相同的思路，只是信號(hào)分解的方法不同；田亞軍等[16]提出了采用AIDNN和U-net兩個(gè)深度神經(jīng)網(wǎng)絡(luò)的地震強(qiáng)反射剝離方法。以上方法都要求正確提取地震子波。如果預(yù)測(cè)的地震子波不準(zhǔn)確，“減去法”會(huì)存在子波殘留，引入虛假的弱信號(hào)。Guo等[17]提出使用反Q濾波方法對(duì)目標(biāo)層的弱反射信號(hào)和相位進(jìn)行補(bǔ)償以增強(qiáng)能量。

近年來(lái)，隨著計(jì)算機(jī)性能的大幅提升，神經(jīng)網(wǎng)絡(luò)算法被廣泛地應(yīng)用到石油勘探領(lǐng)域。一方面，將機(jī)器學(xué)習(xí)應(yīng)用到初至拾取[18-20]和層位拾取[21]，以提高工作效率；另一方面，由于神經(jīng)網(wǎng)絡(luò)強(qiáng)大的擬合能力，避免了求解極其復(fù)雜的非線(xiàn)性映射問(wèn)題，它在地震資料的處理、解釋中被廣泛應(yīng)用，尤其是在地震反演[22-25]中。長(zhǎng)短期記憶(LSTM)循環(huán)神經(jīng)網(wǎng)絡(luò)[26]是一種在時(shí)間序列分析中得到廣泛應(yīng)用的模型，已經(jīng)應(yīng)用于巖相預(yù)測(cè)[27]、速度譜自動(dòng)拾取[28]、測(cè)井曲線(xiàn)重構(gòu)[29]、地震數(shù)據(jù)初至拾取[30]、沉積微相的智能化識(shí)別[31]等方面。

針對(duì)弱信號(hào)增強(qiáng)中“減去法”的局限性，本文提出了一種基于“升弱降強(qiáng)”的新思路。通過(guò)對(duì)反射系數(shù)進(jìn)行分?jǐn)?shù)次冪運(yùn)算，獲得擬反射系數(shù)序列；用原始反射系數(shù)和擬反射系數(shù)計(jì)算合成地震記錄和擬合成地震記錄，生成訓(xùn)練樣本集；然后訓(xùn)練LSTM網(wǎng)絡(luò)，建立合成地震記錄與擬合成地震記錄的映射關(guān)系；最后將該網(wǎng)絡(luò)應(yīng)用于地震數(shù)據(jù)，增強(qiáng)弱地震反射信號(hào)。本文方法不是消除強(qiáng)地震信號(hào)，而是通過(guò)構(gòu)建一種冪次反射系數(shù)模型縮小弱反射系數(shù)與強(qiáng)反射系數(shù)的相對(duì)差異，從而突出弱反射信號(hào)。

1 方法原理

1.1 擬合成地震記錄

在地震剖面上，對(duì)應(yīng)強(qiáng)反射界面和弱反射界面的分別是強(qiáng)同相軸和弱同相軸，反射信號(hào)能量的大小取決于反射系數(shù)的大小及其組合。任意一個(gè)地震反射界面的反射系數(shù)與地震子波褶積可以得到一個(gè)同相軸。當(dāng)界面相距較近時(shí)，弱同相軸被強(qiáng)同相軸掩蓋的程度取決于它們振幅絕對(duì)值之比，而不是振幅的絕對(duì)值。為了減小強(qiáng)、弱同相軸的振幅的比值，本文提出利用分?jǐn)?shù)次冪運(yùn)算對(duì)原始反射系數(shù)序列進(jìn)行重構(gòu)，得到擬反射系數(shù)序列。

設(shè)原始反射系數(shù)序列為r=(r1，r2，…，rN)，擬反射系數(shù)序列為R=(R1，R2，…，RN)，則

(1)

式中：n=1，2，…，N，N為反射系數(shù)個(gè)數(shù)；j1/j2為分?jǐn)?shù)次冪，其中j1、j2為正整數(shù)。

假設(shè)反射系數(shù)|r1|<|r2|，當(dāng)冪指數(shù)j1/j2<1時(shí)，有

(2)

則強(qiáng)、弱反射系數(shù)絕對(duì)值的比值降低，可以減小強(qiáng)、弱反射系數(shù)絕對(duì)值之間的差異。當(dāng)冪指數(shù)j1/j2>1時(shí)，強(qiáng)、弱反射系數(shù)絕對(duì)值的比值增大，相對(duì)增強(qiáng)了強(qiáng)反射系數(shù)。本文應(yīng)用擬反射系數(shù)序列可以對(duì)反射系數(shù)進(jìn)行自由縮放。合成地震記錄s和擬合成地震記錄S可以分別表示為地震子波w與反射系數(shù)序列r、擬反射系數(shù)序列R的褶積，即

(3)

為保持處理前、后數(shù)據(jù)的值域的統(tǒng)一性，對(duì)合成記錄和擬合成地震記錄進(jìn)行標(biāo)準(zhǔn)化處理

(4)

式中：s′、S′分別為標(biāo)準(zhǔn)化后的合成記錄和擬合成記錄；μs、μS分別為s、S的均值；δs、δS分別為s、S的標(biāo)準(zhǔn)差。

1.2 長(zhǎng)短期記憶循環(huán)神經(jīng)網(wǎng)絡(luò)

循環(huán)神經(jīng)網(wǎng)絡(luò)(RNN)允許數(shù)據(jù)在重復(fù)網(wǎng)絡(luò)結(jié)構(gòu)中不斷循環(huán)，可以讓先前步驟中的預(yù)測(cè)結(jié)果加入到當(dāng)前步驟的運(yùn)算中，使RNN具備了長(zhǎng)期記憶功能。因此，不僅當(dāng)前步驟的輸入會(huì)影響RNN的輸出，過(guò)去所有步驟的輸入都會(huì)影響RNN的輸出。

圖1為一個(gè)標(biāo)準(zhǔn)的RNN模型，對(duì)于給定的序列x=(x1，x2，…，xM)(M為序列長(zhǎng)度)，RNN可以通過(guò)迭代

圖1 標(biāo)準(zhǔn)RNN及其展開(kāi)示意圖

ht=tanh(Wxhxt+Whhht+bh)

(5)

給出一個(gè)隱藏狀態(tài)ht。式中：xt為由x分割的第t個(gè)子向量；Wxh、Whh分別為輸入層、隱藏層的權(quán)重系數(shù)矩陣；bh為隱藏層的偏置向量。隱藏層以雙曲正切函數(shù)tanh為激活函數(shù)。

RNN存在梯度消失和梯度爆炸問(wèn)題，因此提出了LSTM網(wǎng)絡(luò)。LSTM網(wǎng)絡(luò)除了增加了一個(gè)tanh函數(shù)層之外，還同時(shí)增加了三個(gè)門(mén)層，即遺忘門(mén)、輸入門(mén)、輸出門(mén)[32](圖2)。

圖2 LSTM網(wǎng)絡(luò)單元結(jié)構(gòu)示意圖

遺忘門(mén)決定上一時(shí)刻的細(xì)胞狀態(tài)Ct-1中哪些被忘記，哪些被保留到當(dāng)前時(shí)刻的細(xì)胞狀態(tài)Ct，遺忘門(mén)的狀態(tài)ft的值域?yàn)?～1，可表示為

ft=σ(Wxfxt+Whfht-1+bf)

(6)

式中：Wxf、Whf分別為遺忘門(mén)中輸入層、隱藏層的權(quán)重系數(shù)矩陣；σ為Sigmoid激活函數(shù)；bf為遺忘門(mén)狀態(tài)f的偏置向量。

(7)

式中：Wxi、Whi分別為輸入門(mén)中輸入層、隱藏層的權(quán)重系數(shù)矩陣；Wxc、WhC分別為候選更新?tīng)顟B(tài)中輸入層、隱藏層的權(quán)重系數(shù)矩陣；bi、bC分別i、C的偏置向量。

當(dāng)前時(shí)刻單元狀態(tài)Ct的計(jì)算式為

(8)

式中“⊙”表示向量中元素按位相乘。

輸出門(mén)用于將當(dāng)前細(xì)胞狀態(tài)和它生成的LSTM網(wǎng)絡(luò)狀態(tài)輸出。輸出門(mén)的狀態(tài)和隱藏狀態(tài)為

(9)

式中：Wxo、Who分別為輸出門(mén)中輸出層、隱藏層的權(quán)重系數(shù)矩陣；bo為輸出門(mén)狀態(tài)o的偏置向量。

1.3 網(wǎng)絡(luò)訓(xùn)練

X=(x1,x2,…,xL)

(10)

(11)

式中：1≤l≤L；1

(12)

式中：Whp和bp分別為輸出單元的權(quán)重系數(shù)矩陣和偏置；t為1～L之間的整數(shù)。與之對(duì)應(yīng)的目標(biāo)輸出為

Y=(y1,y2,…,yL)

(13)

(14)

與xl對(duì)應(yīng)的擬合成地震記錄yl中有N-L+1個(gè)數(shù)據(jù)，預(yù)測(cè)值只有一個(gè)，因此只可以選擇數(shù)據(jù)中的任意1個(gè)樣點(diǎn)進(jìn)行擬合，定義為第k(1≤k≤L-1)個(gè)樣點(diǎn)

(15)

選用均方根誤差作為網(wǎng)絡(luò)的損失函數(shù)，定義為

(16)

地震數(shù)據(jù)服從廣義高斯概率分布，但不同數(shù)據(jù)的值域不同，嚴(yán)重影響LSTM網(wǎng)絡(luò)的輸出結(jié)果，因此將輸入LSTM網(wǎng)絡(luò)的地震數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化

(17)

式中：dn和d′n分別為標(biāo)準(zhǔn)化前、后地震數(shù)據(jù)d和d′的第n個(gè)樣點(diǎn)；μd、δd分別為d的均值和標(biāo)準(zhǔn)差。

2 模型實(shí)驗(yàn)

首先用正演模型(圖3a)驗(yàn)證本文算法的可行性。模型寬度為1000 m，道間距為10 m，包含三層，其中中間層的厚度由10 m線(xiàn)性增大到30 m，各層的速度和密度如表1所示。垂直入射情況下，上層與中間層的反射系數(shù)為0.05，中間層與下層的反射系數(shù)為0.2。本文選用主頻為35 Hz的 Ricker子波[33]與反射系數(shù)褶積合成地震記錄(圖3b)，兩個(gè)擬地震記錄剖面是地震子波分別與1/2次冪、1/3次冪擬反射系數(shù)序列的褶積(圖3c、圖3d)。為了保證次冪運(yùn)算前、后數(shù)據(jù)分布的一致性，需對(duì)地震記錄進(jìn)行標(biāo)準(zhǔn)化處理。

圖3 模型數(shù)據(jù)測(cè)試結(jié)果

表1 模型參數(shù)

相鄰兩道合成地震記錄對(duì)應(yīng)的中間層的厚度差為0.2 m。由圖3b可以看出，弱反射系數(shù)界面對(duì)應(yīng)的同相軸在第75道便無(wú)法識(shí)別。使用模型反射系數(shù)序列的1/2次冪與Ricker子波褶積，其合成記錄剖面在第45道就弱到無(wú)法識(shí)別(圖3c)。使用模型反射系數(shù)序列的1/3次冪與Ricker子波褶積，其合成記錄剖面在第25道才被強(qiáng)反射淹沒(méi)(圖3d)。與原始合成地震記錄(圖3b)相比，使用擬反射系數(shù)序列計(jì)算得到的地震合成記錄(圖3c、圖3d)中弱反射同相軸的能量得到增強(qiáng)，分辨率分別提高了6和10 m。

從圖3b、圖3c的地震剖面中的相同位置選取5道數(shù)據(jù)分別作為輸入數(shù)據(jù)和標(biāo)簽，制作訓(xùn)練集，其中選取的每次輸入的數(shù)據(jù)長(zhǎng)度為40 ms(20個(gè)樣點(diǎn))，將最終訓(xùn)練結(jié)果與標(biāo)簽的均方根誤差作為訓(xùn)練集誤差。用圖3b所示地震剖面作為驗(yàn)證集，將預(yù)測(cè)結(jié)果(圖3e)與對(duì)應(yīng)數(shù)據(jù)(圖3c)的均方根誤差作為驗(yàn)證集誤差。同樣對(duì)圖3b、圖3d所示地震剖面進(jìn)行相同處理，預(yù)測(cè)結(jié)果如圖3f所示。對(duì)比圖3e與圖3c、圖3f與圖3d可見(jiàn)，預(yù)測(cè)結(jié)果與模擬數(shù)據(jù)的誤差極小。

1/2、1/3次冪弱反射增強(qiáng)LSTM網(wǎng)絡(luò)訓(xùn)練集、驗(yàn)證集誤差曲線(xiàn)如圖4所示，可以看出式(15)中選用不同的k值(不同采樣點(diǎn))最終計(jì)算得到的均方根誤差不同，用分割窗口中部對(duì)應(yīng)元素作為標(biāo)簽，擬合誤差較??；1/2、1/3次冪LSTM網(wǎng)絡(luò)的訓(xùn)練集誤差相差不大，但后者的驗(yàn)證集誤差較大。

圖4 模型數(shù)據(jù)1/2、1/3次冪弱反射增強(qiáng)

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

3.1 實(shí)例一

工區(qū)內(nèi)目的層地震數(shù)據(jù)主頻為35 Hz，因此將Ricker子波的主頻也設(shè)定為35 Hz。工區(qū)附近有8口井，用7口井的數(shù)據(jù)制作1/2次冪弱反射增強(qiáng)LSTM網(wǎng)絡(luò)訓(xùn)練集(圖5a)，剩下的一口井的數(shù)據(jù)作為驗(yàn)證集(圖5b)。為了防止過(guò)擬合，設(shè)置當(dāng)誤差小于0.1后停止訓(xùn)練。訓(xùn)練集的均方根誤差為0.09，驗(yàn)證集均方根誤差為0.21。

圖5 測(cè)井?dāng)?shù)據(jù)測(cè)試結(jié)果

圖6 四川盆地應(yīng)用實(shí)例

3.2 實(shí)例二

為了進(jìn)一步展示本文方法的適用性，使用長(zhǎng)慶油田X工區(qū)地震數(shù)據(jù)進(jìn)行測(cè)試。該工區(qū)具有鮮明的河控三角洲前積反射特征(圖7a)。Forest1和Fo-rest2為不同前積期沉積體的分界面。C73為長(zhǎng)73底界對(duì)應(yīng)的反射，為該區(qū)強(qiáng)反射界面。由于本區(qū)相鄰的前積體部分區(qū)域的的彈性參數(shù)差別不大，出現(xiàn)部分區(qū)域反射信號(hào)弱、同相軸不連續(xù)等問(wèn)題。

相對(duì)于原始剖面(圖7a)，本文方法處理后的地震剖面(圖7b)中能量更加均衡，其中強(qiáng)反射信號(hào)減弱，弱反射信號(hào)明顯增強(qiáng)，原本因?yàn)榉瓷淠芰枯^弱而不能顯示的同相軸也突顯出來(lái)，原先不連續(xù)的同相軸更加連續(xù)(紅色箭頭所示)，為后期層位拾取提供了可靠的資料基礎(chǔ)。

圖7 長(zhǎng)慶油田實(shí)際數(shù)據(jù)應(yīng)用

4 結(jié)束語(yǔ)

本文提出的基于“升弱降強(qiáng)”思路的分?jǐn)?shù)次冪弱反射增強(qiáng)方法，縮小了弱反射信號(hào)與強(qiáng)反射信號(hào)的相對(duì)差異，避免了“減去法”中強(qiáng)子波殘余形成的虛假弱反射。本文使用測(cè)井反射系數(shù)合成地震記錄與擬合成地震記錄構(gòu)建訓(xùn)練集，采用LSTM網(wǎng)絡(luò)擬合二者的映射關(guān)系。模型數(shù)據(jù)測(cè)試表明：LSTM網(wǎng)絡(luò)模型可以有效地建立合成地震記錄與擬合成地震記錄之間的映射關(guān)系，并且具有泛化能力，可以直接拓展應(yīng)用于地震數(shù)據(jù)；避免了直接從地震數(shù)據(jù)中提取反射系數(shù)序列，既減小了計(jì)算量，又避免了提取的反射系數(shù)誤差大的問(wèn)題。

對(duì)于在原始地震剖面中小層界面對(duì)應(yīng)的同相軸難以識(shí)別的問(wèn)題，本文方法可以有效增強(qiáng)弱反射，提高地震數(shù)據(jù)對(duì)弱反射界面的識(shí)別能力，為小層劃分提供了一種可行方案，有著廣泛的應(yīng)用前景。