空间变量结构约束下AVO贝叶斯反演

空间变量结构约束下AVO贝叶斯反演摘要引入来源不同的地下先验信息约束地球物理反演，降低反演解的不确定性程度是地球物理反演研究中一个重要的课题。地质变量是在空间（时间）上展布的变量。地质变量除了具有统计意义上的随机性外，还依赖于空间（时间）位置，具有一定的相关性和结构性。本文引入地质统计学里变差函数的思想，建立关于地下空间变量结构的描述作为先验信息，结合贝叶斯推断理论，通过似然函数建立空间变量和地球物理观测的联系，得到AVO反演解估计后验概率密度的解析形式。选取东海某油田一口实际测井数据进行反演试算，应用区间估计和随机模拟的方法对解估计的不确定性程度进行评价，结果表明在空间变量结构约束下，解估计的不确定性程度低于未受到约束的反演解，验证了此文方法的有效性。关键词地质统计学，AVO反演，贝叶斯理论，空间变量结构AVOBayesianinversionconstrainedbyspatialvariablestructureAbstractAnimportantsubjectingeophysicalinversionresearchistodecreasetheuncertaintydegreeofinversesolutionbyintroducingsubsurfaceprioriinformationfromdifferentsourcesasaconstraintofgeophysicalinversion.Geologicalvariableexhibitsspatial(time)variability.Apartfromtherandomnessintermsofstatistics,geologicalvariabledependsonspatial(time)locationaswell,possessingcertaincorrelationandstructuralproperty.Thepaperintroducesthetheoryofvariogramfromgeostaticstocreateadescriptionaboutthestructureofsubsurfacespatialvariableasprioriinformation.CombinedwithBayesianinferencetheory,weestablisharelationshipbetweenspatialvariableandgeophysicalobservationvialikelihoodfunction.Finally,wederivetheanalyticalformofAVOinversesolutionestimation'sposteriorprobabilitydensity.PracticalloggingdatafromawelllocatedinanoilfieldinEastChinaSeaisselectedtoconductpilotinversecalculation.Afterwards,uncertaintydegreeofitssolutionestimationisevaluatedbyperformingintervalestimationandstochasticsimulationmethods.Theconsequencesuggeststhatwithspatialvariablestructureasconstraint,theuncertaintydegreeofsolutionestimationislowerthanitofinversesolutionwhichiswithoutconstraint.Thus,methodproposedinthispaperisprovedtobeeffective.KeywordsGeostatistics,AVOinversion,Bayesiantheory,Spatialvariablestructure0引言众所周知，地球物理反问题是一个多解且不适定的问题。对于地震反演来说，其根本目标是利用观测到的地震数据集，求解地下介质的结构与物理参数[1]。除了地震数据外，如何让引入更多不同来源的信息约束地震反演，降低反问题的多解性，获得更真实描述地下结构与性质的解一直是地震反演的关键问题之一。由于反问题具有多解性，所以在反问题研究中，重要的不是构造出一种广义的解估计，而是要对各种可能的解估计进行评价，这正是由于反问题解的非唯一性所决定的[2]。贝叶斯方法是一个解决上述问题很好的工具，Tarantola（1987,2005）使用关于模型空间的先验信息和似然函数构造了反问题解的后验概率密度，完整的描述了解估计的不确定性问题[3][4]；2001年，Scales和Tenorio关于如何理解在反问题研究中的先验信息与不确定性关系，在其论文中给出了精彩的阐述[5]。Buland（2003,2006）等人应用贝叶斯原理在子波估计、AVO线性反演、时移地震反演方面给出了研究方法和应用实例[6][7][8][9]。2009年，陈建江等人专门论述了AVO反演的不确定性分析[10]。在统计学概念下，将反演的地质变量认为是随机变量，一次观测（如测井数据）均为随机变量的一次实现。经典统计学中，通常假定某随机变量采集的样本是完全随机的或在空间上是完全独立的。但通常地，需要通过反演求得的地质变量在时空域上并不完全独立，而是呈现出一定的空间（时间）相关性。19世纪60年代，法国学者G.Matheron发展起来的地质统计学，考虑到空间变量的结构性，以变差函数、克里金估值、随机模拟为基本手段，是研究空间中某些既显示出随机性又显示出结构性的变量的数学地质方法[11]。1994年，Haas和Dubrule结合地质统计学中的随机模拟思想以及地震数据的反演，提出了现在所谓的地质统计学反演[12]。随机反演（地质统计学反演）要求对地下模型空间产生多个实现，比较常用的方法是MonteCarlo等计算代价高昂的方法[13]，使得在实际应用中受到限制。本文在贝叶斯理论框架下，假设地下介质参数服从对数高斯分布，引入地质统计学变差函数的概念，以地质变量的空间分布结构作为先验对AVO反演进行约束。通过先验的地质变量空间分布和描述地球物理观测和地质变量关系的似然函数，建立了反演解的后验概率密度分布的解析形式。应用随机模拟和区间估计的方法，来定性和定量地评价反演解的不确定性程度。数值算例表明，相比于常规的贝叶斯反演方法，考虑地质变量在空间上的结构性来约束反演，更符合实际情况，不确定性程度降低。1方法原理在一般贝叶斯理论框架下，地球物理反演问题的解以及解的不确定性，以后验概率密度的形式来描述：(|)()(|)(|)()LpLpddmmmddmmm（1）()pm为关于模型m分布的先验概率密度；(|)Ldm为似然函数，是模型空间m向数据空间d的映射的概率表达；分母(|)()Lpddmmm在给定数据情况下为一个常数，也即是概率归一化因子；(|)md是所要求的后验概率密度分布。方便下文讨论，简写(|)Ldm、(|)md为()Lm、()m，忽略掉常数(|)()Lpddmmm，（1）式简写为：()()()Lpmmm（2）接下来，本节将分三个部分分别讨论模型先验分布()pm、似然函数()Lm的建立，以及后验概率密度分布的求取。1.1模型先验概率分布在各向同性介质中，地下介质的弹性性质可由纵波速度、横波速度、密度来完整地描述，同样地也可以由纵波阻抗、横波阻抗、密度来描述。在此本文选择关于地下介质的弹性参数：纵波速度Vp、横波速度Vs、密度来作为反演目标参数。关于模型的先验分布在贝叶斯理论中所起的角色，在理论上一直是概率统计学家争论的焦点[14]。尽管理论上关于模型先验信息的理解还存在争论，但在地球物理实际应用中，在已知地质背景规律或钻井、测井资料指导下反演，已经是大势所趋。关于贝叶斯AVO反演，本文以测井数据作为地下模型的先验信息。但不同于常规的测井约束反演，本文不仅从测井数据里提取关于地下介质的低频信息来约束反演，还进一步挖掘了测井数据里包含的地下介质弹性参数的空间相关性（也即地下空间变量的结构）来进一步约束指导反演。地质统计学里的变差函数是描述空间变量结构的天然数学工具，可以证明在二阶平稳条件下，变差函数与协方差函数等价，在描述空间变量结构方面下文将不再区分协方差函数和变差函数。关于地下空间变量的分布，大量井资料分析表明，弹性介质三参数（纵波速度、横波速度、密度）服从对数高斯分布形式。由以上讨论，本文定义关于空间变量结构的先验分布为：()[ln(),ln(),ln()](,)TtVptVsttNmmmμΣ（3）其中t为地震双程旅行时，mμ为模型期望，mΣ为模型协方差矩阵：(){()}[(),(),()]TVpVstEttttmμm（4）{(),()}(,)ijijCovttttmmΣmmΣ（5）在二阶平稳假设下（5）式可简化为：{(),()}()ijCovttmmΣmmΣ（6）其中12||tt，,,,,,,,,,()()()()()()()()()VpVpVpVsVpVsVpVsVsVsVpVsCCCCCCCCCmΣ（7）关于先验模型均值mμ、协方差函数mΣ矩阵元的计算，是通过测井数据数据的统计分析获得。mμ反映了地质变量的低频信息，mΣ则分反映了地质变量的空间结构信息。由此建立了空间变量的先验信息，可以用来约束反演。注意的是（7）式先验协方差函数的计算稍显复杂，但是它是描述空间变量信息的载体，是进行空间变量结构约束的关键。为了使本文的方法理论简洁明晰，在这里不讨论（7）式的具体计算，其详细的计算方式与讨论见反演算例部分。1.2似然函数在概率统计上的似然函数描述的是估计的模型响应与真实数据的相似程度，由模型向数据的映射，也就是我们称之为的动力系统或正演方程。在一定的假设下，可以建立地下弹性参数映射到地震数据的线性关系。1919年Zoeppritz建立了平面波入射时反射系数、透射系数的理论表达式[15]，但其数学表达式复杂、物理意义不直观，本文采用1980年由Aki和Richard提出来的弱纵波反射系数近似式[16][17]，近似式表达形式与Stolt与Welein(1985)所采用的一样：ln()ln()ln()(,)(,)(,)(,)VpVsVptVsttRtatatatttt（8）其中：21(,)(1tan)2Vpat（9）222()(,)4sin()VsVstatVpt（10）2221()(,)(14sin)2()VstatVpt（11）这里()Vpt、()Vst可以认为是弹性参数的低频背景值，可以通过模型先验均值mμ获得，则系数项Vpa、Vsa、a可以先验求得。（8）式对应的矩阵表达式为：RAm（12）系数矩阵[(,),(,),(,)]VpVsatatatA，m是m关于时间的一阶导数。则地震记录d表示为子波S与反射系数R的褶积加上噪音项e：dSReSAm+e（13）假设噪音服从零均值高斯分布：()Nee0,Σ（14）将褶积运算SA记作算子G，由褶积元算的线性性质可知G为线性算子。由此建立了线性正演方程，将（13）式写为：dSAm+eGme（15）1.3后验概率分布由模型先验概率分布与似然函数，在线性高斯条件下可得后验概率分布为高斯分布[4]：|(,)Nm|dm|dmdμΣ（16）其中：1()(')TTm|dmmmemμμΣGGΣGΣdGμ（17）1()TTm|dmmmemΣΣΣGGΣGΣGΣ（18）其中()mmΣΣ，22()mmΣΣ。式（17）完整地描述了反演解的概率性质，包括了反演解的期望与不确定性。由于后验概