雍 皓 韓 鐸 張俊杰 王俊秋*

(①吉林大學(xué)儀器科學(xué)與電氣工程學(xué)院，吉林長春 130026； ②地球信息探測儀器教育部重點(diǎn)實(shí)驗(yàn)室，吉林長春 130000)

0 引言

在信號(hào)處理領(lǐng)域，通常采用壓縮技術(shù)解決多維海量數(shù)據(jù)的存儲(chǔ)和傳輸問題。Donoho[1]首先提出了壓縮感知(Compressed Sensing，CS)的概念。CS理論認(rèn)為，如果一個(gè)信號(hào)是稀疏的或可壓縮的，那么可利用遠(yuǎn)少于奈奎斯特采樣數(shù)量的測量值無失真地重建[2]。為解決地震勘探采集的地震數(shù)據(jù)因通信帶寬和功耗等限制難以實(shí)時(shí)回傳這一難題，應(yīng)用CS理論將原始地震數(shù)據(jù)壓縮后傳輸，從而減少地震數(shù)據(jù)的傳輸量，提高傳輸效率[3-6]。由CS理論可知：地震數(shù)據(jù)越稀疏，重構(gòu)精度越高； 采用不同的稀疏表示方法對(duì)地震數(shù)據(jù)進(jìn)行稀疏表示，其稀疏度會(huì)有所不同[7-10]。

目前，地震數(shù)據(jù)稀疏表示方法大體上可分為兩類。一類是固定變換域稀疏表示方法，如傅里葉變換、離散余弦變換及小波變換等，均是通過固定公式生成字典。這類方法構(gòu)造字典的過程簡單、快速，但固定基一旦確定，就不能更改，因此稀疏表示效果較差，尤其是當(dāng)野外采集的地震信號(hào)含有噪聲、波形比較復(fù)雜時(shí)。另一類方法是字典學(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)整字典，因此對(duì)野外采集到的地震數(shù)據(jù)有較好的稀疏表示效果。MOD采用交替優(yōu)化方式對(duì)模型進(jìn)行求解，在字典更新階段采用的是全局更新策略，稀疏表示過程中產(chǎn)生的字典性能較低；K-SVD方法采用稀疏編碼和字典更新相互交替完成對(duì)信號(hào)特征的學(xué)習(xí)，得益于其單獨(dú)更新字典原子的策略，對(duì)地震數(shù)據(jù)的稀疏表示效果較好，重構(gòu)精度較高[13-14]。

本文采用一種改進(jìn)的稀疏表示方法[15-19]對(duì)地震數(shù)據(jù)進(jìn)行處理重構(gòu)。該方法采用多層字典學(xué)習(xí)的方案，利用更多的數(shù)據(jù)特征去訓(xùn)練字典，增強(qiáng)稀疏表示方法的特征提取能力； 同時(shí)，在稀疏編碼階段設(shè)置三個(gè)終止條件，使其可以自適應(yīng)地確定不同地震數(shù)據(jù)所選擇的原子數(shù)。與K-SVD方法相比，在重構(gòu)算法相同的前提下，本文方法能夠?yàn)榛贑S的地震數(shù)據(jù)重構(gòu)提供更準(zhǔn)確的稀疏表示，對(duì)地震數(shù)據(jù)有更高的重構(gòu)精度。

1 地震數(shù)據(jù)K-SVD稀疏表示

1.1 地震數(shù)據(jù)K-SVD稀疏表示法

稀疏表示理論作為CS理論的首個(gè)組成成分，為大多數(shù)信號(hào)的稀疏表示提供了新的采樣形式，不同的稀疏表示方法效果不同。作為一種典型的字典學(xué)習(xí)型稀疏表示方法，K-SVD可根據(jù)不同的地震數(shù)據(jù)生成不同的稀疏表示字典，并通過稀疏編碼和字典更新兩個(gè)過程交替完成對(duì)地震數(shù)據(jù)的學(xué)習(xí)過程。K-SVD算法的數(shù)學(xué)表達(dá)式為

i=1,2，…，N

(1)

K-SVD對(duì)地震數(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.2 K-SVD的缺陷

(1)在K-SVD方法中，每個(gè)訓(xùn)練樣本的K是固定的。不同樣本所包含特征的信息量可能不同，對(duì)于一個(gè)包含大量特征的訓(xùn)練樣本，如果K值太小，就沒有足夠的原子表示這些特征，稱為欠擬合； 如果K值過大，多余的原子會(huì)引起表示的過擬合，從而影響稀疏性。

(2)在字典學(xué)習(xí)過程中，不同的目標(biāo)樣本會(huì)產(chǎn)生不同的殘差，進(jìn)而形成客觀特征差異較大的殘差層。在K-SVD方法中，由于不同的殘差層采用同一個(gè)原子表示，因此，一個(gè)原子必須在不止一種類型的客觀特征之間妥協(xié)，降低了CS重構(gòu)的精度。

針對(duì)K-SVD方法的局限性，本文采用AMDL對(duì)地震數(shù)據(jù)進(jìn)行稀疏表示。

2 地震數(shù)據(jù)AMDL稀疏表示方法

AMDL方法已廣泛應(yīng)用于圖像處理、視覺識(shí)別等領(lǐng)域，并獲得了較理想的結(jié)果[22]，但在地震勘探領(lǐng)域中的應(yīng)用至今未見。

AMDL方法首先將字典分成若干子字典，然后分別求出每個(gè)子字典及相應(yīng)的稀疏系數(shù)，最后用它們的乘積之和表示地震數(shù)據(jù)。其稀疏表示模型為

(5)

式中：cmax為字典的最大層數(shù)；Dc為第c層的子字典；Xc為第c層的稀疏系數(shù)矩陣。每個(gè)子字典Dc包含qh個(gè)原子，因此原子總數(shù)q=qhcmax。

圖1給出了字典最大層數(shù)為2時(shí)的地震數(shù)據(jù)AMDL方法流程。AMDL算法在處理地震數(shù)據(jù)時(shí)，首先輸入原始地震數(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)練均分為子字典初始化、稀疏編碼和字典更新三個(gè)階段。

(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]及與其相對(duì)應(yīng)的稀疏矩陣系數(shù)X=[X1,X2,…,Xcmax]。

AMDL方法采用的多層字典學(xué)習(xí)方案除了考慮原始訓(xùn)練樣本包含的特征外，還考慮了不同層的樣本殘差所包含的特征，使其能夠在字典學(xué)習(xí)過程中充分利用不同層次的特征，提高地震數(shù)據(jù)的重構(gòu)精度； 同時(shí)，自適應(yīng)策略可以確定每個(gè)地震數(shù)據(jù)在多層模式下使用原子的數(shù)量，使稀疏度自適應(yīng)于每個(gè)訓(xùn)練樣本，減少地震數(shù)據(jù)在稀疏表示過程中的欠擬合或過擬合問題，提高地震數(shù)據(jù)的稀疏性，進(jìn)而提高重構(gòu)精度。

K-SVD方法的目標(biāo)是盡量減少表示殘差。而AMDL考慮兩個(gè)方面：一是盡可能減少表示殘差； 二是減少地震數(shù)據(jù)稀疏表示中使用的原子數(shù)。換言之，AMDL的主要任務(wù)是使用盡可能少的原子獲得更精確的表示。

3 試驗(yàn)結(jié)果與分析

壓縮比(測量矩陣的行、列數(shù)之比(M/N))是描述對(duì)原始地震數(shù)據(jù)的壓縮程度，壓縮比越大，表示壓縮程度越小，壓縮后的數(shù)據(jù)越多，易于重構(gòu)； 反之，壓縮比越小，表示壓縮程度越高，壓縮后的數(shù)據(jù)越少，不易于重構(gòu)。本試驗(yàn)采用了兩個(gè)地區(qū)的實(shí)測地震數(shù)據(jù)，同時(shí)采取控制變量法，利用K-SVD和AMDL兩種方法在不同壓縮比下分別對(duì)其進(jìn)行稀疏表示，然后采用相同的重構(gòu)算法進(jìn)行重構(gòu)，通過計(jì)算相同壓縮比下的重構(gòu)數(shù)據(jù)信噪比(SNR)和相對(duì)誤差(Error)評(píng)價(jià)兩種方法。

重構(gòu)地震數(shù)據(jù)的信噪比計(jì)算公式為

(7)

重構(gòu)數(shù)據(jù)的相對(duì)誤差計(jì)算公式為

(8)

3.1 遼寧興城地質(zhì)走廊帶實(shí)測地震數(shù)據(jù)試驗(yàn)

選取采樣率為1ms、采樣時(shí)間為0.256s、總計(jì)256道的地震數(shù)據(jù)進(jìn)行試驗(yàn)。首先，利用K-SVD方法和AMDL方法在不同壓縮比下，分別對(duì)地震數(shù)據(jù)進(jìn)行稀疏表示； 然后，利用基追蹤[23-25]重構(gòu)算法進(jìn)行地震數(shù)據(jù)重構(gòu)。試驗(yàn)中，將AMDL方法中的幾個(gè)參數(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時(shí)的重構(gòu)地震數(shù)據(jù)?？梢悦黠@看出，在壓縮比為0.1時(shí)，經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)和原始地震數(shù)據(jù)相似度非常低，重構(gòu)精度遠(yuǎn)低于經(jīng)過AMDL方法稀疏表示后的重構(gòu)精度； 在壓縮比為0.3時(shí)，經(jīng)過K-SVD方法稀疏表示后地震數(shù)據(jù)的重構(gòu)精度有明顯的提升，但仍低于同壓縮比下AMDL方法的重構(gòu)精度； 在壓縮比為0.5時(shí)，經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)和原始地震數(shù)據(jù)相似度較高，重構(gòu)精度與AMDL方法在壓縮比為0.1時(shí)的地震數(shù)據(jù)重構(gòu)精度大致相同。

為了更加直觀地展示地震數(shù)據(jù)的重構(gòu)精度，從原始數(shù)據(jù)和不同壓縮比下的重構(gòu)數(shù)據(jù)中選取第100、125和150道數(shù)據(jù)進(jìn)行波形對(duì)比。圖3a～圖3c分別為原始地震數(shù)據(jù)和在不同壓縮比下經(jīng)過兩種方法稀疏表示后重構(gòu)地震數(shù)據(jù)中第100道、125道和150道數(shù)據(jù)的波形對(duì)比圖，a代表原始地震數(shù)據(jù)，b、d和f分別代表壓縮比為0.1、0.3和0.5時(shí)經(jīng)過K-SVD方法稀疏表示后重構(gòu)地震數(shù)據(jù)，c、e和g分別代表壓縮比為0.1、0.3和0.5時(shí)經(jīng)過AMDL方法稀疏表示后重構(gòu)地震數(shù)據(jù)。圖4a、圖4b分別為地震數(shù)據(jù)經(jīng)過兩種不同方法稀疏表示后，在不同壓縮比下重構(gòu)信噪比和相對(duì)誤差的對(duì)比圖?？梢钥闯?，壓縮比為0.1時(shí)，AMDL方法的重構(gòu)信噪比比K-SVD方法高約10dB，提高了近3.5倍，重構(gòu)誤差比K-SVD方法低約0.8，降低了近80%； 在壓縮比為0.2時(shí)，信噪比提高約20dB，提高了近3.5倍，且重構(gòu)誤差減小約0.53，重構(gòu)誤差降低了近90%； 在其它壓縮比下，信噪比提高了約1～2倍，重構(gòu)誤差降低了80%～90%?？梢?，在相同壓縮比下且使用相同的重構(gòu)算法前提下，AMDL方法對(duì)地震數(shù)據(jù)稀疏表示后重構(gòu)的地震數(shù)據(jù)信噪比更高，相對(duì)誤差更小。

圖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

3.2 甘肅金昌銅鎳礦的實(shí)測地震數(shù)據(jù)試驗(yàn)

為進(jìn)一步證明兩種方法的稀疏表示效果，選取甘肅金昌銅鎳礦的實(shí)測地震數(shù)據(jù)做了相同的對(duì)比試驗(yàn)，本數(shù)據(jù)由500N的輕便小型可控震源產(chǎn)生，且該實(shí)測地震數(shù)據(jù)經(jīng)過預(yù)處理。數(shù)據(jù)的采樣率為1ms，采樣時(shí)間為2s，道數(shù)為72道。同樣利用不同壓縮比下兩種方法分別對(duì)地震數(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時(shí)的重構(gòu)地震數(shù)據(jù)。圖5d～圖5e分別為經(jīng)過K-SVD方法稀疏表示后、壓縮比分別為0.1和0.3時(shí)的重構(gòu)地震數(shù)據(jù)。由圖可見，壓縮比在0.1時(shí)，經(jīng)過AMDL方法稀疏表示后重構(gòu)的地震數(shù)據(jù)具有很高的重構(gòu)精度，可清楚地看出反射層的位置(圖5b箭頭所示)，而經(jīng)過K-SVD方法稀疏表示后重構(gòu)的地震數(shù)據(jù)與原始地震數(shù)據(jù)相似度較低，只能隱約看出反射層的位置(圖5d箭頭所示)。而壓縮比為0.3時(shí)，經(jīng)過K-SVD方法稀疏表示后地震數(shù)據(jù)的重構(gòu)精度有明顯提升，也可以清晰看出反射層位置(圖5e箭頭所示)，且重構(gòu)精度與AMDL方法在壓縮比為0.1時(shí)的地震數(shù)據(jù)重構(gòu)精度幾乎相同。圖6a和圖6b分別為經(jīng)兩種不同方法稀疏表示后、在不同壓縮比下重構(gòu)的地震數(shù)據(jù)信噪比和相對(duì)誤差的對(duì)比，同樣可以看出，壓縮比相同時(shí)，AMDL算法對(duì)地震數(shù)據(jù)稀疏表示后重構(gòu)的地震數(shù)據(jù)信噪比更高、相對(duì)誤差更小。

圖3 遼寧興城實(shí)測地震數(shù)據(jù)不同道波形與兩種方法在不同壓縮比下稀疏表示后重構(gòu)數(shù)據(jù)波形的對(duì)比(a)第100道； (b)第125道； (c)第150道

圖4 遼寧興城實(shí)測地震數(shù)據(jù)在不同壓縮比下兩種方法稀疏表示后的重構(gòu)數(shù)據(jù)的信噪比(a)和相對(duì)誤差(b)對(duì)比

圖5 甘肅金昌實(shí)測地震數(shù)據(jù)與兩種方法在不同壓縮比下稀疏表示后重構(gòu)數(shù)據(jù)的對(duì)比(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ù)不同壓縮比時(shí)兩種方法稀疏表示后的重構(gòu)數(shù)據(jù)的信噪比(a)和相對(duì)誤差(b)對(duì)比

兩組數(shù)據(jù)的試驗(yàn)結(jié)果可以證明，在處理相同的實(shí)測地震數(shù)據(jù)，且使用相同重構(gòu)算法的前提下，相比于K-SVD方法，AMDL方法對(duì)地震數(shù)據(jù)的稀疏表示精度更高、效果更好。

4 結(jié)束語

本文將自適應(yīng)多層字典學(xué)習(xí)(AMDL)方法用于地震數(shù)據(jù)的稀疏表示，旨在字典學(xué)習(xí)過程中利用不同的特征層訓(xùn)練字典，且自適應(yīng)地確定每個(gè)訓(xùn)練樣本使用的原子數(shù)，不僅可利用更多的特征訓(xùn)練字典，而且字典不會(huì)受到其他特征層的影響，同時(shí)避免了稀疏表示過程中出現(xiàn)的欠擬合和過擬合問題。實(shí)測地震數(shù)據(jù)重構(gòu)試驗(yàn)結(jié)果表明，AMDL方法具有較強(qiáng)的稀疏表示能力，相比于K-SVD方法，其在地震數(shù)據(jù)稀疏表示方面更具優(yōu)越性。