Risk Assessment and Source Apportionment of Heavy Metals in Soils from Handan City

: Soil-heavy metals are potentially harmful to the ecosystem and human health. Quantifying heavy metals sources is conducive to pollution control. In this study, 64 surface-soil samples were collected in Handan city. Cr, Mn, Ni, Cu, Zn, Cd and Pb were determined; then, their spatial distribution in the sampling area was drawn by ArcGIS. The pollution index ( PI ) method, geoaccumulation index ( I geo ) method, Nemerow integrated pollution index ( NIPI ) and pollution load index ( PLI ) were used to evaluate the pollution level of heavy metals in surface soil; then, an ecological and health risk assessment of soil-heavy metals was carried out. Combined with the spatial distribution, correlation analysis, cluster analysis, PCA and PMF model, the pollution sources of heavy metals in soil were identiﬁed and apportioned. The results showed that the average content of Cd was nearly ten times that of the background limit, which was the most serious among the studied metals. In terms of non-carcinogenic risk, Cr had the highest value, followed by Pb. In terms of carcinogenic risk, Cd, Cr, and Ni had an acceptable or tolerable risk. Three pollution sources were identiﬁed by cluster analysis and PCA, including trafﬁc sources with Cu, Pb and Cd as main loads, industrial sources with Mn, Cd and Zn as main loads, and natural sources with Cr and Ni as main loads. The PMF model analyzed three main factors: trafﬁc source (17.61%), natural source (28.62%) and industrial source (53.77%). The source categories and the main load elements obtained from the source apportionment results were consistent with the source identiﬁcation results.


Introduction
Soil security is significant for food safety and public health [1]. With the speedy development of the economy, the emission of pollutants by water, atmospheric deposition, and rain leaching quickly increases, which results in rich heavy metal content in soil [2][3][4]. Heavy metals, with toxicity, persistence and biological accumulation, are absorbed by humans through ingestion and skin contact, as well as inhalation, causing a variety of diseases [1,5]. For example, excessive Pb intake will not only harm the nervous system, but also harms the hematopoietic system, leading to anemia [6]. People exposed to Cd for a long time are more likely to suffer from lung cancer and fracture [7]. In addition, excessive Cu, Ni and Cr intake can also have harmful influences on health [8]. Therefore, heavy metal pollution has attracted people's attention [9].
There has been much research on the pollution status, distribution, risk assessment and source of heavy metals in soil [10][11][12]. Some methods have been suggested to evaluate heavy metal pollution level, such as pollution index (PI), geo-accumulation index (I geo ) contamination factor (CF), Nemerow integrated pollution index (NIPI), and pollution load index (PLI) [4,[13][14][15][16]. The PI can be used to assess the contamination degree of a single

Study Region and Soil Sampling
The study region (114 • 24 -114 • 32 E, 36 • 33 -36 • 40 N) covers an area of about 160 km 2 . It is located in Handan City, south of Hebei Province. It is an important national highway transportation hub city in China. The area features a continental monsoon climate with an average annual temperature of 13.5 • C. Handan city is known as the modern "steel capital". Its industry is developed, and it is also an important source of grain, cotton and other crops. The western part of the study area is an industrial zone, with and iron and steel group, fertilizer plant, asphalt-processing plant, refractory plant, metallurgical machinery plant, coking plant, concrete company and other enterprises. Farmland is scattered in the eastern and southwestern parts, and a thermal power plant exists in the northeast. These activities cause heavy metals to enter the soil. The topsoil is comprised of loess and subclade.
Based on the 1.5 km × 1.5 km grid pattern of the whole study region, combined with the terrain distribution of the study area, 64 surface soil samples were collected. The coordinate of each sampling location was located by GPS. The distribution of sampling points is presented in Figure 1. Each sample consisted of five subsamples within 10 m of a specific sampling location. Each sample was packed in 1 kg and brought to the laboratory.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 3 of 18 activities cause heavy metals to enter the soil. The topsoil is comprised of loess and subclade. Based on the 1.5 km × 1.5 km grid pattern of the whole study region, combined with the terrain distribution of the study area, 64 surface soil samples were collected. The coordinate of each sampling location was located by GPS. The distribution of sampling points is presented in Figure 1. Each sample consisted of five subsamples within 10 m of a specific sampling location. Each sample was packed in 1 kg and brought to the laboratory.

Chemical Analysis
The soil samples were oven-dried at 60 °C , and then debris and stones were removed with plastic tweezers [32]. After being ground by agate grinder, the soil sample was sieved through 0.149 mm nylon mesh, then stored in the plastic sample bag. Referring to the research method of Wang [33], about 0.2 g of sample was digested using HCl-HF-HNO3. All solvents were Optima-grade (Fisher Scientific). Then, concentrations of heavy metals in soil were determined by an Element XR inductively coupled plasma mass spectrometer (ICP-MS) (ThermoFisher, Waltham, MA, USA). The relative standard deviations of the soil samples were < 10%. The detection limits of Cr, Mn, Ni, Cu, Zn, Cd and Pb were 0.4, 0.02, 0.3, 0.6, 2.0, 0.03 and 2.0 mg/kg, respectively. To ensure the reliability and quality of the data, standard reference soil (GBW07408) was purchased from the Center of National Reference Materials of China. Recovery values, ranging from 90% to 108%. All observed results were corrected with blanks. The relative standard deviations (RSD) were lower than 10%.

Pollution Evaluation
PI, NIPI, PLI and Igeo [34][35][36] were chosen to assess the pollution level of heavy metals in Handan soil. PI and Igeo are often used to assess the pollution of a single pollutant. NIPI and PLI were used to evaluate the comprehensive pollution of multiple pollutants. The calculation method is as follows:

Chemical Analysis
The soil samples were oven-dried at 60 • C, and then debris and stones were removed with plastic tweezers [32]. After being ground by agate grinder, the soil sample was sieved through 0.149 mm nylon mesh, then stored in the plastic sample bag. Referring to the research method of Wang [33], about 0.2 g of sample was digested using HCl-HF-HNO 3 . All solvents were Optima-grade (Fisher Scientific). Then, concentrations of heavy metals in soil were determined by an Element XR inductively coupled plasma mass spectrometer (ICP-MS) (ThermoFisher, Waltham, MA, USA). The relative standard deviations of the soil samples were < 10%. The detection limits of Cr, Mn, Ni, Cu, Zn, Cd and Pb were 0.4, 0.02, 0.3, 0.6, 2.0, 0.03 and 2.0 mg/kg, respectively. To ensure the reliability and quality of the data, standard reference soil (GBW07408) was purchased from the Center of National Reference Materials of China. Recovery values, ranging from 90% to 108%. All observed results were corrected with blanks. The relative standard deviations (RSD) were lower than 10%.

Pollution Evaluation
PI, NIPI, PLI and I geo [34][35][36] were chosen to assess the pollution level of heavy metals in Handan soil. PI and I geo are often used to assess the pollution of a single pollutant. NIPI and PLI were used to evaluate the comprehensive pollution of multiple pollutants. The calculation method is as follows: where PI i represents pollution index of i th metal. C i and C 0 (mg/kg) stand for the content of i th metal and its background value. The background values refer to the background values of soil elements in China [37]. (PI i ) max represents the maximum PI i . I geoi stands for the geo-accumulation index of i th metal. NIPI represents the Nemerow integrated pollution index. PLI n stands for the pollution load index at the n th sampling site. PLI zone represents the pollution load index across the study area.

Ecological Risk Evaluation Approach
The calculation method of the potential ecological risk index (PERI) is as follows [38]: where RI is the pollutants' total ecological risk, E i represents the ecological risk of i th metal, T i denotes the toxic response factor [34]. The toxic response factor of Cr, Mn, Ni, Cu, Zn, Cd and Pb were 2, 1, 5, 5, 1, 30 and 5, respectively. The classification standard of PI, I geo , PLI, NIPI, and PERI is revealed in Table S1 [8].

Health Risk Evaluation Approach
The caluclation of the health risk is as follows: where ADD ing , ADD dermal , and ADD inh are the mean daily exposure dose to ingestion, dermal and inhalation absorption from soil heavy metals (mg/kg/d) [8]. C stands for the content of metal in soil (mg/kg). A detailed description of the other parameters can be found in Tables S2 and S3. HQ i represents the non-carcinogenic risk of ith heavy metal, while HI stands for the total non-carcinogenic risk [39]. When HI or HQ < 1, the risk is small or negligible. When HI or HQ > 1, the metal may have a non-carcinogenic risk. CR represents carcinogenic risk. Humans suffer unacceptable carcinogenic risk if the risk exceeds 1 × 10 −4 [34]. However, there will be no significant health impact if the risk is lower than 1 × 10 −6 . It is considered that there is an acceptable carcinogenic risk if the risk is between 1 × 10 −6 and 1 × 10 −4 .

Source Identification
Spearman correlation analysis (SCA) was applied to determine the correlation between elements [40]. Elements with P < 0.05 indicated that a prominent correlation exists between elements. The content data of heavy metals were analyzed by cluster analysis (CA) and principal component analysis (PCA) to identify the sources. CA was performed in this study to divide variables into several mutually exclusive clusters. Based on the Ward method, Euclidean distance was utilized to measure the distance between clusters with similar element contents [25]. PCA was adopted to identify the sources of heavy metal pollution in soil by using varimax rotation with Kaiser Normalization [25]. The SCA, CA and PCA were performed in SPSS 19.0.

PMF Model
The PMF model is a multivariable analysis model first proposed by Paatero in 1994 [41,42], which quantitatively analyzes pollution sources. In this study, PMF 5.0 was used for the quantitative analysis of pollution sources. The equation is as follows [43]: where x ij stands for the j th heavy metal content in i number of samples; g ik represents the k th pollution source contribution for i number of samples; f kj refers to the source profile of j th heavy metals for the k th source factor; e ij represents residual value for the j th metal in i samples. The residual matrix e ij is computed by minimizing the goal function Q [12]: where u ij represents the j th metal uncertainty in i number of samples [1], which is calculated as follows [44]: where σ j refers to the relative standard deviation of the content of the j th heavy metal [39]; MDL is the detection limit of the corresponding element.

Spatial Analysis Method
The spatial information of soil materials can be directly displayed by employing geostatistical methods. The spatial distributions of heavy metal concentration and RI were found using the ordinary Kriging interpolation method. All spatial mapping was carried out in ArcGIS 10.2.

Heavy Metals Content in Soil
The results of statistical analysis of heavy-metal soil contents in the study area are displayed in Table 1. The mean concentration of heavy metals followed the order of Mn (588.25 mg/kg), Zn (158.32 mg/kg), Cr (94.4 mg/kg), Ni (47.47 mg/kg), Cu (43.16 mg/kg), Pb (33.27 mg/kg) and Cd (0.92 mg/kg) (from greatest to least). Compared with the background values, the average concentrations of studied metals besides Mn were higher than the background limit of soil in Hebei Province, which indicated that the topsoil's heavy metal content was mainly affected by external and human factors [37]. It is worth mentioning that the average content of Cd was approximately 10 times that of the soil background, indicating that Cd was enriched in the topsoil. According to Table 1, the coefficient of variation (CV) of Cu, Cr, Ni, Cd and Pb was great, which indicated that the distribution of these metals was uneven. In contrast, the CV of Mn and Zn was relatively low, showing a relatively uniform spatial distribution. The K-S test indicated that the concentrations of Zn, Cd and Pb obeyed normal distribution (p > 0.05), while the concentrations of Cr, Mn, Ni and Cu obeyed lognormal distribution.

Spatial Distribution of Soil Heavy Metals
To further characterize the anthropogenic source, we analyzed the spatial distribution of heavy metals in soil using the ordinary Kriging method. In this process, the distribution of Zn, Cd and Pb was directly obtained, while the spatial distribution of Cr, Mn, Ni and Cu was found using inverse logarithmic transformation. The distributions of seven heavy metal concentrations in soils were shown in Figure 2 and the spatial distribution of heavy metals in soil had obvious characteristics in study area area.
The distribution of Mn, Zn and Cd contents was similar, and the high content areas were centered on the industrial area (Figure 2a-c). The content distribution of Cr and Ni was very similar, and the high content areas were in the middle and northeast of the research region (Figure 2d,e). The soil samples from the main roads with a high vehicle density (near South Ring Road, Renmin Road and Lianfang East Road) had high Cu and Pb contents (Figure 2f,g). Some studies found that areas with well-developed transportation had high Pb and Cu concentrations in soil [38].
The spatial distribution map, using the results of the kriging interpolation, offers some information about heavy metals in soil. Southwest Handan has a large concentration of industrial areas. This region has a high number of industrial enterprises, extensive social economic activity, and sizable emissions from steel smelter, coking plant, and traffic. Since the spatial distributions of Zn, Cd and Pb were mainly concentrated in the southwest area of Handan, it was speculated that Zn, Cd and Pb in soil were usually affected by industrial and traffic activities. This is similar to the previously published research in China.

Pollution Index and Geoaccumulation Index
The PI were calculated through Equation (1), as shown in Figure 3a. These heavy metals in Handan city had different levels of pollution compared to the soil background in Hebei. The PI values of the heavy metals in soil were ranged from 0.14 (Cr) to 25.42 (Cd). It should be noted that the PI value of Cd in 95% of soil samples sites was classified as class IV, indicating that the contamination of Cd was serious in almost the entire study area. As shown in Figure 3a, the PI value of Mn in 77% of the sampling points was grouped as class I, and 23% obeyed class II, indicating no pollution in most sampling sites. According to the PI values, Cr, Cu, Ni, Pb and Zn belonged to the pollution II group in most sampling sites.
west area of Handan, it was speculated that Zn, Cd and Pb in soil were usually affected by industrial and traffic activities. This is similar to the previously published research in China. The research carried out by Ma et al. (2016) found a significant positive correlation between the concentrations of Pb and Zn in soil and traffic emissions in Changsha. Industrial activities were also associated with the concentration of Cd [45].

Pollution Index and Geoaccumulation Index
The PI were calculated through Equation (1), as shown in Figure 3a. These heavy metals in Handan city had different levels of pollution compared to the soil background in Hebei. The PI values of the heavy metals in soil were ranged from 0.14 (Cr) to 25.42 (Cd). It should be noted that the PI value of Cd in 95% of soil samples sites was classified as class IV, indicating that the contamination of Cd was serious in almost the entire study area. As shown in Figure 3a, the PI value of Mn in 77% of the sampling points was grouped as class I, and 23% obeyed class II, indicating no pollution in most sampling sites. According to the PI values, Cr, Cu, Ni, Pb and Zn belonged to the pollution II group in most sampling sites.  The Igeo classification of heavy metals is shown in Figure 3b. Compared with the pollution index, the Igeo assessment was divided into seven, more detailed levels. The mean Igeo of the heavy metals in soil were arranged in the following order: Cd (2.47) > Zn (0.31) > Cu (0.01) > Ni (−0.16) > Pb (−0.25) > Cr (−0.48) > Mn (−0.65). Except for Cd, the Igeo of the other six heavy metals ranged from classⅠto class IV, and the Igeo values of most sampling points were in the range of class I-class II, showing that heavy metal pollution in soil was between uncontaminated and moderately contaminated conditions. The Igeo values of Cd in 95% of soil samples belonged to class IV, which showed that the pollution of Cd was serious compared with other metals, ranging from moderate contamination to heavy contamination. The results of Igeo had more non-pollution points than those of the pollution index. Considering the high Cd contamination levels in urban areas, it is necessary to identify the sources and carry out quantitative assessments of risks. The I geo classification of heavy metals is shown in Figure 3b. Compared with the pollution index, the I geo assessment was divided into seven, more detailed levels. The mean I geo of the heavy metals in soil were arranged in the following order: Cd (2.47) > Zn (0.31) > Cu (0.01) > Ni (−0.16) > Pb (−0.25) > Cr (−0.48) > Mn (−0.65). Except for Cd, the I geo of the other six heavy metals ranged from class I to class IV, and the I geo values of most sampling points were in the range of class I-class II, showing that heavy metal pollution in soil was between uncontaminated and moderately contaminated conditions. The I geo values of Cd in 95% of soil samples belonged to class IV, which showed that the pollution of Cd was serious compared with other metals, ranging from moderate contamination to heavy contamination. The results of I geo had more non-pollution points than those of the pollution index. Considering the high Cd contamination levels in urban areas, it is necessary to identify the sources and carry out quantitative assessments of risks.

Ecological Risk
Given that our study area was polluted by multiple heavy metals, we used the quantitative approach developed by Hakanson to evaluate the potential ecological risk of heavy metal pollution in soil. This approach uses both the physicochemical properties and the toxicological profile of the different heavy metals, which can reflect the comprehensive actions of multiple heavy metals and quantitatively divide the potential ecological risk degree. According to the nature and characteristics of heavy metals in soil, the ecological risk index can not only consider the content of pollutants, but also the toxicity of pollutants, which can reflect the comprehensive action of multiple pollutants and quantitatively divide the potential ecological risk degree. The distribution of RI is shown in Figure 4. The ecological risks of Cr, Mn, Ni, Pb and Zn in all sites were class I. In most sampling sites, Cu showed a low ecological risk, and it showed a moderate risk in only 5% of the sampling sites. The risk of Cd was the highest: 50% of the sites were at high risk, and 33% of the sites were at very high risk. This was because Cd had the highest PI and toxic reaction factor [34]. The comprehensive ecological risk of the industrial region and roads with a high vehicle density was higher. The distribution trend of RI was greatly influenced by Cd, indicating that the Cd contained in the soil presented a serious threat to the local ecology. The Er results were similar for the and PI and I geo , with Cd being the most polluted metal. The spatial distribution of Cd in the soil was mostly concentrated in the industrial zone. Some large industrial enterprises were located in the industrial area, such as iron and steel smelters, coking plants, alloy processing plants, etc. This may explain the high concentration of Cd. Therefore, emission-control restrictions on industrial enterprises, industry capacity removal and the promotion of clean energy are currently the major approaches to alleviate heavy-metal pollution in the Handan urban area.
tively divide the potential ecological risk degree. The distribution of RI is shown in Figure  4. The ecological risks of Cr, Mn, Ni, Pb and Zn in all sites were class I. In most sampling sites, Cu showed a low ecological risk, and it showed a moderate risk in only 5% of the sampling sites. The risk of Cd was the highest: 50% of the sites were at high risk, and 33% of the sites were at very high risk. This was because Cd had the highest PI and toxic reaction factor [34]. The comprehensive ecological risk of the industrial region and roads with a high vehicle density was higher. The distribution trend of RI was greatly influenced by Cd, indicating that the Cd contained in the soil presented a serious threat to the local ecology. The Er results were similar for the and PI and Igeo, with Cd being the most polluted metal. The spatial distribution of Cd in the soil was mostly concentrated in the industrial zone. Some large industrial enterprises were located in the industrial area, such as iron and steel smelters, coking plants, alloy processing plants, etc. This may explain the high concentration of Cd. Therefore, emission-control restrictions on industrial enterprises, industry capacity removal and the promotion of clean energy are currently the major approaches to alleviate heavy-metal pollution in the Handan urban area.   Figure 5 presents the classification of NIPI and PLI. The NIPI showed that 89% of the sites had a heavy pollution level, while the pollution load index showed that 56% of the sites had a medium pollution level, and 34% of the sites had a heavy pollution level. The average value of NIPI in the whole study area is at the heavy pollution level, while the pollution load index in the whole study area is at the medium pollution level. The difference in the evaluation results is mainly because the Nemero pollution index highlights the impact of heavy pollutants on the comprehensive evaluation results, while the pollution load index can better reflect the impact of all pollutants on the comprehensive evaluation results. The Cd pollution in this study is serious, at close to 10 times the background value, while heavy metal pollution is very light. The comprehensive pollution level obtained by Nemerow pollution assessment is higher.

Nemerow Integrated Pollution Index and Pollution Load Index
impact of heavy pollutants on the comprehensive evaluation results, while the pollution load index can better reflect the impact of all pollutants on the comprehensive evaluation results. The Cd pollution in this study is serious, at close to 10 times the background value, while heavy metal pollution is very light. The comprehensive pollution level obtained by Nemerow pollution assessment is higher. 8% 3% I Ⅱ Ⅲ Ⅳ Ⅴ

Health Risk Assessment
The results of the health risk assessment are listed in Table 2. For non-carcinogenic risks, the ingestion risk was far greater than that of the dermal and inhalation pathways. The non-carcinogenic risk for children through ingestion, dermal and inhalation was much higher than that for adults, which was consistent with other study conclusions [8,39]. The non-carcinogenic risk of Cr was the highest, followed by Mn, Pb, Ni, Cd, Cu and Zn. The non-carcinogenic risk of seven kinds of soil heavy metals was below 1, which showed that the non-carcinogenic risk of these heavy metals to citizens can be ignored.

Health Risk Assessment
The results of the health risk assessment are listed in Table 2. For non-carcinogenic risks, the ingestion risk was far greater than that of the dermal and inhalation pathways. The non-carcinogenic risk for children through ingestion, dermal and inhalation was much higher than that for adults, which was consistent with other study conclusions [8,39]. The non-carcinogenic risk of Cr was the highest, followed by Mn, Pb, Ni, Cd, Cu and Zn. The non-carcinogenic risk of seven kinds of soil heavy metals was below 1, which showed that the non-carcinogenic risk of these heavy metals to citizens can be ignored. In terms of carcinogenic risk, the ingestion pathway was dominant. Through ingestion and skin contact, the carcinogenic risk in children was slightly greater than that in adults. In contrast, children have a lower carcinogenic risk through the inhalation pathway than adults. The carcinogenic risk of Pb for all people was less than 1 × 10 −6 , indicating that Pb had no significant carcinogenic risk. The carcinogenic risks of Cd, Cr and Ni for children and adults were between 1 × 10 −6 and 1 × 10 −4 , indicating that these heavy metals had acceptable or tolerable carcinogenic risks. Among the heavy metals that can affect soil quality and reduce plant productivity, Cd is one of the most toxic [46]. Long-term Cd exposure leads to different degrees of DNA methylation in the mother and fetus, and Cd disrupts ATP synthesis by reducing the potential of mitochondrial membranes, causing the abnormal apoptosis of germ cells [47].

Correlation Analysis
The correlation coefficient between elements can indicate similarities in their sources [48]. The results of the Spearman correlation analysis (SCA) of seven metals are listed in Table 3. Based on the SCA, there was a markedly positively correlation between Cr and Ni (p < 0.01, R = 0.634). Thus, they possibly had the same sources. Cr and Ni might be come from the parent materials of soil [12]. There was a significant positive correlation between Zn, Cu, Cd and Pb (p < 0.01). Zn, Cu, Cd and Pb were also correlated (p < 0.05), which means that these metals might come from the same source.

Cluster Analysis
Dendrograms were drawn to display the CA results ( Figure 6), which can vividly reflect the distance between the elements and reveal the relationship between the elements. Three clusters were determined: Ni-Cr; Cd-Pd-Cu; Mn-Zn. Cr and Ni were correlated with each other. Cd and Pb clustered with each other and composed another cluster with Cu. Mn and Zn were isolated and joined to the Cd-Pd-Cu cluster. Metals belonging to the same cluster usually have a common source [49]. Cr and Ni were considered to derive from the parent material of soil. Cd, Pd and Cu might originate from anthropogenic sources. According to the cluster analysis, Mn and Zn were closer to the Cd-Pd-Cu group, suggesting that they were likely from anthropogenic sources.

Source Identification by PCA
The sources of heavy metal pollution in soil were identified by PCA. The Kaiser-Mever-Olkin value (0.608) and Bartlett's test (p < 0.001) both showed that the results obtained by PCA were feasible and reasonable. As shown in Table 4 and Figure 7, three factors were identified, and varimax rotation provided a factor loading that corresponded to the principal components. Eigenvalues greater than one were obtained by PCA, accounting for 78.196% of the total variance.

Source Identification by PCA
The sources of heavy metal pollution in soil were identified by PCA. The Kaiser-Mever-Olkin value (0.608) and Bartlett's test (p < 0.001) both showed that the results obtained by PCA were feasible and reasonable. As shown in Table 4 and Figure 7, three factors were identified, and varimax rotation provided a factor loading that corresponded to the principal components. Eigenvalues greater than one were obtained by PCA, accounting for 78.196% of the total variance. Principal component 1 (PC1) explained 26.830% of the total variance. The Cu, Pb and Cd loads were higher, and Zn also had a medium load. PC1 may be the traffic source [12,50]. The gasoline containing Pb is a significant source of Pb in soil, so Pb is often used to identify traffic sources [51,52]. Cu may come from the vehicle's brake system and radiator [53,54]. Cd might come from diesel fuel and lubricating oil leakage [39,42]. Cd could be released from vehicle tires due to the friction between tires and the road surface. Therefore, PC1 was identified as the vehicle emission.
PC2 explained 26.090% of the total variance, and Cd, Mn and Zn were the main loading elements. Mn and Cd may come from coal combustion [39]. The suspended particles from the steel-making process contain Cd and Zn, which leads to their enrichment in the surrounding soil [55]. Therefore, PC2 might be the industrial source. PC3 was mainly loaded by Cr and Ni, explaining 25.276% of the total variance. According to the previous description, PC3 was the natural source. Principal component 1 (PC1) explained 26.830% of the total variance. The Cu, Pb and Cd loads were higher, and Zn also had a medium load. PC1 may be the traffic source [12,50]. The gasoline containing Pb is a significant source of Pb in soil, so Pb is often used to identify traffic sources [51,52]. Cu may come from the vehicle's brake system and radiator [53,54]. Cd might come from diesel fuel and lubricating oil leakage [39,42]. Cd could be released from vehicle tires due to the friction between tires and the road surface. Therefore, PC1 was identified as the vehicle emission.
PC2 explained 26.090% of the total variance, and Cd, Mn and Zn were the main loading elements. Mn and Cd may come from coal combustion [39]. The suspended particles from the steel-making process contain Cd and Zn, which leads to their enrichment in the surrounding soil [55]. Therefore, PC2 might be the industrial source. PC3 was mainly loaded by Cr and Ni, explaining 25.276% of the total variance. According to the previous description, PC3 was the natural source.

Source Apportionment by PMF Model
The PMF model was applied to determine the source categories and quantitatively calculate their contributions of soil heavy metals. Three sources were apportioned by the PMF model, and their factor profiles are presented in Figure 7. Residue analysis showed that most heavy metal values were between −3 and 3. The fitting results of the observed concentrations and predicted concentrations are shown in Figure S1. The R 2 values were greater than 0.51, indicating that the PMF model obtained reliable results and satisfied the research requirement. With the advancement of receptor modeling, together with the Geographic Information System, a visual analysis of the emission source contribution is now possible. When identifying the sources of heavy metals, a study of the spatial relationship between different land-use types and pollution sources is vital. The spatial distributions of source contributions not only provides more detailed heavy-metal source information, it can also be conducive to obtaining reliable source apportionment results.
Factor 1 accounted for 17.61% of the total contribution rate. Pb, Cd and Cu were the main loading elements, accounting for 80.84%, 66.65% and 51.10% of contributions, respectively. The spatial distribution of factor 1 is presented in Figure 8a. Factor 1 has a great influence on the west of the study area, Airport Road, South Ring Road, Fuhe North Street and Lianfang East Road, a source of traffic. According to the previous description, Pb, Cd, and Cu come from a traffic source, so factor 1 is the traffic source.
calculate their contributions of soil heavy metals. Three sources were apportioned by the PMF model, and their factor profiles are presented in Figure 7. Residue analysis showed that most heavy metal values were between −3 and 3. The fitting results of the observed concentrations and predicted concentrations are shown in Figure S1. The R 2 values were greater than 0.51, indicating that the PMF model obtained reliable results and satisfied the research requirement. With the advancement of receptor modeling, together with the Geographic Information System, a visual analysis of the emission source contribution is now possible. When identifying the sources of heavy metals, a study of the spatial relationship between different land-use types and pollution sources is vital. The spatial distributions of source contributions not only provides more detailed heavy-metal source information, it can also be conducive to obtaining reliable source apportionment results.
Factor 1 accounted for 17.61% of the total contribution rate. Pb, Cd and Cu were the main loading elements, accounting for 80.84%, 66.65% and 51.10% of contributions, respectively. The spatial distribution of factor 1 is presented in Figure 8a. Factor 1 has a great influence on the west of the study area, Airport Road, South Ring Road, Fuhe North Street and Lianfang East Road, a source of traffic. According to the previous description, Pb, Cd, and Cu come from a traffic source, so factor 1 is the traffic source. Factor 2 comprised 28.62% of the total variance. The main loading elements were Cr and Ni, with a contribution of more than 70%. The distribution of factor 2 is presented in Figure 8b. The middle of the study area is highly affected by factor 2, and the content has no relationship with the industrial area and traffic roads. The source of Ni is related to the weathering process, and the source of Cr is parent material [29]. Therefore, factor 2 is identified as the natural source.
Factor 3 was the largest of the three factors and made up 53.77% of the total contributions. Mn and Zn had a high load and accounted for 67.58% and 58.49% of the total, respectively. Cd was moderately loaded. The distribution of factor 3 is presented in Figure  8c, and the high-value area is located in the industrial area with the iron and steel group. Therefore, factor 3 can be classified as the industrial source.
Overall, according to the results of spatial distribution, risk assessment and source apportionment, it is noteworthy that Cd was the heaviest pollutant metal in the study area. The spatial distribution of factor 1 and factor 3 with a high Cd concentration loading decreased from southwest to northeast. This indicated that the heavy metal pollution of soil, especially Cd, was typically associated with traffic activities or discharges from steel ammeter and manufacturing industries. As a typical heavy-industry city, the economic structure of Handan is dominated by heavy industry [56]. The Handan urban area is flatter, with less vegetation and bare soil. It is affected by industrial resuspended dust-bound heavy metals that accumulated in the topsoil. Factor 2 comprised 28.62% of the total variance. The main loading elements were Cr and Ni, with a contribution of more than 70%. The distribution of factor 2 is presented in Figure 8b. The middle of the study area is highly affected by factor 2, and the content has no relationship with the industrial area and traffic roads. The source of Ni is related to the weathering process, and the source of Cr is parent material [29]. Therefore, factor 2 is identified as the natural source.
Factor 3 was the largest of the three factors and made up 53.77% of the total contributions. Mn and Zn had a high load and accounted for 67.58% and 58.49% of the total, respectively. Cd was moderately loaded. The distribution of factor 3 is presented in Figure 8c, and the high-value area is located in the industrial area with the iron and steel group. Therefore, factor 3 can be classified as the industrial source.
Overall, according to the results of spatial distribution, risk assessment and source apportionment, it is noteworthy that Cd was the heaviest pollutant metal in the study area. The spatial distribution of factor 1 and factor 3 with a high Cd concentration loading decreased from southwest to northeast. This indicated that the heavy metal pollution of soil, especially Cd, was typically associated with traffic activities or discharges from steel ammeter and manufacturing industries. As a typical heavy-industry city, the economic structure of Handan is dominated by heavy industry [55]. The Handan urban area is flatter, with less vegetation and bare soil. It is affected by industrial resuspended dust-bound heavy metals that accumulated in the topsoil.
Given this, measures should be undertaken to alleviate the contamination. For instance, the capacity of industrials, including iron, steel, chemical, and heavy engineering industries needs to be reduced. The high capacity of industries has seriously hindered the sustainable development of this city. Furthermore, the promotion of pollution reduction, acceleration of economic structure transformation and promotion of ecological environmental protection are also effective ways to alleviate heavy-metal contamination.

Conclusions
Except for Mn, the average contents of other heavy metals exceeded the background values. The average content of Cd was the highest, at about 10 times the background value.
According to the results of the pollution assessment and ecological risk assessment, Cd pollution was the most serious. In terms of health risk, Cr has the highest non-carcinogenic risk, followed by Mn, Pb, Ni, Cd, Cu and Zn. Cd, Cr and Ni had an acceptable or tolerable carcinogenic risk for children and adults, while Pb had no significant carcinogenic risk. For non-carcinogenic risk and carcinogenic risk, the ingestion risk was far greater than other pathways. The non-carcinogenic risk of children far exceeds that of adults.
The distribution of heavy metals showed that the distribution of Cr and Ni was highly similar, and the high-content regions were in the middle and northeast of the research region. The distribution trends of Mn, Zn and Cd were consistent, and the high contents were concentrated near the industrial area. The contents of Cu and Pb were very high on the main roads, with a high vehicle density.
Combined with correlation analysis, cluster analysis and PCA, three pollution sources were identified, which were traffic sources with Cu, Pb and Cd as main loads, industrial sources with Mn, Cd and Zn as main loads, and natural sources with Cr and Ni as main loads. According to the source apportionment by PMF model, there were three factors for heavy metal accumulation in the study area: traffic source (17.61%), natural source (28.62%) and industrial source (53.77%). The source categories and main load elements were consistently obtained from the PCA and PMF models.
This study can not only help to better understand the spatial distribution of heavy metals, but also offer suggestions for heavy metal contamination control in the study area. However, there were certain limitations in this study. Firstly, for the human-health risk assessment, some exposure and toxicological parameters were based on US EPA, which might not be fully applicable, given the location of this study area in China. Secondly, the ecotoxicities of heavy metal had limited bioavailability. The total concentration of heavy metals was used to assess their health risk, which might cause an overestimation of the actual carcinogenic risk.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/app11209615/s1, Figure S1: Scatter plots of predicted and observed concentrations of species using the PMF model, Table S1: Classification standard of PI, Igeo, NIPI, PLI and PERI, Table S2: Abbreviation and reference values for health risk assessment parameters, Table S3: Reference dose (RfD) and slope factor (SF) of each heavy metal.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author. The data are not publicly available due to participants not having consented to the data being shared.