反演原理及公式介紹

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

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

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

4、誤差，實際數(shù)據(jù)將不能精確擬合獲得，故采用最小平方法求解。解的矩陣表示式為 （1-2-5）上式具體計算時可用奇異值分解方法 g=uvt 最后，得 =（gtg）-1gtd=v-1utd （1-2-6）二 約束線性最小平方反演為了得到最合適的解，通?？稍诜匠蘢=gm中加先驗信息，進(jìn)行約束反演。約束方程為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、許多問題都是非線性的，而非線性問題求解通常比較復(fù)雜，這樣就產(chǎn)生這樣一個問題，給定一些非線性問題，而它們又不服從簡單的線性變換，那么能否用通用的方法使我們可以用一些線性反演的方法來估算未知模型參數(shù)，并最終求得問題的解決呢？答案是肯定的。2初始模型和線性化對于非線性問題 di=fi（m1，m2，mp）=fi（m）， i=1，2，n （1-2-10）設(shè)m0為初始模型，則其響應(yīng)為 （1-2-11）現(xiàn)假定f（m）在m0附近是線性的，從而關(guān)于m0的模型響應(yīng)的微小攝動可以用taylor級數(shù)展開為 或簡記為 實際情況要考慮噪聲 d=f（m）+e （1-2-12） 令y=d-f（m0），則有 e=d-=y-ax

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

7、須對mk做x的微小校正。4最速下降（梯度）法初始模型僅在目標(biāo)函數(shù)q的負(fù)梯度方向予以校正，即 （1-2-18） 其中k是合適的常數(shù)，進(jìn)一步推導(dǎo)可得 （1-2-19） 以上方程中以ata-1取代常數(shù)因子2k，將變?yōu)榉匠?-2-16所定義的gauss-newton法，k值決定校正步長。但以上方程并不含有任何逆矩陣，因此較gauss-newton法具備更好的起始收斂特征。最速下降法當(dāng)采用最小平方解法時，其收斂速率將下降，因此不宜在實際反演中應(yīng)用。5對非穩(wěn)定性和非收斂性的補救辦法當(dāng)ata是病態(tài)時，為防止無界解的增大，levenberg（1944）提出了一種阻尼最小平方的方法，該方法可在taylor近似的

8、逐次應(yīng)用過程中，阻滯參數(shù)攝動的絕對值。levenberg建議應(yīng)在ata的主對角線上加一個隨意選取的正的權(quán)因子，并且要顯示出當(dāng)權(quán)因子相等時，q2的剩余和的方向?qū)?shù)為最小。這種想法以后為maequardt（1963，1970）用來開發(fā)了一種非常有用的非線性算法。該技術(shù)稱為嶺回歸（ridge regression）或marquardt-levenberg方法，是地球物理領(lǐng)域最常見的一種反演算法。就其本質(zhì)來講，實際上是gauss-newton法和最速下降法之間的內(nèi)插，一種成功地結(jié)合二者有用特性的混合技術(shù)。五、約束反演：嶺回歸或marquardt-levenberg法1目標(biāo)函數(shù) （1-2-20） 目的：

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

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

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

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

13、方程（1-2-30）與通常的阻尼最小平方或嶺回歸優(yōu)化公式 （1-2-33）的不同之處在于一2m0項，唯一目的是要對參數(shù)增量的變化范圍置一邊界。我們可將約束反演問題定義為求反lagrange函數(shù)的極值問題： （1-2-34） 要搜尋的是最佳擬合數(shù)據(jù)的起始模型的有界攝動。在方程（1-2-29）中以量e取代wtw，我們有： （1-2-35）如果b可以統(tǒng)計地解釋為先驗參數(shù)協(xié)方差矩陣的逆，則上述方程即等效于jackson和matsuura的bayes估計方法，并類似于tarantola和valette的非線性算法。因此，應(yīng)用簡單代數(shù)，我們事實上已經(jīng)導(dǎo)出一種與基于具有先驗數(shù)據(jù)的概率統(tǒng)計處理的數(shù)學(xué)上比較嚴(yán)謹(jǐn)

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

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

16、e-spike inversion)和基于模型的反演(model-based inversion)，后面兩種反演方法通常稱為寬帶反演。遞推反演包括道積分、glog、vlog、seislog、塊狀反演、帶限反演和pivt等；稀疏脈沖反演包括多種實現(xiàn)方法，如模方法、最小熵方法、最大似然方法等；基于模型的反演也包括廣義線性反演（gli）（cooke，1983）、地震巖性模擬（slim）（gelfand，1984）、魯棒的速度反演方法（rovim）（fabre，1989）、寬帶約束反演（bci）（martinez，1988）、parm和jason等。二. 遞推反演方法1 波阻抗遞推公式對于兩層介質(zhì)，反

17、射系數(shù)為： 分別為上下界面的密度，v1、v2為上下界面的速度。當(dāng)?shù)叵聻槎鄬铀浇橘|(zhì)時，任意第i個界面的反射系數(shù)為： （1-3-1） 對應(yīng)的波阻抗為： （1-3-2） 遞推公式： （1-3-3） 如果用經(jīng)過特殊處理的地震剖面記錄道xn代替反射系數(shù)rn，則上式可寫成 （1-3-4）這就是遞推反演的基本公式。2低頻補償?shù)卣鸱瓷湎禂?shù)剖面的頻帶是有限的，它缺失的是高頻成分和低頻成分，對波阻抗反演而言，缺失高頻成分只影響分辨率，而缺失低頻成分就失去了速度的直流分量及速度斜坡的信息。這種波阻抗剖面通常稱為相對波阻抗剖面或剩余波阻抗剖面，剩余波阻抗剖面只能反映波阻抗的相對變化而不能反映波阻抗的真實情況，因此必

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

19、，認(rèn)為x（t）與反射函數(shù)r（t）成比例 r（t）x（t） 則可近似求出任一時刻t的近似波阻抗 （1-3-6） k為比例系數(shù)。 具體實現(xiàn)步驟為： 將地震記錄振幅標(biāo)定到反射系數(shù)數(shù)量級 計算積分道 將積分結(jié)果轉(zhuǎn)換為波阻抗 對轉(zhuǎn)換結(jié)果作帶通濾波得地層相對波阻抗 圖1-1給出了具體的應(yīng)用例子，處理資料為tjh三維工區(qū)一剖面。三 稀疏脈沖反演方法 這種方法假設(shè)地下反射系數(shù)序列是由一系列大的反射系數(shù)疊加在服從高斯分布的小反射系數(shù)背景上構(gòu)成的，主要有：l1范數(shù)反褶積、最小熵方法、最大似然方法等。l1范數(shù)反褶積最早由barrodale于1973年提出，后經(jīng)taylor1979年及oldenburg1983年的研

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

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