2.2 研究方法
本书通过降水统计模型建立、反距离加权法、克里金方法、协克里金方法和改进的协克里金方法进行插值方法优选。
2.2.1 降水量的统计模型
目前对降水的推算主要用以下两种方法:①将影响某一地区降水的要素分为若干修正项,然后进行综合,得出任意点的降水的推算值;②建立降水与其影响因素之间的多元统计模型,用以推算降水值。第二种方法具有计算简便、容易在地理信息系统中实现的优点,因此本文选择第二种方法,即根据实测点的降水,建立降水和位置、海拔、地形因子之间的关系,模拟降水信息的空间变化规律。模型可表示为
式中:P为降水量,mm;h为计算点的高程,m;λ为计算点的纬度;φ为计算点的经度。
2.2.2 反距离加权法
反距离加权法是以插值点与样本点之间的距离为权重的插值方法,插值点越近的样本点赋予的权重越大,其权重贡献与距离成反比。可表示为
式中:Z为插值点估计值;Zi为实测样本值,i=1~n;n为参与计算的实测样本数;Di为插值点与第i个站点间的距离;p为距离的幂,它显著影响内插的结果,幂越高,内插结果越具有平滑的效果,它的选择标准是最小平均绝对误差;当取p=2时,即称作反距离平方加权法。
2.2.3 克里金方法
克里金方法被认为是地学统计最主要的方法之一,是建立在地质统计学基础上的一种插值方法。该方法最早由法国地理学家Matheron和南非矿山工程师Krige提出的,最早用于矿山勘探。这种方法充分吸收了地理统计的思想,认为任何在空间连续变化的属性是非常不规则的,不能用简单的平滑函数进行模拟,只可以用随机表面函数给予恰当描述。克里金方法的关键在于权重系数的确定。该方法在插值过程中根据某种优化准则函数来动态地决定变量的数值,从而使内插函数处于最佳状态。从统计意义上说,克里金法是从变量相关性和变异性出发,在有限区域内对区域化空间变量的取值进行无偏、最优估计的一种统计方法;从插值角度上讲,克里金法以空间结构分析为基础,充分利用了数据空间场的性质,在插值过程中对空间数据求线性最优,可以反映空间场的各向异性,是一种无偏估计的内插方法。该法的最佳适用条件是空间变量存在着空间相关性。
克里金方法的插值公式为
式中:Z(x0)为x0处的估计值;Z(xi)为xi处的观测值;λi为克里金权重系数;n为观测点个数。
权重由克里金方程组式(2.4)来确定
式中:C(xi,yj)为测站样本点之间的协方差;C(xi,xo)为测站样本点与插值点之间的协方差;u为拉格朗日乘子。
插值数据的空间结构特性由半变异函数描述,其表达式为
式中:n为距离为h的采样点对的数目(n对点),采样间隔h也叫滞后;N(h)为被距离区段分割的试验数据对数目,根据试验变异函数的特性,选取适当的理论变异函数模型。
图2.2 半方差图
根据试验半变异函数得到的试验变异函数图,从而确定出合理的变异函数理论模型。对应于h的γ(h)的图被称为“半方差图”。半方差是定量描述区域性变化的第一步,它为空间插值、优化采样提供了有益的信息。为了求得半方差图,必须先得到拟合半方差的理论模型。图2.2表示一个典型的半方差图。
图2.2中的c0为核方差,称为金块效应(nugget effect),是观测误差和距离间隔很小的情况下的空间变化的组合。a称为变程(range),它描述了与空间有关的差异随距离变化的方式。在变程范围内距离越近的点具有更相近的特征,数据点和未知点之间的距离大于变程范围时,表明该数据点与未知点距离太远,对插值没有作用。c称为基台值,它反映了某区域化变量在研究范围内变异的强度。曲线的水平部分称为“梁”。
关于理论半方差函数,假设它是各向同性的。常用的理论半方差函数有:球面模型、指数模型、高斯模型。本文采用球面模型。
2.2.3.1 球面模型
球状模型,其公式为
式中:a为变程;h为延迟;c0+c1为梁。
一般情况下用球面模型拟合效果比较理想。当存在明显的变程和梁,同时核方差也很重要,但数值不太大的情况下,可用球面模型进行半方差拟合。
2.2.3.2 指数模型
当存在明显的变程和梁,而没有渐变的变程,则可用指数模型进行拟合。其公式为
式中a不是变程,该模型的变程为3a。
2.2.3.3 高斯模型
如果核方差和空间变化有关的随机变化相比很小,则可用高斯模型进行拟合。其公式为
式中a也不是变程,该模型的变程为3a。
克里金方法包括普通克里金方法、简单克里金方法、泛克里金方法及协克里金方法等,最常用的是普通克里金插值方法。普通克里金方法认为,当空间变量的结构性成分确定后,剩余的差异变化属于同质变化,不同位置之间的差异仅是距离的函数。
2.2.4 协克里金方法
建立在地质统计学理论基础上的协克里金方法的基本原理与克里金插值方法相同,但它通过考虑一个以上的变量而优化估计。对多个具有空间相关性的空间变量进行估计的克里金方法可以归类为协克里金方法。借助这类方法,可以利用几个空间变量之间的相关性,对其中的一个变量或多个变量进行空间估计,以提高估计的精度和合理性。协克里金是通过建立交叉协方差函数(cross covariance)模型来进行估值的。
设有k个协同区域化变量Z1(x),Z2(x),…,Zk(x),则每一对区域化变量Zk(x)和Zk(x)的交叉变异函数为
虽然协克里金方法的基本原理早已为人们所熟知,但至今它仍是一个不经常使用的方法。其主要原因就在于其符号和算法比较复杂,交叉协方差函数和交叉变异函数的求取比较困难。
2.2.5 改进的协克里金方法
采用协克里金方法,并将高程、坡度、坡向作为影响因素引入到对降水量的空间插值中来。考虑高程、坡度、坡向影响的协克里金估值方法可表示为
式中:Z(x0)为x0点处估计降水量;Z(xi)为第i站的实测降水量;y(xi)为xi点处的高程;n为实测雨量站个数;mh、mp为高程和降雨的全局平均值;ai为协克里金插值综合权重系数,由坡度、坡向、距离权重系数决定。
综合权重系数公式
式中:ai为第i个站点的综合权重;aij为第j个权重对第i个站点的权重;Aij为第i个站点的每j个单个权重因子。
利用这种综合权重计算方法,可定量地表示出不同权重因子对不同降雨量的影响。
2.2.5.1 坡度权重
坡度因子对气象要素的影响是非常大的,不同坡度的地表接受到的太阳辐射相差很大。左大康等发现当坡度在30°以上时,坡度对太阳辐射的影响在所有坡向上都很明显。当相同的太阳高度角时,不同坡度的地表对太阳辐射的吸收、辐射、反射也不同。因此,在气象要素的空间插补中有必要考虑坡度因子对内插结果的影响。研究中,坡度被分为0°~30°、30°~60°、60°~90°3个区间而加以量化,方便了坡度权重的计算,对处在不同的坡度区间的站点赋予不同的整数值,以便于对不同的站点给予不同的权重。某一站点坡度权重的赋值是按以下原则实施的,即:如果站点的坡度和需要插补的DEM格网上的坡度在相同的区间内,那么赋予该站点较高的权重(研究中取0.55),否则赋予较小的权重(研究中取0.2)。考虑到流域地形的复杂情况,该权重系数可根据插值结果对比分析进行优化调整,使插值结果与实测值误差较小。
2.2.5.2 坡向权重
地表坡向也是对气象要素有很大影响的因子之一,尤其对于气温、降水及蒸散发。对于降水,坡向的不同会引起背风坡和迎风坡气流的空间差异,从而形成坡向降水导致降水的空间变化。在本书中,坡向是按顺时针方向把地表分为8个坡向区间(每个区间大小为45°),以正北方向为零计算。采取类似坡度处理的方法,把DEM格网不同的坡向赋予不同的整数值,用1、2、3、4、5、6、7、8分别代表北、东北、东、东南、南、西南、西、西北8个方向,通过坡向区间的划分来表示不同站点对DEM格网气象要素不同的影响力。如果某站点的坡向和该DEM格网处在同一个坡向区间时,赋予该站点较大的权重(本书中取0.3),否则分别给予较小的权重(本书中取0.15和0.007)。该权重系数亦可根据插值结果对比分析进行调整。
2.2.6 模型检验标准
采用交叉验证法来验证插值的效果。在空间插值所得到的降水分布图上可以确定出相应检验站点的模拟值,然后计算检验站点实测值与模拟值的误差,以此来判断估值方法的优劣。采用相对均方差(RMSE)作为评估几种插值方法的插值效果的标准。表达式为
式中:n为检验站点数目;为第i个站点的实测值;为第i个站点模拟值。
相对均方差越小表明误差越小,插值精度高。