Spatial Distribution of Cadmium in Agricultural Soils of Eghlid County , South of Iran

a PhD in Environmental Science, Department of Environment, Faculty of Basic Sciences, Hamedan Branch, Islamic Azad University, Hamedan, Iran b Associate Professor, Department of Environment, Faculty of Basic Sciences, Hamedan Branch, Islamic Azad University, Hamedan, Iran c Assistant Professor, Department of Agriculture, Faculty of Agriculture and Natural Resources, Sayyed Jamaleddin Asadabadi University, Asadabad, Iran

Soil pollution by heavy metals is a significant worldwide problem due to the complex issues caused by this phenomenon (1,2). Irrigation with sewage, application of solid waste as fertilizer, agricultural activities, coal fossil fuel combustion, mining, waste disposal, vehicle exhaustion, and atmospheric deposition can lead to heavy metal contamination of the soil (3)(4)(5). Cadmium (Cd) is one of the most common heavy metals that pollute the environment (6). The Cd in the soil is naturally derived from chemical weathering of parent materials and in many areas, it enters water and soil resources as the result of various human activities, such as waste disposal, mineral extraction, and usage of chemical fertilizers and pesticides (7)(8)(9)(10).
In order to assess the impacts of heavy metals on the environment, first, it is needed to determine their spatial distribution in the soil (11). In the present study, spatial prediction methods (e.g., various Kriging methods) were used to generate spatial distribution maps which revealed the spatial characteristics or patterns of the soils polluted with heavy metals (12,13). Moreover, the adoption of geographic information systems (GIS) that use geostatistical and multivariate analyses is crucial for the assessment of the spatial structure of heavy metals and the physicochemical properties of soil (14).
Geostatistics is a good solution for the generalization of the results obtained from the measured points to other points (15). Geostatistics refers to the methods which make it possible to estimate or simulate the values of a variable at a specific location, based on the theory, regional variables, and the available information (16). Kriging is one of the methods in geostatistics with a high potential for the determination of the spatial distribution of heavy metals in soil (17).
Previous studies have examined the spatial distribution of heavy metal pollution in the surface soil in industrial regions across the globe. For instance, Wang et al. found that the spatial distributions of Cu, Zn, Cr, Cd, As, and also the Hg concentration in the industrial region in Sichuan, China were below the national standards; however, they were slightly higher than the natural baseline values (18).
Spatial distribution of Cd is generally known to be influenced by soil, geology, elevation, climate, vegetation, natural mineralization, and human activity. Such processes influence geochemical variables at different spatial scales, ranging from microscale mineral compositions to macroscale geochemical spheres (19). The zoning and mapping of soil Cd as well as identifying the sites contaminated can largely help the management and prevention of Cd distribution in the environment and living organisms. One of the main problems in the assessment of the contamination of different regions is that it is impossible and costly to draw samples from all points. Jiachun et al. (2006) in a study performed in China found that human intervention in nature was a decisive factor in spatial variations of Cd, Cu, Cr, Pb, and Zn (15). Juang et al. (2001) also reported that Cd was added to the soil by external factors and that 14.46% of the analyzed samples showed moderate to severe signs of Cd contamination. In addition, the spatial analysis indicated that various factories, including metal melting factories in the studied region, were the main sources of Cd that was found in the soil (7). Khodakarami et al. (2012) studied the zoning of Co, Cr, and Ni concentration using GIS and geostatistics and prepared the soil contamination maps for three sub-basins of Hamedan province (20).
Currently, due to the excessive use of chemical fertilizers in the agricultural lands of Eghlid in recent years, the preparation of soil contamination maps seems inevitable, especially in the areas where the potential of contamination is high. Obviously, it is essential to have an understanding of the spatial distribution of heavy metals in the area to help the environmental management and also reduce health-related issues in agricultural products. The present research aimed to investigate the spatial distribution of Cd and prepare the contamination zoning map in the agricultural lands of Eghlid, in Fars province, south of Iran, in order to control soil and crop contamination.

Study area and sample analysis
The soil samples were selected using the systematic random sampling method in the central part of Eghlid. Eghlid is among the northern cities of Fars Province, Iran, that is located between 52° 41′ 53'' eastern longitude and 30° 53′ 42'' northern latitude (21). The present research was performed on 100 soil surface samples that were taken at a 0-15 cm depth. Figure 1 shows the scatter plot of the sampling points. After the samples were airdried, soil aggregation analysis was performed on them using sieve No. 230 (manufactured in Gilson Co.). The particles smaller than 62 μm were separated, powdered in a porcelain mortar, and subsequently homogenized. The Cd was measured by the aqua regia method using (ICP-OES) (Varian 710-ES) (22). Moreover, the soil texture, organic matter, calcium carbonate equivalent, electrical conductivity (EC), and pH were measured in the samples (23).

GIS analysis
A global positioning system device (manufactured in Garmin Ltd.) was used to record the geographic features of the sampling points. Afterward, a GIS was installed, and using Spatial Analyst extension for ArcMap 10.4 GIS (ESRI, Inc., USA), spatial surfaces of Cd and soil properties were obtained using the interpolation method (ordinary Kriging). In the current study, the semi-variance function, as a basic geostatistical tool, and a key function were used to study soil variability. In this regard, the function included some of the major parameters, like nugget (C0), sill (C0+C), nugget effect C0/(C0+C), and range (A0). It can be used to reveal the spatial correlation of soil properties (24). The C0 stands for the change caused by random factors, (C0+C) indicates the total variation of the system, C0/(C0+C) signifies the spatial dependence of soil properties which is an important indicator for the measurement of the spatial variation of the regionalized random variables (25). When C0/(C0+C) < 25%, it means that there is a significant spatial correlation due to the structural factors. Moreover, when 25% < C0/(C0+C) < 75%, there is a medium spatial correlation as the result of both structural and random factors. In addition, when C0/(C0+C) > 75%, the spatial correlation is very weak and is mainly caused by random factors. It should be noted that a nugget effect of 1 indicates that there is a constant variation in all of the scales (24).
The standard sample-based cross-validation method was adopted to evaluate the accuracy of the spatial interpolation. According to the leave-one-out cross-validation method, the accuracy indicators were calculated based on the pairs of estimated-observed data of the sample points, i.e., the mean absolute error (MAE), the Pearson correlation coefficient (r), and the mean bias error (MBE). Strength of the linear relation and mean absolute deviation between the estimated and the observed data were measured using the r and MAE, respectively.
In order to have an accurate spatial interpolation, the value of r must be close to 1 while the value of MAE must be as small as it is possible. Therefore, smaller values of MAE indicated that the estimated values will be closer to the true values (12). The MBE signifies the mean bias between the estimated and true values; the sign of the MBE shows the underestimation (positive value) or overestimation (negative) of the interpolation results (26).

Statistical analysis
In the present study, all the analyses were performed in SPSS software (version 20.0). The mean value and standard deviations of the heavy metals were calculated for each sampling site. Moreover, the normality of the data was assessed by the Kolmogorov-Smirnov (K-S) test. Furthermore, an independent samples t-test was used to compare the concentration of Cd between the agricultural land and the control sites. Besides, Pearson's correlation coefficient was used to calculate the correlations between the contents of the elements among soil samples.
The statistical description of the Cd concentration and some soil properties of the studied soil samples are summarized in Table 1. As shown in Table 1, Cd concentration in the soil ranged from 0.12-7.38 mg kg -1 , with a mean value of 2.8 mg kg -1 . Comparison of the Cd concentration in the agricultural soil samples (0.11-0.25 mg kg -1 with a mean value of 0.16 mg kg -1 ) with pristine intact regions as blank control revealed that the concentration of Cd in control samples was significantly lower than the agricultural soil samples (P<0.01).
Based on the results, the pH rate was within the range of 7.21-8.99 with an average rate of 8.27 which indicated that the soil of this region was neutral and alkaline. As shown in Table 1, the maximum EC in the soil was 954 µs cm -1 . Moreover, the calcium carbonate values ranged from 45.23 to 89.25 which implied high levels of calcium carbonate in the study area. Importance of calcium carbonate is that it increases the soil pH and heavy metals have higher mobility and accessibility in acidic soils. Table 2 summarizes the correlation coefficient between Cd concentration and some soil properties. It must be noted that the Kolmogorov-Smirnov (K-S) test was used to assess the normality of the data distribution before running the correlation test. The descriptive statistics of the studied variables and the p-value of the K-S test are presented in Table 1. Moreover, Figure 2 shows the frequency percentage, normal distribution, Q-Q plot, and changes in the trend of Cd concentration in the studied soil samples.
The results of the K-S test indicated the normal distribution of the data which was also confirmed by the numerical value of the element skewness. Furthermore, the Q-Q plot diagrams showed that the data distribution was consistent with the normal distribution line, indicating the suitability of the data for the spatial analyses. The collected data also indicated the trend of changes in the numerical value of Cd in the samples.

Interpolation
The spatial analysis of the selected soil parameters is a prerequisite for designing effective strategies for the improved management of contaminated soil (27). For that purpose, the selection of the best-fitted variogram model to explore the spatial structure, and ultimately, generate accurate perdition maps at the un-sampled location is required. The variogram models for the soil parameters were estimated based on crossvalidation. In order to investigate the presence or absence of spatial correlation of Cd concentration and soil properties in the study area, according to the ordinary Kriging method, after fitting several semivariogram models, the most appropriate model was selected based on the residual sum of squares (RSS) and R 2 for Cd concentration, organic matter and calcium carbonate percentage, soil texture, pH and EC. The best-fit parameters of the variogram models are presented in Table 3.
Based on the above analyses, the maps of Cd concentration and soil studied properties were prepared based on ordinary Kriging. Figures 3  to 8 show the spatial distribution of Cd concentration, soil texture, organic matter, and calcium carbonate percentage, pH, and EC in the studied area.The contamination factor was used to determine soil contamination with heavy metals. To compute this factor, the content of heavy metals is compared with its natural level to determine soil contamination. In this study, the average concentration of Cd in the intact (untouched or pristine) natural lands Table 2

Figure 2) Frequency percentage and normal distribution (a), Q-Q plot (b), and changes in the trend of cadmium concentration (c) in the analyzed agricultural soil samples collected from Eghlid County
Proper soil management depends on understanding the conditions and distributions of pollutants and soil characteristics in the lands under cultivation. The spatial distribution maps that determine the areas with problems, deficiencies, or environmental contamination could be prepared through soil sampling and analysis. In this case, it is assumed that the soil characteristics measured at the point of interest are also representative of its surrounding sampling points and that the assumption of the spatial variability of soil characteristics is acceptable. Geostatistics and GIS provide the possibility of determining the spatial distribution of soil by having information on numerous soil samples.
In this research, Table 1 shows the statistical description of the Cd concentration and soil properties of the samples taken from the studied area. The coefficient of variation of Cd in the measured points is higher than 50%, indicating that there is a significant change in the Cd concentration in the studied soil samples. The value of this coefficient for Cd in the pristine area as the control area is less than 50 which shows small changes in the Cd concentration in these areas. These results indicate the increasing trend of changes in Cd in the studied area. Moreover, the results of the independent samples t-test (P<0.01) revealed that the agricultural lands and the pristine area had a significant difference regarding the concentration of Cd. These results also indicated that due to intensive cultivations, the use of pesticides and chemical fertilizers in the studied area was high which contains considerable amounts of impurities, such as heavy metals, like As, Cd, Pb, Ni, and Zn (8). This was in line with the results of the study performed by Cheraghi et al. (30). Texture of the soil plays a key role in binding or omitting metals; accordingly, the content of the clay affects the total absorbed heavy metal in the soil. Based on the results of previous studies, approximately, 96% of the content of a number of heavy metals can be detected in the clay fraction (31). The clay in soil has a significant effect on its cation exchange capacity and is considered as one of the sources of negative surface exchange in the soil (32). Besides clay, the other soil characteristics, including pH, organic matter, and carbonate calcium influence the amount of heavy metals found in the soil (33).
Based on Table 2, Cd concentration had a positive significant correlation coefficient with clay percentage and pH. According to the results of several previous laboratory research and studies, there is a negative correlation coefficient between the mobility of heavy metals and pH (34). However, in the present study, Cd concentration had an insignificant relationship with the content of the organic matter, percentage of calcium carbonate, and EC.
A high fraction of clay as a soil substance is a negatively charged colloid that is able to absorb higher amounts of Cd with an electrostatic attraction particularly on the surface of the clay. Clay has a significant effect on the soil chemistry since the surface features differ considerably from the mineral grains that have a larger size. Moreover, the critical limit of heavy metals in the soil has a direct relationship with the amount of clay; accordingly, greater amounts of clay in the soil indicate greater amounts of heavy metals. Sukarjo et al., in their study also found that linear clay content would increase the critical limit of Cd in the soil (32).
Due to the complex nature of the soil and the existence of various sources responsible for the distribution of the heavy metals in the soil, the amount of these contaminants normally undergo dramatic changes in the soil surface, and their distribution pattern generally has a spatial heterogeneity (35). However, identification of the spatial distribution of heavy metals in soil is required for designing the correct management plans to control these dangerous contaminants (36). In the present study, the Kriging method was used for interpolation.
After normalization of the data and selection of an appropriate covariate, geostatistical extension program in GIS was used to perform variogram analyses. The experimental semivariogram was calculated and subsequently, the best model (Table 3), which was selected according to the least value of the residual sum of squares, was fitted with experimental variograms. Analysis of the semi-variogram and evaluation of various methods of interpolation revealed that the semi-variogram of Cd, calcium carbonate, pH, and EC in the studied area and the organic matter followed a linear model and an exponential model, respectively.
Confidence level of all variograms was assessed based on the ratio of nugget variance to sill which is considered a criterion for the classification of the spatial dependence of the studied parameters. Ratio of nugget variance to sill indicated a strong spatial dependence, which was caused by structural factors, among the values of Cd, calcium carbonate, pH, and EC, and a moderate spatial dependence, due to both random and structural factors (37), for organic matters ( Table 3).
The mean square variogram shows the difference between the two variables as a function of the distance between them. Based on the fitted model, the variogram was presented in the mentioned figures. The accuracy of the fitted model was evaluated using R 2 , RSS, MAE, and MBE statistics. It should be noted that if the values of RSS, MAE, and MBE approach zero, it indicates the high accuracy of the estimated value.
The Kriging interpolation maps of the evaluated data of the studied region are shown in figures 3-8. As seen in the Kriging map, most of the area has high concentrations of Cd and is contaminated with it (Table 4). Spatial distribution pattern of soil characteristics revealed that the percentage of the clay in the north and central parts of the studied area was higher than the southeast, while pH and carbonate calcium rates were higher in the northeast and southeast regions, and higher rates of EC and organic matter are found in the northern part of the studied area.
Estimation of the amount of annual chemical fertilizers used in the agricultural lands of Eghlid showed that the central part of the city (the studied area) had the highest consumption rate of chemical fertilizers. Cultivated soils in the central part of the studied area were contaminated with Cd, having concentrations that exceed the standards for the quality of the agricultural soil (Table 4) (38).
Since soil fertility has decreased in the studied area, farmers tend to use more fertilizers. Consequently, such wrong farming practices have resulted in the concentration of metals in the soil (39).
The continuous use of these fertilizers makes the plant absorb and transfer heavy elements to the food chain (34). Martin et al. investigated the amount of heavy metals in the agricultural soils of a basin in Spain, using multivariate analysis and geostatistics and found that chromium and nickel were controlled by parent materials, while Cd, lead, mercury, copper, and zinc were affected by human activities (40). The average concentration of copper in the olive and grape fields was high, mostly owing to the use of chemical fertilizers and soil fertilizing materials in these areas.
In the present study, the contamination factor was used to assess the Cd contamination of soil based on the contamination factor and Hakanson classification (29) (Figure 9). According to the results, the studied area could be classified as an area with high Cd contamination.
The present pioneer study was performed on the spatial distribution of Cd and other soil properties of agricultural land in Eghlid County. Despite the fact that most of the studies area is contaminated with Cd, it seems that this contamination has a significant correlation with pH and clay content in the soil samples. According to the results of the spatial distribution patterns of Cd, it was found that the Cd concentration could be classified in the high contamination level, showing a very high spatial dependency for Cd. Regarding the low Cd concentration in the pristine areas, the results of this study clearly indicate that soil contamination with Cd in the studied area, compared with other factors, can mainly be attributed to the human factors, especially agricultural activities.