Spatial variability of chemical and physical attributes in a soil from La Pampa flat Santafesina

The processes and attributes that influence crop performance vary in space and time. Its magnitude and spatial structure is specific for each lot and its quantification is necessary for the application of site-specific crop management (MSEC). Our objective was to examine the intra-plot spatial variability of chemical and physical attributes of the surface horizon of a typical Argiudol from downtown Santa Fe for the implementation of the MSEC. In a plot of 475 x 100 m, 60 samples were taken distributed in 3 transects of 20 points separated every 25 m and 50 m between transects. In each disturbed sample (0-20 cm) nitrate N (N-NO3), total organic C (TOC), total N (Nt) was determined; Extractable P (P), sand, clay, soil reaction (pH). The bulk density of the soil (Ds) was determined from undisturbed samples from depths 0-10 and 10-20 cm. Descriptive statistics were calculated and the spatial dependence of the attributes was analyzed using empirical variograms and their limits obtained by permutation (envelopes). The coefficients of variation (CV) for TOC, Nt, Sand, Clay, Ds 0-10, Ds 10-20 and pH were low (CV < 15%) while for P and N-NO3 they were moderate (15 < CV < 35%). The spatial dependence of N-NO3 was moderate with a range of 96 m, while the rest of the attributes did not present autocorrelation due to the variation that occurs at distances shorter than those used in the sampling.

Keywords : Argiudols; Geostatistics; Site-specific management.


The processes and properties that regulate crop performance and yield vary in space and time. The assessment of the magnitude and spatial structure of their variation is required in order to apply site-specific crop management (SSCM). Our aim was to assess the intra-field spatial variability of some chemical and physical topsoil properties of a typical Argiudoll from the center of Santa Fe for the use of the SSCM. Sixty soil samples of the top layer (0-20 cm depth) were taken in a 475 x 100 m plot. Sampling was distributed in three transects separated by 50 m with 20 points every 25 m. Disturbed soil samples were used to analyze N of nitrates (NNO 3), total organic C (COT), total N (Nt); available P (P), sand (Sand), clay (Clay) and soil reaction (pH). Soil bulk density (Ds) was determined on non-disturbed soil samples taken at 0-10 and 10-20 cm depths. Summary statistics were calculated and the spatial dependence was determined with sample variograms and envelopes obtained by permutations. Coefficients of variation (CV) for TOC, Nt, Sand, Clay, Ds 0-10, Ds 10-20, and pH were low (CV < 15%) while they were moderate for P and N-NO3 (15 < CV < 35%). Moderate spatial dependence with a range of 96 m was only observed for N-NO3. The other soil properties had no autocorrelation because their variation occurred at distances shorter than 25 m. In this study,

Keywords : Argiudol; geostatistic; Site-specific Management



Soil processes and attributes that influence crop performance vary in space and time (Mulla & Schepers, 1997). Site-specific crop management (MSEC), one of the forms of Precision Agriculture (PA), aims to adjust the application of inputs and agronomic practices to the needs of the soil and crops based on their spatio-temporal variation. (Whelan & McBratney, 2000). Thus, both the magnitude and the spatial structure of the variability within the flock are necessary aspects for the application of the MSEC while the temporal variability tends to make its implementation difficult (Pierce & Nowak, 1999). Pringle et al.
Soil variability is itself scale dependent as it responds to the action of formative factors that act on a continuum of spatial and temporal scales, resulting in nested variation structures (Trangmar et al., 1985). On a regional scale, climate, land use patterns, vegetation type, and relief characteristics are the main determinants of variation (Mallarino & Vittry, 2004). At the plot scale, topography is one of the main factors causing spatial variability through the control of water and sediment distribution (Ceddia et al., 2009). Even in very gently rolling landscapes to widespread plains, the presence of microreliefs determines the formation of contrasting soil complexes a few meters away (Hein et al., 1989). In cultivated soils, management practices such as furrow orientation, nutrient application method, tillage, and compaction can dominate the causes of variability at small scales (Mallarino & Vittry, 2004), altering chemical, physical, and biological attributes of the soil. superficial horizon (Corá et al., 2004).
Since its first applications to soil science in the early 1980s (Campbell, 1978; Burgess & Webster, 1980a and b; Webster & Burgess, 1980; Burgess et al., 1981, Vieira et al., 1981), the Geostatistics has been used for the description and quantification of the spatial variability of the soil through the variography of its attributes (McBratney & Pringle, 1999). Results of these studies, covering various soil and climatic conditions and observation scales, show that soil properties frequently present spatial autocorrelation of different magnitude (Warrick et al., 1986; Tragmar et al., 1985; Whelan, 1998; Mulla & McBratney, 2000 ).
The spatial structure of the variability of the edaphic attributes is specific for each plot (Mallarino, 1996). However, in the absence of previous studies at a site of interest, local references of similar conditions (ie soil type, topography, etc.) are useful to guide sampling (Webster & Oliver, 2007). In sites where there is a description of the spatial continuity of an attribute, it is possible to optimize the sampling intensity based on the desired variability in the predictions, ie kriging variance (McBratney & Webster, 1981).
At the national level, the spatial variability of different soil attributes has been studied on different types of soil and using different sampling schemes (Di Pietro et al., 1985; Correa & Sosa, 1998; Melchiori, 2000; Kemerer & Melchiori, 2004; Zubillaga et al., 2004; Cruzate & Rivero, 2010; Gilli et al., 2010; Rivero et al., 2010), although for the center of the province of Santa Fe the information is scarcer (Alesso & Pilatti, 2008). .
The flat pampa of Santa Fe is characterized by the presence of flat areas with “subnormal” relief and gently undulating areas or with “normal” relief, both with micro relief. Our objective was to examine the spatial variability of chemical and physical attributes of the surface horizon of a production lot in the central area of ​​Santa Fe with normal relief, cultivated in direct sowing with agricultural rotation and application of fertilizers, to verify the feasibility of implementing the site-specific crop management system (MSEC).


The study was carried out in the La Pelada district, Santa Fe province (30° 56′ 17″ S, 60° 55′ 23″ W). The region is characterized by its very gently undulating relief with a predominance of Argiudoles and associated Argialboles; material of loessic or silt-loessic origin and a subhumid-humid climate with annual rainfall that varies between 920 and 1689 mm and a mesothermal regime with an average annual temperature of 19 ºC (De Petre et al., 1977).
A lot was selected for continuous agriculture in direct sowing with a sequence of soybean-corn crops and moderate doses of nitrogenous and phosphorous fertilizers. Within a sector with a general slope of 0.7 %, a systematic sampling was carried out within a 475 x 100 m plot dominated by typical and aquic Argiudols, according to the semi-detailed Soil Map of the province of Santa Fe (INTA, 2009). With a Pentax R 326 EX total station, three parallel transects were drawn 50 m apart with 20 points each (n = 60) every 25 m (Fig. 1). At each site, the height relative to the lowest point (Z, m) was recorded and a disturbed sample was taken from the first 0-20 cm, made up of five subsamples taken in a radius of 2-3 m, and undisturbed samples with cylinders of known volume at 0-10 and 10-20 cm depth.

Figure 1 . Sampling design and altimetry (exaggeration factor in Z axis = 30). In the upper circle it represents the detail of the composite sample taking.
Figure 1 . Sampling design and elevation (Z axis expand factor = 30). The magnified circle shows details on the composite sampling method.

Nitrate N (N-NO 3 , mg kg -1 ) was determined for each disturbed sample by the phenoldisulfonic acid method (Jackson, 1982); Total organic C (TOC, g kg -1 ) estimated from easily oxidizable carbon Walkley & Black (1934) with recovery factor 0.77; Total N (Nt, g kg -1 ) by the Kjeldahl method (Jackson, 1982); extractable P (P, mg kg -1 ) according to Bray & Kurtz (1945); sand between 50 µm and 2 mm (Sand, g kg -1 ) by wet sieving and clay < 2 µm (Clay, g kg -1) by the hydrometer method without pretreatments according to Gee & Or (2002); and soil reaction (pH) by the potentiometric method in soil-water ratio 1:2.5 (Jackson, 1982). The bulk density of the soil (Ds, Mg m -3 ) was determined for the undisturbed samples using the 100 cm -3 cylinder method (Forsythe, 1975).
Univariate summary statistics and exploratory plots were calculated for each attribute using the R statistical package (R Development Core Team, 2011). The presence of outliers was examined using box plots, removing those that were 1.5 times the interquartile range with respect to one of those quartiles. Subsequently, the assumption of normality was verified using the Shapiro-Wilks test and transformations for non-normal distributions were evaluated. Second-order stationarity was verified for each attribute by fitting polynomials based on coordinates, altimetry, or other attributes as covariates. In those cases where the trend explained more than 20% of the total variance (ie
The geostatistical analysis was performed using the geoR package (Ribeiro Jr. & Diggle, 2001). The spatial dependence of the attributes was analyzed using experimental variograms with spacing intervals (lag) of 25 m over an active distance of 200 m with a minimum of 30 pairs per lag (Goovaerts, 1998). The robust estimator of the semivariances ( Equation 1 ) proposed by Cressie & Hawkins (Cressie, 1993) was used:


where: Y(h) is the semivariance, N(h) is the total number of pairs at a distance h, and Z(xi) and Z(xi+h) pairs of observations separated by a distance h.
The minimum and maximum limits of the empirical variograms (envelopes) were calculated by means of permutations of the data on the coordinates under the assumption of spatial independence (Diggle & Ribeiro Jr., 2007). In the cases in which the empirical variogram was located close to or outside the simulated limits, the presence of spatial autocorrelation was corroborated by fitting theoretical variograms by the restricted maximum likelihood method (REML) and subsequent likelihood ratio test ( likelihood ratio test) with respect to a non-spatial model (pure nugget). Spherical semivariance functions (Equation 2 ) and exponential ( Equation 3 ), retaining the model with the lowest Akaike Information Index (AIC).



where: Y(h) is the semivariance between pairs separated by a distance h, c0 is the unstructured variance (nugget), c1 structured variance (partial sill), h is the distance between pairs, and roa is the range of spatial dependence (m).
For the classification of the spatial structure of the attributes, the nugget-plateau relationship (nugget:sill) was used, ie proportion of unexplained variance with respect to the total variance (c0/(c0+c1)), adopting three classes proposed in Cambardella et al. (1994): weak (> 0.75), moderate (0.25-0.75), and strong (< 0.25).
The predictive ability of the model was evaluated using the cross-validation procedure (Webster & Oliver, 2007). From the pairs of observed and predicted values, the square root of the mean square error (RMSE, Equation 4 ), the mean square deviation ratio (MSDR, Equation 5 ) and the G criterion ( Equation 6 ) of goodness of fit were calculated. .




where: is the number of observations used in the cross-validation, z(xi) is the observed value of z at the ith position, Z(xi) is the value estimated by the model for the ith or 2 OK position, variance kriging of the prediction for the point xi , Z mean value of Z in the region.


Table 1 presents the summary statistics of the soil attributes studied after the removal of two outliers of N-NO3, pH and Sand, and a value of Ds 10-20, respectively. The bias and kurtosis coefficients did not differ significantly from those expected for a Normal distribution. In Figure 2The empirical distributions of each attribute are observed together with the normal distribution (1:1 line) and their 95% confidence intervals. For pH and Clay, some observations were located outside the confidence intervals of a Normal distribution verified with the Shapiro-Wilks test (P < 0.05). Given the slight positive asymmetry of the pH, the logarithmic transformation of this variable was evaluated, while for Clay the observations were concentrated around two modes, although not distant from each other (ie 190 and 210 g kg -1 ), possibly due to the analytical technique errors. The transformations explored did not allow the values ​​of these variables to be approximated to a Normal distribution, and we proceeded to use the data without transformation.

Table 1 . Summary statistics of the attributes of the shallow horizon (depth 0-20 cm) of a typical Argiudol from central Santa Fe.
Table 1 . Summary statistics of topsoil properties (0-20 cm depth) of a typical Argiudoll from central Santa Fe.

Figure 2 . Normal probability plots and 95% confidence intervals (95% CI) of the attributes of the surface horizon (depth 0-20 cm) in a typical Argiudol from downtown Santa Fe.
Figure 2 . Normal probability plots with a 95% confidence interval (95% CI) of topsoil properties (0-20 cm depth) of a typical Argiudoll from central Santa Fe.

The variability found for TOC, Nt, Sand, Clay, Ds 0-10, Ds 10-20, and pH was low (ie CV < 15%) and moderate for P and N-NO 3(ie 15 < CV < 35%), according to the classification proposed in Wilding & Drees (1983). These results are within the CV reported in the bibliography. DiPietro et al. (1986) in a micro-scale study obtained lower CVs for P, pH, Nt, TOC and Clay of a typical Argiudol from Pergamino; while Álvarez et al. (2004) found similar results to our study at the plot scale on Argiudols typical of the Pampa Ondulada under direct sowing, except for P whose variability was almost 10 times greater. The increase in the sampling area increases the total variance, although the relationship between this and the different scales do not follow consistent patterns (Wilding & Drees, 1983).
The Nt showed a slight negative trend depending on the relative height (R2 = 0.22), observing the highest values ​​in the lowest positions of the studied area. The empirical variograms for this attribute were calculated after subtracting said trend to satisfy the second-order stationarity assumption.
For most of the attributes, the empirical variograms were located within the limits obtained by simulation under the assumption of spatial independence ( Fig. 3). The absence of spatial correlation observed at this work scale indicates that the resolution used was not appropriate to describe the spatial continuity of these attributes due to the occurrence of spatial variation at distances less than the sampling interval (Webster & Oliver, 2007).

Figure 3 . Empirical variograms (circles) and limits (dashed line) of the attributes of the surface horizon (depth 0-20 cm) in a typical Argiudol from downtown Santa Fe. The semivariance values ​​have been scaled to the sample variance.
Figure 3 . Sample variograms (circles) and envelopes (dashed lines) of topsoil properties (0-20 cm depth) of a Typic argiudoll from central Santa Fe. Variogram values ​​have been scaled to the sample variance.

For N-NO3, Nt and Arena, the empirical variograms in the first distance interval were close to the limits, although only the theoretical model fitted for NNO 3 differed from the hypothesis of spatial independence (P < 0.05). For this attribute, the structured variance represented 74% of the total variance with a spatial dependence range of 96 m ( Fig. 4.a ). Figure 4.b presents the results of the cross-validation procedure. The model obtained allowed the estimation of the N-NO 3 content with an error of ± 2.8 mg kg -1 , which resulted in a 27% increase in accuracy with respect to the estimation of the N-NO 3 content.from the average of the studied area. However, the dispersion diagram shows a marked tendency to overestimate and underestimate the low and high values ​​of N-NO 3 respectively, which translates into a moderate level of correlation. The relationship between the estimation error and the kriging variance show that the fitted model slightly overestimates the latter (MSDR < 1).

Figure 4 . Theoretical variogram model of the N-NO3 content of the superficial horizon (0-20 cm) of a typical Argiudol from downtown Santa Fe: (a) Empirical (circles) and theoretical (line) variogram fitted by REML; and (b) scatterplot of the observed vs. predicted values ​​obtained by the cross-validation procedure.
Figure 4 . Theoretical variogram model of soil N-NO3 concentration (0-20 cm depth) of a typical Argiudoll from central Santa Fe: (a) Sample variogram (circles) and model fitted by REML (line); (b) scatterplot of observed and predicted values ​​from cross-validation procedures.

The scale and design of the sampling, ie number and distribution of the samples, strongly condition the results obtained from the empirical variograms (Diggle & Ribeiro, 2007). These factors should be considered when comparing these results (Stenger et al., 2002). DiPietro et al. (1986), using spacings 10 times smaller than those used in our study, found similar results for Nt and pH, while for P, COT and Clay these authors found high, medium and low degrees of structuring, according to Cambardella et al. (1994), with ranges of spatial dependence of 6.3, 9.8 and 15 m respectively. Alesso & Pilatti (2008) found no spatial autocorrelation for Ds in an Argiudol in downtown Santa Fe with subnormal relief.
Attributes such as granulometry can be controlled by intrinsic variations of soil characteristics, while extrinsic variations such as the application of fertilizers and amendments or tillage can affect the spatial structure at close range of related attributes (Cambardella et al., 1994). Mallarino (1996) found an association between the variability patterns of P and K at short distances and the non-uniform distribution of fertilizers and manures.
The presence of microreliefs even in very gently undulating areas determines the formation of soil complexes in spaces between 15 and 20 m (Hein et al., 1989). This type of variation at distances lower than the resolution used in this study could explain the absence of spatial autocorrelation observed in the Clay and Sand contents. On the other hand, the variations at distances of less than 25 m of the nutrient levels can be attributed to the non-uniform application and removal of fertilizers, while the traffic of machinery can cause variations in the density of the soil.
For the effective application of the MSEC, both the magnitude and the spatial structure of the variability of the attributes that determine the response of the crop are necessary (Pierce & Nowak, 1999). The greater the magnitude of variation, the greater the probability of differentiating input doses, while a strong spatial structure determines smooth and broad patterns of variability that facilitate the application of variable doses (Pringle et al., 2003). In our study, most of the attributes examined presented low levels of variation (ie CV < 20%) that would not justify the differential application of inputs.
The variability observed for P was moderate, but the lack of autocorrelation of the observations at this scale does not allow the application of interpolation techniques, the mean being the best estimator (Webster & Oliver, 2007). An increase in the resolution of the grid would improve the description of spatial continuity at short distances to obtain a reliable map prediction (Kravchenko, 2003).
Finally, the spatial structure observed for NNO 3 allowed to obtain prediction maps, although the accuracy and precision were not adequate. On the other hand, the temporal variability and scant spatial correlation frequently observed in N-NO 3 variograms limit the usefulness of this indicator for the development of maps of application of variable rates of N (Stenger et al., 2002).


The spatial variability of physical and chemical properties of the superficial horizon of a soil from the Santa Fe plain pampas was evaluated. With the exception of N-NO 3 , whose spatial structure was moderate with a dependency range of 96 m, the rest of the evaluated attributes did not present autocorrelation, indicating the occurrence of variation at distances less than the observation scale used. The low variability and the lack of spatial structure did not allow modeling the spatial distribution of the attributes or obtaining soil condition maps that allow guiding the application of the MSEC technique.


To UNL-CONICET for financial support. This work is part of the thesis work of Carlos Agustín Alesso to fulfill the requirements of the Doctorate in Agricultural Sciences of the National University of Córdoba and the CONICET postgraduate scholarship program. To the ranch La Pelada and especially to Ing. Agr. M Lui for his support and good disposition. To Dr. P Cipriotti for his willingness to answer questions.


  1. Alesso, A & MA Pilatti. 2008. Spatial variability of soil density, extractable phosphorus, and exchangeable potassium. XXI Argentine Congress of Soil Science. Potrero de los Funes, Argentina.
  2. Alvarez, R; H Steinbach; B Bauschen & JN Ejalbert. 2004. Spatial variability of soil properties in the rolling pampas: effect on the number of subsamples to take for the diagnosis of fertility. XVI Argentine Congress of Soil Science. Parana, Argentina.
  3. Bray, R & L Kurtz. 1945. Determination of total, organic and available forms of phosphorous in soils. Soils Sci. 59: 29-45. [  Links]
  4. Burgess, TM & R Webster. 1980th. Optimal interpolation and isarithmic mapping of soil properties: I. The semi-variogram and punctual kriging. J. Soil Sci. 31(2): 315-331. [  Links]
  5. Burgess, TM & R Webster. 1980b. Optimal interpolation and isarithmic mapping of soil properties: II. Block kriging. J. Soil Sci. 31(2): 333-341. [  Links]
  6. Burgess, TM; R Webster & AB McBratney. 1981. Optimal interpolation and isarithmic mapping of soil properties. IV. Sampling strategy. J. Soil Sci. 32(4): 643-659. [  Links]
  7. Cambardella, CA; T. B. Moorman; JM Novak; TB Parkin; D. L. Karlen; Turkish RF & AE Konopka. 1994. Field-scale variability of soil properties in central Iowa soils. Soil Sci. Soc. Am. J. 58: 1501-1511. [  Links]
  8. Campbell, JB. 1978. Spatial variation of sand content and pH within single contiguous delineations of two soil mapping units. Soil Sci. Soc. Am. J. 42: 460-464. [  Links]
  9. Ceddia, MB; SR Scallop; HELLO Villela; LS Speck; LHC Anjos & DF Carvalho. 2009. Topography and spatial variability of soil physical properties. Sci. Agr. 66(3): 338-352. [  Links]
  10. Korah, JE; AV Araujo; GT Pereira & JMG Beraldo. 2004. Spatial variability of attributes only for adoption of the precision agriculture system in the sugarcane culture. Rev. Bras. Science. Only 28(6): 1013-1021. [  Links]


Leave a Comment