李婷婷, 段中鈺

(北京信息科技大學(xué) 信息與通信工程學(xué)院，北京 100101)

0 引言

地震數(shù)據(jù)重建是地震數(shù)據(jù)處理中的基礎(chǔ)步驟，通過重建可以增強(qiáng)地震數(shù)據(jù)中的有效信號，抑制不利于后續(xù)數(shù)據(jù)處理和解釋的成分，從而提高地震勘探的準(zhǔn)確性，降低油氣開發(fā)的成本。

近幾年來發(fā)展起來的壓縮感知理論[1-3]為(Compressed Sensing，CS)為地震數(shù)據(jù)重建提供了新的思路，稀疏信號可以在低于奈奎斯特采樣頻率的條件下恢復(fù)出原始數(shù)據(jù)，數(shù)據(jù)的稀疏表示是壓縮感知地震數(shù)據(jù)重建的前提和基礎(chǔ)。地震數(shù)據(jù)的稀疏表示主要有單一變換域的稀疏表示和字典學(xué)習(xí)的稀疏表示。單一變換域的稀疏表示主要是數(shù)據(jù)通過固定正交基進(jìn)行變換(Curvelet變換[4]、Seislet變換[5]、Shearlet變換[6]、小波變換[7]等)，這種方法雖然構(gòu)造簡單, 計算復(fù)雜度低, 但是變換基是固定的, 不能根據(jù)地震數(shù)據(jù)復(fù)雜的細(xì)節(jié)結(jié)構(gòu)自適應(yīng)的調(diào)整。字典學(xué)習(xí)解決了單一變換無法對復(fù)雜地震數(shù)據(jù)進(jìn)行理想稀疏表達(dá)的問題。字典學(xué)習(xí)即用字典中若干個原子來表示信號，每組原子可以獨(dú)立的張成完備的空間，在超完備空間里，信號不是唯一的，這樣的優(yōu)點(diǎn)在于可以通過學(xué)習(xí)訓(xùn)練的方式自適應(yīng)地尋找最佳的基底原子組合，實(shí)現(xiàn)信號的最稀疏表示。學(xué)習(xí)型字典能夠通過樣本地震數(shù)據(jù)不斷地學(xué)習(xí)訓(xùn)練得到含有地震數(shù)據(jù)信息的自適應(yīng)字典，從而刻畫地震數(shù)據(jù)的本質(zhì)特征，得到地震數(shù)據(jù)的最稀疏表示。常用于壓縮感知地震數(shù)據(jù)重建的超完備字典學(xué)習(xí)算法是K-SVD(K-singular value decomposition)[8-9]。周亞同等[10]提出壓縮感知框架下K-SVD字典學(xué)習(xí)對地震數(shù)據(jù)進(jìn)行稀疏表示，再結(jié)合正則化正交匹配追蹤算法有效地實(shí)現(xiàn)了缺失地震數(shù)據(jù)的重建。但是K-SVD字典學(xué)習(xí)算法存在無法得到全局最優(yōu)解并且不能保證收斂的問題，會影響地震數(shù)據(jù)重建精度，且在數(shù)據(jù)量較大的重建任務(wù)下字典訓(xùn)練時間長、迭代較慢。筆者針對以上問題，提出快速字典學(xué)習(xí)算法(Fast Dictionary Learning，F(xiàn)DL),①將稀疏表示目標(biāo)函數(shù)優(yōu)化求解分為兩個子優(yōu)化問題進(jìn)行循環(huán)求解；②在稀疏編碼階段將稀疏約束上限與字典相干性結(jié)合，自適應(yīng)獲取稀疏上限。在此基礎(chǔ)上提出基于壓縮感知理論的FDL重建地震數(shù)據(jù)的模型,使用模擬和大慶油田實(shí)際地震數(shù)據(jù)對本方法進(jìn)行了驗(yàn)證。

1 方法原理

1.1 壓縮感知理論

壓縮感知理論指出對于稀疏信號xi∈Rn(i=1,2,…，N)，可以利用一個與變換基Ψ不相關(guān)的觀測矩陣Φ∈Rm×n(m<

已知原始信號在觀測矩陣的觀測值為：

y=Φx

(1)

用觀測矩陣不相關(guān)的變換基Ψ對原始信號進(jìn)行稀疏表示，得到稀疏矩陣Θ：

x=ΨΘ

(2)

y=ΦΨΘ

(3)

(4)

(5)

1.2 快速字典學(xué)習(xí)算法

字典學(xué)習(xí)稀疏表示的主要思想是利用字典矩陣D稀疏線性表示原始樣本X。用字典學(xué)習(xí)對地震數(shù)據(jù)進(jìn)行稀疏表示可以描述為：

X=DΘ

(6)

其中：X={x1,x2, …,xN}∈Rn×N表示包含N個子地震數(shù)據(jù)塊的輸入矩陣；D∈Rk×K表示含有K個原子dk的字典矩陣；Θ∈RK×N表示與字典D對應(yīng)的稀疏系數(shù)矩陣。

稀疏矩陣和字典矩陣通過求解優(yōu)化問題得到：

(7)

其中，T0為稀疏表示系數(shù)中非零分量數(shù)目的上限。

快速字典學(xué)習(xí)算法稀疏表示目標(biāo)函數(shù)定義如下：

(8)

其中：Θi為矩陣Θ第i列;Dj為矩陣D的第j列。字典矩陣D和稀疏矩陣Θ是非凸的，但在稀疏編碼階段，固定字典D，則對于Θ為凸的；在字典更新階段，固定稀疏矩陣Θ，則對于D為凸的。K-SVD算法稀疏編碼階段通常使用正交匹配追蹤算法，字典更新階段使用奇異值分解，但正交匹配追蹤算法收斂速度較慢，且無法得到全局最優(yōu)解。與K-SVD不同的是，F(xiàn)DL將目標(biāo)函數(shù)分解為兩個子優(yōu)化問題，對于每個子優(yōu)化問題，用塊坐標(biāo)下降法直接求出目標(biāo)函數(shù)的全局最優(yōu)解。因此，F(xiàn)DL算法的稀疏編碼階段優(yōu)化問題為：

(9)

字典更新階段優(yōu)化問題為：

(10)

每次更新Θ的一行或是D的一列，更新速度與字典大小相關(guān)聯(lián)。

為了使字典能夠?qū)Φ卣饠?shù)據(jù)表示更稀疏，將稀疏約束上限與字典的相干性關(guān)聯(lián)，從而獲得自適應(yīng)稀疏約束上限。

引入定理1：

給定字典D∈Rk×K(K>n)，其相干性為μ(D)，假設(shè)x=Dα有稀疏解α，其稀疏度S若滿足:

(11)

則能從字典D中找出最佳的S項(xiàng)原子的線性組合。其中字典的相干性μ(D)∈[0,1]，μ(D)=max|[di,dj]|,i≠j其描述了字典中原子間的最大相似程度，當(dāng)μ(D)值很大時，字典原子間相似程度很強(qiáng)，反之很弱。

定義Lp為每次迭代過程中的稀疏約束上界:

(12)

由定理1可知稀疏約束上限Lp能夠保證稀疏信號被精確重建。使用Lp代替式(9)中的T0，則稀疏編碼階段即為求解如下的優(yōu)化問題:

(13)

2 技術(shù)流程

2.1 基于K-SVD的壓縮感知地震數(shù)據(jù)重建

K-SVD算法地震數(shù)據(jù)重建流程圖如圖1所示。

圖1 K-SVD地震數(shù)據(jù)重建流程圖

K-SVD算法地震數(shù)據(jù)重建：

1)隨機(jī)挑選個K原子形成初始字典，初始化稀疏矩陣為零。

2)稀疏編碼階段。固定字典D，使用正交匹配追蹤算法更新稀疏矩陣Θ。

3)字典更新階段。固定稀疏矩陣Θ，奇異值分解逐列更新字典D。

2.2 基于FDL的壓縮感知地震數(shù)據(jù)重建

FDL算法地震數(shù)據(jù)重建流程圖如圖2所示。

圖2 FDL地震數(shù)據(jù)重建流程圖

FDL算法地震數(shù)據(jù)重建：

1)根據(jù)缺失地震數(shù)據(jù)隨機(jī)挑選K個原子形成初始字典D0，設(shè)置迭代次數(shù)i=0。

2)根據(jù)初始字典D0計算式(12)得到初始化稀疏約束上限L0。

3)稀疏編碼階段。固定字典Di和稀疏約束上限Li，使用塊坐標(biāo)下降法計算式(9)逐行更新得到最優(yōu)的稀疏矩陣Θi+1。

4)字典更新階段。固定稀疏矩陣Θi+1，使用塊坐標(biāo)下降法計算式(10)逐列更新得到本次迭代的最優(yōu)字典矩陣Di+1。

5)用步驟2)、步驟3)、步驟4)得到稀疏約束上限、稀疏矩陣和字典矩陣更新目標(biāo)函數(shù)即式(8)。

6)判斷是否達(dá)到收斂條件，通過計算(f(i)-f(i+1))/f(i)>1×10-5進(jìn)行判斷，若是則返回步驟2)，若否則停止迭代，輸出字典矩陣D和稀疏矩陣Θ。

字典大小對重建精度有影響，但是只能根據(jù)經(jīng)驗(yàn)選取，首先根據(jù)地震數(shù)據(jù)隨機(jī)選擇K個原子，隨后在將字典大小改為K±j,j=1、2…，最終選擇重建信噪比最高，重建性能最好的字典。

3 仿真驗(yàn)證

采用信噪比作為地震數(shù)據(jù)重建性能的評判標(biāo)準(zhǔn)為式(14)。

(14)

3.1 模擬地震數(shù)據(jù)仿真驗(yàn)證

模擬地震數(shù)據(jù)共38道，時間采樣率為0.004 s，如圖3所示，其進(jìn)行隨機(jī)采樣，采樣率為0.3，采樣后數(shù)據(jù)如圖4所示。分別用傳統(tǒng)稀疏表示方法Curvelet、經(jīng)典K-SVD字典學(xué)習(xí)算法以及FDL算法,對缺失地震數(shù)據(jù)稀疏表示然后進(jìn)行重建，重建結(jié)果如圖5～圖7所示。由圖5～圖7可知，Curvelet、K-SVD字典學(xué)習(xí)和FDL算法在壓縮感知的框架下都實(shí)現(xiàn)了缺失地震數(shù)據(jù)的重建，但是從圖5可以看出，部分地震數(shù)據(jù)沒有被完全恢復(fù)，地震波附近出現(xiàn)抖動，且原來未缺失數(shù)據(jù)受到影響(圖中箭頭所示)；從圖6可以看出,K-SVD字典學(xué)習(xí)重建地震數(shù)據(jù)結(jié)果較Curvelet更好，基本實(shí)現(xiàn)了缺失數(shù)據(jù)的重建，但是存在重建振幅與原始數(shù)據(jù)振幅不一致的問題(圖中箭頭所示)；而從圖7可以看出,FDL算法重建后數(shù)據(jù)振幅連續(xù)性良好，重建精度較高，很好地實(shí)現(xiàn)了缺失數(shù)據(jù)的重建。

圖3 原始地震數(shù)據(jù)

圖4 缺失地震數(shù)據(jù)

圖5 Curvelet地震數(shù)據(jù)重建結(jié)果

圖6 K-SVD地震數(shù)據(jù)重建結(jié)果

為了更清楚地對比三種算法，將三種算法的重建結(jié)果與原始地震數(shù)據(jù)進(jìn)行對比，選取缺失較多的25至38道數(shù)據(jù)，對比結(jié)果如圖8所示。對比圖8(a)、圖8(b)、圖8(c)可以看出，Curvelet算法重建振幅與原始數(shù)據(jù)振幅差異較大，K-SVD算法重建振幅差異較小，而FDL算法重建振幅與原始數(shù)據(jù)振幅一致。

圖8 三種算法重建結(jié)果與原始地震數(shù)據(jù)對比

三種算法的重建信噪比與重建時間見表1。結(jié)合表1和圖5～圖7 可以看出，在相同采樣率且K-SVD字典學(xué)習(xí)和FDL字典大小相同的情況下，傳統(tǒng)Curvelet變換的方法重建信噪比最小，部分缺失數(shù)據(jù)未得到恢復(fù)，對缺失地震數(shù)據(jù)的重建效果最差，且重建耗時最長；K-SVD字典學(xué)習(xí)重建信噪比較Curvelet高，重建耗時也相對增加；而FDL地震數(shù)據(jù)重建信噪比最高，對缺失的地震數(shù)據(jù)實(shí)現(xiàn)了恢復(fù)，重建效果最好，且重建耗時最短。

表1 重建方法信噪比及重建時間

三種算法的重建殘差圖對比如圖9所示。三種算法重建模型大部分都能較好的符合原始數(shù)據(jù)，但是相比來說Curvelet和K-SVD算法存在較多異常點(diǎn)，而FDL異常點(diǎn)較少，F(xiàn)DL算法重建性能相比較好。

圖9 三種重建方法殘差圖

表2為字典數(shù)量差異及其重建性能。由表2可知，字典大小的選擇對重建性能有一定的影響，字典越小，字典構(gòu)造時間越短，重建信噪比越??；字典越大，字典構(gòu)造耗時變長，信噪比越大。對于K-SVD字典學(xué)習(xí)和FDL算法，在構(gòu)造相同大小的字典的情況下，F(xiàn)DL算法的重建信噪比更大，字典構(gòu)造時間更短，性能更好。若要得到相同的重建信噪比，K-SVD字典學(xué)習(xí)相較FDL需更長的時間構(gòu)造字典，從而重建所需時間更長。

表2 K-SVD與FDL性能對比

3.2 實(shí)際地震數(shù)據(jù)仿真驗(yàn)證

模擬地震數(shù)據(jù)驗(yàn)證了本文方法的有效性，因此將本文方法應(yīng)用到大慶油田的三維地震數(shù)據(jù)中，選取其中的200道數(shù)據(jù)，200個采樣點(diǎn)，時間采樣率為1 ms，如圖10所示，對該數(shù)據(jù)進(jìn)行隨機(jī)采樣，缺失數(shù)據(jù)如圖11所示。

圖10 原始三維地震數(shù)據(jù)

分別用Curvelet變換、K-SVD算法以及FDL算法對隨機(jī)采樣后的地震數(shù)據(jù)進(jìn)行重建(K-SVD字典學(xué)習(xí)和FDL的字典大小相同)，重建結(jié)果如圖12～圖14所示。對比圖12～圖14，可以看出Curvelet重建地震數(shù)據(jù)效果一般，同相軸模糊不清，而FDL很好地實(shí)現(xiàn)了地震數(shù)據(jù)的重建，同相軸清晰，振幅連續(xù)。為了能更清晰地觀察重建效果，將隨機(jī)缺失較多的50道至60道數(shù)據(jù)與原始地震數(shù)據(jù)進(jìn)行對比(圖15)，對比可以看出，Curvelet變換重建地震數(shù)據(jù)振幅與原始地震數(shù)據(jù)振幅相差較大，K-SVD重建振幅較Curvelet稍好，而FDL重建振幅與原始地震數(shù)據(jù)基本一致。三種算法的重建信噪比與重建時間如表3所示，由表3可知，F(xiàn)DL壓縮感知地震數(shù)據(jù)重建耗時最短。

表3 重建方法信噪比及重建時間

圖12 Curvelet重建地震數(shù)據(jù)結(jié)果

圖13 K-SVD重建地震數(shù)據(jù)結(jié)果

圖14 FDL重建地震數(shù)據(jù)結(jié)果

圖15 三種算法重建結(jié)果與原始地震數(shù)據(jù)對比

4 結(jié)論

筆者提出的基于快速字典學(xué)習(xí)的壓縮感知地震數(shù)據(jù)重建，通過快速字典學(xué)習(xí)對地震數(shù)據(jù)塊進(jìn)行稀疏表示實(shí)現(xiàn)地震數(shù)據(jù)的重建。模擬和大慶油田地震數(shù)據(jù)驗(yàn)證可以得到以下結(jié)論：利用FDL比經(jīng)典K-SVD字典學(xué)習(xí)能夠更好的重建缺失地震數(shù)據(jù)，重建精度更高，重建耗時更短，對地震數(shù)據(jù)中復(fù)雜波形區(qū)域能夠有效的重建。另外，對于含噪聲的地震數(shù)據(jù)如何使用本文提出的方法進(jìn)行去噪，是下一步研究的內(nèi)容。