0 引言

地震反演是获得地下介质弹性参数、物性参数及含流体性质参数的重要方法，是定量描述含油气储层的核心地震技术。根据不适定地球物理反问题求解方式的差异，地震反演分为确定性反演和随机反演两大类。基于最小化原理求解反问题最优解或近似解的方法，通常称为确定性反演，其估计结果是一个相对平滑的唯一解。依据统计学理论的随机反演方法充分考虑了观测数据及反演过程的不确定性，消除了确定性反演中的平滑性，能够得到地下介质随机模型的多个实现，从而获得多个反演结果^[1-6]，在获取反演最优解的同时有效评价反演结果的不确定性。

地下介质普遍存在非均质特征，含油气储层的非均质特征主要表现为储层空间分布范围变化大，因此储层的弹性、物性及含流体性质在空间上是变化而不均匀的^[7-8]。通常使用统计学的随机介质模型描述地下复杂非均匀介质，利用均值、方差、自相关函数等参数表征地下介质空间变化的统计特性^[9-12]。为精细研究地下介质特性及其参数变化，需要借助随机介质模型描述地下介质的小尺度非均匀特性，进而有效刻画地下油气藏的细节信息。相比确定性反演，地震随机反演方法融合了地层格架、测井、地震等多尺度信息，可以获得更丰富的地下地层信息，使反演结果的分辨率更高^[13-15]，可以描述非均质复杂储层。

近几年，随机反演方法在储层预测、流体检测及页岩甜点预测等方面获得了较好的效果。前人^[16-19]利用随机反演方法求取储层弹性参数、物性参数以及流体因子。郭同翠等^[20]利用叠前地质统计学反演方法预测页岩甜点的空间展布。Azevedo等^[21]利用叠前地震随机反演获得纵横波速度、密度及岩相数据。De Figueiredo等^[22-23]基于岩石物理先验模型的联合贝叶斯反演获得了弹性参数、物性参数等储层特征参数。Mohamed等^[24]结合地震属性与叠后地震随机反演获得了砂岩相预测模型。

以贝叶斯推断为基础的地震随机反演方法整合了未知模型参数的先验概率密度分布信息，通过似然函数建立模型参数与观测地震数据间的关系，获得待反演参数的后验概率密度函数^{[22, 24-29]}。在反演过程中，先验信息模型用于描述地下地质体的空间变化特征，对于获得准确的反演结果非常重要，有利于缩小模型参数的求解空间。但常规的地震随机反演方法以测井数据为硬数据，利用变差函数表征地下地质体的空间展布特征^{[13-16, 24-26, 30-31]}。然而，基于测井数据的随机模拟不能描述复杂地下地质体的空间变化特征，从而降低了非均质储层地震定量表征的稳定性和可靠性。

为了充分利用已知地震、测井数据中蕴含的地下地层空间结构信息，本文基于随机介质理论，依据地震和测井数据得到精确表征地层空间变化特征的非均匀介质特征参数，构建非均匀介质特征参数模型作为后续随机反演的先验信息模型，通过优化算法求解后验概率密度函数，形成了基于非均匀介质特征参数的地震随机反演方法。

1 方法原理

本文引入随机介质构建贝叶斯随机反演的先验模型，利用混合型自相关函数表征非均匀介质的空间分布特征，从已知地震、测井数据中估算纵向、横向自相关长度及自相关角度等非均匀介质特征参数，获得反映地下复杂储层小尺度空间变化特征的信息，为后续反演提供先验模型。

1.1 随机介质基本原理

随机介质是利用统计方法描述地下介质不同尺度非均匀性的介质模型，一般可以用均值、标准差和自相关函数等一阶、二阶统计量表征地下介质的空间随机变化特征。随机介质模型可以表征地下复杂介质的非均匀特性，有效反映含油气储层的小尺度变化特征。在二阶平稳假设条件下，随机介质模型定义为^[9-11]

$ \mathit{m}(\mathit{t},\mathit{x}) = {\mathit{m}_0}(\mathit{t},\mathit{x}) + {\mathit{\sigma }_\mathit{m}}(\mathit{t},\mathit{x}) $

(1)

其中

$ {\mathit{\sigma }_\mathit{m}}(\mathit{t},\mathit{x}) = {\mathit{\delta }_m}(\mathit{t},\mathit{x})\mathit{f}(\mathit{t},\mathit{x}) $

(2)

式中：m₀(t, x)为模型的均值，可描述介质的平均特性，表征介质的大尺度非均匀性，x为横坐标，t为时间；σ_m(t, x)描述平均特性上的随机扰动，表征介质的小尺度非均匀性；δ_m(t, x)一般为测井数据的标准差；f(t, x)为空间分布特征，是满足二维自相关函数R(t, x)的均值为0、标准差为1的二维随机序列。因此，可以使用均值为零的平稳空间随机过程表示地下介质的小尺度非均匀性。

1.2 非均匀介质特征参数

混合型自相关函数表征地下介质的小尺度空间分布特征，通过引入粗糙度因子，混合型自相关函数综合了高斯型和指数型的特点。高斯型自相关函数较好地描述了单尺度平滑的非均匀介质。指数型自相关函数描述了随机介质的多尺度、自相似的特性，能更好地表征实际介质的多尺度特性。混合型自相关函数则更灵活地描述地下介质，适用性更强^[9-10]。因此，本文采用混合型自相关函数R(t, x)描述非均匀介质的空间变化特征，即

$ \begin{array}{l} \mathit{R}(\mathit{t},\mathit{x}) = {\rm{exp}}\left\{ { - \left[ {\frac{{{{(\mathit{x}{\rm{cos}}\mathit{\theta } + \mathit{t}{\rm{sin}}\mathit{\theta })}^2}}}{{{\mathit{a}^2}}}} \right.} \right. + \\ \;\;\;\;\;\;\;\;\;\;\;\;\left. {{{\left. {\frac{{{{(\mathit{x}{\rm{sin}}\mathit{\theta } + \mathit{t}{\rm{cos}}\mathit{\theta })}^2}}}{{{\mathit{b}^2}}}} \right]}^{\frac{1}{{1 + \mathit{\eta }}}}}} \right\} \end{array} $

(3)

式中：a、b分别为横向(x方向)、纵向(t方向)自相关长度；θ为自相关角度；η为粗糙度因子，η=0为高斯型自相关函数，η=1为指数型自相关函数。

综合式(1)~式(3)可知，描述地下非均匀介质的特征参数包括均值、标准差、纵向和横向自相关长度及自相关角度等，通过这些参数可以模拟产生满足相应的自相关函数且具有指定均值和标准差的随机介质模型。

1.3 非均匀介质特征参数的估计及建模

为从实际地震数据中获得地下地层空间结构特征参数，首先建立地震记录、地震子波与统计特征参数模型之间的函数关系^[11]。根据随机介质理论可知，在连续介质中，波阻抗Z(t, x)与反射系数r(t, x)的关系为

$ \mathit{r}(\mathit{t},\mathit{x}) = \frac{1}{2}\frac{\partial }{{\partial \mathit{t}}}{\rm{lg}}\left\{ {{\mathit{Z}_{\rm{o}}}(\mathit{t},\mathit{x})[1 + {\mathit{\delta }_\mathit{z}}(\mathit{t},\mathit{x})]} \right\} $

(4)

式中：Z₀(t, x)为表征低频分量的背景波阻抗；δ_z(t, x)为表征高频分量的扰动量。

由式(4)、线性褶积模型以及卷积的微分性质可知，地震数据功率谱S_s(ω, k_x)满足

$ {\left| {{\mathit{S}_{\rm{s}}}(\mathit{\omega },{\mathit{k}_\mathit{x}})} \right|^2} \approx {\left| {{\mathit{S}_{{\mathit{\delta }_{\rm{Z}}}}}(\mathit{\omega },{\mathit{k}_\mathit{x}})} \right|^2} \cdot {\left| {{\mathit{S}_\mathit{\psi }}(\mathit{\omega })} \right|^2} $

(5)

式中：S_δZ(ω, k_x)为随机介质的功率谱，ω为角频率，k_x为波数；S_ψ(ω)为ψ(t)的功率谱，$ \mathit{\psi }(\mathit{t}) = \frac{1}{2}\frac{{\rm{d}}}{{{\rm{d}}\mathit{t}}}\mathit{w}(\mathit{t})$，w(t)为子波。

对S_δZ(ω, k_x)进行傅里叶反变换得到随机介质的自相关函数R(t, x)，进而估算出纵向、横向自相关长度和自相关角度等参数。

利用估算的特征参数求取二维自相关函数及其功率谱，得到非均匀介质的振幅谱，在区间[0，2π)产生均匀分布的随机相位信息。综合振幅谱和相位谱得到非均匀随机介质的频谱函数，由傅里叶反变换得到二维随机序列，对其进行标准化，使其均值为零、标准差为1，再结合由测井数据得到的标准差和均值，从而构建非均匀随机介质模型^[12]。

2 非均匀介质特征参数分析

由于均值表征介质的大尺度平均特性，标准差反映偏离平均值的程度，为更好地突显小尺度非均匀特性，本文给定具体的均值和方差，讨论不同纵向、横向自相关长度以及不同自相关角度对非均匀波阻抗介质模型的影响。

图 1为不同纵向、横向自相关长度及不同自相关角度的非均匀介质模型。由图可见：①自相关函数反映介质在空间的相关程度，a展示了介质在水平方向的相关范围，反映介质的横向尺度；b展示了介质在垂向的相关范围，反映介质的纵向尺度(图 1a、图 1b)。②θ反映随机介质模型的扰动方向(图 1c、图 1d)，通过θ模拟地下介质层位的方向变化。

总之，通过分析不同尺度的非均匀介质参数可知，纵向、横向自相关长度和自相关角度等统计量描述了介质弹性参数的空间扰动，可以刻画介质的空间变化特性。因此，借助非均质介质特征参数构建的模型反映了地下介质的非均质特性，为后续反演提供了更可靠的先验信息。

3 基于非均匀介质参数的随机反演

基于上述方法原理构建表征地下介质非均质特性的先验信息模型，在贝叶斯理论框架下，联合似然函数得到反演目标函数，进而利用优化算法求解目标函数，最终形成基于非均匀介质特征参数的随机反演方法。

3.1 贝叶斯理论

贝叶斯理论融合了待反演参数模型和已知数据的先验信息，通过似然函数转化为后验信息，得到模型参数的概率估计，其后验概率密度函数为^[32]

$ \mathit{p}(\mathit{\boldsymbol{m}}|\mathit{\boldsymbol{d}}) = \frac{{\mathit{p}(\mathit{\boldsymbol{m}})\mathit{p}(\mathit{\boldsymbol{d}}|\mathit{\boldsymbol{m}})}}{{\mathit{p}(\mathit{\boldsymbol{d}})}} \propto \mathit{p}(\mathit{\boldsymbol{m}})\mathit{p}(\mathit{\boldsymbol{d}}|\mathit{\boldsymbol{m}}) $

(6)

式中：p(m | d)为模型参数m的后验概率密度函数，表示观测数据d已知时m的分布规律；p(m)为m的先验概率密度函数，表示d未知时m的分布规律；p(d | m)为似然函数，表示m与d的拟合程度；p(d)为d的分布情况。

3.2 目标函数

本文以构建的非均匀介质模型为地质统计先验模型，并结合待反演参数与地震数据之间的正演关系构建似然函数，进而以贝叶斯理论为基础^[32]，得到表征最终反演结果的后验概率密度函数。为减少反演结果的不确定性，增加了非均匀介质特征参数的约束项。因此，最终建立的目标函数为

$ \begin{array}{l} J = \sum\limits_{k = 1}^N {{{\left( {{s^k} - s_{{\rm{obs}}}^k} \right)}^2}} + \alpha \sum\limits_{k = 1}^N {{{\left( {{R^k} - R_0^k} \right)}^2}} + \\ \;\;\;\;{\alpha _1}\sum\limits_{k = 1}^N {{{\left( {{a^k} - a_{{\rm{obs}}}^k} \right)}^2}} + {\alpha _2}\sum\limits_{k = 1}^N {{{\left( {{b^k} - b_{{\rm{obs}}}^k} \right)}^2}} + \\ \;\;\;\;{\alpha _3}\sum\limits_{k = 1}^N {{{\left( {{\theta ^k} - \theta _{{\rm{obs}}}^k} \right)}^2}} \end{array} $

(7)

式中：s^k为模拟地震记录；s_obs^k为观测地震数据；R^k为反演参数；R₀^k为根据测井数据构建的平滑约束信息；α、α₁、α₂、α₃表示加权系数；a^k、b^k和θ^k分别为反演过程估计的横向、纵向自相关长度和自相关角度；a_obs^k、b_obs^k和θ_obs^k为由估计地震数据得到的特征参数。

3.3 非常快速量子退火优化算法

非常快速量子退火(Very Fast Quantum Annealing，VFQA)算法是量子退火算法(Quantum Annealing，QA)的改进算法，采用依赖于温度的似Cauchy分布产生新的扰动模型，可以加快算法的收敛速度。

为提高随机反演的计算效率，本文采用VFQA算法扰动更新模型数据，扰动方式为

$ \left\{ {\begin{array}{*{20}{l}} {{\mathit{y}_\mathit{i}} = \mathit{T}{\rm{sgn}}(\mathit{u} - 0.5)\left[ {{{\left( {1 + \frac{1}{\mathit{T}}} \right)}^{\left| {2\mathit{u} - 1} \right|}} - 1} \right]}\\ {{\mathit{m}_\mathit{j}} = {\mathit{m}_\mathit{i}} + {\mathit{y}_\mathit{i}}({\mathit{B}_\mathit{i}} - {\mathit{A}_\mathit{i}})} \end{array}} \right. $

(8)

式中：m_j∈ A_i, B_i为修改后模型参数，A_i, B_i ∈ N_x, N_t为模型的搜索空间，j∈ N_x, N_t，N_x×N_t为模型的大小；m_i为修改前模型参数，i∈ N_x, N_t；y_i为由u产生的随机变量，u∈[0, 1]为均匀分布的随机数；T为当前温度。

根据Metropolis接受准则：当模型更新前、后的目标函数之差ΔE < 0时，接受扰动后的模型；当ΔE≥0时，以${\rm{exp}}\left[ {\frac{{ - (\Delta \mathit{E})}}{{\mathit{kL}}}} \right]$的概率接受模型。改进的接受准则为

$ p = \left\{ {\begin{array}{*{20}{l}} {1\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\Delta \mathit{E} < 0}\\ {{\rm{exp}}\left( {\frac{{ - \Delta \mathit{E}}}{{\mathit{kL}}}} \right)\;\;\;\;\;\;\Delta \mathit{E} \ge \mathit{0}} \end{array}} \right. $

(9)

式中k、L均为可调参数。

3.4 基于非均匀介质特征参数的随机反演

由地震数据和测井数据估计非均匀介质参数，基于估计的参数构建随机介质模型作为反演的先验信息；结合地震数据与待反演参数之间的关系构建似然函数；然后在贝叶斯理论框架下得到待反演参数的后验概率分布函数，进而得到反演目标函数；最后利用VFQA优化算法优化目标函数得到最终的反演结果(图 2)。

4 模型试算与实例分析 4.1 模型测试

为验证上述方法的有效性，选用二维波阻抗模型(图 3)测试、分析，该模型具有很明显的非均匀特征。在具体反演过程中，分别在第10、20、30和40道处抽取四口伪井，并给定主频为30Hz的雷克子波，将计算得到的模型的合成地震记录作为实际地震数据进行反演。图 4为估计的b、a及θ。由图可见，估计的三个参数与模型数据吻合较好，如在第250个采样点处，θ(图 4c)与模型数据的角度变化趋势一致，b(图 4a)和a(图 4b)也能反映两个薄层的空间变化特征。

另外，考虑到模型的非均质性，在参数估计时采用加窗处理，由地震数据目标层段的非均质特性确定滑动窗口的尺度，得到的特征参数相当于在窗口范围内对模型进行了平滑处理。因此，估计的参数不能与模型数据完全一致，但已足够反映模型数据的空间趋势变化特征。

基于上述非均匀介质模型的构建方法构建非均匀介质先验信息模型(图 5)，其反映了模型数据(图 3)的空间展布特征。图 6为模型数据的随机反演结果与常规稀疏约束脉冲确定性反演结果。由图可见，随机反演结果的分辨率更高(图 6a)。图 7、图 8分别为第20道和第40道的随机反演结果、常规确定性反演结果与原始模型数据对比。可见，相对于确定性反演结果(图 8)，随机反演结果在幅值和形状上均与模型数据匹配较好(图 7)。图 9为QA算法、VFQA算法接受概率随迭代次数的变化曲线。

由图可见，VFQA算法在迭代约100次时即可收敛(图 9b)，而QA算法则需要迭代近400次才可收敛(图 9a)。

4.2 实例分析

本文选取中国M油田叠后地震数据(图 10)测试、分析，其中间部分反射同相轴变细，局部分辨率较高。图 11为对图 10估计的b、a及θ。由图可见：在图 10的局部分辨率较高位置处，b的中间部分数值偏小(图 11a)，这是由于非均质体的纵向尺度减小所致；a反映了同相轴的横向变化(图 11b)；θ反映了同相轴的倾斜变化，显示地震同相轴先趋于平缓，然后随着道数增加倾角逐渐增大(图 11c)。

首先计算测井数据的标准差等参数，然后由克里金插值建立测井数据插值模型作为均值，进而利用上述建模方法估计特征参数构建地质统计随机介质模型(图 12)。该模型充分融合了已知测井和地震数据蕴含的地下地质信息，可以作为后续反演的先验信息模型。图 13为实际数据随机反演与常规确定性反演结果。由图可见，相对于常规确定性反演结果(图 13b)，基于非均匀介质统计特征参数的随机反演结果的分辨率较高，且横向连续性也较好，更准确地指示了含气储层的位置(图 13a黑色椭圆处)。

5 结论

地震随机反演是描述地下复杂储层的一种非常有效的方法，可以获得高分辨率的储层弹性参数数据体。本文提出了一种基于非均匀介质特征参数的随机反演方法，借助随机介质理论从已知地震和测井数据中提取描述储层非均质特性的纵向、横向自相关长度和自相关角度等统计特征参数，为后续随机反演提供可靠的先验信息模型，实现了波阻抗高分辨率反演。模型试算和实例分析表明，由于充分融合了已知地震和测井数据中蕴含的地下信息，估计的非均匀介质特征参数描述了非均质储层弹性参数的空间扰动特性，反映了储层的空间结构特征，为随机反演提供了可靠的地质统计先验信息，高分辨率随机反演的精度较高。但是，由于需要从地震数据中估计特征参数，地震数据的品质影响反演效果。另外，由于加窗处理使估计的特征参数与地震数据不完全一致，也会影响后续反演精度。因此，需要进一步研究利用地震数据估计特征参数的有效方法。