壓縮感知走進(jìn)地球物理勘探

馬堅(jiān)偉

(哈爾濱工業(yè)大學(xué)地球物理中心/數(shù)學(xué)系,黑龍江哈爾濱150001)

壓縮感知(Compressed Sensing,CS)描述的是可從高度不完備的線性測(cè)量中高精度重構(gòu)未知目標(biāo)。創(chuàng)造性的把L1范數(shù)最小化和隨機(jī)矩陣?yán)碚撚袡C(jī)結(jié)合,可得到稀疏信號(hào)重建的最佳效果。利用壓縮感知重構(gòu)一個(gè)連續(xù)信號(hào),不再與香農(nóng)-奈奎斯特采樣定理所要求的“頻帶”有關(guān),而是與未知信號(hào)的“稀疏度”有關(guān)。CS由DONOHO,CANDS,TAO三位著名數(shù)學(xué)家提出[1-2],并入選了2007年美國(guó)十大科技進(jìn)展,早期的主要幾篇文獻(xiàn)的引用已超過(guò)5×104次。

從數(shù)學(xué)上來(lái)理解,CS就是求解一個(gè)線性方程組Y=AX。Y是已知的觀察數(shù)據(jù),X是要重構(gòu)的目標(biāo)。特別強(qiáng)調(diào)的是,這里重構(gòu)信號(hào)A是一個(gè)非方陣,其矩陣的行數(shù)遠(yuǎn)少于列數(shù)。這相當(dāng)于這個(gè)線性方程組的方程個(gè)數(shù)遠(yuǎn)少于未知數(shù)的個(gè)數(shù)。顯然,這樣的方程組有無(wú)窮多個(gè)解,是一欠定反問(wèn)題。但是CS告訴我們,如果先驗(yàn)知道未知目標(biāo)X有一定的稀疏性或者具有可壓縮性,且A的列之間具有一定的非相關(guān)性(比如隨機(jī)矩陣或勘探中的隨機(jī)接收器布點(diǎn)),那么可以用稀疏約束優(yōu)化的算法X=min‖Y-AX‖2+a‖TX‖0精確求解X。其中,T表示稀疏變換,零范數(shù)表示稀疏約束正則化。

滿(mǎn)足這個(gè)精確重構(gòu)所需要的測(cè)量數(shù)M可以表示為M=Slg(N/S),其中N表示待重構(gòu)信號(hào)的長(zhǎng)度,S表示該信號(hào)的稀疏度。即便在有噪聲的情況下,即Y=AX+Noise,也可以進(jìn)行高精度重構(gòu)。我們補(bǔ)充上面提到的兩個(gè)概念:稀疏和可壓縮性。稀疏是指一個(gè)信號(hào)如果其大部分振幅值為零,只有少數(shù)元素的振幅為非零,則稱(chēng)該信號(hào)為稀疏信號(hào);可壓縮性是指在某種變換域下的稀疏,即該信號(hào)在時(shí)空域下不稀疏,但可以找到一種稀疏變換,如傅里葉變換、小波變換、曲波變換等,使得該信號(hào)經(jīng)過(guò)變換后的系數(shù)只有少數(shù)非零值。我們可以從如下幾個(gè)角度來(lái)理解CS:從數(shù)學(xué)上講,CS本質(zhì)是降維,從低維空間去研究高維空間;從信號(hào)上講,CS本質(zhì)是采樣,從香農(nóng)采樣的頻率相關(guān)到稀疏度相關(guān);從工程上講,CS本質(zhì)是成本,從物理測(cè)量成本轉(zhuǎn)移到數(shù)學(xué)計(jì)算成本,用數(shù)學(xué)計(jì)算來(lái)彌補(bǔ)實(shí)際勘探的采集不足。

從20世紀(jì)40年代開(kāi)始,統(tǒng)治信號(hào)和信息領(lǐng)域的一直是香農(nóng)-奈奎斯特采樣定理,即要想從離散信號(hào)中恢復(fù)連續(xù)信號(hào),其離散信號(hào)的采樣率至少應(yīng)該是此連續(xù)信號(hào)最高頻率的兩倍。可以說(shuō)大部分相關(guān)信號(hào)采集的設(shè)備都是基于香農(nóng)定理設(shè)計(jì)的,如相機(jī)、雷達(dá)、核磁共振等醫(yī)學(xué)成像、通訊、信息等眾多工程領(lǐng)域的設(shè)備。CS理論出現(xiàn)后,其采樣率不再和信號(hào)的頻帶有關(guān),而是和信號(hào)的稀疏度有關(guān),從而打破了傳統(tǒng)理論的束縛。因此很多基于傳統(tǒng)理論設(shè)計(jì)的方法、軟件和設(shè)備都可以得到升級(jí)換代。2007年,萊斯大學(xué)提出的單像素相機(jī)就是一個(gè)很好的例子,不再去追求千萬(wàn)像素的分辨率,而是基于壓縮感知理論用一個(gè)像素的時(shí)間序列成像就可重構(gòu)出高分辨率圖像。

壓縮感知是一個(gè)理論框架,不是單一的技術(shù),其三要素是隨機(jī)測(cè)量、目標(biāo)的稀疏表示和稀疏促進(jìn)的優(yōu)化重構(gòu)算法。地震勘探中,降低野外數(shù)據(jù)采集的成本而又能保證勘探精度是一個(gè)很重要的課題。反過(guò)來(lái)說(shuō),從現(xiàn)有采集的地震數(shù)據(jù)中能否挖掘出更多的寶貴信息來(lái)提高勘探精度？由于CS理論的出現(xiàn),這個(gè)問(wèn)題由香農(nóng)采樣轉(zhuǎn)變到稀疏重構(gòu)的思路上,帶來(lái)了很多變革。CS的成功,在于改變了信號(hào)的采集方式。CS所涉及的技術(shù),在地震勘探的采集、處理、正演模擬、成像反演等方面都會(huì)帶來(lái)變革,但其核心問(wèn)題之一還是地震數(shù)據(jù)和地質(zhì)目標(biāo)的稀疏表示。

1 地震勘探中壓縮感知的最新研究進(jìn)展

學(xué)術(shù)期刊Geophysics已多次組織了與壓縮感知地震勘探相關(guān)的專(zhuān)題,如2010年關(guān)于采樣和波場(chǎng)表示的專(zhuān)題、2015關(guān)于地震數(shù)據(jù)表示的專(zhuān)題。TheLeadingEdge也于2017年8月刊出了壓縮感知和高效采集的專(zhuān)刊。壓縮感知框架已經(jīng)滲入了傳統(tǒng)地震勘探的各個(gè)環(huán)節(jié),例如地震數(shù)據(jù)的采集、壓縮、重構(gòu)、模擬、偏移、反演等,并且在很多新的勘探技術(shù)中也發(fā)揮了作用,例如同時(shí)震源勘探、延時(shí)勘探等。不管應(yīng)用于何種技術(shù),地震數(shù)據(jù)的稀疏表示和隨機(jī)采樣始終是壓縮感知應(yīng)用的重要理論基礎(chǔ)。本章節(jié)首先介紹地震數(shù)據(jù)稀疏表示和隨機(jī)采樣的最新進(jìn)展,然后分別從傳統(tǒng)勘探技術(shù)和新的勘探技術(shù)兩個(gè)方面介紹壓縮感知在地震勘探應(yīng)用中的最新進(jìn)展。

1.1 稀疏變換

地震數(shù)據(jù)的稀疏表示按來(lái)源可分為3類(lèi)。第一類(lèi)是從圖像處理技術(shù)中借用已有的稀疏變換,如曲波變換[3-4]、剪切波變換[5-6]等等。曲波變換與剪切波變換均為多尺度多方向小波變換,適用于表示圖像中的線狀特征,因此也適用于表示地震數(shù)據(jù)的同相軸特征。第二類(lèi)是利用字典學(xué)習(xí)方法從地震數(shù)據(jù)本身自適應(yīng)構(gòu)造稀疏變換,如K-均值奇異值分解(K-SVD)方法[7]、數(shù)據(jù)驅(qū)動(dòng)緊框架方法(data-driven tight frame,DDT)[8]等等。這類(lèi)方法通常在塊尺度下進(jìn)行操作,適用于處理具有復(fù)雜結(jié)構(gòu)的地震數(shù)據(jù)。第三類(lèi)是直接針對(duì)地震數(shù)據(jù)設(shè)計(jì)稀疏變換。物理小波變換[9]分別考慮了地震數(shù)據(jù)在時(shí)間上的小波狀和空間上的拋物狀特征,能夠?qū)ΟB前地震數(shù)據(jù)進(jìn)行稀疏表示。地震小波變換[10-11]考慮了地震同相軸斜率變化的光滑性,實(shí)現(xiàn)了一種針對(duì)地震數(shù)據(jù)的提升小波。非對(duì)稱(chēng)C-子波變換[12]考慮了地震子波的非對(duì)稱(chēng)性,能夠在對(duì)地震數(shù)據(jù)稀疏表示的同時(shí)提取出相關(guān)信息,如走時(shí)、衰減參數(shù)、非對(duì)稱(chēng)性參數(shù)等。

除了傳統(tǒng)意義上的向量的稀疏性,矩陣的低秩性(即奇異值的稀疏)也是一種廣義上的稀疏性,并且結(jié)合隨機(jī)采樣已經(jīng)應(yīng)用于地震數(shù)據(jù)重構(gòu)等[13]。此外,目前還出現(xiàn)了如拓展稀疏以及深度稀疏等新的稀疏性概念。拓展稀疏考慮了能量呈集中分布的情況;深度稀疏來(lái)源于深度學(xué)習(xí)技術(shù),即將稀疏基進(jìn)一步進(jìn)行稀疏表示,實(shí)現(xiàn)一個(gè)遞歸的稀疏表示過(guò)程。

1.2 隨機(jī)采樣

采樣過(guò)程的設(shè)計(jì)需要滿(mǎn)足一定的隨機(jī)性,同時(shí)滿(mǎn)足從高維到低維的采樣。隨機(jī)降維采樣在壓縮感知地震勘探如何體現(xiàn)呢？首先體現(xiàn)在觀測(cè)系統(tǒng)的設(shè)計(jì)上。我們可以采用不滿(mǎn)足奈奎斯特采樣定律的檢波器或者炮點(diǎn)分布來(lái)獲取采樣數(shù)據(jù),以此重構(gòu)完整數(shù)據(jù)。SHAHIDI等[14]引入了一種抖動(dòng)(jitter)隨機(jī)采樣方法,即檢波器隨機(jī)分布在規(guī)則網(wǎng)格所限定的區(qū)域內(nèi)部,可控制相鄰檢波器的最大間距,避免數(shù)據(jù)大片丟失。NAGHIZADEH等[15]總結(jié)了幾種采樣方法,包括規(guī)則采樣、隨機(jī)采樣以及隨機(jī)加規(guī)則采樣,并研究了不同采樣的頻譜特征。NAGHIZADEH[16]還提出了一種雙編織三維地震數(shù)據(jù)采樣方法,這種采樣方法可以配合傅里葉方法進(jìn)行數(shù)據(jù)重構(gòu)。HERRMANN[17]論述了利用隨機(jī)采樣甚至可從更少的數(shù)據(jù)中獲得更多的信息。國(guó)內(nèi)很多學(xué)者也做出了突出貢獻(xiàn),例如陳生昌等[18]提出了幾種隨機(jī)與均勻結(jié)合的高效采集方法,并在數(shù)值試驗(yàn)中取得較理想的效果。

除了觀測(cè)系統(tǒng)的下采樣設(shè)計(jì),還可以利用隨機(jī)同時(shí)震源實(shí)現(xiàn)下采樣。具體來(lái)講,傳統(tǒng)勘探中不同炮是按次序激發(fā),現(xiàn)在可以同時(shí)激發(fā),減少了采樣時(shí)間,但采集中不同炮點(diǎn)的數(shù)據(jù)會(huì)混疊在一起?？梢栽谕瑫r(shí)激發(fā)中添加隨機(jī)性來(lái)解決這類(lèi)問(wèn)題,比如炮點(diǎn)激發(fā)時(shí)間有隨機(jī)延遲(炮編碼[19]),另外可以采用隨機(jī)的激發(fā)波形(波場(chǎng)模擬過(guò)程中可以實(shí)現(xiàn)[20])。這樣在共炮點(diǎn)域具有相關(guān)性的數(shù)據(jù),轉(zhuǎn)換到共檢波點(diǎn)域就是不相關(guān)的,進(jìn)而可以進(jìn)行分離。已有石油公司提出了有效的壓縮感知采集系統(tǒng),并已在實(shí)際勘探中進(jìn)行了應(yīng)用測(cè)試[21-22]。

1.3 模擬、成像、混合震源、延時(shí)勘探

有限差分正演模擬的計(jì)算量取決于模型尺度以及觀測(cè)系統(tǒng)的分布,而不是取決于波場(chǎng)的復(fù)雜度。對(duì)于大尺度細(xì)密度的數(shù)據(jù)生成,即使波場(chǎng)非常簡(jiǎn)單,仍需要大量的計(jì)算資源。HERRMANN等[23]提出了基于壓縮感知的同時(shí)全波形正演模擬,通過(guò)少量的獨(dú)立震源同時(shí)激發(fā),得到混疊波場(chǎng),然后利用壓縮感知恢復(fù)算法得到分離后的波場(chǎng)。NEELAMANI等[20]考慮了波動(dòng)方程格林函數(shù)的稀疏性,在震源同時(shí)激發(fā)時(shí)采用隨機(jī)帶限脈沖波形,利用壓縮感知理論和曲波變換來(lái)恢復(fù)原始波場(chǎng)。

VERA RODRIGUEZ等[24]實(shí)現(xiàn)了基于壓縮感知的微地震偏移成像技術(shù)。LI等[25]利用壓縮感知和隨機(jī)優(yōu)化技術(shù)實(shí)現(xiàn)了地震數(shù)據(jù)的快速全波形反演。作者使用了兩種采樣策略:隨機(jī)震源編碼和隨機(jī)震源位置。LI等[26]和李翔[27]利用稀疏約束的Gauss-Newton方法進(jìn)行基于下采樣數(shù)據(jù)的全波形反演,減弱了下采樣所產(chǎn)生的人為噪聲,改善了反演結(jié)果。王華忠等[28]提出了基于特征波的成像方法,不是使用全波場(chǎng)而是使用壓縮感知框架提取出特征波場(chǎng),保持走時(shí)特征,然后進(jìn)行成像。

為了提高勘探效率,可以采用同時(shí)震源技術(shù)得到混疊的波場(chǎng),然后利用壓縮感知恢復(fù)出沒(méi)有混疊的波場(chǎng)。ABMA等[29]利用傅里葉變換實(shí)現(xiàn)了混疊波場(chǎng)的分離,然后用于偏移成像,得到了與直接使用無(wú)混疊波場(chǎng)相似或更好的結(jié)果。原因是使用波場(chǎng)混疊技術(shù)可以在相同時(shí)間內(nèi)獲取更多炮點(diǎn)位置的波場(chǎng)記錄。LIN等[19]利用復(fù)雜的震源編碼技術(shù)以及曲波變換將不同的混合震源數(shù)據(jù)分離開(kāi)來(lái)。

OGHENEKOHWO等[30]提出了基于壓縮感知的低成本延時(shí)地震勘探技術(shù)。不同時(shí)間的炮點(diǎn)布置在同一區(qū)域,但具體位置每次都不同,不同時(shí)間的勘探可以實(shí)現(xiàn)信息共享。另外,檢波器點(diǎn)是下采樣布置的,因此勘探成本較低。

2 討論與結(jié)論

CS如何在地震勘探開(kāi)發(fā)中發(fā)揮作用？首先是針對(duì)數(shù)據(jù)采集,在CS理論框架下,可以用比傳統(tǒng)方法更少的采集成本實(shí)現(xiàn)同等(甚至更高)精度的反演。這個(gè)“更少”到底是多少,它與前面提到的反演目標(biāo)的稀疏度有關(guān)。這個(gè)稀疏度(包括變換域的稀疏度),對(duì)一般的地震數(shù)據(jù)而言,是有一定的先驗(yàn)知識(shí);對(duì)實(shí)際數(shù)據(jù)的反演精度而言,只能有一個(gè)統(tǒng)計(jì)上的結(jié)論,即高概率重構(gòu)。此外,CS也可以深入挖掘現(xiàn)有勘探采集數(shù)據(jù)的信息,從而提高數(shù)據(jù)重構(gòu)(道加密)和偏移成像的反演精度。

如何才能用好CS？無(wú)非是在實(shí)際應(yīng)用中靈活把握CS的三要素:隨機(jī)采集(包括炮點(diǎn)和檢波器點(diǎn)兩方面的隨機(jī))、目標(biāo)的稀疏表達(dá)(明確具體問(wèn)題的目標(biāo)是高密數(shù)據(jù)還是地下介質(zhì)速度或其它物性參數(shù)？確定目標(biāo)后,再找到或設(shè)計(jì)一個(gè)稀疏變換,進(jìn)行目標(biāo)的稀疏表達(dá))和稀疏約束優(yōu)化重構(gòu)的快速算法。一般來(lái)說(shuō),重構(gòu)更高維的目標(biāo),需要用的采集數(shù)據(jù)(百分比)可更少。因?yàn)楦呔S的目標(biāo),有更多的先驗(yàn)信息可以利用。CS在五維數(shù)據(jù)的插值重構(gòu)比在二維或三維數(shù)據(jù)的重構(gòu)更能發(fā)揮出優(yōu)勢(shì),比如它有更多的斷層信息可利用,所以如何在高維目標(biāo)的稀疏表達(dá)中更好利用結(jié)構(gòu)的幾何信息是一個(gè)很關(guān)鍵的課題。

經(jīng)過(guò)近十年的發(fā)展,壓縮感知理論已開(kāi)始深入到地球物理勘探的方方面面,為降低勘探成本、提高采集效率、提高成像分辨率提供了新的理論支撐。相信在不久的將來(lái),壓縮感知能夠真正應(yīng)用于實(shí)際油氣勘探,加速同時(shí)震源、頁(yè)巖油、深層勘探等非常規(guī)勘探發(fā)展。

壓縮感知可以理解為一種淺層的學(xué)習(xí),結(jié)合當(dāng)前的深度學(xué)習(xí)技術(shù)是未來(lái)發(fā)展的一個(gè)方向。

致謝:感謝李幼銘、張捷、趙邦六、陳生昌、李成博等老師給予的巨大幫助,感謝課題組于四偉博士在此稿件的撰寫(xiě)過(guò)程中做出的很多輔助工作。

[1] DONOHO D.Compressed sensing [J].IEEE Transactions on Information Theory,2006,52(4):1289-1306

[3] MA J,PLONKA G.The curvelet transform[J].IEEE Signal Processing Magazine,2010,27(2):118-133

[5] EASLEY G,LABATE D,LIM W Q.Sparse directional image representations using the discrete shearlet transform[J].Applied and Computational Harmonic Analysis,2008,25(1):25-46

[6] 馮飛,王征,劉成明,等.基于Shearlet變換稀疏約束地震數(shù)據(jù)重建[J].石油物探,2016,55(5):682-691

FENG F,WANG Z,LIU C M,et al.Seismic data reconstruction based on sparse constraint in the shearlet domain[J].Geophysical Prospecting for Petroleum,2016,55(5):682-691

[7] BECKOUCHE S,MA J.Simultaneous dictionary learning and denoising for seismic data[J].Geophysics,2014,79(3):A27-A31

[8] YU S,MA J,ZHANG X,et al.Interpolation and denoising of high-dimensional seismic data by learning a tight frame[J].Geophysics,2015,80(5):V119-V132

[9] ZHANG R,ULRYCH T J.Physical wavelet frame denoising [J].Geophysics,2003,68(1):225-231

[10] FOMEL S,LIU Y.Seislet transform and seislet frame [J].Geophysics,2010,75(3):V25-V38

[11] 劉財(cái),李鵬,劉洋,等.基于seislet變換的反假頻迭代數(shù)據(jù)插值方法[J].地球物理學(xué)報(bào),2013,56(5):1619-1627

LIU C,LI P,LIU Y,et al.Iterative data interpolation beyond aliasing using seislet transform[J].Chinese Journal of Geophysics,2013,56(5):1619-1627

[12] BOβMANN F,MA J.Asymmetric chirplet transform for sparse representation of seismic data[J].Geophysics,2015,80(6):WD89-WD100

[13] YANG Y,MA J,OSHER S.Seismic data reconstruction via matrix completion[J].Inverse Problems and Imaging,2013,7(4):1379-1392

[14] SHAHIDI R,TANG G,MA J,et al.Application of randomized sampling schemes to curvelet-based sparsity-promoting seismic data recovery [J].Geophysical Prospecting,2013,61(5):973-997

[15] NAGHIZADEH M,SACCHI M D.On sampling functions and Fourier reconstruction methods [J].Geophysics,2010,75(6):WB137-WB151

[16] NAGHIZADEH M.Double-weave 3D seismic acquisition,part 1:sampling and sparse Fourier reconstruction[J].Geophysics,2015,80(6):WD143-WD162

[17] HERRMANN F.Randomized sampling and sparsity:getting more information from fewer samples[J].Geophysics,2010,75(6):WB173-WB187

[18] 陳生昌,陳國(guó)新,王漢闖.稀疏性約束的地球物理數(shù)據(jù)高效采集方法初步研究[J].石油物探,2015,54(1):24-35

CHEN S C,CHEN G X,WANG H C.The preliminary study on high efficient acquisition of geophysical data with sparsity constraints [J].Geophysical Prospecting for Petroleum,2015,54(1):24-35

[19] LIN T T Y,HERRMANN F J.Designing simultaneous acquisitions with compressive sensing[J].Expanded Abstracts of 71stEAGE Annual Conference,2009:4-9

[20] NEELAMANI R N,KROHN C E,KREBS J R,et al.Efficient seismic forward modeling using simultaneous random sources and sparsity [J].Geophysics,2010,75(6):WB15-WB27

[21] LI C,MOSHER C,MORLEY L,et al.Joint source deblending and reconstruction for seismic data [J].Expanded Abstracts of 83rdAnnual Internat SEG Mtg,2013:82-87

[22] MOSHER C,LI C,JANISZEWSKI F,et al.Operational deployment of compressive sensing systems for seismic data acquisition [J].The Leading Edge,2017,36(8):661-669

[23] HERRMANN F J,ERLANGGA Y A,Lin T T.Compressive simultaneous full-waveform simulation[J].Geophysics,2009,74(4):A35-A40

[24] VERA RODRIGUEZ I,KAZEMI N.Compressive sensing imaging of microseismic events constrained by the sign-bit[J].Geophysics,2016,81(1):KS1-KS10

[25] LI X,ARAVKIN A Y,VAN LEEUWEN T,et al.Fast randomized full-waveform inversion with compressive sensing[J].Geophysics,2012,77(3):A13-A17

[26] LI X,ESSER E,HERRMANN F J.Modified Gauss-Newton full-waveform inversion explained—Why sparsity-promoting updates do matter [J].Geophysics,2016,81(3):R99-R112

[27] 李翔.基于壓縮感知技術(shù)的全波形反演[J].石油物探,2017,56(1):20-25

LI X.Full-waveform inversion from compressively recovered updates [J].Geophysical Prospecting for Petroleum,2017,56(1):20-25

[28] 王華忠,馮波,王雄文,等.特征波反演成像理論框架[J].石油物探,2017,56(1):38-49

WANG H Z,FENG B,WANG X W,et al.The theoretical framework of characteristic wave inversion imaging [J].Geophysical Prospecting for Petroleum,2017,56(1):38-49

[29] ABMA R,HOWE D,FOSTER M,et al.Independent simultaneous source acquisition and processing [J].Geophysics,2015,80(6):WD37-WD44

[30] OGHENEKOHWO F,WASON H,ESSER E,et al.Low-cost time-lapse seismic with distributed compressive sensing,part 1:exploiting common information among the vintages [J].Geophysics,2017,82(3):P1-P13