PEDOTRANSFER FUNCTION TO ESTIMATE THE SOIL STRUCTURAL “S” INDEX AND SPATIAL VARIABILITY IN AN OXISOL WITHIN A LIVESTOCK FARMING SYSTEM

ABSTRACT

Various attributes are often required to build soil quality indicators. However, determining these attributes is time-consuming and requires several specific devices. Thus, it is desirable to develop indexes that express soil quality based on easily determined attributes. The objective of this study was to investigate the spatial variability of the S index and its correlation with other physical and chemical soil attributes to generate a pedotransfer function for estimating the S index. A georeferenced sampling mesh covering 1.4 ha and employing 71 points was installed. The soil samples were collected in 0.00–0.15 m and 0.15–0.30 m deep layers to determine the physical and chemical attributes. The results indicated that the S index was correlated with porosity, carbon stocks, cation exchange capacity, and particle size fractions. However, macroporosity, microporosity, and sand content were the most suitable attributes for the construction of pedotransfer functions. A principal component analysis indicated that the S index was representative of 9.4% and 11.5% of the total variability in the dataset in the respective soil layers. Spherical semivariogram models showed that the S index was spatially dependent and ranged between 84–188 m. The S index maps estimated by the pedotransfer function resemble the observed S values; therefore, the function can be applied in spatial variability studies.

KEYWORDS
Soil and water management; geostatistics; soil quality index; stepwise procedures

INTRODUCTION

Evaluating soil management by farmers is difficult because obtaining physical attributes of soil is difficult; undisturbed soil samples must be collected and the laboratories used must be enabled in attributes determinations and have trained technicians to interpret the results. In addition, when attribute determination is performed, and one or more have critical levels, those restrictive to plant growth and development result in doubt from farmers about the exact procedures of soil mechanical management. Thus, the scientific community has searched for soil quality indicators that can facilitate decision making by farmers.

As proposed by Dexter (2004 a, b, c), the S index can be calculated from the adjustment coefficients of the soil water retention curve (SWRC). Conceptually, this index represents the slope of the line tangent to the inflection point of the SWRC. Several studies, such as Dexter (2004 a, b, c), Assis Júnior et al. (2016), Magalhães et al. (2018), Rossetti & Centurion (2018), have confirmed the accuracy of the S index using the soil compaction state, which is correlated with the microstructural porosity, bulk density, and soil organic matter content.

Several of the Brazilian Cerrado soils have low or medium clay content and many are considered unsuitable for grain crops cultivation because they have low natural fertility and low organic carbon content, which is responsible for most of a soil's cation exchange capacity (Reichert et al., 2016). These soils are generally used for livestock practices or forest cultivation (silvicultural practices). In addition to the natural soil restriction conditions, the stocking rate and the traffic of animals in the area, which are associated with the deposition of manure, affect the spatial variability in the physical and chemical attributes of soil (Dalchiavon et al., 2017).

Although the S index is widely used as an indicator of the structural quality of soil, no research has reported its spatial dependence on agricultural areas. This is justified by the cost of determining the SWRC at all sampling points in a georeferenced sampling mesh, which is not feasible in commercial agricultural areas. In addition, it may be possible to estimate the S index by pedotransfer functions, as reported by Xu et al. (2017); thus, the index could be determined and mapped using soil attributes that are easier to determine. In Boim et al. (2018), the authors used multiple linear regression (MLR) techniques in a pedotransfer function to estimate the reactive fractions and the availability of Ca, Cr, Cu, Ni, Pb, and Zn in soils cultivated with vegetables in São Paulo, Brazil. Rodríguez-Lado & Martínez-Cortizas (2015) evaluated the performance of MLR, a random forest algorithm, and artificial neural networks for estimating the spatial variability in the soil bulk density measurements using organic carbon and soil particle size fractions.

This study hypothesizes that the S index is spatially dependent and could be estimated by pedotransfer functions using easily determined soil attributes. The objective of the study was to investigate the spatial variability of the S index and its correlation with other physical and chemical attributes of an Oxisol within a livestock farming system. The results can be used to determine which attributes are appropriate for generating a pedotransfer function to estimate the S index.

MATERIAL AND METHODS

The study was conducted at a farm in the municipality of Selvíria, Mato Grosso do Sul, Brazil, at the geographical coordinates 20° 36′ 58.2″ South, 51° 41′ 47.8″ West [WGS 84, EPSG: 4326], at 357 m altitude. The climate of the region, according to the Köppen classification, is defined as Aw (tropical humid), with a rainy summer and a dry winter and an annual average temperature of 25 °C, annual rainfall of 1330 mm and average relative humidity of 66%.

The soil of the site was classified as a Dystrophic Red Yellow Latosol, according to the classification of Embrapa (2018) or an Oxisol, according to the Soil Survey Staff (2014), with a clay content of 130 g kg⁻¹, a silt content of 60 g kg⁻¹ and a sand content of 810 g kg⁻¹ in the 0.00–0.30 m layer, determined by the pipette method. The soil chemical characteristics during the evaluation period were as follows: organic matter content of 17.0 g kg⁻¹, pH in CaCl₂ of 4.9, available P (resin) of 5.8 mg dm⁻³, exchangeable K (resin) of 0.6 mmol_c dm⁻³, exchangeable Ca (resin) of 9.4 mmol_c dm⁻³, exchangeable Mg (resin) of 10.5 mmol_c dm⁻³, Al exchangeable of 0.9 mmol_c dm⁻³, cation exchange capacity of 40.5 mmol_c dm⁻³ and base saturation of 50.1.

The study area was deforested in 1989, when trees of economic interest were removed, and the native vegetation was burned later. The grassland (Urochloa decumbens) was then sowed for livestock farming, with an average stocking rate of three animal unit's (an animal of 450 kg per animal unit) per hectare throughout the year. Adult male animals were used. In 2002, 2008, and 2014, the soil acidity was corrected, and it was fertilized using 1.5 Mg ha⁻¹ of dolomitic limestone and 350 kg ha⁻¹ of NPK 04-14-08 from urea, simple superphosphate, and potassium chloride sources.

A georeferenced sampling mesh with a dimension of 1.4 ha and 71 points was installed in the center of the area (Figure 1). At each point, samples with preserved and deformed structures were collected to determine the physical and chemical attributes of the soil. The preserved structure samples were extracted using 100 cm³ cylindrical rings and the deformed structure samples were extracting with an auger, both from the 0.00–0.15 m and 0.15–0.30 m soil layers.

FIGURE 1
Study site location and representation of the georeferenced sampling mesh, which contained 71 points [WGS 84, EPSG: 4326].

The preserved samples were saturated with water in plastic trays for 48 h and were then subjected to 0, 30, 60, and 100 hPa tension in a sandbox and 300, 600, 1000, 5000, and 15000 hPa tension in pressure chambers with a porous plate (Klute, 1986); these samples were later weighed to determine their water contents. Subsequently, the samples were oven dried at 105 °C for 48 h and reweighed. Macroporosity (MA) was determined to be the difference between the saturated and the 100 hPa soil sample mass; the total porosity (TP) was determined to be the difference between the saturated and dry soil sample mass. Furthermore, the microporosity (MI) was determined to be the difference between the TP and the MA. The bulk density (BD) was calculated as the ratio of the soil dry mass to the volume of the cylindrical rings, according to the methodology proposed by Embrapa (2017).

The SWRC at each sampling point were calculated using volumetric and gravimetric methods, according to the model proposed by Van Genutchen (1980), with Mualem restriction, where m = 1 - 1 / n, according to [eq. (1)], as follows:

(1) $$\text{θ}=\text{θr}+\frac{\left(\text{θs}+\text{θr}\right)}{{\left[{\left(1+\text{α}\Psi \right)}^{\mathrm{n}}\right]}^{\mathrm{m}}}$$

Where:

θ is the soil water content (m³ m⁻³);

θr is the residual water content (m³ m⁻³);

θs is the water content at saturation (m³ m⁻³);

Ψ is the applied tension (hPa), and

α, n, and m are the adjustment parameters of the equation.

The SWRCs of each sampling point are represented in Figure 2.

FIGURE 2
Soil water retention curve in volumetric base at the respective sampling points and (A) 0.00–0.15 m and (B) 0.15–0.30 m layers.

The soil physical quality S index proposed by Dexter (2004 a, b, c), which corresponds to the angle of inclination of the SWRC at the inflection point, was calculated from the adjustment parameters of the SWRC using gravimetric methods, as follows:

(2) $$\mathrm{S}=-\mathrm{n}\left(\text{θs}-\text{θr}\right){\left[1+\frac{1}{\mathrm{m}}\right]}^{-\left(1+\mathrm{m}\right)}$$

Where:

θr is the residual water content (kg kg⁻¹),

θs is the water content at saturation (kg kg⁻¹), and

n and m are the adjustment parameters.

The deformed structure samples were analyzed for their sand, silt, and clay contents using the pipette method, soil organic matter (SOM) using the oxidation dichromate K method, cation exchange capacity (CEC) using the acetate extraction method, and base saturation (V%) by measuring the % of K, Ca, and Mg ions in the CEC, according to Van Raij et al. (1996). The total organic carbon (TOC) and its fractions, the particulate organic carbon (OCp), and the organic carbon associated with the mineral fraction (OCa) were determined according to Cambardella & Elliot (1992). To determine the TOC content, 10 mL of 1N HCl was added to 10 g of soil to remove carbonates (inorganic carbon from limestone). Thereafter, the soil was oven dried at 60 °C for 15 h. The soil was macerated and then analyzed according to the methodology proposed by Van Raij et al. (1996). From the 10 g of sample used to determine the TOC, 5 g were used for particle size separation of the OCa and OCp. These fractions were determined following the same methodology as that used for TOC.

The stocks of TOC (TOC-S), OCp (OCp-S), and OCa (OCa-S) were calculated by multiplying the C (g kg⁻¹) by the soil bulk density (Mg m⁻³) and the depth of the evaluated layer (cm). To compare stocks between equal soil masses, stocks were corrected for the mass of the equivalent soil (Ellert & Bettany, 1995), using the soil bulk density over a reference area (native forest next to the experimental area) as a baseline.

Thereafter, a descriptive analysis was used to calculate the mean, minimum, maximum, and coefficient of variation (CV) values of the data set. The CV was classified according to Wilding & Drees (1983). The hypothesis of normality was verified using the Shapiro and Wilk probability distributions test (p> 0.01) and, when outliers were found, those with values 2.5 times greater than the interquartile range were removed to limit them to 5% of the total observations.

Principal component analysis (PCA) was performed to select attributes that were correlated with the S index in the principal component (PC). PCs with eigenvalues ≥ 1 and explained variance of at least 5% of the total variation of the dataset were examined. To convert the real values of the soil attributes into scores, a sigmoid equation proposed by Mukhopadhyay et al. (2014) was used, with one asymptote tending to 0 and the other tending to 1, as follows:

(3) $$\text{Score}=\frac{1}{\left[1+{\left(\frac{x}{x_0}\right)}^b\right]}$$

Where:

x is the value of the soil attribute (i.e., CEC, TOC-S, clay and sand values);

x₀ is the mean value of the soil attribute, and

b is the slope of the equation.

The slope was 2.5 for the attribute with maximum values as the best condition and −2.5 for the attribute with minimum values as the best condition (i.e., BD and θr).

Each PC explained a certain amount of variation (%) in the total variation of the dataset and this percentage was used to weight the correlated variables in a given PC and to calculate the soil quality index based on PCA (PCA-SQI) (Mukhopadhyay et al., 2014). This index made it possible to verify the representativeness of the S index within the total variability of the dataset.

(4) $$\text{PCA}-\text{SQI}={\displaystyle {\sum }_{\mathrm{i}=1}^{\mathrm{n}}{\mathrm{W}}_{\mathrm{i}}{\text{score}}_{\mathrm{i}}}$$

Where:

W_i is the weighting factor extracted from the PCA, and

score_i is the score of the standardized attribute.

The soil physical and chemical attributes correlated with the S index in the PCs were tested in the construction of a pedotransfer function using MLR techniques with a stepwise procedure (forward and backward) to select the best predictors (equation 5). Of the observations, 70% were used for MLR training and 30% for function validation. The training consisted of adjusting the parameters until the error between the output patterns generated by the function reached the desired minimum value, that is, close to zero.

(5) $$\widehat{\mathrm{S}}={\text{β}}_0+{\text{β}}_1{\mathrm{A}}_1+{\text{β}}_2{\mathrm{A}}_2+\cdots +{\text{β}}_{\mathrm{n}}{\mathrm{A}}_{\mathrm{n}}$$

Where:

Ŝ is the estimated S index;

β₀ is the linear adjustment coefficient;

β_n is the coefficient associated with the attribute, and

A is the selected attribute.

The spatial dependence of the S index was verified by a geostatistical analysis, according to the semivariance calculation, as follows:

(6) $$\text{γ}(\mathrm{h})=\frac{1}{2\mathrm{N}(\mathrm{h})}{\displaystyle {\sum }_{\mathrm{i}=1}^{\mathrm{N}(\mathrm{h})}{\left[\mathrm{Z}\left({\mathrm{X}}_{\mathrm{i}}\right)-\mathrm{Z}\left({\mathrm{X}}_{\mathrm{i}}+\mathrm{h}\right)\right]}^2}$$

Where:

γ (h) is the calculated semivariance;

N(h) is the pair of values measured by the distance h;

Z(X_i) is the value of the attribute Z at point X_i, and

Z (X_i + h) is the value of the attribute Z at point X_i + h. The distance vector (h) is given in m.

Spherical, exponential, and gaussian semivariogram models were tested, and the selection criteria were: a) the number of pairs of points in the first lag > 50, b) a smaller root mean square error (RMSE), c) a higher coefficient of determination (R²), d) a linear coefficient of the cross validation close to zero and an angular coefficient close to one. Using the semivariograms, the data interpolation was performed by the ordinary kriging method to construct spatial distribution maps with equidistant intervals.

Descriptive statistical analyses, the Shapiro and Wilk normality test, PCA, MLR stepwise methods and geostatistics were performed using R software version 3.4.3 (R Development Core Team, 2015) and the “agricolae”, “nortest”, “factoextra”, and “geoR” packages, respectively.

RESULTS AND DISCUSSION

The results of descriptive analysis indicated a higher clay content (148 g kg⁻¹) in the 0.15–0.30 m layer of the Oxisol (Table 1). Therefore, higher levels of sand and silt were observed in the superficial layer, but the variation between layers were of smaller magnitude than were the granulometric variations found within the same layer. In 0.00–0.15 m layer, a clay maximum of 69–175 g kg⁻¹ was observed, and a sand maximum of 773–872 g kg⁻¹ was observed. The CVs were classified as medium for clay and low for sand, in both of the evaluated layers. The texture can modify soil aggregation, porosity, bulk density, the volume of water retained, and its availability to plants (Bonetti et al., 2017). In this case, it is assumed that the texture may be a preponderant factor in the spatial variability of the other soil attributes analyzed, as suggested in the hypothesis of this study.

TABLE 1
Descriptive statistics and the Shapiro and Wilk normality test (S-W) for soil attributes evaluated in (L1) 0.00–0.15 m and (L2) 0.15–0.30 m layers.

Attributes1	Mean		Minimum		Maximum		CV2		S-W3
	L1	L2	L1	L2	L1	L2	L1	L2	L1	L2
S-index	0.086	0.095	0.025	0.031	0.153	0.174	38.7	30.3	0.039ns	0.026ns
α	0.028	0.031	0.013	0.014	0.107	0.062	46.8	30.7	0.000*	0.013ns
n	2.997	2.885	1.382	1.425	5.016	4.263	30.9	23.0	0.074ns	0.016ns
θs (m3 m−3)	0.363	0.369	0.313	0.273	0.442	0.469	6.4	10.8	0.016ns	0.745ns
θ30 (m3 m−3)	0.316	0.293	0.264	0.216	0.380	0.359	8.6	11.2	0.404ns	0.277ns
θ60 (m3 m−3)	0.241	0.204	0.173	0.137	0.311	0.303	16.5	19.8	0.002*	0.080ns
θ100 (m3 m−3)	0.197	0.162	0.155	0.109	0.278	0.261	16.8	22.4	0.000*	0.000*
θ300 (m3 m−3)	0.164	0.131	0.130	0.094	0.235	0.207	12.6	17.5	0.000*	0.000*
θ600 (m3 m−3)	0.157	0.125	0.123	0.092	0.215	0.185	9.7	13.6	0.000*	0.000*
θ1500 (m3 m−3)	0.154	0.122	0.120	0.092	0.196	0.164	7.3	10.5	0.000*	0.036ns
θr (m3 m−3)	0.151	0.119	0.112	0.092	0.169	0.150	5.3	8.3	0.000*	0.154ns
AW (m3 m−3)	0.047	0.042	0.005	0.007	0.123	0.136	67.5	76.1	0.001*	0.000*
MA (m3 m−3)	0.174	0.216	0.080	0.093	0.277	0.283	22.2	18.0	0.862ns	0.117ns
MI (m3 m−3)	0.186	0.153	0.067	0.043	0.262	0.239	17.2	19.6	0.017ns	0.000*
TP (m3 m−3)	0.359	0.368	0.147	0.228	0.442	0.432	9.6	8.3	0.000*	0.000*
Clay (g kg−1)	129.311	149.001	69.208	88.353	175.351	206.607	18.4	17.0	0.178ns	0.815ns
Silt (g kg−1)	46.528	41.515	17.435	14.243	95.131	78.987	40.4	37.5	0.003*	0.000*
Sand (g kg−1)	824.160	809.484	773.869	733.734	872.745	854.564	2.9	2.9	0.267ns	0.010ns
BD (Mg m−3)	1.599	1.613	1.258	1.268	1.795	2.011	5.1	5.7	0.005*	0.000*
SOM (g kg−1)	16.475	10.779	10.000	9.000	23.000	15.000	16.1	10.0	0.104ns	0.000*
TOC-S (Mg m−3)	14.908	12.484	9.049	10.424	20.812	17.373	16.1	10.0	0.104ns	0.000*
OCp-S (Mg m−3)	6.383	6.508	4.113	5.556	13.338	8.605	29.7	11.2	0.000*	0.001*
OCa-S (Mg m−3)	8.426	6.005	2.099	2.996	16.401	11.197	26.5	23.6	0.160ns	0.002*
CEC (mmolc dm−3)	40.122	28.791	31.500	24.500	51.400	39.200	12.4	10.9	0.048ns	0.000*
V (%)	50.075	45.693	27.000	29.000	75.000	62.000	18.4	14.9	0.576ns	0.072ns

• ¹
  θs: Soil water content at saturated point; θ30, θ60, θ100, θ300, θ600 and θ1500, respectively, correspond to the soil water content at 30, 60, 100, 300, 600, 1500 hPa tensions; θr: Residual soil water content; TOC-S: Total organic carbon stock; OCp-S: Particulate organic carbon stock; OCa-S: Organic carbon associated with mineral fraction stock; SOM: soil organic matter content; CEC: Cationic exchange capacity at pH 7.0; V: Base saturation.
• ²
  CV: Corresponds to the coefficient of variation.
• ³
  p-value of the S-W test,
• ^ns
  normal distribution and
• ^*
  non-normal distribution considering p < 0.01.

To characterize the CEC and V%, the criteria established by Sousa & Lobato (2004) for Brazilian Cerrado soils were used. The CEC was low (< 60 mmol_c dm⁻³) in both evaluated layers, and it was 28% lesser in the 0.15–0.30 m layer than that in the upper layer (Table 1). The coefficients of variation were low, indicating a low dispersion of the observed values around the mean. Therefore, the CEC variability in the same layer is low relative to that between layers. However, the V% had a large range of values within the same layer, between 27% and 75% (0.00–0.15 m), and between 29% and 62% (0.15–0.30 m). These results indicate that soil acidity management and fertilization practices should consider the spatial variability of CEC.

The TOC-S were classified as satisfactory for livestock farming system (> 12 Mg ha⁻¹), in both of the evaluated layers (Carvalho et al., 2010) (Table 1). It is assumed that the adopted management system (moderate stocking rate with three animal unit's ha year⁻¹) was favorable for carbon accumulation in the superficial soil layers. The high coefficients of variation and the maximum and minimum amplitude values obtained indicated the existence of sites with values greater than 20 Mg ha⁻¹ TOC-S. The most evolved physical fraction of TOC-S (OCa-S) corresponded to the greatest concentration of carbon found in the 0.00–0.15 m layer, representing 57% of the TOC-S, with an average value of 8.51 Mg ha⁻¹ (Table 1). In the second layer, this condition is altered and the particulate fraction (OCp-S) represents most of the TOC-S, 52%, with 6.37 Mg ha⁻¹. Therefore, it is necessary to conduct local investigations where such discrepant values of OCa-S and OCp-S occurred to verify if these conditions directly reflect the soil water retention and the index S values. The coefficients of variation of TOC S, OCp-S, and OCa-S were classified as medium in the 0.00–0.15 m layer (Table 1), whereas in the 0.15–0.30 m layer, TOC-S and OCp- S were classified as low, and OCa-S was classified as medium.

In contrast to organic carbon stocks, the soil bulk density (BD) did not vary between layers and had a mean of 1.60 kg dm⁻³ (Table 1). However, he highest values, between 1.04 and 2.01 kg dm⁻³, was verified in the second layer evaluated (0.15–0.30 m). Using the coefficients of variation values, both layers were classified as low, but greater variability in the BD values was observed in the second layer.

The soil water contents under applied tensions of 30, 60, 100, 300, 600, and 1500 hPa were higher in the 0.00–0.15 m layer (Table 1). This indicates that this layer is more exposed to soil compaction by animal traffic and has a higher organic matter content. These results are concomitant to the soil water content at the saturation point (θs), which was higher (0.369 m³ m⁻³) in the 0.15–0.30 m layer, indicating the soil pore structure was maintained at depth by root canals formed by the pasture over time (Moraes et al., 2014).

However, the residual soil water content (θr) was 0.151 m³ m⁻³ in the 0.00–0.15 m layer, and this water is considered to be unavailable to the plants (Table 1). Surface soil compaction, caused by animal traffic, may lead to the disassembly of aggregates into microaggregates. These microaggregates have micropores, also known as cryptopores, with diameters < 0.0001 mm, in which water can remain retained with energy > −15000 hPa.

The highest α values were found in the 0.15–0.30 m layer (0.031) and ranged from 0.015 to 0.066 (Table 1). In the second layer, a mean α value of 0.028 was observed. According to Van Genuchten & Nielsen (1985), this SWRC coefficient corresponds to the water capacity drainage and the air intake in the soil pores. This result confirms the maintenance of the soil pores structure in the 0.15–0.30 m layer due to the formation of biopores created by decomposed roots. Therefore, the 0.15–0.30 m layer has better conditions for gas diffusivity in the system relative to the surface soil layer.

The variation of the coefficient n between the 0.00–0.15 and the 0.15–0.30 m layers was small, between 2.989 and 2.852, respectively, as compared to the variations within the same layer, i.e., between 1.392 and 5.255 in the first layer and 1.424 and 4.290 in the second layer (Table 1). According to Magalhães et al. (2018), the coefficient n corresponds to the soil pore size distribution index and together with parameter α, affects the shape of the SWRC, as shown in Figure 2.

The S index had high coefficients of variation, > 35%, which explains the high amplitudes of the maximum and minimum values, between 0.025 and 0.174, in both of the evaluated layers (Table 1). In addition, S index values below 0.035 indicated “soil of low structural quality” at specific sites in the area.

The Shapiro & Wilk normality test indicated a normal distribution for the S index, n, θs, θ30, MA, MI, clay, sand, SOM, TOC-S, OCa-S, V% and CEC in the 0.00–0.15 m layer, and the S index, α, n, θs, θ30, θ60, θ1500, θr, MA, clay, and V% in the 0.15–0.30 m layer. The other attributes had non-normal distributions, and the test indicated that the mean and the observed variance were not statistically equal to the values estimated from the statistical model.

Seven PCs were separated in the PCA analysis, with eigenvalues between 8.978 and 1.086, which explained 87.2% of the total variability of the dataset (Table 2). The first two components (PC1 and PC2) retained more than 52% of this variance and correlated positively with the S index, n, MA, sand, and SOM and negatively with the θ30, θ60, θ100, θ300, θ600, θ1500, θr, available water (AW), MI, clay, TOC-S, OCa-S, and CEC. The correlations obtained for soil water retention, carbon stocks (TOC-S and OCa-S), and soil granulometry (clay and sand contents) may be important covariates for estimating the S index. According to Rawls et al. (2003), the content and composition of organic matter affect both the soil structure and its adsorption properties; therefore, water retention can be altered by changes in the organic matter that occur due to climate and management changes.

TABLE 2
Principal component analysis of the physical and chemical soil attributes in the 0.00–0.15 and 0.15–0.20 m layers.

PCA	PC 1	PC 2	PC 3	PC 4	PC 5	PC 6	PC 7
Eigenvalue	8,978	3,568	2,519	2,080	1,479	1,226	1,086
Explained variance (%)	37,407	14,868	10,497	8,667	6,164	5,110	4,526
Accumulated variance (%)	37,407	52,275	62,772	71,439	77,603	82,713	87,239
Attributes 1	Correlation coefficient 2
S-index	0,666 *	0,561 *	-0,134	-0,355	0,063	-0,055	-0,086
α	0,075	-0,090	-0,439	0,716 *	0,063	0,323	0,151
n	0,523 *	0,621 *	0,176	-0,394	-0,006	-0,008	-0,039
θs	0,016	0,091	-0,866 *	-0,018	0,145	-0,143	-0,243
θ30	-0,529 *	0,137	-0,230	-0,687 *	0,050	-0,265	-0,151
θ60	-0,844 *	-0,199	-0,037	-0,414	0,006	-0,158	-0,035
θ100	-0,930 *	-0,283	-0,091	-0,158	-0,017	0,001	0,023
θ300	-0,938 *	-0,074	-0,163	0,011	-0,018	0,196	0,091
θ600	-0,925 *	0,107	-0,174	0,005	-0,010	0,225	0,105
θ1500	-0,880 *	0,297	-0,177	-0,022	-0,002	0,220	0,101
θr	-0,770 *	0,506 *	-0,167	-0,075	0,004	0,171	0,061
AW	-0,696 *	-0,630 *	-0,016	-0,150	-0,024	-0,096	-0,007
MA	0,774 *	0,147	-0,511 *	0,162	-0,007	-0,162	-0,095
MI	-0,830 *	-0,143	-0,114	-0,082	-0,121	-0,155	0,273
TP	0,154	0,045	-0,807 *	0,130	-0,138	-0,382	0,164
Clay	0,003	-0,788 *	-0,024	0,067	0,293	-0,121	-0,090
Silt	-0,214	0,145	0,008	0,087	-0,577 *	0,171	-0,625 *
Sand	0,148	0,738 *	0,020	-0,132	0,094	0,009	0,537 *
BD	0,032	-0,191	0,345	0,252	-0,303	-0,569 *	0,253
SOM	-0,699 *	0,526 *	0,171	0,303	0,216	-0,190	-0,096
TOC-S	-0,594 *	0,390	0,184	0,434	0,285	-0,307	-0,186
OCp-S	-0,002	-0,081	0,097	0,016	0,847 *	-0,020	-0,161
OCa-S	-0,556 *	0,411	0,136	0,375	-0,253	-0,257	-0,094
CEC	-0,680 *	0,430	0,126	0,211	0,137	-0,072	-0,176

• ¹
  θs: Soil water content at saturated point; θ30, θ60, θ100, θ300, θ600 and θ1500, respectively, correspond to the soil water content at 30, 60, 100, 300, 600, 1500 hPa tensions; θr: Residual soil water content; TOC-S: Total organic carbon stock; OCp-S: Particulate organic carbon stock; OCa-S: Organic carbon associated with mineral fraction stock; SOM: soil organic matter content; CEC: Cationic exchange capacity at pH 7.0; V: Base saturation.
• ²
  CV: Corresponds to the coefficient of variation.
• ³
  p-value of the S-W test,
• ^ns
  normal distribution and
• ^*
  non-normal distribution considering p < 0.01. ^2*Correlations considered in the principal component's interpretation. Underlined values indicate the attributes chosen for SQI calculation.

The other principal components PC3, PC4, PC5, PC6 and PC7, represented the θs conditions and were positively influenced by MA and TP; the α conditions were negatively influenced by θ30; and the OCp-S accumulation was negatively influenced by the silt content (Table 2).

Using the PCA-SQI calculation, the S-index, soil water retention, α, TOC-S, and OCa-S were the attributes that retained the largest portion of the total variability of the dataset (Figure 3). However, the representativity of the S index was only 9.4% in the 0.00–0.15 m layer and 11.5% in the 0.15–0.30 m layer. In this case, the PCA-SQI indicator suggests that the joint evaluation of SWR, TOC-S, and OCa-S is a more appropriate approach for ascertaining soil quality in the 0.00–0.15 m layer. In the 0.15–0.30 m layer, the expression of the S index increased due to the lower representativity of water retention and carbon stocks, which confirms the results presented earlier in the descriptive analysis.

FIGURE 3
Soil attributes contribution to the total variability of the data set.

Using the PCA-SQI values, the soil quality is higher in the 0.00–0.15 m layer (0.45), which may be explained by the higher water retention and OCa-S and TOC-S values. In the equation used to calculate the PCA-SQI, these attributes were tied to the highest adjustment coefficient (0.382) (Figure 3). This explains the greater preponderance of the total value. Similarly, Rawls et al. (2003) reported the effects of the organic carbon content on soil water retention. The authors found that in sandy soils with a low carbon content, the initial carbon input led to a significant increase in the water retention. However, when this contribution occurred in clayey soils, the soil water retention decreased. If the organic carbon content was high, increasing it would result in an increase in the water retention, independent of the textural class (as observed in the present study). Several studies, such as those by Naderi-Boldaji & Keller (2016) and Van Lier (2014), emphasize that soil quality indices express the complexity of the system, and that the use of a single indicator (such as the S index) should be viewed with great caution along with other indicators.

The stepwise procedures suggested that the MA, TP, and sand contents are appropriate attributes for estimating the S index because the MLR was significant at 0.01 probability, with an R² of 0.57, a RMSE below 0.2 and an Akaike's information criterion (AIC) of −711.3 (Table 3). The 100 samples (70% of the dataset) used for MLR training were sufficient to fit model under the study conditions. The MLR adjustment coefficients were all significant at 5% probability, with low standard deviations. In addition, the correlation coefficient obtained (R = 0.75) indicated the strong relationship between the selected attributes and the S index.

TABLE 3
Statistical summary of the proposed multiple linear regression used to estimate the S index as a function of the correlated attributes in the PCs of the principal component analysis.

Summary1	Ŝ = a + b (MA) + c (MI) + d (Sand)
Number of samples for training	100
R	0,75
R2	0,57
RMSE	0,174
F	38,17
p-value	8,815e-16**
AIC	-711,3
Summary of coefficients 2	a	b	c	d
Coefficient value	-0,2876	0,2045	-0,2827	4,72e-04
Standard error	0,0806	0,0816	0,1126	0,0001
T	-3,5680	2,5050	-2,5100	5,6240
p-value	5,92e-04**	1,41e-02*	1,39e-02*	2,27e-07**

• ¹
  RMSE: Root mean squarer error; F: Fisher test calculated F-value; AIC: Akaike criteria information.
• ²
  t: Student test calculated t-value.
• ^{* and **}
  indicates significance at 0.05 and 0.01 probability, respectively.

The estimates of the S index in the validation dataset (42 samples) reaffirm the quality of the proposed function, which had an R² of 0.61 and an R of 0.78 between the estimated and observed values, which were close to a 1:1 line (Figure 4). The model was significant at 0.01 probability. The limitations of the pedotransfer function were exposed in the estimates of the minimum values (< 0.065); the function tended to overestimate the minimum values and underestimate the maximum values (> 0.120). However, for S values between 0.065 and 0.120, the model estimates are consistent with the observed values.

FIGURE 4
Estimated S-Index values by MLR as a function of observed S-index values of the validation dataset (42 samples).

The semivariance calculated for the S index in the 0.00–0.15 m and 0.15–0.30 m layers were adjusted to the theoretical models of spherical semivariograms, respectively (Figure 5A, C). The quality of the semivariogram adjustments was confirmed by the RMSE values (< 2.84e-4) and the moderate correlation (R > 0.30; R < 0.50) in the cross-validation between the observed and estimated values. The degree spatial dependence (DSD) in the 0.00–0.15 m layer was 51.1% and was 50.1% in the layer 0.15–0.30 m, indicating a moderate contribution (C) of the semivariance in the establishment of the threshold (C₀ + C) and a high nugget effect (C₀) of the semivariogram (Cambardella et al., 1994). The semivariance was better adjusted to the spherical and gaussian models for the S index estimated by the pedotransfer function, in the 0.00–0.15 m and 0.15–0.30 m layers, respectively (Figure 5B, D). The R² was above 0.80, the RMSE close to zero, and the DSD value was 65.3 e 57%, respectively, indicating lower dispersion of values relative to the semivariogram estimates and a low variation in the S values at the sampling points, relative to the adjustments obtained for the S index determined by the SWRC. The correlations obtained between the observed and estimated values by the semivariograms in the cross-validation were satisfactory, with R > 0.50 and the linear and angular coefficients close to zero and one, respectively.

FIGURE 5
Experimental and fitted theorical semivariograms for the S-index observed values in the 0.00–0.15 m (A) and 0.15–0.30 m (C) layers, and the S-index values estimated by MLR in the 0.00–0.15 m (B) and 0.15–0.30 m (D) layers. Model adjustment parameters: spherical (Sph) and gaussian (Gau) model; nugget (C₀); scale (C); range (A₀); coefficient of determination (R²); root mean square error (RMSE). Cross-validation parameters: linear coefficient (y₀); angular coefficient (b) and correlation coefficient (R).

Based on the kriging maps, it was possible to observe the effects of the semivariograms adjustments on the reach and the establishment of the threshold (Figure 6). The S index maps of the 0.00–0.15 m layer presented a lower spatial continuity of the class's values, characteristic of spherical semivariogram models (Figures 6A, B). This is consistent with what has been discussed previously; the animal component attributed spatial variability to the soil attributes in surface layer, which occurs most effectively above 0.05 m (Moraes et al., 2014). Under a high stocking rate, which is common in an extensive livestock system, the distribution of animal weight occurs over a small area, delimited by its hull, and generates sites with less organic material and high BD values. As previously discussed, related to the quality of the pedotransfer function, it was evident in the estimated S index map that the regions delimited by the minimum (0.070 to 0.075) and maximum (0.100 to 0.105) values are smaller than the observed S index map, exposing the limitations of MLR. However, the S maps for the 0.15–0.30 m layer show that delimitations of the classes of the minimum values between 0.080 and 0.085 contrasted in several places in the right quadrant of the maps (Figure 4A). However, the values represented by the mean and maximum values were satisfactorily estimated in this layer.

FIGURE 6
Observed and estimated S index kriging maps in the (A) 0.00–0.15 and (B) 0.15–0.30 m layers.

The lack of similarity between the S index and the other indicators used in soil quality classification, such as the PCA-SQI, shows that the proposed critical limits for S, such as 0.035, are generally not valid and do not apply to all conditions. To overcome this problem, it is necessary to know the spatial variability of the area, both in terms of the crop productivity and the soil attributes. In addition, the use of the S index should be considered to be part of a minimum set of soil quality assessment indicators, such as carbon stocks, particle size distribution, bulk density, etc., which are much more easily determined and consistent. Research efforts should focus on assessing soil quality for assessing land degradation from a more complex point of view or using a more integrated approach. In general, the estimates of the S index for the pedotransfer function and the spatial distribution maps obtained are similar to the observed S index values and maps, and the function can be applied in spatial variability studies.

CONCLUSIONS

The S index presented spatial dependence, and its range varied between 84 and 188 m in the spherical semivariogram models.

The S index can be applied for spatial variability studies of soil attributes and for precision agriculture to assist decisions regarding physical soil management.

Among the attributes that were most easily determined, only macroporosity, microporosity, and sand contents were appropriate for constructing a pedotransfer function.

The S index maps estimated by the pedotransfer function resemble the observed S values; thereby proving the applicability of the function in spatial variability studies.

• Area Editor: Fernando França da Cunha

REFERENCES

• Assis Júnior RN, Mota JCA, Freire AG, Alencar TLD (2016) Pore network of an Inceptisol under different uses and relativized S-index as an indicator of soil physical quality. Pesquisa Agropecuária Brasileira 51(9):1575-1583. DOI: http://dx.doi.org/10.1590/s0100-204x2016000900057
  » http://dx.doi.org/10.1590/s0100-204x2016000900057
• Boim AG, Rodrigues SM, Santos-Araújo SN dos, Pereira E, Alleoni LR (2018) Pedotransfer functions of potentially toxic elements in tropical soils cultivated with vegetable crops. Environmental Science and Pollution Research 25(13):1-11. DOI: http://doi.org/10.1007/s11356-018-1348-0
  » http://doi.org/10.1007/s11356-018-1348-0
• Bonetti AJ, Anghinoni I, Moraes MT de, Fink JR (2017) Resilience of soils with different texture, mineralogy and organic matter under long-term conservation systems. Soil and Tillage Research 174(12):104-112. DOI: https://doi.org/10.1016/j.still.2017.06.008
  » https://doi.org/10.1016/j.still.2017.06.008
• Cambardella CA, Moorman TB, Parkin TB, Karlen DL, Novak JM, Turco RF, Konopka AE (1994) Field-scale variability of soil properties in central Iowa soils. Soil science society of America journal 58(5):1501-1511.
• Cambardella CA, Elliott ET (1992) Particulate soil organic-matter changes across a grassland cultivation sequence. Soil science society of America journal 56(3): 777-783. DOI: https://doi.org/10.2136/sssaj1992.03615995005600030017x
  » https://doi.org/10.2136/sssaj1992.03615995005600030017x
• Carvalho JLN, Raucci GS, Cerri CEP, Bernoux M, Feigl BJ, Wruck FJ, Cerri CC (2010) Impact of pasture, agriculture and crop-livestock systems on soil C stocks in Brazil. Soil and Tillage Research 110(1):175-186. DOI: https://doi.org/10.1016/j.still.2010.07.011
  » https://doi.org/10.1016/j.still.2010.07.011
• Dalchiavon FC, Montanari R, Andreotti M (2017) Production and quality of Urochloa decumbens (stapf) r.d. webster forage co-related to the physical and chemical properties of the soil. Revista Ceres 64(3):315-326. DOI: http://dx.doi.org/10.1590/0034-737x201764030013
  » http://dx.doi.org/10.1590/0034-737x201764030013
• Dexter AR (2004a) Soil physical quality: Part I. Theory, effects of soil texture, density, and organic matter, and effects on root growth. Geoderma 120(3-4):201-214. DOI: https://doi.org/10.1016/j.geoderma.2003.09.004
  » https://doi.org/10.1016/j.geoderma.2003.09.004
• Dexter AR (2004b) Soil physical quality: Part II. Friability, tillage, tilth and hard-setting. Geoderma 120(3-4):215-225. DOI: https://doi.org/10.1016/j.geoderma.2003.09.005
  » https://doi.org/10.1016/j.geoderma.2003.09.005
• Dexter AR (2004c) Soil physical quality: Part III: Unsaturated hydraulic conductivity and general conclusions about S-theory. Geoderma 120(3-4):227-239. DOI: https://doi.org/10.1016/j.geoderma.2003.09.006
  » https://doi.org/10.1016/j.geoderma.2003.09.006
• Ellert BH, Bettany JR (1995) Calculation of organic matter and nutrients stored in soils under contrasting management regimes. Canadian Journal of Soil Science 75(4):529-538. DOI: https://doi.org/10.4141/cjss95-075
  » https://doi.org/10.4141/cjss95-075
• Embrapa - Empresa Brasileira de Pesquisa Agropecuária (2017) Manual de métodos de análise de solo. Rio de Janeiro, Centro Nacional de Pesquisa de Solos, 575p.
• Embrapa - Empresa Brasileira de Pesquisa Agropecuária (2018) Brazilian Soil Classification System. Brasília, Embrapa Solos, 5 ed. 531p.
• Klute A (1986) Water retention: laboratory methods. In: Klute A. Methods of soil analysis. Madison, American Society of Agronomy, p635-662.
• Magalhães WDA, Freddi ODS, Wruck FJ, Petter FA, Tavanti RF (2018) Soil water retention curve and s index as soil physical quality indicators for integrated production systems. Engenharia Agrícola 38(1):64-73. DOI: http://dx.doi.org/10.1590/1809-4430-eng.agric.v38n1p64-73/2018
  » http://dx.doi.org/10.1590/1809-4430-eng.agric.v38n1p64-73/2018
• Moraes A, Carvalho PCF, Anghinoni I, Lustosa SBC, Andrade SEVG, Kunrath TR (2014) Integrated crop–livestock systems in the Brazilian subtropics. European Journal of Agronomy 57(7):4-9. DOI: https://doi.org/10.1016/j.eja.2013.10.004
  » https://doi.org/10.1016/j.eja.2013.10.004
• Mukhopadhyay S, Maiti SK, Masto RE (2014) Development of mine soil quality index (MSQI) for evaluation of reclamation success: a chronosequence study. Ecological Engineering 71(10):10-20. DOI: https://doi.org/10.1016/j.ecoleng.2014.07.001
  » https://doi.org/10.1016/j.ecoleng.2014.07.001
• Naderi-Boldaji M, Keller T (2016) Degree of soil compactness is highly correlated with the soil physical quality index S. Soil and Tillage Research 159(6):41-46. DOI: https://doi.org/10.1016/j.still.2016.01.010
  » https://doi.org/10.1016/j.still.2016.01.010
• R Development Core Team (2015) A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0. URL: http://www.R-project.org
  » http://www.R-project.org
• Rawls WJ, Pachepsky YA, Ritchie JC, Sobecki TM, Bloodworth H (2003) Effect of soil organic carbon on soil water retention. Geoderma 116(1-2):61-76. DOI: https://doi.org/10.1016/S0016-7061(03)00094-6
  » https://doi.org/10.1016/S0016-7061(03)00094-6
• Reichert JM, Rosa VT da, Vogelmann ES, Rosa DP da, Horn R, Reinert DJ, Sattler A, Denardin JE (2016) Conceptual framework for capacity and intensity physical soil properties affected by short and long-term (14 years) continuous no-tillage and controlled traffic. Soil and Tillage Research 158(4):123-136. DOI: https://doi.org/10.1016/j.still.2015.11.010
  » https://doi.org/10.1016/j.still.2015.11.010
• Rodríguez-Lado L, Martínez-Cortizas A (2015) Modelling and mapping organic carbon content of topsoils in an Atlantic area of southwestern Europe (Galicia, NW-Spain). Geoderma 245(4):65-73. DOI: https://doi.org/10.1016/j.geoderma.2015.01.015
  » https://doi.org/10.1016/j.geoderma.2015.01.015
• Rossetti KDV, Centurion JF (2018) Use of S-index as a structural quality indicator for compacted Latosols cultivated with maize. Revista Caatinga 31(2):455-465. DOI: http://dx.doi.org/10.1590/1983-21252018v31n222rc
  » http://dx.doi.org/10.1590/1983-21252018v31n222rc
• Soil Survey Staff (2014) Keys to soil taxonomy. Washington, United States Department of Agriculture, 12 ed. 360p.
• Sousa DMG, Lobato E (2004) Cerrado: correção do solo e adubação. 2 ed. Brasília, Embrapa Informações Tecnológicas, 2 ed. 416p.
• Van Genuchten MT (1980) A closed-form equation for predicting the hydraulic conductivity of unsaturated soils 1. Soil science society of America journal 44(5):892-898. DOI: https://doi.org/10.2136/sssaj1980.03615995004400050002x
  » https://doi.org/10.2136/sssaj1980.03615995004400050002x
• Van Genuchten MT, Nielsen DR (1985) On describing and predicting the hydraulic properties. Annales Geophysicae 3:615-628.
• Van Lier J (2014) Revisiting the S-index for soil physical quality and its use in Brazil. Revista Brasileira de Ciência do Solo 38(1):1-10. DOI: http://dx.doi.org/10.1590/S0100-06832014000100001
  » http://dx.doi.org/10.1590/S0100-06832014000100001
• Van Raij B, Cantarela H, Quaggio JA, Furlani AMC (1996) Recomendação de adubação e calagem para o estado de São Paulo. Campinas, Instituto Agronômico e Fundação IAC, 285p.
• Wilding LP, Drees LR (1983) Spatial variability and pedology. In: Wilding LP, Smeck NE, Hall GF. Pedogenesis and Soil Taxonomy: 1. Concepts and Interactions. Amsterdam, Elsevier, p83-116.
• Xu C, Xu X, Liu M, Yang J, Zhang Y, Li Z (2017) Developing pedotransfer functions to estimate the S-index for indicating soil quality. Ecological indicators 83(12): 338-345. DOI: http://dx.doi.org/10.1016/j.ecolind.2017.08.01
  » http://dx.doi.org/10.1016/j.ecolind.2017.08.01

Publication Dates

• Publication in this collection
  17 Feb 2020
• Date of issue
  Jan-Feb 2020

History

• Received
  08 Apr 2019
• Accepted
  18 Nov 2019