Optimal interpolation and isarithmic mapping of soil properties

The principle of optimal estimation using regionalized variable theory is extended from that of a single soil property to situations where there are two or more spatially interdependent ones. Auto and cross semi-variograms express the spatial relations among the properties concerned. They can be estimated from data and can then be used to interpolate the values of a variable by co-kriging from measurements of it plus data on one or more other properties that have been more intensively sampled. The technique of co-kriging is described and illustrated by a case study of the particle size distribution at Woburn experimental farm. There was a strong co-regionalization with common anisotropy between topsoil silt, subsoil silt and subsoil sand. This allowed topsoil silt to be estimated and mapped by co-kriging more precisely than by kriging from data on topsoil silt alone. When the auto and cross semi-variograms for a set of variables are known in advance or estimated from reconnaissance they can be used to plan an optimal sampling scheme. The main variable is sampled on a rectangular grid with finer grids for subsidiary variables. The maximum kriging variances are calculated for a range of sample spacings and relative sampling intensities. Those that match the maximum tolerable variance are potentially useful. The optimum scheme is the one that achieves the desired precision for least cost. For Woburn it is shown that measuring a main variable would need to cost at least 5 times that of a subsidiary variable to make a design for co-kriging economically sound. Such differences are unlikely for particle size fractions. Nevertheless there are many other instances in soil research where there are large differences in cost. If there is also a strong co-regionalization then savings should be possible by designing a sampling scheme that takes advantage of co-kriging.     

A comparison of parametric and non-parametric methods for modelling a coregionalization

Two methods for modelling a coregionalization were compared, the traditional parametric linear model of coregionalization (LMCR) and a non-parametric method based on a Fourier transform of the empirical (cross-) correlogram maps. The methods were compared in terms of how well they fit the experimental correlograms, and the prediction quality of their co-kriged estimates. Three datasets were compared, each representing different situations where we might use co-kriging.We found that both methods were somewhat restricted in how well they could represent the experimental correlograms because of the constraint that any coregionalization model must be positive-definite. There was little to distinguish between both methods in terms of how well the models fitted the raw correlogram data.The cokriged estimates from both methods were very similar in terms of their accuracy however the kriging variances from the LMCR were a better reflection of the prediction error. The non-parametric modelling is substantially faster than modelling the LMCR so if the only interest is in obtaining cokriged estimates then it should seriously be considered. In cases where the kriging variances are of interest then the LMCR should be used.     

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.     

A comparison of some robust estimators of the variogram for use in soil survey

The standard estimator of the variogram is sensitive to outlying data, a few of which can cause overestimation of the variogram. This will result in incorrect variances when estimating the value of a soil property by kriging or when designing a sampling grid to map the property to a required precision. Several robust estimators of the variogram, based on location and scale estimation, have been proposed as improvements. They seem to be suitable for analysis of soil data in circumstances where the standard estimator is likely to be affected by outliers. Robust estimators are based on assumptions about the distribution of the data which will not always hold and which need not be made in kriging or in estimating the variogram by the standard estimator. The estimators are reviewed. Simulation studies show that the robust estimators vary in their susceptibility to moderate skew in the underlying distribution, but that the effects of outliers are generally greater. The estimators are applied to some soil data, and the resulting variograms used for ordinary kriging at sites in a separate validation data set. In most cases the variograms derived from the standard estimator gave kriging variances which appeared to overestimate the mean squared error of prediction (MSEP). Kriging with variograms based on robust estimators sometimes gave kriging variances which underestimated the MSEP or did not differ significantly from it. Estimates of kriging variance and the MSEP derived from the validation data were generally close to estimates from cross‐validation on the prediction set used to derive the variograms. This indicates that variogram models derived from different estimators could be compared by cross‐validation.     

Sampling and bulking strategies for estimating soil properties in small regions

Estimates of mean values of soil properties within small rectangular blocks of land can be obtained by kriging provided the semi-variogram is known. This paper describes optimal rectangular grid sampling configurations whereby estimation variances can be minimized. For linear semi-variograms square blocks are best estimated by sampling at the nodes of a centrally placed grid with its interval equal to the block side divided by the square root of the sample size. For spherical semi-variograms the same configuration is almost optimal. The estimation variance of a bulked sample can be identical with that of a kriged estimate where the semi-variogram is linear and equal portions of soil are taken from each node on the optimally configured grid and provided the soil property is additive. For spherical semi-variograms the above is approximately true. Comparisons with estimates that take no account of known spatial dependence show that the true variances can be much less than those apparent using classical theory, and the necessary sampling effort much less. Within block-variances are often needed for planning, and an appendix gives two-dimensional auxiliary functions from which they can be calculated for linear and spherical semi-variograms.     

OPTIMAL INTERPOLATION AND ISARITHMIC MAPPING OF SOIL PROPERTIES III CHANGING DRIFT AND UNIVERSAL KRIGING

Universal kriging is a form of interpolation that takes account of local trends in data when minimizing the error associated with estimation. The presence of such trends, or drifts as they are known, is identified qualitatively, and their form found quantitatively by structural analysis, which simultaneously estimates semi-variances of the differences between the drift and actual data. The resulting semi-variograms are then used for the interpolation. The method was applied to measurements of electrical resistivity made in the soil at 1 m intervals at Bekesbourne, Kent. Analysis showed that the data could be adequately represented as a series of linear drifts over distances of 4 m to 8 m and with negligible nugget variance. Semi-variances of residuals from the drift were computed, and used to krige missing values and so complete an isarithmic map of the site. The method is by no means universally applicable in soil survey, mainly because of the large nugget variances usually encountered. These effectively prevent any distinction between constant and changing drift. They arise in part because measurements are made on small widely separated volumes of soil. Universal kriging is likely to be profitable only where measurements are made on contiguous volumes of soil or after substantial bulking.     

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.     

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.     

Fuzzy classification methods for determining land suitability from soil profile observations and topography

Because conventional Boolean retrieval of soil survey data and logical models for assessing land suitability treat both spatial units and attribute value ranges as exactly specifiable quantities, they ignore the continuous nature of soil and landscape variation and uncertainties in measurement which can result in the misclassification of sites that just fail to match strictly defined requirements. This paper uses fuzzy classification to determine land suitability from (i) multivariate point observations of soil attributes, (ii) topographically controlled site drainage conditions, and (iii) minimum contiguous areas, and compares the results obtained with conventional Boolean methods. The methods are illustrated using data from the Alberta Agricultural Department experimental farm at Lacombe in Alberta, Canada. Data on site elevation and soil chemical and physical properties measured at 154 soil profiles were interpolated by ordinary block kriging to 15 m × 15 m cells on a 50 × 50 grid. The soil property data for each cell were classified by Boolean and fuzzy methods. The digital elevation model created by interpolating the elevation data was used to determine the surface drainage network and map it in terms of the numbers of cells draining through each cell on the grid. This map was reclassified to yield Boolean and fuzzy maps of surface wetness which were then intersected with the soil profile classes. The resulting classification maps were examined for contiguity to locate areas where a block of minimum size (45m × 45m) could be located successfully. In this study Boolean methods reject larger numbers of cells than fuzzy classification, and select cells that are insufficiently contiguous to meet the aims of the land classification. Fuzzy methods produce contiguous areas and reject less information at all stages of the analyses than Boolean methods. They are much better than Boolean methods for classification of continuous variation, such as the results of the drainage analysis.     

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.     

Bayesian estimation of survival from mark-recapture data

An understanding of survival patterns is a fundamental component of animal population biology. Mark-recapture models are often used in the estimation of animal survival rates. Maximum likelihood estimation, via either analytic solution or numerical approximation, has typically been used for inference in these models throughout the literature. In this article, a Bayesian approach is outlined and an easily applicable implementation via Markov chain Monte Carlo is described. The method is illustrated using 13 years of mark-recapture data for fulmar petrels on an island in Orkney. Point estimates of survival are similar to the maximum likelihood estimates (MLEs), but the posterior variances are smaller than the corresponding asymptotic variances of the MLEs. The Bayesian approach yields point estimates of 0.9328 for the average annual survival probability and 14.37 years for the expected lifetime of the fulmar petrels. A simple modification that accounts for missing data is also described. The approach is easier to apply than augmentation methods in this case, and simulations indicate that the performance of the estimators is not significantly diminished by the missing data.     

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.     

Geostatistical modelling and mapping of nematode-based soil ecological quality indices in a polluted nature reserve

Nematodes are indicators of soil quality and soil health. Knowledge of the relationships between nematode-based soil quality indices and environmental properties is beneficial for assessing environmental threats on soil biota. This study evaluated the spatial distribution of nematode-based soil quality indices in a 23-ha heavy metal-polluted nature reserve using geostatistical methods. We expected that a selection of abiotic soil properties (pH and moisture, clay, organic matter, cadmium (Cd), and zinc (Zn) contents) could explain a significant portion of the spatial variation of the indices and that regression kriging could more accurately model their spatial distribution than ordinary kriging. A stratified simple random sampling scheme was used to select 80 locations where soil samples were taken to extract nematodes and derive the indices. The area had a distinct gradient in soil properties with Cd and Zn content ranging from 0.07 to 68.9 and 5.3 to 1329 mg kg^-1, respectively. Linear regression models were fitted to describe the relationships between the indices and soil properties. By also modelling the spatial correlation structure of regression residuals using spherical semivariograms, regression kriging was used to produce maps of the indices. The regression models explained between 21% and 44% of the total original variance in the indices. Soil pH was a significant explanatory variable in almost all cases, while heavy metal conent had a remarkably low effect. In some cases, the regression residuals had spatial structure. Independent validation indicated that in all cases, regression kriging performed slightly better because of having lower values of the root mean square prediction error and a mean prediction error closer to zero than ordinary kriging. This study showed the importance of soil properties in explaining the spatial distribution of biological soil quality indices in ecological risk assessment.     

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.     

基于不同地表曲面模型预测土壤有机碳含量

Local terrain attributes,which are derived directly from the digital elevation model,have been widely applied in digital soil mapping.This study aimed to evaluate the mapping accuracy of soil organic carbon (SOC) concentration in 2 zones of the Heihe River in China,by combining prediction methods with local terrain attributes derived from different polynomial models.The prediction accuracy was used as a benchmark for those who may be more concerned with how accurately the variability of soil properties is modeled in practice,rather than how morphometric variables and their geomorphologic interpretations are understood and calculated.In this study,2 neighborhood types (square and circular) and 6 representative algorithms (Evans-Young,Horn,Zevenbergen-Thorne,Shary,Shi,and Florinsky algorithms) were applied.In general,35 combinations of first-and second-order derivatives were produced as candidate predictors for soil mapping using two mapping methods (i.e.,kriging with an external drift and geographically weighted regression).The results showed that appropriate local terrain attribute algorithms could better capture the spatial variation of SOC concentration in a region where soil properties are strongly influenced by the topography.Among the different combinations of first-and second-order derivatives used,there was a best combination with a more accurate estimate.For different prediction methods,the relative improvement in the two zones varied between 0.30％ and 9.68％.The SOC maps resulting from the higher-order algorithms (Zevenbergen-Thorne and Florinsky) yielded less interpolation errors.Therefore,it was concluded that the performance of predictive methods,which incorporated auxiliary variables,could be improved by attempting different terrain analysis algorithms.     

Spatio‐temporal change of soil organic matter content of Jiangsu Province,China, based on digital soil maps

Estimation of spatio‐temporal change of soil is needed for various purposes. Commonly used methods for the estimation have some shortcomings. To estimate spatio‐temporal change of soil organic matter (SOM) in Jiangsu province, China, this study explored benefits of digital soil maps (DSM) by handling mapping uncertainty using stochastic simulation. First, SOM maps on different dates, the 1980s and 2006–2007, were constructed using robust geostatistical methods. Then, sequential Gaussian simulation (SGS) was used to generate 500 realizations of SOM in the area for the two dates. Finally, E‐type (i.e. conditional mean) temporal change of SOM and its associated uncertainty, probability and confidence interval were computed. Results showed that SOM increased in 70% of Jiangsu province and decreased in the remaining 30% during the past decades. As a whole, SOM increased by 0.22% on average. Spatial variance of SOM diminished, but the major spatial pattern was retained. The maps of probability and confidence intervals for SOM change gave more detailed information and credibility about this change. Comparatively, variance of spatio‐temporal change of SOM derived using SGS was much smaller than sum of separate kriging variances for the two dates, because of lower mapping variances derived using SGS. This suggests an advantage of the method based on digital soil maps with uncertainty dealt with using SGS for deriving spatio‐temporal change in soil.     

Estimating variograms of soil properties by the method-of-moments and maximum likelihood

Variograms of soil properties are usually obtained by estimating the variogram for distinct lag classes by the method‐of‐moments and fitting an appropriate model to the estimates. An alternative is to fit a model by maximum likelihood to data on the assumption that they are a realization of a multivariate Gaussian process. This paper compares the two using both simulation and real data. The method‐of‐moments and maximum likelihood were used to estimate the variograms of data simulated from stationary Gaussian processes. In one example, where the simulated field was sampled at different intensities, maximum likelihood estimation was consistently more efficient than the method‐of‐moments, but this result was not general and the relative performance of the methods depends on the form of the variogram. Where the nugget variance was relatively small and the correlation range of the data was large the method‐of‐moments was at an advantage and likewise in the presence of data from a contaminating distribution. When fields were simulated with positive skew this affected the results of both the method‐of‐moments and maximum likelihood. The two methods were used to estimate variograms from actual metal concentrations in topsoil in the Swiss Jura, and the variograms were used for kriging. Both estimators were susceptible to sampling problems which resulted in over‐ or underestimation of the variance of three of the metals by kriging. For four other metals the results for kriging using the variogram obtained by maximum likelihood were consistently closer to the theoretical expectation than the results for kriging with the variogram obtained by the method‐of‐moments, although the differences between the results using the two approaches were not significantly different from each other or from expectation. Soil scientists should use both procedures in their analysis and compare the results.     

A comparison of kriging, cubic splines and classification for predicting soil properties from sample information

The clay content of the topsoil in two regions of contrasting physiography was predicted from sample data using four different procedures. The predictors were the means of mapped classes, the usual kriging estimator, a cubic spline interpolator and a kriging estimator within classes using a pooled within-class variogram. The performances of the procedures were evaluated and compared. In the first region, Sandford St Martin on Jurassic sediments where there were some abrupt changes in soil, the classification predicted best within those classes bounded by sharp change. Elsewhere the usual kriging performed somewhat better, and kriging within classes was still more precise. In the second region, Yenne on the alluvial plain of the Rhone where the soil varied gradually, kriging performed better than classification, though a small improvement resulted from combining kriging with classification. Both prediction by class means and kriging attempt to minimize the estimation variance, and their mean prediction variances were close to the theoretical values overall. Spline interpolation is more empirical, and though it followed the abrupt changes better than kriging, it fluctuated excessively elsewhere, and its overall performance was poorer than that of kriging.