Spatial variability of eight hazardous heavy metals in agricultural topsoils: a case study from Xuzhou area, Jiangsu, China

ABSTRACT

Understanding the concentration of agricultural topsoil heavy metals and their spatial variability is crucial in terms of keeping a close track of the quality of regional environment. In this study, a total of 139 sites were sampled from the agricultural topsoils in Xuzhou area, Jiangsu, China. The contents of eight hazardous heavy metals (Cu, Zn, Ni, Pb, Cr, Hg, As, Cd) were measured for each soil sample. We used semivariogram to study the spatial patterns of these heavy metals, and obtained their spatial distribution maps in Xuzhou area. The descriptive statistics and the correlation between each heavy metal were examined. Our data show that the concentrations of all eight hazardous heavy metals do not exceed the threshold values that are recommended in the national soil quality standard (GB15618–2018). These new results provide a useful basis for the policy design regarding the agriculture, industry, topsoil contamination treatment, and regional environment in this region.

KEYWORDS:

• Agricultural topsoils
• heavy metals
• spatial distribution
• spatial variability
• Xuzhou
• China

1. Introduction

The quality of the agricultural topsoils is critical for policy design in order to use the land properly, control the topsoil contamination level, and establish an environmentally friendly ecosystem in the urban area. Heavy metals in agricultural topsoils are very useful in order to reflect the contamination level of regional environment, and the most commonly studied hazardous heavy metals are Cu, Zn, Ni, Pb, Cr, Hg, Ad, and Cd [1–5]. To evaluate the agricultural topsoil quality in the right manner, it is important to consider not only the concentration but also the spatial distribution of the heavy metals [6,7]. One of the difficulties associated with the provenance analysis of the heavy metals in the topsoils is that often time there is more than one source of pollution, such as emission from vehicles, pollution from coal burning, and wasters generated by factories, etc. Moreover, any contamination of agricultural topsoils could potentially lead to groundwater contamination, which, in turn, contaminates an even larger range of the topsoils since the heavy metals are highly mobile [3, 8,9]. Therefore, topsoil heavy metal contents are of great importance regarding the safety of the farmland products, and the agricultural economy [10,11].

Geostatistical analysis has been routinely used in the spatial mapping of topsoil heavy metals [12,13,14,15,16]. For instance, the semivariogram can reveal the spatial variations of the sample values in the background of the regional values. The Kriging interpolation can reflect the spatial correlation in the function of distance or direction, which provides an explanation for the spatial variations on a regional scale. However, it is worthwhile to point out that even though the geostatistics has been successfully applied to many studies, it is still debated how the sampling resolution could affect the final conclusion of the studies [17].

The industry of Xuzhou area heavily concentrates on coal-mining. However, the spatial distribution of the agricultural topsoil heavy metals has not gained much attention [18]. In this study, the spatial distribution and variability of eight hazardous heavy metals (Cu, Zn, Ni, Pb, Cr, Hg, As, Cd) was studied using techniques including geostatistics and geographical information system (GIS) in order to investigate the spatial distribution of these heavy metals and hence, offer valuable suggestions for the regional soil quality management.

2. Materials and methods

2.1. Study area overview

This study focuses on the agricultural topsoils in Xuzhou area (GPS: 34.2°N, 117.3°E), which lies in the northwestern corner of the Jiangsu Province, China (Figure 1). Xuzhou is well-known as one of the major coal-mining cities in China, with a total area of 11,258 km². The landform of Xuzhou is overall flat but becomes hilly in the central and east areas. The climate is characterized by the temperate continental monsoon, and the mean annual precipitation and temperature are 900 mm and 14°C, respectively.

Figure 1. (a) The map of China, with the red star showing the location of the Xuzhou area. (b) Sample locations in the Xuzhou area. (c) An example of the distribution of samples with the background of a satellite image

2.2. Chemical analyses

Sampling sites were designed in order to be spatially representative for the agricultural topsoils in Xuzhou area, and were also distributed in a way to help employ geostatistical analysis (Figure 1). A total of 139 samples were collected from the topsoils within a depth of 0–15 cm, of which 89 samples were collected from 19. Sample locations are shown in Figure 1. At each sampling site, 3 subsamples of surface soil were collected and mixed thoroughly to obtain an approximate 1 kg soil sample. The soil specimens were carefully cleaned by a wooden spatula in order to remove any potential contamination by metal, plant remains, and roots, etc. Foreign rock pieces were also picked out. Then, samples were air-dried at room temperature and sieved using a 150 µm size sieve.

The sieved soil samples were digested by HCl–HNO₃–HClO₄–HF and the concentrations of Cu, Zn, Ni, Pb, Cr, Cd were determined by an Agilent 7700 Series inductively coupled plasma mass spectrometry (ICP–MS). Hg and As concentrations were measured by an atomic fluorescence spectrometer (AFS-8330, Jitian Analytical Instrument Co., Beijing, China) after the digestion of HNO₃–HCl for Hg analysis and HCl–CH₄N₂S for As. The standard reference materials (GSS–3 and GSS–13) were used in all analyses to control the accuracy of measurement with recoveries of 88.4–115.2% for 8 heavy metals. Raw data of the heavy metal concentrations can be found in supplementary Table S1.

2.3. Statistical analyses

Conventional statistical methods and histograms are only descriptive, which are unable to take the spatial distribution of the data into consideration. The four commonly used linear geostatistical methods are, specifically, the variances of estimation, dispersion, regularization, and kriging. However, these methods only incorporate the structural function of the random function, i.e. covariance or variogram. On the opposite, the geostatistical analysis usually starts with constructing a model in order to characterize the spatial structure of the regionalized variables, such as the heavy metals in the topsoils. The prerequisite to the application of geostatistics is that the presence of a spatial structure where observations close to each other are more similar than those that are far apart (spatial autocorrelation) [20]. The average degree of dissimilarity between a sampled value and a nearby sample value can be determined by experimental semivariogram [21]. Therefore, it can give autocorrelation at various distances.

Geostatistics uses the method called semivariogram to obtain the spatial variability of a regionalized variable, and provides the input parameters for the spatial interpolation of kriging [22]. The value of the experimental semivariogram for a separating distance of h (referred to as the lag distance) is half of the average squared difference between the value at Z(x_i) and at Z(x_i+h),

(1) $\gamma \left(h\right)=\frac{1}{2}Var\left[Z\left(x_i+h\right)-Z\left(x_i\right)\right]$(1)

(2) $\gamma \left(h\right)=\frac{1}{2N\left(h\right)}\sum _{1=0}^{N\left(h\right)}{\left[Z\left(x_i+h\right)-Z\left(x_i\right)\right]}^2$(2)

where N(h) is the total number of paired sample points separated by the lag distance h. For the sampling at an irregular interval, the distances between the sample pairs are unlikely to be exactly equal to h. Hence, the lag distance h is often represented by a distance band. Z(x_i) is the measured sample value at the point i, and Z(x_i + h) is measured sample value at point i + h (Equations (1)(1) $\gamma \left(h\right)=\frac{1}{2}Var\left[Z\left(x_i+h\right)-Z\left(x_i\right)\right]$(1) and (2(2) $\gamma \left(h\right)=\frac{1}{2N\left(h\right)}\sum _{1=0}^{N\left(h\right)}{\left[Z\left(x_i+h\right)-Z\left(x_i\right)\right]}^2$(2) )). The variogram model is chosen out of a set of functions (e.g. linear, gaussian, spherical, exponential) that could best describe the spatial relationships [23]. We determine a model as appropriate if the chosen model can best-fit the shape of the curve of the experimental variogram to the shape of the curve of the mathematical function [23]. The fitted model provides information about the spatial structure as well as the variogram parameters (range, sill, and nugget) for kriging interpolation. By using the appropriate variogram model, the distance-dependent coefficients can also be estimated and graphically interpreted [23]. Raw data were analyzed with the help of various software packages. The descriptive statistical parameters were calculated using SPSS statistics software. The geostatistical analyses and the probability calculation were carried out using VarioWin software. Semivariograms were made using GS+ software. The maps of eight heavy metal elements spatial distributing pattern were generated using ArcGIS software by kriging interpolation.

3. Results and discussion

3.1. Agricultural topsoil heavy metal concentrations and correlations

Raw data were evaluated and the descriptive statistical parameters of eight topsoil heavy metals are summarized in Table 1. Based on the results, none of eight heavy metals pass the Kolmogorov–Smirnov (K–S) test of normality [24] Table 1. However, in order to perform geostatistical analysis, the normal distribution is required. Therefore, the data were transformed to follow a lognormal distribution first.

Table 1. Summary of statistical parameters of topsoil heavy metals (unit: mg/kg)

Element	Cu	Zn	Ni	Pb	Cr	Hg	As	Cd
Minimum	2.74	37.22	15.45	4.98	15.12	0.01	3.73	0.00
Maximum	245.9	398.30	57.12	125.46	229.65	0.20	18.80	0.74
Median	18.40	64.41	28.12	21.43	57.45	0.09	8.80	0.10
Mean	22.38	71.29	30.15	24.65	57.60	0.08	9.39	0.13
S.D.	21.73	36.62	7.60	13.01	20.56	0.05	2.85	0.11
C.V.	0.97	0.51	0.25	0.53	0.36	0.59	0.30	0.83
Skewness	8.81	6.23	1.19	4.46	4.38	0.04	1.13	3.23
Kurtosis	90.23	51.84	1.46	28.69	35.45	1.26	1.98	12.80
K–S test	0.24	0.20	0.13	0.18	0.14	0.13	0.13	0.20
K–S p	0.00	0.00	0.023	0.00	0.00	0.02	0.012	0.00
n	127	127	127	137	137	137	137	135
Distribution	lognormal	lognormal	lognormal	lognormal	lognormal	lognormal	lognormal	lognormal
Xuzhou background	12.61	91.10	32.80	16.30	55.50	0.03	11.20	0.29
China Mean	20.00	69.96	25.04	18.80	57.77	0.08	6.16	0.15
GB15618–2018	100	300	190	170	250	1.0	20	0.6

p < 0.05, S.D. = standard deviation, C.V. = coefficient of variation.

Results show that all eight heavy metals exhibit a large range of values within the study area. For example, the maximum values of each heavy metal are usually a couple of, sometimes even more than ten times of the minimum values (Table 1). The median and mean values are close to each other. However, when compared with the regional background or the mean values of China, several heavy metals show high concentrations in the topsoils. Specifically, Cu, Zn, Ni, Pb, Hg, and As show mean values higher than China mean, which reflects that in the Xuzhou region, these heavy metals are potentially pollution sources, and can be used for the topsoil metal contamination assessment. In contrast, Cr and Cd have values lower than China mean, indicating these two metals are not the major concern in terms of the topsoil contamination in Xuzhou area. Considering that Xuzhou is a city famous for its production of coal, as well as its large amount of heavy industries, these heavy metals could be originated from a combination of these activates. It is also possible that the high contents of heavy metals in topsoils are caused by anthropogenic activity, such as the use of fertilizers and pesticides, vehicle emissions, industrial wastes, and pollution from coal burning especially in the winter. Based on the Standard of Chinese Environmental Quality for Soils (GB 15618–2018), none of the concentrations of these topsoil heavy metal is higher than the threshold values recommended for cultivated areas (Table 1), indicating that soil environmental quality in the study area was of little threat in terms of both human health, environment, and agricultural activities [25].

All eight heavy metals show large values of their standard deviation (S.D.). The coefficient variation (C.V.) values of eight heavy metals in the study range from 0.25 to 0.97, suggesting that they vary moderately across the study area (Table 1). For instance, Cu has the highest C.V. value of 0.97, indicating that Cu is very likely to be influenced by extrinsic factors such as human activities. Cd is similar to Cu. On the opposite, Ni, Cr, and As have low C.V. values, meaning that these three metals are affected by intrinsic factors and have minimal variation. Their contents are nearly constant across the study area (Table 1). The rest of the heavy metals such as Zn, Pb, and Hg have C.V. values around 0.5, indicating that they have both extrinsic and intrinsic effects (Table 1).

We evaluated how these heavy metals correlate. The Pearson correlation coefficients and their significance levels between all the variables are presented in Table 2. We observed the strongest correlations between Zn and Pb (Table 2). Zn is also strongly correlated with Cd. The correlations between As and Ni and between As and Hg are also strong (Table 2). It is noted that lots of these heavy metals are chalcophile elements, meaning they prefer to bond with sulfur, and they are closely related to each other. Therefore, it is reasonable to see a correlation between them. On the contrary, Hg shows no or little correlation with Zn, As, and Cr (Table 2), suggesting that the Hg might have a distinct soil parent source compared to other heavy metals.

Table 2. Pearson correlation coefficients of topsoil heavy metals

	Cu	Zn	Ni	Cr	Hg	As	Cd	Pb
Cu	1
Zn	0.531**	1
Ni	0.308**	0.291**	1
Cr	0.215*	0.451**	0.475**	1
Hg	0.145	−0.038	0.316**	−0.035	1
As	0.241**	0.361**	0.584**	0.456**	0.038	1
Cd	0.212*	0.576**	0.249**	0.325**	–0.149	0.262**	1
Pb	0.524**	0.709**	0.389**	0.347**	–0.079	0.230**	0.298**	1

* Correlation is significant at the 0.05 level; ** Correlation is significant at the 0.01 level.

3.2. Spatial structure and distribution of topsoil heavy metals

The semivariograms could help us understand the level of continuity of the heavy metals spatially, as well as their ranges of spatial dependence. Using the GS+ software, we modeled the eight heavy metals in the topsoils, and obtained their semivariograms parameters such as nuggets, sills, and ranges. The statistics of the semivariograms are summarized in Table 3. Results show that different models are needed in order to best fit different heavy metals. For examples, Cu, Pb, and Hg are fitted well with an exponential model, while Cd and Zn are better fitted using a linear model. For As and Ni, we used a spherical and a gaussian model, respectively (Table 3; Figure 2). The application of different semivariogram models indicates the fundamental differences in the spatial structure of each heavy metal.

Table 3. Semivariogram parameters and models of topsoil heavy metals

Heavy medal	Model	Nugget (C0)	Sill (C+ C0)	Nugget/Sill C0/(C+ C0)	Range (km)
Hg	Exponential	0.247	0.591	41.8%	0.807
Cd	Linear	0.327	0.339	96.4%	1.0859
As	Spherical	0.0001	0.0902	0.1%	0.0902
Cr	Gaussian	0.0262	0.104	25.1%	0.0849
Cu	Exponential	0.0227	0.192	11.8%	0.099
Ni	Gaussian	0.0001	0.0552	0.2%	0.0901
Pb	Exponential	0.013	0.147	8.8%	0.204
Zn	Linear	0.0827	0.0827	100%	1.0844

Figure 2. Semivariograms of topsoil heavy metals

Specifically, Hg and Cd have relatively larger nugget values compared with other heavy metals, suggesting they carry a strong variance within a shorter distance. On the opposite, the nugget values of As and Ni are small, indicating that they do not vary significantly in a short distance (Table 3; Figure 2). Sill values show the highest level of variance. In our dataset, Hg and Cd, spatially, vary the most, while metals like As, Ni, and Zn comparably do not show a high level of variance (Figure 2). In terms of the range value, for Cd and Zn which are fitted by a linear model, their variances do not change with distance, and they have the largest range values, which means that the concentration of these metals in the topsoils might be dependent on the parent soils. Hg also has a large range value, meaning that it is more likely controlled by the original soil material. Instead, As, Cr, Cu, Ni, and Pb achieve their sill values within relatively small range values (Figure 2), indicating their concentration might be affected by a hotspot that yields a radiating effect spatially.

We could use C₀/(C+ C₀) ratio to evaluate the spatial dependence of soil properties. Conventionally, if the C₀/(C+ C₀) ratio is above 75%, the variable only has weak spatial dependence. A C₀/(C+ C₀) value between 25% and 75% indicates a moderate spatial dependence, and a strong spatial dependence is expected when C₀/(C+ C₀) ratio is below 25% [26]. Weak spatial dependence usually indicates extrinsic dominated factors such as anthropogenic activities, while strong spatial dependence is controlled by intrinsic dominated factors like parent soil materials etc. [17, 20, 21]. The C₀/(C+ C₀) ratios of Hg and Cr are 41.8% and 25.1%, respectively (Table 3), which can be categorized as moderate spatial dependence, indicating that the anthropogenic factors changed their spatial correlation extrinsically through either fertilization or industrial production, or both. Cd and Zn have the largest C₀/(C+ C₀) ratios, consistent with their linear model fit. The C₀/(C+ C₀) ratios of As, Cu, Ni, and Pb are below 25% (Table 3), suggesting that these heavy metals have strong spatial dependence and are likely affected by intrinsic factors such as parent soil materials.

ArcGIS is commonly used to interpret spatial variability of heavy metals and to monitor the environment. Results of the semivariogram (e.g. range, sill, and nugget) were utilized to calculate sample weighing factors for spatial interpolation. A simple point kriging procedure was used in which the nearest 8 sampling sites and a maximum searching distance equal to the range distance of the variable were considered. Importantly, the mapping of the conditional probabilities of the properties of the topsoils is critical for policy making and agricultural and environmental management. The threshold values of this property would serve as the basis for the policy design, such as outlining safe or hazardous areas and identifying zones that are either suitable for crop growth or those that must be properly treated. In Figure 3, the spatial distribution maps of eight heavy metals are shown, which are obtained from the semivariograms from the study area. Several heavy metals, like Ni, Pb, and Cu show similar geographic distributions, where the high contents are in the eastern Xuzhou area, and low contents are in the central and western of the area (Figure 3). Zn, Cr, and Cd show near evenly distributed patterns, with similar values across the entire study area (Figure 3). However, As shows a different pattern from all other heavy metals, with high content in the western and central Xuzhou area, and low content only in the eastern corner of the area (Figure 3). Hg is a hazardous metal; the higher soil Hg concentrations are found in almost the entire region except for the southwest margin (Figure 3). In general, the input of Cr and Ni to the topsoils from human activities such as fertilizers and limestone used in constructions is lower than the background value of the soil. Cu and Zn are commonly used in pesticides. Pb usually comes from manures, vehicle exhausts, and coal combustion. Especially important and relevant of the Xuzhou area is the production and use of coal, which reaches the maximum in the cold seasons. Coal burning could lead to high content of Hg, Zn, and Pb from air fallout [27]. Cu, Zn, Pb, Cd, Hg, and As mainly result from a combination of anthropic activities such as industrial, agricultural, and transport. The sites with high-content of heavy metals are mostly in coincidence with the industry locations. All eight heavy metals that contain less than the threshold values, therefore, can be regarded as safe for crop growth in the Xuzhou area.

Figure 3. The spatial distribution maps of the eight heavy metal contents. unit: mg/kg

4. Conclusion

Using the topsoil samplings of the Xuzhou area, a total of eight hazardous heavy metal concentrations were measured. Results show that all heavy metals have a lognormal distribution. There are strong correlations between Zn and Pb, as well as Cd. A significant correlation was observed between As and Ni and between As and Hg, suggesting that they might have the same provenance, i.e. the same topsoil parent materials. The distribution and concentration of these heavy metals are likely influenced, to different degrees, by anthropic and industrial activities. Results from this study provide clear maps showing the spatial variability of eight heavy metals in the Xuzhou area. Our study is valuable for monitoring the quality of topsoils, which has useful implications for policy design in terms of topsoil management in the future. Importantly, all eight heavy metal contents are less than the national threshold values (GB15618–2018). Although these heavy metals are not yet a threat to the environment or to human health at present, they should be closely monitored. The variation of the heavy metal concentrations in the agricultural topsoils bear long-term meanings for improving agriculture and the ecosystem in the region.

Acknowledgments

This study is supported by Grants KC17155 and NYKJ8067. We thank the constructive comments from anonymous reviewers.

Disclosure statement

There are no relevant financial or non-financial competing interests to report.

Supplementary material

Supplemental data for this article can be accessed here.

Additional information

Funding

This work was supported by the Grants KC17155 and NYKJ8067 [Grants KC17155 and NYKJ8067].