Estimation and implications of instrumental drift, random measurement error and nugget variance of soil attributes—a case study for soil pH

We propose a general model for soil pH measurement that includes instrumental drift, random measurement error, and random and correlated spatial variation. Methods for estimating these four components are described in detail. For soil pH in water, instrumental drift, random measurement error and random spatial variation (nugget effect) were greater than the corresponding quantities for soil pH in CaCl₂. For both pH measurements, instrumental drift was quite marked. Measurement error and nugget effect were of a similar size. A modified kriging method is presented that takes into account the four-component model proposed here. It is concluded that, for measuring soil chemical attributes, grid layouts should be supplemented by additional sites for the estimation of short-range variation, that laboratory sampling designs should include controls, and that field measurements should be adjusted for instrumental drift prior to being used for spatial contouring or kriging.     

污灌区土壤重金属空间结构与分布特征

该文探讨了地统计学插值模型应用于土壤重金属污染评价的适用条件,分析了北野场污灌区土壤重金属的正态分布特征和主导分布趋势,提出了不同重金属因子适宜的地统计插值模型。结果表明:土壤重金属空间变异系数处于12%～37%之间,整体变异性不大;Pb、Cd等因子的块金效应分别为0.90、0.87,空间相关性较弱,受人为随机因素的影响较大,As、Cr、Zn、Cu值等因子的块金效应分别为0.52、0.51、0.51和0.46,空间相关性中等,受人为随机因素和空间结构因素的共同作用,Hg和Ni的块金效应分别为0.253和0.06,空间相关性较强,受空间结构性因素的影响较大,可能受原生地质的影响较大。污灌区重金属含量存在增加趋势,灌区土壤重金属含量自北向南总体呈现较少趋势,这与灌区上下游灌溉保证率不同有关,表明长期污水灌溉可导致土壤重金属富集,但与全国其他污灌区相比,北野场污灌区土壤重金属污染相对较轻,应加强污水处理利用避免土壤重金属污染。     

OPTIMAL INTERPOLATION AND ISAR1THMIC MAPPING OF SOIL PROPERTIES

Kriging is a means of spatial prediction that can be used for soil properties. It is a form of weighted local averaging. It is optimal in the sense that it provides estimates of values at unrecorded places without bias and with minimum and known variance. Isarithmic maps made by kriging are alternatives to conventional soil maps where properties can be measured at close spacings. Kriging depends on first computing an accurate semi-variogram, which measures the nature of spatial dependence for the property. Estimates of semi-variance are then used to determine the weights applied to the data when computing the averages, and are presented in the kriging equations. The method is applied to three sets of data from detailed soil surveys in Central Wales and Norfolk. Sodium content at Plas Gogerddan was shown to vary isotropically with a linear semi-variogram. Simple punctual kriging produced a map with intricate isarithms and fairly large estimation variance, attributed to a large nugget effect. Sloniness on the same land varied anisotropically with a linear semi-variogram. and again the estimation error of punctual kriging was fairly large. At Hole Farm. Norfolk, the thickness of cover loam varied isotropically, but with a spherical semi-variogram. Its parameters were estimated and used to krige point values and produce a map showing substantial short-range variation.     

利用数字高程模型改进高山灰岩坑土壤pH值预测

Among spatial interpolation techniques,geostatistics is generally preferred because it takes into account the spatial correlation between neighbouring observations in order to predict attribute values at unsampled locations.A doline of approximately 15 000 m 2 at 1 900 m above sea level (North Italy) was selected as the study area to estimate a digital elevation model (DEM) using geostatistics,to provide a realistic distribution of the errors and to demonstrate whether using widely available secondary data provided more accurate estimates of soil pH than those obtained by univariate kriging.Elevation was measured at 467 randomly distributed points that were converted into a regular DEM using ordinary kriging.Further,110 pits were located using spatial simulated annealing (SSA) method.The interpolation techniques were multi-linear regression analysis (MLR),ordinary kriging (OK),regression kriging (RK),kriging with external drift (KED) and multi-collocated ordinary cokriging (CKmc).A cross-validation test was used to assess the prediction performances of the different algorithms and then evaluate which methods performed best.RK and KED yielded better results than the more complex CKmc and OK.The choice of the most appropriate interpolation method accounting for redundant auxiliary information was strongly conditioned by site specific situations.     

On spatial prediction of soil properties in the presence of a spatial trend: the empirical best linear unbiased predictor (E-BLUP) with REML

Geostatistical estimates of a soil property by kriging are equivalent to the best linear unbiased predictions (BLUPs). Universal kriging is BLUP with a fixed‐effect model that is some linear function of spatial co‐ordinates, or more generally a linear function of some other secondary predictor variable when it is called kriging with external drift. A problem in universal kriging is to find a spatial variance model for the random variation, since empirical variograms estimated from the data by method‐of‐moments will be affected by both the random variation and that variation represented by the fixed effects. The geostatistical model of spatial variation is a special case of the linear mixed model where our data are modelled as the additive combination of fixed effects (e.g. the unknown mean, coefficients of a trend model), random effects (the spatially dependent random variation in the geostatistical context) and independent random error (nugget variation in geostatistics). Statisticians use residual maximum likelihood (REML) to estimate variance parameters, i.e. to obtain the variogram in a geostatistical context. REML estimates are consistent (they converge in probability to the parameters that are estimated) with less bias than both maximum likelihood estimates and method‐of‐moment estimates obtained from residuals of a fitted trend. If the estimate of the random effects variance model is inserted into the BLUP we have the empirical BLUP or E‐BLUP. Despite representing the state of the art for prediction from a linear mixed model in statistics, the REML–E‐BLUP has not been widely used in soil science, and in most studies reported in the soils literature the variogram is estimated with methods that are seriously biased if the fixed‐effect structure is more complex than just an unknown constant mean (ordinary kriging). In this paper we describe the REML–E‐BLUP and illustrate the method with some data on soil water content that exhibit a pronounced spatial trend.     

Spatial prediction of soil properties using EBLUP with the Matérn covariance function

Spatial prediction with the presence of spatially dense ancillary variables has attracted research in pedometrics. While soil survey and analysis of soil properties are still expensive and time consuming, the secondary data can be made available on a dense grid for the whole area of interest. The main aim of using the ancillary data is to enhance prediction of soil properties by making use of the ancillary variables as covariates. Methods that can be used for this purpose are kriging with external drift, cokriging, regression kriging, and REML-EBLUP (Residual Maximum Likelihood-Empirical Best Linear Unbiased Predictor). Regression kriging is a sub-optimal method that has been utilised extensively because it is easy to use and has been shown empirically to perform as well as other methods. A statically sound method is REML-EBLUP. This paper examines the use of REML-EBLUP in combination with the Matérn covariance function for spatial prediction of soil properties. Methods for estimating parameters of the Matérn variogram using REML, and prediction with EBLUP are described. The prediction capability of REML-EBLUP, regression kriging, and ordinary kriging is compared for four datasets. Results show that although REML-EBLUP generally improves the prediction, the improvement is small compared with regression kriging. Thus, for practical applications regression kriging appears to be a robust method. REML-EBLUP is useful when the trend is strong, and the number of observations is small (< 200). We concluded that improvement in the prediction of soil properties does not rely on more sophisticated statistical methods, but rather on gathering more useful and higher quality data.     

Optimal interpolation and isarithmic mapping of soil properties: I The semi‐variogram and punctual kriging

Kriging is a means of spatial prediction that can be used for soil properties. It is a form of weighted local averaging. It is optimal in the sense that it provides estimates of values at unrecorded places without bias and with minimum and known variance. Isarithmic maps made by kriging are alternatives to conventional soil maps where properties can be measured at close spacings. Kriging depends on first computing an accurate semi‐variogram, which measures the nature of spatial dependence for the property. Estimates of semi‐variance are then used to determine the weights applied to the data when computing the averages, and are presented in the kriging equations. The method is applied to three sets of data from detailed soil surveys in Central Wales and Norfolk. Sodium content at Plas Gogerddan was shown to vary isotropically with a linear semi‐variogram. Ordinary punctual kriging produced a map with intricate isarithms and fairly large estimation variance, attributed to a large nugget effect. Stoniness on the same land varied anisotropically with a linear semi‐variogram, and again the estimation error of punctual kriging was fairly large. At Hole Farm, Norfolk, the thickness of cover loam varied isotropically, but with a spherical semi‐variogram. Its parameters were estimated and used to krige point values and produce a map showing substantial short‐range variation.     

Multivariate conditional stochastic simulation of soil heterotrophic respiration at plot scale

Soil heterotrophic respiration fluxes at plot scale exhibit substantial spatial and temporal variability. Within this study secondary information was used to spatially predict heterotrophic respiration. Chamber-based measurements of heterotrophic respiration fluxes were repeated for 15 measurement campaigns within a bare 13 × 14 m² soil plot. Soil water contents and temperatures were measured simultaneously with the same spatial and temporal resolution. Further, we used measurements of soil organic carbon content and apparent electrical conductivity as well as the prior measurement of the target variable. The previous variables were used as co-variates in a stepwise multiple linear regression analysis to spatially predict bare soil respiration. In particular the prior measurement of the target variable, the soil water content and the apparent electrical conductivity, showed a certain, even though limited, predictive power. In the first step we applied external drift kriging and regression kriging to determine the improvement of using co-variates in an estimation procedure in comparison to ordinary kriging. The improvement using co-variates ranged between 40 and 1% for a single measurement campaign. The difference in improving the prediction of respiration fluxes between external drift kriging and regression kriging was marginal. In a second step we applied sequential Gaussian simulations conditioned with external drift kriging to generate more realistic spatial patterns of heterotrophic respiration at plot scale. Compared to the estimation approaches the conditional stochastic simulations revealed a significantly improved reproduction of the probability density function and the semivariogram of the original point data.     

Regionalization of organic pollutants in Bavarian soils: The performance of Indicator Kriging

The evaluation of the impact of additional soil pollutants has to be contrasted against the naturally occurring pollutant concentration, i.e., the background concentration. Because background concentrations have to represent areal entities, point information has to be extrapolated into the area using interpolation methods. Thus, the accuracy of the interpolation method is crucial for the correct designation of background values to the areas. For the area of Bavaria (SE Germany), the actual background values of organic and inorganic soil pollutants were derived from >337,000 data from 5000 horizons based upon 1134 soil profiles. Background values were determined for the different soil depth compartments (O layer, topsoil, subsoil, and parent material) and land uses (agriculture, forestry). For interpolation between the nodes, Indicator Kriging was applied. The kriged total area was subdivided into 6 subareas of different background concentrations using percentile thresholds. To derive representative background concentrations, the reliable segregation of the total area into subareas and, thus, a robust interpolation method is a prerequisite. In this study, the robustness of the applied Indicator Kriging should be tested. Influences of data transformations and different kriging methods upon the demarcation of subareas should be investigated for the organic sum‐parameter EPA‐PAH. Neither a data transformation nor the comparison with Ordinary Kriging yielded significant deviations in the assigned subareas. Furthermore, cross‐validation as well as addition of synthetic noise was used to check the susceptibility of the method to artifacts and changes in the data set. After random splitting of the original data set into 4 subsets and re‐arrangement to 6 half‐sets, subsequent Indicator Kriging produced 6 results with mainly identical subarea configurations. Cross‐validation, i.e., comparison of points from the kriging surface (validation data set) with the calibration data set, yielded considerable residuals between estimates and measurements. Based on these normally distributed residuals, random numbers with identical statistical moments were generated and used as measurement errors for another kriging run. This synthetic noise was added to the corresponding result based on the calibration half‐set. The resulting subareas changed only slightly for the most polluted region, but considerably for the other regions. The chosen interpolation method provides sufficient stability to demarcate the relevant areas with elevated pollution in Bavaria. For other areas, its stability is less clear. Here, additional soil samples are required.     

Selection of a Geostatistical Method to Interpolate Soil Properties of the State Crop Testing Fields using Attributes of a Digital Terrain Model

The three most common techniques to interpolate soil properties at a field scale—ordinary kriging (OK), regression kriging with multiple linear regression drift model (RK + MLR), and regression kriging with principal component regression drift model (RK + PCR)—were examined. The results of the performed study were compiled into an algorithm of choosing the most appropriate soil mapping technique. Relief attributes were used as the auxiliary variables. When spatial dependence of a target variable was strong, the OK method showed more accurate interpolation results, and the inclusion of the auxiliary data resulted in an insignificant improvement in prediction accuracy. According to the algorithm, the RK + PCR method effectively eliminates multicollinearity of explanatory variables. However, if the number of predictors is less than ten, the probability of multicollinearity is reduced, and application of the PCR becomes irrational. In that case, the multiple linear regression should be used instead.     

OPTIMAL INTERPOLATION AND ISARITHMIC MAPPING OF SOIL PROPERTIES

Soil properties mapped in two intensive surveys had large nugget variances, leading to large estimation variances and erratic isarithms when mapped by punctual kriging. It is likely that both surveyors and survey clients are interested in average values of soil properties over areas rather than point values, and such values can be obtained by block kriging. Estimation variances are very much smaller, and maps of sodium and stone content at Plas Gogerddan, Central Wales, kriged over blocks 920m², and thickness of cover loam at Hole Farm, Norfolk, kriged over blocks of 400m², are much smoother than the punctually kriged maps. The map of Hole Farm has a distinct and meaningful regional pattern.     

Disaggregation of legacy soil data using area to point kriging for mapping soil organic carbon at the regional scale

Legacy data in the form of soil maps, which often have typical property measurements associated with each polygon, can be an important source of information for digital soil mapping (DSM). Methods of disaggregating such information and using it for quantitative estimation of soil properties by methods such as regression kriging (RK) are needed. Several disaggregation processes have been investigated; preferred methods include those which include consideration of scorpan factors and those which are mass preserving (pycnophylactic) making transitions between different scales of investigation more theoretically sound. Area to point kriging (AtoP kriging) is pycnophylactic and here we investigate its merits for disaggregating legacy data from soil polygon maps. Area to point regression kriging (AtoP RK) which incorporates ancillary data into the disaggregation process was also applied. The AtoP kriging and AtoP RK approaches do not involve collection of new soil measurements and are compared with disaggregation by simple rasterization. Of the disaggregation methods investigated, AtoP RK gave the most accurate predictions of soil organic carbon (SOC) concentrations (smaller mean absolute errors (MAEs) of cross-validation) for disaggregation of soil polygon data across the whole of Northern Ireland.Legacy soil polygon data disaggregated by AtoP kriging and simple rasterization were used in a RK framework for estimating soil organic carbon (SOC) concentrations across the whole of Northern Ireland, using soil sample data from the Tellus survey of Northern Ireland and with other covariates (altitude and airborne radiometric potassium). This allowed direct comparison with previous analysis of the Tellus survey data. Incorporating the legacy data, whether from simple rasterization of the polygons or AtoP kriging, substantially reduced the MAEs of RK compared with previous analyses of the Tellus data. However, using legacy data disaggregated by AtoP kriging in RK resulted in a greater reduction in MAEs. A jack-knife procedure was also performed to determine a suitable number of additional soil samples that would need to be collected for RK of SOC for the whole of Northern Ireland depending on the availability of ancillary data. We recommend i) if only legacy soil polygon map data are available, they should be disaggregated using AtoP kriging, ii) if ancillary data are also available legacy data should be disaggregated using AtoP RK and iii) if new soil measurements are available in addition to ancillary and legacy soil map data, the legacy soil map data should be first disaggregated using AtoP kriging and these data used along with ancillary data as the fixed effects for RK of the new soil measurements.     

Neural network ensemble residual kriging application for spatial variability of soil properties

High quality, agricultural nutrient distribution maps are necessary for precision management, but depend on initial soil sample analyses and interpolation techniques. To examine the methodologies for and explore the capability of interpolating soil properties based on neural network ensemble residual kriging, a silage field at Hayes, Northern Ireland, UK, was selected for this study with all samples being split into independent training and validation data sets. The training data set, comprised of five soil properties: soil pH, soil available P, soil available K, soil available Mg and soil available S,was modeled for spatial variability using 1) neural network ensemble residual kriging, 2) neural network ensemble and 3) kriging with their accuracies being estimated by means of the validation data sets. Ordinary kriging of the residuals provided accurate local estimates, while final estimates were produced as a sum of the artificial neural network (ANN) ensemble estimates and the ordinary kriging estimates of the residuals. Compared to kriging and neural network ensemble,the neural network ensemble residual kriging achieved better or similar accuracy for predicting and estimating contour maps. Thus, the results demonstrated that ANN ensemble residual kriging was an efficient alternative to the conventional geo-statistical models that were usually used for interpolation of a data set in the soil science area.     

A coherent geostatistical approach for combining choropleth map and field data in the spatial interpolation of soil properties

Information available for mapping continuous soil attributes often includes point field data and choropleth maps (e.g. soil or geological maps) that model the spatial distribution of soil attributes as the juxtaposition of polygons (areas) with constant values. This paper presents two approaches to incorporate both point and areal data in the spatial interpolation of continuous soil attributes. In the first instance, area-to-point kriging is used to map the variability within soil units while ensuring the coherence of the prediction so that the average of disaggregated estimates is equal to the original areal datum. The resulting estimates are then used as local means in residual kriging. The second approach proceeds in one step and capitalizes on: 1) a general formulation of kriging that allows the combination of both point and areal data through the use of area-to-area, area-to-point, and point-to-point covariances in the kriging system, 2) the availability of GIS to discretize polygons of irregular shape and size, and 3) knowledge of the point-support variogram model that can be inferred directly from point measurements, thereby eliminating the need for deconvolution procedures. The two approaches are illustrated using the geological map and heavy metal concentrations recorded in the topsoil of the Swiss Jura. Sensitivity analysis indicates that the new procedures improve prediction over ordinary kriging and traditional residual kriging based on the assumption that the local mean is constant within each mapping unit.     

县域土壤质量数字制图方法比较

土壤质量研究几乎涵盖土壤研究的所有领域,土壤质量制图理论与方法是土壤质量研究的一项重要研究内容。该研究以北京市密云县为研究区,基于土壤质量评价最小数据集和指数和法计算的土壤质量指数,探究了在地学模型支持下区域土壤质量数字制图方法。研究设计了5种区域土壤质量数字制图方法,并比较了不同方法的空间数字制图精度。结果显示,目前广泛使用的基于参评指标空间插值结果的土壤质量数字制图方法精度最低、工序较繁琐,且无法反映研究区景观高度异质的特点;而基于计算后的土壤质量指数（soil quality index,SQI）,借助于地统计学方法的土壤质量数字制图方法相对比较科学合理,其中又以基于计算后的SQI和回归克里格法预测效果最好,均方根误差最小,仅为0.01897,相对于基于参评指标空间插值结果的土壤质量数字制图方法,精度相对提高率最大,达到50%以上。综合考虑空间制图精度、工序的繁简程度,在该研究设计的5种方法中基于计算的SQI和回归克里格法最佳,该法避免了地统计插值在景观高度异质区的应用局限性,预测结果与实际最为相符。     

The use of region partitioning to improve the representation of geo statistically mapped soil attributes

An attempt to improve the representation of a geo statistically mapped soil attribute, clay content of the surface soil, through partitioning of the study area into two new regions was made. A topographic boundary divided the study area into hill and plain regions. Possible global non-stationarity or non-stationarity within the two newly defined regions was dealt with through the use of intrinsic random functions (IRF) of order k. Cross-validation of generalized covariance functions suggested that ordinary kriging might also have been appropriate. Exponential variogram models were subsequently fitted to the experimental variograms for each region. IRF-k block kriging and ordinary block kriging were then used as the primary methods of estimation. Both IRF-k and ordinary kriging performed badly in the vicinity of the topographic boundary when global models were used. This discontinuity was removed, at the expense of the introduction of some additional edge effects, when the hill and plain regions were kriged using models appropriate to each region. Independent zero-order generalized covariance functions with nugget and linear terms and exponential variogram models produced similar representations of clay content within each region, when used with their respective estimators. Splitting the region resulted in a 6% reduction in mean absolute deviation and a 14% reduction in mean squared deviation of predicted clay contents compared with a global model.     

Multiscale sources of spatial variation in soil. I. The application of fractal concepts to nested levels of soil variation

Soil variation has often been considered to be composed of‘functional’ or ‘systematic’ variation that can be explained, and random variation (‘noise’) that is unresolved. The distinction between systematic variation and noise is entirely scale dependent because increasing the scale of observation almost always reveals structure in the noise. The white noise concept of a normally distributed random function must be replaced to take into account the nested, autocorrelated and scale-dependent nature of unresolved variations. Fractals are a means of studying these phenomena. The Hausdorff-Besicovitch dimension D is introduced as a measure of the relative balance between long- and short-range sources of variation; D can be estimated from the slope of a double logarithmic plot of the semivariogram. The family of Brownian linear fractals is introduced as the model of ideal stochastic fractals. Data from published and unpublished soil studies are examined and compared with other environmental data and simulated fractional Brownian series. The soil data are fractals because increasing the scale of observation continues to reveal more and more detail. But soil does not vary exactly as a Brownian fractal because its variation is controlled by many independent processes that can cause abrupt transitions or local second order stationarity. Estimates of D values show that soil data usually have a much higher proportion of short-range variation than landform or ground water surfaces. The practical implication is that interpolation of soil property values based on observations from single 10 cm auger observations will be unsatisfactory and that some method of bulking or block kriging should be used whenever longrange variations need to be mapped.     

A practical two‐step approach for mixed model‐based kriging,with an application to the prediction of soil organic carbon concentration

Soil scientists often use prediction models to obtain values at unsampled locations. The spatial variation in the soil is best captured by using the empirical best linear unbiased predictor (EBLUP) based on a restricted maximum likelihood (REML) approach that efficiently exploits available data on both mean trends and correlation structures. We proposed a practical two‐step implementation of the REML approach for model‐based kriging, exemplified by predicting soil organic carbon (SOC) concentrations in mineral soils in Estonia from the large‐scale digital soil map information and a previously established prediction model. The prediction model was a linear mixed model with soil type, physical clay content (particle size < 0.01 mm) and A‐horizon thickness as fixed effects and site, transect, plot, year, year‐transect random intercepts and site‐specific random slopes for clay content. We used only the site‐specific intercept EBLUPs for estimating spatial correlation parameters as they described most of the variation in the random effects (86.8%). Fitting an exponential correlation model to these EBLUPs resulted in an estimated range of 10.5 km and the estimated proportion of the variance from the nugget effect was 0.23. The results of a simulation study showed a downwards bias that decreased with sample size. The results were validated through an external dataset, resulting in root mean square errors (RMSE) of 1.06 and 1.07% for the two‐step approach for kriging and the model with only fixed effects (no kriging), respectively. These results indicate that using the two‐step approach for kriging may improve prediction.     

Area-to-Point Kriging of Soil Phosphorus Composite Samples

Mapping of phosphorus (P) is based on sampling and laboratory analysis. Although laboratory analysis is costly, the number of samples is restricted in practice. In zone sampling, areas of the field are used to composite samples from sets of sampling points to reduce efforts. This study introduces area-to-point (ATP) kriging for downscaling composite samples with different sizes and shapes of the sampling areas. ATP kriging makes use of the coordinates of the sampling points of the composite samples. The applicability was tested on a simulated data set as well as on a spatially dense sample set of P measurements. Validation shows that ATP kriging outperforms point kriging with centroid interpolation. The root mean square error (RMSE) is reduced from 39.5 to 33.5 mg kg^?1. ATP kriging predictions were better at retaining the P value of the sampling area. The smoothing effect of interpolation and the aggregation effect of compositing the samples were reduced.