1、第一章 反演理論第一節 基本概念一 反演和正演1 反演反演是一個很廣的概念，根據地震波場、地球自由振蕩、交變電磁場、重力場以及熱學等地球物理觀測數據去推測地球內部的結構形態及物質成分，來定量計算各種有關的物理參數，這些都可以歸結為反演問題。在地震勘探中，反演的一個重要應用就是由地震記錄得到波阻抗。有反演，還有正演。要正確理解反演問題，還要知道正演的概念。2 正演正演和反演相反，它是對一個假設的地質模型，給定某些參數（如速度、層數、厚度）用理論關系式（數學模型）推導出某種可測量的量（如地震波）。在地震勘探中，正演的一個重要應用就是整理合成地震記錄。3 例子考慮地球內部的溫度分布，假定地球內部的溫

2、度隨深度線性增加，其關系式可表示成： T(z)=a+bz 正演：給定a和b，求不同深度z的對應溫度T(z) 反演：已經在不同點z測得T(z)，求a和b。二反演問題描述和公式表達的幾個重要問題1應用哪種參數化方式離散的還是連續的？2地球物理數據的性質是什么？觀測中的誤差是什么？3問題能不能作為數學問題提出，如果能夠，它是不是適定的？4對問題有無物理約束？5能獲得什么類型的解，達到什么精度？要求得到近似解、解的范圍、還是精確解？6問題是線性的還是非線性的？7問題是欠定的、超定的、還是適定的？8什么是問題的最好解法？9解的置信界限是什么？能否用其它方法來評價？第二節 反演的數學基礎一解超定線性反問題

3、 1簡單線性回歸 可利用最小平方法確定參數a、b使誤差的平方和最小。 （1-2-1） 擬合公式為： （1-2-2） 該方法的公式原來只適用于解超定問題，但同樣適用于欠定問題，當我們有多個參數時，稱為多元回歸，在地球物理領域廣泛采用這種方法。此過程用矩陣形式表示，則稱為廣義最小平方法矩陣方演。2非約束最小平方法反演廣義矩陣方法由前面討論可知，參數估計的最小平方方法用矩陣公式表示，所得到的算法等價于一個或多個模型參數的一個或多個數據集反演，步驟為：問題定義矩陣公式最小平方解 線性問題采用廣義矩陣形式 d=Gm （1-2-3） 對于精確的數據模型，參數m為 m=G-1d （1-2-4） 但是由于試驗

4、誤差，實際數據將不能精確擬合獲得，故采用最小平方法求解。解的矩陣表示式為 （1-2-5）上式具體計算時可用奇異值分解方法 G=UVT 最后，得 =（GTG）-1GTd=V-1UTd （1-2-6）二 約束線性最小平方反演為了得到最合適的解，通常可在方程d=Gm中加先驗信息，進行約束反演。約束方程為Dm=h （1-2-7） D一般為只有對角線有值的矩陣，我們希望朝著偏置使得最小。 =（d-Gmd-Gm）+2（Dm-hDm-h） （1-2-8） 如果D是單位矩陣，可以得到約束解 =（GTG+2I）（GTd+2h） （1-2-9） 式中，稱為Lagrange乘子。三解非線性反演問題1思路在實際工作中

5、許多問題都是非線性的，而非線性問題求解通常比較復雜，這樣就產生這樣一個問題，給定一些非線性問題，而它們又不服從簡單的線性變換，那么能否用通用的方法使我們可以用一些線性反演的方法來估算未知模型參數，并最終求得問題的解決呢？答案是肯定的。2初始模型和線性化對于非線性問題 di=fi（m1，m2，mp）=fi（m）， i=1，2，n （1-2-10）設m0為初始模型，則其響應為 （1-2-11）現假定f（m）在m0附近是線性的，從而關于m0的模型響應的微小攝動可以用Taylor級數展開為 或簡記為 實際情況要考慮噪聲 d=f（m）+e （1-2-12） 令y=d-f（m0），則有 e=d-=y-Ax

6、 （1-2-13） e=y-Ax這樣，非線性問題轉化成線性問題，我們可以用線性的方法求出問題的解。四、無約束非線性反演1問題的公式化目標函數： q=eTe=(d-f(m)T(d-f(m) （1-2-14） 利用前述結果，上式改寫為 q=eTe=(y-Ax)T(y-Ax) （1-2-15）2問題的解法：Gauss-Newton法對參數攝動的最小平方解 （1-2-16） 將攝動（x=m）應用于起始模型m0，迭代公式如下： (1-2-17) 其中mk為Jacobian矩陣A的賦值。3Gauss-Newton法的局限性當ATA病態（本征值很小或近于0）時，計算的解會大到令人難以置信。因此在實踐當中，必

7、須對mk做x的微小校正。4最速下降（梯度）法初始模型僅在目標函數q的負梯度方向予以校正，即 （1-2-18） 其中k是合適的常數，進一步推導可得 （1-2-19） 以上方程中以ATA-1取代常數因子2k，將變為方程1-2-16所定義的Gauss-Newton法，k值決定校正步長。但以上方程并不含有任何逆矩陣，因此較Gauss-Newton法具備更好的起始收斂特征。最速下降法當采用最小平方解法時，其收斂速率將下降，因此不宜在實際反演中應用。5對非穩定性和非收斂性的補救辦法當ATA是病態時，為防止無界解的增大，Levenberg（1944）提出了一種阻尼最小平方的方法，該方法可在Taylor近似的

8、逐次應用過程中，阻滯參數攝動的絕對值。Levenberg建議應在ATA的主對角線上加一個隨意選取的正的權因子，并且要顯示出當權因子相等時，q2的剩余和的方向導數為最小。這種想法以后為Maequardt（1963，1970）用來開發了一種非常有用的非線性算法。該技術稱為嶺回歸（Ridge Regression）或Marquardt-Levenberg方法，是地球物理領域最常見的一種反演算法。就其本質來講，實際上是Gauss-Newton法和最速下降法之間的內插，一種成功地結合二者有用特性的混合技術。五、約束反演：嶺回歸或Marquardt-Levenberg法1目標函數 （1-2-20） 目的：

9、誤差和攝動量均取極小。其中攝動量是新增的約束條件，從本質上講，嶺回歸法實際上是約束非線性最小平方法。是Lagrange乘子，可認為是阻尼因子。如果賦值近于0，則其解近似于Gauss-Newton解。2問題的求解求解方法與非約束最小平方法相同，最終的解為： （1-2-21） 而后可將解xr用于迭代過程 （1-2-22） 其中A是k+1次迭代對mk求的值 （1-2-23）嶺回歸法實際上是最速下降法和Gauss-Newton法二者相結合的混合技術。當初始模型與問題的解相差甚遠時，最速下降法起主要作用；而當接近于最終解時，最小平方法起主要作用。六非線性偏置估計對一組既不完整又不準確的數據進行解釋時，通

10、常比較明智的做法是尋找一個和先驗數據相一致的模型，這些先驗數據可以是先前的地球物理研究數據，地質數據、測井數據，這些附加的先驗信息可以幫助我們從不準確的實際數據得出的所有的解中求出最可信的一個，附有先驗信息的反演問題可在一個統一的偏置估計框架內進行討論。此方法強調實際過程的簡單有效，為清楚起見，在此種方法中將初始模型和先驗信息加以區別。1理論基礎偏置估計的理論很簡單，其基本原理類似于約束線性最小平方反演方法。特別的是除起始（或初始）模型m0外引入了先驗信息h。同時，用對角線加權矩陣W=-1I來比例數據方程，使求解過程穩定。2應用先驗信息的非線性反演為設有p個參數，h為先驗數據，Dm=h形式的約

11、束方程可表示為 （1-2-24） 為使相鄰物理參數之間的差異降至最小平滑度，需采取TwoneyTikhonoy平滑度措施。 （1-2-25）我們的目的是要使m偏向于h，不妨將問題簡單陳述為：給定一組有限的不準確的觀測數據，在所有等效解中求其真解（考慮數據和模型誤差）并使之與觀測數據相吻合，且滿足模型參數的可靠估計。從數學意義來講，上述問題就等效于對預測誤差eTe和最終解與特定約束的偏差極小 （1-2-26） 如果f(m)是連續的并且可微，則可用Taylor定理將其相對于初始模型m0 展開，從而給出方程(1-2-26)的線性近似（1-2-27） 令B=T，展開上式，并將偏微分置0，最后得偏置解為

12、 (1-2-28) 迭代公式 （1-2-29）如果先驗信息有疑義（或不可信），那就需要將約束置為，即h=0，0，0T，而且所有的元素均置為相等的常數（0），這樣所有的參數都具有相等的權重。在這種情況下，可以方便地由一單值未定乘子所取代。這樣就有參數校正解 （1-2-30） 其迭代公式 （1-2-31） 因為D=1，這里2I用以控制求解的步長，而2mk有助于減小其向零矢量h的位置，我們可以將這種方法稱為平滑度約束反演或最小偏置算法。3與標準方法的關系在偏置估計中，如果0，那么上述所有迭代估計公式均會簡化為改進了的加權經典最小平方化式 （1-2-32） 偏置估計方法的穩定性和有效性主要取決于和D。

13、方程（1-2-30）與通常的阻尼最小平方或嶺回歸優化公式 （1-2-33）的不同之處在于一2m0項，唯一目的是要對參數增量的變化范圍置一邊界。我們可將約束反演問題定義為求反Lagrange函數的極值問題： （1-2-34） 要搜尋的是最佳擬合數據的起始模型的有界攝動。在方程（1-2-29）中以量E取代WTW，我們有： （1-2-35）如果B可以統計地解釋為先驗參數協方差矩陣的逆，則上述方程即等效于Jackson和Matsuura的Bayes估計方法，并類似于Tarantola和Valette的非線性算法。因此，應用簡單代數，我們事實上已經導出一種與基于具有先驗數據的概率統計處理的數學上比較嚴謹

14、的非線性反演法相類似的方法，但是應該注意到，Tarantola和Valette的里程碑方法中的反演理論和先驗信息的使用均與我們的方法不盡相同。我們的主要興趣在于迫使最終解盡可能與那些先驗參數估計相一致，因此方程（1-2-35）右端最后一項不為0，因為在實際情況下，已知的先驗參數估計很少。在Tarantola和Valette算法中，h即是實際起始模型m0。在這種情況下，正如Pous等人指出的那樣，方程（1-2-35）的最后一項在第一次迭代中應為0。我們將h和m0分開，這點和Jackson和Matsuura的作法大致相當。但在Jackson和Matsuura的算法中，方程（1-2-35）括號中的量

15、乘了一個適當選擇的因子b(0<b<1)。如此說來，我們的簡單方法或許更為通用。第三節 地震勘探中的反演方法一 地震反演的分類地震反演通常分成疊前和疊后反演兩大類：疊前反演應用較少，較成熟的是AVO反演；疊后反演的大量應用是波阻抗反演，這是當前地震資料處理的重要結果。從反演方法上可將地震反演劃分為基于波動理論的波動方程反演和基于Robinson褶積模型反演兩大類，在實際工作中主要是基于Robinson模型的反演。我們通常所說的地震資料波阻抗反演指的是基于Robinson褶積模型的疊后地震資料反演，目前常見的有遞推反演（Recursive inversion），稀疏脈沖反演(Spars

16、e-spike inversion)和基于模型的反演(Model-based inversion)，后面兩種反演方法通常稱為寬帶反演。遞推反演包括道積分、GLOG、VLOG、SEISLOG、塊狀反演、帶限反演和PIVT等；稀疏脈沖反演包括多種實現方法，如模方法、最小熵方法、最大似然方法等；基于模型的反演也包括廣義線性反演（GLI）（Cooke，1983）、地震巖性模擬（SLIM）（Gelfand，1984）、魯棒的速度反演方法（ROVIM）（Fabre，1989）、寬帶約束反演（BCI）（Martinez，1988）、PARM和Jason等。二. 遞推反演方法1 波阻抗遞推公式對于兩層介質，反

17、射系數為： 分別為上下界面的密度，V1、V2為上下界面的速度。當地下為多層水平介質時，任意第i個界面的反射系數為： （1-3-1） 對應的波阻抗為： （1-3-2） 遞推公式： （1-3-3） 如果用經過特殊處理的地震剖面記錄道Xn代替反射系數Rn，則上式可寫成 （1-3-4）這就是遞推反演的基本公式。2低頻補償地震反射系數剖面的頻帶是有限的，它缺失的是高頻成分和低頻成分，對波阻抗反演而言，缺失高頻成分只影響分辨率，而缺失低頻成分就失去了速度的直流分量及速度斜坡的信息。這種波阻抗剖面通常稱為相對波阻抗剖面或剩余波阻抗剖面，剩余波阻抗剖面只能反映波阻抗的相對變化而不能反映波阻抗的真實情況，因此必

18、須在剩余波阻抗剖面上，再加上合理的低頻成分，進行低頻補償。（1） 利用聲波測井資料補償低頻這是最常用的方法。（2） 利用疊加速度補償低頻 但疊加速度補償因為實際三維速度場精度有限，會出現低頻缺口，造成聲測井曲線的值偏高和偏低的振蕩。（3） 利用地質模型補償低頻這種方法比較費時。低頻缺口在波阻抗反演中是常見的，有時也是比較嚴重的問題。所幸其橫向上速度相對變化通常是正確的，仍然能確定目的層段上有意義的巖性變化。3. 一個簡單應用道積分 圖1-1 道積分剖面 該方法不做低頻補償，得到的是相對波阻抗。 用連續時間函數表示（1-3-4）式 （1-3-5）如果經反褶積處理后的地震道x（t）的脈沖寬度足夠小

19、，認為x（t）與反射函數r（t）成比例 r（t）x（t） 則可近似求出任一時刻t的近似波阻抗 （1-3-6） K為比例系數。 具體實現步驟為： 將地震記錄振幅標定到反射系數數量級 計算積分道 將積分結果轉換為波阻抗 對轉換結果作帶通濾波得地層相對波阻抗 圖1-1給出了具體的應用例子，處理資料為TJH三維工區一剖面。三 稀疏脈沖反演方法 這種方法假設地下反射系數序列是由一系列大的反射系數疊加在服從高斯分布的小反射系數背景上構成的，主要有：L1范數反褶積、最小熵方法、最大似然方法等。L1范數反褶積最早由Barrodale于1973年提出，后經Taylor1979年及Oldenburg1983年的研

20、究，改進成為一種獨特的反褶積方法，它的特點是對子波的各種相位特性都有較好的適應性。常規脈沖反褶積及預測反褶積都要假定子波是最小相位的，并且反射系數是白噪。在這兩個條件下，反褶積的求解運算工作只有在最小二乘的意義下（與期輸出波形均方根誤差最小），才能得到一組Toeplitz方程組，才能用萊文森遞推法快速求解反褶積因子。誤差的最小平方就被稱為其范數為2。而L1范數是不做平方的判斷，而用誤差的絕對值之和作為標準，故稱其范數為L1范數。1985年王承曙等又提出Lp范數反褶積，即其判斷的范數可以不是2，也不是1，而是一個任意的正整數p。由于采用了L1范數，帶來的好處是對子波的相位特性放寬了限制，但是在計

21、算中沒有脈沖反褶積那樣簡單了。它一般是從線性規劃的理論出發，求解一組超定方程組的最優解。在求解過程中必須反復迭代，或者化為一組非線性方程組，用非線性規劃方法迭代求解。最小熵方法由R.A.Wiggins首先提出，它以方差模為判斷準則，信號的規則性達到最大，熵為最小。Wiggins的方差模的定義如下： (1-3-7)式中是地震道數據在某個時窗內的第個數值。當波形很突出時，達到很大振幅值，于是達到極大值，此時認為效果達到最佳。此方法對子波的相位特性不做約束，而且在一定程度上可以把混合相位子波向零相位靠攏。但該方法假設反射系數是稀疏的，只有當有少數大的脈沖存在時效果才很好，所以在具有亮點強波的剖面上，往往得到較好的反演效果。稀疏脈沖反演方法的輸出為矩形波阻抗曲線形式，地層邊界清晰，對厚層碳酸鹽巖地區較為合適。然而其致命的弱點是要求反射系數是稀疏的，而實際上大多數地震道的反射系數是稠密的。四 基于模型的反演1 流程框圖模型為基礎的方法，或簡稱模型法，首先構造一個地質模型，并將其與