Published in Soil Sci. Soc. Am. J. 68:2031-2041 (2004).
© 2004 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
Division S-8—Nutrient Management & Soil & Plant Analysis
Site-Specific Soil Fertility Management
A Model for Map Quality
T. G. Mueller,
N. B. Pusuluri,
K. K. Mathias,
P. L. Cornelius and
R. I. Barnhisel
Dep. of Agronomy, N-122 Agronomy Science North, Univ. of Kentucky, Lexington, KY 40502
* Corresponding author (mueller{at}uky.edu)
 |
ABSTRACT
|
|---|
The performance of site-specific fertility management (SSFM) systems depends on the quality of soil property maps used to develop variable-rate fertilizer recommendations. Map quality assessment, however, may be too expensive for routine site-specific soil sampling. The objectives of this study were (i) to evaluate the quality of soil property maps created with ordinary kriging for five fields in Kentucky, and (ii) to develop a model describing the relationship between map quality and statistical properties of data. Five fields across Kentucky were sampled on 30.5-m grids and samples were analyzed for pH, buffer pH (bpH), P, K, Ca, and Mg. For each field, four 61.0 and nine 91.5-m data subsets were extracted from the 30.5-m grid. Semivariograms could only be adequately modeled for the 30.5- and 61.0-m grid datasets. Therefore, only these data sets were interpolated with ordinary kriging. Map quality was evaluated with an independent data set. Multiple stepwise regression was used to model map quality using data from several Kentucky fields and from a previously published Michigan study. Prediction efficiency (PE) was a function of the relative structural variability, range of spatial correlation, and grid increment (R2 = 0.82). The range of spatial correlation was the major factor controlling map quality within the range of variation studied. This model may potentially be a useful tool for the development of sampling designs for site-specific management.
Abbreviations: bpH, buffer pH • PE, prediction efficiency • SSFM, site-specific fertility management • VIDS, validation with an independent data set
|
INTRODUCTION
|
|---|
SITE-SPECIFIC FERTILITY MANAGEMENT is based on the premise that the quantity of fertilizer and lime producing the maximum economic crop response varies spatially and temporally within agricultural fields in ways that may be adequately predicted and managed. While many factors influence economic crop response, recommendations for P and K fertilizer and lime applications are based primarily on soil test values. Therefore, the quality of soil property maps is fundamental to site-specific P, K, and lime management (Sawyer, 1994; Pierce and Nowak, 1999). Inaccurate or imprecise soil property maps may explain why some SSFM studies have produced poor results (Mueller et al., 2001). Consequently, map quality evaluation is critical for assessing or predicting the performance of site-specific P, K, and lime management.
Map quality is usually determined by comparing mapped and observed values with various quantitative and qualitative analytical techniques. Plots of predicted vs. measured values should always be visually examined to assess prediction quality. If the maps are of high quality, the scatter of data will adhere closely to the 1:1 line. Residual (i.e., predicted minus measured) maps or semivariograms of the residuals may be examined to determine whether and how prediction errors vary spatially. Some of the quantitative measures include map precision, which is the standard deviation of the residuals. Map accuracy is the square of bias, and bias is the average of the residuals. The MSE is the sum of precision and accuracy. While these measures are informative, they are scaled. Therefore, errors cannot be meaningfully compared between variables and, in some cases, across locations or time.
Some measures of map quality can be compared across variables, locations, and time. Correlations between predicted and measured values have been used in this way; however, they only assess the regression relationship between predicted and measured and not the deviation of predicted and measured values from a 1:1 line. The PE, a reduction-in-error index, is the standard measure of the quality of spatial predictions (Agterberg, 1984; Gotway et al., 1996; Kravchenko and Bullock, 1999; Mueller et al., 2001; Kravchenko 2003).
Map quality assessment is impacted by the methods used to calculate predicted and measured values. Two general techniques are used. Cross-validation with replacement is a rapid, inexpensive procedure for comparing predicted and measured values. It unfortunately does not adequately describe spatial prediction errors in many situations (Isaaks and Srivastava, 1989; Mueller et al., 2001). Consequently, some have urged that cross-validation results be interpreted with caution (Isaaks and Srivastava, 1989; Cressie, 1993; Goovaerts, 1997; Mueller et al., 2001). Validation with an independent data set (VIDS) is a superior and more dependable method for measuring residuals. Although this procedure has sometimes been referred to as jackknife analysis (e.g., Deutsch and Journel, 1998; Kravchenko, 2003), we refer to it as VIDS because the standard jackknife approach is a delete-one method (Hinkley, 1983). Unfortunately, VIDS is expensive and, therefore, may not be practical for routine site-specific soil sampling. Alternative approaches for assessing map quality are needed.
Currently, there are no existing statistical models for predicting the quality of soil property maps. However, studies have shown that the degree (Kravchenko, 2003) and range (Mueller, 1998) of spatial autocorrelation and sampling intensity (Wollenhaupt et al., 1994; Gotway et al., 1996; Mueller et al., 2001, Kravchenko, 2003) impact map quality. Other factors such as sampling design (e.g., orientation of samples in the field, subsampling, and compositing procedures), laboratory procedures (e.g., laboratory error and quality control), and mapping procedures (e.g., interpolation method and parameterization, smoothing, contour frequency) may also impact map quality (i.e., PE).
The development of models for map quality could be very beneficial for many areas of study if predictions are sufficiently accurate. In agriculture, they might have utility as tools for assessing and implementing site-specific management systems. There were two objectives in this study. The first, to evaluate the quality of soil property maps created with ordinary kriging for five fields in Kentucky. The second, to use the results of these evaluations to develop a mathematical model describing the relationship between map quality and statistical properties of the mapped variable.
|
MATERIALS AND METHODS
|
|---|
To achieve a broad inference space, five fields were chosen from various agriculturally important and diverse physiographic regions in the Commonwealth of Kentucky: the Mississippian Embayment [Calloway County location (36°32'3''N, 88°27'7''W], Mississippian Plateaus [Hardin County location (37°31'7''N, 85°56'13''W)], Bluegrass Region [Shelby (38°1'32''N, 85°27'40''W) and Nelson County (37°39'45''N 85°30'11''W) locations], and Western Kentucky Coal Fields [Hopkins County location (37°26'16''N, 87°26'5''W)]. First-order soil surveys were created for the Hardin and Calloway County locations (Mueller et al., 2004). Digital second-order soil surveys were obtained online for Hopkins (Soil Conservation Service, 1977) and Shelby (Soil Conservation Service, 1980) County from the National SSURGO database (Natural Resource Conservation Service, 2003). The study area in the second-order Nelson County soil survey (Soil Conservation Service, 1971) was digitized. Then it was rectified and georeferenced with ArcGIS (Redlands, CA).
All of the fields were in no-till or reduced-till cropping systems immediately before sampling. The fields had been in 2-yr corn (Zea mays L.)–soybean [Glycine max (L.) Merr.] or corn–wheat (Triticum aestivum L.)–double-crop-soybean rotations for 10 to 20 yr with the exception of the Hardin County location that had been in pasture for more than 10 yr until 1 yr before sampling.
Within and across locations, there were considerable differences in drainage class, erosion and slope phases, and taxonomy (Table 1). The Calloway and Hopkins soils formed in loess or silty alluvium. The soils from the Hardin and Shelby locations were mainly derived from limestone residuum overlain by loess, except for small areas occupied by the Lindside and Nolin soil series, which developed in mixed alluvium. The soils of the Nelson location developed either in limestone, shale, sandstone, siltstone, or some combination of these parent materials.
Soil samples were obtained from these fields at points on 30.5-m regular grids (Fig. 1)
. Soil samples were not collected from areas that were out of production at the time of sampling including grass waterways (Calloway County location), sinkholes (Hardin County location), and gullies (Hopkins, Nelson, and Shelby County locations). The total number of sampling points (n) per field was as follows: Calloway 163, Hardin 104, Hopkins 265, Nelson 248, and Shelby 588. Additional validation data sets (Calloway 59, Hardin 53, Hopkins 71, Nelson 65, and Shelby 70) were obtained randomly within regular grids that varied in size between fields. Our rationale for varying the size of these grids was to obtain a minimum of 60 validation points, a number that was found to be adequate in a similar study (i.e., Mueller et al., 2001). This goal was nearly achieved for most of the locations. At each grid point, five subsamples (one at the grid point and four within a 7-m radius) were obtained using a 2.1-cm-diam. core to a depth of 20 cm, and these samples were combined to form a composite at each grid or sample point. Soils were air dried at 25°C and ground to pass a 2-mm sieve. Standard soil analyses were conducted by the soil testing laboratory at the Department of Regulatory Services, University of Kentucky. Analyses included pH (1:1 soil-to-water mixture), bpH [Shoemaker-Mc-Lean-Pratt (SMP) buffer] P, K, Ca, and Mg (Mehlich III extractable) (University of Kentucky, 2003). Subsets were extracted from the original 30.5-m grid dataset (i.e., four 61.0-m and nine 91.5-m grids). For each location, an additional dataset was created by combining the validation data values and the 30.5-m grid data values. These datasets will be referred to as the FULL datasets.
View larger version (21K):
[in this window]
[in a new window]
|
Fig. 1. Sampling design. The 30.5-m grids are indicated with crosses and the validation points are indicated with circles.
|
|
Semivariance values were calculated with SAS (Cary, NC) procedure VARIOGRAM for each of the datasets. Additionally, semivariogram values were calculated for several Michigan datasets from Mueller et al. (2001) to fairly compare the Michigan and Kentucky map statistics. The Michigan datasets included one 30.5-m grid, four 61.0-m grids, one 100-m grids, one validation dataset, and one FULL dataset. Soil test results (i.e., soil pH, bpH, P, K, Ca, and Mg) were analyzed from each dataset. Semivariogram values were determined for each 3.05-m lag. This small lag interval was used because it improved the resolution of the semivariogram behavior at smaller lag distances (i.e., <30.5 m for the FULL dataset) and some of the major lag distances (e.g., grid increment, the square root of two times the grid increment, and two times the grid increment). Unfortunately, the use of small lag intervals also resulted in considerable noise for larger lags. To remove this noise, lag intervals with less than a threshold number of observed pairs were ignored. The thresholds were reduced for smaller lag distances (i.e., <30.5 m) and increased for datasets with more observations.
Semivariograms quality for the 91.5- (five Kentucky datasets) and 100-m (Michigan dataset) grids were generally very poor. Therefore, they were not modeled. Exponential models were fit based on visual assessment. Parameters included the nugget, partial sill, and range of spatial correlation. The total sill is the sum of the nugget and partial sill. The range of spatial correlation refers to the practical range. For exponential models, the practical range is the lag distance at which the semivariance is equal to 95% of the total sill. The relative structural variability is the ratio of the partial to total sill. For more detailed descriptions of these parameters, see Isaaks and Srivastava (1989) and Schabenberger and Pierce (2002). For each of the six FULL datasets (including the Michigan data from Mueller et al., 2001), semivariogram values for individual pairs were calculated with Variowin (Golden Software, Golden, CO). These values were used to assess differences between modeled and observed semivariogram values for individual pairs.
Interpolated grids (1 by 1 m) were calculated with ordinary kriging (Surfer, Golden Software, Golden, CO) for the 30.5- and 61.0-m grid datasets. The VIDS was conducted by obtaining predicted values for each measured value in the validation datasets using bilinear interpolation. Errors associated with bilinear interpolation were very small, as indicated by preliminary analyses and other studies (Mueller et al., 2001).
The predicted and measured values were used to calculate various measures of map quality. The MSE was the sum of accuracy and precision. It was defined in Eq. [1].
 | [1] |
where vi was the difference between predicted value and observed value at location si, i = 1, ..., nv, and nv was the number of values in the check data set. The RMSE was defined as the square root of the MSE. Accuracy (the square of bias) was defined in Eq. [2]
 | [2] |
and precision (the variance of the residuals) was defined in Eq. [3]
 | [3] |
where 
represents the mean of the residuals. An accurate map will have no bias. In conversation, the expression great accuracy and precision implies small errors. However, mathematically, accuracy and precision have been defined as errors in Eq. [2] and [3]. To avoid confusion, accuracy (Eq. [2]) and precision (Eq. [3]) will be referred to as errors of accuracy and errors of precision in this paper. The PE, a reduction-in-error index, was defined in Eq. [4].
 | [4] |
where MSEfield average was the MSE obtained by using the field average values (obtained from the 30.5-m grid datasets) as estimate for all validation data points. This index has was refered to as goodness by Agterberg (1984).
Stepwise regression (Proc REG, SELECTION = Stepwise; SLENTRY = 0.5; SLSTAY = 0.15) was used to model PE. Diagnostics for assessing multicolinearity (i.e., intercept adjusted variance-inflation factors and condition indices statistics) were calculated by using the VIF and COLLINOINT model options. The INFLUENCE option was used to calculate the RStudent totalistic for identifying potential outliers. The PRESS option was used to perform leave-one-out validation (sometimes referred to as cross-validation) to the uncertainty of the regression model.
Because there were multiple observations for some of the grid increments (e.g., four 61.0-m grids), the mean values of these observations were used to develop the model. If individual observations had been used, the 61.0-m grids to have a substantially greater influence on the estimates of model parameters than the 30.5-m grids. Candidate regressor variables for entry in the model included the sampling interval, the coefficient of variation, and geostatistical parameters. The natural logs, squares, cubes, and fourth powers of these variables (with the exception of sampling intensity) were also used as candidate variables. To avoid multicolinearity problems, only a variable (e.g., CV) or one of its transformations (e.g., natural log of CV, the square, cube, or fourth power of CV) was allowed to remain in the model. Since the goal of this analysis was to relate map quality to the underlying statistical properties of the spatial processes, estimates of the spatial variability statistics were obtained from the FULL datasets and the CV values were obtained from the 30.5-m grid datasets.
Estimated and measured PE values were used to calculate prediction MSE (Eq. [1]), accuracy (Eq. [2]), and precision (Eq. [3]) statistics for the multiple regression model. These comparisons were performed using the statistical parameters and PE values from the following: (i) the data that was generate the model (i.e., individual 30.5-m grids and mean values for the 61.0-m grid), (ii) the raw observations model before summarization (i.e., the individual 30.5- and 61.0-m grid observations), and (iii) published values reported in a map quality simulation study (Kravchenko, 2003).
The simulated values used from the Kravchenko study used conditioning data from three soil test variables measured in an Illinois study. The CV values of the three variables varied considerably (i.e., 12, 40, and 67%). For each of the three variables, conditional simulations were generated with three different nugget-to-sill ratios (i.e., 0.1, 0.3, and 0.6) and with a range of spatial correlation of 97 m (spherical semivariogram model parameters). So that results could be compared across studies, exponential model parameters were fit to the specified spherical models. Our results were compared with reported map quality values (Tables 1, 2, and 3, under the columns entitled "Kriging with spatial structure unknown" in Kravchenko, 2003) which were the mean results from 100 simulations. We only considered reported results for grid increments within (i.e., 40, 50, and 60 m) and nearly within (i.e., 30 m) the range of variability of the data that was used to generate the map quality regression model.
|
RESULTS AND DISCUSSION
|
|---|
Semivariogram Analyses and Map Quality Assessments
The semivariogram model parameters for the FULL data sets (Table 1) were in-line with those reported in other studies (e.g., McBratney and Pringle, 1997; Cambardella and Karlen, 1999; Mueller et al., 2001). In some situations the semivariance values did not reach a plateau and, therefore, intrinsic stationarity was assumed. These cases were indicated in Table 1. Unfortunately, the parameters for these models had little utility for agronomic interpretation because an infinite combination of range and sill values could have been used to fit to the empirical semivariograms. Therefore, the variables that exhibited the intrinsic behavior were only interpolated. They were not, however, included in the development of the multiple regression model describing PE variability.
Prediction efficiency values ranged between –62 and 83%. Other published studies have generally reported values within this range (i.e., Gotway et al., 1996; Kravchenko and Bullock, 1999; Mueller et al., 2001; Kravchenko, 2003). Differences between studies may not only have been due to differences in spatial structure and sampling intensity but also to differences in laboratory analysis and other sampling methods. For example, in this study, five subsamples were collected within a 7.0-m radius and composited. In two Michigan studies, five subsamples were collected within a radius of 1.5 m (i.e., Mueller et al., 2001; Mueller and Pierce, 2003). The impact of sampling methods on map quality deserves further consideration.
Prediction efficiency values for the 30.5-m grids were fairly large (i.e., between 70 and 83%) for approximately 13% of the variables (i.e., Mg at Hardin County field, pH and bpH at Hopkins County field, and P for Shelby County field) and marginal to poor (i.e., between –21 and 30%) for the 23% of the variables. Managing at the 30.5-m grid-scale would not have been practical for the SSFM of corn–soybean–wheat rotations. This scale, however, could potentially have been realistic for some high-value crops.
The performance of kriging using 61.0-m grids was of greater practical interest because site-specific management may have been economically feasible at this scale if maps were of adequate quality. Surprisingly, prediction efficiencies were relatively large in some cases (e.g., pH and bpH for Hopkins County, and P for the Shelby County data sets; Table 2).
Because there were 150 plots of predicted vs. measured (five locations, five prediction datasets per location, and six variables), only plots were created for pH, P, and K from the 61.0-m grids for the Hopkins (Fig. 2)
and Shelby County (Fig. 3)
fields. These were chosen because pH, P, and K were considered important management variables, management was potentially feasible at this scale, and PE values were relatively large for these fields.
View larger version (27K):
[in this window]
[in a new window]
|
Fig. 2. Plots of predicted vs. measured for pH, P, and K for the four 61.0-m grids (i.e., G61a, G61b, G61c, G61d) at the Hopkins County location.
|
|
View larger version (29K):
[in this window]
[in a new window]
|
Fig. 3. Plots of predicted vs. measured for pH, P, and K for the four 61.0-m grids (i.e., G61a, G61b, G61c, G61d) at the Shelby County location. PE = prediction efficiency.
|
|
Unfortunately, there were not predictable patterns associated with specific variables. For example, the clusters of points for pH in the Hopkins County field (Fig. 2) and P in the Shelby County field (Fig. 3) more closely adhered to the 1:1 line than for P and K in the Hopkins County field (Fig. 2) and pH and K in Shelby County field (Fig. 3). The inconsistent behavior of specific soil test results can also be observed in Table 3.
For the 61.0-m grids, as average prediction quality increased, PE became less variable. For example, prediction efficiencies were stable across the various 61.0-m grid datasets for pH in Hopkins County (Fig. 2) and P in Shelby County (Fig. 3.). However, they were not as stable for P and K values at Hopkins County location (Fig. 2) and pH and K values at the Shelby County location (Fig. 3), which were of lower average PE. A mapping scale that is adequate for management for a particular variable should not depend substantially on the orientation of the grid in the field.
Map Quality Model Development
During preliminary regression analysis, the initial model had an R2 value of 0.67 and included four significant factors (i.e., the natural log of the relative structural variability, the range of spatial correlation, sample intensity, and the natural log of the coefficient of variation). The RStudent statistic for one of the observations in the dataset was 2.7. Because RStudent values > 2.0 are indicative of potential outliers (SAS Institute, 2001; Schabenberger and Pierce, 2002), the observation was temporally deleted. During subsequent and successive iterations, an inordinate number of pH and bpH observations were temporally removed based to this criterion compared with other variables. While this approach for identifying potential outliers is an accepted tool for deleting outliers, it usually results in the removal of too many data observations (Schabenberger and Pierce, 2002). Therefore, the RStudent statistic was used only as a diagnostic tool. The preliminary regression analyses indicated potential problems with some of the observations, particularly for pH and bpH.
Plots of difference between modeled and actual semivariance values vs. lag distance revealed spatial anomalies for some variables (i.e., pH and bpH from the Calloway, Hardin, Hopkins, and Nelson locations) (e.g., Fig. 4a)
. These patterns either did not occur or occurred for P, K, Ca, and Mg at all five locations (e.g., Fig. 4b). The trends were minimal for pH and bpH at the Shelby location (e.g., 4c). These perceived spatial trends may have been due to trends in the chemical analysis of the data. The soil samples from the Shelby location were sent in a random order to the soil testing laboratory. However, the samples from the other locations were sent to the laboratory in an order that corresponded to the grid sample numbers, and hence spatial location in the field. It appears that laboratory assessments of pH were likely influenced by the preceding pH measurement during analysis. Because the analyses were performed in the order in which the samples occurred in the field, the correlation between sequential analyses distorted measurements of spatial autocorrelation. This distortion may have explained the large RStudent statistics for pH and bpH during the initial regression analyses. Therefore, pH and bpH data from the Calloway, Hardin, Hopkins, and Nelson datasets were not used to develop the final model. Clearly, randomizing samples before laboratory analysis is of critical importance, particularly when the spatial structure of the data is in question. However, the actual impact of these errors on map quality is not well understood. Note, for example, that PE values for pH and bpH at the Hopkins County location were large despite the problems with the semivariogram data.
View larger version (19K):
[in this window]
[in a new window]
|
Fig. 4. Model semivariance minus measured semivariance of individual semivariograms of individual pairs for (a) pH in Hardin County, (b) P in Hardin County, and (c) pH in Shelby County. PE = prediction efficiency.
|
|
To offset the substantial loss of degrees of freedom in the regression analysis caused by the removal of the pH and bpH observations, the Michigan dataset from the Mueller et al. (2001) study was used for model development rather than for model validation as had been previously planned. However, pH and bpH in the Michigan dataset exhibited the same spatial trends that occurred for pH and bpH for four of the five Kentucky locations. Therefore, these values were not used in the development of the final regression model. Clearly, problems with the errors associated with pH and bpH were not laboratory specific. It may be possible to avoid some of these problems by adjusting the settings of the pH meter (i.e., increasing the stability thresholds before a final pH reading is given). However, this may not be practical for commercial laboratory analyses.
The final multiple regression model (Table 4) explained 82% of the variability in observed PE (Table 3). The data used to develop the model included observations from a Michigan study (Mueller et al., 2001) and excluded observations for which the semivariograms did not reach a clear plateau in the analysis of the FULL data set (Table 2). Additionally, pH and bpH data from the Calloway, Hardin, Hopkins, Michigan, and Nelson locations were also excluded. Multicolinearity was not problematic, as indicated by small intercept adjusted variance inflation factors (Table 3) and condition indices. Regression parameters can be affected when variance inflation factors (Gunst and Mason, 1980) and condition indices (Belsley et al., 1980; SAS Institute, 2001) are >10.
The model (Table 4) was consistent with the work of others who have demonstrated or reported that the degree of spatial autocorrelation (Kravchenko, 2003), the range of spatial correlation (Mueller, 1998), and sampling intensity (Wollenhaupt et al., 1994; Gotway et al., 1996; Mueller et al., 2001, Kravchenko, 2003) impacted map quality. The CV or its transformations were not significant factors in the final model, consistent with the findings of Kravchenko (2003) but inconsistent with the findings of Gotway et al. (1996), Kravchenko and Bullock (1999), and Mueller et al. (2001). The R2 value (i.e., 0.82) was surprisingly large considering the diversity of soils and that this was a simple three-factor model. Clearly, the factors controlling map quality were not greatly affected by differences in soil parent materials.
Because the model (Table 4) included a nonlinear term (i.e., the cube of the range of spatial correlation), it was difficult to evaluate the relative importance of the three factors in the model based on the coefficients. To visually compare these factors, the partial contributions to the predictions (i.e., regression coefficient x repressor variable value) were plotted against the untransformed values of the regressor variables (Fig. 5)
. The y axes were adjusted to have ranges of variation (i.e., 150%) between minimum and maximum prediction efficiencies. This allowed the impacts of the three factors on map quality to be compared fairly. The ranges of the x axes were based on the maximum and minimum values of the observations used to generate the model. Within the range of the data used to generate the model, the impact of the geostatistical properties of the data (i.e., the range of spatial correlation and the relative structural variability) was substantially greater than the impact of sampling intensity. Of these variables, unfortunately, only sample intensity could be readily managed to improve map quality.
View larger version (24K):
[in this window]
[in a new window]
|
Fig. 5. A visual representation of the contributions of the factors in the regression model (Table 4) to the final estimates. The partial contributions to the predictions (i.e., regression coefficient x repressor variable value) were plotted against the untransformed values of the regressor variables.
|
|
The excellent correspondence between predicted and measured values for the data used to develop the model was not surprising (Fig. 6a)
considering the large model R2 value (i.e., 0.82). Whenever the same data that is used to develop a model is also used to evaluate it, the RMSE statistic may underrepresent the prediction errors that would have been encountered by those who would use the model for prediction. Therefore, leave-one-out cross-validation was conducted (i.e., see the SAS PRESS options described in the methods section). The RMSE value with this approach was 11.1%, which was similar to the model RMSE (i.e., 10.8%). This indicated that the model parameters were fairly insensitive to the individual data values used to develop the model. This is a desirable property for a regression model.
View larger version (20K):
[in this window]
[in a new window]
|
Fig. 6. Regression model (Table 4) evaluation. Modeled vs. measured prediction efficiency for (a) the data that was used to generate the model, (b) the original data with the exception of pH and buffer pH values from the Calloway, Hopkins, Hardin, Nelson, and Michigan locations, the original data from the pH and buffer pH values from the Calloway, Hopkins, Hardin, Nelson, and Michigan locations, and (d) reported results from Kravchenko (2003). RMSEs and errors of accuracy (ACC) and precision (PREC) listed. Accuracy and precision have been defined Eq. [2] and [3], respectively.
|
|
We next assessed the uncertainty associated with individual observations. Recall that the data used to develop the model included the mean prediction efficiencies of the four 61.0-m grids. As expected, errors of prediction were substantially larger when the model was used to predict the values of individual observations (Fig. 6b). Nevertheless, the data adhered fairly well to the 1:1 line on the plots of predicted vs. measured. However, errors of prediction for the pH and bpH (Fig. 6c) were fairly large. The large accuracy errors (i.e., bias = 17.9) occurred because the model underestimated prediction quality for pH and bpH. Predictions were off, presumably because one or both of the geostatistical parameters (i.e., relative structural variability or the range of spatial correlation) were poor estimates because spatial trends were associated with the order of sampling and analysis as previously described (e.g., Fig. 4a).
Prediction efficiency was also underestimated by the model (Fig. 6d) when it was applied to the Kravchenko (2003) simulation study data. This was likely because of methodological differences between the studies. Specifically, Kravchenko (2003) used the validation data sets to optimize the number of points used in ordinary kriging. Therefore, it would have been suppressing if prediction values had not been greater in that study. Whether this accounted for the entire error associated with accuracy (i.e., 319.1) is unknown. Known differences in semivariogram modeling (i.e., exponential vs. spherical) were taken into account (see the Materials and Methods section) and therefore likely did not impact bias. The error of precision for the Kravchenko (2003) data (i.e., 77.4% squared, indicated in Fig. 6d) was even smaller than the error of precision obtained for the data that was used to develop the regression model (i.e., 119% squared, indicated in Fig. 6a). The error of precision was small because the reported prediction efficiencies from the Kravchenko study were the mean prediction efficiencies of 100 observations. The small error of precision and linear behavior of predicted and measured values suggested that the model would likely have performed very well if identical procedures had been used for interpolation and spatial variability assessment.
|
CONCLUSIONS
|
|---|
The model could be used to aid with site-specific sampling design. The first step would be to conduct an initial survey of soil properties potentially along one or more transects. The second step would be to use the model, statistical properties of the transect data, and desired sampling increment to predict map quality. Prediction uncertainty could be taken into account with our reported model RMSE values. Plots of predicted vs. measured and PE values from this and other studies could be used to determine whether estimated prediction efficiencies would be adequate for management. If not, then sampling intensity could be adjusted to improve prediction quality, or alternative sampling methods could also be considered (e.g., zone sampling).
|
ACKNOWLEDGMENTS
|
|---|
The authors acknowledge funding from the USDA (Special Grant Contract No. 99-34408-7561). In addition, we express our appreciation to Mike Ellis, the Luck brothers, Rick Murdock, the Peterson brothers, and Charlie Stuecker for providing access to their farms to conduct this research. We would especially like to thank Danna Reid, Frank Sikora, and other staff of the University of Kentucky soil testing lab for conducting the laboratory analyses for this study. Finally, we gratefully thank Meng Fengxuan for the preparation of the data files used for the analyses in this publication.
Received for publication July 17, 2003.
|
REFERENCES
|
|---|
- Agterberg, F.P. 1984. Trend surface analysis. p. 174–171. In G.L. Gaile and C.J. Willmott (ed.) spatial statistics and models. Reidel, Dordrecht, the Netherlands.
- Belsley, D.A., E. Kuh, and R.E. Welsch. 1980. Regression diagnostics. John Wiley & Sons, New York.
- Cambardella, C.A., and D.L. Karlen. 1999. Spatial analysis of soil fertility parameters. Precis. Agric. 1:5–14.
- Cressie, N. 1993. Statistics for spatial data. John Wiley & Sons, New York.
- Deutsch, C.V., and A.G. Journel. 1998. GSLIB: Geostatistical software library and user's guide. 2nd ed. Oxford Univ. Press, New York.
- Goovaerts, P. 1997. Geostatistics for natural resources evaluation. Oxford Univ. Press, New York.
- Gotway, C.A., R.B. Ferguson, G.W. Hergert, and T.A. Peterson. 1996. Comparison of kriging and inverse-distance methods for mapping soil parameters. Soil Sci. Soc. Am. J. 60:1237–1247.[Abstract/Free Full Text]
- Gunst, R.F., and R.L. Mason. 1980. Regression analysis and its application: A data-oriented approach. Marcel Dekker, New York.
- Hinkley, D. 1983. Jackknife methods. p. 280–287. In S. Kotz and N.L. Johnson (ed.) Encyclopedia of statistical sciences: Vol IV. John Wiley & Sons, New York.
- Isaaks, E., and R. Srivastava. 1989. An introduction to applied geostatistics. Oxford Univ. Press, New York.
- Kravchenko, A.N. 2003. Influence of spatial structure on accuracy of interpolation methods. Soil Sci. Soc. Am. J. 67:1564–1571.[Abstract/Free Full Text]
- Kravchenko, A., and D.G. Bullock. 1999. A comparative study of interpolation methods for mapping soil properties. Agron. J. 91:393–400.[Abstract/Free Full Text]
- McBratney, A.B., and M.J. Pringle. 1997. Spatial variability in soil—Implications for precision agriculture. In J. Stafford (ed.) Precision agriculture. BIOS Scientific Publ., Oxford, UK.
- Mueller, T.G. 1998. Accuracy of soil property maps for site-specific management. Ph.D. Diss. Michigan State Univ., East Lansing, MI (Diss Abstr. 99–22353. Diss. Abstr. Int. 60B:0901).
- Mueller, T.G., N.J. Hartsock, A.D. Karathanasis, and P.L. Cornelius. 2004. Interpreting soil electrical conductivity and terrain attribute variability with soil surveys. Precis. Agric. (In press).
- Mueller, T.G., and F.J. Pierce. 2003. Enhancing spatial estimates of soil carbon with terrain attributes. Soil Sci. Soc. Am. J. 67:258–267.[Abstract/Free Full Text]
- Mueller, T.G., F.J. Pierce, O. Schabenberger, and D.D. Warncke. 2001. Map quality for site-specific management. Soil Sci. Soc. Am. J. 65:1547–1558.[Abstract/Free Full Text]
- Natural Resource Conservation Service. 2003. National SSURGO database data access [Online]. Available at http://www.ftw.nrcs.usda.gov/ssur_data.html [verified 29 July 2004]. USDA and National Cartography and Geospatial Center, Fort Worth, TX.
- Pierce, F.J., and P. Nowak. 1999. Aspects of precision agriculture. Adv. Agron. 87:1–85.
- SAS Institute. 2001. SAS OnlineDoc: Version 8. SAS Inst., Cary, NC.
- Sawyer, J.E. 1994. Concepts of variable rate technology with considerations for fertilizer application. J. Prod. Agric. 7(2):195–201.
- Schabenberger, O., and F.J. Pierce. 2002. Contemporary statistical models for the plant and soil sciences. CRC Press, Boca Raton, FL.
- Soil Conservation Service. 1971. Soil survey of Nelson County, Kentucky. USDA, NRCS, Fort Worth, TX.
- Soil Conservation Service. 1977. Soil survey of Hopkins County, Kentucky. USDA, NRCS, Fort Worth, TX.
- Soil Conservation Service. 1980. Soil survey of Shelby County, Kentucky. USDA, NRCS, Fort Worth, TX.
- University of Kentucky. 2003. Soil testing methods [Online]. Available at http://soils.rs.uky.edu/newpage1.htm [updated 20 May 2004; verified 29 July 2004] Division of regulatory services, Lexington, KY.
- Wollenhaupt, N.C., R.P. Wolkowski, and M.K. Clayton. 1994. Mapping soil test phosphorus and potassium for variable-rate fertilizer application. J. Prod. Agric. 7:441–448.
TOP
ABSTRACT
INTRODUCTION
