雍 皓 韓 鐸 張俊杰 王俊秋*
(①吉林大學(xué)儀器科學(xué)與電氣工程學(xué)院,吉林長春 130026; ②地球信息探測儀器教育部重點(diǎn)實(shí)驗(yàn)室,吉林長春 130000)
在信號處理領(lǐng)域,通常采用壓縮技術(shù)解決多維海量數(shù)據(jù)的存儲和傳輸問題。Donoho[1]首先提出了壓縮感知(Compressed Sensing,CS)的概念。CS理論認(rèn)為,如果一個信號是稀疏的或可壓縮的,那么可利用遠(yuǎn)少于奈奎斯特采樣數(shù)量的測量值無失真地重建[2]。為解決地震勘探采集的地震數(shù)據(jù)因通信帶寬和功耗等限制難以實(shí)時回傳這一難題,應(yīng)用CS理論將原始地震數(shù)據(jù)壓縮后傳輸,從而減少地震數(shù)據(jù)的傳輸量,提高傳輸效率[3-6]。由CS理論可知:地震數(shù)據(jù)越稀疏,重構(gòu)精度越高; 采用不同的稀疏表示方法對地震數(shù)據(jù)進(jìn)行稀疏表示,其稀疏度會有所不同[7-10]。
目前,地震數(shù)據(jù)稀疏表示方法大體上可分為兩類。一類是固定變換域稀疏表示方法,如傅里葉變換、離散余弦變換及小波變換等,均是通過固定公式生成字典。這類方法構(gòu)造字典的過程簡單、快速,但固定基一旦確定,就不能更改,因此稀疏表示效果較差,尤其是當(dāng)野外采集的地震信號含有噪聲、波形比較復(fù)雜時。另一類方法是字典學(xué)習(xí)型稀疏表示法,如最優(yōu)方向法(Method of Optimal Directions,MOD)和K—奇異值分解(K-Singular Value Decomposition,K-SVD)方法[11-12]等,均能學(xué)習(xí)地震數(shù)據(jù)的特點(diǎn),然后根據(jù)這些特征自適應(yīng)調(diào)整字典,因此對野外采集到的地震數(shù)據(jù)有較好的稀疏表示效果。MOD采用交替優(yōu)化方式對模型進(jìn)行求解,在字典更新階段采用的是全局更新策略,稀疏表示過程中產(chǎn)生的字典性能較低;K-SVD方法采用稀疏編碼和字典更新相互交替完成對信號特征的學(xué)習(xí),得益于其單獨(dú)更新字典原子的策略,對地震數(shù)據(jù)的稀疏表示效果較好,重構(gòu)精度較高[13-14]。
本文采用一種改進(jìn)的稀疏表示方法[15-19]對地震數(shù)據(jù)進(jìn)行處理重構(gòu)。該方法采用多層字典學(xué)習(xí)的方案,利用更多的數(shù)據(jù)特征去訓(xùn)練字典,增強(qiáng)稀疏表示方法的特征提取能力; 同時,在稀疏編碼階段設(shè)置三個終止條件,使其可以自適應(yīng)地確定不同地震數(shù)據(jù)所選擇的原子數(shù)。與K-SVD方法相比,在重構(gòu)算法相同的前提下,本文方法能夠?yàn)榛贑S的地震數(shù)據(jù)重構(gòu)提供更準(zhǔn)確的稀疏表示,對地震數(shù)據(jù)有更高的重構(gòu)精度。
稀疏表示理論作為CS理論的首個組成成分,為大多數(shù)信號的稀疏表示提供了新的采樣形式,不同的稀疏表示方法效果不同。作為一種典型的字典學(xué)習(xí)型稀疏表示方法,K-SVD可根據(jù)不同的地震數(shù)據(jù)生成不同的稀疏表示字典,并通過稀疏編碼和字典更新兩個過程交替完成對地震數(shù)據(jù)的學(xué)習(xí)過程。K-SVD算法的數(shù)學(xué)表達(dá)式為
i=1,2,…,N
(1)
K-SVD對地震數(shù)據(jù)進(jìn)行稀疏表示通過以下三步進(jìn)行。
(1)初始化,即從原始地震數(shù)據(jù)中選取少部分?jǐn)?shù)據(jù)作為字典的初始值。
(2)稀疏編碼,即固定初始字典D,采用正交匹配追蹤算法(OMP)[20-21]求解稀疏矩陣X。求解過程可表示為
(2)
(3)字典更新,即固定稀疏系數(shù)矩陣X,更新字典D。其過程為
(3)
Ej,R=UΣVT
(4)
重復(fù)步驟(2)~步驟(3),直至收斂到指定的誤差或達(dá)到指定的迭代次數(shù)。
(1)在K-SVD方法中,每個訓(xùn)練樣本的K是固定的。不同樣本所包含特征的信息量可能不同,對于一個包含大量特征的訓(xùn)練樣本,如果K值太小,就沒有足夠的原子表示這些特征,稱為欠擬合; 如果K值過大,多余的原子會引起表示的過擬合,從而影響稀疏性。
(2)在字典學(xué)習(xí)過程中,不同的目標(biāo)樣本會產(chǎn)生不同的殘差,進(jìn)而形成客觀特征差異較大的殘差層。在K-SVD方法中,由于不同的殘差層采用同一個原子表示,因此,一個原子必須在不止一種類型的客觀特征之間妥協(xié),降低了CS重構(gòu)的精度。
針對K-SVD方法的局限性,本文采用AMDL對地震數(shù)據(jù)進(jìn)行稀疏表示。
AMDL方法已廣泛應(yīng)用于圖像處理、視覺識別等領(lǐng)域,并獲得了較理想的結(jié)果[22],但在地震勘探領(lǐng)域中的應(yīng)用至今未見。
AMDL方法首先將字典分成若干子字典,然后分別求出每個子字典及相應(yīng)的稀疏系數(shù),最后用它們的乘積之和表示地震數(shù)據(jù)。其稀疏表示模型為
(5)
式中:cmax為字典的最大層數(shù);Dc為第c層的子字典;Xc為第c層的稀疏系數(shù)矩陣。每個子字典Dc包含qh個原子,因此原子總數(shù)q=qhcmax。
圖1給出了字典最大層數(shù)為2時的地震數(shù)據(jù)AMDL方法流程。AMDL算法在處理地震數(shù)據(jù)時,首先輸入原始地震數(shù)據(jù)Y、層數(shù)cmax、每層最大稀疏水平qmax、每層最大迭代循環(huán)數(shù)tmax、系數(shù)ε以及閾值λ。第一層直接利用原始訓(xùn)練樣本Y1=Y訓(xùn)練子字典D1,第一層的殘差r1=Y1-D1X1,將其作為第二層子字典的訓(xùn)練樣本Y2=r1來訓(xùn)練子字典D2。
圖1 AMDL方法地震數(shù)據(jù)稀疏表示流程
與K-SVD類似,AMDL方法中每層子字典的訓(xùn)練均分為子字典初始化、稀疏編碼和字典更新三個階段。
(1)初始化階段。從地震數(shù)據(jù)Y中選取少部分?jǐn)?shù)據(jù)作為第一層子字典的訓(xùn)練樣本Y1,并將Y1作為子字典D1的初始值; 每層交替執(zhí)行t次稀疏編碼和字典更新,直至收斂到指定的誤差,或達(dá)到指定的迭代次數(shù)。
(2)稀疏編碼階段。首先,固定初始化子字典Dc,求解稀疏系數(shù)矩陣Xc,即
(6)
(3)字典更新階段。固定稀疏系數(shù)矩陣Xc,采用K-SVD的更新策略進(jìn)行Dc的原子更新。
每層重復(fù)執(zhí)行上述階段,直至遍歷設(shè)定的層數(shù)cmax,終止計(jì)算過程,最后輸出更新后的字典D=[D1,D2,…,Dcmax]及與其相對應(yīng)的稀疏矩陣系數(shù)X=[X1,X2,…,Xcmax]。
AMDL方法采用的多層字典學(xué)習(xí)方案除了考慮原始訓(xùn)練樣本包含的特征外,還考慮了不同層的樣本殘差所包含的特征,使其能夠在字典學(xué)習(xí)過程中充分利用不同層次的特征,提高地震數(shù)據(jù)的重構(gòu)精度; 同時,自適應(yīng)策略可以確定每個地震數(shù)據(jù)在多層模式下使用原子的數(shù)量,使稀疏度自適應(yīng)于每個訓(xùn)練樣本,減少地震數(shù)據(jù)在稀疏表示過程中的欠擬合或過擬合問題,提高地震數(shù)據(jù)的稀疏性,進(jìn)而提高重構(gòu)精度。
K-SVD方法的目標(biāo)是盡量減少表示殘差。而AMDL考慮兩個方面:一是盡可能減少表示殘差; 二是減少地震數(shù)據(jù)稀疏表示中使用的原子數(shù)。換言之,AMDL的主要任務(wù)是使用盡可能少的原子獲得更精確的表示。
壓縮比(測量矩陣的行、列數(shù)之比(M/N))是描述對原始地震數(shù)據(jù)的壓縮程度,壓縮比越大,表示壓縮程度越小,壓縮后的數(shù)據(jù)越多,易于重構(gòu); 反之,壓縮比越小,表示壓縮程度越高,壓縮后的數(shù)據(jù)越少,不易于重構(gòu)。本試驗(yàn)采用了兩個地區(qū)的實(shí)測地震數(shù)據(jù),同時采取控制變量法,利用K-SVD和AMDL兩種方法在不同壓縮比下分別對其進(jìn)行稀疏表示,然后采用相同的重構(gòu)算法進(jìn)行重構(gòu),通過計(jì)算相同壓縮比下的重構(gòu)數(shù)據(jù)信噪比(SNR)和相對誤差(Error)評價兩種方法。
重構(gòu)地震數(shù)據(jù)的信噪比計(jì)算公式為
(7)
重構(gòu)數(shù)據(jù)的相對誤差計(jì)算公式為
(8)
選取采樣率為1ms、采樣時間為0.256s、總計(jì)256道的地震數(shù)據(jù)進(jìn)行試驗(yàn)。首先,利用K-SVD方法和AMDL方法在不同壓縮比下,分別對地震數(shù)據(jù)進(jìn)行稀疏表示; 然后,利用基追蹤[23-25]重構(gòu)算法進(jìn)行地震數(shù)據(jù)重構(gòu)。試驗(yàn)中,將AMDL方法中的幾個參數(shù)設(shè)定為:層數(shù)cmax=4,每層最大稀疏水平qh=10,每層最大迭代循環(huán)數(shù)tmax=3,系數(shù)ε=0.5,閾值λ=0.0001。
圖2a為原始地震數(shù)據(jù),圖2b~圖2d和圖2e~圖2g分別為經(jīng)過AMDL方法和K-SVD方法稀疏表示后,壓縮比為0.1、0.3和0.5時的重構(gòu)地震數(shù)據(jù)。可以明顯看出,在壓縮比為0.1時,經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)和原始地震數(shù)據(jù)相似度非常低,重構(gòu)精度遠(yuǎn)低于經(jīng)過AMDL方法稀疏表示后的重構(gòu)精度; 在壓縮比為0.3時,經(jīng)過K-SVD方法稀疏表示后地震數(shù)據(jù)的重構(gòu)精度有明顯的提升,但仍低于同壓縮比下AMDL方法的重構(gòu)精度; 在壓縮比為0.5時,經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)和原始地震數(shù)據(jù)相似度較高,重構(gòu)精度與AMDL方法在壓縮比為0.1時的地震數(shù)據(jù)重構(gòu)精度大致相同。
為了更加直觀地展示地震數(shù)據(jù)的重構(gòu)精度,從原始數(shù)據(jù)和不同壓縮比下的重構(gòu)數(shù)據(jù)中選取第100、125和150道數(shù)據(jù)進(jìn)行波形對比。圖3a~圖3c分別為原始地震數(shù)據(jù)和在不同壓縮比下經(jīng)過兩種方法稀疏表示后重構(gòu)地震數(shù)據(jù)中第100道、125道和150道數(shù)據(jù)的波形對比圖,a代表原始地震數(shù)據(jù),b、d和f分別代表壓縮比為0.1、0.3和0.5時經(jīng)過K-SVD方法稀疏表示后重構(gòu)地震數(shù)據(jù),c、e和g分別代表壓縮比為0.1、0.3和0.5時經(jīng)過AMDL方法稀疏表示后重構(gòu)地震數(shù)據(jù)。圖4a、圖4b分別為地震數(shù)據(jù)經(jīng)過兩種不同方法稀疏表示后,在不同壓縮比下重構(gòu)信噪比和相對誤差的對比圖。可以看出,壓縮比為0.1時,AMDL方法的重構(gòu)信噪比比K-SVD方法高約10dB,提高了近3.5倍,重構(gòu)誤差比K-SVD方法低約0.8,降低了近80%; 在壓縮比為0.2時,信噪比提高約20dB,提高了近3.5倍,且重構(gòu)誤差減小約0.53,重構(gòu)誤差降低了近90%; 在其它壓縮比下,信噪比提高了約1~2倍,重構(gòu)誤差降低了80%~90%。可見,在相同壓縮比下且使用相同的重構(gòu)算法前提下,AMDL方法對地震數(shù)據(jù)稀疏表示后重構(gòu)的地震數(shù)據(jù)信噪比更高,相對誤差更小。
圖2 遼寧興城實(shí)測地震數(shù)據(jù)及兩種方法在不同壓縮比下稀疏表示后重構(gòu)的地震數(shù)據(jù)(a)原始地震數(shù)據(jù); (b)AMDL法,壓縮比為0.1; (c)AMDL法,壓縮比為0.3; (d)AMDL法,壓縮比為0.5;(e)K-SVD,壓縮比為0.1; (f)K-SVD,壓縮比為0.3; (g)K-SVD,壓縮比為0.5
為進(jìn)一步證明兩種方法的稀疏表示效果,選取甘肅金昌銅鎳礦的實(shí)測地震數(shù)據(jù)做了相同的對比試驗(yàn),本數(shù)據(jù)由500N的輕便小型可控震源產(chǎn)生,且該實(shí)測地震數(shù)據(jù)經(jīng)過預(yù)處理。數(shù)據(jù)的采樣率為1ms,采樣時間為2s,道數(shù)為72道。同樣利用不同壓縮比下兩種方法分別對地震數(shù)據(jù)進(jìn)行稀疏表示。然后利用基追蹤重構(gòu)算法進(jìn)行數(shù)據(jù)重構(gòu)。AMDL方法的參數(shù)與3.1相同。
圖5a為原始地震數(shù)據(jù)。圖5b~圖5c分別為經(jīng)AMDL方法稀疏表示后壓縮比分別為0.1和0.3時的重構(gòu)地震數(shù)據(jù)。圖5d~圖5e分別為經(jīng)過K-SVD方法稀疏表示后、壓縮比分別為0.1和0.3時的重構(gòu)地震數(shù)據(jù)。由圖可見,壓縮比在0.1時,經(jīng)過AMDL方法稀疏表示后重構(gòu)的地震數(shù)據(jù)具有很高的重構(gòu)精度,可清楚地看出反射層的位置(圖5b箭頭所示),而經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)與原始地震數(shù)據(jù)相似度較低,只能隱約看出反射層的位置(圖5d箭頭所示)。而壓縮比為0.3時,經(jīng)過K-SVD方法稀疏表示后地震數(shù)據(jù)的重構(gòu)精度有明顯提升,也可以清晰看出反射層位置(圖5e箭頭所示),且重構(gòu)精度與AMDL方法在壓縮比為0.1時的地震數(shù)據(jù)重構(gòu)精度幾乎相同。圖6a和圖6b分別為經(jīng)兩種不同方法稀疏表示后、在不同壓縮比下重構(gòu)的地震數(shù)據(jù)信噪比和相對誤差的對比,同樣可以看出,壓縮比相同時,AMDL算法對地震數(shù)據(jù)稀疏表示后重構(gòu)的地震數(shù)據(jù)信噪比更高、相對誤差更小。
圖3 遼寧興城實(shí)測地震數(shù)據(jù)不同道波形與兩種方法在不同壓縮比下稀疏表示后重構(gòu)數(shù)據(jù)波形的對比(a)第100道; (b)第125道; (c)第150道
圖4 遼寧興城實(shí)測地震數(shù)據(jù)在不同壓縮比下兩種方法稀疏表示后的重構(gòu)數(shù)據(jù)的信噪比(a)和相對誤差(b)對比
圖5 甘肅金昌實(shí)測地震數(shù)據(jù)與兩種方法在不同壓縮比下稀疏表示后重構(gòu)數(shù)據(jù)的對比(a)原始地震數(shù)據(jù); (b)AMDL法,壓縮比為0.1; (c)AMDL法,壓縮比為0.3; (d)K-SVD,壓縮比為0.1; (e)K-SVD,壓縮比為0.3
圖6 甘肅金昌實(shí)測數(shù)據(jù)不同壓縮比時兩種方法稀疏表示后的重構(gòu)數(shù)據(jù)的信噪比(a)和相對誤差(b)對比
兩組數(shù)據(jù)的試驗(yàn)結(jié)果可以證明,在處理相同的實(shí)測地震數(shù)據(jù),且使用相同重構(gòu)算法的前提下,相比于K-SVD方法,AMDL方法對地震數(shù)據(jù)的稀疏表示精度更高、效果更好。
本文將自適應(yīng)多層字典學(xué)習(xí)(AMDL)方法用于地震數(shù)據(jù)的稀疏表示,旨在字典學(xué)習(xí)過程中利用不同的特征層訓(xùn)練字典,且自適應(yīng)地確定每個訓(xùn)練樣本使用的原子數(shù),不僅可利用更多的特征訓(xùn)練字典,而且字典不會受到其他特征層的影響,同時避免了稀疏表示過程中出現(xiàn)的欠擬合和過擬合問題。實(shí)測地震數(shù)據(jù)重構(gòu)試驗(yàn)結(jié)果表明,AMDL方法具有較強(qiáng)的稀疏表示能力,相比于K-SVD方法,其在地震數(shù)據(jù)稀疏表示方面更具優(yōu)越性。