THE USE OF MULTIVARIATE STATISTICAL ANALYSIS AND SOIL QUALITY INDICES AS TOOLS TO BE INCLUDED IN REGIONAL MANAGEMENT PLANS. A CASE STUDY FROM THE MASHHAD PLAIN, IRAN

Several methods have been used to model reality and explain soil pedogenesis and evolution. However, there is a lack of information about which soil properties truly condition soil quality indicators and indices particularly at the pedon scale and at different soil depths to be used in land management planning. Thus, the main goals of this research were: i) to assess differences in soil properties (particle size, saturation point, bulk density, soil organic carbon, pH and electrical conductivity) at different soil depths (0-30 and 30-60 cm); ii) to check their statistical correlation with soil quality indicators (CEC, total N, Olsen-P, available K, exchangeable Na, calcium carbonate equivalent, Fe, Mn, Zn, and Cu); and, iii) to elaborate a soil quality index and maps for each soil layer. To achieve this, forty-eight soil samples were analysed in the laboratory and subjected to statistical analyses by ANOVA, Spearman Rank coefficients and Principal Component Analyses. Finally, a soil quality index was developed based on indicators of sensitivity. The study was conducted in a semiarid catchment in northeast Iran with irrigated farming and well-documented land degradation issues. We found that: i) organic carbon and bulk density were not similar in the topsoil and subsoil; ii) calcium carbonate and sand content conditioned organic carbon content and bulk density; iii) organic carbon showed the highest correlations with soil quality indicators; iv) particle size conditioned cation-exchange capacity; and, v) heavy metals such as Mn and Cu were highly correlated with organic carbon due to non-suitable agricultural practices. Based on the communality analysis to map of soil quality, CEC, Mn, Zn, and Cu had the highest weights (≥0.11) at both depths, coinciding with the same level of relevance in the multivariate analysis. Exchangeable Na, CaCO 3, and Fe had the lowest weights (≤0.1) and N, P, and K had intermediate weights (0.10.11). In general, the map of the soil quality index shows a lower soil quality in the subsoil increment than in the topsoil. Rodrigo-Comino et al. 688 Cuadernos de Investigación Geográfica 45 (2), 2019, pp. 687-708 El uso de análisis estadístico multivariante y de índices de calidad del suelo como herramientas que deben incluirse en los planes de gestión regional. El caso de la llanura de Mashhad, Irán RESUMEN. Muchos métodos permiten modelizar la realidad para tratar de explicar la formación y evolución del suelo. Sin embargo, a una escala suelo y a diferentes profundidades, no existe una gran cantidad de estudios que también mencionen los indicadores de calidad para ser incluidos en la ordenación del territorio. Por lo tanto, los principales objetivos de esta investigación son: i) evaluar las diferencias entre propiedades del suelo (tamaño de partículas, punto de saturación, densidad aparente, carbono orgánico, pH y conductividad eléctrica) a diferentes profundidades (0-30 cm y 30-60 cm); ii) su posible correlación estadística con indicadores de calidad de suelo (capacidad de intercambio catiónico, N total, Olsen-P, K disponible, Na intercambiable, carbonato cálcico equivalente, Fe, Mn, Zn y Cu); y, iii) elaborar un índice de calidad del suelo con su respectivo mapa para cada capa estudiada. Para alcanzar dicho objetivo, 48 muestras fueron recogidas y analizadas en el laboratorio, comparadas con diferentes métodos estadísticos (ANOVA, correlación de Spearman y análisis de componentes principales) y utilizadas para la creación de un índice de calidad del suelo en una cuenca de un ambiente semiárido al Noreste de Irán, caracterizada por una agricultura por irrigación y documentados procesos de degradación. Las principales conclusiones alcanzadas fueron: i) el carbono orgánico y la densidad aparente fueron diferentes para cada profundidad; ii) el carbonato cálcico y el contenido en arenas estaba altamente correlacionado con el carbono orgánico y la densidad aparente; iii) el contenido en carbono orgánico mostró la mejor correlación con los indicadores de calidad del suelo; iv) el tamaño de las partículas del suelo condicionó la capacidad de intercambio catiónico; y, v) los metales pesados como el Mn y el Cu estuvieron fuertemente correlacionados con el carbono orgánico debido, posiblemente, a las malas prácticas agrícolas. Basados en el análisis de comunalidad del índice de calidad del suelo, para ambas profundidades, la capacidad de intercambio catiónico, Mn, Zn y Cu obtuvieron los mayores pesos, coincidiendo con su importancia dentro del análisis multivariante realizado anteriormente. Por el contrario, el Na intercambiable, CaCo 3 y el Fe obtuvieron los valores más bajos. N, P y K alcanzaron valores intermedios. En general, el mapa de calidad del suelo mostró valores menores en el horizonte sub-superficial que en la capa superficial.


Introduction
Land degradation affects soil quality in many semiarid areas and is one of the most important concerns of humankind (Armenise et al., 2013;Qi et al., 2009). Conservation strategies designed by policymakers and farmers should be structured to protect soil quality and fertility (Calleja-Cervantes et al., 2015). Specifically, soil erosion and contamination are two of the most important causes of soil quality degradation (García-Ruiz et al., 2015;. In natural semiarid landscapes, high temperatures and evapotranspiration rates and extreme rains concentrated in a few events are the most important factors driving soil quality (Zeraatpishe and Khormali, 2012). The generation of soil crusts (Singer and Shainberg, 2004) or the loss of fine particles such as silt or clays (Bruand and Tessier, 2000;Torres-Sallan et al., 2017) are clear signals of decreasing soil quality. In anthropogenic areas such as cultivated valleys, researchers also have to pay attention to other variables such as the intensive use of organic amendments (Fernández-Calviño et al., 2017), irrigation systems (García-Garizábal et al., 2017) and land use changes such as grazing or abandonment (Gutzler et al., 2015;Khaledian et al., 2017a). The establishment of clear soil degradation signals such as soil compaction, increase in heavy metals, loss of soil fertility and topsoil depletion can allow us to quantify degradation (Mwango et al., 2015;Pulido et al., 2017).
Several recent studies have focused on detecting soil parameters that can be used as soil quality indicators in topsoil (Liu et al., 2017). Results have demonstrated that some well-known chemical and physical properties have significant correlations. For example, soil texture has a high correlation with cation-exchange capacity (CEC), aggregate stability and soil water content (Khaledian et al., 2017a;Mamedov et al., 2017); soil organic carbon (SOC) is also well-correlated with aggregate stability, soil nutrients and biodiversity (Acín-Carrera et al., 2013;Wang et al., 2016); and bulk density (BD) with electrical conductivity (EC), clay content and soil porosity (Nasri et al., 2015). Heavy metals such as Fe, Mn, Zn and Cu are widely used as soil quality indicators because of their effect on soil fertility and corresponding plant responses (Pezeshki and DeLaune, 2012;Sparrow and Uren, 2014;Tao et al., 2007).
However, research into soil properties and soil quality indicators is not common in semiarid landscapes in developing countries with a combination of traditional land uses such as grazing and irrigated farming. In Iran, research has focused on the development of models and pedotransfer functions to clarify these kinds of relationships (Hosseini et al., 2017;Keshavarzi et al., 2018;Shiri et al., 2017). In addition, there is a lack of studies focused on the variation in and interpretation of soil properties at different depths and how they could influence soil quality indicators. Thus, the first goal of this research was to evaluate differences in soil properties at different soil depths and their statistical correlations with soil quality indicators. Moreover, the research aimed to discover which group of factors can explain the variance and weight of soil properties and quality indicators at the pedon scale in a semi-arid catchment. Thus, we developed a soil quality index and soil quality maps. To achieve this goal, a semiarid catchment was selected in northeast Iran.

Study area
The study area is in Khorasan-e-Razavi Province, northeast Iran, specifically in the Mashhad Plain between latitudes 35°59′ N to 37°04′ N and longitudes 58°22′ E to 60°07′ E (Fig. 1). Soil samples were collected from sites with a mean elevation of 1200 m a.s.l. The highest altitudes in nearby areas are up to 2700 m a.s.l. The geological units are shown in Figure 2. Volcanic rocks such as andesites are located in the northwest. Along the west and north parts, limestones and dolomites are the most common rock units. Red marls, gypsiferous marls, sandstones and conglomerates are found in between. The second most extensive geological unit is the fluvial valley characterized by conglomerates and piedmont. Granites and metamorphic rocks such as schists, phyllites, slates and meta-sands can be found in the southern part of the study area. The surrounding landscapes are composed of pediment plains, upper terraces, gravelly colluvial fans, and plateaus. The climate is semi-arid, with rainfall averaging about 222.1 mm yr -1 and temperatures about 15.8°C (Keshavarzi et al., 2016). The soils are classified as Calcaric Cambisols, Gypsic Regosols, Calcaric Regosols and Calcaric Fluvisols (WRB Working Group, 2014) with horizons that are high in calcium carbonate content and low in soil organic carbon. The main land uses in the study area are also shown in Figure 2. Irrigated farming and a mixture of different crops are located along the valley. Rangelands are the most common land use in the mountainous area. Dryland farming practices are carried out in the eastern part of the study area and at specific points along the water courses. Low to very low forest areas are scarce.

Soil sampling and laboratory analysis
Forty-eight representative soil profiles were selected using a stratified random sampling technique, similar to the procedures explained by Arnab (2017) and Christofides (2005). Soil samples were taken from two different depth intervals, 0-30 cm and 30-60 cm. We decided to follow this separation instead of soil horizons because our goal was to focus on soil quality in the topsoil, or the part of the soil influenced by crop roots and water infiltration. Large plant materials and pebbles in the samples were separated by hand and discarded. Soil samples were air-dried, homogenized, sieved through 2 mm mesh and stored in polyethene bags at ambient temperature. Soil organic carbon (SOC) content was determined by the Walkley and Black method (Walkley and Black, 1934) with dichromate extraction and titrimetric quantization. Soil particle size (clay < 0.002 mm, silt=0.002-0.05 mm, and sand=0.05-2 mm) was analyzed by the hydrometer method (Gee and Bauder, 1986). Soil water content based on the degree of saturation (SP) was also quantified following Blake and Hartge (1986) and Gardner (1986). The clod method (Blake and Hartge, 1986) was used to determine bulk density (BD) with three replications per sample site and depth interval. Electrical conductivity (EC) and pH values were measured in saturated paste extract using a digital EC-pH meter.
The exchangeable sodium (Na) was measured using the method of Harron et al. (1983). CEC was analyzed according to Bower et al. (1952). Available phosphorus (P) and available potassium (K) were measured using standard procedures (Sparks et al., 1996). The calcium carbonate equivalent (CCE) was obtained by the back-titration technique (Nelson, 1982). Nitric acid-extractable forms of Zn, Cu, Mn and Fe were determined after extraction with 4 M HNO 3 (1:8 mass ratio of soil to final mixture) at 80 º C for 16 h. After cooling, samples were filtered through coarse acid washed cellulose filters. Ten g of soil was extracted with 20 ml of 0.005 M DTPA (diethylenetriaminepentaacetic acid) + 0.01 M TEA + 0.01 M CaCl 2 , adjusted to pH 7.3 and shaken for 2 h to measure Cu, Mn, Fe, and Zn DTPA forms (Lindsay and Norvell, 1978). In both cases, metal concentrations were determined using an atomic absorption system (Shimadzu model AA.670).

Multivariate analysis
Soil parameters were divided into two different groups: physical and chemical soil properties such as soil particle size, SP, BD, OC, pH, and EC) and soil quality indicators such as soil nutrients and heavy metals. We included CEC in the soil quality indicator group because it is also considered an indicative property of minerals and nutrients in relation to soil texture and organic matter (Sulieman et al., 2018). The statistical parameters average, median, standard deviations and outliers of soil particle size (clay, silt, and sand), SP, BD, OC, pH, and EC were calculated and depicted in box plots to describe variability in the soil profile and the catchment (Fig. 3). We applied the One-Way ANOVA to test statistical differences between soil properties at the different depths when the properties were normality distributed as determined by the Shapiro-Wilk test and variances were equal. The Kruskal-Wallis One Way Analysis of Variance on Ranks test was used to determine statistical differences between soil properties when the Shapiro-Wilk normality test failed but variances were equal. The Tukey test was used to test statistical differences between soil properties when results were not normally distributed and variances were not equal. SigmaPlot v. 12.0 (Systat Software, Inc.) was used for these statistical analyses. Since not all the data followed a normal distribution, a Spearman rank coefficient was conducted to observe the correlation between variables.
To assess the soil properties that influenced soil quality indicators, a principal component analysis (PCA) was performed. Interpreting their influence on each factor, we were able to distinguish how important each factor was at the pedon and catchment scales at different soil depths. We selected a total of 6 components. Prior to conducting the PCA test, variables were centred and normalized in order to assign equal weights to all the variables and not misleading directions. To verify these assumptions, we applied a Kaiser-Mayer-Olkin (KMO) test of sphericity and calculated the coefficients and determinants. A PCA with an orthogonal rotation method (Varimax) and a correlation matrix, using factors with eigenvalues >1 was used.

Soil quality index and mapping
The development of a soil quality index and map usually includes three steps. First, appropriate indicators must be chosen. Second, the indicators are scored (0-1) in order to make a classification. A linear scoring method was used since the indicators are expressed with different numerical scales. Finally, the obtained indicator scores are combined into a unique and specific index, which can be mapped. The index in this study was based on an indicator sensitivity, where three different functions were assigned: 1) a "more is better" function; 2) a "less is better" function; and, 3) an "optimal range" function (Andrews et al., 2002;Mukherjee and Lal, 2014;Qi et al., 2009). To obtain the soil quality index (SQI), soil quality indicators have to be previously weighted. We assigned the weight value of each indicator using the communality value through the factor analysis that we elaborated previously. The equation used to develop the SQI was (Equation 1): (1) where SQI is the soil quality index; W i is the weight of each indicator; and, N i represents the score of each indicator. In the current research, ten different soil quality indicators (CEC, total N, Olsen-P, available K, exchangeable Na, calcium carbonate equivalent, Fe, Mn, Zn, and Cu) were used in the soil quality assessment and mapping. The soil quality map was generated using ordinary kriging in ArcMap 10.5 (ESRI, USA).

Soil properties and soil quality indicators at different depths
Soil properties such as particle size, SP, BD, OC, pH values and EC at different depths (0-30 and 30-60 cm) were depicted in box plots to observe the average, median, maximum and minimum values (Fig. 3). Statistical differences between soil properties at different depths were also calculated (Table 1). At 0-30 cm, particle size was characterized by a loam texture (32.1%, 49.6% and 18.3% for sand, silt, and clay, respectively). Similar average values were also observed at the 30-60 cm depth. SP (from 36 to 38%), pH (8.1), and EC (from 1.9 to 2.3 dS m -1 ) values were also similar at the two depths and did not show statistically significant differences. OC and BD also did show significant differences (Table 1) at the 0-30 and 30-60 cm depths. Average OC levels were low (0.6% and 0.4% at the 0-30 and 30-60 cm depths, respectively). Average BD values were 1.39 g/cm -3 in the topsoil and 1.55 g/cm -3 in the subsoil.  The values and ranges of soil quality indicators (CEC, total N, Olsen-P, available K and exchangeable Na, CCE, Fe, Mn, Zn and Cu) in the box plots (Fig. 4) and statistical comparisons are shown in Table 1. CEC did not show statistical differences with depth. Soil nutrients (total N, Olsen-P, and available K) showed statistical differences at the different depths, with the highest values in the 0-30 cm layer. CCE and exchangeable Na were not statistically different at the two depths; however, values for both were lower in the 0-30 cm depth than the 30-60 cm depth. Fe and Cu contents in the 0-30 and 30-60 cm depth intervals were also not significantly different. On the other hand, Mn and Zn showed significant statistical differences with both of them having higher contents in the 0-30 cm layer.

Correlations among variables at different soil depths
Observing the relationships between selected soil properties and soil quality indicators allowed us to distinguish which parameters played important roles in the PCA. As almost none of the variables followed a normal distribution, a Spearman rank coefficient was conducted (Sp rc ). Positive and negative results higher than 0.45 were considered to be indicative of correlation.
Results of the Spearman analysis for the 0 to 30 cm interval are shown in Table 2. The highest correlations were found between sand and silt contents and CEC. OC had the highest number of correlations with soil quality indicators (four), those being total N (Sp rc > 0.80), Olsen-P (Sp rc > 0.55), Mn (Sp rc > 0.57), and Cu (Sp rc > 0.57). Other relevant correlations between soil properties and soil indicators were: EC with exchangeable Na, total N with Olsen-P and Mn with Cu. At the 30-60 cm depth sand and silt again registered high correlations with CEC, and OC showed several correlations with the quality indicators total N, Fe, Mn and Cu (Table 3). Other strong correlations were found between EC and exchangeable Na, Fe and Mn and Mn and Cu (Table 3).

Principal component analysis
The total percentage of variance was added after conducting the principal component analysis (PCA) ( Table 4). It was noted that six factors were able to explain a large amount of the variability (82.1%) at the 0-30 cm depth; three components accounted for almost half (49.5%). For the 30 to 60 cm depth, six components were able to explain 80.1% of the variability in the data; the first three components explained a total of 51.6%. At the 0-30 cm depth interval, the first component (accounting for 20.1% of the variability in soil analysis) was mainly determined by parameters related to soil particle sizes such as sand, clay, BD, OC and SP (Table 5). The second component illustrated the positive correlation between silt and CCE, and the negative correlation with sand (-0.78). The third component was composed of heavy metals such as Mn and Cu with CCE. The fourth component showed the positive correlations between EC and exchangeable Na and a negative correlation with pH values. The fifth component grouped positive correlations of OC and soil nutrients (total N and Olsen-P). The last component included increases in available K and decreases in Zn. Grey colours represent the selected variables.

Soil quality indices
The results of the paired-samples using the T-test for both the calculated soil quality indexes and the trend of soil quality index for topsoil and subsoil data are shown in Table 7 and Figure 5. The results showed that there was a significant difference (< 0.01), with average soil quality index values of 35.76 and 31.62 in the 0-30 and 30-60 cm depth, respectively.  To generate the soil quality index for the different soil depths, each soil quality indicator was assessed and weighted (Tables 8 and 9). Based on the communality analysis carried out before, CEC, Mn, Zn, and Cu had the highest weights (≥ 0.11) for both depths while exchangeable Na, CaCO 3 , and Fe showed the lowest weights (≤0.1). N, P, K, and Cu had intermediate weights (0.1-0.11).  Two maps of the soil quality indexes were generated (Fig. 6). Coefficients near zero indicate very low soil quality while results close to 100 represent areas with high quality as indicated by our selected soil quality indicator criteria. Higher soil quality was observed over larger areas in the 0-30 cm interval than at the 30-60 cm interval. The highest 0-30 cm soil quality index (>48) was observed in the northwest part of the study area, while the lowest index (<28) was in the middle part. A drastic decrease in soil quality was noted from the middle to the southern parts of the study area at 30-60 cm. Only the northwest and northeast areas can be considered to be at an intermediate level of soil quality.

Discussion
The identification of soil properties that control soil quality indicators is considered vital to plan appropriate soil management (Calleja-Cervantes et al., 2015). In semiarid areas with deep valleys that are cultivated and irrigated, such as in northeast Iran, appropriate soil and water management techniques are needed as important tools to conserve soil fertility and productivity (Samani et al., 2018). Assessing the relationships between soil properties and soil quality indicators, as was done in this research, will allow us to refine our interpretations of pedogenesis and soil horizon evolution due to a range of events from simple alterations of soil properties to degradation processes (Brimhall et al., 1991;Qadir et al., 2013;Sulieman et al., 2016).
In our study area, OC and BD were the soil properties with the highest variability at the pedon scale. In semiarid areas with strong human impacts due to cultivation, livestock or abandonment, the distribution of soil properties changes across the landscape and with depth (Acín-Carrera et al., 2013;Pulido et al., 2017). In most natural soils, OC is higher in the topsoil layers and decreases with depth, but in agricultural areas, the spatiotemporal variability can be very high and can affect plant vigour (Novara et al., 2018). Therefore, human impacts can be considered a factor driving these soil property changes. In the Mashhad Plain, the higher OC seen in the topsoil layer may be due to the use of amendments for the irrigated fields (García-Garizábal et al., 2017). Increases in BD have been shown to reduce levels of soil organic matter (Brevik and Fenton, 2012), because it reduces soil aeration, root penetration, and makes it more difficult for organisms to move and conduct typical activities (Guo et al., 2016). Calcium carbonate content and textural changes were driving factors in the BD differences. The PCA demonstrated a high positive correlation between BD, clay content, and sand. We concluded that soil particle size and carbonate content contributed to significant statistical differences between soil depths. It is also well known that BD increases with increasing depth, even without changes in texture or carbonate content, because the weight of the overlying soil compresses the deeper, underlying soil (Al-Shammary et al., 2018).
Cation exchange capacity is a useful soil quality indicator that represents the soil's ability to hold positively charged ions. Correlation analyses, PCA and the soil quality index showed that at 0-30 and 30-60 cm there was a strong correlation between CEC and soil particle size, specifically, a negative relationship with sand and positive correlation with silt and high soil quality. A strong correlation with clay would typically be expected (Mamedov et al., 2017;Sulieman et al., 2018). However, Khaledian et al. (2017a) also found that CEC was not strongly correlated with clay in Spain when they modelled CEC for different countries. Khaledian et al. (2017a) hypothesized that the lack of correlation between clay and CEC in the Spanish soils might be due to the long-term intensive agricultural management that had led to the creation of highly degraded soils. That situation is highly similar to our experimental area. However, in our study area, the physiographic position on the alluvial plain with conglomerates and fine materials (Fig. 2) possibly means the soils developed in a thick sediment-dominated environment belonging to the Quaternary period, and this also may have affected this correlation (Keshavarzi et al., 2016(Keshavarzi et al., , 2018. Organic carbon had the highest number of correlated soil quality indicators; two of them were heavy metals (Mn and Cu). Pezeshki and DeLaune (2012) concluded that Mn is an essential plant micronutrient that is well correlated with OC and total N and P. They also observed that a decrease in soil tillage or increase in soil humidity could affect their levels. However, Sparrow and Uren (2014) indicated that high levels of Mn could damage crops. Thus, appropriate management of Mn is essential, with organic amendments often producing better results than chemical products (Tao et al., 2007). We also observed that increases in some soil nutrients such as Olsen-P and available K or exchangeable Na corresponded to increased levels of Fe, Cu, and Zn. High concentrations of heavy metals in soils are considered dangerous for plant, animal, and human health (Rinot et al., 2018) and increasing their levels is not typically recommended (Chen et al., 2015;Fernández-Calviño et al., 2017), although Fe, Cu, and Zn are all essential human nutrients and do need to be present in small amounts (Steffan et al., 2018). We observed a positive correlation at the surface (0-30 cm) between organic matter, Fe, Mn and Cu; other researchers have also reported positive correlations between organic matter and metals (Gao and Chen, 2012;Kargar et al., 2015;Khaledian et al., 2017b). Further studies in this area could be focused on the relationships between organic matter, heavy metals and exchangeable Na, and soil fertility by applying pedotransfer functions or fractal analysis (Hosseini et al., 2017;Keshavarzi et al., 2018).
Our maps show that soil quality was higher in the 0-30 cm layer than at 30-60 cm. We observed that the areas of lowest soil quality extended from the middle part of the study area to the south, coinciding with flat and cultivated areas. Soil can serve as a filter that protects groundwater from pollutants, but if not well managed, the utility of soil as a filter can be overwhelmed or eliminated (Keesstra et al., 2012). In fact, in some cases, soil can actually become a source of water contamination (Steffan et al., 2018). At this level, possible nature-based solutions must be evaluated, by using concepts such as connectivity (López-Vicente et al., 2013) that allow us to think of suitable management as a guide to evaluate control strategies. To achieve this goal, new advances in mapping the spatial distribution of soil properties and land use changes can contribute to improvements in soil management (Kharazmi et al., 2018). This coincides with the main goal of our research, which was to use knowledge of the interrelationships between soil parameters to allow predictions of their spatial variability at the pedon scale.

Conclusions
We demonstrated that the statistical assessment of soil properties and soil quality indicators can be a useful tool to detect the factors that condition soil evolution and fertility. In the Mashhad Plain of northeast Iran, we found significant differences in select soil properties between 0-30 and 30-60 cm depth intervals. The multivariate analysis demonstrated that OC, bulk density, and CEC showed the highest correlations with the other soil physicochemical properties, showing their potential as soil characterization predictors. A soil quality index was developed using CEC, Mn, Zn, and Cu as soil quality indicators. Soil quality maps showed that heavy metals were key factors causing soil degradation in the 30-60 cm depth interval. These kinds of relationships are influenced by several factors, which in turn has the potential to affect soil and water quality and fertility. This must be considered in land management plans. Therefore, one direction of future research should focus on detecting the spatial variation of potential contaminants across the Mashhad Plain, paying particular attention to hotspots with high levels of soil contamination.