Spatial distribution of potentially toxic element and the potential ecological risk of wolfberry-soil system in genuine producing area of Ningxia，China


 The potentially toxic element (PTE) content and physical and chemical properties of 37 wolfberry and root soil samples obtained from Zhongning county and Shapotou district (Ningxia Hui Autonomous Region) were analyzed, and a potential ecological risk assessment was performed. In addition, multiple linear regression and geostatistical methods were used to predict and analyze the PTE concentrations and potential ecological risks of both soil and wolfberries. The PTE content in the soil was within the standard limits, and only one wolfberry sample contained Cd exceeding the standard limit (by 2.7%). The other potentially toxic elements (PTEs) did not exceed the standard limits. In the regression fitting, the models obtained for Cu, Zn, Pb, and Ni contents in the fruits showed a good correlation (p < 0.05), the PTE content of the soil showed an overall trend of high concentrations in the central part of the studied region. The ecological risk of Hg in the soil is rated as high, whereas the risks of the remaining PTEs were mostly medium or low. Further, the predicted PTE content of the wolfberries was generally high in the western region and low in the eastern region, although there were slight variations for different PTEs in some areas.

crop for both medicines and food, and knowledge about their nutritional value has gradually increased [15] . Specifically, wolfberries are rich in vitamins, proteins, minerals, polyphenols, carotenoids, amino acids, and other nutrients, making them an important raw material for the development of nutritional health food products [16] . However, there have been few geographical analyses and predictive studies of the spatial distribution of wolfberries, and there have only been a few studies concerning the PTE pollution and evaluation of the contamination of wolfberries [17][18] . This is crucial because, if the PTE content of the soil used for cultivation is high, wolfberry contamination may occur.
Further, although the PTE content of wolfberries is typically reported to be lower than the Chinese standard limits for Chinese medicinal materials, further investigation is required if cultivation and commercialization are to be increased [18] .
In this study, the Ni, Cu, Zn, As, Cd, Pb, Cr, and Hg contents of wolfberry fruits and the soil around the roots was investigated. First, the PTE concentration of the soil and potential ecological risks were evaluated, and a spatial prediction of the PTE content was carried out. In addition, because of the discontinuous nature of wolfberry plantations, we used a novel method to predict the contamination of wolfberry fruits via soil contamination, which was assessed using spatial prediction analysis based on geostatistical principles. Specifically, we evaluated the soil system around the roots of the wolfberry plants and explored the correlation and relationship between the soil PTE content and its physicochemical properties (pH, organic matter content, and cation exchange capacity) and the PTE content of wolfberry fruits.

Analysis on the characteristics of PTE content in soil and wolfberries
As shown in Table 1, the average contents of Ni, Cu, Zn, As, Cd, Pb, Cr, and Hg in the soil were 40.95, 5.27, 42.46, 8.36, 0.026, 23, 52.89, and 0.069 mg/kg, respectively.
The contents of Ni, Pb, and Hg were slightly higher than the soil background value [19] , whereas the contents of Cu, Zn, As, Cd, and Cr were lower than the background value.
Crucially, all of these values are lower than the standard limits for these metals [20] . The enrichment ratios were in the order of Hg > Ni > Pb > Cr > Cd > Cu = Zn = As = 0, and the enrichment ratios of Hg, Ni, Pb, Cr, and Cd were 97.3%, 64.86%, 51.35%, 35.14%, and 2.7%. On the basis of the enrichment ratios, the enrichment proportion of Ni, Pb, and Hg is greater than 1, which indicates that these species are accumulated in the soil in this area. The order of the coefficients of variation for these eight elements was Cd > Hg > Pb > Cr > Cu > Ni > As > Zn, ranging from 14.63% to 92.72%.
According to Wilding's classification of soil variability based on coefficients of variation, that of Cd, which is greater than 90%, indicates great dispersion and strong variability, i.e., high spatial variability [21] , possibly because this element has been disturbed by human activities for a long time. In contrast, the coefficient of variation of Zn is less than 15%, and its degree of dispersion is small, indicating that its spatial variability is weak. Those of the other elements were moderately variable.
The results for the wolfberry samples in the study area are shown in Table 2 Thus, the Pb contents in the samples are lower than the standard limit. The Cd content in one sample exceeded the standard by 2.7%. The Chinese Pharmacopoeia [22] and "Green standards of medicinal plants and preparations for foreign trade and economy" [23] stipulate the following limits: As ≤ 2 mg/kg, Hg ≤ 0.2 mg/kg, Pb ≤ 5 mg/kg, Cd ≤ 0.3 mg/kg, and Cu ≤ 20 mg/kg. Thus, the Cu, As, and Pb in the wolfberry samples did not exceed the standards. In China, except for those of As, Hg, Pb, Cd, and Cu, national standard limits for the other PTEs Chinese herbal medicines have not established [24] .
The mean BCF values for Ni, Cu, Zn, As, Cd, Pb, and Cr in the wolfberries were 0.023, 1.762, 0.475, 0.026, 5.821, 0.017, and 0.056, respectively. The mean BCF of Cd (5.821) was the highest, and that of Pb (0.017) was the lowest. In particular, the BCF values for Cu and Cd are greater than 1, indicating that these elements were highly enriched in the wolfberries. As an organism at a high trophic level, wolfberries contain Cu and Cd concentrations higher than those in the environment, and the contents increase with longer growth and development time.

Correlation evaluation of the wolfberry-soil system
The Pearson correlation coefficient and statistical significance are shown in Supplementary Table S1 online. There was a significant positive correlation between the Cu content of the wolfberries and the As content in the soil (p = 0.022), indicating that the wolfberry plants absorb Cu from the soil (p < 0.1). The Pb content in the wolfberries was negatively correlated with Ni and Zn contents in the soil (p = 0.06, p = 7 0.093, respectively), which indicates that these metals may play an inhibitory role in the absorption of Pb from the soil. There were no significant correlations between the other PTEs. On adding the physicochemical properties of the root soil (pH, organic matter content, and cation exchange capacity) to the model, the pH was found to be negatively correlated with the Ni, Cu, and Cr contents but positively correlated with As, Cd, and Pb in the wolfberry fruit. The organic matter content and cation exchange capacity in the soil were positively correlated with the PTE content of the wolfberries; in particular, the Cu and Zn contents were positively correlated with the cation exchange capacity (p = 0.006 and p = 0.036, respectively), which indicates that a high cation exchange capacity promotes the absorption of Cu and Zn by wolfberries.
In this study, a multiple linear regression equation with the PTE content of the wolfberry fruit as the dependent variable and the soil indexes as independent variables was obtained. As shown in Supplementary Table S2 online, when the Cu, Zn, Pb, and Ni contents of the wolfberry fruits were used as dependent variables (p < 0.05), the R value for the fits for Cu and Pb were high, indicating a good fit. When the As and Cr contents in wolfberry fruits was used as dependent variables to fit the equation, the pvalue slightly exceeded 0.05 (p = 0.053). When the Cd content in wolfberry fruits was used as the dependent variable to fit the equation and the level () was set to 0.15, no independent variable satisfactorily fit the model. However, only when  was set to 0.2, the soil pH showed correlation (for these parameters, p ranged from 0.186 to > 0.05).
In the Pearson correlation analysis, the correlation between the independent and dependent variables was significant (p = 0.186), indicating that the multivariate linear 8 relationship was weak.

Spatial distribution of soil PTEs and potential ecological risk
Supplementary Table S3  In contrast, the p-value of As is 0.396 and the Z-score is less than 1, indicating a weak spatial correlation. Therefore, the As content in the soil is random.
The results of the semivariogram analysis of the seven PTEs in the soil of the wolfberry root system in the study area are shown in Table S4. The semivariogram nuggets ( 0 ) of the seven PTEs are all greater than 0. 0 /( + 0 ) is the ratio of the spatial variation caused by the random part to the total spatial variability. The ratios of the seven detected PTEs in the root soil samples ranged from 0.012 to 0.276, and that of As was higher than 0.25, indicating that the other six PTEs had a strong spatial correlation. The ratios for Zn, Pb, and Cr are all less than 0.1, which is consistent with 9 the spatial autocorrelation conclusion given by the Moran's I analysis above, suggesting that that analytical method is reliable. The variogram model type for Ni, Cu, As, and Hg soil contents was Gaussian, whereas the rest were spherical. The R 2 values of each model are between 0.381 and 0.879. Thus, these fit well with the corresponding theoretical model, which helps the Kriging interpolation analysis. The semivariogram ranges of the seven PTEs are between 6870 and 99731.49, which is larger than the sampling distance, thus meeting the requirements for spatial analysis.
Comparing two interpolation graphs (GS+ 10 in Fig. 1 and ArcGIS 10.0 in Fig. 2), we can prove the applicability of the semivariogram fitted by GS+ 10 in ArcGIS 10.0. The 3D inner illustration ( Fig. 1) prepared by GS+ 10 shows intuitively the trend for each PTE within the distribution range. Because the plane rectangular coordinate system must be converted when entering coordinates in GS+ 10, the interpolation diagram shows the opposite direction to convention.
In the areal distribution ( Fig. 2) based on the GCS-Krasovsky-1940 coordinate system and D-Krasovsky-1940 datum plane, the trends for each PTE in the distribution range are consistent with the 3D inner illustration produced by GS+ 10, which proves the applicability and accuracy of the semivariogram parameters obtained in GS+ 10 for Kriging interpolation and extrapolation in ArcGIS 10.0. Generally speaking, the eastwest trend is obvious, whereas the north-south trend is weak. The western region comprises sandy land and mountainous areas, and the contents of Zn, As, and Pb in the soil are relatively low, whereas the contents of Cr, Hg, and Ni are relatively high. In the middle cultivated area, all elements show peaks, indicating that the cultivated land area has a higher PTE content than those of the mountains and sandy land, most likely because of cultivation, sewage discharge, and fertilizer use, which may all contain PTEs.
The highest Cu content is close to the lowest value, indicating that Cu migration in soil is not significant. In contrast, the As content is high, suggesting that attention is required in this area. There is also a relatively high As content in the cultivated land. The eastern and southeastern regions are plain regions and have high values of Ni and Cr, although the Cu content is also high in this region. Hg , there are 1 serious risk level, 9 high risk levels, 24 high risk levels, 2 medium risk levels, and 1 low risk level. According to the overall risk index, there were no points at which total was rated as heavy or serious risk, although 1 point was rated as high risk, 18 points as medium risk level, and 18 points as low risk. The contribution of each PTE to total were Ni = 3.33%, Cu = 0.74%, Zn = 0.45%, As = 4.38%, Pb = 3.48%, Cr = 1.10%, and Hg = 86.51%. Except for Hg, the contribution rate of other PTEs is less than 5%, and these six PTEs show no potential ecological risks.

Supplementary
Supplementary After comparing the two interpolation results, the error in the stable distribution was found to be smaller, so it was used for further analysis. Using the GCS-Krasovsky-1940 coordinate system and D-Krasovsky-1940 reference plane, the Kriging interpolation method was used for interpolation prediction analysis for the risks and total risks of seven PTEs. The results are shown in Fig. 3. Six PTEs (Ni, Cu, Zn, As, Pb, and Cr) were assessed as low risk in the whole region, whereas Hg is a high risk, and there are serious risks in the eastern region. On the basis of the Etotal distribution, the rating of the total risk has declined compared to previous years, and most areas are medium risk, some areas in the middle and southeast regions are low risk, and only a few areas are high risk.  Supplementary Table S7 online, the proportion of seven PTE elements in the wolfberries is between 0.030 and 0.361, and the Cu, As, Pb, and Cd contents in the wolfberries were higher than 0.25, which indicates that these six PTE contents are 12 moderately spatially correlated. The Zn, Ni, and Cr contents range from 0 to 0.25, which indicates that there is a strong spatial correlation between the three elements in wolfberries. The variogram model type for Pb, Cd, and Cr contents in wolfberries was Gaussian, that for Cu content was exponential, and those of the other metals were spherical. The In contrast, the Cr distribution is smooth, and the first 20% of the concentration range (low concentration) is located east of the middle region, with no high-to-low west-to-13 east trend. The concentrations of As, Cd, Ni, and Zn in the wolfberries fall in the first 20% of the concentration range, whereas, for Cr and Cu, this concentration range is distributed over a smaller area. With regard to the PTEs limits for Chinese herbal medicines discussed in Section 3.1, the spatially predicted concentrations of all PTEs did not exceed the standard limits.

Discussion
Overall, in the soil, Ni, Pb, and Hg have high enrichment ratios and moderate spatial variability, which indicates that these three elements are greatly influenced by agricultural fertilization and irrigation (Table 1). However, all PTE contents are below the standard limits and, thus, present no risk. In a study by Li [25] , the average contents of PTEs in the root soil of the study area were 22.8, 0.225, and 12.2 mg/kg for Pb, Cd, and As, respectively. Thus, the current values represent decreases of 33% for Pb and 90% for Cd compared to those 7 years ago. In a study by Bai [26] , the samples exceeding the background values of Pb, Ni, Zn, Cu, Cr, and Cd in Ningxia soil accounted for more than 83% of the total samples, which is higher than the results obtained in this study, and the coefficient of variation for Zn was slightly different.
Except for Cd, none of the other PTE contents in the wolfberries exceeded the standard limits, and only one sample contained Cd exceeding the standard by 0.05 mg/kg. Therefore, the quality of wolfberries in the study area is high, and they can be safely consumed. The test results provide a reference for the safe production of wolfberries based on standard limits for PTE elements. Qi et al. and Geng et al. studied the PTE content in wolfberries from different areas, and only Cd (35.56%) was found 14 to be present above the standard limit in Qinghai Province [18,20] . In a study of wolfberry plantations in western China by Xiao et al., it was also concluded that farmland was polluted during years of cultivation, causing, in particular, cadmium pollution, suggesting that wolfberries are subject to excessive Cd exposure, and the mechanism of Cd adsorption by wolfberries should be studied further and considered during planting [27] .
For Cd and Cu, the BCF was greater than 1, and their concentrations in the wolfberry fruits are far greater than those in the soil, indicating that absorption and accumulation occur. Petukhov [28] found that the BCF value of Cu in Chinese herbal medicine is generally greater than 1; for example, the BCF for chamomile is 4.67, whereas that for Cd is also higher than 1 (as found in this study). However, Garg et al. [29] determined BCF ranges for Cu and Cd in different crops of 0.011-0.148 and 0.006-0.179, respectively. Possibly, the enrichment of Cu and Cd in Chinese herbal medicines is greater than those in these crops, although this requires further research. Qi et al. [18] found that, compared with those of other PTEs, the enrichment coefficient of Cd in wolfberries was also larger, ranging from 0.318 to 0.903, and the average enrichment coefficients of Ni, As, and Pb were less than 0.1, which is consistent with the results of this study.
In a correlation study of the PTE content of wolfberries and the PTE content and physicochemical soil properties, Qi et al. [18] concluded that, except for Cd in wolfberries and Pb in root soil, there was no clear correlation between PTEs in other wolfberry-soil systems, which is consistent with the results of this study. Wang et al.
15 [30] found that the adsorption of PTE ions in soil increases with increase in pH, resulting in a decrease in active PTE ions in soil, which affects the absorption of these metals by the wolfberry plants. In this study, this was also observed; however, there was a negative correlation between pH and Cr content, but the correlation was very weak and requires further study. Geng et al. [31] reported that soil organic matter affects the absorption of PTEs (especially As) because it contains a large number of active functional groups, and Chen et al. [12] successfully fitted a linear regression equation for the As content of five crops and that of the root soil. However, good fits were only obtained for corn and red jujube. Possibly, the analysis is limited to these crop species or the study involved too few sampling points. Adams et al. [32] reported that the Cd content and soil pH are the main factors affecting the Cd content in wheat and barley using multiple regression analysis. At the same time, McBride [33] obtained an equation by fitting that showed that the soil Cd content and pH significantly affected the Cd content of lettuce, beet, and maize leaves. Zhu et al. [34] also obtained a regression equation for the Cd content of Sanchi roots as a dependent variable and Zn concentration and pH value of the soil as independent variables, achieving r = 0.765 at p < 0.01. However, although many previous studies have proven the applicability and universality of this method, the multivariate linear relationship between Cd content of wolfberries and the soil PTE content and physicochemical properties was not strong, although those of other metals were. This may be a result of the limited sample size. However, if this method were used for the analysis of widely planted crops, more samples could be obtained, and the results could be improved. In addition, the determination of more soil physicochemical properties and an increase in the number and relevant types of independent variables could increase model accuracy.
In this study, the dispersion of As in the soil was found to be significant, whereas Chen et al. [12] found that the As content in soil in Xinjiang and other places in China showed significant spatial autocorrelation. Therefore, we cannot exclude a spatial autocorrelation in the As content in soil. Possibly, this was not observed in this study because of the small sample size or inappropriate sample point selection. The spatial variability of soil As content was assessed to be moderately correlated, which may be caused by the superposition of structural factors (such as topography, soil type, parent material, climate) and random factors (such as fertilization, cultivation, planting system, land use intensity, and other human activities), which is similar to the trend obtained by Moran's I analysis. The nugget of each PTE indicates that there are positive basal effects caused by sampling error, short-distance variation, and random variations in the soil of the study area [35] . The research area used by Bai et al. [26] overlaps with our research area. Because of the high density of research points, the overall trend of overlapping regions is the same, but different results were obtained, which is mainly manifested in the hierarchy and variability of the interpolation map. For example, compared with our Ni distribution map, that of Bai et al. has four content levels, whereas, in the same area, our map only contains two content levels. Considering that the analysis of more common crops would enable the analysis of more root soil, this hierarchy gap could be avoided. Geng [20] evaluated the PTE content of the soil in a wolfberry plantation in Qinghai. Using the potential ecological risk assessment method, 38.9% of the sampling sites were assessed to be medium risk for Ni, whereas most of the other PTEs were rated as low risk. A high Hg content in the soil may be a major feature of wolfberry plantations. To date, studies have not investigated the potential ecological risks of Hg in the soil in this study area, so the Hg risk grade presented here will enable the development of risk management strategies at the government level. In this study, we innovatively combined the risk model with Kriging interpolation to visualize the risks, thus providing a reference to monitor changes in soil PTE content and improvements in soil health.
The Chinese government plans to increase the area used for cultivating wolfberries in the Ningxia Hui Autonomous Region, making the need to forecast and analyze the PTE content and ecological risk of wolfberries in the study area crucial.
Because the planting of wolfberries is not spatially continuous, which is not consistent with the requirements for data in geostatistical spatial interpolation, but the soil indicators are continuous, we attempted to use the predicted PTE contents of wolfberries after linear fitting of the soil data for interpolation of the study area.
Although the results contained errors, this novel concept for the interpolation of PTE concentrations and the health risk of crops is promising. In previous studies, the health risks of tomatoes and Chinese rice were assessed, and this method can also be combined with geospatial analysis to visualize the health risk presented by different crop cultivation areas [36][37] . The resulting interpolation map can provide a reference for the site selection of wolfberry cultivation sites, as well as legislation for PTE standards for wolfberries. 18

PTE determination
The Ningxia Hui Autonomous Region is the main production area of Chinese wolfberries and the origin of this plant; specifically, Zhongwei is a key planting area.
Shapotou District and Zhongning County in Zhongwei are located from 104°17' E to 106°7' E and 36°6' N to 37°50' N [38] . As a result of its suitable latitude, long sunshine period, large diurnal temperature range, and abundant water resources, this region has been developed into a high-quality planting area [30][31][32][33][34][35][36][37][38][39][40] . Therefore, these two areas were used in this study (Fig. 5). Using the spot placement method [41] , 37 wolfberry samples ( Fig. 5) and surface soil (from the surface to 20-cm depth, approximately 300 g) of the soil around the wolfberry root system was collected. Supplementary Table S8 online, and the methods of determination, basis, units, and detection limits for PTEs, as well as the pH, organic matter content, and cation exchange capacity of the soil, are listed in Supplementary Table S9 online. The collected soil samples were spread flat to allow the removal of foreign objects such as stones and animal and plant residues in the soil samples. The samples were allowed to dry in the air and were not exposed to direct sunlight. Then, they were crushed with wooden sticks, any impurities were removed, mixed, and crushed. The samples were selected using the quartering method. Each sample was passed through a 10-mesh (2 mm) nylon sieve, and all the sieved samples were placed on a colorless polyethylene film, fully stirred and mixed; then, two parts were selected following the quartering method. One part was transferred to the sample 19 inventory, and the other part was finely ground. The finely ground sample was then divided into two parts by the quartering method, and the sample was further ground and passed through a 100-mesh (0.149 mm) sieve. These samples were used for analysis.

The instrumentation used for PTE analysis is listed in
Each wolfberry weighed about 0.2 g averagely, and each sample was accurately weighed into a digestion tube. Then, 10 mL of mixed acid (HNO3:HClO4 = 4:1) was placed in a 50-mL digestion tube, stoppered, sealed, and left overnight [42] . For each sample, a parallel experiment was carried out, and standard addition and reagent blank samples were prepared. The next day, the plug was removed, and 4 mL of 30% H2O2 was added to each digestion tube. The tube was placed on a plate in a graphite furnace to catch the acid and heated. The graphite furnace was programmed to increase from room temperature to 60 °C, from 60 to 120 °C, and, then, to 180 °C. Subsequently, the sample was held at this temperature for 4-6 h. After the acid had been completely removed, the digestion tube was cooled to room temperature, and the sample plug was added. Before the metal determination, ultrapure water was added to the digestion tube several times in small amounts, and the tube was shaken to dissolve the deposited crystals, transferred to a 10-mL sample tube, and, again, shaken well.

Bioconcentration factor of wolfberries
To compare the ability of different PTEs to migrate from the soil into the wolfberries, the bioconcentration factor (BCF) was used [29] . The BCF was calculated using equation (1). If BCF ≤ 1.00, a plant absorbs but does not accumulate metals.
However, if BCF > 1.00, a plant accumulates metals [28] . Thus, this is an indicator that characterizes the degree of concentration or enrichment of a chemical substance by an 20 organism.

=
(1) Here, is the concentration of different PTEs in wolfberry fruit (mg/kg), and is the concentration of PTEs in the soil (mg/kg).

Multiple linear regression
Multiple linear regression was used to analyze the linear relationship between different dependent and multiple independent variables. F-tests were used to assess model quality based on the selected independent variables; because the number of samples was small, the α value was set to 0.10-0.15 [43] . In this study, the PTE contents of the wolfberries were used as dependent variables, and the PTEs in the soil around the root and related physical and chemical properties (pH, organic matter, and cation exchange capacity) were introduced as independent variables. A stepwise regression method was used to perform multiple linear regression. The general form of the equation is given by equation (2).
Here, Y is the predicted dependent variable (PTE content in wolfberry samples), 0 to are partial regression coefficients, and 1 to are different independent variables (for example, soil PTE content and other physical and chemical properties).

Normality tests and conversion
As in statistics, in geostatistics, we assume that a sample obeys the normal distribution [35] . The normality tests revealed that, of those of the eight PTEs studied, the concentration distributions of Ni, Zn, and As in the soil were normal, and the p-values of the other five PTEs were all less than 0.05. In ArcGIS 10.0, three methods are provided for the normal transformation of the original data: log, Box-Cox, and Arcsin.
As shown in Supplementary Table S10

Spatial autocorrelation analysis
In this study, ArcGIS (version 10.0) was used to analyze the spatial autocorrelation of the concentration of seven detected PTEs in the soil. Spatial autocorrelation refers to the potential interdependence between the observed data of certain variables in the same partition. Moran's I, a global aggregation characteristic, is the most common statistic used to characterize spatial autocorrelation analysis. After the variance has been normalized, the composite index lies between -1 and 1. When Moran's I > 0, there is a positive spatial correlation (clustering), and larger I values indicate a stronger spatial autocorrelation; when Moran's I < 0, there is a negative spatial correlation (dispersion), and smaller indices indicate greater spatial differences.
When Moran's I = 0, the space is distributed randomly [44] . 22 Geostatistics is based on the theory of regionalized variables, for which a function of the variation (a semivariogram) is the basic tool. This can be used to describe the structural characteristics and randomness of the key elements of the soil or the method of spatial correlation and dependence [12] . The general form of a semivariogram is given by equation (3).

Soil potential ecological risk index method
To evaluate the PTE pollution in soil, we used the potential ecological risk index proposed by Hakanson [45] . This method is uses principles of sedimentology, PTE behavior, and environmental behavior to evaluate the pollution of soil or sediments by PTEs [46] . In addition to the PTE content of the soil, the risk index also considers factors, such as multi-element synergy, toxicity, pollutant concentration, and environmental 23 sensitivity to PTE pollution, and adopts a comparable and equivalent attribute index classification method for evaluation [47] . The potential ecological risk index is calculated using equations (4) and (5).
Here, total is the potential ecological risk index, is the individual coefficient of potential ecological risk of metal element i, is the toxicity response coefficient of the metal element (Supplementary Table S11 Table   S11 online) [19] .

Data Availability
All data generated or analyzed during this study are included in this published article.
a u indicates that it is not detected. b -indicates that the limit standard is not specified in China. whereas the trends for Ni, Zn, Pb, and As are high in the north and low in the south. In contrast, uniquely, Hg has a high distribution in the south and low distribution in the north. The highest value for 34 Ni appeared in the west wave crest, and the lowest value was located in the trough in the southeast. The Cu peak appears on the east side of the wave crest, whereas the lowest point is located in the west surrounded by other low values. The peak and valley values for Zn both appear in the middle region, the distance between the extreme points is close, and the Zn content changes significantly. The highest value for As appeared in the northeast, and the lowest value appeared in the west. The strong fluctuation in Pb content, whose distribution shows multiple peaks. Concerning the Cr distribution, there is a peak in the east region, a low value in the middle region, and a small peak in the middle lowcontent region. The Hg content did not show any sudden increases or decreases, but a peak appeared in the western region, and the lowest point appeared north of the central region. 38 Figure 5 Geomorphology of the study area and sampling sites distribution Note: The southern part of the study area is mountainous area, the northwest part is desert area, and the dark green area along the Yellow River in the middle part is plain cultivated area. ArcGIS Desktop 10 software with the authorization number of EFL564098460