1、空間變量結構約束下AVO貝葉斯反演摘要引入來源不同的地下先驗信息約束地球物理反演，降低反演解的不確定性程度是地球物理反演研究中一個重要的課題。地質變量是在空間(時間)上展布的變量。地質變量除了具有統計意義上的隨機性外，還依賴于空間(時間)位置,具有一定的相關性和結構性。本文引入地質統計學里變差函數的思想，建立關于地下空間變量結構的描述作為先驗信息，結合貝葉斯推斷理論，通過似然函數建立空間變量和地球物理觀測的聯系，得到AVO反演解估計后驗概率密度的解析形式。選取東海某油田一口實際測井數據進行反演試算應用區間估計和隨機模擬的方法對解估計的不確定性程度進行評價，結果表明在空間變量結構約束下，解估計的

2、不確定性程度低于未受到約束的反演解，驗證了此文方法的有效性。關鍵詞地質統計學，AVO反演，貝葉斯理論，空間變量結構AVOBayesianinversionconstrainedbyspatialvariablestructureAbstractAnimportantsubjectingeophysicalinversionresearchistodecreasetheuncertaintydegreeofinversesolutionbyintroducingsubsurfaceprioriinformationfromdifferentsourcesasaconstraintofgeophy

3、sicalinversion.Geologicalvariableexhibitsspatial(time)variability.Apartfromtherandomnessintermsofstatistics,geologicalvariabledependsonspatial(time)locationaswell,possessingcertaincorrelationandstructuralproperty.Thepaperintroducesthetheoryofvariogramfromgeostaticstocreateadescriptionaboutthestructu

4、reofsubsurfacespatialvariableasprioriinformation.CombinedwithBayesianinferencetheory,weestablisharelationshipbetweenspatialvariableandgeophysicalobservationvialikelihoodfunction.Finally,wederivetheanalyticalformofAVOinversesolutionestimationsposteriorprobabilitydensity.Practicalloggingdatafromawelll

5、ocatedinanoilfieldinEastChinaSeaisselectedtoconductpilotinversecalculation.Afterwards,uncertaintydegreeofitssolutionestimationisevaluatedbyperformingintervalestimationandstochasticsimulationmethods.Theconsequencesuggeststhatwithspatialvariablestructureasconstraint,theuncertaintydegreeofsolutionestim

6、ationislowerthanitofinversesolutionwhichiswithoutconstraint.Thus,methodproposedinthispaperisprovedtobeeffective.KeywordsGeostatistics,AVOinversion,Bayesiantheory,Spatialvariablestructure0引言眾所周知，地球物理反問題是一個多解且不適定的問題。對于地震反演來說，其根本目標是利用觀測到的地震數據集，求解地下介質的結構與物理參數1。除了地震數據外，如何讓引入更多不同來源的信息約束地震反演，降低反問題的多解性，獲得更真

7、實描述地下結構與性質的解一直是地震反演的關鍵問題之一。由于反問題具有多解性，所以在反問題研究中，重要的不是構造出一種廣義的解估計，而是要對各種可能的解估計進行評價，這正是由于反問題解的非唯一性所決定的。貝葉斯方法是一個解決上述問題很好的工具，Tarantola(1987,2005)使用關于模型空間的先驗信息和似然函數構造了反問題解的后驗概率密度，完整的描述了解估計的不確定性問題；2001年，Scales和Tenorio關于如何理解在反問題研究中的先驗信息與不確定性關系，在其論文中給出了精彩的闡述肌Buland(2003,2006)等人應用貝葉斯原理在子波估計、AVO線性反演、時移地震反演方面給

8、出了研究方法和應用實例。2009年，陳建江等人專門論述了AVO反演的不確定性分析a】。在統計學概念下，將反演的地質變量認為是隨機變量，一次觀測(如測井數據)均為隨機變量的一次實現。經典統計學中，通常假定某隨機變量采集的樣本是完全隨機的或在空間上是完全獨立的。但通常地，需要通過反演求得的地質變量在時空域上并不完全獨立，而是呈現出一定的空間(時間)相關性。19世紀60年代，法國學者G.Matheron發展起來的地質統計學，考慮到空間變量的結構性，以變差函數、克里金估值、隨機模擬為基本手段，是研究空間中某些既顯示出隨機性又顯示出結構性的變量的數學地質方法11。1994年，Haas和Dubrule結合

9、地質統計學中的隨機模擬思想以及地震數據的反演，提出了現在所謂的地質統計學反演12。隨機反演(地質統計學反演)要求對地下模型空間產生多個實現，比較常用的方法是MonteCarlo等計算代價高昂的方法13，使得在實際應用中受到限制。本文在貝葉斯理論框架下，假設地下介質參數服從對數高斯分布，引入地質統計學變差函數的概念，以地質變量的空間分布結構作為先驗對AVO反演進行約束。通過先驗的地質變量空間分布和描述地球物理觀測和地質變量關系的似然函數，建立了反演解的后驗概率密度分布的解析形式。應用隨機模擬和區間估計的方法，來定性和定量地評價反演解的不確定性程度。數值算例表明，相比于常規的貝葉斯反演方法，考慮地

10、質變量在空間上的結構性來約束反演，更符合實際情況，不確定性程度降低。1方法原理在一般貝葉斯理論框架下，地球物理反演問題的解以及解的不確定性，以后驗概率密度的形式來描述：b(mId)L(dIm)p(m)(1)JL(dIm)p(m)dmp(m)為關于模型m分布的先驗概率密度；L(dIm)為似然函數，是模型空間m向數據空間d的映射的概率表達；分母JL(dIm)p(m)dm在給定數據情況下為一個常數，也即是概率歸一化因子；(mId)是所要求的后驗概率密度分布。方便下文討論，簡寫L(dIm)、1)式簡寫為：似然函數L(m)的建立，以(mId)為L(m)、b(m)，忽略掉常數JL(dIm)p(m)dm,b

11、(m)二L(m)p(m)(2)接下來，本節將分三個部分分別討論模型先驗分布p(m)、及后驗概率密度分布的求取。1.1模型先驗概率分布在各向同性介質中，地下介質的彈性性質可由縱波速度、橫波速度、密度來完整地描述,同樣地也可以由縱波阻抗、橫波阻抗、密度來描述。在此本文選擇關于地下介質的彈性參數:縱波速度Vp、橫波速度Vs、密度P來作為反演目標參數。關于模型的先驗分布在貝葉斯理論中所起的角色，在理論上一直是概率統計學家爭論的焦點14。盡管理論上關于模型先驗信息的理解還存在爭論，但在地球物理實際應用中，在已知地質背景規律或鉆井、測井資料指導下反演，已經是大勢所趨。關于貝葉斯AVO反演，本文以測井數據作

12、為地下模型的先驗信息。但不同于常規的測井約束反演，本文不僅從測井數據里提取關于地下介質的低頻信息來約束反演，還進一步挖掘了測井數據里包含的地下介質彈性參數的空間相關性(也即地下空間變量的結構)來進一步約束指導反演。地質統計學里的變差函數是描述空間變量結構的天然數學工具，可以證明在二階平穩條件下，變差函數與協方差函數等價，在描述空間變量結構方面下文將不再區分協方差函數和變差函數。關于地下空間變量的分布，大量井資料分析表明，彈性介質三參數(縱波速度、橫波速度、密度)服從對數高斯分布形式。由以上討論，本文定義關于空間變量結構的先驗分布為：TOC o 1-5 h z HYPERLINK l bookm

13、ark22 o Current Document m(t)二lnVp(t),lnVs(t),lnp(t)tNg,E)(3)mm其中t為地震雙程旅行時，g為模型期望，E為模型協方差矩陣：mm HYPERLINK l bookmark24 o Current Document g(t)=Em(t)=卩(t),卩(t),卩(t)t(4)mVpVsE二Covm(t),m(t)二E(t,t)(5)mijmij在二階平穩假設下(5)式可簡化為:E二Co譏m(t),m(t)二E(t)(6)mijm其中t=|ttI，12(C(t)Vp,VpC(t)Vs,VpIC(t)VpC(t)Vp,VsC(t)Vs,VsC

14、(t),VsC(t)Vp,PC(t)Vs,C(t)丿,7)關于先驗模型均值gm、協方差函數Em矩陣元的計算，是通過測井數據數據的統計分析獲得。gm反映了地質變量的低頻信息，Em則分反映了地質變量的空間結構信息。由此建立了空間變量的先驗信息，可以用來約束反演。注意的是(7)式先驗協方差函數的計算稍顯復雜，但是它是描述空間變量信息的載體,是進行空間變量結構約束的關鍵。為了使本文的方法理論簡潔明晰，在這里不討論(7)式的具體計算，其詳細的計算方式與討論見反演算例部分。1.2似然函數在概率統計上的似然函數描述的是估計的模型響應與真實數據的相似程度，由模型向數據的映射，也就是我們稱之為的動力系統或正演方

15、程。在一定的假設下，可以建立地下彈性參數映射到地震數據的線性關系。1919年Zoeppritz建立了平面波入射時反射系數、透射系數的理論表達式15，但其數學表達式復雜、物理意義不直觀，本文采用1980年由Aki和Richard提出來的弱縱波反射系數近似式16，近似式表達形式與Stolt與Welein(1985)所采用的一樣：dInVp(t)dInVs(t)dInp(t)TOC o 1-5 h zR(t,0)=a(t,0)+a(t,0)+a(t,0)(8)vpdtVsdtPdt其中：1a(t,0)=(1+tan20)(9)Vp2a(t,0)=-4Vs(t)sin20(10)VsVp(t)2a(t

16、,0)-1(1-4Vs(t尸sin20)(11)p2Vp(t)2這里Vp(t)、Vs(t)可以認為是彈性參數的低頻背景值，可以通過模型先驗均值A獲得，則m系數項a、a、a可以先驗求得。(8)式對應的矩陣表達式為：VpVsPR-Am(12)系數矩陣A-a(t,0),a(t,0),a(t,0)，m是m關于時間的一階導數。VpVsP則地震記錄d表示為子波S與反射系數R的褶積加上噪音項e:d-SR+e-SAm+e(13)假設噪音服從零均值高斯分布：eN(0,E)(14)e將褶積運算SA記作算子G，由褶積元算的線性性質可知G為線性算子。由此建立了線性正演方程，將(13)式寫為：d-SAm+e-Gm+e(

17、15)1.3后驗概率分布由模型先驗概率分布與似然函數，在線性高斯條件下可得后驗概率分布為高斯分布4：TOC o 1-5 h zmIdN(p,E)(16)m|dm|d其中：g=p+EGt(GEGt+E)-i(d-Gp)(17)m|dmmmemE=E-EGt(GEGt+E)-iGE(18)m|dmmmem6力2其中，二學E(t)，=-E(t)。m6Tmm6T2m式(17)完整地描述了反演解的概率性質，包括了反演解的期望與不確定性。由于后驗概率分布是高斯分布，反問題的極大后驗概率解也就是后驗概率的條件期望MAPm=g。這里取得的反演解是真實彈性參數的對數ln(Vp)、ln(Vs)、ln(Den)，m

18、id必須經過指數變換才可得到Vp、Vs、Den。2數值反演算例為了驗證本文方法的可行性，選擇東海某油田一口包含縱波速度、橫波速度、密度測井數據的勘探井進行試驗。測井數據經過標定，深度域變換到時間域(雙程旅行時)，井數據選取時窗為3500ms3700ms，采樣率1ms。井曲線及其選用的低頻背景見圖1。圖1縱波速度、橫波速度、密度(紅色實線)；模型低頻背景(黑色虛線)Fig.1P-wavevelocity,S-wavevelocity,density(redsolidline);Modellow-frequencybackground(blackdotline)由井曲線通過Aki-Richard縱

19、波近似反射系數公式計算得到的反射系數與主頻50Hz雷克子波褶積獲得角道集地震記錄，圖2所示。記錄時長與井數據相同3500ms3700ms,采樣率1ms。計算得到的合成記錄作為反演時地震數據的輸入。3500352535503575600SE姑36503675370037253750合成地震記錄5101520253035Angle圖2部分角度合成地震記錄Fig.2Partialanglesyntheticseismogram對測井數據Vp、Vs、Den取對數，采用正態概率紙來檢驗對數彈性參數是否符合高斯分布。由圖3可得，認為認為測井數據的對數彈性參數的高斯分布假設是可以接受的。關于模型先驗信息，模

20、型均值取圖1所示低頻背景的對數，二階平穩假設下認為模型參數在每個采樣點上的對數方差為常數b2（t）=0.0162，b2（t）二0.0117，b2（t）=0.000403。關于模VpVsp型參數的空間結構，也即先驗協方差函數計算見下一部分。99:99979-H-5999o.n.IJ._u_uoonun-9075502597985o999999o.go.gooDuo552o.o0.100.05.020.010.0030.001ln(Vp)(m/s)0.100.05.20.010.003.1ln(Vs)(m/s)0.003.0010.100.05.20.017.87.847.88ln(Den)(Kg

21、/ms)97985nu999999.93_UO._uon-5_u760_u52o圖3對數彈性參數在正態概率紙上的分布Fig.3Distributionoflogarithmelasticparametersonnormalprobabilitypaper2.1先驗空間變量結構計算地質統計學里，關于空間變量的結構是利用變差函數工具來進行描述的，變差函數反映了變量自身在空間上的結構性以及變量之間在空間上的依賴關系I。在本文中的空間變量就是需要反演的彈性參數，縱波速度Vp、橫波速度Vs、密度Den。一般地，在應用地震數據反演前，我們已經通過鉆井獲得了關于井眼處地下介質的一些信息（如巖性測井、物性測井

22、、電性測井）。通過井數據的統計分析，可以計算地下介質彈性參數的變差函數。將井上計算到的變差函數引入到反演中來進行空間變量結構的約束，更符合地下地質客觀實際。由于變差函數需要滿足正定性條件，通常的方式是運用一些已知的、滿足正定條件的變差函數的線性組合來建立變差函數模型。一些被證明有效的、常用的、基本的變差函數包括有球形模型、指數模型、高斯模型、DeWijsian模型等，具體的表達式及應用見參考文獻11。圖4是三種常用的變差函數模型的形態曲線。圖4三種常用變差函數模型的形態曲線Fig.4Curvesofthreecommonvariogramfunctionmodel型型型狀數斯球指高可以證明，在

23、二階平穩下，變差函數與協方差函數之和為一常數，所以關于(7)式協方差矩陣元的計算可以通過變差函數的概念來求取。在二階平穩假設下，變差函數與協方差函數在反應空間結構性方面具有等價的意義，在本文接下來，為了與上述反演理論一致，將使用協方差函數這個術語。在井上計算ln(Vp)、ln(Vs)、ln(Den)的自協方差函數以及相互之間的互協方差函數，通過與理論協方差函數的套合，本文選擇高斯模型作為理論協方差函數，協方差函數的高斯模型的理論表達式見式(19),套合參數與效果見表1與圖5。C(h)=Ce-(h/a)2(19)其中C為空間變量方差，即分離距離為0時的協方差，a為變程。值得注意的是，在選擇理論模

24、型以及求取理論模型參數時，最小二乘法的曲線擬合一般來說不是一種最佳的，一般來說甚至不是合適的選擇變差函數模型的方法，必須通過人為的觀察和套合來產生理論模型。在套合過程中，應將更多的注意力集中在較小分離距離的協方差函數值。在這里我們選擇前15個分離距離來套和協方差函數。U.U20.0150.010.0060C12468101214ln(Vp)自協方差函數圖5協方差函數計算套和Fig.5Nestedanalysisofcovariancefunctionln(Vp)自協方差ln(Vs)自協方差ln(Den)自協方差ln(Vp).ln(Vs)互協方差ln(Vp).ln(Den)互協方差ln(Vs).

25、ln(Den)互協方差方差C0.01620.01174.01e-40.01120.00170.0016變程a8.205811.00716.295711.12875.94267.0830(ms)表1套和協方差函數參數Table1.Covariancefunctionparametersofnestedanalysis2.2反演結果不確定性分析這里主要討論在空間變量結構約束下與沒有空間結構約束情況下的反演結果比較。在這里選用極大化后驗概率解作為反演解，也即后驗驗概率的數學期望。圖7是空間結構約束下的反演，圖6是不受空間結構約束的常規反演。從圖上可知，對于兩種反演結果的精度來說,兩種反演方法都取得了

26、很好的解估計。受空間結構約束與否，似乎并不影響反演的準確度。正如前文的討論可知，衡量一個解估計的好壞，除了反演精度之外，另一重要的方面是反演結果的不確定性程度。對于不確定性程度評估，此文應用了兩種方法，一種方法是是對后驗解空間進行隨機模擬，定性地探索后驗解估計空間的不確定性程度；另一種方法是，應用統計推斷上置信區間的概念來進行分析，這是一種定量化的不確定性分析方法。juoO圖6不受空間結構約束Vp、Vs、Den反演結果（黑色點線為先驗低頻背景值，紅色階梯線為實際井數據，藍色線為反演結果）Fig.6InverseresultsofVp,Vs,Denwithoutspatialstructural

27、constraint(Blackdotlinereferstopriorilow-frequencygroundvalue,redstep-linereferstopracticalwelldataandbluelinereferstoinverseresult)圖7空間結構約束Vp、Vs、Den反演結果黑色點線為先驗低頻背景值，紅色階梯線為實際井數據，藍色線為反演結果）Fig.7InverseresultsofVp,Vs,Denwithspatialstructuralconstraint(Blackdotlinereferstopriorilow-frequencygroundvalue,

28、redstep-linereferstopracticalwelldataandbluelinereferstoinverseresult)2.2.1隨機模擬不確定性分析事實上，由（16）式得到的是關于解空間的概率表達，它給出了每一個可能的解估計的概率。由（17）式和（18）式計算的期望和方差，唯一地決定了關于解空間后驗高斯分布的所有性質。由于具有后驗分布的解析表達，所以可以應用隨機模擬這種方法來部分地探索解空間，給出定性的關于解空間的不確定性評價。在這里后驗解空間進行50次隨機模擬實現，實現結果如圖8和圖9。相比于在空間結構約束下的隨機模擬實現，不受約束的解波動性更強，變化范圍較大，說明解的

29、不確定性程度更強。而受約束的解，則很好地限制在離真實井數據較近的范圍內波動。由這50次的隨機模擬結果，可以定性地說明此文的方法確實一定程度上減小了解的不確定性程度。35003525355035753600362536503675370037253750Den24002500260027002800Eg/mJ圖8不受空間結構約束解的隨機模擬實現灰色線為隨機模擬解，紅色階梯線為實際井數據）Fig.8Stochasticsimulationrealizationwithoutspatialstructuralconstraint(Greylinereferstostochasticsimulatio

30、nsolution;redstep-linerepresentspracticalwelldata)juFl-3VPjo535253525350535U53飛5.375530500256663336/3ju375237050025666333jo37523720002500300035004000nVsjo51324002500260027002800Kg/mJ5o5553v5Eh3jlj6352635063jo37635237圖9空間結構約束下解的隨機模擬實現灰色線為隨機模擬解，紅色階梯線為實際井數據）Fig.8Stochasticsimulationrealizationconstrain

31、edbyspatialstructure(Greylinereferstostochasticsimulationsolution;redstep-linerepresentspracticalwelldata)2.2.2統計推斷不確定性分析統計推斷上研究解估計的可靠性程度，主要內容包括點估計和區間估計。本文采用區間估計里置信區間的概念來說明解的不確定性程度。具體的表述方式是，在給定同樣解估計、同等置信水平下，應用解估計所處的置信區間的大小來表達表不確定性程度,置信區間越小，表示解的波動性越小，解的不確定性程度越小。對于高斯分布來說，在給定置信水平為1-Q的情況下，估計量的置信區間為：X土zb

32、(20)0/2其中X為估計量，b為估計量標準差，z為0分為點。這里取0.05，置信水平01-0二0.95，則z=1.96。（21）式變為X士1.96b。將式（18）后驗期望作為解估計02量，式（19）可以計算解估計在每個采樣點上的后驗方差，所以可以得到Vp、Vs、Den的后驗置信區間。圖10與圖11給出了不受空間結構約束和在空間結構約束下的先驗置信區間和后驗置信區間。圖中明顯看到,在空間結構約束下的置信區間相比于沒受約束的要減小很多,從這個意義上講認為受空間結構約束的獲得的反演結果是一個更好的解估計，因為解估計被很好地限制了在兩條窄的綠線之間，波動性更小，可靠程度更高。一個量化的解估計評價見表

33、2,列出了兩種解估計對Vp、Vs、Den后驗置信區間相對于先驗置信區間的平均減小程度。350035035753600367636503675370040006000m/s35003526355035763600362636503G7637003253750Vs200030004000m/sDen35003625.35503600362536503675;70037-25-37502400-.2500260Q-2700圖10不受空間結構約束的解估計先驗、后驗置信區間（黑色點線為先驗置信區間，紅色階梯線為實際井數據，綠色線為后驗置信區間）Fig.10Prioriandposteriorconfid

34、enceintervalofsolutionestimationwithoutspatialstructuralconstraint(Blackdotlinereferstoprioriconfidenceinterval,redstep-linereferstopracticalwelldataandgreenlinereferstoposteriorconfidenceinterval)ju53525350355.3050025666333(SUI)auIIHju37_L523735ju35%.50L!zl3ju63o563ju372537200030004U00nVsjo53525350

35、53下-5.3JLI_til352635063ju376237圖11空間結構約束下解估計先驗、后驗置信區間(黑色點線為先驗置信區間，紅色階梯線為實際井數據，綠色線為后驗置信區間)Fig.11(BlackPrioriandposteriorconfidenceintervalofsolutionestimationwithspatialstructuralconstraintdotlinereferstoprioriconfidenceinterval,redstep-linereferstopracticalwelldataandgreenlinereferstoposteriorconfid

36、enceinterval)VpVsDen未約束28.5126.6517.19空間結構約束66.8163.0135.47表2解估計后驗置信區間平均減小百分比（%）Table2.Averagedecreasepercentofsolutionestimationposteriorconfidenceinterval3結論文章假設地下介質彈性參數的對數符合高斯分布情況，在貝葉斯理論框架下，建立先驗信息和似然函數的概率表達，得到了反演解估計的理論表達式，該理論表達式由一個后驗概率密度函數的形式給出，完整地描述了反演解的位置信息(期望)和不確定性程度(方差)。在建立先驗信息時，引入地質統計學里變差函數的

37、概念，以協方差函數的形式建立反演參數的空間結構信息，以此約束地震反演。通過某實際工區的單井反演算例，以置信區間的概念和隨機模擬的方法評估解的不確定程度，驗證了在空間結構約束下，相比于未約束的反演解估計不確定性程度要低。進一步研究的建議：本文只是通過單井的數據建立了垂向上的空間變量結構信息，如何引入橫向上甚至全三維方向上空間變量結構信息約束地震反演，是下一步要研究且富有意義的課題。參考文獻(References)1Cooke,D.A.,andW.A.Schneider,Generalizedlinearinversionofreflectionseismicdata:Geophysics,198

38、3,48(6),665676,楊文采.地球物理反演的理論與方法.北京:地質出版社,1997,1-13.YangWencai,Theoryandmethodofgeophysicalinversion,Beijing:ChinaLandPress,1997,1-13Tarantola,A.,Inverseproblemtheory:Methodsfordatafittingandmodelparameterestimation,1987:Elsevier,ScientificPubl.Co.Inc.Tarantola,A.,Inverseproblemtheoryandmethodsformod

39、elparame-terestimation,2005:SocietyforIndustrialandAppliedMathematics.JohnA.ScalesandLuis,Priorinformationanduncertaintyininverseproblems.Geophysics,2001,66(2):389-397Buland,A.,andH.Omre,BayesianlinearizedAVOinversion:Geophysics,2003a,68(1),185-198.Buland,A.,andH.Omre,Bayesianwaveletestimationfromseismicandwe