何 鹏,张 媛,高文波,蹇东南,林正雨
(1.四川省农业科学院农业信息与农村经济研究所,成都 610066;2.四川省农业科学院大数据中心,成都 610066; 3.成都理工大学旅游与城乡规划学院,四川成都 610059;4.四川农业大学资源学院,成都 611130)
随着生态学、农业科学、资源环境科学等相关学科的发展,气象信息在区域和全球尺度变化研究中发挥着重要作用,特别是各类生态系统模型、作物模型都需要空间化的气象要素,并且所要求的时间和空间分辨率也越来越高。例如,在全球变化研究中,MT-CLIM、FOREST-BGC等景观、区域、全球尺度的生态系统模型,都需要空间化的气温、降水、太阳辐射等环境因子作为输入参数[1]。因此,利用现有的地面气象观测资料和相关技术,开发空间化的气候数据产品成为近年来生态学、农业区划、资源和环境科学研究的重点任务之一[2-5]。
然而有限的气象站点仅代表一定区域内气象要素的分布情况,其观测获取的只是局部、离散、有限的空间点数据,对于没有气象站点分布的广大区域,只能通过以点代面或者空间内插等方法推算求得。国内外学者已经在气象要素推算方面做了大量研究。国外,Eleanor等人引入地形,运用多元回归法对气温和降水进行空间化处理[6];Marquinez等人引入高程、坡度、坡向、离海岸线的距离和离相对西边的距离等因子,利用多元回归方法和GIS技术拟合降水空间分布[7];20世纪90年代以来,大量区域性乃至全球性栅格气象数据库不断建立和完善,如基于PRISM建立的美国、加拿大、中国、蒙古、欧洲等国家和地区的空间气象数据信息系统[8];基于Anusplin插值方法的澳大利亚、南非等国家的气象数据信息系统[9,10];基于DAYMET的美国生物气象数据信息系统[11]。国内,20世纪80年代,以傅抱璞、沈国权等为代表的学者通过数理统计模拟方法对气温、降水等要素进行推算,获得了非气象站点所在区域的气象要素空间分布情况[12-14];20世纪90年代以来,随着GIS技术在国内的流行,于贵瑞、陈晓峰等利用GIS软件并结合地形等因素,对区域气象要素空间分布进行了深入研究[3-5,15,16],廖顺宝等、蔡福等还对气温空间化的插值方法、精度进行研究[1,17,18],廖玉芳等以湖南省为例,进行了精细化农业气候区划研究,建立相关的气候空间化数据集[19,20]。
影响气温空间分布的因素很多,主要以海拔高度和地形条件最为显著。一般情况下,随着海拔高度的增加,气温下降,但其变化速率因山地性质和气候条件而不同。在对气温素进行空间插值时,由于海拔高度的影响非常显著,海拔高度的误差对气温空间插值精度的影响不能忽略[21-24]。四川省处于我国东部季风区与西南青藏高寒区的交接地带,气候垂直变化十分显著,加之地形复杂,导致气候的区域差异很大。江蕾采用全国633个站点数据,主要考虑地貌拟合插值模型,生成1km空间分辨率的气温数据,该研究对于四川省这一地貌复杂区的气温空间插值误差研究不够[25]。为了对比和分析海拔高度对四川省多年平均气温空间分布的影响,文章利用5种插值方法,同时考虑数据源密度(3种密度)和空间化栅格大小(10个空间化栅格尺度)对四川省144个气象观测站点的多年平均气温进行了空间插值研究,并利用预留检验站点法对结果进行评估。
气温数据来源于中国气象数据网(http://data.cma.cn/)申请下载的“中国地面累年值月值数据集(1981—2010年)”,包括气温、降水、气压、湿度、风向和风速共五大类12小类的年均气象数据,数据制作时间2012年,共包括四川省145个国家级地面站点(图1)。参考四川省气候区划(自然资源区划),并考虑原有站点分布疏密现状,采用分层随机抽样法,全省预留检验站13个(约10%),剔除数据缺失较多的站点1个,实际131个站点参与插值运算。
气象站点密度设定。数据源密度变化分析以预留的13个检验站点为基础,仍然采用分层随机抽样方法,并考虑原有气象站点分布的疏密度,131个站点按照1/3的比例递减,得到低密度(43个站点),中密度(87个站点),高密度(131个站点)共3个级别。
图1 四川省气象站点分布
数字高程模型(DEM)来源为国际科学数据服务平台(http://datamirror.csdb.cn/)申请下载的ASTER GDEM第一版本(V1)数据,数据时间2009年,数据投影为UTM/WGS84,空间分辨率30m,空间范围为四川省全域。
空间化栅格大小设定。空间化尺度变化分析结合四川省幅员面积和已有研究文献,依次为100、250、500、1 000、2 000、3 000、4 000、5 000、7 500、10 000m共10个尺度等级。
①反距离权重法(IDW)和加入海拔影响因子的反距离权重法(IDWE)。反距离权重法是对采样点进行线性加权决定输出的栅格值,加权值与距离成反比,采样点离输出栅格越远,它对输出栅格的影响越小。反距离权重法的表达式为:
(1)
考虑了海拔高度影响后,由于气温随着海拔高度增高而降低,进行气温空间插值时,海拔高度是附加的重要变量,其计算公式为
(2)
式(2)中,Hk为第k个点的海拔高度,β为气温随海拔高度的递减率,其余参数与IDW同。
②普通克立格法(OK)和加入海拔影响因子的普通克立格法(OKE)。普通克立格法又称局部最优线性无偏估计法,建立在变异函数理论结构分析基础之上。所谓线性是指估计值是样本值的线性组合,即加权线性平均;无偏是指理论上估计值的平均值等于实际样本值的平均值,即估计的平均误差为0,最优是指估计的误差平方差最小。其计算公式为:
(3)
式(3)中,Z(ti,j)为待插值的栅格点要素值;Z(tk)为第k个点的要素值;λk为权重系数,其和等于1。
加入海拔影响因子的普通克立格法基于普通克立格法并充分考虑了海拔因子的影响,其计算公式为:
(4)
③多元回归分析法(MRM)。多元回归分析模型是应用数理统计方法,建立气象要素与影响气象要素空间插值因子的多元回归模型。该文是基于经度、纬度和海拔高度的多年平均气温多元回归模型,并计算残差。
多元回归分析模型表达式为:
T(λ、φ、h)=b0+b1λ+b2φ+b3h+ε
(5)
式(5)中,b0为常数,b1,b2,b3为待定系数,λ为经度,φ为纬度,h为海拔高度,ε为残差。
表1 不同密度气象站点海拔高度、经纬度与平均气温的相关系数
利用地理信息系统软件ArcGIS 10.1的数据管理模块(Data Management Tools)、转换模块(Conversion Tools)、地统计模块(Geostatistical Analyst Tools)和空间分析模块(Spatial Analyst Tools)进行数据预处理和插值,并使用预留站点法检验插值精度。属性数据的处理采用SPSS 20。主要工作如下:
(1)将气象站点的经度、纬度和多年平均气温在ArcGIS 10.1中生成点文件,分别使用两种常规空间插值方法(逆距离权重法和普通克立格法)生成多年平均气温栅格图。
(2)参考游松财[21]提出的气温垂直递减法。首先,在考虑了海拔高度的影响下,根据气温垂直递减率和各气象站的海拔高度将各气象站的年平均气温订正到虚拟海平面上;其次,运用两种常规插值方法(逆距离权重法和普通克立格法)对虚拟海平面上的月平均气温进行空间插值,得到虚拟海平面上的气温空间分布;最后,在GIS环境下,利用DEM数据和气温垂直递减率进行栅格图层的代数运算得到平均气温的空间分布。
(3)利用气象站点的海拔高度、经度、纬度和多年平均气温数据,在SPSS 20中进行多元回归分析,拟合出回归模型。将四川省DEM数据在ArcGIS 10.1中生成value值分别为经度(x)和纬度(y)的栅格图层,使用栅格计算器根据回归模型生成多年平均气温栅格图。
在使用气象站点数据进行气候指标插值结果的精度检验时,交叉检验(cross-validation)是最常用的一种方法。该方法假设某一个气象站点的气温未知,用其他所有气象站点来估算该气象站点的气温,通过计算所有气象站点的实际气温值与估算值之间的误差来评判空间插值方法的精度。对不同的插值方法,交叉验证可以准确的验证不同插值方法之间的相对精度,使用时一般用平均绝对误差(MAE)和均方根误差(RMSE)作为检验精度的标准。MAE反映样本数据估值的总体误差或精度水平,RMSE则反映样本数据的估值灵敏度和极值。计算公式分别为:
(6)
(7)
采用1km×1km栅格尺度,分别采用IDW、IDWE、OK、OKE和MRM 5种方法,以144个气象站点进行多年平均气温数据的空间化(图2),并对144个站点逐一做交叉验证,统计平均绝对误差和均方根误差(表2),其中IDWE精度最高,OKE精度次之,其中未引入海拔因子的IDW、OK模型,精度明显低于有海拔因子参与的模型。
表2 不同模型空间化误差 ℃
图2 5种方法144个站点空间化结果
采用1km×1km栅格尺度,分别采用IDW、IDWE、OK、OKE和MRM 5种方法,按照低、中、高3个气象站点密度等级进行多年平均气温数据的空间化,并统计出13个检验站点的MAE、RSME(表3)。(1)应用IDW和OK插值模型,采用中等密度站点对年平均气温插值精度最好,高密度站点插值精度次之,但两者差异不显著低密度站点插值精度最低,3种数据源密度之间的插值精度差异较大;(2)应用IDWE、OKE、MRM插值模型,随着站点密度的增大,多年平均气温数据空间化误差总体上都保持减小的趋势,3种数据源密度插值精度之间的差异较小,采用高密度站点插值年平均气温误差最小;(3)3种插值模型中,IDWE方法插值精度最好,低、中、高3种密度的年平均气温误差分别为0.74℃、076℃、0.67℃;OKE模型次之,其MAE值分别为0.87℃、0.82℃、0.80℃。
表3 不同站点密度插值精度检验结果 ℃
图3 5种模型方法下10种栅格尺度插值精度的MAE和RSME
针对131个气象站点(即高密度),分别采用上文中5种方法,按照10个尺度等级进行多年平均气温数据的空间化,并统计出13个检验站点的MAE、RSME值,结果见图4。在空间化尺度由100m增大到10km的过程中:(1)IDW和OK两种插值模型随着空间化栅格尺度的增大,其空间化误差变化不大。应用IDW和OK插值模型,平均气温误差分别为0.92~0.99℃和1.26~1.34℃;(2)IDWE、OKE模型下,随着空间化尺度的增大,其空间化误差逐渐增大。空间化尺度小于2km时,其平均绝对误差小于1℃。(3)应用MRM插值模型下,总体来说,随着空间化尺度的增大,其空间化误差逐渐增大。除采用1km栅格尺度外,其余栅格尺度平均气温误差均超过1℃;(4)总体来说,加入海拔因子的IDWE、OKE、MRM插值模型比IDW和OK插值模型精度更高,其中以IDWE插值精度最好,OK模型精度最差。
上述分析表明,5种插值方法中,以IDWE插值精度最好。该部分采用插值最好的IDWE模型,分别采用高、中、低3种站点密度,按照10个尺度等级进行多年平均气温数据的空间化,定量统计13个检验站点的平均绝对误差(表4)。(1)随着空间化栅格尺度的增加,空间化误差逐渐增加;栅格尺度小于等于2 000m时,3种站点密度的平均气温误差小于1℃。其中高密度站点插值精度最好;(2)同样的空间化栅格尺度下,高密度站点的空间化误差最小。高、中、低3种密度之间的空间化误差差距不明显。100m×100m尺度下3种密度的MAE在0.64~0.76℃,10km×10km尺度下3种密度MAE值在1.37~1.69℃。
表4 不同数据源密度和不同空间化栅格大小的空间化误差分布 ℃
游松财等研究表明[21],由于国家气象局提供的气象站点的经纬坐标只精确到分,其空间误差大约±30″(约1km),因此在气温内插时考虑到海拔高度的影响,并且使用DEM数据就必须考虑两者的海拔高度差异[13]。利用四川省144个气象站点的经度及纬度资料,从ASTER GDEM中提取相应的海拔高度,发现有20个站点(约占气象站点总数的13.9%)具有100m以上的差异,其中有9个站点(约占气象站点总数的6.3%)具有200m以上的差异,已有的研究表明全国基础气象站点(约635个)中仅有8%左右站点海拔高度与DEM海拔高度有100m以上的差距,说明我国西南地区由于海拔较高,地貌类型复杂,给气象数据空间化增加了不确定因素。该部分采用IDWE、OKE、MRM 3种插值方法,分析了气象站海拔与DEM海拔插值之间的气温绝对误差。气象站海拔与DEM海拔插值之间的气温绝对误差的空间分布见图4。
表5 海拔高度差异及其温度差异对比
图4 IDWE(a)、OKE(b)、MRM(c)3种插值模型下气象站海拔与DEM海拔插值的气温差
表6 站点海拔与站点DEM值参与插值的各海拔高度其温差 ℃
图4表明,采用IDWE、OKE、MRM 3种插值模型,采取气象站海拔数据插值的气温(Tms)与DEM海拔数据插值的气温(TDEM)均存在差异,3种插值模型下的气温差空间分布相似。小部分区域TDEM小于Tms,主要位于四川的东部区域;大部分区域TDEM大于Tms;IDWE、OKE、MRM模型下,TDEM与Tms之间的气温差分别在-0.98~2.77℃、-0.98~2.69℃、-1.06~3.05℃范围内。3种插值模型下,气温差异最大的地方主要位于四川的西北部、中北部以及东南部。以栅格面积统计,3种插值模型空间化结果中,气温差绝对值大于1℃的区域分别占全省面积的0.76%、1.13%、1.52%,其中气温差绝对值大于1.5℃的区域仅分别占全省面积的0.21%、0.26%、0.31%,所以两种海拔数据插值结果的差异性并不显著。
按500m×500m的栅格密度分别以气象站海拔和DEM海拔进行气温插值,将海拔分为10个级别,分别统计不同方法下各海拔范围内两者之间的气温差(表6)。结果表明:(1)采用3种方法对不同海拔插值的气温存在差异,在0~500m范围内两者之间差异较小,气温差小于0.1℃。(2)海拔大于1 500m时,气温差普遍大于0.3℃;海拔在2 000~2 500m范围内时,气温差平均值为0.709℃;海拔在2 500~3 000m范围内时,气温差平均值为0.706℃;海拔大于5 000m时,气温差平均值为0.637℃。
图5 四川省农业气候区划
以四川省地貌区划和气候区划为参照,将全省分为4个区[26],见图5。采用131个站点,5种方法分别对各区进行插值后拼接,与直接对全省131个站点插值的结果进行对比,仍然以13个校验站点进行分析,其中IDWE表现仍为最优(图6),其MAE、RSME分别为0.44℃、0.62℃,按不同地貌区划进行分区插值,各种方法的精度明显优于全域(图7),且MAE与RSME值的差距缩小。
图6 采用5种方法分区与不分区插值精度比较
图7 IDWE模型下不分区与分区插值结果图(a.不分区插值;b.分区插值;c.已有的四川省多年平均气温分布数据)
对比分区与不分区IDWE模型下的插值结果图,不分区拟合结果年均温的值域为-9.64~21.76℃,分区拟合结果年均温的值域为-17.19~22.01℃,结合已有的四川省多年平均气温统计数据,可以看出分区拟合后气温的值域跨度更符合真实气温分布状况。
空间化模型的选择对空间化精度有较大影响。通过精度分析,以空间自相关为主导的传统空间化模型,如IDW和OK,虽然理论和方法都比较成熟,但对真实气温环境的模拟能力较低,其相对误差较大。而IDWE、OKE和MRM等方法较全面的考虑了海拔、经纬度等因素的影响,对于海拔变化较大地区的气温有更好的表现。
研究结果表明,数据源密度(站点数量)对空间化结果有一定影响,但不显著。随着站点数量的逐渐增加,所有空间化方法的插值精度均有所提高。加入海拔因子的IDWE、OKE、MRM等模型除整体精度高于不考虑海拔因子的方法外,随着站点密度的增大,空间化误差总体减小,精度逐步提高,但3种数据源密度插值精度之间的差异较小,表明即使减少气象站点数据,加入海拔因子的插值模型也能较好地模拟气温。其原因应该是空间自相关模型是直接对站点观测值进行插值,没有考虑其他影响因素,因此站点数量越多即样本数越多,插值精度自然越高。而考虑了海拔、经纬度等因素的异相关模型,方法核心是建立气温因子与其他环境因子的关系,样本数大小反倒对插值精度影响不大。由于一般科研人员能共享的地面气象观测站点数较少,这一验证结果也对我国地面气象站点观测结果插值精度提供了有利的佐证。
栅格大小对空间化精度有一定的影响。一般空间化模型的栅格大小有经纬度和地理空间距离两种,两种方法的空间化尺度一般从100km(1°)到30m(1″)不等。该次研究随着空间化尺度从100m到10km的逐级变化,对于只考虑空间自相关性的模型,如IDW和OK,因未考虑海拔因素,随栅格尺度的变化,误差变化很小。但是,考虑了海拔影响的模型,如IDWE、OKE和MRM,随着空间化栅格单元的增大,表现出较显著的误差增大、空间化精度减小的趋势。究其原因,是因为考虑海拔因素的情况下,栅格尺度越大,栅格之间的方差缩小,栅格面趋于平滑连续,模拟真实情况的能力下降,又因垂直递减率在0~1之间,在垂直递减率的作用下这种拟合能力下降的效应就被放大。
采用气象局提供的海拔数据和提取DEM的海拔数据进行插值,不同范围内均存在误差,尤其是四川的西北部、中北部以及东南部的气温差异较大。海拔在2 000~3 000m和大于5 000m时,气温差异达到0.6℃以上。因此,在地形复杂的地区,对气温插值时要适当考虑两种不同海拔数据源带来的差异。但从四川省插值结果栅格统计面积来看,误差绝对值大于1℃的区域占全省面积仅为1%左右,在相对宏观的数据分析时可以不考虑两种海拔数据源对空间化误差的影响。
在地貌类型多样,气候条件复杂的地区,采用分区拟合能降低空间化误差,类似的研究思路在其他社会经济数据空间化时多有体现[27]。区域地貌、气候的差异导致地区与地区间存在明显气压差异的,其垂直递减率亦不相同,这是导致分区拟合优于全域拟合的重要原因。由于历史原因,川西高原和四川盆周山区属于气象观测资料稀疏区,虽然分区拟合能提高空间化精度,但对川西地区提升有限。未来需要加密自动观测站网,在地质灾害易发区、输油管线和输电通道等重要基础设施布局区、高原牧区完善气象灾害监测能力,提升中小尺度天气变化监测能力,发挥分区拟合等不同空间化方法的优势提升插值精度。
根据上述分析,为减少空间化误差,空间化过程中宜采用空间异相关模型,如考虑了海拔因子影响的空间自相关模型,或者同时考虑了海拔、经纬度、距海岸线距离等因素的多元回归模型并用残差校正;栅格单元应符合气象观测数据所能代表的范围,一般宜采用500m(15")至2 000m(1′)的栅格单元;气象站点应具有一定的代表性,覆盖各种地貌类型,空间分布上相对均匀,在气象观测资料稀疏区加密自动观测站点。由于观测资料本身受到时空代表性有限、分辨率不足、时空覆盖范围有限及不同观测资料之间的系统偏差等问题的限制,以观测资料为基础进行多元资料的融合分析及利用数值模式进行资料的同化分析与再分析,是获得质量可靠、空间覆盖完整、分辨率高的气象数据产品的普遍做法[28]。对于四川省累年年平均气温、累年月平均气温等数据空间化产品生产加工,宜采用IDWE、OKE、MRM和Anusplin 等空间异相关模型,栅格单元的大小控制在2km以内,可以保障产品的平均绝对误差低于1℃。