Modelling spatial uncertainty of soil erodibility factor using joint stochastic simulation

Soil erosion varies greatly over space and is commonly estimated using the revised universal soil loss equation (RUSLE). Neglecting information about estimation uncertainty, however, may lead to improper decision‐making. One geostatistical approach to spatial analysis is joint stochastic simulation, which draws alternative, equally probable, joint realizations of a regionalized variable. Differences between the realizations provide a measure of spatial uncertainty and allow us to carry out an error propagation analysis. The objective of this paper was to assess spatial uncertainty of a soil erodibility factor (K) model resulting from the uncertainties in the input parameters (texture and organic matter). The 500 km² study area was located in central‐eastern Sardinia (Italy) and 152 samples were collected. A Monte Carlo analysis was performed where spatial cross‐correlation information through joint turning bands simulation was incorporated. A linear coregionalization model was fitted to all direct and cross‐variograms of the input variables, which included three different structures: a nugget effect, a spherical structure with a shorter range (3500 m) and a spherical structure with a longer range (10 000 m). The K factor was then estimated for each set of the 500 joint realizations of the input variables, and the ensemble of the model outputs was used to infer the soil erodibility probability distribution function. This approach permitted delineation of the areas characterized by greater uncertainty, to improve supplementary sampling strategies and K value predictions. Copyright © 2007 John Wiley & Sons, Ltd.     

A stochastic-geometric model of soil variation

Stochastic models of soil variation are used in geostatistical analysis, but in general they bear no relation to our mechanistic understanding of the processes in soil that cause its properties to vary spatially. It is proposed that we require a suitable stochastic model in which space is partitioned into discrete domains as a first step towards random spatial models that incorporate our understanding of processes in soil. Even though the soil is essentially continuous in its spatial variation, there are components of soil variation (e.g. differences between parent materials) which are discontinuous. This paper shows how variogram models can be derived directly from the Poisson Voronoi Tessellation (PVT), a stochastic-geometric partition of d -dimensional space. The PVT variogram models, for d = 2 and 3, were fitted to variograms estimated from data over disparate scales, including computerized tomographic images of soil aggregates (pixels of a few tens of micrometres long) and the land systems of Swaziland. In all cases, PVT variogram models fitted better than the conventional geostatistical ones. The good performance of PVT variogram models at these disparate scales encourages further work on tessellation models for soil variation. In principle such models could incorporate information on underlying factors of soil formation such as the spatial distribution of individual plants, the origin and growth of microbial colonies, spatial processes in soil chemistry (such as reaction–diffusion processes) and geometrical information on boundaries between geological strata or contrasting plant communities. PVT models may therefore be one component of a random model of soil variation which reflects our understanding of soil-forming processes, and so have a stronger scientific basis than the models that are now in standard use.     

Estimation of Linear Models of Coregionalization by Residual Maximum Likelihood

Observations of ancillary soil properties spatially correlated to a soil property of interest may be used to increase the precision and reduce the sampling costs of a geostatistical survey. The relationship between such coregionalized properties must be expressed as a linear model of coregionalization but the conventional estimator of the linear model of coregionalization is biased unless the mean value of each property is constant across the study region. However, the mean value of a soil property may vary according to a spatial trend or a deterministic relationship with other factors which vary within the study region. We therefore propose that a linear mixed model should be fitted to coregionalized soil properties by residual maximum likelihood. This approach simultaneously fits spatial trends or deterministic relationships and random effects to the observations with minimum bias. We implement a residual maximum likelihood estimator for coregionalized properties and suggest a criterion to decide what order of spatial trend and which deterministic relationships should be included in the model. The effectiveness of the estimator is proved upon simulated data and upon observations of zinc and cadmium concentrations from the Swiss Jura.     

Scale-dependent correlation between topsoil copper and cobalt concentrations in Scotland

The soil of south-east Scotland is locally deficient in copper and cobalt. Measurements from nearly 3000 fields for which the soil association is known were analysed to study the coregionalization of the two elements and to assess the influence of parent material on the metals' concentrations. The experimental auto- and cross-variograms revealed distinct local (1.5 km) and regional (20 km) scales of spatial variation. A combination of indicator variograms of the soil associations had the same spatial structures, suggesting that parent material influences the concentrations of the metals. The coregionalization between copper and cobalt was modelled as a linear combination of three spatial structures. The resulting structural correlation coefficients showed the two elements to be fairly strongly positively correlated at the regional scale. Kriging allowed determination and mapping of each spatial component; these maps were then compared with the spatial distribution of soil associations in the region. An analysis of variance was performed before and after filtering out the nugget and short-range spatial components. Classification by soil association (parent material) accounted for a large proportion of the variance at the regional scale, suggesting that the parent material contributes substantially more to the trace element content of the soil than had been thought earlier.     

Factor Score Mapping of Soil Trace Element Contents for the Shenzhen Area

The abundances of copper, lead, zinc, nickel, chromium, cobalt, mercury, vanadium, and manganese in eighty three surface soil samples collected from Shenzhen area were determined. The correlation among element contents and the factors affecting the contents were studied using principle component analysis. The factor scores of the first two components were analyzed as regionalized variables using variogram analysis and Kriging. The interpolated factor scores were then mapped to show the common features in spatial distribution of a set of elements with similar geochemical behavior. It was demonstrated that the spatial distribution patterns of parent material in the area are the primary factor governing the spatial variation of the elements.     

Coregionalization of trace metals in the soil in the Swiss Jura

The provenance of trace metals in soil, whether from the parent material or from pollution, is rarely known with certainty, and the metals' history must usually be pieced together from fragmentary statistical information. This is particularly true in the Swiss Jura where the concentrations of several heavy metals around La Chaux de Fonds exceed the statutory recommended thresholds for safety. The topsoil of the 14.5-km² region was sampled on a square grid at 250-m intervals with additional nesting with distances of 100 m, 40 m, 16 m and 6 m. The concentrations of seven potentially toxic metals, namely Cd, Co, Cr, Cu, Ni, Pb and Zn, were measured. Their coregionalization could be represented by a linear model consisting of a nugget component plus two spherical structures with ranges of 0.2 km and 1.3 km. The short-range component dominated the variograms of Cd, Cr, Cu and Pb; the long-range component dominated those of Co and Ni; the variogram of Zn combined the two in approximately equal proportions. The coregionalization matrices contain moderate correlation among the nugget and the short-range components, notably between Cu and Pb, between Cd and Zn, and between Cr, Ni and Zn. The strongest correlations are at the long range between Co, Cr and Ni, and to a somewhat smaller degree between Zn and Co. Analysis of variance showed Co and Ni to be related to geology, and to the Argovian formation in particular. The indicator variogram of this formation has also a short-range component. The analysis also showed Cr and Cu to be related to land use (in different ways). Copper and Pb are strongly correlated to one another and distinct from the five other metals. The long-range structure is almost certainly a geological effect, whereas the one of short range probably results from both the geology and human activities.     

Multiscale spatial variability of CO2 emissions and correlations with physico-chemical soil properties

Spatial variability of greenhouse gas (GHG) emissions from agricultural lands is not well known although it has a great impact on the accuracy of GHG budget.The objectives of this study were to assess the spatial variability of CO₂ emission fluxes (CO_2-flux) and correlate these emissions with soil physico-chemical properties at two spatial scales and at different depths using a new geostatistical approach (coregionalization analysis with a drift, CRAD) that performs multiscale spatial analysis.Two agricultural sites with sandy and loamy soils were instrumented at 108 geo-referred sampling points and at two depths during spring 2007 where soil surface CO_2-flux and soil physico-chemical parameters were measured. The CO_2-flux presented spatial patterns characterized by different scales (i.e., non-spatial, small spatial and large spatial scale components), each describing a different fraction of its variability. About a quarter of CO_2-flux variability at the first site and one fifth at the other site was attributed to the non-spatial component. Strongest correlations were obtained between CO_2-flux and soil temperature, water saturation (S_w), elevation, electrical conductivity, soil bulk density, and the C/N ratio, but with differences between sites. Correlations were much stronger at large scale. Analyzing correlations between CO_2-flux and soil properties without discriminating for scales can miss important scale-dependent processes controlling soil gas emissions. Scales at which these processes vary should therefore be taken into account.     

Changes in variance and correlation of soil properties with scale and location: analysis using an adapted maximal overlap discrete wavelet transform

The magnitude of variation in soil properties can change from place to place, and this lack of stationarity can preclude conventional geostatistical and spectral analysis. In contrast, wavelets and their scaling functions, which take non‐zero values only over short intervals and are therefore local, enable us to handle such variation. Wavelets can be used to analyse scale‐dependence and spatial changes in the correlation of two variables where the linear model of coregionalization is inadmissible. We have adapted wavelet methods to analyse soil properties with non‐stationary variation and covariation in fairly small sets of data, such as we can expect in soil survey, and we have applied them to measurements of pH and the contents of clay and calcium carbonate on a 3‐km transect in Central England. Places on the transect where significant changes in the variance of the soil properties occur were identified. The scale‐dependence of the correlations of soil properties was investigated by calculating wavelet correlations for each spatial scale. We identified where the covariance of the properties appeared to change and then computed the wavelet correlations on each side of the change point and compared them. The correlation of topsoil and subsoil clay content was found to be uniform along the transect at one important scale, although there were significant changes in the variance. In contrast, carbonate content and pH of the topsoil were correlated only in parts of the transect.     

Assessment of Spatial Distribution of Selected Soil Properties using Geospatial Statistical Tools

To gain additional knowledge and better understand forest soil management on a small scale, geostatistical analytical tools were employed to examine the spatial distribution in dry aggregate mean weight diameter (MWD) and other selected soil properties and to assess the possible relationships between MWD and other soil properties. Selected properties of forest soils collected along a 300-m transact in the Nimbia Forest Reserve of Nigeria exhibited moderate to high variability in distribution with sodium ion displaying the greatest variability [coefficient of variation (CV, 91.2%)] and principal component analysis revealed the exchange complex cluster as influencing total variation of field soil properties. The autocorrelation function showed significant spatial correlation from 1 lag in soil organic carbon up to 17 lags (51 m) in soil moisture content (θ). The spherical and Gaussian semivariogram models described the spatial structure of most soil properties; however, for clay, cation exchange capacity (CEC), and soil organic carbon (SOC), an exponential model analyzed their spatial dependence.     

Robust estimation of the pseudo cross‐variogram for cokriging soil properties

The pseudo cross‐variogram can be used for cokriging two or more soil properties when few or none of the sampling locations have values recorded for all of them. The usual estimator of the pseudo cross‐variogram is susceptible to the effects of extreme data (outliers). This will lead to overestimation of the error variance of predictions obtained by cokriging. A solution to this problem is to use robust estimators of the pseudo cross‐variogram, and three such estimators are proposed in this paper. The robust estimators were demonstrated on simulated data in the presence of different numbers of outlying data drawn from different contaminating distributions. The robust estimators were less sensitive to the outliers than the non‐robust one, but they had larger variances. Outliers tend to obscure the spatial structure of the cross‐correlation of the simulated variables as described by the non‐robust estimator. The several estimators of the pseudo cross‐variogram were applied to a multitemporal data set on soil water content. Since these were obtained non‐destructively, direct measurements of temporal change can be made. A prediction subset of the data was subsampled as if obtained by destructive analysis and the remainder used for validation. Estimators of the auto‐variogram and pseudo cross‐variogram were applied to the prediction data, then used to predict the change in water content at the validation sites by cokriging. The estimation variances of these predictions were best calculated with a robustly estimated model of coregionalization, although the validation set was too small to conclude that the non‐robust estimators were unsuitable in this instance.     

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.     

Assessing heavy metal contamination in soils of the Zagreb region (Northwest Croatia) using multivariate geostatistics

The assessment of soil contamination and location of pollution sources represent a crucial issue in soil remediation. Topsoil samples were collected in the Zagreb area (Northwest Croatia) and the total contents of trace and major elements were determined. A multivariate geostatistical analysis was used to estimate soil chemical composition variability. Factorial Kriging Analysis (FKA) was used to investigate the scale-dependent correlation structure of some variables by modelling co-regionalization of ten chemical variables, co-kriging specific factors and mapping them. The FKA provided two regionalized factors at different spatial scales of variability: the first factor at shorter range for Zn, Pb, Cd, Cu and Ni indicated different sources of anthropogenic contamination, whereas Ca (mainly loading on the longer range factor) was related to the lithology and parent material composition. The methodology used has proved to be a useful tool to separate geological and anthropogenic causes of variation in soil heavy metal content and to identify common pollution sources.     

Two robust estimators of the cross‐variogram for multivariate geostatistical analysis of soil properties

If we wish to describe the coregionalization of two or more soil properties for estimation by cokriging then we must estimate and model their auto‐ and cross‐variogram(s). The conventional estimates of these variograms, obtained by the method‐of‐moments, are unduly affected by outlying data which inflate the variograms and so also the estimates of the error variance of cokriging predictions. Robust estimators are less affected. Robust estimators of the auto‐variogram and the pseudo cross‐variogram have previously been proposed and used successfully, but the multivariate problem of estimating the cross‐variogram robustly has not yet been tackled. Two robust estimators of the cross‐variogram are proposed. These use covariance estimators with good robustness properties. The robust estimators of the cross‐variogram proved more resistant to outliers than did the method‐of‐moments estimator when applied to simulated fields which were then contaminated. Organic carbon and water content of the soil was measured at 256 sites on a transect and the method‐of‐moments estimator, and the two robust estimators, were used to estimate the auto‐variograms and cross‐variogram from a prediction subset of 156 sites. The data on organic carbon included a few outliers. The method‐of‐moments estimator returned larger values of the auto‐ and cross‐variograms than did either robust estimator. The organic carbon content at the 100 validation sites on the transect was estimated by cokriging from the prediction data plus a set of variograms fitted to the method‐of‐moments estimates and two sets of variograms fitted to the robust estimates. The ratio of the actual squared prediction error to the cokriging estimate of the error variance was computed at each validation site. These results showed that cokriging using variograms obtained by the method‐of‐moments estimator overestimated the error variance of the predictions. By contrast, cokriging with the robustly estimated variograms gave reliable estimates of the error variance of the predictions.     

土壤空间变异研究中的半方差问题

简要回顾了土壤空间变异的研究。根据地质统计学理论和多年从事土壤空间变异研究的经验,对土壤空间变异研究的关键问题——半方差函数的基本假设、取样、模型选取及模型的检验进行了讨论,并对确定半方差函数模型应注意的问题提出建议。在保证取样样本容量的前提下,检查测定数据是否服从内蕴假设;注意提高每一个估算值的置信水平;尽量选择安全型模型作为半方差函数模型;对确定的半方差模型进行统计检验。由此可以求得较为客观合理的半方差模型。     

黄土高原六道沟小流域坡面土壤入渗特性的空间变异研究

水蚀风蚀交错带是黄土高原水土流失最为严重的地区,研究其不同地类土壤的入渗特性及其空间变异规律,有助于揭示黄土高原土壤侵蚀过程和生态环境的脆弱性以及流域水文模型精度的提高。本研究选取在水蚀风蚀交错带的强烈侵蚀中心六道沟小流域进行,采用双环定水头入渗法,在六道沟流域内一个长375m的完整天然坡面上网格布点,应用传统统计学方法和地统计学方法对稳定入渗率和前30min累积入渗量等土壤入渗特性重要参数的空间变异结构进行了研究,结果表明:(1)土壤稳定入渗率和前30min累积入渗量的变异系数分别为0.480和0.404,在坡面上的变异程度均呈现中等变异性;(2)两个入渗特性参数的试验半方差函数与理论半方差函数均拟合较好,自相关特征距离分别为126m和226m;块金值均大于0,表明在采样间隔范围内可能存在更小尺度的空间变异,要进一步研究在采样间隔内是否具有更小尺度的空间相关特征,可以增加采样密度来分析。     

典型喀斯特林地土壤养分空间变异的影响因素

为了探明喀斯特森林生态系统土壤养分空间异质性的成因及其对养分生物地球化学过程的指示意义,该研究以广西木论国家级自然保护区典型喀斯特峰丛洼地为研究对象,利用地统计学和经典统计方法分析了土壤养分的空间变异特征,并探讨了其主要影响因子。结果表明,研究区土壤有机碳(SOC)和全氮(TN)的块金值/基台值较大,分别为49.9%和28.6%,表现为中等程度的空间自相关,全磷(TP)和全钾(TK)的块金值/基台值较小,分别为10.4%和2.9%,表现为强烈的空间自相关,说明随机因素对TP和TK的影响相对较小;逐步回归分析表明,各环境因子对TK的方差解释最大,对SOC的方差解释最小。其中,土壤交换性Ca2+离子和凋落物中N含量是SOC和TN的主要控制因素,随着交换性Ca2+和凋落物中N含量升高,土壤SOC和TN积累增加;TP的控制因素比较单一,仅受凋落物中P含量影响。TK的影响因素比较复杂,除主要受交换性Ca2+控制外,凋落物N:P比、海拔高度和黏粒含量也有显著影响。     

应用土壤质地预测干旱区葡萄园土壤饱和导水率空间分布

田间表层土壤饱和导水率的空间变异性是影响灌溉水分入渗和土壤水分再分布的主要因素之一,研究土壤饱和导水率的空间变化规律,有助于定量估计土壤水分的空间分布和设计农田的精准灌溉管理制度。为了探究应用其他土壤性质如质地、容重、有机质预测土壤饱和导水率空间分布的可行性,试验在7.6 hm2的葡萄园内,采用均匀网格25 m×25 m与随机取样相结合的方式,测定了表层(0～10 cm)土壤饱和导水率、粘粒、粉粒、砂粒、容重和有机质含量,借助经典统计学和地统计学,分析了表层土壤饱和导水率的空间分布规律、与土壤属性的空间相关性,并对普通克里格法、回归法和回归克里格法预测土壤饱和导水率空间分布的结果进行了对比。结果表明:1)土壤饱和导水率具有较强的变异性,平均值为1.64 cm/d,变异系数为1.17;2)表层土壤饱和导水率60%的空间变化是由随机性或小于取样尺度的空间变异造成;3)土壤饱和导水率与粘粒、粉粒、砂粒和有机质含量具有一定空间相关性,而与土壤容重几乎没有空间相关性;4)在中值区以土壤属性辅助的回归克里格法对土壤饱和导水率的预测精度较好,在低值和高值区其与普通克里格法表现类似。研究结果将为更好地描述土壤饱和导水率空间变异结构及更准确地预测其空间分布提供参考。     

黄土高原小流域土壤有机碳空间变异性研究

Soil organic carbon (SOC) has great impacts on global warming, land degradation and food security. Classic statistical and geostatistical methods were used to characterize and compare the spatial heterogeneity of SOC and related factors, such as topography, soil type and land use, in the Liudaogou watershed on the Loess Plateau of North China. SOC concentrations followed a log-normal distribution with an arithmetic and geometric means of 23.4 and 21.3 g kg^-1, respectively, were moderately variable (CV = 75.9%), and demonstrated a moderate spatial dependence according to the nugget ratio (34.7%). The experimental variogram of SOC was best-fitted by a spherical model, after the spatial outliers had been detected and subsequently eliminated. Lower SOC concentrations were associated with higher elevations. Warp soils and farmland had the highest SOC concentrations, while aeolian sand soil and shrublands had the lowest SOC values. The geostatistical characteristics of SOC for the different soil and land use types were different. These patterns were closely related to the spatial structure of topography, and soil and land use types.