Bayesian approach for sample size determination, illustrated with Soil Health Card data of Andhra Pradesh (India).

A crucial decision in designing a spatial sample for soil survey is the number of sampling locations required to answer, with sufficient accuracy and precision, the questions posed by decision makers at different levels of geographic aggregation. In the Indian Soil Health Card (SHC) scheme, many thousands of locations are sampled per district. In this paper the SHC data are used to estimate the mean of a soil property within a defined study area, e.g., a district, or the areal fraction of the study area where some condition is satisfied, e.g., exceedence of a critical level. The central question is whether this large sample size is needed for this aim. The sample size required for a given maximum length of a confidence interval can be computed with formulas from classical sampling theory, using a prior estimate of the variance of the property of interest within the study area. Similarly, for the areal fraction a prior estimate of this fraction is required. In practice we are uncertain about these prior estimates, and our uncertainty is not accounted for in classical sample size determination (SSD). This deficiency can be overcome with a Bayesian approach, in which the prior estimate of the variance or areal fraction is replaced by a prior distribution. Once new data from the sample are available, this prior distribution is updated to a posterior distribution using Bayes' rule. The apparent problem with a Bayesian approach prior to a sampling campaign is that the data are not yet available. This dilemma can be solved by computing, for a given sample size, the predictive distribution of the data, given a prior distribution on the population and design parameter. Thus we do not have a single vector with data values, but a finite or infinite set of possible data vectors. As a consequence, we have as many posterior distribution functions as we have data vectors. This leads to a probability distribution of lengths or coverages of Bayesian credible intervals, from which various criteria for SSD can be derived. Besides the fully Bayesian approach, a mixed Bayesian-likelihood approach for SSD is available. This is of interest when, after the data have been collected, we prefer to estimate the mean from these data only, using the frequentist approach, ignoring the prior distribution. The fully Bayesian and mixed Bayesian-likelihood approach are illustrated for estimating the mean of log-transformed Zn and the areal fraction with Zn-deficiency, defined as Zn concentration <0.9 mg kg -1, in the thirteen districts of Andhra Pradesh state. The SHC data from 2015-2017 are used to derive prior distributions. For all districts the Bayesian and mixed Bayesian-likelihood sample sizes are much smaller than the current sample sizes. The hyperparameters of the prior distributions have a strong effect on the sample sizes. We discuss methods to deal with this. Even at the mandal (sub-district) level the sample size can almost always be reduced substantially. Clearly SHC over-sampled, and here we show how to reduce the effort while still providing information required for decision-making. R scripts for SSD are provided as supplementary material.

Propagation of uncertainties in soil and pesticide properties to pesticide leaching.

In the new Dutch decision tree for the evaluation of pesticide leaching to groundwater, spatially distributed soil data are used by the GeoPEARL model to calculate the 90th percentile of the spatial cumulative distribution function of the leaching concentration in the area of potential usage (SP90). Until now it was not known to what extent uncertainties in soil and pesticide properties propagate to spatially aggregated parameters like the SP90. A study was performed to quantify the uncertainties in soil and pesticide properties and to analyze their contribution to the uncertainty in SP90. First, uncertainties in the soil and pesticide properties were quantified. Next, a regular grid sample of points covering the whole of the agricultural area in the Netherlands was randomly selected. At the grid nodes, realizations from the probability distributions of the uncertain inputs were generated and used as input to a Monte Carlo uncertainty propagation analysis. The analysis showed that the uncertainty concerning the SP90 is 10 times smaller than the uncertainty about the leaching concentration at individual point locations. The parameters that contribute most to the uncertainty about the SP90 are, however, the same as the parameters that contribute most to uncertainty about the leaching concentration at individual point locations (e.g., the transformation half-life in soil and the coefficient of sorption on organic matter). Taking uncertainties in soil and pesticide properties into account further leads to a systematic increase of the predicted SP90. The important implication for pesticide regulation is that the leaching concentration is systematically underestimated when these uncertainties are ignored.

Impact of model uncertainty on soil quality standards for cadmium in rice paddy fields.

At present, soil quality standards used for agriculture do not consider the influence of pH and CEC on the uptake of pollutants by crops. A database with 750 selected paired samples of cadmium (Cd) in soil and paddy rice was used to calibrate soil to plant transfer models using the soil metal content, pH, and CEC or soil Cd and Zn extracted by 0.01 M CaCl2 as explanatory variables. The models were validated against a set of 2300 data points not used in the calibration. These models were then used inversely to derive soil quality standards for Japonica and Indica rice cultivars based on the food quality standards for rice. To account for model uncertainty, strict soil quality standards were derived considering a maximum probability that rice exceeds the food quality standard equal to 10 or 5%. Model derived soil standards based on Aqua Regia ranged from less than 0.3 mg kg⁻¹ for Indica at pH 4.5 to more than 6 mg kg⁻¹ for Japonica-type cultivars in clay soils at pH 7. Based on the CaCl2 extract, standards ranged from 0.03 mg kg⁻¹ Cd for Indica cultivars to 0.1 mg kg⁻¹ Cd for Japonica cultivars. For both Japonica and Indica-type cultivars, the soil quality standards must be reduced by a factor of 2 to 3 to obtain the strict standards. The strong impact of pH and CEC on soil quality standards implies that it is essential to correct for soil type when deriving national or local standards. Validation on the remaining 2300 samples indicated that both types of models were able to accurately predict (> 92%) whether rice grown on a specific soil will meet the food quality standard used in Taiwan.

Impacts of sampling design and estimation methods on nutrient leaching of intensively monitored forest plots in the Netherlands.

Element fluxes through forest ecosystems are generally based on measurements of concentrations in soil solution at regular time intervals at plot locations sampled in a regular grid. Here we present spatially averaged annual element leaching fluxes in three Dutch forest monitoring plots using a new sampling strategy in which both sampling locations and sampling times are selected by probability sampling. Locations were selected by stratified random sampling with compact geographical blocks of equal surface area as strata. In each sampling round, six composite soil solution samples were collected, consisting of five aliquots, one per stratum. The plot-mean concentration was estimated by linear regression, so that the bias due to one or more strata being not represented in the composite samples is eliminated. The sampling times were selected in such a way that the cumulative precipitation surplus of the time interval between two consecutive sampling times was constant, using an estimated precipitation surplus averaged over the past 30 years. The spatially averaged annual leaching flux was estimated by using the modeled daily water flux as an ancillary variable. An important advantage of the new method is that the uncertainty in the estimated annual leaching fluxes due to spatial and temporal variation and resulting sampling errors can be quantified. Results of this new method were compared with the reference approach in which daily leaching fluxes were calculated by multiplying daily interpolated element concentrations with daily water fluxes and then aggregated to a year. Results show that the annual fluxes calculated with the reference method for the period 2003-2005, including all plots, elements and depths, lies only in 53% of the cases within the range of the average +/-2 times the standard error of the new method. Despite the differences in results, both methods indicate comparable N retention and strong Al mobilization in all plots, with Al leaching being nearly equal to the leaching of SO(4) and NO(3) with fluxes expressed in mol(c) ha(-1) yr(-1). This illustrates that Al release, which is the clearest signal of soil acidification, is mainly due to the external input of SO(4) and NO(3).

Validation of regression models for nitrate concentrations in the upper groundwater in sandy soils.

For Dutch sandy regions, linear regression models have been developed that predict nitrate concentrations in the upper groundwater on the basis of residual nitrate contents in the soil in autumn. The objective of our study was to validate these regression models for one particular sandy region dominated by dairy farming. No data from this area were used for calibrating the regression models. The model was validated by additional probability sampling. This sample was used to estimate errors in 1) the predicted areal fractions where the EU standard of 50 mg l(-1) is exceeded for farms with low N surpluses (ALT) and farms with higher N surpluses (REF); 2) predicted cumulative frequency distributions of nitrate concentration for both groups of farms. Both the errors in the predicted areal fractions as well as the errors in the predicted cumulative frequency distributions indicate that the regression models are invalid for the sandy soils of this study area.

National baseline survey of soil quality in The Netherlands.

To determine background values of the 252 chemical compounds listed in Dutch legislation, a survey was designed with the aim of estimating percentiles of the cumulative frequency distributions and areal fractions exceeding the former, legislative reference values. Sampling locations were selected by probability sampling, so that in estimating the target quantities no model assumptions on the spatial variation were required, and valid estimates could be obtained by design-based inference. Strata in random sampling were formed by overlaying an aggregated soil map and land use map. For most of the heavy metals the areal fraction with concentrations in the topsoil (0-10 cm) exceeding the reference value was smaller than the allowable maximum of 5%. For these compounds a background value was determined smaller than the reference value. Exceptions were V, Co, Ba and Cu, for which a background value was defined (slightly) larger than the reference value.

Can we gain precision by sampling with probabilities proportional to size in surveying recent landscape changes in the Netherlands?

Seventy-two squares of 100 ha were selected by stratified random sampling with probabilities proportional to size (pps) to survey landscape changes in the period 1996-2003. The area of the plots times the urbanization pressure was used as a size measure. The central question of this study is whether the sampling with probabilities proportional to size leads to gain in precision compared to equal probability sampling. On average 1.03 isolated buildings per 100 ha have been built, while 0.90 buildings per 100 ha have been removed, leading to a net change of 0.13 building per 100 ha. The area with unspoiled natural relief has been reduced by 2.3 ha per 100 ha, and the length of linear relicts by 137 m per 100 ha. On average 74 m of linear green elements have been planted per 100 ha, while 106 m have been removed, leading to a net change of -31 m per 100 ha. For the state variables 'unspoiled natural relief', 'linear relicts', 'removed linear green elements', and 'new-removed linear green elements' there is a gain in precision due to the pps-sampling. For the remaining state variables there is no gain or even a loss of precision ('new buildings', 'removed buildings', 'new-removed buildings', 'new linear green elements'). Therefore, if many state variables must be monitored or when interest is not only in the change but also in the current totals, we recommend to keep things simple, and to select plots with equal probability.

Uncertainty and sensitivity analysis of spatial predictions of heavy metals in wheat.

Heavy metals seriously threaten the health of human beings when they enter the food chain. Therefore, policymakers require precise predictions of heavy metal concentrations in agricultural crops. In this paper we quantify the uncertainty of regression predictions of Cd and Pb in wheat (Triticum aestivum L.) and the contributions to the uncertainties in these predictions associated with inputs to the regression model. For each node of the 500- x 500-m grid covering the arable soils in The Netherlands, a latin hypercube sample size of 1000 is constructed from the uncertainty distributions of the explanatory variables (pH, soil organic matter [SOM], and heavy metal concentration in soil), the regression coefficients, and the random term of the regression model. This sample is used as input for the regression model to obtain 1000 values from the uncertainty distributions of the log(Cd) and log(Pb) concentration in wheat. There were no nodes where the recent EU quality standards for Cd and Pb (0.2 mg kg(-1) fresh wt.) in wheat were almost certain to be exceeded. For most nodes with clay soils, the quality standard for Cd in wheat almost certainly will not be exceeded; for Pb this is much less certain. The uncertainty in the Cd concentration in soil contributes most to the uncertainty in the predicted Cd concentrations in wheat (36% on the average), followed by the random term of the regression model (23%). For Pb the contribution of the random term is by far the largest (52%).

A method to combine non-probability sample data with probability sample data in estimating spatial means of environmental variables.

In estimating spatial means of environmental variables of a region from data collected by convenience or purposive sampling, validity of the results can be ensured by collecting additional data through probability sampling. The precision of the pi estimator that uses the probability sample can be increased by interpolating the values at the nonprobability sample points to the probability sample points, and using these interpolated values as an auxiliary variable in the difference or regression estimator. These estimators are (approximately) unbiased, even when the nonprobability sample is severely biased such as in preferential samples. The gain in precision compared to the pi estimator in combination with Simple Random Sampling is controlled by the correlation between the target variable and interpolated variable. This correlation is determined by the size (density) and spatial coverage of the nonprobability sample, and the spatial continuity of the target variable. In a case study the average ratio of the variances of the simple regression estimator and pi estimator was 0.68 for preferential samples of size 150 with moderate spatial clustering, and 0.80 for preferential samples of similar size with strong spatial clustering. In the latter case the simple regression estimator was substantially more precise than the simple difference estimator.

Mapping the probability of exceeding critical thresholds for cadmium concentrations in soils in The Netherlands.

The probability of exceeding critical thresholds of Cd concentrations in the soil was mapped at a national scale. The critical thresholds in soil were based on food quality criteria for Cd in crops or in organs of cattle (Bos taurus), and were calculated by inverting a regression model for the Cd concentration in the crop, with the Cd concentration in soil, soil organic matter (SOM) content, clay content, and pH as predictors. The probability of exceeding the critical threshold for Cd in soil per node of a 500- x 500-m grid was approximated by Monte Carlo simulation, using the estimated cumulative distribution functions (cdf) of SOM, clay, pH, and Cd as input. The cdfs were estimated by simple indicator kriging with local prior means. For SOM, clay, and pH, detailed maps of soil type and land use were used to define subregions with assumed constant local means of the indicators (a priori distributions). The cdfs were sampled by Latin hypercube sampling. We accounted for correlation between the actual and critical Cd concentrations in soil by drawing Cd values from cdfs conditional on SOM and clay. The estimated probability for grassland is negligible, even in areas with high Cd concentrations in soil, and for maize (Zea mays L.) land the probability is almost everywhere smaller than 5%. For arable soils, however, these probabilities commonly are larger than 5% when sugar beet (Beta vulgaris L.) or wheat (Triticum aestivum L.) is taken as a reference crop, and locally exceed 50%.