基于数据同化的地下水埋深插值研究

马 欢 岳德鹏 YANG Di 于 强 张启斌 黄 元

(1.北京林业大学精准林业北京市重点实验室， 北京 100083； 2.佛罗里达大学地理系， 盖恩斯维尔 32611)

基于数据同化的地下水埋深插值研究

马 欢1岳德鹏1YANG Di2于 强1张启斌1黄 元1

(1.北京林业大学精准林业北京市重点实验室， 北京 100083； 2.佛罗里达大学地理系， 盖恩斯维尔 32611)

以西北干旱区典型县域磴口县为研究区，以2015年8月份40个地下水采样点的样品数据为基础，引入集合卡尔曼滤波(EnKF)数据同化将其优化作为主变量，以蒸散发量反演结果以及归一化植被指数(NDVI)数据为协变量，进行协同克里金插值，同时与未采用同化的协同克里金插值结果以及经同化采用普通克里金插值结果进行交叉验证。结果表明：三者在较大空间尺度上对地下水埋深空间分布趋势的模拟基本一致，南部沙漠地区整体较高，在空间分布上表现为明显的地理规律性。同化后的数据进行协同克里金插值的结果改善最显著，平均误差、均方根误差、平均标准误差均优于未同化插值结果，其中平均误差仅为0.270 5 m。与普通克里金插值方法相比，协同克里金插值考虑蒸散发与NDVI的协同作用，精度明显提高，平均误差减小0.409 7 m，均方根误差减小0.078 4 m，平均标准误差减小1.016 7 m。

地下水埋深； 干旱区； 数据同化； 协同克里金法； 插值

引言

地下水是能够影响地表生态环境系统的一个重要影响因子。地下水的变化往往会引起整个生态环境系统天然平衡状态的失衡和破坏。过量开采地下水，造成地下水水位下降，在世界范围内总的来说会诱发4方面的环境地质问题：区域性地下水降落漏斗、地面沉降、海水入侵和水质恶化等。在干旱地区浅层地下水水位大幅度下降，疏干了原有的沼泽湿地，原有的绿洲会变成沙漠，原有的景观会遭到破坏[1-2]。因此，对西北干旱地区而言，地下水作为水环境和生态环境状况的重要指标，实时掌握地下水分布变化有着重要的意义。

空间插值方法是研究地下水属性的空间变异特征和分布的有效工具之一，并已广泛地应用在地下水资源领域[3-5]。由于通过分布在各处的监测站点只能获取有限的点尺度上地下水埋深信息，空间上是间断的，不便于在区域面尺度上开展细节研究。利用空间插值方法对已有站点数据进行处理可以获得空间连续的地下水埋深数据，但插值结果受制于能获取资料站点的密度，在站点分布稀疏的区域，地下水数据具有不确定性、离散性高，很难得到精确的插值效果。在多数研究中，往往不重视进行空间数据分析、插值方法的筛选和插值参数的优化。不适合的插值方法和插值参数会造成严重扭曲，使数据具有的某种空间变化趋势难以辨识[6]。

插值方法的比较与参数优化已成为空间插值技术在各领域应用研究的热点，实际工程中许多变量之间存在一定的相关性，其中一些有较丰富的资料，另一些变量的资料比较少，有的变量很容易测量，有的变量则难以测量或者费时费力。能否利用资料多、容易测量的变量资料来提高资料少、难以测量的变量估计精度，协同克里金法就可以解决此类问题，利用多变量进行高精度估值[7-9]。

地下水埋深的空间分布特征由于存在取样数据之间自相关和互相关的特性，未知点的埋深属性可通过特定的方法内插或外推，经典的数理统计方法难以确定周围的数据样本点以及相应的插值系数。通常情况下，地下水埋深与蒸散发和NDVI数据在统计量及物理意义上存在某种联系[10-11]，由于蒸散发和NDVI这2种协同变量的获取相对简单，所以可以利用大范围蒸散发数据、NDVI及部分范围的地下水埋深数据，利用协同克里金方法估计大范围的地下水埋深。

本文以西北干旱区典型城市磴口县为研究区，以2015年8月份的地下水埋深监测数据为基础数据，采用数据同化技术中的集合卡尔曼滤波算法模拟修正地下水埋深数据，以此为主变量，以蒸散发反演结果以及NDVI作为协变量，进行小样本地下水埋深协同克里金插值，分析地下水的空间分布规律，以期为干旱区水资源均衡管理和分水政策的制定提供依据。

1 材料与方法

1.1 研究区概况

磴口县位于内蒙古巴彦淖尔市西南部(106°9′～107°10′E、40°9′～40°57′N)，地处内蒙古河套平原源头，黄河中上游，背靠狼山山脉，西邻乌兰布和沙漠。全县东西长约92 km，南北宽约65 km，总面积约为4 167 km2。境内海拔高度为1 030～2 046 m，整个地形除山区外东南高西北低，总体渐向鄂尔多斯高原倾斜。磴口县属温带大陆性季风气候，历年平均风速3.0 m/s，瞬间最大风速28 m/s，多年平均降水量143.9 mm，多年平均蒸发量2 327 mm，多年平均气温7.6℃。黄河流经磴口县52 km，年径流量310亿m3，共有水域面积2 406万m2。河套地区地下水埋深0.5～3 m，沙区地下水埋深3～10 m，山前洪积扇地下水埋深3～30 m，地表水资源较丰富，主要得益于黄河来水，但巴彦淖尔市总体水资源紧缺，时空分布差异大。

1.2 数据来源与处理

地下水埋深样点共40个，分别采集于磴口县水务局、中国林科院沙林中心、内蒙古巴彦淖尔市黄灌局设立的监测井，如图1所示。这些井多数从2006年开始观测，有较为完整的资料。本文选取2015年8月份地下水埋深数据作为空间插值的原始数据。选取磴口县2015年8月份Landsat TM影像，利用ENVI5.1软件对影像进行波段合成、图像增强和几何校正处理，通过建立解译标志与实地调查得到磴口县2015年土地利用现状图以及NDVI分布图。利用SEBAL模型反演出磴口县2015年8月份蒸散发数据。其他地表温度、水文、植被及社会经济数据来自巴彦淖尔市磴口县水土保持局、国土资源局、统计局等。选用ArcGlS 10.2作为数据统计分析和插值平台，建立地理信息数据库，将数字化底图和监测井地下水埋深数据导入数据库。

图1 磴口县地下水监测井位置分布Fig.1 Location of groundwater monitoring wells in Dengkou County

1.3 co-Kriging模型

协同克里金法(co-Kriging)建立在协同区域化变量理论基础之上，利用多个区域化变量之间的互相关性，通过建立交叉协方差函数和交叉变异函数模型，用易于观测和控制的变量对不易观测的变量进行局部估计[12]。在本研究中以不易观测的地下水变量筛选蒸散发以及NDVI作为协同变量，建立交叉变异函数，与普通克里金法相比，可更有效改进地下水数据的估计精度和采样效率。

取其中一个主要变量作为主变量，即目标变量，本文为地下水埋深数据。同时选取另外的一个或者多个变量作为次要变量，本文选取相关性较高的蒸散发数据以及NDVI数据。主变量的自相关性和其他类型变量之间的交叉相关性都有助于对结果做出更好的预测。这主要因为它不只是对各变量做自相关性预测，还包括对所有变量之间交叉相关性的估计[13]。co-Kriging插值模型步骤为：

(1)

(2)

式中ai、bj——协同克里金权重系数

(1)无偏性条件

(3)

式中μu——蒸散发的平均值μv——NDVI的平均值

(2)最优性条件

在满足无偏条件下，协同克里金估计方差为

(4)

(5)

式中μ1、μ2——拉格朗日乘子

对F求偏导数并令其为零，得协同克里金线性方程组

(6)

整理此协同克里金线性方程组为

(7)

式(7)是一个n+m+2阶线性方程组，解该方程组可得协同克里金权重系数ai、bj，再代入协同克里金线性估计式(2)，得到协同克里金线性无偏最优估计量。

根据协同克里金方程组，地下水协同克里金方差为

(8)

若有多个变量，则求解λak的协同克里金方程组为

(9)

(10)

1.4EnKF-co-Kriging模型

集合卡尔曼滤波(Ensemble Kalman filter，EnKF)通过引入观测数据利用集合卡尔曼滤波同化算法提高模型模拟精度，首先通过扰动生成模型预报的状态变量集合，通过卡尔曼滤波进行分析，同时考虑模型预测误差和观测值的误差，利用卡尔曼增益(Kalman lain)对模拟结果进行校正，以每一时刻状态变量集合的平均值为该时刻状态变量的最优估计[14-15]。通过向前积分，计算不同时刻的概率密度函数所对应的统计特性(如均值和协方差)。

本文基于EnKF方法的数据同化步骤为：

(1)初始化。给定Q个符合高斯分布的状态变量，在本文中为各测站点的地下水实测数据。

(11)

式中q——模型状态变量的个数Q——集合的个数l——时间刻度

(12)

本研究将M()引入协同克里金内插半方差函数模型，通过交叉验证进行筛选，可以更加精确的模拟地下水埋深。

(3)计算l+1时刻的卡尔曼增益矩阵Ll+1。

(13)

其中

(vl+1～N(0,R))

(14)

(15)

(16)

式中Yl+1——l+1时刻的状态变量的观测值vl+1——期望为0、方差为R的高斯白噪声

(5)判断是否到结束时刻，即目标年份。如果未到结束时刻，那么返回步骤(2)，否则结束。

1.5 精度验证

采用交叉验证方法可以对插值方法的空间插值精度进行比较，在上述研究基础上，以站点实测值为检验点，采用交叉验证的方法对比3种插值方法的插值精度。交叉检验的方法很多，包括相关系数、均值、相对误差、均方根误差等。本文将平均误差(Mean error)、均方根误差(Root mean square error)、平均标准误差(Mean standard error)作为评估插值效果的标准[16-19]。

2 结果与分析

2.1 变量空间变异

图2 蒸散发、NDVI与地下水埋深相关性分析Fig.2 Correlation analysis of evapotranspiration, NDVI and groundwater depth

地下水埋深是地下水环境要素的重要指标，能较好地反映地下水的情况。为了合理利用地下水，需对地下水埋深及其影响因素有较好的认识。植被分布情况、环境和土地利用类型等因素都会对地下水埋深产生影响。由于气候变化和人类活动导致地下水水位频繁、快速的波动，进而诱发依赖地下水的植被蒸腾量变化，是全球半干旱地区广泛存在的一种现象。在依赖地下水植被广泛分布的地区，植被蒸腾同蒸散发一样，是地下水排泄不可忽略的水均衡要素。同时，蒸散发量也是判断植被水分胁迫响应及植被对地下水依赖程度的重要指标。许多学者对此进行了研究，发现半干旱区植被覆盖状况及蒸散发对地下水变化的响应作用较大[20-22]，蒸散发数据利用SEBAL模型的反演结果与NDVI均以像元为单位，因此在整个研究区域有大量数据点作为协变量，且均匀分布于整个研究区域。

将40个地下水埋深取样点的数据与蒸散发数据、NDVI数据进行相关性分析(图2)，地下水埋深与蒸散发两者拟合的决定系数R2=0.810 6，地下水埋深与NDVI两者拟合的系数R2=0.722 5，可以看出，蒸散发、NDVI与地下水埋深均具有强相关性，研究区地下水埋深与蒸散发以及NDVI在统计量及物理意义上存在某种联系。对三者进行单样本柯尔莫戈洛夫-斯米尔诺夫(Kolmogorov-Smirnov)检验，地下水埋深与蒸散发、NDVI的检验均满足正态分布[23]。基于以上分析，可进一步利用协同克里金理论，将蒸散发数据以及NDVI数据作为本次研究的协变量，对地下水埋深进行空间分析。

利用蒸散发数据、NDVI作为本次研究的协变量，利用基于EnKF同化的40个取样点的地下水埋深优化数据作为主变量，对地下水埋深进行空间插值研究，利用较容易获取的大范围数据补充地下水埋深数据的不足。变异函数应用于分析区域化变量的变异特征及结构性状，但由于样本数量有限，通常变异函数值是由有限的实测样本值计算所得。因此，当定量地描述变量在整个区域的特征时，必须通过变异函数进行推断[24]。由于实际工作中区域化变量的变异通常很复杂，它可能在不同方向上有不同的变异性，或者在同一方向上包含着不同尺度多层次的变异性，因此在计算出实验变异函数后，常采用多种理论模型进行拟合[25]。

为了增加方法的准确性，本文采用11种半方差函数模型对地下水埋深插值的结果进行实证对比，精度评价见表1，发现稳定模型是研究区地下水埋深插值的最优模型，说明不同的指标应根据参数设置比较才能确定最优的协同克里金插值。

表1 不同半方差函数模型空间内插结果的精度评价

由于地下水埋深研究是对中尺度范围的空间内插进行分析，因此，平均误差和均方根误差是主要的精度评价指标和模型选择的依据[26]。在ArcGIS软件中的交叉验证是利用参加空间插值的数据组进行的验证，包含5个评价指标：平均误差反映估计值误差的大小，其值越接近0越好；均方根误差、平均标准误差、标准化平均误差三者反映插值模型的反演灵敏度、极值效应，其值越小越好；标准化均方根误差越接近1则越好。综合考虑平均标准误差、标准化平均误差、标准化均方根误差，交叉验证分析数据表明，地下水埋深的半方差函数模型在选择稳定模型进行时内插精度最高，其经过拟合后的稳定半方差数学方程表示为：

地下水埋深-地下水埋深模型

h=0.001 114 3C0+1.114 3S(0.136 44,2)

地下水埋深-蒸散发模型

h=-0.418 61S(0.136 44,2)

蒸散发-蒸散发模型

h=0.001 045 3C0+1.045 3S(0.136 44,2)

式中h——半方差函数值S——稳定模型半方差数学方程C0——块金值，反映最小抽样尺度以下变量的变异性及测量误差

模型结果见图3。

图3 稳定模型模拟协方差分析Fig.3 Covariance analysis of stable model

2.2 地下水埋深插值

图4 2015年8月磴口县地下水埋深空间分布图(单位：m)Fig.4 Spatial distributions of groundwater depth of Dengkou County in August 2015

利用2015年8月1日—8月30日的地下水埋深数据作为同化的初始状态变量，同化间隔为24 h，并只在表层加误差(包括初始误差、模型误差、观测误差)，其他各层误差均设为0。利用EnKF融合协同克里金插值模型，得到优化的地下水埋深数据，以此作为插值主变量，以蒸散发数据与NDVI数据为协变量，对磴口县2015年地下水埋深数据进行空间插值，利用Jenks分类方法重分成8类得到2015年地下水埋深空间分布结果(图4a)。同时分别将未利用数据同化的地下水数据采用协同克里金插值模型的模拟结果(图4b)、经同化采用普通克里金插值模型的模拟结果(图4c)，与协同克里金插值模型结果(图4a)进行对比，并从物理分布规律和统计学误差方面进行分析。

各插值结果可以直观地反映该地区地下水因子的整体连续变化情况，总体上来看，三者在较大空间尺度上对地下水埋深空间分布趋势的模拟基本一致，大范围地下水埋深数据均比较符合一般意义上的生态埋深标准范围(0～10 m)，河套地区地下水埋深0.5～3 m，沙区地下水埋深3～10 m。地下水埋深最高值区域出现在南部沙漠地区，其值均大于5 m；最低值区域出现在东部黄河流域附近区域，其值范围为0.5～2 m，流域附近水资源较为丰富，地下水埋深较浅；东南部磴口县城区地下水埋深较大，表明随着城市化发展，人类活动用水较多，导致地下水埋深较深。从图4中可以看出3种方法插值结果的主要差别在磴口县西北部，EnKF-co-Kriging模型的结果中西北部的地下水埋深最大，最高值达到了10.83 m，西北部山前监测点较少，由于EnKF-co-Kriging模型考虑NDVI因素，山前植被覆盖率较高，蒸散发作用较强，使得其插值结果中地下水埋深较深。通过小样本插值后可更加了解地下水埋深分布范围和分布面积，对比离散型监测数据可直观地反映区域地下水分布变化规律。

从较小的空间尺度观察插值结果，发现EnKF-co-Kriging生成的插值结果图更为平滑，反映出该方法具有较高的耐抗性，能使等值线的平滑度与精确度之间达到较好的平衡，利于总体趋势与局部趋势的良好表现。而另外2种方法产生的等值线会围绕埋深极大值或极小值采样点发生大曲率的弯曲，甚至闭合形成小等值线圈，插值结果较为破碎，说明其对极端值敏感，耐抗性低，对总体趋势的反映产生较大干扰。

2.3 交叉验证

为了比较上述3种方法的插值精度，需要分别计算地下水埋深40个样点插值估算值和实测值之间的平均误差、均方根误差和平均标准误差，结果见表2。

表2 3种插值方法交叉验证结果

从表2中可以看出，在平均误差和均方根误差方面，EnKF-co-Kriging方法均小于co-Kriging插值方法，在进行数据同化后，2种插值方法的平均误差和均方根误差都有所减小，可见数据同化方法可以大大提高插值精度；同时，EnKF-co-Kriging方法的平均误差和均方根误差均小于EnKF-Kriging插值方法，说明蒸散发、NDVI的协同作用对地下水埋深影响较大，对地下水埋深数据的插值精度来说，考虑蒸散发和NDVI的协同作用进行插值在精度上要优于普通的克里金插值法。因此，基于EnKF的协同克里金方法整体要优于另外2种插值方法得到的结果。

综上所述，在采取一系列插值优化措施后，3种插值方法均能在较大空间尺度上，比较客观地模拟出地下水埋深的空间分布趋势，在较小空间尺度上，EnKF-co-Kriging对局部趋势的模拟最为准确和合理。结合交叉验证的结果综合评定，EnKF-co-Kriging是最适合研究区地下水埋深的插值方法。

3 结论

(1) 磴口县地下水埋深数据利用集合卡尔曼滤波进行同化，将同化的修正值作为协同克里金插值的主变量，依据蒸散发数据以及NDVI数据的协同作用进行协同克里金插值，得到地下水的空间分布变化。同化后的数据进行协同克里金插值的模拟结果改善最显著，平均误差小于未同化协同克里金插值0.057 9 m，均方根误差减小0.048 1 m，平均标准误差减小0.461 6 m。

(2) 地下水埋深的空间分布特征由于存在取样数据之间自相关和互相关的特性，利用地下水埋深与蒸散发以及NDVI在统计量及物理意义上存在的某种联系，将蒸散发数据及NDVI数据作为协变量进行协同克里金插值，可以弥补主变量(地下水埋深)样本点少、变异方向不明显等缺点，与普通克里金插值方法相比，协同克里金插值的精度明显提高，平均误差减小0.409 7 m，均方根误差减小0.078 4 m，平均标准误差减小1.016 7 m。

(3) 运用协同克里金插值方法对地下水数据分析的结果表明，三者在较大空间尺度上对地下水埋深空间分布趋势的模拟基本一致，大范围地下水埋深数据均比较符合一般意义上的生态埋深标准范围(0～10 m)，空间分布规律与实地情况基本一致，南部沙漠地区整体较高，在空间分布上表现为明显的地理规律性。

1 郑丽萍. 地下水资源的保护与管理[J].农业与技术，2010(3)：58-61. ZHENG Liping. Protection and management of groundwater resources[J].Agriculture and Technology,2010(3):58-61.(in Chinese)

2 樊自立，陈亚宁，李和平，等. 中国西北干旱区生态地下水埋深适宜深度的确定[J].干旱区资源与环境，2008,22(2)：1-5. FAN Zili, CHEN Yaning, LI Heping, et al. Determination of suitable ecological groundwater depth in arid areas in northwest part of China[J]. Journal of Arid Land Resources and Environment,2008,22(2):1-5.(in Chinese)

3 蒋庆. 地下水时空变化及监测网多目标优化研究[D].武汉：华中科技大学，2008. JIANG Qing. Research on spatio-temporal changes of the groundwater and monitoring network multi-objective optimization[D]. Wuhan:Huazhong University of Science and Technology,2008.(in Chinese)

4 王春梅，孟庆岩，占玉林，等. 西北干旱区农田春小麦蒸散量的空间插值方法[J/OL].农业机械学报，2014，45(11)：172-179.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20141127&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2014.11.027. WANG Chunmei, MENG Qingyan, ZHAN Yulin, et al. Spatial interpolation analysis of spring wheat evapotranspiration measurements in Northwest China[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2014,45(11):172-179. (in Chinese)

5 何红艳，郭志华，肖文发. 降水空间插值技术的研究进展[J].生态学杂志，2005,24(10)：1187-1191. HE Hongyan, GUO Zhihua, XIAO Wenfa. Review on spatial interpolation techniques of rainfall[J]. Chinese Journal of Ecology,2005,24(10):1187-1191.(in Chinese)

6 文雯，周宝同，汪亚峰，等. 基于辅助环境变量的土壤有机碳空间插值——以黄土丘陵区小流域为例[J].生态学报，2013,33(19)：6389-6397. WEN Wen, ZHOU Baotong, WANG Yafeng, et al. Soil organic carbon interpolation based on auxiliary environmental covariates:a case study at small watershed scale in loess hilly region[J]. Acta Ecologica Sinica,2013,33(19):6389-6397.(in Chinese)

7 王长虹，朱合华，钱七虎. 克里金算法与多重分形理论在岩土参数随机场分析中的应用[J].岩土力学，2014，35(增刊2)：386-392. WANG Changhong, ZHU Hehua, QIAN Qihu. Application of Kriging methods and multi-fractal theory to estimate of geotechnical parameters spatial distribution[J].Rock and Soil Mechanics, 2014,35(Supp.2):386-392.(in Chinese)

8 DAI H Y, REN L Y, WANG M, et al. Water distribution extracted from mining subsidence area using Kriging interpolation algorithm[J]. Transactions of Nonferrous Metals Society of China,2011，21(Supp.3):723-726.

9 PANG Su, LI Tingxuan, WANG Yongdong, et al. Spatial interpolation and sample size optimization for soil copper(Cu) investigation in cropland soil at county scale using CoKriging[J]. Agricultural Sciences in China,2009,8(11):1369-1377.

10 郝兴明，陈亚宁，李卫红，等. 塔里木河中下游荒漠河岸林植被对地下水埋深变化的响应[J].地理学报，2008，63(11)：1123-1130. HAO Xingming, CHEN Yaning, LI Weihong, et al. Response of desert riparian forest vegetation to groundwater depth changes in the middle and lower Tarim River[J]. Acta Geographica Sinica,2008,63(11):1123-1130.(in Chinese)

11 赵文智，刘鹄. 荒漠区植被对地下水埋深响应研究进展[J].生态学报，2006，26(8)：2702-2708. ZHAO Wenzhi, LIU Hu. Recent advances in desert vegetation response to groundwater table changes[J]. Acta Ecologica Sinica,2006,26(8):2702-2708. (in Chinese)

12 于正军，董冬冬，宋维琪，等. 相带控制下协克里金方法孔隙度预测[J].地球物理学进展，2012，27(4)：1581-1587. YU Zhengjun, DONG Dongdong, SONG Weiqi, et al. Porosity prediction with co-Kriging method controlled by sedimentary fades[J]. Progress in Geophysics,2012,27(4):1581-1587.(in Chinese)

13 俞超锋. 基于线性矩的协克里金插值在降雨频率分析中的应用[C]∥中国水利学会2013学术年会论文集——S5河口治理与保护，2013：1773-1780.

14 刘昭，周艳莲，居为民，等. 基于集合卡尔曼滤波同化方法的农田土壤水分模拟[J].应用生态学报，2011，22(11)：2943-2953. LIU Zhao, ZHOU Yanlian, JU Weimin, et al.Simulation of cropland soil moisture based on an ensemble Kalman filter[J]. Chinese Journal of Applied Ecology,2011,22(11):2943-2953. (in Chinese)

15 于强，岳德鹏，YANG Di，等. 基于EnKF-MCRP模型的生态用地扩张模拟研究[J/OL]. 农业机械学报，2016，47(9)：285-293.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20160939&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2016.09.039. YU Qiang, YUE Depeng, YANG Di, et al. Simulation on ecological land use expansion based on EnKF-MCRP model[J/OL]. Transactions of the Chinese Society for Agricultural Machinery,2016,47(9):285-293. (in Chinese)

16 田雷，董德明，魏强，等. 3种空间插值方法在道路尘中Pb监测数据统计处理中的应用比较[J].吉林大学学报:理学版，2011，49(5)：964-968. TIAN Lei, DONG Deming, WEI Qiang, et al. Comparison of three spatial interpolation methods statistical processing of lead monitoring data in road dust[J].Journal of Jilin University:Science Edition, 2011,49(5):964-968. (in Chinese)

17 郑小波，罗宇翔，于飞，等.西南复杂山地农业气候要素空间插值方法比较[J].中国农业气象，2008，29(4)：458-462. ZHENG Xiaobo, LUO Yuxiang, YU Fei, et al. Comparisons of spatial interpolation methods for agro-climinate factors in complex mountain areas of Southwest China[J]. Chinese Journal of Agrometeorology,2008,29(4):458-462. (in Chinese)

18 李巍，范文义，毛学刚，等. 降雨量空间插值方法比较研究[J].安徽农业科学，2014，42(12)：3667-3669. LI Wei, FAN Wenyi, MAO Xuegang, et al. A comparative study of spatial interpolation methods for precipitation[J]. Journal of Anhui Agricultural Sciences,2014,42(12):3667-3669. (in Chinese)

19 云安萍，鞠正山，胡克林，等. 基于距离反比法的土壤盐分三维空间插值研究[J/OL].农业机械学报，2015，46(12)：148-156.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20151221&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2015.12.021. YUN Anping, JU Zhengshan, HU Kelin, et al. Three limensional spatial interpolation of soil salinity based on inverse distance weighting method[J/OL]. Transactions of the Chinese Society for Agricultural Machinery,2015,46(12):148-156. (in Chinese)

20 王利书，悦琳琳，唐泽军，等. 气候变化和农业发展对石羊河流域地下水位的影响[J/OL].农业机械学报，2014，45(1)：121-128.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20140120&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2014.01.020. WANG Lishu, YUE Linlin, TANG Zejun, et al. Influence of climate change and agricultural development on groundwater level in Shiyang River Basin[J/OL].Transactions of the Chinese Society for Agricultural Machinery,2014,45(1):121-128. (in Chinese)

21 黄金廷. 半干旱区蒸散发对地下水变化响应机制研究[D].西安：长安大学，2013. HUANG Jinting. The responses of evapotranpiration to the groundwater changes in the semi-arid area[D]. Xi’an: Chang’an University,2013. (in Chinese)

22 马玉蕾，王德，刘俊民，等. 黄河三角洲典型植被与地下水埋深和土壤盐分的关系[J].应用生态学报，2013，24(9)：2423-2430. MA Yulei, WANG De, LIU Junmin, et al. Relationships between typical vegetations, soil salinity and groundwater depth in the Yellow River Delta of China [J]. Chinese Journal of Applied Ecology,2013,24(9):2423-2430. (in Chinese)

23 徐驰，曾文治，黄介生，等. 基于高光谱与协同克里金的土壤耕作层含水率反演[J].农业工程学报，2014，30(13)：94-103. XU Chi, ZENG Wenzhi, HUANG Jiesheng,et al. Modelling soil water content in plow layer using co-Kriging interpolation based on hyperspectral data[J].Transactions of the CSAE,2014,30(13):94-103. (in Chinese)

24 叶回春，黄珊瑜，张世文，等. 土壤有机碳空间变异性对采样密度的响应研究[J/OL].农业机械学报，2014，45(12)：215-223.http:∥www.j-csam.org/jcsam/ch/reader/view_abstract.aspx?flag=1&file_no=20141232&journal_id=jcsam.DOI:10.6041/j.issn.1000-1298.2014.12.032. YE Huichun, HUANG Shanyu, ZHANG Shiwen, et al. Spatial variability response of soil organic carbon to sampling density change[J/OL]. Transactions of the Chinese Society for Agricultural Machinery,2014,45(12):215-223. (in Chinese)

25 张朔，鲁学军，赵英俊，等. 空间变异函数等效应椭圆套合方法及其应用[J].地球信息科学学报，2009，11(3)：342-348.ZHANG Shuo, LU Xuejun, ZHAO Yingjun, et al. Equbalent ellihtic nested structure and its application in variogram nesting[J]. Journal of Geo-information Science,2009,11(3):342-348. (in Chinese)

26 王文，刘永伟，寇小华，等. 基于集合卡尔曼滤波和HYDRUS-1D模型的土壤剖面含水量同化试验[J].水利学报，2012，43(11)：1302-1311. WANG Wen, LIU Yongwei, KOU Xiaohua, et al. EnKF and HYDRUS-1D based data assimilation experiments for improving soil moisture profile prediction[J]. Journal of Hydraulic Engineering,2012,43(11):1302-1311. (in Chinese)

Interpolation of Groundwater Depth Based on Data Assimilation

MA Huan1YUE Depeng1YANG Di2YU Qiang1ZHANG Qibin1HUANG Yuan1

(1.BeijingKeyLaboratoryofPrecisionForestry,BeijingForestryUniversity,Beijing100083,China2.DepartmentofGeography,UniversityofFlorida,GainesvilleFL32611,USA)

Groundwater monitoring is limited by practical conditions, and only limited monitoring results can be obtained when it is observed. As a kind of geostatistical interpolation method, cooperative Kriging (co-Kriging) method can effectively represent the transformation of discrete point-like information to planar continuous information. Dengkou County, a typical county in the arid region of Northwest China, was selected as the study area. The sampled data from 40 groundwater sampling sites in 2015 was selected as the main variable. And this data optimized by EnKF was used as the basic data of co-Kriging interpolation. The evapotranspiration results and NDVI data were selected as the covariates. Co-Kriging interpolation was carried out by using the sampled data from 40 groundwater sampling sites in August, 2015, as the main variable, which were optimized by EnKF, and the evapotranspiration results and NDVI data were used as the covariates. Meanwhile, the results of co-Kriging interpolation without using EnKF model and Kriging interpolation optimized by EnKF model were used to verify the accuracy. The results showed that the spatial distribution trend of groundwater depth was basically the same at large scale, the value in the southern desert region was higher, and the spatial distribution showed obvious geography regularity. The most significant improvement was achieved with EnKF model. Based on this improvement, the mean error, root mean square error and mean standard error were all better than those without assimilation, with the mean error of 0.270 5 m. Compared with the ordinary Kriging interpolation method, co-Kriging model took the synergistic effect of evapotranspiration and NDVI into consideration, and the precision was obviously improved. The mean error was decreased by 0.409 7 m, the root mean square error was decreased by 0.078 4 m and the mean standard error was decreased by 1.016 7 m. This study can provide a scientific basis for spatial visualization simulation and reasonable management of water resources in arid areas.

groundwater depth; arid region; data assimilation; co-Kriging; interpolation

10.6041/j.issn.1000-1298.2017.04.027

2017-01-11

2017-02-11

国家自然科学基金项目(41371189)和“十二五”国家科技支撑计划项目(2012BAD16B00)

马欢(1992—)，女，博士生，主要从事3S技术在生态环境中的应用研究，E-mail： xiaohuan27@163.com

岳德鹏(1963—)，男，教授，博士生导师，主要从事景观生态学和土地评价研究，E-mail： yuedepeng@126.com

K903

A

1000-1298(2017)04-0206-09