Printer Friendly

Towards cost-effective estimation of soil carbon stocks at the field scale.


The implementation of carbon-trading schemes or market-based instruments (MBI) for soil carbon requires accurate estimates of the status and change of soil carbon stocks. Such schemes often perform soil carbon budgeting at the farm scale. However, the field tends to be the natural management unit of a farm, and therefore, the effects of management practices on soil carbon stocks should be quantified at the field scale. Soil carbon stocks can vary greatly within a field according to locally variable factors such as plant growth, soil texture, drainage conditions, and landform. This means that a large number of measurements and, hence, substantial financial resources might be required to produce estimates of the field mean with acceptable precision. Therefore, the cost-effective sampling of carbon at the field scale is a pressing problem.

There are many ways in which the efficiency of soil carbon estimation at the field scale can be managed. Allen et al. (2010) suggest that classical, design-based sampling schemes and analyses are more efficient than model-based or geostatistical methods. Various types of design-based sample schemes are described by de Gruijter et al. (2006). Bricklemyer et al. (2005) note that by bulking soil cores within a small area, the short-scale variation of carbon that is not relevant to management practices can be removed. For this reason, the standard sampling unit for soil organic carbon in Australia is a 25 m by 25 m square (McKenzie et al. 2000; Bowman et al. 2009). Remote sensing and landscape digital data have been used to predict soil carbon levels at the field scale. Simbahan et al. (2006) used slope, surface reflectance, elevation, and electrical conductivity to stratify the field for sampling. This reduced sampling costs without a reduction in accuracy. Miklos et al. (2010) took a similar approach to estimate soil carbon levels at the farm scale, but they also utilised radiometric surveys to stratify the study area. Davis et al. (2004) noted that soil organic carbon was most variable in the top 20 cm of the soil profile and that it was correlated with soil texture, structure, gamma radiometric data, climate, and management history. Brus et al. (1999) considered how the costs of sampling could be managed. Such costs will vary according to the number of cores extracted, the distance travelled between sample sites, and the number of laboratory analyses. These costs can be reduced through clustering of sample sites and bulking.

The aims of this paper are to determine the most efficient strategy to estimate soil carbon stocks at the field scale and to quantify the associated costs. We consider these aims in the context of an Australian example of an MBI for soil carbon (Lachlan Catchment Management Authority 2011; Murphy et al. 2013). For this pilot scheme, it was decided that the field-mean soil carbon stocks needed to be estimated with a standard error <2 t/ha. Sanderman et al. (2010) suggested that this is the magnitude of the change in carbon stocks caused by various standard changes in management practices. The sample size required to achieve this level of accuracy will depend on the variability of carbon within the particular field being considered. Previous field-scale surveys of soil carbon have observed very different amounts of within-field variability. Knowles and Singh (2003) recorded coefficients of variation (CVs) of 30-40% for carbon percentage and 9-15% for bulk density in a Vertosol. Shukla et al. (2004) recorded CVs of 15-35% for carbon and <15% for bulk density. Much smaller carbon CVs, were recorded by Goidts et al. (2009) (8-15%), Conant and Paustian (2002) (13-24%), and Kravchenko et al. (2006) (11-28%). Hence, the amount of sampling required is likely to vary between different fields. Domburg et al. (1994) describe how the variogram must be known to estimate the expected errors from design-based surveys. Therefore, we conducted a detailed survey of soil carbon stocks within a 68-ha field in New South Wales and estimated its variogram. We used this to compare the relative merits of different sample designs to estimate the carbon stocks within the field. Some of these designs include local clusters of sample sites referred to as quadrats. We also consider whether digital elevation models (DEMS), terrain analyses, electrical conductivity surveys, and gamma radiometric surveys can be used to stratify the field and improve the sampling efficiency.

Statistical background

Design-based sampling and analysis

Classical or design-based statistical methods use probabilistic sampling of a spatial variable to estimate its mean over a study region and to determine the uncertainty of that estimate. Design-based methods are often favoured over geostatistical techniques in this context because they do not require that a model of variation is fitted to the observed data. Therefore, there is no need to confirm the validity of a model before inferences can be drawn, and fewer observations are generally required. The most basic design-based strategy is simple random sampling where every location within a study region is equally likely to be sampled. If the n measured values are denoted [z.sub.i] = 1,..., n, then the mean is estimated by:

[[??].sub.SI] = 1/n [n.summation over (i=1)] [z.sub.i] (1)

the sampling variance of this mean is estimated by:


and the standard error of the mean is estimated by calculating the square root of [??]. The efficiency of simple random sampling can be improved upon if the study area is divided into L strata and [n.sub.i] sampling locations are randomly selected from each stratum l = 1, ..., L. This approach assumes that it is possible to select strata within which the variation is less than the global variation. Strata can be selected according to the value of a covariate. Alternatively, Brus et al. (1999) advocated spatial stratification where the study area is divided into equal, compact subregions by a clustering algorithm (Walvoort et al. 2010). Here the assumption is that the variable is spatially correlated so that the variation within each subregion is smaller than the global variation. The global mean for stratified random sampling is estimated by the weighted average of the within-stratum means:

[[??].sub.ST] = [L.summation over (l=1)][w.sub.l][[??].sub.l] (3)

where [[??].sub.l] is the mean for stratum l estimated by simple random sampling and the corresponding weight [w.sub.t] is equal to the relative area of this stratum. The variance of this estimate of the mean is in turn estimated by:


In the spatial random-sampling strategy of Brus et al. (1999), the strata have equal area and the [w.sub.t] are all equal to 1/L. The cost effectiveness of design-based strategies can be improved by including local clusters of points. This approach assumes that that the costs associated with extracting clustered soil cores are relatively small, and therefore, it is an efficient way to account for short-range variation. In Australia, cluster designs based on a 25 m by 25 m square are specified (McKenzie et al. 2000; Bowman et al. 2009). If the variable being measured is additive, then the clustered soil cores can be bulked and the cluster mean measured directly. Then these cluster means can be treated as individual observations and the global mean and standard error estimated using the formulae that are appropriate for the design used to select cluster locations.

The expressions for the estimation variances of the global means (Eqns 2 and 4) can only be calculated once the variable has been measured. Therefore, it is not possible to use these expressions to decide how many observations are required in a particular design. However, the expected estimation variances for a particular variable and design-based strategy can be determined if the variable is assumed to be a realisation of a specified random function. If the random function is assumed to be second-order stationary, then it is described by its variogram.

[gamma](h) 1/2 E[[{z(x) - z(x + h)}.sup.2]] (5)

Here, z(x) is the value of the variable at location x and h is a lag vector separating two observations. The variogram describes how the correlation between observations of a variable decreases with separation. In general, the variogram can vary according to the length and direction of lag h, but often the variogram is assumed to vary only according to the length, which we denote h. In the next section we describe how the variogram can be estimated from observed data. Such data are unlikely to be available before conducting design-based sampling of a variable within a particular study region. However, one can approximate the variogram for such a variable based on previous surveys of the same property across a similar study region. In this paper we conduct a model-based survey of carbon stocks with the intention of estimating the variogram. Then we use the estimated variogram to predict the errors that can be expected from different design-based strategies.

Domburg et al. (1994) fully describe how variograms can be used to predict the errors that will occur from design-based strategies. For example, the sampling variance from simple random sampling can be predicted by:


where [bar.[gamma]] denotes the expected semi-variance between any two sites within the study area. This expected semi-variance can be estimated through simulation. Multiple realisations of pairs of observations are simulated at different locations within the study region and the mean semi-variance between these pairs is calculated. We describe how realisations of a variable can be simulated from their variogram in the next section. The same approach can be used if the observations are actually mean values for local clusters by simulating pairs of clusters. The corresponding formula for stratified random sampling is:


where [[bar.[gamma]].sub.i] denotes the mean semi-variance between all pairs of points in stratum 1. The variograms can differ between strata.

Estimation of variograms and simulation of random functions

We fitted nested nugget and exponential variogram models to the observed data by maximum likelihood (Webster and Oliver 2007). The nested nugget and exponential model for h > 0 is written as:

[gamma](h) = [c.sub.0] + [c.sub.1]{1-exp(-h/r)} (8)

If h = 0 then [gamma](h)= 0. The parameters of this model are the nugget variance [c.sub.0], the sill variance of the exponential component [c.sub.1], and a spatial parameter r. Since the model approaches [c.sub.0] + [c.sub.l] asymptotically, use of the exponential model implies that the variable is spatially correlated for all h and the range of correlation is infinite. The effective range of a model is defined as the lag distance at which 95% of the spatial correlation has been accounted for. For the exponential model, this is approximately 3r. The maximum likelihood estimator selects the model parameter values which maximise an analytical expression that the observed data would have arisen from the proposed model. This approach assumes that the variable is Gaussian. If this assumption is not plausible because the variable is highly skewed (skewness >1), then a log transform can be applied to the data before model estimation.

Once a variogram model has been estimated, it can be used to predict maps of the variable by kriging or to simulate realisations of the variable at specified sites. In this paper we simulate realisations of carbon stocks using the LU approach (Webster and Oliver 2007). This method applies a decomposition to the covariance matrix of the variable at the specified sites. The result is a large number of realisations of the variable where the correlation between each observation is consistent with the fitted variogram model.


The site

The study area was a cropping field in central-western New South Wales on the lower plains of the Macquarie River. The field is within the Carrabear-Macquarie meander plain and the main soils are Red Chromosols and Red Kandosols (Duncan et al. 2008) or Red Luvisols (WRB 2006), with small areas of more sandy soils with old levees. It is an inactive alluvial plain with alluvial sediments in the range 15 000-150 000 years old, but there is the possibility of some more recent alluvial activity in some of the lower slope areas. The land management history of the site is relatively simple. It was cleared in the early 1960s, when it was initially under a wheat-pasture rotation. In 1994 it was converted to a continuous cropping paddock with the rotation wheat-canola-wheat-legume crop-wheat. The legume crops included lupins and chickpeas. Initially from 1994 the cropping practice was minimum tillage, but from 2007 the cropping practice has been zero-till using a disc seeder and stubble retention, with controlled traffic using designated tramlines. Yields of wheat have averaged 2-2.5 t/ha between 2005 and 2008.

Soil sampling

A 100-m regular grid was laid over the whole field and aligned with the existing vehicle tramlines. Although the tramlines were used to set up the grid, the samples were taken well clear of the vehicle tracks, at least 5 m away. This provided 75 primary sampling sites within the field, to which were added 15 companion sites. The companion sites were randomly allocated to 15 of the primary sites; at eight of these points the sample was taken 10m from the primary point and the remaining seven were taken 1 m from the primary point. The purpose of the companion points was to enable some measurement of short-scale variability. Sampling points were located using a decimetre-accurate Trimble DGPS (Trimble Navigation Ltd, USA) with OmniSTAR HP correction.

A vehicle-mounted hydraulic corer was used to sample the soil to a depth of 0.5 m at each site. The soil cores were cut into depth increments 0-0.1, 0.1-0.2, 0.2 0.3, and 0.3 0.5m for laboratory analysis. At the time of sampling, the soil was in a moderately moist condition to ensure that coherent, intact cores were taken with little disturbance. The depth of the hole from which the core was taken was measured and compared with the length of the core to ensure that no compaction was occurring and that the core reflected the bulk density of the undisturbed soil. Samples were placed in labelled plastic bags for transport to the laboratory. The samples were weighed and a subsample was taken to measure the moisture content for estimating bulk density. The soil was air-dried at 40[degrees]C for 48 h and passed through a 2-mm mesh sieve. The soil was further ground to <53 [micro]m. Soil carbon percentages were determined with a vario MAX CNS analyser (Elementar, Germany). Bulk densities (t/[m.sup.3]) were estimated with the moisture content of the subsample, wet weight, and the volume of the core. Carbon storage (t/ha) over each depth increment was calculated from the percentage soil carbon and the bulk density. Since the soils were alluvial, no gravel was detected.

Sensor data

Soil apparent electrical conductivity measurements ([EC.sub.a], mS/m) were acquired using an EM38-DD electromagnetic induction instrument (Geonics Ltd, Canada). The depths of investigation of the vertical mode (ECv) and horizontal mode (ECh) were 1.5 m and 0.75 m, respectively.

A gamma radiometric survey had been conducted by the Geosciences Australia and State and Territory Geological Surveys. This provided maps of the ground level concentrations of potassium (K), uranium (U), and thorium (Th) based on airborne gamma-ray spectrometry (Minty et al. 2009). The gamma radiometric maps had 100-m pixels and this project used K (%K levelled and merged), equivalent-U (ppm U), equivalent-Th (ppm Th), and the total dose (Tdose).

A multispectral imagery SPOT5 satellite image was acquired on 29 February 2008 with a pixel size of 10m from the Department of Environment, Climate Change and Water, NSW. The four bands of SPOT5 are BI (green), B2 (red), B3 (near infrared), and B4 (short-wave infrared). Green normalised difference vegetation index (NDVI) and photosynthetic vigour ratio (PVR) were calculated from the SPOT imagery. The NDVl was defined as (B3 B2)/(B3 + B2) and PVR as B1/B2.

An ASTER digital elevation model (m) was obtained from NASA's Land Processes Distributed Active Archive Centre (LP DAAC). The DEM used for this research has a pixel size of 30m. The terrain attributes were processed in SAGA GIS software (Bock et al. 2007). The potential usefulness of the terrain attributes is restricted by the low relief of the paddock, which is on an alluvial plain.

Field and landform classification

The field was subjectively classified into different landforms (Fig. 1). Use was made of the published information in Duncan et al. (2008) and the local knowledge of authors and advisors of the influence of landforms on soil types. Obvious features such as the ridge line, depression, and lower slope area were recognised. Although the approach is subjective, many features in the landscape within a field are readily apparent. Such classifications are potentially useful for stratifying fields for efficient sampling of soil carbon. These landform units within the study area included:

(i) the main ridge (MR), an area of higher ground and an obvious former terrace or levee;

(ii) the road ridge (RR), an area of lower ground than the main ridge that runs along the access road;

(iii) southern flat (SF), an area of lower ground south of the main ridge;

(iv) lower slope (LS), an area of lower slope which borders the edge of the paddock;

(v) depression (D), an internal depression within the paddock that the corer showed to have iron and magnesium nodules at the 20-30 cm layer, and so generally has poorer drainage than the rest of the paddock.

To give an indication of the characteristics of the field, the means and standard deviations for the sensor data and the soil carbon values were calculated for each of the units of the field and landform classification.

Geostatistical analysis and assessment of different sampling strategies

Nested nugget and exponential variogram models were fitted by maximum likelihood to the observed carbon storage (t/ha) over depths 0-10, 10-20, 20-30, 30-50, 0-30, and 0-50 cm. These variograms were to be used to determine the expected errors for different design-based strategies following the approach of Domburg et al. (1994). One of the design-based strategies considered in this paper stratifies the sampling according to a landform classification. From Eqn 7, we see that in this situation a different variogram is required for each stratum. There were insufficient data to fit independent variograms for each of these strata. We therefore assumed that all of the strata had a common spatial parameter r and a common nugget to total sill ratio [c.sub.o]/ ([c.sub.0] + [c.sub.1]. We allowed the variance to vary between the strata. This model was fitted by maximum likelihood.

Design-based strategies

Two tests of different design-based strategies were conducted. The first considered the use of the different quadrat designs. The quadrats considered were based on a radial design with cores located in eight compass directions from a central point (Murphy et al. 2013). The cores were bulked before laboratory analysis. The smallest quadrat was referred to as the 20-m design. The distance of the core from the central point was randomly selected to be 10, 15, or 20m (Fig. 2). The same design was scaled to a 40-m quadrat and a 100-m quadrat. Strategies without local clustering were also considered.

The quadrats were positioned randomly in the field. The approach of Domburg et al. (1994) and the variograms fitted to the field observations were used to determine the expected standard errors upon using different numbers of quadrats to estimate the carbon stocks from 0-10 and 0-30 cm. The standard errors for the different quadrat designs were plotted against (i) the number of quadrats and (ii) the total costs of conducting the design-based strategy.

The second test compared the performance upon positioning quadrats by simple random sampling or stratified random sampling for 20-, 40-, and 100-m radius designs. The most cost-effective quadrat design from the first test was used. Three different types of stratification were considered. These were the spatial stratification suggested by Bins et al. (1999), stratification according to landform unit, and stratification according to a sensor output (the yield map from 2008), which was among the most strongly correlated with carbon stocks. For the landform stratification, each stratum contained at least one quadrat. Other quadrats were allocated to strata according to their relative area. The sensor data were used to divide the field into three, equally sized strata. These were sampled randomly. Again, the expected standard errors were determined according to Domburg et al. (1994) and were plotted against both the number of quadrats and the complete cost of the survey.

Estimation of costs of sampling

The costs of several potential sampling protocols or methodologies were determined. The full set of potential options for sampling are summarised in Table 1, where some of the basic questions that need to be addressed when sampling soil carbon at the field scale are presented. The potential options were selected to give different estimates of soil carbon levels with differing accuracies or uncertainty and requiring different levels of cost.

The estimation of costs was based on AUD$ in November 2011. The estimation of costs was similar to that used by Brus et al. (1999). The total cost was the sum of the staff fieldwork costs, equipment hire, and laboratory analyses as below:

[C.sub.tot] = [C.sub.fld] + [C.sub.eqp] + [C.sub.lab]

The cost of the field work and the cost of equipment are proportional to t, the time in hours (h) taken to conduct the survey and the hourly rates are denoted [C.sub.fld] and [C.sub.eqp]. Sampling for soil carbon to 30 cm requires a soil-coring machine, and [C.sub.eqp] was estimated to be $62.5/h in November 2011. Sampling using a soil-coring machine requires a trained operator and an assistant, which includes a Senior Technical Officer at $72/h and a Technical Officer at $60/h. The cost of the GPS unit was considered negligible. Hence, [C.sub.fld] was equal to $132, the sum of the labour rates.

The time taken to conduct the survey is the sum of the time taken to select the location of the random samples for the cores or the quadrats using a GIS or other system ([t.sub.pln]); the time to travel to these locations ([t.sub.acc]); the time to set up the randomised sampling pattern within the quadrats ([t.sub.qdr]); the time to take the cores, cut up the cores, and complete the sampling including bulking ([t.sub.sam]); and the time taken to travel the distance to and from base to the field site ([t.sub.dist]). We have not included [t.sub.pln] in our calculations, and while this would add to the cost, it is generally considered to be a small cost compared with the field, equipment, and laboratory costs. The time spent travelling from the office to the field ([t.sub.dist]) was not included since this is a fixed cost that is independent of the design-based strategy and is highly variable depending on the location of the field relative to the office.

The assumed [t.sub.acc] is a function of the number of quadrats to be visited. We assumed that it would generally take 5 rain to locate and travel to quadrat sites. However, once one quadrat had been located, one nearby could be reached in 2 min. If the number of quadrats was <20, then they would be sparsely distributed across the field and all of the travel times would be 5 min. Any quadrats beyond the initial 20 could be found in 2 rain. It was also assumed that [t.sub.acc] and all of the other sampling costs were independent of the type of stratification used in the sample design.

The assumed values of [t.sub.qdr] for each quadrat depended upon the quadrat design. There are two possible methods of laying out the quadrats:

(i) The standard grid pattern as recommended in Australian monitoring protocols. This takes ~45 min to set up a grid of 100 cells, of which 10 are selected at random for sampling (Bowman et al 2009).

(ii) The radial sampling system used in our sampling tests (Fig. 2) as used in the Lachlan project (Lachlan Catchment Management Authority 2011). This method takes ~20 min to set up for quadrats with radius 20 or 40 m, but this time is estimated to double if the presence of vegetation prevents the use of tape measures in setting up quadrats with 100-m radius.

Costs are estimated only for the radial sampling system using eight cores per quadrat. The sampling time, [t.sub.sam], can be difficult to predict since it depends on conditions. For the 30 cm core, low and high times of 5 and 10 rain per core are suggested. A value of 7.5 min per core was used in our calculations. For the 10 cm core, a time of 0.5 min per core was assumed since this is a more rapid process. The effectiveness and time taken for the sampling is dependent on the moisture content and condition of the soil. For this paper, the assumption was made that a sampling program would not be undertaken unless the soil was moist and in a good condition for sampling.

The laboratory costs ([C.sub.lab]) were calculated using the relationship:

[C.sub.lab] = [C.sub.soc] [n.sub.soc] + [] []

where [C.sub.soc] is the cost of analysing for soil carbon per sample, [n.sub.soc] is the number of samples analysed for soil carbon, [] is the cost of one measurement for bulk density, and [] is the number of bulk density measurements. These costs were estimated at $20 per sample for [C.sub.soc] based on costs from local laboratories and $50 per measurement for [] based on recent sets of measurements, including the use of personnel time and facilities. While bulk density measurements do not require any sophisticated equipment, the measurement can be labour-intensive. We assume that labour accounts for all of the laboratory costs of measuring bulk density, whereas 50% of the costs of measuring carbon are associated with labour and 50% with capital equipment.

A spreadsheet to calculate the cost of each of the different sampling strategies is available from the authors on request.


Summary statistics and characteristics of the field

Summary statistics for the field are shown in Table 2. The carbon stocks in each layer were not substantially skewed, and therefore, these variables have not been log-transformed. Average carbon stocks of 21.29 and 53.09t/ha were measured in the 0-10 and 0-30cm layers, respectively. The CVs for carbon (%) for the different layers were all >30%, whereas CVs for bulk density were all <10%.

A summary of the largest correlations between the soil carbon variables and sensor data for the field is given in Table 3. We note that since the sampling was not conducted at random, it is not possible to confirm the significance of these correlations. The spatial distributions of the sensor data are shown in Fig. 3. Carbon stocks are highly correlated with carbon percentage, indicating a minor role of bulk density in determining the variability in soil carbon stock. The various sensor data are, at best, only weakly correlated with the soil carbon variables. Apparent EC has a weak relationship with soil carbon in the 0-10cm layer and in the 20-30cm layer. The radiometric data also have weak correlations with carbon in these two layers. There is probably some pedological cause for these relationships, such as parent material and texture effects. The multispeetral SPOT imagery has a weak relationship with soil carbon in the 0-10 cm layer but a stronger relationship to soil carbon in the 20-30cm layer. Total yield (or more specifically the yield in 2008) is the property most strongly correlated with the 0-30cm carbon stocks. While these relationships are possibly associated with plant growth, again they are too weak to be of sufficient value in linear models for the spatial prediction of soil carbon. The soil carbon stocks in individual layers all have strong relationships to cumulative soil carbon stores, as expected.

Table 4 shows that the field and landform classification does not have a substantial impact on soil carbon stores in the 0-30 cm layer, and has a small effect in the 0-10 cm layer. However, there are some weak relationships between some units and sensor data. For example, the depression area is distinguished by yields, the EM data, and the radiometric data. The eastern section of the field appears to have a different radiometric signal to the western section. All of these correlations are fairly weak and we expect there to be little potential for the field and landform classification to improve the sampling design.

Geostatistical analysis of field-scale survey

There are stark differences between the estimated variograms of carbon stocks over the different depth increments (Fig. 4, left). The 30-50 cm layer has the longest effective range of close to 600 m. This might be because carbon variation in this layer is largely unaffected by management and plant growth as seen closer to the surface. The 0-10 and 20-30 cm layers have similar variograms with effective ranges of ~150 m. In contrast, the variogram of the 10-20 cm layer is close to pure nugget. We hypothesise that this is because the 10-20cm layer shows no spatial relationships to either plant growth or the pedological influences of the site. This is in contrast to the 0-10 cm layer, which shows some relationships to the plant growth parameters, and the deeper layers (20-30 and 30-50 cm), which show relationships with the pedological parameters of the site (see Table 3). The variogram for the 0-30 cm layer has a relatively small effective range of ~100 m and a substantial nugget to total variance ratio of ~0.5

Various features are evident in the corresponding maps of carbon stocks (Fig. 4, right). In the west of the field there is an area with low carbon stocks in the 0-10 cm layer but large stocks in the 30-50 cm layer. This maybe a consequence of cut and fill during erosion and deposition events associated with the breakout stream which flows on the south-west edge of the paddock. It is possible that surface soil materials have become buried adjacent to this stream. There are areas of low stocks at 31-50 cm in the north-east of the field. The maps of the 10-20, 20-30, and 0-30 cm layers are more uniform, reflecting the large nugget effects and small effective ranges of the fitted variogram models.

Tests of cost effectiveness of design-based sampling strategies

Bulking the soils within quadrats reduces the standard error for a given number of laboratory analyses of soil carbon and bulk density compared with the use of single cores (Fig. 5). This is not unexpected given the substantial short-scale variation evident in the estimated variograms. When the standard errors are plotted against costs, the advantages of bulking in quadrats are slightly diminished. This is because of the substantial time required to set up a quadrat and to drill multiple cores. The balance of laboratory costs to field costs shifts as the sampling methodology is shifted from single cores to bulking based on quadrats. This can be seen in Table 5 where laboratory costs make up 80-90% of the costs under the randomised single cores but <50% of the costs when cores are bulked in quadrats. The quadrats with 40-m radius are more cost-effective than those with 20-m radius since they account for more of the short-range variation without substantially increasing costs. Quadrats with 100-m radius are even more cost-effective, even though they take longer to set out.

Positioning the quadrats according to a spatially stratified design (Brus et al. 1999) leads to smaller standard errors than simple random sampling (Fig. 6). Stratification according to landform unit actually increases the standard errors. This is because there was no significant relationship between landform classes and carbon stocks, and because the need to sample within all of the landform units led to a relatively non-uniform design. The performance of the stratification based on 2008 yield data was similar to simple random sampling.

These results indicate that to obtain an estimate of the soil carbon store for the 0-30 cm layer with a standard error <2 t/ha would cost ~$2500. To obtain an estimate of the soil carbon store for the 0-10 cm layer with a standard error of 1 t/ha would cost $1400.


Our study has found that the cost of estimating the 0-30 cm carbon store across the field with a standard error of 2 t/ha is ~AU$2500. If this finding were repeated across New South Wales then a soil carbon monitoring scheme would lead to considerable expense. However, there are several reasons why this result might not be representative of all fields. The CV of 30% for carbon recorded in our survey is larger than most of the CV values recorded in the previous surveys listed in the Introduction. This large variability is probably a consequence of the alluvial origin of the soils and the influence of soils of different age and texture as well as the influence of land management. Kravchenko et al. (2006) showed that the spatial variability of soil carbon is highest for soils under no-till cropping systems such as those in the study field. Furthermore, the relationships between carbon stocks and the various sensor data are weaker than we might expect, and the spatial correlation in the 0-30 cm carbon stocks is particular weak. The weak relationships between landform units and carbon stocks might reflect the complex pedology of the field and that the strongest differences in pedology occur below 30 cm. Pedological relationships might be stronger on bedrock-derived soils or in fields wither stronger contrasts in parent materials and landforms (e.g. Miklos et al. 2010; Qin et al. 2012; Simbahan et al. 2006; Huang et al. 2007).

Given the weak spatial correlation and weak correlation with auxiliary information, there was little potential for stratification to improve the efficiency of sampling. We expect there will be more potential on other fields, and it is important that further surveys are conducted to see if this is the case. Some benefits of bulking within quadrats were observed. This reflects the substantial short-scale variation of carbon. In terms of cost effectiveness, we might have expected to see a larger difference due to the reduced laboratory costs, but some of this benefit was negated by the time required to set out each quadrat. The cost of drilling to 30 or 40 cm is substantial and large savings can be made if carbon is only measured to 10 cm. These costs can be reduced further if a hand core is used to sample to 10 cm in good sampling conditions (moist soil) and the cost of hiring a corer is eliminated (see Table 5). However, it is likely that carbon stored at depth can represent a large proportion of the sequestration potential of the soil.

Future studies need to consider whether the quadrats can be laid out less precisely and more quickly, and whether other methods can make field sampling more efficient. For example, it might be practical to bulk the soils from multiple quadrats and reduce laboratory costs. Since bulk density is more expensive to measure and less variable than carbon, there is potential to measure bulk density less frequently than carbon and further reduce costs. This idea is supported by Holmes et al. (2011), who concluded that bulk density has only a small effect on the variation in soil carbon store, and confirms the need for more investigations into the possibility of minimising the proportion of costs allocated to measuring bulk density when estimating soil carbon stores. The soil in this study had a low capacity to shrink and swell, as did the Tenosols, Chromosols, and Kandosols investigated by Holmes et al. (2011). Further investigation is required to determine whether bulk density sampling costs might similarly be reduced for soils with significant swelling potential.

This study considered only how sampling should be conducted to measure the status of carbon stocks. The longterm aim of monitoring schemes is to measure the change in carbon stocks. However, we will only be able to assess the suitability of sample designs to perform this task once we have re-surveyed carbon so that we are able to estimate the variogram of carbon change. Goidts et al. (2009) considered this problem and suggest that the repeated sampling of the same quadrats is a statistically efficient way in which change can be monitored. Lark (2009) recommended the use of paired sites to detect changes in soil carbon levels over time.

In this study, we determined the intensity of sampling required within our field by estimating variograms of carbon and bulk density. It is not practical to estimate variograms in every field. However, it is hoped that future work will demonstrate that the within-field variability of carbon stocks, the form of variograms, and hence the sampling requirements are strongly related to known factors such as soil type, parent material, and land management history.


We have set out a logical framework to determine the cost effectiveness of different design-based sampling strategies for the estimation of carbon stocks. On the study field, the cost was AU$2500 to estimate 0-30 cm carbon stocks with a standard error of 2 t/ha. The study field is of alluvial origin and has complex pedology. Therefore, it might be more variable than most fields, and an accurate estimation of its carbon stocks more expensive. Other paddocks on different soil types and parent materials and with different land management histories may have different patterns of variation and, hence, different sampling requirements. It is important to repeat this procedure on other fields to determine how representative these costs are.

More information is required on the nature of carbon variability across the landscape to give informed advice regarding the use of quadrats and their optimal size and the general potential of stratifications based on landform units or sensor data to reduce sampling costs. We hope that similar variograms will be observed within soil mapping units or soil types and that these factors can be used to guide sampling requirements. 10.1071/SR12119

Received 10 May 2012, accepted 1 January 2013, published online 5 February 2013


We thanks Angus O'Brien for allowing us to sample his paddock so intensely; Michael Short for EC measurements in the field; Professor Alex McBratney, Brett Whelan, and Budiman Minasny for advice on planning the field work and analysis; David Duncan for advice on paddock selection, sampling, and coring; Ian Packer for cost estimates on bulk density measurements; Aaron Simmons and Elizabeth Warden for advice on sampling costs in the field; Andrew Rawson, Warwick Badgery, Jason Crean, and lan Packer for general discussion and comments on the manuscript; Kate Lorimer-Ward, manager of the CAMBI Project, and the NSW Catchment Action program for support. Ben Marchant's work was funded by Biotechnology and Biological Sciences Research Council (BBSRC). He is grateful to the University of Sydney for the award of an International Visiting Research Fellowship.


Allen DE, Pringle MJ, Page KL, Dalai RC (2010) A review of sampling designs for the measurement of soil organic carbon in Australian grazing lands. The Rangeland Journal 32, 227-246. doi: 10.1071/RJ09043

Bock M, Bohner J, Conrad O, Kothe R, Ringeler A (2007) Methods for creating functional soil databases and applying digital soil mapping with SAGA GIS. In 'Status and prospect of soil information in south-eastern Europe: soil databases, projects and applications'. EUR 22646. (Eds T Hengl, P Panagos, A Jones, G Toth) pp. 149-162. (EN Scientific and Technical Research Series, Office for Official Publications of the European Communities: Luxemburg)

Bowman G, Chapman GA, Murphy BW, Wilson BR, Jenkins BR, Koen T, Gray JM, Morand DT, Atkinson G, Murphy C, Murrell A, Milford HB (2009) 'Protocols for soil condition and land capability monitoring.' NSW Natural Resources Monitoring, Evaluation and Reporting Program. (Department of Environment, Climate Change and Water NSW: Sydney) Available at: soils/SoilsProtocols.pdf

Bricklemyer RS, Miller PR, Paustian K, Keck T, Nielsen GA, Antle JM (2005) Soil organic carbon variability and sampling optimization in Montana dryland wheat fields. Journal of Soil and Water Conservation 60, 42-51.

Brus DJ, Spatjens LEEM, Gruijter JJ (1999) A sampling scheme for estimating the mean extractable phosphorus concentration of fields for environmental regulation. Geoderma 89, 129-148. doi: 10.1016/S00167061(98)00123- 2

Conant RT, Paustian K (2002) Spatial variability of soil organic carbon in grasslands: implications for detecting change at different scales. Environmental Pollution 116, S127-S135. doi:10.1016/S0269-7491 (01)00265-2

Davis AA, Stolt MH, Compton JE (2004) Spatial distribution of soil carbon in southern New England hardwood forest landscapes. Soil Science Society of America Journal 68, 895-903. doi:10.2136/sssaj2004.0895

de Gruijter J, Brus D, Bierkens M, Knotters M (2006) 'Sampling for natural resource monitoring.' (Springer-Verlag: Berlin)

Domburg P, de Gruijter JJ, Brus DJ (1994) A structured approach to designing soil sampling schemes with prediction of sampling error from variograms. Geoderma 62, 151-164. doi:10.1016/0016-7061(94) 90033-7

Duncan D, Wooldridge A, Brennan N, Agar B, Murphy B, Welch A, Andersson K, Kellett J, Lawrie J, Kew G, Andersson K (2008) Soil Information Package, Nyngan 1:250 000 Sheet. Central West Catchment Management Authority and NSW Department of Environment and Climate Change.

Goidts E, Van Wesemael B, Crucifix M (2009) Magnitude and sources of uncertainties in soil organic carbon (SOC) stock assessments at various scales. European Journal of Soil Science 60, 723-739. doi:10.1111/ j. 1365-2389.2009.01157.x

Holmes KW, Wherret A, Keating A, Murphy DV (2011) Meeting bulk density sampling requirements efficiently to estimate soil carbon stocks. Soil Research 49, 680-695. doi:10.1071/SRl 1161

Huang X, Senthilkumar S, Kravchenko A, Thelen K, Qi J (2007) Total carbon mapping in glacial till soils using near-infrared spectroscopy, Landsat imagery and topographical information. Geoderma 141, 34-42. doi: 10.1016/j.geoderma.2007.04.023

Knowles TA, Singh B (2003) Carbon storage in cotton soils of northern NSW. Australian Journal of Soil Research 41, 889-903. doi:10.1071/ SR02023

Kravchenko AN, Robertson GP, Bullock DG (2006) Management practice effects on surface total carbon: Differences in spatial variability patterns. Agronomy Journal 98, 1559-1568.

Lachlan Catchment Management Authority (2011) Soil Carbon Pilot Project. Lachlan Catchment Management Authority, NSW, Australia. Available at: Pilot/Fact_Sheets/FactSheet_1_Soil_Carbon_Pilot_Overview.pdf, and Sheets/Fact_Sheet_3_Soil_Carbon_calculations_for the_pilot.pdf

Lark RM (2009) Estimating the regional mean status and change of soil properties: two distinct objectives for soil survey. European Journal of Soil Science 60, 748-756. doi: 10.1111/j. 1365-2389.2009.01156.x

McKenzie N, Ryan P, Fogarty P, Wood J (2000) Sampling, measurement and analytical protocols for carbon estimation in soil, litter and coarse woody debris. National Carbon Accounting System Technical Report No. 14, Australian Greenhouse Office, Canberra.

Miklos M, Short MG, McBratney AB, Minasny B (2010) Mapping and comparing the distribution of soil carbon under cropping and grazing management practices in Narrabri, north west NSW. Australian Journal of Soil Research 48, 248-257. doi: 10.1071/SR09111

Minty B, Franklin R, Milligan P, Richardson M, Wilford J (2009) The radiometric map of Australia. Exploration Geophysics 40, 325-333. doi: 10.1071/EG09025

Murphy B, Badgery W, Simmons A, Rawson A, Warden E, Andersson K (2013) Soil testing protocols at the paddock scale for contracts and audits. Market Based Instruments for Soil Carbon (CAMBI). Version 5.3. New South Wales Department of Primary Industries. Orange, Australia.

Qin C-Z, Zhu A-X, Qiu W-L, Lu Y-J, Li B-L, Pei T (2012) Mapping soil organic matter in small low-relief catchments using fuzzy slope position information. Geoderma 171-172, 64-74. doi: 10.1016/j.geoderma.2011. 06.006

Sanderman J, Farquhason R, Baldock J (2010) Soil carbon sequestration potential: A review for Australian agriculture. Report for the Australian Department of Climate Change and Energy and Efficiency. CSIRO, Division of Land and Water, Urrbrae, S. Aust.

Shukla MK, Slater BK, Lal R, Cepuder P (2004) Spatial variability of soil properties and potential management classification of a chernozemic field in Lower Austria. Soil Science 169, 852 860. doi:10.1097/ 00010694-200412000-00004

Simbahan GC, Dobermann A, Goovaerts P, Ping JL, Haddix ML (2006) Fine resolution mapping of soil carbon based on multivariate secondary data. Geoderma 132, 471 489. doi:10.1016/j.geoderma.2005.07.001

Walvoort DJJ, Brus DJ, de Gruijter JJ (2010) An R package for spatial coverage sampling and random sampling from compact geographical strata by k-means. Computers & Geosciences 36, 1261-1267. doi: 10.1016/j.cageo.2010.04.005

Webster R, Oliver MA (2007) 'Geostatistics for environmental scientists.' (Wiley: Chichester, UK)

WRB (2006) 'World reference base for soil resources 2006.' 2nd edn. IUSS Working Group. World Soil Resources Reports No. 103. (FAO: Rome)

K. Singh (A,E), B. W. Murphy (B), and B. P. Marchant (C,D)

(A) Faculty of Agriculture and Environment, University of Sydney, Sydney, NSW 2006, Australia.

(B) Office of Environment and Heritage, Cowra, NSW 2794, Australia.

(C) Rothamsted Research, Harpenden, Hertfordshire AL5 2JQ, UK.

(D) Current address: British Geological Survey, Keyworth, Nottingham, NG12 5GG, UK.

(E) Corresponding author. Email:

Table 1. Sampling strategies for estimating soil carbon levels
at the paddock scale Some of the decisions that affect the
accuracy of the estimate and the cost of the sampling strategy
are included

Is the paddock                                    Yes
  stratified based on     This can be done by observing features in the
  observable soil,        field, interviewing the landholder for land
  landform and land       use history, and using available natural
  use features?           resource information manually, or using
                          digital soil mapping  methodologies
                          (as demonstrated by Simbahan et al. 2006,
                          Huang et al. 2007, Qin et al. 2012)

Is the paddock                                    Yes
  stratified based on
  a spatial
  statistical design
  (e.g. Brus et al.
How is the location       Completely              All cores contained
  of cores in the           randomised (A)          within quadrats.
  paddock determined?       within each             Quadrats located
                            stratification          at randomised
                            unit                    points within

Number of cores or        5-100                   3, 5, 10, 20, 30
Radius of the             Not applicable          20, 40, 100 m
  quadrat (m)
Sample design within      Not applicable          Radial (rapid)
  the quadrats?                                     or grid (B)

Number of cores within    Not applicable          3,5,8,10
  each quadrat?

Bulking of cores          Not applicable          Bulked or not
  within each quadrat?                              bulked

Depth of cores?           10 or 30 cm             10 or 30 cm

Is the paddock                                    No
  stratified based on
  observable soil,
  landform and land
  use features?

Is the paddock                                    No
  stratified based on
  a spatial
  statistical design
  (e.g. Brus et al.
How is the location       Completely              All cores contained
  of cores in the           randomised (A)          within quadrats.
  paddock determined?       within the whole        Quadrats located
                            paddock                 at randomised
                                                    points within the
                                                    whole paddock

Number of cores or        5100                    3, 5, 10, 20, 30
Radius of the             Not applicable          20m, 40m, 100m
  quadrat (m)
Sample design within      Not applicable          Radial (rapid) or
  the quadrats?                                     grid (B) (slower)

Number of cores within    Not applicable          3, 5, 8, 10
  each quadrat?

Bulking of cores          Not applicable          Bulked or not
  within each quadrat?                              bulked

Depth of cores?           10 or 30 cm             10 or 30 cm

(A) Tbere are several options for selecting randomised points within
paddocks or stratification units that can affect the accuracy of
estimates of soil carbon and the costs of making the measurements
(see Brus et al. 1999).

(B) Not tested in this paper.

Table 2. Summary statistics of the field-scale survey

Depth                    Mean              Coefficient of
(cm)                                        variation (%)
               Soil organic carbon (%)

0-10                     1.46                   30.8
10 -20                   1.08                   44.1
20-30                    0.89                   44.9
30-50                    0.76                   56.5

               Bulk density (t/[m.sup.3])

0-10                     1.46                    9.6
10-20                    1.60                    5.6
20-30                    1.66                    6.0
30-50                    1.68                    6.0

               Soil carbon store (t/ha)

0-10                     21.29                  30.9
10-20                    16.99                  35.8
20-30                    14.81                  43.3
0-30                     53.09                  23.2

Table 3. Informative correlations between sensor data and carbon
measurements at the 100 sites across the field

SOC, Soil organic carbon; ECh, ECv; electrical conductivity
horizontal and vertical modes; K, potassium; Tdose, total
dose; Th, thorium; Wl, wetness index; SPOT band 3, satellite
image near infrared; NDVI, normalised difference vegetation
index; PVR, photosynthetic vigour ratio

                          SOC%                    SOC store
                  0-10 cm    0-10 cm    0-20 cm    0-30 cm    0-50 cm

SOC%               1.000      0.936      0.763      0.579      0.087
0-10 cm
ECh               -0.246     -0.315     -0.185     -0.047      0.270
ECv               -0.255     -0.330      0.195      0.062      0.254
K                 -0.217     -0.244     -0.086      0.058      0.344
Tdose             -0.195      0.223     -0.096      0.024      0.297
Th                -0.179     -0.219     -0.105      0.011      0.271
W1                 0.099      0.200      0.258      0.249      0.165
SPOT Band 3        0.158      0.175      0.060      0.036     -0.148
NDVI               0.149      0.174      0.101      0.056     -0.012
PVR               -0.250     -0.164     -0.167      0.113      0.259
Sum yield         -0.012     -0.020      0.004      0.145      0.162
Yield 2008         0.042     -0.001      0.055      0.190      0.187
Bulk density      -0.155      0.182      0.109      0.118     -0.060
0-10 cm

Table 4. Mean of sensor data and soil carbon data within each
landform classification unit

ECv, Electrical conductivity vertical mode; SPOT band 3, satellite
image near infrared; Tdose, total dose; K, potassium; Th, thorium

                        Eastern pedological section of paddock

                     Main           South flat        Depression

                           Apparent electrical conductivity

ECv                   48.27             50.39             13.38

                             Yield (t/ha) and plant growth

Yield 2005            2.91              2.65              2.02
Yield 2008            2.83              2.48              2.28
SPOT Band 3          74.56             83.88             92.14

                            Radiometric spectrophotometry

Tdose (nGy/h)        73.93             73.67             68.96
K (%)                 1.09              1.10              1.03
Th (ppm)              6.98              6.36              5.74

                            Soil carbon store (t/ha)

0-10 cm               22.61             24.88             24.01
0-30 cm               55.14             58.62             54.69

                  Western pedological section of paddock

                      Road ridge            Lower slope

ECv                     78.40                  58.06

Yield 2005               2.76                   2.54
Yield 2008               2.72                   2.18
SPOT Band 3             79.88                  83.86

Tdose (nGy/h)           87.23                  83.61
K (%)                    1.29                   1.25
Th (ppm)                 9.20                   8.61

0-10 cm                 19.71                  18.64
0-30 cm                 55.13                  49.61

Table 5. Total costs for different sampling strategies and
percentage allocated to different procedures

For sampling strategy: nSI, n individual cores whose locations have
been selected by simple random sampling; nQ(m), n quadrats randomly
located with m cores bulked within the quadrat. Travel refers to
travel between sites within the paddock and finding the sites using a
GPS; setting up sites refers to laying out quadrats and finding sites
within quadrats. SOC, Soil organic carbon; the cost of laboratory
measurement of SOC was assumed to be 50% equipment and 50% labour
preparation costs; BD, bulk density

Sampling                          Field (%)
strategy       Labour: travel     Labour:      Hire of        Sum
               and setting up     coring       corer
                                                            0-30 cm
100SI                2.4            4.7          3.3          10.5
50SI                 3.0            4.6          3.6          11.2
10S1                 4.5            4.5          4.3          13.3
1OQ(10)             12.1            24.3         17.2         53.7
1OQ(5)              12.4            15.5         13.2         41.0
                                                            0-10 cm
100S1                7.1            1.4          4.0          12.6
50S1                 8.6            1.3          4.7          14.6
IOSI                12.5            1.3          6.5          20.3
1OQ(10)             32.8            6.6          18.7         58.1
1OQ(5)              30.8            3.8          16.4         51.0

Sampling                         Laboratory (%)
strategy        SOC: capital          SOC:            BD       Sum
                and equipment        labour         (100%

100SI               12.8              12.8           63.9      89.5
50SI                12.7              12.7           63.4      88.7
10S1                12.4              12.4           61.9      86.6
1OQ(10)              6.6               6.6           33.1      46.3
1OQ(5)               8.5               8.5           42.1      59.0

100S1               12.5              12.5           62.5      87.4
50S1                12.2              12.2           61.0      85.4
IOSI                11.4              11.4           56.9      79.7
1OQ(10)              6.0               6.0           29.9      41.9
1OQ(5)               7.0               7.0           35.0      49.0

Sampling        Total
strategy         cost

100SI           23464
50SI            11829
10S1             2424
1OQ(10)          4531
1OQ(5)           3559

100S1            8005
50S1             4100
IOSI             878
1OQ(10)          1673
1OQ(5)           1429
COPYRIGHT 2012 CSIRO Publishing
No portion of this article can be reproduced without the express written permission from the copyright holder.
Copyright 2012 Gale, Cengage Learning. All rights reserved.

Article Details
Printer friendly Cite/link Email Feedback
Author:Singh, K.; Murphy, B.W.; Marchant, B.P.
Publication:Soil Research
Article Type:Report
Geographic Code:8AUNS
Date:Nov 1, 2012
Previous Article:Soil organic carbon content and storage of raised field wetlands in different functional zones of a typical shallow freshwater lake, China.
Next Article:Particulate organic matter in soil under different management systems in the Brazilian Cerrado.

Terms of use | Privacy policy | Copyright © 2021 Farlex, Inc. | Feedback | For webmasters |