, 刘代志 2 1.清华大学 自动化系, 智能技术与系统国家重点实验室, 北京 100084;
2.火箭军工程大学 907教研室, 西安 710025
收稿日期: 2015-12-01
基金项目: 国家自然科学基金资助项目(41374154)
作者简介: 陈鼎新(1986-), 男, 博士研究生。
通讯作者: 陆文凯, 教授,E-mail: lwkmf@tsinghua.edu.cn
摘要：时空Kriging算法的核心, 是将变差函数的概念扩展到时空域, 变差函数的构建过程基础是计算时间切片和空间切片的向量距离。该文讨论了向量距离对构建时空变差函数的影响, 提出了空间距离加角度差异的向量距离模型。以地磁场观测数据作为对象, 分别用L1范数、L2范数和新距离模型对数据进行分析, 比较3种距离定义下的时空Kriging插值性能。结果表明: 加入了角度信息的向量距离, 能够更有效地表征数据, 提高时空Kriging的插值精度。
关键词： 地磁场 时空Kriging 变差函数 向量角度
Vector distance direction information for spatio-temporal Kriging
CHEN Dingxin1,2, LU Wenkai1
, LIU Daizhi21.State Key Laboratory of Intelligent Technology and Systems, Department of Automation, Tsinghua University, Beijing 100084, China;
2.Staff Room 907, PLA Rocket Force University of Engineering, Xi'an 710025, China
Abstract:The spatio-temporal Kriging method can be significantly improved by extending the variogram definition to the space-time domain. The key step in constructing the spatio-temporal variogram is to calculate the vector distances between the time slices and the space slices. This study analyzes the influence of the vector distance on the spatio-temporal variogram construction and presents a vector distance model that includes both the magnitude and the direction information. The algorithm was evaluated using magnetic field data with the evaluations based on the L1 norm and the L2 norm. The results show that the model with the additional direction information in the vector distance, more effectively represented the data characteristics which improved the spatio-temporal Kriging interpolation accuracy.
Key words: geomagnetic fieldspatio-temporal Krigingvariogramvector direction
时空Kriging是Kriging方法在时间—空间联合域的扩展，其利用数据量的时间相关性进行插值分析，对于研究连续变化的物理量具有优势。目前，时空Kriging在空气质量监测[1, 2]、 下雨量建模[3, 4]，风力数据插值[5]、 地下水位分析[6, 7]等领域已取得成功应用。Gething等[8]甚至用其分析门诊疟疾病人等连续变化的社会科学数据。
1 时空Kriging原理时空Kriging通过解Kriging方程组[1, 3]
| $\left\{ \begin{matrix} \sum\limits_{j=1}^{ N}{{{\lambda }_{j}}{{\gamma }_{s.t}}({{h}_{s0}},{{h}_{t0}}),} \\ \sum\limits_{i=1}^{N}{{{\lambda }_{i}}=1, } \\\end{matrix} \right.$ | (1) |
| $z_{st}^{*}(s,t)=\sum\limits_{l({{r}_{s}},{{r}_{t}})}{{{\lambda }_{i}}{{z}_{st}}({{s}_{i}},{{t}_{i}})}.$ | (2) |
| ${{\hat{\gamma }}_{s,t}}({{h}_{s}},{{h}_{t}})=\frac{1}{|L({{r}_{s}},{{r}_{t}})|}\times {{\sum\limits_{L({{r}_{s}},{{r}_{t}})}{\left[ Z(s+{{h}_{s}},t+{{h}_{t}})-Z(s,t) \right]}}^{2}}.$ | (3) |
| $\gamma (x,h)=\frac{1}{2}Var\left[ Z(x)-Z(x+h) \right]=\frac{1}{2}E{{\left[ Z(x)-Z(x+h) \right]}^{2}},$ | (4) |
| $\left\{ \begin{matrix} {{{\hat{\gamma }}}_{s,t}}({{h}_{s}})=\frac{1}{2N}{{\sum\limits_{i=1}^{N}{\left[ Z({{x}_{i}})-Z({{x}_{i}}+{{h}_{s}},) \right]}}^{2}} \\ {{{\hat{\gamma }}}_{s,t}}({{h}_{t}})=\frac{1}{2M}{{\sum\limits_{j=1}^{M}{\left[ Z({{t}_{j}})-Z({{t}_{j}}+{{h}_{t}}.) \right]}}^{2}} \\\end{matrix} \right.$ | (5) |
| $\left\{ \begin{matrix} {{\gamma }_{s}}({{{\dot{h}}}_{s}})=\frac{1}{2N}\sum\limits_{i=1}^{N}{d{{(Z({{x}_{i}}),Z({{x}_{i}}+{{h}_{s}}))}^{2}},} \\ {{\gamma }_{t}}({{{\dot{h}}}_{t}})=\frac{1}{2N}\sum\limits_{j=1}^{N}{d{{(Z({{x}_{j}}),Z({{x}_{j}}+{{h}_{t}}))}^{2}}.} \\\end{matrix} \right.$ | (6) |
2 向量距离对于向量X=(x1,x2,…,xn)，范数定义为[13]
| ${{\left\| X \right\|}_{p}}={{\left( {{({{x}_{1}})}^{p}}+{{\left( {{x}_{2}} \right)}^{p}}+\cdots +{{\left( {{x}_{n}} \right)}^{p}} \right)}^{\frac{1}{p}}}$ |
| ${{\left\| X \right\|}_{1}}=\left| {{x}_{1}} \right|+\left| {{x}_{2}} \right|\cdots \left| {{x}_{n}} \right|,$ | (7) |
| ${{\left\| X \right\|}_{2}}=\sqrt{{{\left( {{x}_{1}} \right)}^{2}}+{{\left( {{x}_{2}} \right)}^{2}}+\cdots {{\left( {{x}_{n}} \right)}^{2}}}.$ | (8) |
| $\left\{ \begin{matrix} {{\gamma }_{s}}\left( {{{\dot{h}}}_{s}} \right)=\frac{1}{2N}\sum\limits_{i=1}^{N}{\left\| Z\left( {{x}_{i}} \right)-Z\left( {{x}_{i}}+{{h}_{s}} \right) \right\|_{1}^{2},} \\ {{\gamma }_{t}}\left( {{{\dot{h}}}_{t}} \right)=\frac{1}{2M}\sum\limits_{j=1}^{M}{\left\| Z\left( {{t}_{j}} \right)-Z\left( {{t}_{j}}+{{h}_{t}} \right) \right\|_{1}^{2}.} \\\end{matrix} \right.$ | (9) |
| $\left\{ \begin{matrix} {{\gamma }_{s}}\left( {{{\dot{h}}}_{s}} \right)=\frac{1}{2N}\sum\limits_{i=1}^{N}{\left\| Z\left( {{x}_{i}} \right)-Z\left( {{x}_{i}}+{{h}_{s}} \right) \right\|_{2}^{2},} \\ {{\gamma }_{t}}\left( {{{\dot{h}}}_{t}} \right)=\frac{1}{2M}\sum\limits_{j=1}^{M}{\left\| Z\left( {{t}_{j}} \right)-Z\left( {{t}_{j}}+{{h}_{t}} \right) \right\|_{2}^{2}.} \\\end{matrix} \right.$ | (10) |
| $\begin{align} & D\left( {{Z}_{1}},{{Z}_{2}} \right)=p{{\left\| {{Z}_{1}}-{{Z}_{2}} \right\|}_{2}}+\left( 1-p \right)\left[ {{\cos }^{-1}}\left( \frac{Z_{1}^{T}{{Z}_{2}}}{{{\left\| {{Z}_{1}} \right\|}_{2}}{{\left\| {{Z}_{2}} \right\|}_{2}}} \right) \right], \\ & 0\le p\le 1, \\ \end{align}$ | (11) |
将式(11)代入式(6)，得到新向量距离定义下的条件变差函数
| $\left\{ \begin{matrix} {{\gamma }_{s}}\left( {{{\dot{h}}}_{s}} \right)=\frac{1}{2N}\sum\limits_{i=1}^{N}{D{{\left( Z\left( {{x}_{i}} \right),Z\left( {{x}_{i}}+{{h}_{s}} \right) \right)}^{2}}}, \\ {{\gamma }_{t}}\left( {{{\dot{h}}}_{t}} \right)=\frac{1}{2M}\sum\limits_{i=1}^{M}{D{{\left( Z\left( {{t}_{j}} \right),Z\left( {{t}_{j}}+{{h}_{t}} \right) \right)}^{2}}.} \\\end{matrix} \right.$ | (12) |
3 实验结果及分析地磁台站信号是连续的时间序列，可以将时空Kriging应用于地磁场的插值过程。将上述的向量距离应用到时空变差函数的计算中，并讨论结果。
3.1 交叉验证实验实验采用的数据，是中国地磁台网等32个观测台站从2002年1月1日到2002年9月30日的地磁场日均值序列。台站位置如图1所示，五角星分别表示32个观测台站的地理位置，覆盖了经度87.2°E—126.6°E，纬度19.0°N—49.6°N的范围。对32个观测台站进行交叉验证实验。
 |
| 图 1 台站分布图 |
| 图选项 |
实验中，计算某一时刻t的值，用其前20 d的观测量建立数据库，即用区间[t-20,t-1]内的所有观测值进行分析。性能的检验标准是误差统计量MAE (mean absolute error)[8]，
| $MAE = {1 \over n}\sum\limits_{i = 1}^n {\left| {{Z^ * }\left( {{s_i},t} \right) - Z\left( {s,t} \right)} \right|} .$ |
| $MSE = {1 \over n}\sum\limits_{i = 1}^n {{{\left( {{Z^ * }\left( {{s_i},t} \right) - Z\left( {{s_i},t} \right)} \right)}^2}} .$ |
表 1 条件变差函数的拟合参数
| 距离定义 | 域 | 函数 | 块金值 | 基台值 | 量程 | |
| L1 | 空间 | sphere | 0.048798 | 0.277273 | 37.41952 | |
| 时间 | Gauss | 0.00921 | 0.488384 | 0.148264 | ||
| L2 | 空间 | sphere | 0.083601 | 0.415782 | 35.5823 | |
| 时间 | Gauss | 0.01052 | 0.603239 | 0.152459 | ||
| 新定义 | 空间 | sphere | 0.03772 | 0.230226 | 29.11301 | |
| 时间 | Gauss | 0.00196 | 0.341643 | 0.115696 |
表选项
本文采用Product-Sum模型，利用条件分布函数来构建时空域变差函数[17-18]。3种定义下，求得时空变差函数，并分别代入式(1)中，可以解得Kriging方程组，进而算出时空中任意点的地磁场值。
从图3可以看出，相对于L1范数、 L2范数，新定义的距离使实验结果性能相对稳定，MAE集中在1左右，MSE集中在[0,5]的范围内。
将图3中各时刻MAE和MSE进行时间域的统计分析，结果如表2所示。L2范数下的误差均值比L1范数小，新定义的距离下误差均值比L2范数小，且MAE和MSE两项的方差都是最小。 3种定义下，L2范数优于L1范数。新定义距离的结果误差最小，性能最稳定，具有明显优势。角度信息的加入，使得时空变差函数在拟合过程中充分考虑了向量信息，从而提高了拟合精度，直至提升时空Kriging插值的总体精度。
表 2 交叉验证结果的时间统计量
| 距离定义 | MAE的均值 | MSE的均值 | MAE的方差 | MSE的方差 |
| L1 | 1.049 596 15 | 2.343 945 44 | 0.088 487 45 | 4.036 574 39 |
| L2 | 1.044 532 17 | 2.287 708 79 | 0.093 078 87 | 3.765 384 67 |
| 新定义 | 1.029 387 35 | 2.104 110 77 | 0.073 844 09 | 2.822 751 15 |
表选项
 |
| 图 2 不同距离定义下的条件变差函数拟合结果 |
| 图选项 |
 |
| 图 3 每个时间切片处的结果分析 |
| 图选项 |
4 结 论本文介绍了时空Kriging插值过程中变差函数的时空域扩展，并从向量距离的概念出发，讨论了不同距离定义对变差函数拟合的影响。将角度信息引入向量距离，定义了新的向量距离形式。对地磁场数据的交叉验证结果表明，引入向量角度的概念后，向量距离能够更好地反映各时间切片(或空间切片)间的差异，从而更好地拟合变差函数，提高插值的精度。
参考文献
| [1] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Iaco S D, Myers D E, Posa D. The linear coregionalization model and the product-sum space-time variogram[J].Mathematical Geology,2003, 35(35): 25–38. |
| [2] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Carla N, Amilcar S. Geostatistical space-time simulation model for air quality prediction[J].Environmetrics,2005, 16(4): 393–404. |
| [3] | Journal of Central South University(Science and Technology), 41(2):649-654.--> Li S, Shu H, Xu Z. Spatial-temporal statistics and analysis of rainfall in Jilin Province[C]//International Workshop on Computer Science for Environmental Engineering and Ecoinformatics, 2011, Kunming, China. Berlin:Springer-Verlag, 2011:255-261. |
| [4] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Tang Y. Comparison of semivariogram models for Kriging monthly rainfall in eastern China[J].Journal of Zhejiang University Science,2002, 3(5): 584–590. |
| [5] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Tilmann G. Nonseparable, stationary covariance functions for space-time data[J].Journal of the American Statistical Association,2002, 97(458): 590–600. |
| [6] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Cichota R, Hurtado A L B, Lier Q D J V. Spatio-temporal variability of soil water tension in a tropical soil in Brazil[J].Geoderma,2006, 133(s 3-4): 231–243. |
| [7] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Lark R M, Bellamy P H, Rawlins B G. Spatio-temporal variability of some metal concentrations in the soil of eastern England, and implications for soil monitoring[J].Geoderma,2006, 133(s 3-4): 363–379. |
| [8] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Gething P W, Akinson P M, Noor A M. A local space-time Kriging approach applied to a national outpatient malaria data set[J].Computers & Geosciences,2007, 33(10-50): 1337–1350. |
| [9] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Noel C. Spatial prediction and ordinary Kriging[J].Mathematical Geology,1988, 20(4): 405–421. |
| [10] | Journal of Central South University(Science and Technology), 41(2):649-654.-->De C L, Myers D, Posa D. Estimating and modeling space-time correlation structures[J].Stat Probab Lett,2001, 51(1): 9–14. |
| [11] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Iaco S D, Myers D E, Posa D. Nonseparable space-time covariance models:Some parametric families[J].Mathematical Geology,2002, 34(1): 23–42. |
| [12] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Siders I V, Gabella, Erdin R, et al. Real-time radar-rain-gauge merging using spatiio-temporal co-Kriging with external drift in the alpine terrain of Switzerland[J].Journal of the Royal Meteorological Society,2014, 140(680): 1097–1111. |
| [13] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Kalivas J H. Overview of two-norm (L2) and one-norm (L1) Tikhonov regularization variants for full wavelength or sparse spectral multivariate calibration models or maintenance[J].Journal of Chemometrics,2012, 26(6): 218–230. |
| [14] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Jin L H, Li D H, Song Y M. Combining vector ordering and spatial information for color image interpolation[J].Image and Vision Computing,2009, 27(4): 410–416. |
| [15] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Astola J, Haavisto J, Neuvo Y. Vector median filters[J].Proceedings of the IEEE,1990, 78(4): 678–689. |
| [16] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Trahanias P E, Karakos D G, Venetsanopoulos A N. Directional processing of color images:Theory and experimental results[J].IEEE Transactions on Image Processing,1996, 5(6): 868–880. |
| [17] | Journal of Central South University(Science and Technology), 41(2):649-654.-->De C L, Myers D, Posa D. Estimating and modeling space-time correlation structures[J].Stat Probab Lett,2001, 51(1): 9–14. |
| [18] | Journal of Central South University(Science and Technology), 41(2):649-654.-->Cesare L D, Myers D E, Posa D. Product-sum covariance for space-time modeling:An environmental application[J].Environmetrics,2001, 12(1): 11–23. |
| [19] | Journal of Central South University(Science and Technology), 41(2):649-654.-->李莎, 舒红, 董林. 基于时空变异函数的Kriging插值及实现[J].计算机工程与应用,2011, 47(23): 25–38. LI Sha, SHU Hong, DONG Lin. Research and realization of Kriging interpolation based on spatial-temporal variogram[J]. Computer Engineering and Applications, 2011,47(23):25-26. (in Chinese) |
