土壤属性空间预测精度与不确定性分析——有限最大似然法和高程辅助变量的应用

1、摘 要土壤属性空间变化的准确预测是土壤学及环境科学等相关学科的一个共同主题。在全球变化、资源与环境、生态多样性、食品安全以及人口与耕地等一系列问题的相关研究过程中，对土壤属性空间分布信息，不论在信息的数量方面还是在准确性方面，都提出了越来越高的要求。本文从准确性和不确定性两个方面，对有限最大似然法（REML）和高程辅助变量在土壤属性空间预测中的应用进行了探讨和研究，主要内容与结论如下：（1）以北京市大兴研究样区作为研究样区，以土壤有机质、土壤含水量、土壤速效钾和土壤有效锰四个土壤属性作为目标变量，通过选择不同的样本点数，对比分析REML法和传统的矩量法（MoM）在计算变异函数准确性方面的表现。

2、结果表明：在样本点数从150个逐渐减少到50个过程中，当样本点数小于70个后，所有目标变量的预测精度都开始明显的下降。在样本点数比较少的情况下，REML法估计的变异函数比MoM法估计的变异函数更准确一些。 （2）以土壤有机质和土壤有效锰作为目标变量，通过选择不同的空间尺度，分析空间尺度对土壤属性变异函数准确性的影响。结果表明：空间尺度对于认识土壤目标变量的空间变异特征具有明显的影响；REML法在提高空间预测精度方面效果是否明显与研究区的空间尺度也紧紧相关，本文中，REML法估计的变异函数在20km和30km两个尺度范围内比MoM法变异函数的预测精度高。（3）以平谷区为研究样区，以高程作为辅助变

3、量，以土壤有机质、土壤速效钾、土壤有效磷、土壤有效铁等为目标变量，通过利用普通克里格法（OK）、协克里格法（CK）和回归克里格法（RK）三种方法，对比分析高程数据是否可以用来提高上述土壤属性的空间预测精度。结果表明：高程作为辅助变量，可以用来提高土壤属性的空间预测精度，但这种方法并不适合所有土壤属性，本文中土壤有效铜、土壤有效铁和土壤有效锰三种微量元素的预测精度没能够被提高；利用高程对土壤属性进行空间预测时，在选用最适宜的方法之前，应该对土壤变量的空间结构、土壤属性的全局趋势、土壤属性与高程之间的线性相关关系、结构相关关系等进行仔细分析。（4）利用高程、坡度和地形湿度指数等地形指标组合成不同的

4、外部趋势模型，通过RK法和有限最大似然法为基础的线性无偏最优估计法（BLUP-REML）对平谷区土壤有机质的空间分布进行预测。结果表明：地形指标高程、坡度和地形湿度指数等可以作为辅助变量用来提高土壤有机质的空间预测精度，尤其是通过不同地形指标之间的组合，能进一步的提高预测精度；在样本点数较少的前提下，在外部趋势模型对土壤有机质空间变化解释程度偏低的情况下， BLUP-REML法在预测的稳定性方面比RK法表现得更好。（5）将高程辅助数据与序贯高斯协模拟方法相结合，对密云样区和平谷样区的土壤有机质的空间分布进行模拟，通过对模拟结果的准确性、局部不确定性和空间不确定性进行评价。结果表明，在条件高斯模

5、拟的计算环境下，可以利用高程作为辅助变量来提高土壤有机质空间预测精度和空间不确定性模拟的准确性，但没能提高局部不确定性模拟的准确性。关键词：土壤属性，空间预测精度，不确定性，地统计方法，高程，有限最大残差法AbstractSpatial prediction of soil properties has become a common topic in soil science research. There are increasing demands for the information of spatial distribution of soil characteristics, b

6、oth in the amount of information and in terms of accuracy, for the study of the issues related to changes in the global resources and the environment, biodiversity, food security and population and arable land. In this paper, we discussed the issues of enhancing the accuracy of spatial prediction of

7、 soil from the use of geostatistics methods and auxiliary variables. Main contents and conclusions are as follows:Soil organic matter, soil water content, soil available potassium and manganese have been examined here from Beijing Daxing district with different sample sizes to compare the performanc

8、e of the MoM and REML. The results shows the size of data set have a marked impact on the accuracy of variogram. When the sample is less than 70 points, the prediction accuracy of target variables are beginning to decline significantly. The REML variograms are more precise than the MoM variograms in

9、 the cases of small size of sample.To assess the impact of spatial scale on variogram, the performance of MoM variogram and REML variogram were examined with the data of soil organic matter and soil available manganese from different scale plots. The results showed that the spatial scale had a obvio

10、us effect on the accuracy of variogram computed by both MoM and REML. In this paper, with soil organic matter and soil manganese as the target variable, REML variograms are more accurate than MoM variograms with the 20km scale and 30km scale plots.This paper compares the accuracy of predictions of t

11、he soil organic matter, soil available potassium, phosphorus, iron and other soil variables in Pinggu District from ordinary kriging, cokriging and regression kriging. For the last two, intensive elevation data were used as auxiliary variable. The results showed that elevation can be used to improve

12、 the accuracy of spatial prediction of soil properties, but this method is not suitable for all soil variables. The results of factor analysis suggested that linear relationship, structure relationship between elevation and soil properties and the global trend of soil properties should be carefully

13、checked before the selection of the most appropriate method.The performance of the empirical best linear unbiased predictor with residual maximum likelihood (BLUP-REML) with that of regression kriging (RK) were compared for predicting soil organic matter (SOM) with the presence of different external

14、 drifts. Terrain attributes (elevation, slope and topographic wetness index) calculated from DEM was used as external drift variable. The results suggest that topographical data can further improve the accuracy of the spatial predictions of SOM by using RK and REML. Although the improvement with REM

15、L-EBLUP is not noticeable compared with RK, the prediction capability of REML-EBLUP is more stable than RK especially when a smaller proportion of the variation in target variable is accounted for by a trend model.The sequential Gaussian simulation (SGS) and sequential Gaussian co-simulation (SGCS)

16、algorithms were compared with respect to the accuracy of predictions as well as to the uncertainty inherent in the spatial prediction of soil organic matter. The SGS algorithm accounted for only the SOM data. The SGCS accounted for both SOM data and intensive elevation data. The root mean square err

17、ors revealed that the more accurate simulations were those accounting for intensive elevation information by the SGCS method for the two areas compared with SGS. As regards modeling local uncertainty, SGCS performed better at modeling prediction uncertainty than SGS. In addition, the results of asse

18、ssing the standard deviation confirmed that the exhaustive elevation data could be used to reduce the spatial uncertainty of SOM by SGCS compared with the SGS algorithm.Keywords: soil properties, accuracy, uncertainty, geostatistics, residual maximum likelihood, elevation 93目 录第一章 绪论11.1选题背景与研究意义11.

19、2 研究现状综述11.3主要研究内容61.4 技术路线7第二章 样本点数对变异函数准确性的影响92.1 引言92.2 研究目的92.3 理论与方法92.4 研究区概况与数据获取112.5 结果分析132.6本章小结29第三章 空间尺度对变异函数的影响313.1 前言313.2 研究目的313.3 研究区介绍与数据获取313.4 结果分析333.5 本章小结40第四章 利用高程辅助提高土壤属性空间预测精度424.1 引言424.2 研究目的424.3 理论与方法434.4 研究区概况与数据获取454.5 结果分析474.6 本章小结56第五章 利用不同外部趋势变量预测土壤有机质空间分布575.1

20、 引言575.2 研究目的575.3 理论与方法575.4 研究区概况与数据获取595.5 结果分析615.6 本章小结70第六章 高程辅助变量在不确定性分析中的应用726.1 引言726.2 研究目的726.3 理论与方法726.4 研究区概况与数据获取756.5结果分析766.6 本章小结82第七章 主要结论与研究展望837.1 主要结论837.2 研究展望84参考文献85致 谢90个人简历91中国农业大学博士学位论文 第一章 绪论第一章 绪论1.1选题背景与研究意义土壤作为地球表层系统自然地理环境的重要组成部分，在全球气候变化、资源与环境、生态多样性、食品安全以及人口与耕地等一系列问题中

21、起着十分重要的作用。例如，土壤产生与排放CO2、CH4和N2O等温室气体的过程，是陆地生态系统碳氮循环的一个重要过程，是土壤碳氮库的主要输出途径，对全球气候变化，尤其是全球变暖起着重要的影响作用。随着上述热点问题的逐步深入研究，对相关的土壤物理、化学和生物过程的理解和认识显得越来越重要，因此对土壤属性空间变化的信息，不论在信息量方面还是在信息的准确性方面，都提出了越来越高的要求。例如，在区域土壤水分和溶质运移动力学模拟研究中，需要区域内每个具体点位的土壤饱和与非饱和导水率和土壤水分特征曲线等土壤水力学参数作为模型运算的基础；在现代精准农业中，为了实现田间的点位管理，就必须对田间具体点位的土壤养

22、分信息有一个清楚的认识。然而，土壤作为母质、地貌、气候、生物和时间等成土因素相互作用的产物，土壤属性在空间和时间尺度上表现出复杂的差异性。由于土壤的空间变异特性，在目前相关土壤研究工作中经常面临的一个主要问题是缺少必要的、准确的和足够的土壤属性空间分布信息。田间调查采样及随后的室内分析是目前唯一能够获取土壤准确信息的方法。土壤的田间调查一方面需要大量的人力、物力和时间，另一方面获取的土壤信息却十分有限，不能满足各方面对土壤信息的需求。为了有效的解决土壤信息的供求矛盾，同时也为了客观真实的反映区域内总的土壤变化特征，就需要通过一定的方法利用田间调查的准确土壤信息对未采样点的土壤属性值进行预测。在

23、提高土壤属性空间预测精度方面，目前地统计学中备受关注的两条主要途径是：（1）利用有限最大似然法估计变异函数；（2）利用在空间上与土壤相关的变量来辅助预测。本文将从样本点数、空间尺度以及外部趋势度等方面对上述两条途径的应用进行探讨。研究结果将有有助于对土壤属性空间预测精度的预测方法体系有更全面的认识，对正确认识土壤空间变化提供更完善的理论基础。这对于当前我国开展的测土施肥和耕地质量动态监测工作具有重要的实践意义。1.2 研究现状综述为了获取更准确的土壤属性预测值，土壤学家在土壤空间预测方面一直努力寻求更好的预测方法，一些先进的预测方法不断地被引进到土壤学中。例如，Burgess and Webs

24、ter （1980）将普通克里格方法应用到土壤属性的空间预测；Lark and Webster（1999）将小波分析方法应用到土壤学中进行尺度相关的分析；Christakos（1998）首先用贝叶斯最大熵方法解决土壤属性时空变化问题。在这些方法中，地统计法方法得到了广泛地关注和研究。地统计学起源于20 世纪60 年代，早期主要应用于采矿业和石油勘探中。地统计学的研究对象是空间变量（也叫区域化变量）。空间变量的特点是变量在空间上不一定是完全随机的或完全独立的，即在空间上可能是相互联系的。地统计学将空间变量在空间某一位置的测量值看作是满足固有假设条件的随机函数在这个位置上的一次实现值。由于不同位置

25、的样本值之间可能是相互联系的，就需要分析样本的空间位置是否含有必要的信息，来揭示变量性质随距离的相关关系，即变量的空间变异结构(张仁铎, 2005)。有关地统计学的基本理论在此不再累述，具体内容参考相关文献（Goovaerts，1997；Stein，1999；王政权，1999；Webseter and Oliver，2001；Wachernagel，2003；张仁铎，2005；）。如今，越来越多涉及到空间分析的学科求助于地统计学的研究工具，地统计学已经被广泛用于土壤学、生态学、环境科学、地理学等诸多领域的研究中。自从Burgess and Webster（1980）将地统计学引入土壤学后，地统

26、计学已经在土壤学的各个领域得到广泛地研究与应用。例如，土壤养分方面（胡克林等，1999； 张世熔等，2003；许文强等，2006；王淑英等，2008）、土壤有机质方面（张世熔等，2002；黄元仿等，2004）、土壤质地方面（高峻等，2003）、土壤污染方面（杨劲松等，2006； 张红艳等，2006；钟晓兰等，2007）、土壤盐分方面（胡克林，2000；杨劲松等, 2006）、土壤水分方面（李笑吟等, 2006; 杨劲松等, 2006）等。上述文献的侧重点都是利用克里格方法对各土壤属性的空间分布进行插值预测，然后对土壤属性的预测结果进行评价。在应用克里格法对土壤属性进行预测的过程中，还需要对克里

27、格方法的适用范围、克里格法预测能力表现以及克里格法预测结果存在的误差等问题有所了解和认识。相关文献对这些问题进行了探讨和研究。下面从变异函数准确性、地统计混合模型（geostatistics hybrid method）结合辅助变量的应用、预测结果不确定性分析三个方面进行阐述。1.2.1关于变异函数准确性的研究地统计方法以区域化变量理论为基础，将土壤特性的空间变化作为一个随机函数的实现。用来表征随机变量空间结构的变异函数是地统计学的核心。传统的估计变异函数的方法是Matheron（1965）提出的矩量法（method of moment, MoM）。MoM法估计变异函数包括两个阶段：第一，利用

28、样本数据对不同分离距离计算样本变异函数；第二，对样本变异函数进行拟合得到一个连续的变异函数。为了准确地估计变异函数，用来计算样本变异函数值的数据量应足够大。Webster and Oliver（1992）根据Monte Carlo方法决定的置信限显示至少需要100150个样本数据，利用MoM法估计的变异函数才具有可靠性。但是，在实际应用中，经常出现样本数据比较少的情况。因此，需要一种方法能够利用比较少的样本数据估计出比较可靠的变异函数。Pardo-Iguzquiza（1998）提出用最大似然法（maximum likelihood, ML）估计变异函数是另外一种可以用于土壤学的方法。ML是一种

29、估计空间变量协方差参数的方法。Mardia（1984）首先在回归分析中用ML估计随机变量的空间协方差模型。随后，Kitanidis（1987）将ML法明确地用于地统计分析。Lark（2000）分别通过模拟数据和田间实测数据对MoM法和ML法进行了比较，模拟数据的结果显示ML法在采样点数为60的时候的结果与MoM法在采样点数为90到120点的结果一样，但是田间实测数据的结果却显示两种方法之间没有明显区别。相关文献(Cressie, 1993; Kitanidis and Lane, 1985) 指出由于ML法是同时估计趋势参数和协方差参数，导致了ML法估计的方差参数是存在偏差的。Patterso

30、n and Thompson（1971）提出的REML法克服了ML法的缺点。REML首先对原始数据进行线性合并为广义增量（generalized increments），将趋势去除，然后再通过最大似然函数估计协方差参数。REML在统计理论上是最合理的一种方法，最近已经在计量土壤学中引起了关注。Lark and Cullis（2004）详细的描述了如何利用REML法估计土壤属性的协方差参数，并建议在土壤属性空间预测过程中使用REML法。Kerry and Oliver（2007a）通过四个不同研究区域评价了样本点数量对REML法在估计变异函数方面的影响，结果显示，当样本点数小于100时，利用RE

31、ML法变异函数的预测结果比利用MoM法变异函数的预测结果要准确。Minasny and McBratney（2007）根据他们的研究结果指出，当样本点个数小于200时，就可以考虑使用REML法分析变异函数分析。1.2.2地统计混合模型方法结合辅助变量在土壤属性空间预测中的应用在地统计学中，土壤属性变量被作为一种区域化变量来分析，土壤属性的空间变异性可由两部分来代表，即空间确定性部分和随机部分，也称趋势成分和随机成分(张仁铎, 2005)。在普通克里金法中，假定确定性部分在空间上常量，主要是估计随机部分，在实际应用中，这些假定可能会限制克里金法的应用。例如，由于地形地貌的影响，土壤特征可能在研究

32、区内表现出明显的趋势，在这种情况下，确定性部分在空间上不再是常量，也就不能满足普通克里格法对数据的要求。泛克里格法（universal kriging：UK）是一种用来分析具有空间趋势变量的最佳线性估值方法。这种方法通常假定土壤变量的空间确定性部分（称为空间趋势）能够用空间坐标的多项式函数来表示。这样就可以首先通过空间坐标预测未采样点的土壤趋势部分，再利用普通克里格法预测未采样点上的随机部分，最后将空间趋势部分和随机部分进行叠加得到未采样点的土壤预测值。但是，用空间坐标来预测土壤变量的空间趋势也存在一定的局限性。外部趋势克里格法（kriging with external drift： KED

33、）和回归克里格（regression kriging：RK）作为泛克里格法的有效的推广形式，可以利用另外一种辅助变量的函数来表示土壤变量的空间趋势部分，弥补了UK法的缺陷。通常要求辅助变量具有以下特点，（1）获取比较容易或者测量费用相当廉价；（2）在研究区内信息量特别丰富；（3）与目标土壤变量存在一定的相关关系。RK法和KED法在数学推导公式上是相同的，两者的区别在于：在KED法中用来描述趋势的线性函数的系数在局部搜寻邻域内被看作常数，他们在克里格预测过程中被隐含地估计；而在RK法中，估计趋势线性函数的系数是独立与克里格法过程的。另外一种可以利用辅助变量的方法是多元克里格方法中的协克里格法（c

34、okriging：CK），CK法是普通克里格的拓展，考虑了空间两个或多个变量的空间交叉相关关系。有关UK法、KED法、RK法和CK法的基础理论在此不再详述，具体参考相关文献（Odeh et al., 1995；Hengl et al., 2004，2007；Goovaerts，1997；Webseter and Oliver，2001；Wachernagel，2003）。关于上述几种地统计混合模型方法的预测能力和适用范围，许多文献通过不同的土壤变量和辅助变量对这几种方法进行了比较。例如，Meul and Van Meirvenne（2003）通过不同类型的非平稳土壤质地数据和辅助变量高程，对O

35、K、CK、UK和RK法进行了比较，结果显示CK法和UK法可以很好的解释局部非平稳性。Simbahan et al（2006）利用表面反射率、土壤电导率和相对高程等多个外部辅助变量，通过CK、KED和RK三种方法，对土壤有机碳进行空间预测，结果表明，在这几种地统计混合方法，RK法是一种应用灵活且预测精度较高的方法。Baxter and Oliver（2005）以高程为辅助变量，通过CK、RK和KED三种方法，对两个研究区的土壤矿物氮和潜在有效氮进行分析，结果表明KED法和CK法都具有比较好的预测结果。从上述文献的结果可以看出，对于不同的土壤变量，不同的研究区域和不同的研究尺度，没有一种最好的预测

36、方法，每种方法都有一定特点和适用范围。下面简要论述了每种方法的应用条件及其不足之处。CK法在理论上是一种比较完善的多元克里格法。CK法要求的条件为：（1）目标变量与辅助变量都满足固有假定条件；（2）目标变量与辅助变量在空间上是相联系的，即在同一位置上，目标变量的取值与辅助变量的取值是相关的（Frogbrook and Oliver, 2001）；（3）目标变量的样本点比较少，而辅助变量的采样点则比较多。CK法在实际应用中存在两个方面的不足：（1）CK法要求目标变量的变异函数、辅助变量的变异函数以及二者之间的协变异函数是由共同的基本变异函数模型线性组合而成（即：协区域化线性模型），在协区域化线性

37、模型的拟合过程中受到的限制条件比较多（Wachernagel, 2003）；（2）CK法在预测过程中还存在计算量大、计算时间长的问题(Goovaerts, 1997)。因此，在应用中CK法并不是一种实用的方法。KED法的应用条件是：（1）辅助变量的空间趋势与目标变量的空间趋势是相关的；（2）目标变量去除趋势后的残差能够用地统计方法模拟；（3）目标变量与辅助变量是线性相关的，辅助变量必须在所有目标变量的采样点上是已知的，并且在预测格网的所有点上也是已知的。KED法的缺陷是：由于趋势系数的估计是在克里格法估值过程中同时进行的，当辅助变量与目标变量的关系不是简单的线性关系时，KED法的应用受到限制（

38、Deutsch and Journel, 1998）。RK法的应用条件与KED法相同，但是RK法明确地把趋势估计与残差插值分离，这样就可以利用任意复杂的回归形式，比较简单实用，因此RK法是一种应用最广泛的方法（Baxter and Oliver, 2005; Finke et al., 2004; Herbst et al., 2006; Simbahan et al., 2006），而且许多文献（Hong et al., 2005; Odeh et al., 1995; Simbahan et al., 2006）也报道了RK法的预测精度要高于OK法、CK法和RK法等方法的预测精度。但是，C

39、ressie（1993）和Lark et al.（2006）指出，RK法也存在两方面的缺陷：（1）趋势一般通过普通最小二乘法（OLS）来估计，虽然这个估计是无偏差的，但是这个估计没有产生最小方差估计；（2）利用回归残差估计的变异函数在理论上是有偏差的，而且偏差随着分离距离的加长而逐渐变大。RK法的缺陷导致两方面的结果，一方面是这种偏差将导致随机变量总的变化被低估，另一方面是拟合的变异函数不能正确的表达随机变量的空间相关性。鉴于RK法存在上述的缺陷，Lark et al.（2006）建议采用统计上更稳健的一种方法，即以有限最大似然法为基础的线性无偏最优估计方法（BLUP-REML），来预测具有空

40、间趋势的土壤变量。Lark (2006)指出BLUP-REML法能够无偏的估计残差的变异函数，能够得到最小方差的趋势估计，还能够得到未采样点上估计的趋势成分和随机成分的和，而且还知道预测方差。Minasny and McBratney（2007）利用了不同研究区的不同土壤属性数据对该方法进行了检验，结果显示当样本点个数小于200时，BLUP-REML法是一种有效的方法。然而，目前有关BLUP-REML法在土壤属性的空间预测方面的具体案例还比较少，对BLUP-REML法的预测能力及影响因素还需要进一步检验与评价。随着科学技术的发展，在空间上与土壤相关的辅助变量越来越容易获取，这些辅助变量主要包括

41、数值高程模型（DEM）、遥感影像数据和电磁波（EMI）扫描数据等。利用这些信息量丰富的数据来辅助提高土壤变量的空间预测精度，在计量土壤学土壤学中已经引起极大的关注和研究。例如，Odeh et al.（1995）利用高程数据通过RK法提高了上层土壤沙粒和下层土壤粘粒的空间预测精度。许红卫 等（2007）研究了土壤养分(土壤全N、全P、有机质、速效N、速效K) 与水稻各生长期水稻冠层光谱的关系, 并将光谱指数作为协因子, 进行土壤养分的Cokriging 插值研究，结果显示，利用光谱指数作为协因子可以提高土壤养分的空间预测精度。Huang et al.（2007）利用近红外波谱测量数据和LandS

42、at ETM影像数据，结合地形信息对土壤有机碳的空间预测进行了分析，结果证实用近红外波谱测量数据或LandSat ETM影像数据结合地形数据对提高土壤有机碳预测精度是一种有效的方法。Liu et al.（2008）利用电导率（ECa）、遥感高光谱反射率和DEM数据对土壤渗透力（soil drainage）的空间分布进行了预测，结果表明，单独利用地形指标作为辅助变量不能提高土壤渗透力的预测精度，但将地形指标与电导率或遥感数据结合后，就可以有效的提高土壤渗透力的预测精度。从上述可以看出，利用信息量丰富的辅助变量，可以有效地提高土壤属性的空间预测精度。在与土壤相关的辅助变量中，地形指标变量是一类主要

43、的辅助变量，相关文献也证实了地形指标与土壤属性之间明显的相关关系。More（1993）通过用地形指标作为解释变量，发现坡度的倾角、地形湿度指数可以解释土壤有机质、pH值、有效磷和粘粒含量等土壤变量空间变化的50%；Iqbal et al.（2005）通过在逐步线性回归分析，地形变量之间的组合能够介绍土壤空间变化的10%62%。可以看出，地形在土壤空间分布起着重要的影响作用。1.2.3土壤属性空间预测结果不确定性分析利用土壤采样点数据，采用克里格等方法对未采样点的土壤属性值进行预测，已经被深入研究和应用。然而，克里格估计值是有误差的，预测结果具有一定的平滑效应，即在土壤属性实际值比较低的点却有一

44、个高的估计值，或者在真实值高的地方又低估了。土壤管理者如果以克里格预测结果作为施肥决策依据，就可能导致错误的施肥决策。例如，在土壤养分含量缺乏的地点应该增加施肥却没有进行，在不需要施肥的地方却又增加施肥，导致资源浪费和其他环境问题。上述这种决策风险是不可避免的。如果这些风险能够被量化，管理者将能够对这些风险进行更有效地管理。这就需要我们对未采样点土壤属性值进行不确定性分析，估计未采样点土壤属性值超过或低于某一预设阈值的概率。目前，这种不确定性分析在决策过程和模型评价方面等方面受到越来越多的关注。地统计学中指示克里格法在解决这些问题方面比较突出，经常被用于土壤重金属污染、土壤养分分布等的相关研究

45、中。对于某一具体位置，指示克里格法用条件累计分布函数估计某一土壤变量取值大于给定阈值的概率。李笑吟等（2006）以普通克里格法和指示克立格法为基本方法,研究土壤水分空间分布及土壤水分有效性水平在空间上的概率分布。徐英等（2006）运用非参数地质统计学的多元指示克立格法, 结合单变量指示克立格法, 对黄河河套灌区的盐渍化土地上两个比较关键时期的土壤水分和盐分进行空间分布的分析, 给出了同一时期土壤水盐和不同时期土壤盐分满足一定标准的综合概率图。然而，指示克里格法以克里格估值作为基础，使这种方法不可避免的表现出克里格插值方法的缺点。克里格法估计预测结果的不确定性主要是通过计算预测值和相应的误差方差

46、，然后对二者进行计算构造高斯置信区间(Cinnirella et al., 2005)。这一过程存在三方面的不足，（1）克里格法的预测结果存在平滑效应，把目标变量空间变化的局部细节给平滑掉了，通常表现为把应该为高值的地方给估低了，把应该为低值的地方估高了；（2）克里格预测误差方差与实际预测值没有关系，它只依赖于采样数据的空间分布和变异函数模型；（3）计算置信区间时，需要假设预测误差为正态分布，然而土壤变量和其它环境变量经常呈偏态分布，这一假设通常很难满足(Bourennane et al., 2007)。由此可以看出，克里格法用估计值和估计误差构造置信区间并不是一个满意的解决方案。条件模拟方法

47、可以克服克里格法在估计不确定性方面的缺陷。条件模拟的主要思想是通过产生一组等概率的土壤属性空间分布图，然后根据所模拟的这些结果图间的差别来评价预测值的不确定性，而且这种方法还可以很好地再现原始变量的空间结构。条件模拟还可以同时考虑几个位置的空间不确定性。这一点也是指示克里格法的缺陷。所以条件随机模拟方法在土壤及其他学科中的应用受到越来越多的关注。例如，李保国（2002）等利用序贯高斯模拟对农田土壤表层饱和导水率进行了模拟，结果表明条件模拟与实测数据具有相同的统计分布特征及相同的空间结构。Zhao et al.（2005）利用条件随机指示模拟方法对河北省土壤有机质碳密度的空间预测结果的不确定性进

48、行了定量地评价。史舟（2007）以杭嘉湖平原为样区，利用普通克立格法和序贯高斯模拟方法对土壤Cu 的空间分布进行估值和模拟，并利用序贯指示模拟进行不确定性评价。结果表明由普通克里格法得到的土壤Cu的空间分布整体比较连续，具有明显的平滑效应，估值结果数据的分布频率趋于平缓；序贯高斯模拟结果整体分布相对离散，突出了原始数据分布的波动性，其模拟结果数据的分布频率相对集中。总的来看，在进行空间预测和不确定性分析时随机模拟方法比克里格方法显示出更多的优点，因此在土壤学科中得到越来越广泛地应用，尤其是对于研究那些具有空间不确定性且会对环境带来不良影响的风险性变量，更具有实际意义。1.3主要研究内容从上述研

49、究现状可以看出，目前提高土壤属性空间预测精度的基本途径主要有两方面：利用更完善的预测方法和利用信息量更丰富的辅助变量。在预测方法方面，将REML法估计变异函数与地统计方法结合引入到土壤属性空间预测领域时间并不长，相关研究还比较少。还有许多问题需要进一步的探讨与研究。同样，在利用辅助变量提高土壤空间预测方面也存在一些问题需要回答。为此，本文从预测方法和利用辅助变量两个方面对以下主要内容进行研究和阐述。（）评价样本点数对变异函数准确性的影响。变异函数是对目标土壤变量空间变异结构的一种表现形式。一方面不同的样本点数影响变异函数的估计结果，另一方面，目标变量本身的内在特性应该对变异函数具有更重要的影响

50、作用。本文将以同一研究区的土壤有机质、土壤含水量、土壤速效钾和土壤有效锰四个土壤属性作为研究目标变量，通过选择不同的样本点数，对比分析MoM法和REML法在计算变异函数准确性方面的表现。（2）评价空间尺度对计算变异函数的影响。土壤的空间变化存在于不同的空间尺度上，而不同空间尺度上具有不同的主要控制因素和过程。当我们用变异函数来描述土壤变量空间变异结构的时候，同样需要考虑土壤变化所处的空间尺度。本文以土壤有机质和土壤有效锰为目标变量，在保持样本点个数不变的情况下，通过选择不同的空间尺度样区，对比分析MoM法和REML法在计算变异函数准确性方面的表现。（3）利用高程作为辅助变量提高土壤属性空间预测

51、精度。本文将以高程数据作为辅助变量，以土壤有机质、土壤速效钾、土壤有效磷、土壤有效铁等为目标变量，通过利用OK法、CK法和RK法三种方法，对比分析高程数据是否可以用来提高上述土壤目标变量的空间预测精度。探讨土壤养分变量与高程数据之间应该满足哪些条件，就可以明确地利用高程数据来辅助提高土壤变量的空间预测精度。（4）评价RK法和BLUP-REML法在不同外部趋势度下的表现。地形在土壤的空间分布中起着重要的作用，本文将利用高程、坡度和地形湿度指数等地形指标作为辅助变量，通过组合成不同的外部趋势模型，利用RK法和BLUP-REML法来预测土壤有机质的空间分布。比较RK法和BLUP-REML法在不同程度

52、的外部趋势情况下预测能力的表现，评价高程、坡度和地形湿度指数等地形指标作为辅助变量在土壤有机质空间预测中的作用。（5）高程辅助变量在模拟不确定性方面的应用。本文将高程辅助数据与序贯高斯协模拟方法相结合，对土壤有机质的空间分布进行模拟，对预测结果的准确性、局部不确定性模拟的准确性和空间不确定性模拟的准确性进行评价。探讨在条件模拟计算环境下高程是否有助于提高土壤有机质的空间预测精度，是否有助于提高预测结果不确定性模拟的准确性。1.4 技术路线本文将按照图1-1所示的技术路线对上述主要内容进行研究：不同土壤属性变异函数准确性比较分析不同土壤属性同一辅助变量同一土壤属性不同辅助变量MoM法OK法REM

53、L法CK法RK法BLUP-REML法序贯高斯模拟序贯高斯协模拟不同样本点数不同空间尺度利用辅助变量比较分析局部不确定性对比分析辅助变量土壤属性空间预测精度分析预测结果不确定性分析参考方法空间不确定性对比分析土壤属性图1-1 本研究技术路线总体框图Fig. 1-1 Global framework for the research中国农业大学博士学位论文第二章 样本点数对变异函数准确性的影响第二章 样本点数对变异函数准确性的影响2.1 引言土壤属性的变异函数作为其空间变异结构的一种表现形式，一方面受土壤属性本身的内在特征决定，另一方面受样本点数和变异函数估计方法的影响。关于对样本点数最小需求这一

54、问题，不同的文献有不同的结论。Journel and Huijbregts (1980)指出一般必须大于30个样本点，而Webster and Oliver（1992）则建议至少需要100150个样本数据MoM法估计的变异函数才具有可靠性。REML法在理论上是一种更合理的计算变异函数的方法，Kerry and Oliver（2007a）的研究结果表明，当样本点数小于100时利用REML法变异函数的预测结果比利用MoM法变异函数的预测结果要准确，而Minasny and McBratney（2007）指出当样本点个数小于200时，就可以考虑使用REML法来代替MoM法分析变异函数分析。然而，对于

55、REML法来说，至少需要多少个样本点数其变异函数才具有可靠性，上述文献也没有明确结论。由于在土壤学中REML法还是一种比较新的方法，相关研究还比较少，仍然需要更多的土壤属性研究实例对REML法的适用范围和预测能力进行评价和研究。2.2 研究目的本章将以同一研究区的土壤有机质、土壤含水量、土壤速效钾和土壤有效锰四个土壤属性作为研究目标变量，通过选择不同的样本点数，对比分析MoM法和REML法在计算变异函数准确性方面的表现。本章的目的主要有以下三个方面：（1）进一步评价土壤样本点数对变异函数准确性的影响；（2）比较不同土壤变量随样本点数变化其变异函数准确性的变化方式；（3）比较MoM法和REML法

56、计算的变异函数的准确性随着样本点数变化的具体区别。2.3 理论与方法地统计分析的核心工具是变异函数。变异函数用来表征随机变量的空间变异结构。用来计算变异函数的方法有多种。目前最常用的估计变量函数的方法是MoM法。另外一种最新引进地统计的参数化估计方法是REML法。下面分别介绍两种方法的主要步骤及相关数学推导公式，有关这两种方法的详细理论参考相关文献(Pardo-Iguzquiza, 1997; Webster and Oliver, 2001)。2.3.1 MoM法估计变异函数变异函数与一定分离距离的两个随机变量取值的平方差相关。实际应用中，通过一定范围分离距离（例如1020km）内所有数据对

57、的平方差的平均值与各对样本数之间的距离的平均值来计算样本变异函数。MoM法首先根据样本数据计算样本变异函数值。设是土壤变量在个位置的样本数据值系列，在固有假设条件下，即： （2-1）样本的变异函数通过下式来计算： （2-2）式中代表样本变异函数，代表变异函数模型，是分离距离，是在之间用来计算样本变异函数值的样本的对数。根据目标变量的采样点数据利用公式（2-2）计算不同分离距离的样本变异函数值，然后利用不同分离距离的样本变异函数值，通过加权最小二乘法来拟合适当的变异函数模型。最常用的变异函数理论模型有球状模型（spherical）、指数模型（exponentail）和高斯模型（gaussian）等。有关这些模型详细介绍参考相关文献(Goovaerts, 1997; 张仁铎, 2005) 。变异函数模型的参数就可以直接参与到克里格计算中。在选择分离距离是有两种情况：（1）对于不规则的空间数据，需要以实际分离距离为中心定义一个分离距离等级，这个等级包括一定范围内的距离和方向。分离距离等级的选择对样本变异函数的影响至关重要，如果分离距离等级太短，变异函数比较嘈杂，如果分离距离太长，变异函数太光滑；（2）对于格网数据，格网间距通常被用来作为分离距离等级。2.3.2 REML法估计变异函数Patterson and Thompson (1971)提出的有限最大似然法