楊 帆,王長鵬,張春霞,張講社,熊 登

1.長安大學(xué)理學(xué)院,西安 710064

2.西安交通大學(xué)數(shù)學(xué)與統(tǒng)計學(xué)院,西安 710049

3.東方地球物理公司物探技術(shù)研究中心,河北 涿州 072751

0 引言

地震數(shù)據(jù)重建是地震研究中的一個重要課題,其目的是重建規(guī)則采樣或不規(guī)則采樣的缺失痕跡[1]。重建的質(zhì)量直接影響到后續(xù)的處理步驟,如全波形反演[2-3]、偏移矢量傾斜[4]等。人們提出許多不同的方法用于地震數(shù)據(jù)重建,包括促進(jìn)稀疏性的L1范數(shù)最小化方法、基于預(yù)測濾波器的方法、基于壓縮感知的方法等,其中:促進(jìn)稀疏性的L1范數(shù)最小化方法要求變換域中的數(shù)據(jù)通過變換[5-9]稀疏地表示,其局限性是過于注重固定變換中的稀疏表示,所以缺乏對復(fù)雜數(shù)據(jù)的適應(yīng)能力;基于預(yù)測濾波器的方法[10]利用了信號在f-x域(頻率-空間域)和f-k域(頻率-波數(shù)域)的線性可預(yù)測性,但是需要等距地震信號數(shù)據(jù);基于壓縮感知的方法[11-12]利用地震數(shù)據(jù)的稀疏特征重建地震數(shù)據(jù)。近年流行的神經(jīng)網(wǎng)絡(luò)方法也被應(yīng)用于地震數(shù)據(jù)重建[13-17],然而為了獲得精確的結(jié)果,它需要大量的數(shù)據(jù)和時間訓(xùn)練網(wǎng)絡(luò),以學(xué)習(xí)輸入和輸出之間的非線性映射。

低秩矩陣補全方法對數(shù)據(jù)適應(yīng)性強,重建效率高,不需要大量的樣本來提取數(shù)據(jù)的特征,且對硬件要求友好,已在地震數(shù)據(jù)重建方面得到廣泛的應(yīng)用。由于秩函數(shù)的非凸性和不連續(xù)性,該模型一般是非確定性多項式時間難題(non-deterministic polynomial-time hard, NP-hard),而核范數(shù)作為秩函數(shù)的凸松弛,在矩陣補全問題中取代了秩函數(shù)。這個優(yōu)化問題可以用奇異值閾值(singular value threshold, SVT)[18]法解決,該算法對特定矩陣的奇異值使用迭代軟閾值操作,并使用交替迭代方法求解優(yōu)化問題,可以有效重建地震數(shù)據(jù)。低秩矩陣補全問題的另一個模型是核范數(shù)正則化線性最小二乘問題,加速近端梯度(accelerated proximal gradient, APG)[19]算法是基于壓縮傳感的快速迭代收縮閾值算法,可以有效求解該問題。Yang等[20]首次將核范數(shù)與地震數(shù)據(jù)重建聯(lián)系起來,并通過紋理塊算子[21]對地震數(shù)據(jù)矩陣進(jìn)行預(yù)處理。Zhang等[22]提出基于最大化最小化框架的非凸對數(shù)和函數(shù)算法(nonconvex log-sum function-based majorization-minimization framework algorithm, LSMM),該算法提出對數(shù)和函數(shù)代替核范數(shù)來逼近矩陣的秩函數(shù)。由于對數(shù)和函數(shù)的非凸性,采用majorization-minimization算法進(jìn)行求解。Zhang等[23]利用截斷核范數(shù)正則化重建了地震數(shù)據(jù),并用交替方向乘子法(alternating direction method of multipliers,ADMM)進(jìn)行求解。Wang等[24]首次將加權(quán)核范數(shù)最小化(weight nuclear norm minimization, WNNM)模型應(yīng)用于地震數(shù)據(jù)重建,將L1范數(shù)的迭代加權(quán)原理擴(kuò)展到核范數(shù)最小化問題中[25-26],取得了較好的效果。然而,權(quán)重向量的增加導(dǎo)致一些奇異值在迭代過程中丟失,重建精度也會隨之下降。

本文提出基于聯(lián)合加速近端梯度和對數(shù)加權(quán)核范數(shù)最小化(joint accelerated proximal gradient and log-weighted nuclear norm minimization, APG-LWNNM)的地震數(shù)據(jù)重建方法。首先通過紋理塊算子對原始地震數(shù)據(jù)進(jìn)行低秩預(yù)處理;然后使用APG算法對低秩地震數(shù)據(jù)進(jìn)行初步重建,提出APG-LWNNM算法,使該算法解決優(yōu)化問題并重建缺失數(shù)據(jù),以在迭代過程中減少矩陣奇異值損失,提高重建精度;最后用合成地震數(shù)據(jù)和真實地震數(shù)據(jù)對比實驗證明APG-LWNNM算法的有效性。

1 紋理塊預(yù)處理

對于給定的原始數(shù)據(jù)矩陣X∈Rm×n和紋理塊尺寸p,

可以將X分解為mn/p2個子塊,標(biāo)記為Bj∈Rp×p(j∈[1,2,…,mn/p2])。對于分塊矩陣:

以B1為例,可以重新組合為

則紋理塊矩陣PT(X)可以定義為

PT(X)=[b1,b2,…,bmn/p2]∈Rp2×(mn/p2)。

(1)

式中,PT：Rm×n→Rp2×(mn/p2)表示紋理塊算子。相應(yīng)地,PT′:Rp2×(mn/p2)→Rm×n表示紋理塊逆算子。

紋理塊矩陣PT(X)具有低秩結(jié)構(gòu),并且缺失的數(shù)據(jù)會隨機分散在整個矩陣中,避免了全零列的出現(xiàn)。這與低秩矩陣補全理論[18-19]一致。

2 低秩矩陣補全

低秩矩陣補全問題的模型為

(2)

式中:Μ∈Rm×n,為要恢復(fù)的矩陣;rank(X)為矩陣X的秩。低秩矩陣補全問題的基本理論是通過問題(2),求得最優(yōu)的低秩矩陣X*,使得X*在已知項上的元素與M足夠接近,在未知項上的元素即為恢復(fù)的缺失值[27]。

2.1 核范數(shù)最小化

由于秩函數(shù)的非凸性和不連續(xù)性,直接求解問題(2)是NP-hard的。核范數(shù)是矩陣單位球上rank(X)的凸近似,常用來近似rank(X),則問題(2)可轉(zhuǎn)化為

(3)

由于核范數(shù)不能對rank(X)產(chǎn)生很好的近似,Gaiffas等[25]將迭代加權(quán)原理應(yīng)用到了核范數(shù)中,提出了加權(quán)核范數(shù)來逼近rank(X):

(4)

2.2 核范數(shù)正則化線性最小二乘

低秩矩陣補全問題的另一個模型是核范數(shù)正則化線性最小二乘問題。該正則化問題是無約束非光滑凸優(yōu)化問題的一種特殊情況,其目標(biāo)函數(shù)是核范數(shù)和具有Lipschitz連續(xù)梯度的凸光滑函數(shù)之和,即

(5)

式中:λ>0,為超參數(shù);‖·‖2為矩陣的譜范數(shù)。APG算法[19]可以有效地求解核范數(shù)正則化線性最小二乘問題,它結(jié)合了梯度下降法和近端算子,在每次迭代中,通過計算當(dāng)前矩陣的梯度和近端算子來更新矩陣的近似解,其中,近端算子可以有效地對矩陣進(jìn)行低秩約束,從而得到更加合理的解。

3 APG-LWNNM算法

3.1 LWNNM算法

LWNNM算法目標(biāo)函數(shù)為

(6)

將X的奇異值分解記為X=UΣVT,其中:

U=[u1,u2,…,ur];

V=[v1,v2…,vr];

Σ=diag(σ1,σ2,…,σr)。

由

(7)

(式中,〈·,·〉為矩陣的內(nèi)積運算),此時固定M,只考慮X,結(jié)合式(7),則式(6)的目標(biāo)函數(shù)可轉(zhuǎn)化為

(8)

將M的奇異值分解記為M=U′Σ′V′T,其中:

式(8)可轉(zhuǎn)換為

(9)

(10)

固定j,則目標(biāo)函數(shù)為

(11)

此時式(11)易求解。當(dāng)

時,目標(biāo)函數(shù)取得最小值。

綜上所述,式(6)可通過對數(shù)加權(quán)軟閾值算子求得唯一解:

(12)

(13)

(14)

其中,

(15)

由此可知,對數(shù)權(quán)重向量的引入減少了矩陣奇異值的損失,保留了更多的信息量。

3.2 APG-LWNNM算法

本文提出的APG-LWNNM算法基于低秩矩陣補全理論,主要分為四階段進(jìn)行。第一階段用紋理塊算子對輸入數(shù)據(jù)進(jìn)行低秩預(yù)處理;第二階段用APG算法求得初步重建結(jié)果;第三階段用本文所提LWNNM算法解決優(yōu)化問題并重建地震數(shù)據(jù);最后用紋理塊逆算子進(jìn)行逆處理。

具體流程如圖1所示。

圖1 APG-LWNNM算法流程圖

aPG-LWNNM算法如下。

輸入:缺失地震數(shù)據(jù)M,k,τ,λ,tol,λtar;

計算:1)紋理塊矩陣PT(Μ);

2)通過APG算法得到X;

賦值:ω=σ(X);

Whileλ>λtar

θ=∞;

Whileθ>tol

Xold=X;

PΩ(M)];

θ=‖X-Xold‖/‖Xold‖;

end

賦值:ω=σ(X),λ=k×λ;

End

輸出:紋理塊逆矩陣PT′(X)。

其中:k=0.02;tol=10-3;λ=1;λtar=10-4·‖PT(M)‖;τ<10-6,否則容易求得局部最優(yōu)解。具體實驗參數(shù)依輸入地震數(shù)據(jù)而定。

4 實驗與分析

本文以信噪比和重建誤差為評價指標(biāo)。信噪比的概念來源于電子通信,代表信號和噪聲的有效功率之比。在插值方法中,將原始數(shù)據(jù)和重建結(jié)果之間的差距視為噪聲,噪聲越大說明重建精度越低。信噪比本質(zhì)上為相對誤差,重建誤差本質(zhì)上為絕對誤差。信噪比表示為

(16)

重建誤差表示為

e=‖X-X*‖2。

(17)

式中:X為原始地震數(shù)據(jù);X*為重建結(jié)果。信噪比數(shù)值越大、重建誤差數(shù)值越小,表明重建越度越高。

4.1 消融實驗

選取一組地震數(shù)據(jù)共300道,300個采樣點,采樣間隔為0.004 s,并對該數(shù)據(jù)進(jìn)行隨機缺失,缺失率為20%～80%。分別使用APG-LWNNM算法和LWNNM算法進(jìn)行數(shù)據(jù)重建,重建結(jié)果的信噪比和重建誤差見表1。

表1 APG-LWNNM算法與LWNNM算法的重建性能比較

由表1可見,在同組數(shù)據(jù)相同缺失率條件下,APG-LWNNM算法的重建結(jié)果信噪比更高,重建誤差更低,即重建精度更高。因此,使用APG算法求得初步重建結(jié)果是必要的。

4.2 實驗結(jié)果與分析

為評估APG-LWNNM算法的有效性,本節(jié)將在合成地震數(shù)據(jù)和真實地震數(shù)據(jù)上進(jìn)行實驗,并與現(xiàn)有主流的低秩地震數(shù)據(jù)重建算法SVT、WNNM以及LSMM進(jìn)行比較,按原算法進(jìn)行參數(shù)設(shè)置。

合成地震數(shù)據(jù)(圖2a)共260道,300個采樣點,采樣間隔為0.004 s。對隨機缺失40%地震道數(shù)據(jù)(圖2b),分別用SVT、WNNM、LSMM和本文所提APG-LWNNM等4種算法進(jìn)行重建?？梢钥闯?SVT算法(圖2c)對部分缺失數(shù)據(jù)未完全重建;WNNM算法(圖2d)和LSMM算法(圖2e)重建結(jié)果較好,但對數(shù)據(jù)細(xì)節(jié)的恢復(fù)不夠完備,仍有部分?jǐn)?shù)據(jù)無法恢復(fù);APG-LWNNM算法(圖2f)對該缺失數(shù)據(jù)能夠較好地恢復(fù),數(shù)據(jù)細(xì)節(jié)恢復(fù)更完備,生成的紋理更加豐富。

上述4種算法在該數(shù)據(jù)上重建結(jié)果的信噪比和重建誤差見表2。其中APG-LWNNM算法重建結(jié)果的信噪比最高,重建誤差最低。

表2 合成地震數(shù)據(jù)重建性能比較

為了更加直觀地對重建結(jié)果進(jìn)行比較,對第162道地震數(shù)據(jù)重建結(jié)果進(jìn)行單道比較。SVT、WNNM、LSMM和APG-LWNNM等4種算法重建結(jié)果單道比較結(jié)果(圖3a—d)以及相應(yīng)的殘差(圖3e—h)表明,APG-LWNNM算法的重建結(jié)果更好。

本文在Mobil Avo Viking Graben Line 12疊前真實數(shù)據(jù)集和Netherlands F3真實數(shù)據(jù)集上評估APG-LWNNM方法。

Mobil Avo Viking Graben Line 12疊前地震數(shù)據(jù)集是一個二維野外數(shù)據(jù)集,包含北海的1 001個炮點集。每個炮點集有120道,均勻分布,距離為25 m。每道地震數(shù)據(jù)有1 500個時間采樣點,采樣間隔為0.004 s。用于實驗的數(shù)據(jù)大小為120×300。

Netherlands F3地震數(shù)據(jù)集是一個位于荷蘭近海的小型海洋野外地震數(shù)據(jù)集,包含601個炮點集。每個炮點集有951道,每道地震數(shù)據(jù)有463個時間采樣點,采樣間隔為0.004 s。用于實驗的數(shù)據(jù)大小為200×180。

Mobil Avo Viking Graben Line 12和Netherlands F3地震數(shù)據(jù)集下載地址為https://wiki.seg.org/wiki/,所有數(shù)據(jù)均可為公眾使用。

首先對Mobil Avo Viking Graben Line 12疊前地震數(shù)據(jù)(圖4a)進(jìn)行30%地震道隨機缺失(圖4b),分別使用SVT、WNNM、LSMM以及APG-LWNNM等4種算法對缺失數(shù)據(jù)進(jìn)行重建,可以看出:SVT算法(圖4c)不能很好地重建缺失數(shù)據(jù),紋理模糊;WNNM算法(圖4d)和LSMM算法(圖4e)重建結(jié)果一般,紋理不夠豐富;而APG-LWNNM算法(圖4f)重建結(jié)果生成的紋理較豐富。特別是在矩形框內(nèi),APG-LWNNM算法重建結(jié)果更符合原始數(shù)據(jù)。

4種算法在Mobil Avo Viking Graben Line 12疊前地震數(shù)據(jù)上重建結(jié)果的信噪比和重建誤差見表3,其中本文所提APG-LWNNM算法的信噪比最高,重建誤差最低。

為了直觀地比較比較重建結(jié)果,接下來用評價指標(biāo)對單個地震道作進(jìn)一步分析。將SVT、WNNM、LSMM和APG-LWNNM等4種算法重建結(jié)果第2道進(jìn)行單道比較,由重建結(jié)果(圖5a—d)及相應(yīng)的殘差(圖5e—h)明顯看出,APG-LWNNM算法的重建結(jié)果更接近于原始地震數(shù)據(jù),與原始數(shù)據(jù)殘差最小,表明APG-LWNNM算法的重建精度優(yōu)于其他3種算法。

接著選取Netherlands F3地震數(shù)據(jù)集(圖6a)共200道,180個采樣點,采樣間隔為0.004 s。隨機缺失比例為60%(圖6b),分別使用SVT、WNNM、LSMM以及APG-LWNNM等4種算法對缺失數(shù)據(jù)(圖6b)進(jìn)行重建(圖6c—f), 其信噪比和重建誤差見表4,第59道地震數(shù)據(jù)單道比較及相應(yīng)的殘差見圖7。

表4 Netherlands F3地震數(shù)據(jù)重建性能比較

a. 原始真實地震數(shù)據(jù);b. 隨機缺失60%真實地震數(shù)據(jù);c. SVT算法重建結(jié)果;d. WNNM算法重建結(jié)果;e. LSMM算法重建結(jié)果;f. APG-LWNNM算法重建結(jié)果。

a. SVT算法重建結(jié)果;b. WNNM算法重建結(jié)果;c. LSMM算法重建結(jié)果;d. APG-LWNNM算法重建結(jié)果;e. SVT算法殘差;f. WNNM算法殘差;g. LSMM算法殘差;h. APG-LWNNM算法殘差。

結(jié)合表4、圖6和圖7可以分析得到:SVT算法重建結(jié)果仍有部分地震道殘缺,且單道圖與原始數(shù)據(jù)有明顯誤差,欠理想;WNNM算法和LSMM算法結(jié)果較好,但其信噪比較低,重建誤差較大;APG-LWNNM算法能夠較完整地重建該缺失數(shù)據(jù),且其信噪比最高,重建誤差最低。由第59道地震數(shù)據(jù)單道比較可知,APG-LWNNM算法的重建結(jié)果(圖7d)與原始數(shù)據(jù)更接近,振幅差值(圖7h)最小,進(jìn)一步表明APG-LWNNM方法的重建精度更高。

5 結(jié)論與建議

本文提出基于APG-LWNNM的地震數(shù)據(jù)重建方法,主要原理為低秩矩陣補全理論,得到如下結(jié)論:

1)APG-LWNNM算法重建結(jié)果精度高。在迭代過程中,對數(shù)權(quán)重向量的引入減少了數(shù)據(jù)矩陣奇異值的損失,盡可能保留了矩陣的信息量。APG-LWNNM算法的重建結(jié)果生成了更豐富的紋理,并且獲得了最高的信噪比和最低的重建誤差,在定性和定量分析上均有提升。

2)APG-LWNNM算法重建效率高。與深度學(xué)習(xí)方法相比,低秩矩陣補全方法無需學(xué)習(xí)地震數(shù)據(jù)特征擬合網(wǎng)絡(luò)參數(shù),可以加快數(shù)據(jù)重建速度,對硬件要求更友好。

3)APG-LWNNM算法對地震數(shù)據(jù)適應(yīng)性強。與深度學(xué)習(xí)方法相比,低秩矩陣補全方法可以對缺失地震數(shù)據(jù)直接重建,不需要大量的數(shù)據(jù)樣本訓(xùn)練多個隱藏層的網(wǎng)絡(luò)。

低秩性是低秩矩陣補全理論的前提,并且奇異值表示了矩陣的信息量。因此,后期研究方向是如何更高效地低秩預(yù)處理地震數(shù)據(jù)和如何在迭代過程中更多地保留矩陣奇異值。