Developing scoring functions to assess soil quality at a regional scale in rangelands of SW Spain

The drawing of maps of soil quality at a large scale is increasingly being more useful to land planners and stakeholders. Nevertheless, it involves different methodological steps from the description of soil profiles in the field until the regional mapping of integrative soil quality index (IQI) values. The development of proper scoring functions is a paramount task for the calculation of these IQI values since every parameter needs to be standardized accordingly and weighting factors are usually estimated by multivariate techniques. The main goal of this study was to map soil quality in the Spanish region of Extremadura (commonly known by its rangelands called dehesas). To do that, i) we gathered information from 194 soil profiles described throughout the region, ii) we calculated the weighting factors of ten meaningful parameters used as indicators by using multivariate techniques (Principal Component Analysis, PCA; and Analytic Hierarchy Process, AHP), and iii) we developed standard scoring functions (SSFs) that represent the singularity of every variable (less is better, more is better). We established upper and lower limits for standardizing the values of each indicator properly. Regarding weighting factors, soil texture was highlighted by the PCA and nutrients by the AHP. Once IQI values were calculated, two regional maps of soil quality were drawn by using interpolation methods (ordinary kriging). The IQI maps showed remarkable spatial differences in soil quality presumably induced by land management. We conclude this methodology could be useful and we encourage other colleagues to test its effectiveness in places where soil data are available.


INTRODUCTION
Soil quality can be defined as the capacity that a soil has to function properly within its ecosystem limits, i.e., providing food for people and animals, storing and cleaning water, keeping good habitat conditions for several species, etc. (Bünemann et al., 2018). The generalized interest of scientists for measuring soil quality emerged in the early 1990s after many decades of research focused on the estimation of erosion rates as the way of evaluation of the land management practices effects (Andrews et al., 2004). Nevertheless, the assessment of soil quality has been a controversial topic that has always required a wide knowledge of the system object of study (Liu et al., 2017).
The use of specific indicators is the most consistent methodology known so far for assessing soil quality (Laudicina et al., 2017). In addition, Zornoza et al. (2015) remarked the importance of adopting soil quality indicators of different natures (physical, chemical, and biological origins) to achieve a holistic point of view. These indicators must be selected and undoubtedly interpreted, and after that, the most important task is the integration in a comprehensive index (Andrews et al., 2002). This indexing process requires indicators that have been previously scored and combined (Shiri et al., 2017). Therefore, a soil quality index should return a single value from a sum or multiplication (using weights or not) from different indicators scored (Mukherjee and Lal, 2014). One of the most known indexes is the integrative quality index (IQI) proposed by Qi et al. (2009).
The weighting method is a permanent challenging question since the variety of criteria can be used (Lima et al., 2013). The weights of each variable are usually obtained from statistical-based procedures (e.g., multivariate analysis). Principal Component Analysis (PCA), for instance, has been frequently used in the process of obtaining a minimum data set (MDS), i.e., to reduce the total number of indicators (Stefanoski et al., 2016). Meanwhile, the Analytic Hierarchy Process (AHP) has been already utilized for many authors (Ying et al., 2007) to determine the weight in the evaluation of environmental quality at a regional scale.
The standardization of each value is also a controversial topic since the statistical heterogeneity of the parameters that are usually considered as indicators. The most common procedure is to convert their real values into values within a range from 0 to 1 (0 = lowest quality and 1 = highest quality) by using commonly known methods such as min-max scaling and/or z-score standardization. Pulido et al. (2017), for instance, classified real values in five classes following a rationale based on a Likert scale.
The development of standard scoring functions (SSFs) allows, therefore, obtain more accurate standardized values since they are based on the previous knowledge of every variable (e.g., bulk density is better when its values decrease meanwhile soil organic carbon is better when its values increase). They have been widely used for assessing soil quality in agricultural soils of many parts of the world (Li et al., 2018). Although their calculation is based on consistent methods largely validated at a regional scale, at least in agricultural soils, there are no published studies (up to our knowledge) on rangelands, where millions of people depend directly on keeping soils in good conditions that can guarantee enough pasture production for feeding their livestock. Pulido et al. (2017) assessed soil quality at a farm-scale in rangelands of Extremadura Region (Spain). Nevertheless, they did not use both weighting factors and scoring functions for calculations. In addition, there is not a regional map that allows the assessment of soil quality at a large scale. For these reasons, the main goal of this study was to draw maps of soil quality for Extremadura. To do that, we calculated weighting factors by using multivariate techniques (PCA and AHP) and to develop SSFs for every indicator.
As indicators, we selected ten soil properties (sampling at the layer 0.00-0.30 m): sand, silt, clay, bulk density (BD), total nitrogen, available phosphorous (P), potassium (K), pH(H 2 O), soil organic carbon (SOC), and cation exchange capacity (CEC). We hypothesize the AHP generates more realistic maps than PCA since we can select the most influencing factors from our knowledge on the spatial variability of the soils of the study area.

Study area
This study was carried out in the Autonomous Region of Extremadura (SW Spain) (Figure 1). Extremadura covers 41,633 km 2 (8.2 % of Spain) where the 50-60 % out of its territory is occupied by rangelands and grasslands, the remaining land is agricultural (Guadiana River) and mountainous areas (Spanish Central System at the North). The rangelands are located at an elevation ranging from 300 to 700 m a.s.l. dominated by undulated terrains (2-8 % of inclination) over old erosion surfaces of shales, schists, greywackes, and granites. The climate is Mediterranean with dry and warm summers from June to September. Soils are sandy loam-textured, shallow, and poor in nutrients the main soil types present in the area are: Leptosols and (Epileptic) Cambisols (IUSS Working Group WRB, 2015).
Extremadura is inhabited by more than one million people (2 % out of whole Spain), recording a population density of ≈25 people km -2 , spread out by two provinces (Cáceres and Badajoz) and 388 municipalities. Extremadura's population has still a noticeably rural character since more than 25 % of people are living in places with a population lower than 10,000 inhabitants. Thirteen percent of active people work in the primary sector (mainly agriculture and livestock husbandry) providing 6 % of the regional gross domestic product (GDP). So, although rangelands do not suppose a significant contribution in terms of GDP, they are particularly interesting in terms of livestock production (735 tons of meat per year), biodiversity conservation (habitat of special protection by the European Union), and rural idiosyncrasy.

Soil profile description
The soil profile dataset is composed of 194 soil profile descriptions: 47 out of them were conducted during the fieldwork of the Ph.D. thesis of the first author (Pulido Fernández, 2014), 116 during former research projects between 2003 and 2006, and 31 descriptions were taken from the Extremenian Catalog of Soils.

Work scheme and rationale
Figure 2 shows a conceptual scheme in which the rationale of this study can be understood. The core of this study is the Qi et al. (2009) equation for calculating the integrative quality index (IQI) and the final target is the drawing of maps of soil quality at a regional scale. Since we selected ten indicators from our previous knowledge, and we gathered information from about 194 soil profiles described throughout Extremadura, the first step was the calculation of weighting factors for each indicator. For doing that, we used multivariate Rev Bras Cienc Solo 2020;44:e0200090 techniques such as Principal Component Analysis (PCA) and Analytic Hierarchy Process (AHP). Secondly, we developed standard scoring functions (SSFs) to standardize every indicator. We distinguished two types of functions: less is better (L) and more is better (M). Finally, two IQIs were calculated for the 194 sampling points and these values were interpolated (ordinary kriging) to obtain two regional maps of soil quality.

Weighting factors
We utilized two sets of weighting factors, one obtained from the PCA and the other one from the AHP. Both analyses were made by using IBM SPSS 22.0 and MO Excel 2013 software packages, respectively.
The PCA weighting factors were assigned from the commonality of each indicator calculated according to equation 1.

Weighting factor Ix = (Commonality Ix)/(Sum of every communalities)
Eq. 1 The weighting factors of the AHP set were obtained following these methodological steps (Saaty, 1990

Standard scoring functions
We applied standard scoring functions (SSFs) (Mukherjee and Lal, 2014) in which the calculation of the upper and lower limits are needed. Since there is no information about the upper and lower thresholds in the specific rangeland areas of Extremadura, the minimum and maximum observational values of the variables in the region were considered as lower and upper thresholds, respectively. Regarding the fact that the range of soil reaction can affect soil quality, the values of the range had a score equal to 1 and increasing the distance from this range (more or less), the score decreased. The final results were classified as "less is better" and "more is better".

Soil quality mapping
Once every indicator was both standardized (range: 0-1) by the scoring functions and its two weighting factors (PCA and AHP) were assigned, the two IQI (Qi et al., 2009) values (PCA and AHP) of each one of the 194 sampling points were calculated. These values were interpolated by using ordinary kriging, embedded in the Geostatistical Wizard Tool by ArcGIS 10.5 (ESRI), and the resulting raster was clipped to the boundaries of the region of Extremadura. Finally, the two maps generated were interpreted to decide which one is more reliable.

Descriptive statistics of the whole dataset
The total dataset is comprised of ten soil quality indicators reporting information about soil texture, acidity, cations and nutrients contents, and bulk density. Table 1 summarizes the descriptive statistics of the 194 soil profiles at a sampling layer from 0.00-0.30 m. Soil texture is characterized by 14.1, 38, and 47.9 % of clay, silt, and sand contents, respectively. It is important to highlight that extreme maximum values can also be found (45.0 % clay, 62.7 % silt, and 82.1 % sand), being clay and silt contents the most variable.
The mean pH(H 2 O) values were close to 5.6, with maximum values of 8.2 and a minimum of 4.1, giving us information about the acidity of the rangeland's soils. The average value of CEC was 10.65 cmol c kg -1 , reaching maximum values of 31.54 cmol c kg -1 and minimum ones of 3.33 cmol c kg -1 . Total N reaches mean contents of 0.14 %, increasing in some places till 0.49 %. The mean Olsen-P value was 17.08 mg kg -1 and available K 0.49 mg kg -1 . These values are much related to the low SOC (1.9 %) found in the dehesas.

Correlations among variables
The Pearson linear correlation coefficients of the soil quality indicators are shown in figure 3. Regarding soil texture variables (clay, silt, and sand contents) silt and sand showed a high negative correlation (r = -0.85). This fact is due to the dominance of loam soils with slight variations in texture, depending on a higher content of sand (sandy loam) or silt (silty loam). Clay, by its part, was positively correlated with the cation exchange capacity (r = 0.70). Soil compaction (expressed as bulk density) was negatively correlated with almost all variables. Chemical and biological properties were logically correlated in a positive way among them and negatively with bulk density and the sand content.
The CEC seemed to be the variable that explains the best the existing relationships between every indicator. It was positively correlated to the pH (r = 0.43) and with the parameters usually utilized as fertilizers: N (r = 0.36), P (r = 0.41), and K (r = 0.40). Contrariwise, it was negatively correlated to bulk density (r = -0.29) that was also negatively correlated to N (r = -0.47). Remarkable are also the role played by P and SOC as indicators of fertility. They were positively correlated between them (r = 0.37) and kept positive correlations with K (K-P: 0.53; K-SOC: 0.34) and N (N-P: 0.29; N-SOC: 0.42).

Principal component analysis
The PCA sorted the ten variables studied according their importance. The CEC and clay content showed the highest predictive power (0.779 and 0.684, respectively) while the lowest values were recorded with SOC and bulk density (0.365 and 0.317, respectively). Nutrients content (N, P, and K) showed importance lower than soil texture, pH, and CEC. However, their variables showed similar power values of prediction (Figure 4).    (Burt, 2004); pH (1:2.5 soil/water ratio); cation exchange capacity (CEC) (MAPA, 1982); total nitrogen according to Kjeldahl (1883); phosphorous according to Olsen et al. (1954); available potassium according to MAPA (1982); soil organic carbon according to Walkley and Black (1934); and bulk density according to Coile (1936). SD: standard deviation; CV: coefficient of variation.

Standard scoring functions, weight assignment and integrative quality index
Each indicator was transformed and normalized in a range from 0 to 1 and their functions varied according to their logical behavior, i.e., properties that improve if their values  increase ("more is better" pattern, M) and those that improve when their values decrease ("low is better" pattern, L). The detailed equations are shown in the table 3, including the upper and lower limits within the system to function properly.
The weight of each indicator is calculated by its communality obtained by two ways: PCA and AHP analysis (table 4). In the PCA, soil texture, pH, CEC, and available K obtain a weight higher than 0.100, meanwhile, in the AHP, the highest weights are obtained by the available K and Olsen-P.
Then, in table 5, we show the descriptive statistics of the integrated quality index in the 194 sites according to its weights returned by PCA and AHP analysis, respectively. The PCA returned an IQI with mean values higher than AHP but this last shows a wider range of qualities (0.154-0.705). Figures 5 and 6 show the spatial distribution of IQI throughout Extremadura using PCA and AHP, respectively. In the first map (PCA), the IQI values ranged from 0.203 to 0.682, meanwhile in the second one (AHP), the IQI values ranged from 0.154 to 0.705. This higher data range of AHP could be influenced by the higher weights obtained for the contents of phosphorus (0.220) and potassium (0.201).    In both maps, the province of Badajoz (southern half of the territory) evidences a higher quality in its soils, particularly using AHP procedure. The soils with lower values (green colors) correspond to the northern part of the region (mountainous areas in which data of soil properties are scarce) and eastern part of the province of Badajoz dominated by treeless grasslands. Remarkable differences have been observed in the province of Cáceres (northern half of the region) in which rangeland soils are mostly considered as of medium quality (predominance of yellow instead of red color).

DISCUSSION
This study can be interpreted as an up-scaling (from farm to the regional scale) of previous studies carried out by Pulido et al. (2017). The role played by potassium and clay (in sandy/silty loam soils) looks keys in the understanding on soil quality in rangelands of SW Spain. From a methodological point of view, we propose an index for assessing soil quality more sensitive (weight vs. average values), integrative (more indicators), and easy to use (there are much information free on Internet). The pros and constraints of using different analysis (PCA vs. AHP) in the scoring process have been also addressed.
The analysis performed in this study evidenced a great importance of some couples of variables such as CEC, clay (PCA), phosphorous, and potassium (AHP) on soil quality. These findings were in consonance with two indicators (CEC and K) also proposed by Pulido et al. (2017) in their research at a farm scale. Nevertheless, it keeps alive the never ended discussion about the origin of a good soil. That is, it shows a high quality because belong to areas with deep soils rich in clay and nutrients (e.g., on alluvial deposits) or its quality has been improved by land management practices such as a long-term fertilization (animal excreta, improved pastures, etc.) or avoiding soil erosion.
Phosphorous is considered by a majority of regional farmers and extensionists (personal communications in several meetings) as the key nutrient for fertilizing natural pastures  due to their rapid effects on pasture production (Barbosa et al., 2009;Recena et al., 2017). Potassium has been hardly ever considered in any soil indicators-based system in rangeland areas. Pulido et al. (2017) selected it due to its strong correlations with pasture production as a good explaining factor of its spatial variability regardless rainfall amounts. Potassium is commonly known for increasing the tolerance of plants against droughts and frosts (Behera and Shukla, 2015;Blanchet et al., 2017), which are recurrent phenomena in inner Mediterranean regions such as Extremadura (Schnabel, 1997).
The role of clay on soil quality in this land system can be controversial. It is obvious that clay contributes to a higher retention of cations (Sulieman et al., 2018) and some nutrients (Hosseini et al., 2017) due to its surface negative charge, although a high content of clay and silt can increase soil bulk density and consequently induce soil compaction (Bienes et al., 2016). Anyway, the affecting power of clay has been neglected by Pulido et al. (2017) due to the low range of soil texture studied at a farm scale (mostly sandy-loam texture).
Soils with high percentages of sand usually show high bulk density values because their soil density is naturally higher, and it does not mean soil compaction (Dexter, 2004). In addition, sandy soils in Extremadura are logically heavier because quartz is abundant in their mineralogy.
The positive correlation between clay and CEC has been also confirmed by other authors in arid and semi-arid environments (Koganti et al., 2017). It can be an indication of a positive effect of clay on fertility retaining cations and nutrients (although not for soil organic carbon and nitrogen) and reducing the acidity of soils (Silva et al., 2018). Nevertheless, the range of clay content has not been so large in this study.
The PCA also confirmed the logical negative effect of bulk density, considered as land degradation indicator by Pulido et al. (2017). These latter only found evidences of the negative effect of bulk density on pasture production (Pulido et al., 2018) although it still remains unclear if it is a consequence of a reduction in the fertility or in the capacity of storing water. Authors such as Montagu et al. (2001), Anderson et al. (2006), Kraaijvanger et al. (2015), and Bogunovic et al. (2017), among others, have already found evidences of a reduction in fertility induced by soil compaction.
In rangelands, soil compaction generally reduces fertility by two interconnected processes: animal trampling and the removal of tree and pasture litter accumulated on the soil surface (Pulido-Fernández et al., 2013). Under circumstances of overgrazing animal trampling exerts a pressure on soil comparable to heavy machinery that makes decrease its porosity and increases its penetration resistance limiting grasses root activity (Greenwood and McKenzie, 2001;Bilotta et al., 2007). It also increases runoff and soil erosion as well as the washing of nutrients (Bogunovic et al., 2020).
Although the lithology dominant in Extremenian rangelands are siliceous rocks (shales and granites) and the dominant land management is the livestock husbandry, we consider the necessity of combining in further research data taken at a regional scale with studies based on local knowledge from a farm scale. In addition, the values returned in this IQI should be compared with topographic factors (Manandhar and Odeh, 2014), soil maps (Stahr and Clemens, 2017;Nabiollahi et al., 2018), and information provided by the Spanish National Inventory of Soil Erosion, among other sources, to detect regional patterns of soil quality or particular influences. This IQI should be also validated in other areas with different dominant land uses such as vineyards, olive groves, or chestnut tree orchards.

CONCLUSIONS
We drew two regional maps (Extremadura, SW Spain) of soil quality by interpolating 194 punctual values of an integrated quality index (IQI) calculated both from two sets of weighting factors (PCA and AHP) and two types (less is better, much is better) of standard scoring functions (SSFs). Methodologically speaking, it supposed the transformation of 194 soil profile descriptions into valuable maps of soil quality. The most reliable map seemed to be that obtained from the AHP method since the weighting of indicators was more focused on the content of nutrients, which it is one of the limiting factors that explains the geographical distribution of rangelands. This study represents a methodological advance regarding previous studies aimed at a farm scale but further research is still needed in terms of spatialization of soil properties considering also topographical and managerial factors.