一種類RNN的改進ISTA稀疏脈沖反褶積

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

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

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

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

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

1 類RNN的改進ISTA原理

1.1 ISTA

稀疏脈沖反褶積可以表達為:

(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的準確確定,但是在實際應(yīng)用中,很難準確確定地震子波,所以傳統(tǒng)的ISTA在反褶積應(yīng)用過程中因為地震子波的問題而存在較大的局限性。

1.2 類RNN的改進ISTA原理

RNN是一種節(jié)點定向連接成環(huán)的人工神經(jīng)網(wǎng)絡(luò),其本質(zhì)特征是在處理單元之間既有內(nèi)部的反饋連接又有前饋連接,可以在網(wǎng)絡(luò)運行中對內(nèi)部參數(shù)不斷進行優(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)值符號)連接,即上一步處理的隱藏層(Nt-1)將會影響到下一步處理的隱藏層(Nt)。RNN結(jié)構(gòu)包含了RNN正向計算(黑色細箭頭方向)和反向誤差修正(紅色細箭頭方向)兩個核心內(nèi)容,其中M和U為計算權(quán)值,V為激活函數(shù)。本節(jié)提出的類RNN的改進ISTA同樣包含了正向計算和反向誤差修正兩個重要環(huán)節(jié)。

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

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

圖3 類RNN的改進ISTA

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

2 改進算法正反向計算方法

2.1 正向算法

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

(4)

式中:Et為中間變量。

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

(5)

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

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

2.2 反向誤差修正算法

圖4 正向算法流程

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

圖5 改進算法的反向誤差修正流程

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

圖14 頻譜對比

5 結(jié)論

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

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