Study area
This study was carried out in a selected area of 225 ha in the Arıklı soil series located in the Çukurova region of southern Turkey, which is one of the most fertile plains used for intensive agricultural production (Celik et al. 2011) (Fig. 1). The study coordination is between 4,081,768 and 4,083,268 N and 709,876–711,376 E (WGS-84, UTM-m, 36 Zone; 9-12 m above sea level).
Fig. 1
Arıklı soil series is the most common series with 24% (56,083 ha) in the region. The soil has a clay texture and slope of less than 1, formed on the alluvial deposits of old river terraces (Dinç et al. 1995), and classified as Entic Chromoxerert (Soil Survey Staff 1999). The prevailing climate of the study area is Mediterranean with the summers hot and dry and the winters rainy and mild. According to the long term (1929–2019) meteorological data (TSMS 2020), mean annual temperature, total annual evapotranspiration and total annual rainfall in the study area are 19.2 °C, 1557 mm and 639 mm, respectively. About 75% of the annual rainfall falls in the winter and spring months. Soil temperature and moisture regimes of the study area are thermic and xeric, respectively.
In the study area, agricultural production based on arable crops has been carried out for long period of time (more than 50 years). Crop production in arable lands at the time of soil sampling was majorly corn (Zea mays L.), soybean (Glycine max. L.), cotton (Gossypium hirsutum L.), melon (Cucumis melo L.) and watermelon (Citrullus lanantus). Soils in the arable lands are conventionally tilled with a moldboard plow (Celik et al. 2011). In addition, residues are removed from the land by burning in most places.
The proliferation of citrus species suitable for cultivation in the plain in recent years, increased the area that were converted to horticulture lands, namely to mandarin (Citrus reticulata L.), grapefruit (Citrus paradisi L.), lemon (Citrus lemon L.) and pomegranate (Punica granatum L.) that have been grown for 17 years on average in the plain with the oldest since 1995 and the youngest since 2015.
Soil sampling and analyses
The study was carried out on a land section of 1.5 km x 1.5 km width within the Arıklı soil series (Fig. 1). This section was divided into square grids of 250 m x 250 m for soil sampling. Soil samples were collected at 49 main points, approximately at the corner point of each grid. In addition, intermediate sampling lines as transect were placed between the corner point of 6 grids in order to more accurately predict changes in soil properties at distances shorter than 250 meters. On these lines, a total of 36 soil samples were collected at 25 m, 75 m and 150 m distances in the south-north and east-west directions, and 6 samples were collected from these points at a distance of 175 m in the north-east, north-west, south-east and south-west directions. In total, 42 intermediate sampling points along a transect added up to 91 soil samples. Soil samples were collected at a depth of 0-20 cm, which is the preferred depth in many studies where the effects of different land use were studied (Wang et al. 2012; Wu et al. 2019). In the study, disturbed and undisturbed soil samples were collected at each point in September 2018. All samples were carefully transported to the laboratory. After the removal of roots and coarse particles, soil samples were prepared for physical, chemical, biological and fertility analysis. All analyses were performed with three soil replicates. The list of soil properties and analytical protocols selected are presented in Table 1.
Table 1
Indicator selection and SQI
SQI is generally determined with four steps: (i) selection of indicators, (ii) scoring of selected indicators, (iii) determination of indicator weights using FAHP, and (iv) integrating the weighted indicator scores into a comparative SQI using Weighted Linear Combination (Gozukara et al. 2022). In this study, physical, chemical, biological and fertility indicators were carefully selected by soil experts. While deciding on the indicators with the opinion of experts, based on the recommendation of the experts, climate, soil properties and land use were considered.
AS, BD, AWC and WFPS used in this study are common physical properties used in soil quality related studies (Cherubin et al. 2016; Budak et al. 2018; Dengiz 2020; Çelik et al. 2021). Blanco Canqui and Ruis (2018) indicated that changes in soil physical properties affect all ecosystem services important to life that soil provides such as food, feed, fiber and fuel production, etc. Therefore, soil physical properties were the indicators considered in assessment of soil quality.
Soil chemical characteristics such as soil pH, EC, SAR, CaCO3, nutrients were preferred by researchers in many studies due to their effect on nutrient availability, root growth and fertility (Karaca et al. 2021; Mahajan et al. 2021; Gozukara et al. 2022; Kılıc et al. 2022). Intensive use of chemical fertilizers, especially nitrogen fertilizers, for agricultural production in the study area (İbrikçi et al. 2015) highlights soil pH as the core soil properties responsible for nutrient availabilities in the area (Zurovec et al. 2021). As excessive and long term chemical fertilization may cause salinity problems in the soil, EC was taken into account as one of the parameters to be monitored. SAR, an important determinant of sodium problem, is also one of the indicators that should be considered in study areas where furrow irrigation methods are used.
OC is an important property affecting soil physical, chemical, biological and fertility properties and processes associated with them. OC improves soil structure (Ozlu and Kumar 2018), enhances microbial activity (Koçak and Cenkseven 2021), soil functions and soil quality (Çelik et al. 2021), and increases crop productivity (Lal 2004) and thus is preferred in many studies as an important indicator of SQI (Raiesi and Kabiri 2016; Liu et al. 2018; Acir et al. 2022). In addition, MBC, BGEA and qCO2, which are significantly related to organic carbon (Bastida et al. 2008), were preferred in soil quality studies as they are properties that respond quickly to change (Bhaduri et al. 2017).
After selection of indicators, standard scoring functions were used to convert them to unitless scores in the second step (Table 2). In the literature, there are three scoring functions, “more is better”, “less is better” and “mid-point is optimum”, according to where the best soil functionality is considered as low, high or medium values or ranges (Liebig et al. 2001). In this study, AS, AWC, OC, MBC, BGEA, K, Fe, Cu, Mn and Zn due to their roles in soil fertility were scored with “more is better” function while BD, EC, SAR, CaCO3 and qCO2 due to their indicators of compaction, salinity and poor quality were scored with “less is better” function. When scoring pH, P and WFPS, “mid-point optimum” function is preferred (Andrews et al. 2004; Wienhold et al. 2009). In the scoring of pH, the most desirable limits for plant growth are between 6.0 and 7.0 (Marschner 2011). In this study, as all pH values are above 7.0, “less is better” function was used for soil pH. Bai et al. (2013) reported that the critical Olsen-P levels for optimum crop yield ranged from 10.9 mg kg-1 to 21.4 mg kg-1. Therefore, P values within the range were scored as 1. “Less is better” function was used for P values above the range and “more is better” function was used for P values below the range. Linn and Doran (1984) stated that WFPS is around 60% for optimum microbial activity. In this study, this value was accepted as the critical value. “Less is better” function was used for WFPS > 60% and “more is better” function was used for WFPS < 60%.
Table 2
Fig. 2
Table 3
In the third step, all indicators were grouped into four sub-criteria: physical, chemical, biological and fertility (Fig. 2). After that, the scored indicators in the groups and groups in the SQI were weighted with FAHP. AHP is one of the multi-purpose decision making methods which considers the opinion of a group or individual as a priority evaluating qualitative and quantitative variables together (Saaty 1980). AHP has been deployed in solving many real life decision making problems (Dağdeviren 2007). However, AHP does not necessarily fully reflect the complexity of the non-numeric nature of the pairwise comparisons due to the single numeric values taken from the 1–9 scale (Pedrycz and Song 2014). AHP is also criticized for its inadequacy to handle situations of uncertainty and indecision (Deng 1999). In order to enhance the effectiveness of AHP in decision making process, the method was integrated with fuzzy logic and FAHP concept emerged. Van Laarhoven and Pedrycz (1983) proposed the first triangular method adopted in FAHP. Chang (1996) further developed it for pairwise comparison of triangular fuzzy numbers (Table 3) and has been used in many studies in recent years (Turan et al. 2020; Kılıc et al. 2022; Kaya et al. 2022). To combine fuzzy logic and AHP, different methods have been proposed in the literature (Buckley 1985; Chang 1996; Deng 1999). In this study, the extent analysis method developed by Chang (1996) was used to integrate both methods. The main component of the equation used in FAHP is given in Eq. (1).
where l ≤ m ≤ u, l and u are the lower and upper values of M respectively, and m is taken as the modal value. l, m, and u are the the triangular fuzzy number. M is the set of the defined elements. if all the triangular fuzzy numbers are equal (l = m = u) it is discounted as a fuzzy number. After determining the relative importance levels of the indicators, the selected indicators were combined into an overall SQI using Eq. (2).
where SQI is the soil quality index, W is the relative importance of the indicator with FAHP and S is the indicator score.
Descriptive statistics and spatial variability of SQI
Descriptive statistics were conducted using the SPSS 21.0 statistical software (SPSS, Inc., Chicago, IL). Spatial variability of SQI was determined using geostatistical techniques. Modeling of experimental semivariogram was performed with GS+ (version 7) (Gamma Design Software). The best fit model for SQI was determined with RSS (Residual Sum of Squares) and r2 values of the semivariogram model. Models with RSS values close to zero and r2 values close to 1.0 were selected as the best fitted model (Yang et al. 2011). The most appropriate parameters (range, sill, nugget, lag number) obtained by semi-variogram models were used to obtain the spatial distribution map of SQI. The kriging map of SQI was created using the geostatistical extension of ArcGIS 10.2. SQI was mapped considering very low (<40%), low (40-55%), medium (55-70%), high (70-85%) and very high (>85%) grades specified by Gugino et al. (2009).