潘樹林,閆 柯,楊海飛,蔣從元,秦子雨

(1.西南石油大學(xué)地球科學(xué)與技術(shù)學(xué)院,四川成都610500;2.長慶油田長北作業(yè)分公司,陜西榆林719000;3.四川職業(yè)技術(shù)學(xué)院電子電氣工程系,四川遂寧629000)

稀疏脈沖反褶積可利用有限帶寬地震記錄反演得到地下反射系數(shù)序列,是一種常用的地震資料高分辨率處理方法[1]。該方法在傳統(tǒng)的最小二乘反褶積的基礎(chǔ)上,加一個(gè)稀疏先驗(yàn)約束項(xiàng),通過求解帶約束的目標(biāo)函數(shù),進(jìn)而得到高分辨率的處理結(jié)果。求解該目標(biāo)函數(shù)的方法主要有以下兩類。第一類為匹配追蹤類貪婪迭代算法[2-5],直接使用L0范數(shù)作為稀疏約束項(xiàng),通過枚舉匹配迭代的方法求解未知參數(shù)。該類方法絕對收斂,但運(yùn)算量巨大。第二類為CHEN等[6]提出的基追蹤類凸優(yōu)化算法,此類方法利用L1范數(shù)代替L0范數(shù),以便使用線性規(guī)劃來求解。該類方法僅在保證期望足夠稀疏時(shí),才能獲得理想的反演結(jié)果。為了進(jìn)一步改善基追蹤類方法對數(shù)據(jù)稀疏度的要求,張繁昌等[7]、曹靜杰等[8]通過引入地質(zhì)模型約束和非凸Lp范數(shù)正則化方法對其進(jìn)行了改進(jìn)。基追蹤類方法結(jié)構(gòu)簡單,便于實(shí)現(xiàn),計(jì)算效率高,是當(dāng)前較主流的求解方法。應(yīng)用該類方法時(shí),要得到好的反褶積效果,地震子波必須準(zhǔn)確,但在實(shí)際生產(chǎn)中,難以從地震數(shù)據(jù)中提取出非常準(zhǔn)確的地震子波,因而制約了稀疏脈沖反褶積方法的應(yīng)用效果。

針對上述問題,國內(nèi)外學(xué)者進(jìn)行了很多研究。第一類方法為DONOHO等[9]采用過完備字典的方法,通過將地震子波擴(kuò)充為包含多個(gè)頻率成分的過完備字典,消除了單個(gè)頻率成分子波帶來的缺點(diǎn)。這類算法最大的特點(diǎn)在于利用多個(gè)頻率成分的過完備子波字典代替單個(gè)子波字典,求解結(jié)果的精度取決于子波字典的準(zhǔn)確性。第二類方法為字典學(xué)習(xí)的方法[10-12],通過學(xué)習(xí)修改的方式來更新地震子波,但是該類方法只是利用誤差來交替修改反射系數(shù)與地震子波,而不是直接通過誤差反傳來修改,因此計(jì)算誤差較大。

迭代閾值收縮算法(iterative shrinkage-thresholding algorithm,ISTA)是基追蹤類方法中一種重要的求解算法,其收斂性受輸入子波的影響。循環(huán)神經(jīng)網(wǎng)絡(luò)(recurrent neural network,RNN)是一種遞歸神經(jīng)網(wǎng)絡(luò),可以根據(jù)輸入數(shù)據(jù)進(jìn)行自主學(xué)習(xí),改善處理結(jié)果,在很多領(lǐng)域取得了較好的應(yīng)用效果。本文通過對比ISTA結(jié)構(gòu)和RNN網(wǎng)絡(luò)結(jié)構(gòu),結(jié)合兩種方法提出了一種類RNN的改進(jìn)ISTA稀疏脈沖反褶積方法,該算法是將ISTA轉(zhuǎn)換為類RNN網(wǎng)絡(luò)結(jié)構(gòu),并引入RNN中BPTT的誤差反傳的思想修改地震子波,克服了子波不準(zhǔn)確造成的反演不收斂問題。最后用理論模型數(shù)據(jù)和吐哈盆地雁木西地區(qū)的實(shí)際資料對本文方法進(jìn)行了應(yīng)用測試。

1 類RNN的改進(jìn)ISTA原理

1.1 ISTA

稀疏脈沖反褶積可以表達(dá)為:

(1)

其中,y表示地震記錄向量,W表示由地震子波組成的矩陣,x表示待求解反射系數(shù)向量,λ為正則化參數(shù)。

BECK等[13]提出了ISTA,用于求解具有L2-L1結(jié)構(gòu)的復(fù)合函數(shù)。利用迭代軟閾值的方法求解(1)式,可得:

(2)

(3)

ISTA迭代流程如圖1所示。

圖1 ISTA迭代流程示意

顯然,ISTA求解反射系數(shù)的關(guān)鍵在于地震子波矩陣We的準(zhǔn)確確定,但是在實(shí)際應(yīng)用中,很難準(zhǔn)確確定地震子波,所以傳統(tǒng)的ISTA在反褶積應(yīng)用過程中因?yàn)榈卣鹱硬ǖ膯栴}而存在較大的局限性。

1.2 類RNN的改進(jìn)ISTA原理

RNN是一種節(jié)點(diǎn)定向連接成環(huán)的人工神經(jīng)網(wǎng)絡(luò),其本質(zhì)特征是在處理單元之間既有內(nèi)部的反饋連接又有前饋連接,可以在網(wǎng)絡(luò)運(yùn)行中對內(nèi)部參數(shù)不斷進(jìn)行優(yōu)化,這是一類適用于處理序列數(shù)據(jù)的神經(jīng)網(wǎng)絡(luò),它與基礎(chǔ)的神經(jīng)網(wǎng)絡(luò)最大的不同就是在層之間的神經(jīng)元也建立了權(quán)值連接[14-20]。

圖2為RNN基本結(jié)構(gòu)及展開形式。左側(cè)為RNN基本結(jié)構(gòu),垂向黑色粗箭頭表示“循環(huán)”,這意味著RNN的循環(huán)體現(xiàn)在隱藏層;右側(cè)為RNN展開結(jié)構(gòu),隱藏層各神經(jīng)元之間通過不同權(quán)值(見小箭頭上的權(quán)值符號(hào))連接,即上一步處理的隱藏層(Nt-1)將會(huì)影響到下一步處理的隱藏層(Nt)。RNN結(jié)構(gòu)包含了RNN正向計(jì)算(黑色細(xì)箭頭方向)和反向誤差修正(紅色細(xì)箭頭方向)兩個(gè)核心內(nèi)容,其中M和U為計(jì)算權(quán)值,V為激活函數(shù)。本節(jié)提出的類RNN的改進(jìn)ISTA同樣包含了正向計(jì)算和反向誤差修正兩個(gè)重要環(huán)節(jié)。

對比圖1和圖2可以發(fā)現(xiàn),ISTA與RNN在結(jié)構(gòu)上有著相似之處。將ISTA中的每一次迭代看作一個(gè)時(shí)間層,輸入向量中的每一個(gè)值看作是一個(gè)神經(jīng)元,軟閾值函數(shù)hθ等價(jià)于激活函數(shù)V,那么ISTA便可以看作一個(gè)類似的RNN網(wǎng)絡(luò),如圖3所示。

圖2 RNN基本結(jié)構(gòu)及展開形式

圖3 類RNN的改進(jìn)ISTA

由圖3可知,矩陣S為循環(huán)層參數(shù),即為時(shí)間方向參數(shù),We為各層神經(jīng)元之間的參數(shù),即為層方向參數(shù)。

2 改進(jìn)算法正反向計(jì)算方法

2.1 正向算法

由圖2可知,RNN的正向計(jì)算可表示為:

(4)

式中:Et為中間變量。

公式(3)可以改寫為如下形式:

(5)

其中,公式(5)中的y相當(dāng)于公式(4)中的輸入xt,S相當(dāng)于隱藏層權(quán)值M,xk相當(dāng)于公式(4)中的隱藏層參數(shù)Nt,hθ相當(dāng)于激活函數(shù)。

由此可以給出改進(jìn)ISTA的正向計(jì)算流程(圖4),其計(jì)算步驟與常規(guī)的ISTA一致。

2.2 反向誤差修正算法

圖4 正向算法流程

圖5中,xr為訓(xùn)練時(shí)給定的反射系數(shù),xk為訓(xùn)練過程中計(jì)算獲得的反射系數(shù),h′為軟閾值求導(dǎo),Hk,Hk+1為梯度下降過程中的控制函數(shù)。由圖5可知,需要先通過迭代逐步反向求解各層誤差,然后求和,得到關(guān)于初始輸入的誤差δB、δθ、δWe,進(jìn)一步利用圖5流程最后一步的梯度下降方法更新這3個(gè)參數(shù)。

圖5 改進(jìn)算法的反向誤差修正流程

結(jié)合正向算法與反向算法,便可以根據(jù)訓(xùn)練數(shù)據(jù)多次迭代修改與地震子波有關(guān)的參數(shù),得到更好的收斂效果,更加有利于后續(xù)的反褶積計(jì)算。

3 理論模型數(shù)據(jù)試算

首先合成了一個(gè)不含噪聲的地震單道數(shù)據(jù),使用50Hz零相位Ricker子波與第1和第2個(gè)反射系數(shù)進(jìn)行褶積,40Hz零相位Ricker子波與第3和第4個(gè)反射系數(shù)進(jìn)行褶積,30Hz零相位Ricker子波與第5和第6個(gè)反射系數(shù)進(jìn)行褶積,最終得到如圖6所示的理論地震道,該數(shù)據(jù)包含512個(gè)采樣點(diǎn),6個(gè)反射位置。

圖6 理論模型a 地震子波; b 反射系數(shù); c 合成地震記錄

利用圖6中的模型數(shù)據(jù),令初始輸入的子波矩陣為零相位40Hz的Ricker子波,結(jié)合前文所述的正向算法與反向算法,對初始子波進(jìn)行修正,經(jīng)過本文算法優(yōu)化后的子波如圖7所示,從圖7可以看出,雖然輸入僅有40Hz的子波,但是經(jīng)過算法自我修正后,可以得到與圖6a多個(gè)子波接近的結(jié)果,且經(jīng)過算法修正后獲得的子波集合與期望子波集合基本吻合。使用優(yōu)化后的子波矩陣進(jìn)行反褶積處理,得到如圖8所示的結(jié)果。給定同樣的初始子波(40Hz的Ricker子波),進(jìn)行傳統(tǒng)ISTA反褶積方法處理,結(jié)果如圖8中藍(lán)線所示。

由圖8可見,由于傳統(tǒng)ISTA不具備子波修正能力,當(dāng)采用40Hz子波作為輸入子波時(shí),算法收斂后計(jì)算結(jié)果與期望輸出存在明顯差異,反演結(jié)果出現(xiàn)了許多虛假信息。而類RNN的改進(jìn)ISTA在計(jì)算過程中不斷優(yōu)化子波集,最終在不同的位置獲得不同的最優(yōu)子波。在此基礎(chǔ)上計(jì)算的結(jié)果明顯優(yōu)于傳統(tǒng)ISTA。這說明類RNN的改進(jìn)ISTA能夠通過反向修正子波集合,消除初始輸入子波不準(zhǔn)確對反褶積帶來的影響。

為了說明本文方法的抗噪性,在圖6c合成地震記錄的基礎(chǔ)上加入隨機(jī)噪聲,得到一個(gè)如圖9所示的低信噪比地震記錄,傳統(tǒng)ISTA與本文方法的反褶積結(jié)果如圖10所示。

圖7 原始子波與優(yōu)化后的子波

圖8 傳統(tǒng)ISTA和改進(jìn)ISTA反褶積處理的結(jié)果

圖9 含噪聲地震記錄(信噪比=-4.55dB)

圖10 含噪聲地震記錄反褶積結(jié)果

由圖10可以看出,傳統(tǒng)的ISTA在地震資料信噪比較低時(shí),難以獲得準(zhǔn)確的反褶積結(jié)果,而采用改進(jìn)ISTA可獲得較為可靠的反褶積結(jié)果,說明該方法具有良好的抗噪能力。

4 實(shí)際地震資料應(yīng)用

圖11為吐哈盆地雁木西地區(qū)T3井區(qū)Line 299線主要目的層段局部地震剖面。分別采用傳統(tǒng)ISTA和改進(jìn)ISTA對其進(jìn)行反褶積處理,結(jié)果如圖12所示。

對比圖11和圖12可以看出,采用傳統(tǒng)ISTA進(jìn)行反褶積處理后的剖面整體同相軸變多,變細(xì),且以井旁地震道分辨率提高最為明顯,但是距離井位置越遠(yuǎn)處理效果越不明顯(圖12a)。相比之下,采用改進(jìn)ISTA的反褶積剖面中,不論是井旁道還是距離井位置較遠(yuǎn)的地震道,分辨率均得到了顯著提升,且與測井合成記錄的吻合程度大大提高,說明采用改進(jìn)ISTA的反褶積處理結(jié)果可靠。從圖11和圖12中井旁紅色方框中局部放大的地震剖面(圖13)可以看出,采用傳統(tǒng)ISTA處理后,原剖面中部分無法分辨的同相軸得到了有效分離,但是分辨率仍然偏低,而采用改進(jìn)ISTA處理后的剖面,不僅使原來無法分離的同相軸得到了有效分離,而且同相軸連續(xù)性增強(qiáng),與井資料對比良好,分辨率明顯提高。

圖11 雁木西地區(qū)T3井Line 299線主要目的層段局部地震剖面

圖12 采用傳統(tǒng)ISTA(a)和改進(jìn)ISTA(b)進(jìn)行反褶積處理的結(jié)果(Line299)

圖13 采用不同方法反褶積結(jié)果剖面細(xì)節(jié)對比(Line299)a 原始剖面;b 傳統(tǒng)ISTA;c 改進(jìn)ISTA

分別對圖13中的3個(gè)剖面進(jìn)行頻譜分析,結(jié)果如圖14所示。顯然,基于改進(jìn)ISTA的反褶積處理結(jié)果有效頻帶寬度從原始剖面的10～40Hz拓展為5～50Hz左右,明顯優(yōu)于傳統(tǒng)ISTA。

圖14 頻譜對比

5 結(jié)論

1) 與傳統(tǒng)的基于ISTA的稀疏脈沖反褶積方法相比,采用類RNN的改進(jìn)ISTA實(shí)施反褶積處理,可以獲得更加穩(wěn)定、可靠的反褶積處理結(jié)果,且方法的抗噪性更佳。

2) 基于類RNN的改進(jìn)ISTA的稀疏脈沖反褶積方法可極大提高地震資料分辨率。實(shí)際資料的處理結(jié)果表明,采用該方法可使原始地震資料的有效頻帶拓展約1.5倍,而且除拓展高頻外,還可以極大地保護(hù)低頻有效信息,處理結(jié)果更有利于進(jìn)一步儲(chǔ)層精細(xì)描述或反演處理。