喻思羽
(长江大学地球科学学院 湖北 武汉 430100)
地质统计学提出于20 世纪50~60 年代,起初用于南非金矿储层预测。地质统计方法论是一种建立在空间自相关性分析基础上的空间插值方法体系。使用地质统计学方法可以得到最优无偏估计量的插值模型,与传统的各类样条函数插值法和反距离加权插值法相比,地质统计模型不仅能有效表达空间变量的分布特点,而且可以量化评估结果误差及不确定性,因此地质统计方法广泛应用于地质、石油、固体矿产等地下资源勘查领域,而且在环境保护、卫星遥感、地下水资源、农业、现代商业等其他领域也取得了丰硕成果。
为了让学生更容易理解变差函数概念,本文深入浅出地从传统统计学过渡到变差函数,着重讲解其物理意义及其特点,以多个实例演示不同空间结构模型的变差函数形态,让学生掌握变差函数的本质。本文首先从传统变量与区域化变量的异同点入手,由传统的双变量相关性引入到空间单变量的自相关性,然后再由空间自相关性过渡到变差函数的概念及其物理意义,最后讲解理论变差函数模型的参数意义和不同空间结构的变差函数曲线特征,为学生使用变差函数分析地质特征打下基础。
地质统计学的研究对象称为区域化变量,这种变量与传统统计变量的不同之处在于,区域化变量是一种分布在空间区域中的变量,并且具有一定结构性和随机性。区域化变量在现实中经常遇见,例如气温、降雨量、储层地质体的物性参数等。空间数据采样具有局限性——某个位置只能采样一次,无法通过对某个位置进行大量采样获取该点的统计量(均值、方差等)。地质统计学重点集中于区域化变量的空间自相关性。以例子简单阐述传统变量的相关性以及与区域化变量的自相关性的区别。假如有一个小卖部,我们记录N 天的雪糕销量与当天气温的数值,得到了如图1(p59)的一组数据,并建立销量与气温两个变量的映射关系,对数据进行可视化(图2,p59),其中图2[a]是对销量与气温的归一化曲线图,每天统计的两个变量大致具有规律:当一个变量增加时,另外一个变量也增加,反之亦然,反映了两个变量具有一定的相关性。为了定量计算两个变量的相关性,建立如图2[b]的变量映射散点图,该图中每个点的X 轴代表气温变量(归一化数据),Y 轴代表销量变量(归一化数据),将每个点投影到45 度对角线上,可以统计得到当前日期的两个变量相关性值,很明显当点偏离对角线距离越远,说明该点将拉低两个变量的相关性,反之将提高相关性。
图1 某小卖部的雪糕销量(元)与当天的气温(摄氏度)
图2 气温(度)与雪糕销量(元)的相关程度
协方差是衡量两个变量的相关性常用指标,协方差计算公式如下:
基于协方差对“气温与销量”的关系分析可知两个变量在变化过程中是同方向变化,还是反方向变化,以及同向或反向变化的程度。例如气温与销量的关系,气温升高,销量也随之增加,说明这两个变量是同向变化的,则协方差为正,否则两个变量的协方差为负。
以上是传统变量的相关性分析方法,在地质统计学中,通常研究的对象称之为区域化变量,属于随机场的范畴。为了便于描述,下面以一维随机过程为例讲解如何计算区域化变量的相关性,进而引出变差函数的概念。如图3 所示是一个一维随机过程,随机变量Z随着时间t 随机波动。为了计算时间间隔为△t 的两个时间点的变量Z(t)与Z(t+△t)之间的相关性,需要从随机过程f 中提取若干个间隔△t 的点对。如图3 所示,t1 和t1+△t 属于一个点对,t2和t2+△t 属于一个点对。注意这里计算内容称为单个变量的自相关性,即变量Z随着时间t 变化的自相关性,而不是两个变量的相关性,此时如果把横轴的时间改为空间位置或者距离,计算结果就变换为空间自相关性了。如图4所示,得到间距△t 的点对映射。此时将变量Z(t)与Z(t+△t)作为两个不同的变量,基于传统变量的协方差计算公式即统计间距△t 的点对相关性值,即变量t 的自相关性。
图3 随机过程中间距为△x
图4 某小卖部的雪糕销量与当天的气温
传统地质统计学是以空间两点相关性,以变差函数为工具完成空间相关性分析,进而进行未采样区域的预测。1965 年法国统计学家Matheron 提出基于距估计的变差函数,定义为在相距h(称为滞后距,属于矢量)的区域化变量Z(x)与Z(x+h)增量的方差之半,计算公式定义为
其中x 是采样点的位置,变差函数反映了区域化变量的空间结构性,用于评价区域化变量在某个方向某个距离的变化程度。实验变差函数是指应用观测值计算的变差函数,计算公式如下
其中是空间两点的位置矢量,N()为点对的总数,和代表两个点,是区域化变量。
为了在教学中使学生易于理解变差函数,以下面实例深入浅出地解释变差函数的物理意义。如图5 所示,在研究区域内有多个采样数据点(实心点表示采样点,空白位置表示未采样点)。点与点之间为规则网格分布,两个点之间的水平或垂直间距为100 米。
图5 某区域内的区域化变量分布
现计算滞后距等于1(东西方向、间距100 米)的空间两点相关性。如图6[a]所示,首先从图5 中提取滞后距为1 的点对,将所有的点对按照起点值作为横轴、终止点值作为纵轴投在图6[b]的散点图上。根据公式2 计算图6[b]中全部点到45 度对角线的距离,统计得到滞后距等于h1的变量Z 增量的方差,即算出该滞后距的变差函数值。地质统计学理论强调距离相近的观测值比距离较远的观测值更加相似,即方差较小(Matheron,1963),反之随着距离增加,区域化变量的空间自相关性逐渐减小。通过计算不同滞后距的变差函数值,建立变差函数值与滞后距h 的函数关系(图7),就能算出评价空间相关性的一组特征定量化参数,常用参数包括变差函数的块金值、变程和基台值。由图7 可以发现,随着滞后距的增加,区域化变量Z 的增量方差首先快速增加,然后逐渐收敛,表明变差函数值不会无限增加。
图6 间距h=100 的点对样本映射关系
图7 实验变差函数
变差函数曲线能准确表征出区域化变量的结构性变化规律,例如地质结构渐变、突变、周期性以及变化幅度都能敏感地投射到变差函数的曲线波动特征上。以两组实例展示模型结构的实验变差函数特征。图8(p61)中有3个模型,模型中有两种相类型,用编码0 和1 代表,图8[a1]中以背景相编码为0,目标相编码为1 的垂直高度为10 个网格单元,图8[b1]中以背景相编码为0,目标相编码为1的垂直高度为30 个网格单元,图8[c1]中以背景相编码为1,目标相编码为0 的垂直高度为30 个网格单元。图8[a2]、图8[b2]和图8[c2]分别是图8[a1]、图8[b1]和图8[c1]的垂直方向的实验变差函数,滞后距取值范围为0~80,同时进行理论变差函数拟合。观察模型与对应变差函数曲线特征可以发现,图8[a1]的实验变差函数在滞后距约为10 网格单元的位置达到第一个转折点,这与目标相的垂向高度吻合。分析原因为:由实验变差函数公式(公式3)可知,该模型的变差函数值大小取决于垂直方向上相距h 的点对差是否等于1,如果点对的差等于1,则能增加变差函数公式的分子项,如果点对的差为0,不会增加变差函数分子项。随着h 的增加,如果分子项没有增加,则将降低变差函数最终结果。因此实验变差函数曲线的几个拐点对应着实验变差函数分子项出现变化的位置。因此图8[a2]中h 达到40 时纵轴值开始下滑,反映了h 大约40 之后,点对的差值等于1 的比例开始快速减少。在图8[b2]中,曲线的第一个转折点对应的h 值为30,与图8[b1]的编码为1 的相高度一致。值得注意的是,变差函数的大小只与某个滞后距的点对数量(实验变差函数的分母项)和点对差的绝对值(实验变差函数的分子项)有关。因此如果点对数量相同,点对差绝对值之和越大,则变差函数越大,互换背景相和目标相的编码对变差函数值没有本质影响,例如图8[b1]和图8[c1]中的两个模型,它们的相编码进行互换,对应的实验变差函数曲线特征(图8[b2]与图8[c2])完全一样。
图8 三个模型与实验变差函数,其中a2 是a1 的垂直方向变差函数,b2 是b1 的垂直方向变差函数,c2 是c1 的垂直方面变差函数
本文针对石油地矿类高校中地质统计学的关键概念——变差函数讲授过程较为枯燥和公式化的问题,将变差函数的基本思想和物理意义作为重点,通过多个实例逐步将变差函数内涵展现给学生。此外还阐述了变差函数在地质研究中的作用特点。本文重点辅助地质统计学的入门学习,对变差函数的各种高级复杂特点不做探讨。