EN中文

地球信息科学学报 >

2018 , Vol. 20 >Issue 4: 505 - 514

DOI: https://doi.org/10.12082/dqxxkx.2018.170592

地球信息科学理论与方法

面向稀疏散布数据集的时空Kriging优化

  • 杨明远 , * ,
  • 刘海砚 ,
  • 季晓林 ,
  • 郭文月 ,
  • 陈思雯
展开
  • 信息工程大学,郑州 450001

作者简介：杨明远（1992-）,男,硕士生,研究方向为数字地图制图。E-mail:

收稿日期: 2017-12-07

  要求修回日期: 2018-02-04

  网络出版日期: 2018-04-20

基金资助

国家自然科学基金项目（41501446）

地理信息工程国家重点实验室开放基金资助项目（SKLGIE2015-M-4-3）

收起

An Improved Spatio-temporal Kriging Algorithm to Sparse Scattered Dataset

  • YANG Mingyuan , * ,
  • LIU Haiyan ,
  • JI Xiaolin ,
  • GUO Wenyue ,
  • CHEN Siwen
Expand
  • Information Engineering University, Zhengzhou 450001, China
*Corresponding author: YANG Mingyuan, E-mail:

Received date: 2017-12-07

  Request revised date: 2018-02-04

  Online published: 2018-04-20

Supported by

National Natural Science Foundation of China, No.41501446

Open Research Fund Program of State Key Laboratory of Geoinformation Engineering, No.SKLGIE2015-M-4-3.

Copyright

《地球信息科学学报》编辑部 所有
Fold

摘要

时空Kriging法通过将变异函数向时空域进行扩展得到时空变异函数,有效地利用时空邻近的采样点综合进行插值,由于时空稀疏散布数据集具有单一时刻下样本点数量少以及时空分布不规律的特点,难以满足使用时空Kriging插值法的基本条件,导致插值精度不高,据此本文提出了优化方法：通过多时段叠置拟合空间变异函数的方法,综合利用时空邻域内的采样点以解决单一时刻下空间邻域内数量不足情况;控制时间变异对空间变异函数拟合的误差影响;采用积合式模型构建时空变异函数进行插值。最后使用Argo海温数据进行插值实验,在相同条件下与时空Kriging法以及时空权重法的交叉验证结果对比得出,该方法在保证拟合所需采样点数量要求的同时,有效削减了一般时空Kriging法中时间变异对空间变异函数拟合结果的干扰,插值结果的绝对误差均值从0.5降低至0.2以内,稳定性进一步增强,改善了时空Kriging法在稀疏散布数据条件下精度上的不足。

关键词： 克里金; 时空插值; 变异函数; 空间探针; 时间切片

本文引用格式

杨明远 , 刘海砚 , 季晓林 , 郭文月 , 陈思雯 . 面向稀疏散布数据集的时空Kriging优化[J]. 地球信息科学学报, 2018 , 20(4) : 505 -514 . DOI: 10.12082/dqxxkx.2018.170592

Abstract

Spatio-temporal Kriging is efficiently used in interpolation by adjacent sampling point in space-time. The core is to extend the spatial variogram into space time. Because sparse scattered dataset is lack of sampling point in single time slice and the distribution of point is non-uniform, we propose an improved method of spatio-temporal Kriging against low precision. First, the trend surface of the sample is obtained by using cubic polynomial. The sample data is decomposed into trend terms and residual terms, because original trending sample data can't satisfy the stationary assumption required for Kriging interpolation. Then, the time variogram is fitted with less data coming from a nearby stationary station. The sampling position of the stationary station is constant and the sampling frequency is consistent. It is suitable for fitting the variogram because its observation sequence is longer, although few in quantity. Meanwhile,instead of fitting method by dataset in all-time, we adopt the strategy of multi-period overlap fitting to obtain a more reasonable spatial variogram. The length of the time segment is selected according to the degree of time variation. The variable values of the sampling point in each sub-period are calculated and overlaid to fit spatial variogram. In such way, the spatio-temporal variogram is constructed based on product-sum model, which is used to estimate the variable value in space and time. In final phrase, interpolation is performed using the spatio-temporal weights solved by Kriging equations. To verify the effectiveness of the proposed method, a comparison with the existing interpolation method is made by sea temperature data of Argo buoys from China Argo real-time data center and moored buoys from Pacific Marine Environmental Laboratory. From the comparison of the cross-verification results of interpolation, we judge the accuracy and stability of the method based on MAE and MSE. Compared to general spatio-temporal kriging and spatio-temporal weight interpolation, the proposed method is increased by 69.5% and 38.9% respectively in accuracy, and increased by 61.9% and 48.9% respectively in stability. The proposed method is improved based on spatio-temporal Kriging, considering the structural characteristics and spatial and temporal variation characteristics of sparse scattered datasets. Spatio-temporal variogram is more scientific and practical, which is constructed through the proposed method. Interpolation precision and stability are also improved significantly.

Key words： Kriging; spatio-temporal interpolation; variogram; space probes; time slices

1 引言

当前地球科学研究热点中,大量无法通过卫星遥感获取的信息仍需借助其他实地观测设备进行采样,由于受到经济成本、自然环境的影响,采样结果的空间分布往往呈现出离散、不规则的特征,不足以满足实际应用的要求,需要采用合适的空间插值方法对采样数据进行插值以获得连续、均匀、利于分析的网格化资料。目前主要的插值方法包括：距离反比加权法、样条函数法、趋势面法和Kriging法。
克里金（Kriging）插值法又称自协方差最优插值法,是一种基于地统计学的插值方法。该方法以区域化变量理论为基础,以变异函数模型为主要工具,在对采样点的空间自相关分析基础上,依据自然现象的空间变异规律进行插值,从而可以得到无偏最优估计量。随着对诸多自然领域（如气象、海洋、土壤）随机过程的不断探索,人们发现除了分析变量的空间变异规律以外,还应充分考虑其时间变异规律[1,2],由此时空Kriging法应运而生。
时空Kriging法通过将变异函数向时空域进行扩展得到时空变异函数,有效地利用时空邻近的采样点综合进行插值,弥补了单一时刻空间点数量不足带来的插值精度问题,实现了时空离散点到连续体的映射。通常时空理论变异模型采用分离型和非分离型2类[3,4]。Mitchell等[5]都提出了分离模型,该模型一般表示为纯空间和纯时间变异函数的乘积,其特点表现在创建过程简单,计算效率较高,但无法达到较为理想的实验条件,损失了区域化变量的时空交互性。Cressie等[6,7]提出了非分离模型,该模型是将纯空间和纯时间变异函数模型通过乘积、线性组合、积分变换而来。李莎等[8]、胡丹桂[9]等对气温变化分别通过积合式和积分式对时空变异函数建模进行插值;王建民等[10]和陈鼎新等[11]率先将时空Kriging法引入变形监测和地磁场2类数据精度要求较高的领域并针对数据特征加以改进,实现插值的精度提升;文献[12]、[13]中对近期普遍关注的PM2.5预测和分析问题用时空Kriging法进行了研究,并取得较好的结果。
在上述应用场景中,观测站采样在空间维上位置相对固定且分布均匀,在时间维上采样时间同步、频率一致,样本数量也能得到保证,样本数据集在时空域中呈现密集稳定的特征。但对于一些移动观测设备（如高空气球、大洋浮标等）的采样过程则难以达到上述理想条件,该类采样设备的观测位置不固定、时刻不同步、频率不一致,导致样本数据集整体呈现稀疏散布的特征。图1为2类时空数据集对比示意图。在应用时空Kriging法时,由于单一时刻下的样本点数量少且分布不均造成的空间变异函数拟合误差较大,以及由于位置不断移动造成的时间变异函数无法获得,最终导致整体时空变异函数难以构建或误差较大。
Fig. 1 Comparison of two kinds of dataset

图1 2类数据对比示意图

本文重点研究时空Kriging法在稀疏散布的时空数据集条件下的改进方法,提出多时段叠置拟合方法。该方法采用时间间隔较短的子时段划分来计算空间变异函数值的方式,以降低时间变异对空间变异函数拟合结果的影响,进一步提升时空Kriging法的插值精度。通过以Argo海洋浮标温度剖面数据集为例进行插值实验,并与时空Kriging插值法[14]和海温重构常用的时空权重插值法[15]的交叉验证结果进行对比分析。

2 时空Kriging算法

Kriging插值法由南非矿产工程师Kriging首先提出,用于寻找金矿。该方法建立在变异函数理论及结构分析基础上,以区域化变量作为研究对象,在有限区域内对区域化变量进行线性无偏最优估计[16]。时空Kriging是在普通Kriging法的基础上,由空间域扩展到时空域,通过建立时空变异函数提供结构信息,对待估点进行插值。

2.1 时空变异函数

变异函数作为Kriging插值的基本工具,能够反映区域化变量的空间变化特征,特别是透过随机性反映区域化变量的结构性。设区域化变量 Z x ∈ S 二阶平稳的条件下,定义空间变异函数的公式为：
γ h = 1 2 Var Z x - Z x + h = 1 2 E Z x - Z x + h 2 （1）
式中： x 表示样本点在空间域内的坐标; h 为空间距离。
扩展至时空域后,设时空区域化变量 Z s , t ∈ S × T 满足二阶平稳假设, S 表示空间域, T 表示时间域,时空变异函数可以被定义为：
γ h s , h t = 1 2 Var Z s + h s , t + h t - Z s , t = 1 2 E Z s + h s , t + h t - Z s , t 2 （2）
式中： s , t 表示样本点在时空域中的坐标; h s 和 h t 分别为空间和时间距离。
为了定量的描述整个区域化变量的特征,需要通过计算得到的实验变异函数值,并选择合适的理论模型来分别拟合出最优的时间和空间变异函数。拟合方法通常采用最小二乘拟合[17]和加权回归法拟合[18]。利用分离型或非分离型模型构建时空变异函数,拟合得到时间和空间变异函数模型。因此,解决时空散布数据的Kriging插值问题的关键在于如何获得合理的时空变异函数。

2.2 时空Kriging插值

通过引入时空维的相关性和随机性构建有效的时空变异函数,再将空间Kriging插值公式扩展到时空域得到时空Kriging插值公式。
Z * s , t = ∑ i = 1 n λ i Z s i , t i （3）
式中： s i , t i 表示某一样本点 i 在时空域中的坐标。
待求权重系数 λ i i = 1 , 2 , ⋯ , n 需满足 Z * s , t 为 Z s i , t i 的无偏估计量,且估计方差最小。在二阶平稳假设条件下,通过时空Kriging方程组计算可得：
∑ j = 1 n λ j γ s i , t i , s j , t j - μ = γ s i , t i , s , t ∑ i = 1 n λ i = 1 （4）

3 时空Kriging优化方法

时空Kriging插值法的关键在于时空变异函数的拟合结果,而这主要取决于纯空间和纯时间变异函数的拟合。针对稀疏散布的时空数据存在的空间和时间变异函数难以拟合的现状,本文对该方法进行了改进：定义 Z t 为时空范围 S × T （ s ∈ S ） 内空间点 s 处的空间探针,即 Z s = Z s , t 1 , Z s , t 2 , ⋯ , Z s , t n ;定义 Z t 为时空范围 S × T （ t ∈ T ） 内 t 时刻的时间切片,即 Z t = Z ( s 1 , t ) , Z ( s 2 , t ) , ⋯ , Z ( s n , t ) ,如图2所示。
Fig. 2 Space probe and time slice

图2 空间探针与时间切片

3.1 时间变异函数拟合

时空稀疏散布数据通常由移动采样设备采样生成,造成同一位置的观测值时间序列无法得到,进一步导致时间变异函数难以构建。而固定站具有采样位置不变且频率一致的特点,尽管站点少但便于获得长时间观测序列。因此,本文借助附近少量 k 个固定观测站的空间探针分别进行时间维实验变异函数值的估计。某固定观测站 s 处的时间实验变异函数值计算公式为：
γ t * s , h t = 1 2 n ∑ i = 1 n Z s , t i - Z s , t i + h t 2 （5）
由于固定站数量较少,估值点位置处的时间变异情况与各固定站的变异情况受空间距离影响不同,因此本文根据计算结果分别拟合出 k 个固定观测站位置处的时间变异函数 γ t s j , h t ( j = 1 , 2 , … , k ) ,将每个时间变异函数的块金指数 C 0 、拱高 C 、变程 a ,通过距离反比加权估计得到估值点处的时间变异函数 γ t ( h t ) 。 计算公式如下：
γ t h t = 1 k ∑ j = 1 k λ j γ t s j , h t （6）
式中： λ j 为固定站到估值点的距离权重系数。

3.2 多时段叠置的空间变异函数拟合

实际工作中,采样点对数目不可能无限多,通常要求变程以内各距离上点对数目不应少于20对,以保证变程范围内实验变异函数值能准确地反映区域化变量的空间变异性。移动观测设备的采样时刻和更新频率不一致,造成单一时间切片内的样本数量过少,拟合得到的空间变异函数不具有代表性,无法描述当前区域的空间变异特征。而采用全时段内的采样点拟合空间变异函数则无法排除时间变异带来的误差影响,因此本文采用多时段叠置的方法进行变异函数的拟合。单个子时段空间实验变异函数值计算公式为：
γ s * h s = 1 2 mn ∑ j = 0 m ∑ i = 1 n Z s i , t - Z s i + h s , t ± j 2 （7）
式（7）将原有求单一时间片的空间实验变异函数值沿着时间维向两侧做适度延展（延展的间隔 m 依据时间变异函数确定,通常不大于3）,采用间隔较小的子时段内的样本点分别进行空间实验变异函数值的计算,将各子时段内的实验变异函数值计算结果叠置后共同进行拟合得到 γ s ( h s ) 。较短的子时段内,采样点之间实验变异函数值很大程度地排除了时间变异的干扰,多个子时段的计算结果叠置后保证了拟合空间变异函数所需的数量,同时避免了全时段内样本点之间的时间变异程度过大对变异函数拟合精度的影响。

3.3 时空变异函数构建与插值

将3.1节和3.2节中得到空间和时间变异函数采用积合式模型构建时空变异函数,积合式是一种常见的非分离模型,其构建方法如下：
γ h s , h t = k 1 C t 0 + k 2 γ s h s + k 1 C s 0 + k 3 γ t h t - k 1 γ s h s γ t h t （8）
k 1 = C s 0 + C t 0 - C st 0 , 0 / C s 0 C t 0 k 2 = C st 0 , 0 - C t 0 / C s 0 k 3 = C st 0 , 0 - C s 0 / C t 0 （9）
式中： γ h s , h t 为时空变异函数; C st 0 , 0 、 C s 0 、 C t 0 对应时空、空间和时间变异函数的基台值。
将时空变异函数的计算结果出代入到Kriging方程组式（4）求解时空采样点的时空权重 λ i ,最终通过时空Kriging插值公式（3）计算得到插值结果。

4 实验与分析

本文采用实际的Argo海洋观测浮标数据对改进后的时空Kriging算法的效果进行检验,实验分为时空变异函数建模和优化时空Kriging插值精度对比2部分。

4.1 实验数据

Argo（Array for real-time geostrophic oceanography）计划通过投放浮标以获得立体、广域、实时、多要素综合集成的水文剖面数据。浮标随海流自由漂浮,每个Argo浮标隔几天获取一条0~2000 m深度范围内的温盐剖面数据,由于投放批次不同,导致获取信息的频率不一致和获取时间不同步,其数据集是一种典型的时空稀疏散布数据。为了验证本文方法的有效性和稳定性,选择不同区域不同季节的Argo温度数据进行实验。如图3所示实验区A选择西太平洋（140°~160° E,10°~30° N）,时间范围为2015年6-7月。实验区B选择西太平洋（150°~170° E,0°~20° N）,时间范围为2015年1-2月。实验数据集来自中国Argo实时资料中心。
Fig. 3 Test region

图3 试验区域

TAO（Tropical Atmosphere Ocean）锚式浮标同为海洋观测浮标,具有位置固定的特点,能够锚定位置并获得完整的时间观测序列,但规模远不及Argo浮标,在实验中作为辅助数据用于获得时间变异函数。

4.2 数据预处理

（1）垂直方向的温度数据标准化
Argo浮标数据和TAO浮标数据在垂直方向均是离散、不统一的,需要将所有剖面数据的观测值插值到统一标准层。由于浅层的海水受风浪和日照影响变化幅度较大,而深层的海水则相对稳定,为减少数据量又不丢失细节因此采取不等距间隔的设置方式。采用Akima法可将离散的数据点拟合形成曲线,相比于3次样条函数拟合的曲线更为光滑、自然。拟合插值结果如图4所示。
Fig. 4 Vertical interpolation results of Akima

图4 垂直方向Akima插值结果

（2）趋势项剔除
区域化变量的二阶平稳是影响Kriging插值精度的重要因素,而海温的时空分布较不平稳。本文使用的Argo剖面数据在空间上采样点之间跨度很大,同一天的平均分辨率仅为3°。图5为20 m水深的海温散点图,可以看出水下温度在空间上分布呈现明显的低纬度地区较高的趋势,这影响了样本分布的平稳性,进一步导致空间变异函数的拟合过程易产生孔穴效应。本文在插值前采用局部趋势面分析法,对估值点周围的时空数据进行趋势面拟合。将空间样本数据分解为趋势量和剩余量,对剩余量进行Krging插值后,将插值结果与趋势量相加得到最终结果。该类方法也称为残差Krging法。
Fig. 5 Trend elimination of temperature

图5 趋势面剔除

4.3 结果与分析

4.3.1 实验1：时空变异函数建模
时空变异函数的构建是时空Kriging插值的关键步骤,建模的结果将直接影响到最终的插值精度。本节将采用本文提出的多时段叠置拟合和一般时空Kriging的全时段拟合2种方式来构建时空变异函数,并对比分析2种方式构建的变异函数在结构上存在的差异。实验分为时间变异函数拟合、空间变异函数拟合,以及利用积合式模型对时空变异函数进行建模3部分。可选择的用于拟合的理论模型包括球状模型、指数模型、高斯模型,实践证明95%以上的变异函数散点图都可以用球状模型拟合,实验中也采用该模型。
（1）时间变异函数拟合
TAO锚式浮标具有固定位置,且观测周期始终为1 d,便于用于拟合时间变异函数。首先,搜索待估点附近的锚式浮标数据,时间范围为 t - 20 , t + 20 （ t 为估值日期,单位为天）。其次,以范围内的第一个TAO浮标 s 1 为例,利用式（5）计算其实验变异函数值（半方差值）得到变异函数云图（图6（a））,并对各时间间隔下的实验变异函数值求平均后采用球状模型进行拟合,拟合后的时间变异函数结果（图6（b））。最后,将空间范围内全部锚式浮标对应的时间变异函数 γ t s i , h t 参数加权平均得到最终的时间变异函数式（（10））。
γ t h t = 0.06 h t = 0 0.81 × 3 h t 2 × 22 - h t 3 2 × 22 3 + 0.06 0 < h t < 22 0.87 h t > 22 （10）
Fig. 6 Fitting of time variogram

图6 时间变异函数拟合

（2）空间变异函数拟合
采用4.2节得到的剔除趋势后Argo自由浮标数据,样本集为空间范围为 s - 10 , s + 10 （ s 为估值点位置,单位为°）,时间范围为 t - 20 , t + 20 （ t 为估值点位置,单位为d）。实验分别用单一子时段、全时段和多时段叠置3种方法对空间实验变异函数值计算,每个区间内的点取平均值便于拟合空间变异函数,并对比拟合结果。其中,基于球面模型的地理空间距离计算采用Haversine公式,结果如下：① 单一子时段内计算得到的实验变异函数值如 图7（a）所示;② 目前普遍采用的方法,忽略时间变异影响直接对整个时段计算得到的空间实验变异函数值如图7（b）所示,该结果较单时段的结果有了较大改善,基本能够实现拟合;③ 利用本文提出的式（7）对多个子时段的实验变异函数值叠置后的结果如图7（c）所示;图7（d）为全时段和多时段的实验变异函数值拟合曲线对比。采用2类方法计算,经最小二乘法拟合得到的空间变异函数分别如下：
γ multi - period h s = 0.03 h s = 0 0.43 × 3 h t 2 × 1250 - h t 3 2 × 1250 3 + 0.03 0 < h s < 1250 0.46 h s > 1250 （11）
γ all - time h s = 0.13 h s = 0 0.30 × 3 h t 2 × 1350 - h t 3 2 × 1350 3 + 0.13 0 < h s < 1350 0.43 h s > 1350 （12）
Fig. 7 Fitting of space variogram

图7 空间变异函数拟合

（3）时空变异函数建模
拟合得到时间变异函数基台值 C t 0 = 0.81 ,2种方式拟合得到的空间变异函数基台值, C s 0 分别为0.46和0.43;时空基台值取时空样本变异函数最大值 C st 0 , 0 = 1 。由于2种方法得到的空间变异函数基台值相差不大,代入式（11）求解积合式模型的3个参数基本保持一致： k 1 = 0.71 , k 2 = 0.42 , k 3 = 0.67 。根据式（10）基于积合式模型构建的2类时空变异函数如下：
γ multi - p eriod h s , h t = γ multi - period h s + γ t h t - 0.67 × γ multi - period h s γ t h t （13）
γ all - time h s , h t = γ all - time h s + γ t h t - 0.67 × γ all - time h s γ t h t （14）
本文借助锚式浮标解决了Argo数据条件下区域内时间变异函数无法拟合的现状。在空间变异函数拟合部分,从图7（d）中可看出2种方法拟合的空间变异函数的变程和基台值相差不大,表明描述的海温区域化变量空间自相关范围大小和变化幅度大小是一致的。但采用全时段拟合的空间变异函数块金值相比于本文提出的多时段叠置拟合的结果更大,分析其主要原因在于全时段拟合的方法没有排除时间变异的干扰,误将时间变异成分理解成为了块金效应所产生的随机性,因此该方法会导致时空变异函数在近距离内产生较大误差,从而降低插值精度。而本文提出的多时段叠置拟合方法,在每个子时段内的时间变异程度极低,排除了时间变异干扰后采样点之间的空间变异情况,更纯粹、更符合实际情况;同时,各时段的实验变异函数值的叠置后保证了拟合过程所需样本点充足,因而拟合出的空间变异函数更加合理,进一步影响构建的时空变异函数。
4.3.2 实验2：优化时空Kriging插值与精度对比
本节采用上述模型对实验区域进行插值实验,实验过程采用“留一法”交叉验证：在2个实验区的时间范围内,逐日进行查询区域内的采样点,并每次预留一个观测点作为估值点进行插值,流程如图8所示。将插值结果于实际观测结果进行了比对得到实际误差。其中,实验区A的时间范围为2015年6月6日至2015年7月23日,区域内采样点数量为479个;试验区B的时间范围为2015年1月6日至2015年2月22日,区域内采样点数量为353个。本文截取部分日期的插值结果,如图9所示。
Fig. 8 Interpolation process of proposed Spatio-temporal Kriging

图8 本文Kriging插值流程图

Fig. 9 Error of interpolation in partial date

图9 部分日期插值结果误差

将本文方法分别与一般时空Kriging法以及时空权重插值进行横向对比,并对算法可靠性进行验证。同样采用趋势面拟合后的残差进行插值,其区别为：① 本文的时空Kriging优化方法采用的是实验1中多时段叠置拟合得到的时空变异函数参数;② 一般时空Kriging法采用实验1中全时段拟合得到的参数构建的时空变异函数进行插值;③ 时空权重插值法主要涉及：空间影响半径和事件影响范围2个参数,参考文献分别设置为500 km和7 d。
每日插值结果的平均误差对比如图10、11所示,以预测值与实测值之间的平均绝对误差（Mean Absolute Error,MAE）以及均方差（Mean Square Error,MSE）作为衡量精度依据（表1）。不难看出,在面向时空稀疏散布的Argo海温数据的插值实验中,Levy方法和一般时空Kriging方法的MAE分别在[0.5, 0.6]和[0.4, 0.5]且MSE较大,本文提出的优化算法MAE保持在0.2左右,精度和稳定性较前2种方法均有所提升。主要原因在于：① 3种方法均通过时空邻近点分配权重进行插值,时空权重法通过时空距离分配权重进行插值,仅考虑了估值点与邻近采样点之间的相互关系,存在精度不足的问题;② Kriging法建立在对采样点的时空自相关分析基础之上,依据自然现象的时空变异规律进行插值。由于考虑了预测点与采样点之间的相互关系,同时也考虑了每个采样点之间的时空关系,在插值精度上较时空权重插值具有一定优势;③ 本文的研究对象为时空稀疏散布数据,面临的主要困难在于单一时刻下空间邻域内的采样点数量不足,导致空间变异函数难以获得,需要补充时空邻域中的采样点拟合空间变异函数,一般时空Kriging法在未经优化的情况下,忽视了长时段内采样点之间的时间变异对空间变异函数值计算的干扰,导致最终拟合的空间变异函数表现出随机性过大。本文的时空Kriging优化方法通过子时段划分计算空间变异函数值,叠置拟合构建变异函数的方式有效削减了时间变异的影响,使得最终的插值精度进一步提高。
Fig. 10 Comparison of mean error in different day(area A)

图10 每日插值的平均误差对比（实验区A）

Fig. 11 Comparison of mean error in different day(area B)

图11 每日插值的平均误差对比（实验区B）

Tab. 1 Comparison of error statistics

表1 误差统计结果对比

误差统计量 本文方法 时空Kriging 时空权重法
实验区A 实验区B 实验区A 实验区B 实验区A 实验区B
MAE 0.160 0.204 0.540 0.603 0.407 0.473
MSE 0.067 0.102 0.176 0.301 0.131 0.168

5 结论

本文针对稀疏散布的时空数据集的结构特点,分析现阶段插值方法中存在问题,充分考虑时空变异特征,并基于时空Kriging法进行优化。首先对数据本身的趋势项进行剔除并验证其平稳假设条件,引入TAO辅助数据对时间变异函数进行解算,采用多时段叠置拟合空间变异函数的方法削弱了时间变异的干扰,基于积合式模型构建时空变异函数用于插值计算,通过交叉验证的方式对方法进行检验,并与其它插值方法进行了对比分析。结果表明,该方法构建的时空变异函数更适用于稀疏散布型的时空数据的Kriging插值,插值结果的误差降低,稳定性提高。自适应性不强需要人工干预是Kriging法的短板,后续研究中,将进一步挖掘时空位置与变异函数拟合模型选择的关联关系,以提高变异函数拟合部分的自动化程度,提升整体计算效率。

The authors have declared that no competing interests exist.

参考文献

原文顺序 | 文献年度倒序 | 文中引用次数倒序
[1]
李莎,舒红,徐正全.东北三省月降水量的时空Kriging插值研究[J].水文,2011,31(3):31-35.

[ Li S, Shu H, Xu Z Q.Study on spatial-temporal Kriging interpolation of monthly precipitation in three provinces of northeast China[J]. Journal of China Hydrology, 2011,31(3):31-35. ]

[2]
郑兴文,舒红,胡泓达.气温数据的时空变异结构分析[J].地理与地理信息科学,2015,31(6):56-61.时空地统计是传统地统计学结合时间维的扩展,其行业应用非常广阔,然而目前时空地统计模型的计算实现研究较少,时空变异函数是时空地统计的核心,研究时空变异函数理论模型的实现是时空地统计技术突破的关键。该文探讨了时空变异函数的积和模型理论,讨论了时空变异函数的两种实现方法和数据时空平稳化的预处理方法。以东北三省月平均气温数据为实验数据,检验了两种时空变异函数构造算法的效果,分析获知东北三省月平均气温时空结构性显著,空间结构性比时间结构性更为显著。基于现有实验数据,发现在空间48km内和时间3个月内东北月均值气温是相关的,超出这个范围则相关性逐步消失。该文工作提供了一条实验可行的时空地统计计算的技术路线。

DOI

[ Zheng X W, Shu H, Hu H D.Analysis of spatial-temporal variation of temperature data[J]. Geography and Geo-Information Science, 2015,31(6):56-61. ]

[3]
Cesare L D, Myers D E, Posa D.Estimating and modeling space-time correlation structures[J]. Statistics and Probability Letters, 2001,51(1):9-14.In this paper, a class of product–sum covariance models has been introduced for estimating and modeling space–time correlation structures. It is shown how the coefficients of this class of models are related to the global sill and “partial” spatial and temporal sills; moreover, some constraints on these sills have been given in order to assure positive definiteness of the product–sum covariance model. A brief comparative study with some other classes of spatial–temporal covariance models has been pointed out.

DOI

[4]
Ma C.Families of spatio-temporal stationary covariance models[J]. Journal of Statistical Planning and Inference, 2003,116(2):489-501.This paper provides simple methods for constructing new families of spatio-temporal stationary covariance models from purely spatial (or purely temporal) stationary covariance models. As the application of the methods developed, we introduce spatio-temporal stationary covariance models with the Gaussian and related spatial margins, and develop the Heine family and the Whittle–Matérn family of spatio-temporal stationary covariance models.

DOI

[5]
Mitchell M W, Genton M G, Gumpertz M L.Testing for separability of space-time covariances[J]. Environmetrics, 2005,16(8):819-831.Separable space-time covariance models are often used for modeling in environmental sciences because of their computational benefits. Unfortunately, there are few formal statistical tests for separability. We adapt a likelihood ratio test based on multivariate repeated measures to the spatio-temporal context. We apply this test to an environmental monitoring data set. Copyright 0008 2005 John Wiley &amp; Sons, Ltd.

DOI

[6]
Gneiting T.Nonseparable, stationary covariance functions for space-time data[J]. Journal of the American Statistical Association, 2002,97(458):590-600.Geostatistical approaches to spatiotemporal prediction in environmental science, climatology, meteorology, and related fields rely on appropriate covariance models. This article proposes general classes of nonseparable, stationary covariance functions for spatiotemporal random processes. The constructions are directly in the space090009time domain and do not depend on closed-form Fourier inversions. The model parameters can be associated with the data''s spatial and temporal structures, respectively; and a covariance model with a readily interpretable space090009time interaction parameter is fitted to wind data from Ireland.

DOI

[7]
Noel Cressie, Hsin-Cheng Huang.Classes of nonseparable, spatio-temporal stationary covariance functions[J]. Journal of the American Statistical Association, 1999,94(448):1330-1340.

DOI

[8]
李莎,舒红,徐正全.利用时空Kriging进行气温插值研究[J].武汉大学学报·信息科学版,2012,37(2):237-241.以黑龙江省近37a的月均气温数据为研究对象,介绍了一类积分式时空协方差（变异）函数模型进行时空Kriging插值。针对月均气温呈现出的明显季节变化,对各站点的气温进行去季节项处理,并在此基础上建立时空变异函数。将空间维的普通Kriging插值扩展至时空维,同时考虑空间和时间相关性对研究变量进行时空估计,并将估计结果与空间Kriging插值效果进行了比较。结果表明,时空插值效果理想,插值精度较空间Kriging更高。

[ Li S, Shu H, Xu Z Q.Interpolation of temperature based on spatial-temporal Kriging[J]. Geomatics and Information Science of Wuhan University, 2012,37(2):237-241. ]

[9]
胡丹桂,舒红,胡泓达.时空CoKriging的变异函数建模[J].华中师范大学学报(自科版),2015,49(4):596-602.在对地观测中,所研究的地学变量不仅具有时间、空间特征,还受其他变量的影响,采用多元时空相关数据,可以提高时空估值的精度.时空CoKriging是多元时空插值中一种常用的方法,建立时空变异函数和协变异函数是时空CoKriging插值的关键一步.以东北三省为试验区,利用该地区气象站观测数据的月平均空气相对湿度为主变量,引入同时间同位置的月平均空气温度作为协变量,对空气相对湿度和空气温度进行时空变异函数和时空协变异函数建模.实验结果表明,采用和度量时空模型的时空变异函数的随机性空间结构建模的实际拟合效果较好.

[ Hu D G, Shu H, Hu H D.Variogram modeling in space-time CoKriging[J]. Journal of Central China Normal University(Natural Sciences), 2015,49(4):596-602. ]

[10]
王建民,张锦,邓增兵,等.时空Kriging插值在边坡变形监测中的应用[J].煤炭学报,2014,39(5):874-879.lt;p>边坡形变监测中的监测点布设受限现场条件，布设位置可能不尽合理且数量有限，监测区域数据处理时需要进行时空插值。将经典的Kriging空间插值进行时空扩展，利用时空一类积和式方法建立了时空变异函数模型；分析了参数与变形的之间的关系；给出了时空插值的计算方法。以平朔露天矿边坡监测为实验对象，将整体时空Kriging插值与单一的Kriging空间插值结果进行对比分析，结果表明，同时兼顾了时间与空间关联性的时空Kriging插值精度有所提高。将监测区进行网格化，初步建立了监测区域的三维动态形变场模型，从整体上分析了边坡的稳定性。</p>

DOI

[ Wang J M, Zhang J, Deng Z B.Slope deformation analyses with space-time Kriging interpolation method[J]. Journal of China Coal Society, 2014,41(5):874-879. ]

[11]
陈鼎新,刘代志,曾小牛.时空Kriging算法在区域地磁场插值中的应用及改进[J].地球物理学报,2016,59(5):1743-1752.

[ Chen D X, Liu D Z, Zeng X N. Application and improvement of spatial temporal Kriging in geomagnetic field interpolation, 2016,59(5):1743-1752. ]

[12]
徐文,黄泽纯,张倩宁.基于时空模型的PM2.5预测与插值[J].江苏师范大学学报(自然科学版),2016,34(3):70-75.

[ Xu W, Huang Z C, Zhang Q N.Prediction and interpolation of PM2.5 based on space-time model[J]. Journal of Jiangsu Normal University(Natural Sciences Edition), 2016,34(3):70-75. ]

[13]
梅杨,杨勇,李浩.时空理论变异函数模型及其精度影响[J].测绘科学,2017,42(6):1-5.针对已有的时空模型研究深度不够的问题,该文主要研究两类时空理论变异函数模型,即分离和非分离模型,并分析其对时空预测精度的影响。首先,阐述了两类模型中较为流行的6种时空理论子模型及其参数物理意义;最后以2014年山东省PM2.5日均浓度为例,给出各模型理论表达式,并在Matlab中对各模型运用遗传算法进行参数拟合和时空预测精度验证。实验结果显示,6种模型的整体预测精度顺序与模型拟合精度顺序相当,说明时空理论模型同空间理论模型一样,其模型拟合精度与后续预测精度成正相关。

DOI

[ Mei Y, Yang Y, Li H.Study of spatio-temporal theory model and its influence on the spatio-temporal prediction accuracy[J]. Science of Surveying and Mapping, 2017,42(6):1-5. ]

[14]
杨胜龙,马军杰,伍玉梅.基于Kriging方法Argo数据重构太平洋温度场研究[J].海洋渔业,2008,30(1):13-18.

[ Yang S L, Ma J J, Wu Y M.Study on the reconstruction of Pacific temperature arena with Argo data based on Kriging methods[J]. Marine Fisheries, 2008,30(1):13-18. ]

[15]
Zeng L, Gad L.Space and time aliasing structure in monthly mean polar-orbiting satellite data[J]. Journal of Geophysical Research Atmospheres, 1995,100(D3):5133-5142.Monthly mean wind fields from the European Remote Sensing Satellite (ERS1) scatterometer are presented. A banded structure which resembles the satellite subtrack is clearly and consistently apparent in the isotachs as well as the u and v components of the routinely produced fields. The structure also appears in the means of data from other polar-orbiting satellites and instruments. An experiment is designed to trace the cause of the banded structure. The European Centre for Medium-Range Weather Forecasts gridded surface wind analyses are used as a control set. These analyses are also sampled with the ERS1 temporal-spatial sampling pattern to form a simulated scatterometer wind set. Both sets are used to create monthly averages. The banded structures appear in the monthly mean simulated data but do not appear in the control set. It is concluded that the source of the banded structure lies in the spatial and temporal sampling of the polar-orbiting satellite which results in undersampling. The problem involves multiple timescales and space scales, oversampling and undersampling in space, aliasing in the time and space domains, and preferentially sampled variability. It is shown that commonly used spatial smoothers (or filters), while producing visually pleasing results, also significantly bias the true mean. A three-dimensional spatial-temporal interpolator is designed and used to determine the mean field. It is found to produce satisfactory monthly means from both simulated and real ERS1 data. The implications to climate studies involving polar-orbiting satellite data are discussed.

DOI

[16]
Oliver M A, Webster R.Kriging: A method of interpolation for geographical information systems[J]. International Journal of Geographical Information Science, 1990,4(3):313-332.Geographical information systems could be improved by adding procedures for geostatistical spatial analysis to existing facilities. Most traditional methods of interpolation are based on mathematical as distinct from stochastic models of spatial variation. Spatially distributed data behave more like random variables, however, and regionalized variable theory provides a set of stochastic methods for analysing them. Kriging is the method of interpolation deriving from regionalized variable theory. It depends on expressing spatial variation of the property in terms of the variogram, and it minimizes the prediction errors which are themselves estimated. We describe the procedures and the way we link them using standard operating systems. We illustrate them using examples from case studies, one involving the mapping and control of soil salinity in the Jordan Valley of Israel, the other in semi-arid Botswana where the herbaceous cover was estimated and mapped from aerial photographic survey.

DOI

[17]
田垅,刘宗田.最小二乘法分段直线拟合[J].计算机科学,2012,39(b06):482-484.

[ Tian L, Liu Z T.Least-squares method piecewise linear fitting[J]. Computer Science, 2012,39(b06):482-484. ]

[18]
王珂靖,蔡红艳,杨小唤.多元统计回归及地理加权回归方法在多尺度人口空间化研究中的应用[J].地理科学进展,2016,35(12):1494-1505.对统计型人口数据进行格网形式的空间化可更直观地展示人口的空间分布,但不同的人口空间化建模方法和不同的格网尺度在表达人口空间化结果方面存在差异.本文在人口特征分区的基础上,引入DMSP/OLS夜间灯光对城镇用地进行再分类,采用多元统计回归和地理加权回归方法（GWR）,开展人口统计数据空间化多尺度模型研究,生成1km、5km和10km等3个尺度的2010年安徽省人口空间数据,并对3个尺度下2个模型结果进行精度评价与比较.结果表明：人口空间数据精度不仅与建模所用方法关系密切,还受到建模格网尺度大小的影响.基于多元统计回归方法的模型估计人口数与实际人口的平均相对误差值随着尺度的增加而降低,而基于GWR方法获得的人口空间数据误差值随着尺度的增加而升高.整体来看,基于GWR方法的1km研究尺度的人口空间数据平均相对误差最低（22.31％）.区域地形地貌条件与人口空间数据误差有较强的关联,地貌类型复杂的山区人口空间数据误差较大.

DOI

[ Wang K J, Cai H Y, Yang H Y, Yang X H.Multiple scale spatialization of demographic data with multi-factor linear regression and geographically weighted regression models[J]. Progress in Geography, 2016,35(12):1494-1505. ]

Options
文章导航

/

〈 〉