自適應步長FISTA算法稀疏脈沖反褶積

潘樹林 閆 柯 李凌云 蔣從元 石林光

(①西南石油大學地球科學與技術學院，四川成都 610500；②中國石油西南油氣田分公司,四川成都 610051；③勝利油田分公司物探研究院，山東東營 257000；④四川職業技術學院電子電氣工程系,四川遂寧 629000)

0 問題的提出

稀疏脈沖反褶積是目前較為常用的提高地震資料分辨率的方法，其目的是利用帶限的地震記錄反演地下反射系數，從而得到高分辨率的地震數據[1-2]。稀疏脈沖反褶積的原理是在傳統的最小二乘反褶積的基礎上，加上一個稀疏先驗約束項，通過求解L2范數與L1范數之和形式的目標函數得到高分辨率的反演結果[3-4]。

基追蹤[5]方法利用Nesterov[6-8]提出的近端梯度的概念，使用二次函數近似表示目標函數，最終利用梯度法求解。該方法具有結構簡單、實現容易、計算效率高等特點，是主流求解方法。

FISTA算法是Beck等[9]提出的基追蹤方法，用于求解具有復合結構函數的最值

(1)

式中：F(x)為具有復合結構的目標函數，x為反射系數向量；f(x)為具有Lipschitz連續梯度的凸函數；φ(x)為在可行域Q上不可微的凸函數。根據Nesterov[6-8]提出的近端梯度法，式(1)對應的近端算子可以寫作

(2)

式中x,y∈R，y為x的近似值。由式(2)可以看出，Nesterov[6-8]采用二次函數近似式(1)的f(x)，通過求解式(2)迭代求解式(1)。式(2)的常數L即為內部梯度的步長，它決定了FISTA算法的收斂度。通常無法準確給出L，FISTA算法采用線性搜索方法尋找最佳L[10]。然而，線性搜索只能使L朝著增大的方向搜索，嚴重影響了FISTA算法的收斂性。線性搜索算法的主要缺點為：①當L的初始估計值L0大于實際值時，整個迭代過程的L偏大，無法快速收斂[11]；②當函數f(x)的局部曲率在初始迭代點附近較大，但在最優點附近減小時，如果L一直偏大，很可能導致“周波跳躍”而無法收斂。

針對上述問題前人進行了研究[12-16]，但效果不好。為此，本文提出了一種改進的自適應步長FISTA方法，在沒有外部參數調整的情況下仍能夠求解，并利用實驗驗證了方法的可靠性。

1 FISTA方法原理

時間域矩陣形式的稀疏脈沖反褶積的目標函數為[17]

(3)

式中：S為地震記錄向量；W為由地震子波組成的矩陣；λ為正則化參數變量。式(3)前半部分為最小二乘約束項，代表通過褶積模型生成的地震記錄與實際地震記錄的相似程度，后半部分為數據的稀疏先驗項。

根據前文所述，式(3)可以利用FISTA算法求解，其求解過程如圖1所示。由于采用線性搜索的方法，常數L只能沿著增大的方向變化，因此L的初始值決定了算法的收斂度，從而造成FISTA算法的不可控性。

圖1 FISTA計算流程圖

2 自適應步長FISTA方法原理

當常數L不變時，輔助序列tk+1與L無關，可以直接由下式確定[18]

(4)

為了能夠滿足理論上的收斂性，構造一個新的輔助序列，該輔助序列不僅只由tk和tk+1決定，還與前、后兩次迭代的常數Lk和Lk+1有關，即

(5)

由此，結合FISTA算法和式(5)可以得到一種自適應步長的改進算法(下文簡稱改進算法)，算法流程如圖2所示。改進算法在不同模型下均可取得滿意的計算效果。通過對比FISTA算法流程(圖1)和改進算法流程(圖2)可知，在每次迭代中，當L不滿足F[pLk(yk)]≤QLk[pLk(yk),yk]時： ①圖1通過Lk=nLk-1線性增大常數L，而在下一次迭代中直接使用Lk。②圖2則在每次迭代中，仍通過Lk=bLk-1線性增大常數L，但是在下一次迭代前，先通過Lk=aLk-1收縮L，然后利用收縮后的Lk進行后續迭代；與此同時，在每一次迭代過程中，利用前、后兩次迭代常數L的信息更新輔助序列tk，使其在理論上可收斂。

圖2 自適應步長FISTA方法流程圖

其中01。a為收縮因子，在每一次迭代前收縮常數L，b為放大因子，在每一次迭代中放大常數L。參數a，b的取值對計算結果影響較小，經實驗得a=0.7、b=2.0

3 理論模型分析

為了驗證改進算法的稀疏脈沖反褶積效果，首先合成一個理論模型(圖3)，分別利用FISTA算法和改進算法對理論模型進行稀疏脈沖反褶積。分兩種情況進行計算：①初始輸入常數L值過大，L=10L0(圖4)；②初始輸入常數L值過小，L=0.3L0(圖5)。整個過程中取a=0.7、b=2.0。

由圖4可見：當初始L過大時，改進算法反褶積效果優于傳統基于FISTA算法(圖4a)；當初始L過大時，傳統方法無法調整L的大小，改進算法能夠自適應地調整L(圖4b)；在相同迭代次數的情況下，由于常數L的調整使改進算法收斂度遠遠高于傳統方法(圖4c)。由圖5可見：當初始L過小時，改進算法與傳統方法計算過程大致相同，因此改進算法的計算效果與傳統算法的差異不大(圖5a)；改進算法在相同迭代次數下，收斂度略高于傳統算法(圖5b、圖5c)。

采用帶隨機噪聲的數據進一步測試算法效果。在理論合成地震記錄(圖3c)中加入隨機噪聲，得到了一個低信噪比理論模型(圖6)，對該資料分別使用過大和過小L值進行反演，獲得不同方法反褶積結果(圖7)?？梢?，基于傳統FISTA的稀疏脈沖反褶積方法在低信噪比情況下無法得到準確的反褶積結果，由改進算法能得到較準確的反褶積結果。上述分析進一步證明了改進方法具有更好的收斂性，能在不同信噪比下得到理想的反演結果，較常規FISTA算法具有更好的抗噪能力。

圖3 理論模型

圖4 初始L過大時反褶積結果

圖5 初始L過小時反褶積結果

圖6 含噪聲理論模型

圖7 不同方法反褶積結果(SNR=3.85dB)

4 實際數據處理

選取吐哈盆地雁木西北部地震資料測試改進算法的適用性[19-20]。圖8為雁木西T3井區Line 299測線地震剖面，可見地震分辨率較低，無法滿足精細解釋的需求。為此，分別采用基于傳統FISTA算法的稀疏脈沖反褶積以及改進算法對圖8數據進行高分辨率處理。由處理結果可見： 在相同的時間段內，改進算法反褶積結果(圖9b)的同相軸數目較FISTA算法(圖9a)多；經過基于FISTA算法的稀疏脈沖反褶積以及改進算法處理的地震資料的分辨率都得到提升，但是后者的頻帶更寬(圖10)、細節刻畫更清楚(圖11c)，利于高精度構造解釋。

圖8 雁木西T3井區Line 299測線地震剖面

圖9 反褶積剖面

圖10 圖9數據的頻譜

5 結束語

本文提出了一種基于自適應步長FISTA算法的稀疏脈沖反褶積方法，該方法在FISTA算法的基礎上，通過在每一次迭代之前適當減小常數L，然后利用線性搜索的方式尋找最優的常數L，以達到自適應調整L的目的。為了使算法達到理論收斂，通過結合前、后兩次的L，對傳統FISTA算法的輔助序列進行修改，最終使整套算法在理論上得以收斂。理論模型與實際地震資料的處理、分析結果表明，改進方法具有更好的收斂性，能在不同信噪比下得到理想的反演結果，較常規FISTA算法具有更強的抗噪能力。

圖11 圖9的局部放大