Spatial Analysis of Chemical and Textural Soil Attributes in a Multistrata Agroforestry System

Abstract

The objective was to evaluate the spatial distribution of chemical and textural soil variables in a multistrata agroforestry system. A total of 73 georeferenced soil samples were collected at depths of 10-20 cm and 20-40 cm. The studied parameters were: pH _H2O , potential acidity (H+Al), calcium (Ca ²⁺ ), magnesium (Mg ²⁺ ), aluminum (Al ³⁺ ), sodium (Na ⁺ ), potassium (K ⁺ ), phosphorus (P), organic carbon (C _org ), cation exchange capacity (T-value), base saturation (V-value), total clay, total sand, and silt. Principal component analysis (PCA) was performed in R software using the FactoMineR and Factoextra packages. For variables with spatial dependence, ordinary kriging was performed using the best-fitted model. For variables without spatial dependence, inverse distance weighted (IDW) interpolation was applied (power = 2). The spherical model was the best fit for chemical attributes. IDW interpolation accurately mapped the textural attributes. It was concluded that geostatistics enabled a detailed analysis of chemical and textural attributes.

Keywords:
Spatial models; kriging; soil conservation; sustainability; agroecology

1. INTRODUCTION

The pursuit of food production systems capable of ensuring people’s food security and nutrition is acknowledged as fundamental strategy to fight poverty. This goal is part of the Sustainable Development Goals (SDGs) agenda followed by 193 countries (Abhiyan, ( 2017 )). Accordingly, agroforestry systems (AFS) stand out for their high stability and flexibility, which can be achieved by combining tree species to shrubs, herbaceous plants, and/or animals distributed in space and time. This approach enables diversifying food production over the year (Steenbock & Vezzani, ( 2013 )).

AFSs have the advantage of promoting significant biomass production, mainly when they are implemented in tropical and subtropical areas featured by highly weathered and leached soils (Higuchi et al., ( 1998 ); Santos et al., ( 2004 )). This ability arises from an integrated action between plants and organisms that mediate cyclic organic matter deposition, decomposition, mineralization and humification processes taking place in the soil (Andrade et al., ( 1999 ); Vale Júnior et al., ( 2011 ); Pezarico et al., ( 2013 )).

According to Uzêda et al. (( 2017 )), AFSs also provide soil protection through canopy cover and litterfall, since they reduce direct soil exposure to solar radiation and rainfall impact on it. Moreover, different root systems found in AFSs help improving soil quality (Canuto, ( 2017 )).

Systematic organic matter input through management practices implemented in AFSs has positive influence on integrated and complex chemical, physical and biological processes (Araújo & Melo, ( 2010 )). It happens because this input stimulates the soil microbial community to decompose the produced residues, and this decomposition plays fundamental role in nutrient cycling, humic substances’ synthesis, as well as in compounds’ aggregation and degradation processes (Iwata et al., ( 2012 ); SEN et al., ( 2020 )).

Therefore, it is essential to both monitor and gather qualitative and quantitative information at the time to plan, implement and monitor AFSs (Araújo et al., ( 2019 )). Geostatistics is an efficient tool used for such a purpose, mainly in spatial analysis applied to soil physical and chemical attributes (Sibaldelli et al., ( 2015 )). Precise mapping of soil variables capable of influencing AFS components’ growth and production processes stands out among several advantages associated with the application of geostatistical methods in AFS management procedures. For this reason, it can help decision-making focused on specific and rational soil management in these systems (Silva et al., ( 2016 )).

The hypothesis of this study is that the chemical and textural soil attributes present spatial continuity structure in agroforestry systems and enable spatial estimates aimed at optimizing resource management. Thus, the aim of the current study was to assess the spatial distribution of soil chemical and textural variables in a multistrata agroforestry system.

2. MATERIALS AND METHODS

2.1. Study site and data collection

The investigated multistrata agroforestry system started being implemented in December 2014, within a 2,000-m² area, at the experimental field of the Agronomy Institute of Federal Rural University of Rio de Janeiro (UFRRJ), Seropédica County-RJ ( Figure 1 ). Terrain in the study site ranges from gently undulating to undulating, and it presents slopes ranging from 3% to 8%. Climate in this region is classified as Aw, humid tropical, which is featured by dry winter and rainy summer (Alvares et al., ( 2014 )).

Figure 1.
Study site located at Federal Rural University of Rio de Janeiro, Seropédica County, RJ.

Plant species’ selection for the agroforestry system was based on information about their biological function, adaptation to soil and weather conditions, tolerance to relevant diseases observed in the investigated region, as well as about diversification of products used for subsistence and income generation purposes in agroecological family farming systems. Soil preparation encompassed plowing and harrowing, as well as dolomitic limestone incorporation sixty days before planting the cassava stakes (Saracura variety) in alternation to pigeon pea ( Cajanus cajan cv. Fava Larga), in rows.

Persistent and perennial plant species were introduced in the system after cassava harvest. These species comprised banana ( Musa paradisiaca cv. BRS Princesa), coffee ( Coffea canephora var. Robusta) and peach palm ( Bactris gasipaes ) planted in alternating rows. They were interspersed with dense rows planted with gliricidia ( Gliricidia sepium ), flemingia ( Flemingia macrophilla ), vinhático ( Phathymenia reticulata ), Brazilian peppertree ( Schinus terebinthifolius ) and annatto ( Bixa orellana ). Guapuruvu ( Schizolobium parahyba ) and embaúba ( Cecropia angustifolia ) were planted between coffee plants, and pineapple (cv. Pérola) was planted between banana plants. Annual crops, such as sweet potato ( Ipomoea batatas ), peanut ( Arachis hypogae ), arrowroot ( Marantha arundinacea ), taro ( Colocasia esculenta ) and cassava (var. Saracura), were grown in the alleys between perennial species. The soil in the study site was classified as belonging to the planosol group, and its texture class was featured as sandy loam 40-cm down into the ground.

A total of 73 georeferenced soil samples were collected at depths of 10-20 cm and 20-40 cm. The studied parameters were: pH _H2O , potential acidity (H+Al), calcium (Ca ²⁺ ), magnesium (Mg ²⁺ ), aluminum (Al ³⁺ ), sodium (Na ⁺ ), potassium (K ⁺ ), phosphorus (P), organic carbon (C _org ), cation exchange capacity (T-value), base saturation (V-value), total clay, total sand, and silt. Descriptive analysis was applied to each variable in order to find their mean values, standard deviations and coefficients of variation. Moreover, data normality was assessed through Shapiro-Wilk test, at 95% confidence level. Furthermore, Principal Component Analysis (PCA) was carried out based on using FactoMineR (Lê et al., ( 2008 )) and Factoextra (Kassambara & Mundt, ( 2017 )) statistical packages in R software (R Core Team, ( 2015 )) to assess both the association and similarity of the analyzed variables.

2.2. Spatial analysis

Spatial continuity was assessed in the first geostatistical analysis phase. It was done by using the experimental semivariogram generated by the estimator applied to the semivariance ( Equation 1 ) of the soil chemical and textural variables. Three theoretical models were fitted through the Maximum Likelihood Method, in association with the GeoR package (Ribeiro Júnior & Diggle, ( 2001 )), in R software (R Core Team, ( 2015 )). The exponential ( Equation 2 ), spherical ( Equation 3 ) and Gaussian ( Equation 4 ) models were the adopted theoretical models.

\begin{equation}\widehat{\gamma}(h) = \frac{1}{2N(h)}\sum_{i = 1}^{N(h)}\left\lbrack Z\left( x_{i} \right) - Z(x_{i} + h) \right\rbrack^{2}\qquad{(1)}\end{equation}

\begin{equation}\widehat{\gamma}(h) = C_{0} + C\left\lbrack 1 - e^{\left( - \frac{h}{a} \right)} \right\rbrack\qquad{(2)}\end{equation}

\begin{equation}\widehat{\gamma}(h) = C_{0} + C\left\lbrack 1.5\frac{h}{a} - 0.5\left( \frac{h}{a} \right)^{3} \right\rbrack if:\ h\ <\ a\qquad{(3)}\end{equation}

\begin{equation}\widehat{\gamma}(h) = C_{0} + C\ if:\ h\ \geq\ a\end{equation}

\begin{equation}\widehat{\gamma}(h) = C_{0} + C\left\lbrack 1 - e^{\left( - \frac{h}{a} \right)^{2}} \right\rbrack\qquad{(4)}\end{equation}

Wherein, \begin{equation}\widehat{\gamma}(h)\end{equation}: semivariance; \begin{equation}Z\left( x_{i} \right)\end{equation} regionalized variable value at point \begin{equation}x\end{equation} ; \begin{equation}Z\left( x_{i} + h \right)\end{equation} variable value at point \begin{equation}x + h\end{equation} ; \begin{equation}N(h)\end{equation}: number of pairs separated by distance \begin{equation}h\end{equation} ; \begin{equation}C_{0}\end{equation}: nugget effect; \begin{equation}C\end{equation}: contribution; \begin{equation}a\end{equation}: range.

Criteria set to assess the quality of fits comprised Akaike Information Criterion (AIC) ( Equation 5 ) and Spatial Dependence Index (SDI) ( Equation 6 ). According to Cambardella et al. (( 1994 )), SDI can be classified as weak, when values are lower than 0.25; moderate, when values range from 0.25 to 0.75; and strong, when values are higher than 0.75. Other criteria defined after data cross-validation were also used, namely: Reduced Mean Error ( Equation 7 ) and Standard Deviation of Reduced Errors ( Equation 8 ). The main criterion adopted for model selection purposes was the lowest AIC value; models presenting difference in AIC values lower than 2 were statistically similar (Burnham & Anderson, ( 2002 )). However, in this case, selection criteria would prioritize the lowest values observed for absolute error (AR), reduced mean error ( \begin{equation}\overline{ER}\end{equation} ) and standard deviation of reduced errors (S _er ).

\begin{equation}AIC = - 2\ ln(L)\ + 2K\qquad{(5)}\end{equation}

\begin{equation}SDI = \frac{C}{(C_{0} + C)}\qquad{(6)}\end{equation}

\begin{equation}\overline{ER} = \frac{1}{n}\ \ \sum_{i = 1}^{n}{\ \frac{z\left( x_{i0} \right) - \widehat{z}\left( x_{i0} \right)}{\sigma\left( x_{i0} \right)}}\qquad{(7)}\end{equation}

\begin{equation}S_{er} = \sqrt{\frac{1}{n}\ \ \sum_{i = 1}^{n}{\left\{ \frac{z\left( x_{i0} \right) - \widehat{z}\left( x_{i0} \right)}{\sigma\left( x_{i0} \right)} \right\}^{2}\ }}\qquad{(8)}\end{equation}

Wherein, \begin{equation}L\end{equation}: maximum likelihood function of the model; K: number of parameters in the model; \begin{equation}C:\end{equation} contribution; \begin{equation}C_{0}:\end{equation} nugget effect; n: number of observations; \begin{equation}z\left( x_{i0} \right)\end{equation}: observed value; \begin{equation}\widehat{z}\left( x_{i0} \right):\end{equation} estimated value; \begin{equation}\sigma\left( x_{i0} \right):\end{equation} standard deviation of kriging; \begin{equation}i:\end{equation} sample point ranging from 1 to n.

Interpolation was applied to variables of interest presenting spatial dependence. It was done by using ordinary kriging ( Equation 9 ) in the best selected model. This technique was applied to generate variables-of-interest estimates at non-sampled points. Variables that did not show spatial dependence were subjected to inverse distance weighted interpolation (IDW) at power = 2 ( Equation 10 ).

\begin{equation}\widehat{Z}\left( x_{0} \right) = \sum_{i = 1}^{n}{\lambda_{i}Z(x_{i})}\qquad{(9)}\end{equation}

\begin{equation}\widehat{Z}\left( x_{0} \right)\frac{\sum_{i = 1}^{n}\left( \frac{1}{d_{i}^{2}}Z(x_{i}) \right)}{\sum_{i = 1}^{n}\left( \frac{1}{d_{i}^{2}} \right)}\qquad{(10)}\end{equation}

Wherein, \begin{equation}\widehat{Z}\left( x_{0} \right)\end{equation}: estimate at point \begin{equation}{\ x}_{0}\end{equation} ; \begin{equation}Z\left( x_{i} \right):\end{equation} observed value at point \begin{equation}x_{i}\end{equation} ; n: number of sample points; \begin{equation}\lambda_{i}:\end{equation} kriging-point weight; h: distance between points; p: power.

3. RESULTS AND DISCUSSION

Mean values recorded for all assessed chemical attributes have decreased from the 0-20 cm layer to the 20-40 cm layer ( Table 1 ). With respect to textural attributes, clay and silt recorded the highest means in the 20-40 cm soil layer; sand content, in its turn, recorded lower values. Another noteworthy factor lies on the reduced variability observed for the analyzed chemical attributes in the deepest soil layer, except for organic carbon, which recorded the highest heterogeneity level in the 20-40 cm layer. In addition, sand, clay and silt contents recorded increased variability in the deepest soil layer. Total cation exchange capacity (T) in the 0-20 cm layer was the only variable unable to show normal distribution in Shapiro-Wilk test, at 95% confidence level.

PCA application enabled obtaining two principal components for each of the assessed soil layers, based on the total variance recorded for the investigated chemical and textural attributes ( Figure 2 ). Accordingly, 31.4% (Dim1) and 21% (Dim2) were observed for the 0-20 cm soil layer, whereas 30.6% (Dim1) and 16.4% (Dim2) were orecorded for the 20-40 cm layer. Textural attributes presented low contribution to principal components at the 0-20 cm layer, whereas chemical attributes presented a more homogeneous contribution in the 20-40 cm layer.

PCA results ( Figure 2 ) pointed out high correlation among variables, such as organic carbon, pH in water, and sum of bases in both soil layers. According to Hongyu et al. (( 2015 )), the slope between vectors has shown correlation among variables. In other words, sharper angles meant stronger correlation, whereas angles closer to 90 degrees indicated lack of correlation. These findings were supported by maps generated based on kriging ( Figure 3 ), according to which, the highest values observed for these variables were distributed in similar areas. Seyedmohammadi & Matinfar (( 2018 )) observed that soil CEC was highly correlated to Dim1 (r = 0.81, P < 0.01) deriving from soil physical and chemical properties. The aforementioned authors used Dim1 as auxiliary variable to predict soil CEC based on the cokriging method. According to them, cokriging application provided more accurate results than kriging. This study has evidenced that Dim1 (deriving from soil physical and chemical properties) application as auxiliary variable to predict soil CEC through cokriging was effective in improving prediction accuracy. They concluded that the principal component presenting high correlation to the primary variable can be used as auxiliary variable to increase prediction accuracy.

Table 1.
Descriptive statistics of soil chemical and textural attributes in the agroforestry system (0-20 cm and 20-40 cm soil layers).

Variable	Unit	0-20 cm				20-40 cm
		x	sd	CV%	p	x	sd	CV%	p
H+Al	cmolc dm - ³	4.60	2.56	55.71	3.43E-04	4.52	1.48	32.80	1.20E-02
Al	cmolc dm - ³	0.33	1.39	421.84	6.43E-18	0.22	0.14	65.96	4.30E-06
S	cmolc dm - ³	4.24	2.33	55.05	9.99E-09	3.20	1.39	43.44	1.04E-09
T	cmolc dm - ³	8.84	3.00	33.98	1.15E-01	7.72	1.97	25.56	1.58E-01*
V	%	50.10	21.66	43.24	4.98E-03	41.55	11.60	27.93	2.70E-02
pH H20	1:2.5	5.19	0.61	11.69	9.98E-04	4.84	0.50	10.42	5.98E-07
Corg	%	1.23	0.67	54.50	1.50E-07	0.94	0.61	64.78	1.10E-12
Sand	%	74.03	4.09	5.52	4.15E-03	72.22	5.59	7.74	9.01E-09
Clay	%	11.03	4.34	39.38	2.89E-04	11.27	5.23	46.43	7.58E-03
Silt	%	15.01	5.49	36.57	8.20E-02	16.59	7.38	44.51	7.16E-08

x: mean; sd: standard deviation; CV%: coefficient of variation (expressed in percentage); p: p-value in Shapiro-Wilk normality test; H+Al: hydrogen + aluminum; Al: aluminum; S: sum of bases; T: total cation exchange capacity; V: base saturation indices; pH _H20: potential of hydrogen in water; Corg: organic carbon; Sand: sand content; Clay: clay content; Silt: silt content.

Figure 2.
Principal Component Analysis applied to soil chemical and physical attributes in the 0-20 cm and 20-40 cm soil layers, in agroforestry system.

Based on the assessment parameters appplied to all geostatistical models, most chemical features have shown strong or moderate spatial dependence and SDI values higher than 25% ( Table 2 ). Therefore, one can infer that soil chemical features are spatially structured and that spatial estimates can be performed based on ordinary kriging techniques. All assessed attributes presented strong spatial dependence (SDI > 75%) at the 0-20 cm soil layer. Aluminum, sum of bases, and pH showed strong spatial dependence (SDI > 75%) at the 20-40 cm layer, whereas the other variables recorded moderate spatial dependence (SDI ranging from 25% to 75%). Based on the study conducted by Silva et al. (( 2016 )) in an agroforestry system, in Seropédica County, Rio de Janeiro State, Mg ²⁺ and base saturation were the only chemical attributes that did not show spatial continuity structure, and it has evidenced the potential of geostatistical techniques to be used in soil features’ spatial estimates. Pure nugget effect was observed in the variographic analysis conducted during the aforementioned study.

The best fitting parameters were found based on using the spherical model for most variables, mostly in the 0-20 cm soil later. This model was followed by the exponential one and, in few cases, by the Gaussian model, in the 20-40 cm layer. Other studies focused on investigating geostatistical analysis application in agroforestry systems have also evidenced that these models are the best to model semivariance in the assessed attributes (Silva et al., ( 2016 )). Semivariogram models selected for each variable, spatial dependence index and the model fitting statistics are shown in Tables 2 and 3.

If one takes into consideration the models selected for each variable, the overall mean ranges from 28.65 m, in the 0-20 cm soil layer to 15.56 m, in the 20-40 cm layer. This finding pointed out the mean distance needed for a given variable to present spatial correlation. Thus, this parameter is essential to help determining spacing between sampling points for soil variables in agroforestry systems. Aluminum recorded the lowest range (4.21 m) in the 0-20 cm soil layer and the second lowest range (5.12 m) in the 20-40 cm layer. This finding suggested that aluminum concentrations in the soil were similar within a 5m radius.

Variable H+Al presented the largest range in both assessed layers (47.29 m, in the 0-20 cm layer; and 38.53 m, in the 20-40 cm layer). Therefore, sampling points set approximately 40 meters apart from each other were enough to detect spatial variations in the concentrations of this specific attribute.

It was not possible applying geostatistical analysis to textural attributes like sand, clay and silt either because the semivariogram structure did not allow model fitting or because, whenever it did, the models showed weak spatial dependence. This finding was confirmed by high nugget effect values, according to which, variations in these attributes resulted from random factors rather than from spatial dependence contributions. Consequently, the spatialization process was conducted through the inverse distance weighted method.

Table 2.
Geostatistical parameters set for soil chemical attributes in the 0-20 cm soil layer.

Variable	Model	SDI	Co	C	A	AR	ER	Sre	AIC
H+Al	Spherical	100.00	0.000	7.616	47.294	-1.52E-03	0.000	0.001	245.49
	Exponential	100.00	0.000	7.629	29.778	-4.22E-03	-0.001	0.008	246.94
	Gaussian	87,34	1.172	8.084	27.309	3.61E-03	-0.004	0.031	256.13
Al	Spherical	100.00	0.000	1.913	4.210	-1.29E-02	0.253	2.166	262.53
	Exponential	20.58	1.520	0.394	0.841	-4.23E-03	0.254	2.169	262.53
	Gaussian	13.37	1.657	0.256	2.386	-9.51E-03	0.258	2.204	262.53
S	Spherical	95.64	0.304	6.676	19.862	1.59E-02	0.012	0.101	305.29
	Exponential	98.38	0.198	12.037	23.068	2.44E-02	0.014	0.121	303.64
	Gaussian	71.43	1.532	3.830	1.100	1.65E-08	0.017	0.147	337.75
T	Spherical	94.47	0.752	12.854	35.517	4.48E-03	0.000	0.003	324.08
	Exponential	99.42	0.070	12.065	15.303	1.28E-02	0.004	0.034	325.07
	Gaussian	77.04	3.561	11.953	26.803	1.55E-02	0,006	0.048	331.05
V	Spherical	93.41	38.135	540.372	57.190	9.90E-02	-0.005	0.043	577.81
	Exponential	95.41	31.626	657.131	40.000	1.15E-01	-0.004	0.032	579.01
	Gaussian	84.64	73.698	406.067	13.961	1.49E-01	0.002	0.021	579.68
pH	Spherical	100.00	0.000	0.384	21.760	-1.15E-03	0.001	0.010	80.82
	Exponential	100.00	0.000	0.971	40.297	7.67E-04	0.003	0.030	75.86
	Gaussian	86.97	0.074	0.492	13.807	2.02E.-03	0.008	0.070	79.91
Corg	Spherical	88.61	0.057	0.446	13.424	7.39E-03	0.017	0.14	140.65
	Exponential	100.00	0	0.523	6.447	7.54E-03	0.017	0.15	140.41
	Gaussian	53.91	0.297	0.347	22.774	8.07E-03	0.020	0.17	143.75

SDI: spatial dependence index; Co: nugget effect; C: contribution; A: range; EA: absolute error; ER: mean standardized error; Ser: standard error of the standardized errors; AIC: Akaike information criterion; H+Al: hydrogen and aluminum; Al: aluminum; S: sum of bases; T: total cation exchange capacity; V: base saturation index; pH: pH in water; Corg: organic carbon; NA: model could not be fitted; NV: cross-validation could not be performed for the model.

Table 3.
Geostatistical parameters for soil chemical attributes in the 20-40 cm soil layer.

Variable	Model	SDI	Co	C	A	AR	ER	Sre	AIC
H+Al	Spherical	36.68	1.506	0.872	38.531	5.31E-03	0.010	0.085	263.69
	Exponential	100.00	0.000	2.174	3.444	8.44E-03	0.010	0.087	264.46
	Gaussian	35.54	1.677	0.925	24.611	5.78E-03	0.011	0.090	264.31
Al	Spherical	100.00	0.000	0.020	5.125	2.95E-04	0.016	0.137	-70.73
	Exponential	100.00	0.000	0.020	1.101	9.10E-04	0.022	0.187	-70.43
	Gaussian	100.00	0.000	0.020	2.398	5.01E-04	0.018	0.152	-70.60
S	Spherical	56.65	1.249	1.633	63.726	-2.90E-02	0.011	0.090	255.57
	Exponential	100.00	0.000	2.394	6.111	2.64E-03	0.036	0.304	253.75
	Gaussian	100.00	0.000	2.207	4.951	1.30E-02	0.037	0.314	251.55
T	Spherical	51.05	2.003	2.088	17.645	6.41E-03	0.009	0.076	305.61
	Exponential	100.00	0.000	3.976	4.098	1.80E-02	0.013	0.112	304.42
	Gaussian	100.00	0.000	3.954	4.573	2.92E-02	0.015	0.128	300.92
V	Spherical	52.38	67.044	73.752	17.753	-2.82E-02	0.004	0.033	563.21
	Exponential	54.38	78.686	93.814	26.159	-2.30E-02	0.003	0.022	561.48
	Gaussian	100.00	0.000	136.820	4.469	1.16E-02	0.007	0.060	561.24
pH	Spherical	86.89	0.034	0.226	9.731	1,10E-03	0.007	0.059	107.34
	Exponential	100.00	0.000	0.259	3.603	1.95E-03	0.010	0.085	108.30
	Gaussian	74.26	0.067	0.193	4.794	NV	NV	NV	107.68
Corg	Spherical	57.94	0.165	0.228	15.112	8.34E-04	0.048	0.410	134.57
	Exponential	100.00	0.000	0.380	3.719	2.53E-03	0.050	0.428	135.42
	Gaussian	47.10	0.201	0.179	6.757	-5.53E-03	0.035	0.298	134.71

SDI: spatial dependence index; Co: nugget effect; C: contribution; A: range; EA: absolute error; ER: mean standardized error; Ser: standard error of the standardized errors; AIC: Akaike information criterion; H+Al: hydrogen and aluminum; Al: aluminum; S: sum of bases; T: total cation exchange capacity; V: base saturation index; pH: pH in water; Corg: organic carbon; NA: model could not be fitted; NV: cross-validation could not be performed for the model.

Spatial estimate maps were generated through ordinary kriging, based on the best models selected for each soil chemical feature in both assessed soil layers ( Figure 3 ). All chemical features presented concentration variations between the 0 - 20 cm and 20 - 40 cm soil layers, to varying extents.

The smallest variations were observed for aluminum content ( Figure 3 E), aluminum saturation index ( Figure 3 I) and potassium content ( Figure 3 N). These features mostly maintained lower classes in both assessed soil layers. Attributes, such as sodium ( Figure 3 A), hydrogen and aluminum ( Figure 3 D), total cation exchange capacity ( Figure 3 G), base saturation index ( Figure 3 H), sodium saturation index ( Figure 3 J), pH in water ( Figure 3 K), organic carbon ( Figure 3 L) and phosphorus ( Figure 3 M), recorded decreasing concentrations as soil layers got deeper. This spatial pattern was observed in the same region where the study site is located in - higher value classes were recorded in the 0-20 cm soil layer and lower value classes were observed in the 20-40 cm layer. Calcium content presented the most significant spatial distribution similarity between the two assessed soil layers.

Silva et al. (( 2016 )) observed negative spatial correlation between pH and Al ³⁺ , organic carbon and cations, phosphorus and total clay, and silt and sand. This finding allows inferring that ordinary kriging can be applied to spatialize soil chemical and textural attributes in agroforestry systems.

Navidi & Seyedmohammadi (( 2022 )) have analyzed spatial variability in soil CEC, total nitrogen, available potassium and phosphorus contents based on the Ordinary Kriging method. NRMSE (normalized root mean square error), NMAE (normalized mean absolute error) and R² (coefficient of determination) values observed by the aforementioned authors pointed out the maps’ good accuracy. According to them, spatial distribution maps can be used as appropriate basis to achieve more accurate and specific nutritional management in agricultural areas, as well as to help protecting the environment by preventing the contamination of underground water resources.

Soil texture attributes did not show spatial dependence; therefore, they were spatialized through the inverse distance weighted method. Plotted maps have evidenced higher heterogeneity in the distribution of sand, clay and silt concentrations in both assessed soil layers ( Figure 4 ). Clay prevailed in the 0-20 cm soil layer, whereas silt prevailed in the 20-40 cm layer.

Results reported by Seyedmohammadi et al. (( 2019 )) have evidenced accurate spatial distribution of soil texture. According to the aforementioned authors, understanding the spatial variability of soil properties, such as texture, can be an important tool in land-use planning processes aimed at mitigating potential soil losses.

Figure 3.
Spatial distribution of soil chemical attributes recorded for the 0-20 cm and 20-40 cm soil layers in the investigated agroforestry system, based on ordinary kriging interpolation. Wherein A: sodium (cmol _c .dm ^-3 ); B: calcium (cmol _c .dm ^-3 ); C: magnesium (cmol _c .dm ^-3 ); D: hydrogen + aluminum (cmol _c .dm ^-3 ); E: aluminum (cmol _c .dm ^-3 ); F: sum of bases (cmol _c .dm ^-3 ); G: total cation exchange capacity (cmol _c .dm ^-3 ); H: base saturation index (%); I: aluminum saturation index (%); J: sodium saturation index (%); K: pH in water (1:2.5); L: organic carbon (%); M: phosphorus (mg.kg ^-1 ); and N: potassium (mg.kg ^-1 ).

Figure 4.
Spatial distribution of soil texture attributes in the 0-20 cm and 20-40 cm soil layers, in the investigated agroforestry system through inverse distance weighted interpolation. A represents sand content (%), B represents clay content (%), and C represents silt content (%).

4. CONCLUSIONS

The investigated soil chemical attributes were spatially structured, and it enabled performing estimates across the study site through geostatistical methods. On the other hand, textural did not present spatial dependence; thus, they should be spatialized through deterministic methods, such as inverse distance weighted interpolation.

Overall, the concentrations of chemical and textural attributes decreased as soil layers got deeper, except for silt content, which presented the highest strata classes in the 20-40 cm layer.

Textural attributes presented low contribution in the 0-20 cm soil layer, whereas chemical attributes showed a more homogeneous contribution in the 20-40 cm layer.

REFERENCES

• Abhiyan WNT. Sustainable Development Goals: Agenda 2030. India 2017: A Civil Society Report. [cited 2017]. Available from: < https://www.kractivist.org/wp-content/uploads/2017/07/Civil-society-Report-on-SDGs.pdf >.
  » https://www.kractivist.org/wp-content/uploads/2017/07/Civil-society-Report-on-SDGs.pdf
• Alvares CA, Stape JL, Sentelhas PC, Gonçalves JLM, Sparovek G. Köppen’s climate classification map for Brazil. Meteorologische Zeitschrif, 2014; 22(6): 711-728.
• Andrade AG, Caballero SSU, Faria SM de. Ciclagem de nutrientes em ecossistemas florestais. Rio de Janeiro: Embrapa Solos; 1999.
• Araújo ASF, Melo WJ. Soil microbial biomass in organic farming system. Ciência Rural, 2010; 40(11): 2419-2426.
• Araújo EJG, Lima DO, Curto RA, Silva EV, Arthur Júnior JC. Space pattern of forest species and its relationship with agricultural factors agroflorestal system. Floresta, 2019; 49(2), 335-344.
• Burnham KP, Anderson DR. Model selection and multimodel inference: a practical information-theoretic approach. 2 ed. New York: Springer; 2002.
• Canuto JC. Sistemas Agroflorestais: experiências e reflexões. In: Vaz, P. Agroflorestas, clareiras e sustentabilidade. Brasília, DF: Embrapa; 2017.
• Cambardella CA, Moorman TB, Novak JM, Parkin TB, Karlen DL, Turco RF et al. Field-Scale Variability of Soil Properties in Central Iowa Soils. Soil Science Society of America Journal, 1994; 58(5): 1501-1511.
• Higuchi N, Santos J, Ribeiro RJ, Minette L, Biot Y. Biomassa da parte aérea da vegetação da floresta tropical úmida de terra-firme da Amazônia Brasileira. Acta Amazonica, 1998; 28(2): 153-166.
• Hongyu K, Sandanielo VLM, Oliveira Junior GJ. Análise de Componentes Principais: resumo teórico, aplicação e interpretação. Engineering and Science, 2015; 1(5): 83- 90.
• Iwata BF, Leite LFC, Araújo ASF, Nunes LAPL, Gehring C, Campos LP. Sistemas agroflorestais e seus efeitos sobre os atributos químicos em argissolo vermelho-amarelo do cerrado piauiense. Revista Brasileira de Engenharia Agrícola e Ambiental, 2012; 16(7): 730- 738.
• Kassambara A, Mundt F. Extract and Visualize the Results of Multivariate Data Analyses. [cited 2017]. Available from: < https://CRAN.R-project.org/package=factoextra >
  » https://CRAN.R-project.org/package=factoextra
• Lê S, Josse J, Husson F. FactoMineR: An R Package for Multivariate Analysis. Journal of Statistical Software, 2008; 25(1): 1-18.
• Navidi M N, Seyedmohammadi J. Mapping and spatial analysis of soil chemical effective properties to manage precise nutrition and environment protection. International Journal of Environmental Analytical Chemistry, 2022; 102(8): 1948-1961.
• Pezarico CR, Vitorino ACT, Mercante FM, Daniel O. Indicadores de qualidade do solo em sistemas agroflorestais. Revista de Ciências Agrárias, 2013; 56(1): 40-47.
• R Core Team. R: A language and environment for statistical computing.Vienna, Austria. R Foundation for Statistical Computing, 2015.
• Ribeiro Júnior PJ, Diggle, PJ. GeoR: a package for geostatistical analysis. RNews, 2001; 1(2): 15-18.
• Santos SEM, Miranda IS, Tourinho MM. Estimativa de biomassa de sistemas agroflorestais das várzeas do rio Juba, Cametá, Pará. Acta Amazonica, 2004; 36(1): 01-08.
• Sen DOU, Shan J, Song X, CAO R, Wu M, Li C, Guan S. Are humic substances soil microbial residues or unique synthesized compounds? A perspective on their distinctiveness. Pedosphere, 2020; 30(2):159-167.
• Seyedmohammadi J, Matinfar H R. Statistical and geostatistical techniques for geospatial modeling of soil cation exchange capacity. Communications in Soil Science and Plant Analysis, 2018; 49(18): 2301-2314.
• Seyedmohammadi J, Navidi M N, Esmaeelnejad L. Geospatial modeling of surface soil texture of agricultural land using fuzzy logic, geostatistics and GIS techniques. Communications in Soil Science and Plant Analysis, 2019; 50(12): 1452-1464.
• Sibaldelli RNR, Carvalho JFC, Oliveira MCN. Uso de Geoestatística no Estudo da Variabilidade Espacial da Capacidade de Troca de Cátions do Solo. Global Science and Technology, 2015, 8:141-156.
• Silva CSD, Pereira MG, Delgado RC, Silva EV. Spatialization of soil chemical and physical attributes in an agroforestry system, Seropédica, Brazil. Cerne, 2016; 22(4): 407-414.
• Steenbock W, Vezzani FM. Agrofloresta: aprendendo a produzir com a natureza. 1 ed. Curitiba: Bambual; 2013.
• Uzêda MC, Tavares PD, Rocha FI, Alves RC. Paisagens agrícolas multifuncionais: Intensificação ecológica e segurança alimentar. 1 ed. Brasília: Embrapa; 2017.
• Vale Júnior JF, Souza MIL, Nascimento PPRR, Cruz DLS. Solos da Amazônia: etnopedologia e desenvolvimento sustentável. Revista Agroambiente On-line, 2011; 5(2): 158-165.

Edited by

Associate editor:

Fernando Gomes http://orcid.org/0000-0003-0363-4888

Publication Dates

• Publication in this collection
  09 Sept 2024
• Date of issue
  2024

History

• Received
  25 June 2024
• Accepted
  16 Aug 2024