Spatial Distribution and Source Apportionment of Soil Heavy Metals in Pearl River Delta, China

: The contents of ten heavy metals (HMs) (Cu, Pb, Zn, Cd, Ni, Cr, Hg, As, Co, and Mn) in 80 surface soil samples (0–20 cm) were investigated in the Pearl River Delta (PRD), Guangdong Province, China. The average contents of Cu, Pb, Zn, Cd, Ni, Cr, Hg, As, Co, and Mn were 16.45, 40.20, 45.10, 0.09, 12.93, 47.93, 0.13, 14.44, 5.68, and 199.66 mg/kg, respectively. The soil quality was generally good, though slightly higher levels (1.17, 1.61, 1.67, and 1.62 times) of soil Pb, Cd, Hg, and As contents were observed compared with the soil background values. The spatial distribution of soil HM pollution in the PRD showed that 36% of sample sites were evaluated as sites without soil pollution, 32% as sites with slight pollution, 20% as sites with nearly moderate pollution, 9% as sites with moderate pollution, and 3% as site with serious pollution. Source apportionment analysis showed that the source of 64.33% of soil HMs in the PRD could be explained by natural and industrial sources, 24.80% by transportation, and 10.87% by agricultural activities.


Introduction
To achieve sustainable development, it is necessary to collect basic environmental information for resources such as soil and to investigate the contamination status of heavy metals (HMs) in rapidly developing economic regions. Soil HMs have received long-term attention due to the persistent threat they pose to regional ecosystems [1]. The excessive accumulation of HMs in the soil at high rates in the long-term is an irreversible process that reduces soil quality and ultimately damages agricultural lands [2]. Natural factors, such as type of geological formation, can lead to soil HM pollutions [3]. Agricultural and industrial activities have intensified the input of HMs to soil ecosystems [4]. HMs, such as Cu, Zn, Ni, Cd, Pb, Cr, Hg, and As, not only affect the quality of crops, but also accumulate in the human body through food-chain accumulation and harm public health, resulting in most countries listing them as priority pollutants to control [5][6][7]. The abovementioned eight HMs have been listed as priority pollutants by the United Environmental Protection agency, based on a series of hazards. However, previous studied have found that inadequate Co in the human body can lead to anemia and anorexia, and excessive Co intake can cause pancreatic failure, cardiomyopathy, different types of cancer, hypothyroidism, and bone marrow proliferation [8]. The excessive intake of Mn in childhood can cause a neurotoxic syndrome that impedes dopamine balance and behavior control [9]. It had also been reported that Ni, Cr, and Mn were important factors leading to endemic diseases in the

Study Area
The studied sites included municipal and rural districts in the PRD with a total area of 5.5 × 10 4 km 2 . It consisted of nine cities, including Guangzhou, Shenzhen, Dongguan, Foshan, Zhuhai, Zhaoqing, Huizhou, Zhongshan, and Jiangmen ( Figure 2). Previous studies have pointed out that the PRD has gradually become the region experience the most rapid change of the ecological environment and the most serious resource depletion and environmental pollution in China, especially as regards soil HMs pollution [21], which could increase the risk of HM contamination in food chains [22]. The climate is a typical subtropical monsoon climate, with an average annual temperature of 21.9 °C and an average annual rainfall of 2200 mm. According to the statistics of the soil census office of Guangdong province, the dominant soil types can be divided into seven categorizations: paddy soil, red soil, accumulated soil, lateritic soil, yellow soil, lime soil, and salinized soil. The main types of parent rocks are granite, quaternary soil, sand shale, and metamorphic rock. The climate is hot and humid, with leaching and sedimentation ongoing features within the soil of the PRD. As a result, the soil profile is continuously depleted of alkali metals, alkali earth metals, and soluble salts, enriched with hydrogen ions and ironaluminum oxides. Under such conditions, the soil is acidic in the PRD [23].

Study Area
The studied sites included municipal and rural districts in the PRD with a total area of 5.5 × 10 4 km 2 . It consisted of nine cities, including Guangzhou, Shenzhen, Dongguan, Foshan, Zhuhai, Zhaoqing, Huizhou, Zhongshan, and Jiangmen ( Figure 2). Previous studies have pointed out that the PRD has gradually become the region experience the most rapid change of the ecological environment and the most serious resource depletion and environmental pollution in China, especially as regards soil HMs pollution [21], which could increase the risk of HM contamination in food chains [22]. The climate is a typical subtropical monsoon climate, with an average annual temperature of 21.9 • C and an average annual rainfall of 2200 mm. According to the statistics of the soil census office of Guangdong province, the dominant soil types can be divided into seven categorizations: paddy soil, red soil, accumulated soil, lateritic soil, yellow soil, lime soil, and salinized soil. The main types of parent rocks are granite, quaternary soil, sand shale, and metamorphic rock. The climate is hot and humid, with leaching and sedimentation ongoing features within the soil of the PRD. As a result, the soil profile is continuously depleted of alkali metals, alkali earth metals, and soluble salts, enriched with hydrogen ions and ironaluminum oxides. Under such conditions, the soil is acidic in the PRD [23].

Soil Sample Treatment and Analysis
According to the topographic and hydrogeological conditions in the PRD, a total of 80 sampling sites (0-20 cm) were collected, with a sampling density of about 1 site per 30 km 2 and an area of about 20 m 2 for each randomly selected sampling site. GPS was used to locate the geographic locations of the sampling sites ( Figure 2). Two surface soil samples were collected at each sampling site for geochemical analysis (Cu, Pb, Zn, Cd, Ni, Cr, Hg, As, Co, and Mn). Mixed soil samples after collection were naturally air-dried in a lab, rocks, plant branches and leaves, grass roots and other sundries were then removed, and the samples were properly mashed with wooden sticks. Then, an appropriate amount of the soil samples was taken by the quartering method and grounded with a ceramic mortar. Finally, soil samples were placed into a sealed bottle after being passed through a 100-mesh sieve.

Soil Sample Treatment and Analysis
According to the topographic and hydrogeological conditions in the PRD, a total of 80 sampling sites (0-20 cm) were collected, with a sampling density of about 1 site per 30 km 2 and an area of about 20 m 2 for each randomly selected sampling site. GPS was used to locate the geographic locations of the sampling sites ( Figure 2). Two surface soil samples were collected at each sampling site for geochemical analysis (Cu, Pb, Zn, Cd, Ni, Cr, Hg, As, Co, and Mn). Mixed soil samples after collection were naturally air-dried in a lab, rocks, plant branches and leaves, grass roots and other sundries were then removed, and the samples were properly mashed with wooden sticks. Then, an appropriate amount of the soil samples was taken by the quartering method and grounded with a ceramic mortar. Finally, soil samples were placed into a sealed bottle after being passed through a 100mesh sieve.
The soil pH values were measured using a pH meter (ST3100, OHAUS Instruments Co., Changzhou, China) in 1:2.5 (w/v) deionized water. Soil Cu, Pb, Zn, Cd, Ni, Cr, Co, and Mn was detected using an Inductively Coupled Plasma-Atomic Emission Spectrometer (ICP-AES, ICPE-9800 Shimadzu Co., Kyoto, Japan) after digesting with HNO3-HF-HClO4. The elements Hg and As were determined by atomic fluorescence spectrometry (AFS230E Haiguang Analytical Instrument Co., Beijing, China). Each sample was analyzed in duplicate and all the reagents were of analytical grade. Standard quality assurance and quality control measures were followed during the analysis of the HMs. The recovery rates of the heavy metals ranged from 84 to 120%, and the results of the experiments are consistent with quality control standards (Table S1).

Evaluation of Soil HM Pollution
The geo-accumulation index ( ) is a widely used evaluation index for soil HM pollution, which can reflect not only the influence of a single HM in a specific environment, but also the composite influence of multiple HMs. According to the , the pollution degree rating is divided into seven levels (Table S2) [24]. The specific calculation equation used, Equation (1), is shown below: where, is the actual concentration of an element ; represents the geochemical background value of an element , taking the geochemical background values of element in the Guangdong province as the fiducial values [25]. is the coefficient set for changes in background values that may be caused by diagenesis (generally 1.5). The soil pH values were measured using a pH meter (ST3100, OHAUS Instruments Co., Changzhou, China) in 1:2.5 (w/v) deionized water. Soil Cu, Pb, Zn, Cd, Ni, Cr, Co, and Mn was detected using an Inductively Coupled Plasma-Atomic Emission Spectrometer (ICP-AES, ICPE-9800 Shimadzu Co., Kyoto, Japan) after digesting with HNO 3 -HF-HClO 4 . The elements Hg and As were determined by atomic fluorescence spectrometry (AFS230E Haiguang Analytical Instrument Co., Beijing, China). Each sample was analyzed in duplicate and all the reagents were of analytical grade. Standard quality assurance and quality control measures were followed during the analysis of the HMs. The recovery rates of the heavy metals ranged from 84 to 120%, and the results of the experiments are consistent with quality control standards (Table S1).

Evaluation of Soil HM Pollution
The geo-accumulation index (I geo ) is a widely used evaluation index for soil HM pollution, which can reflect not only the influence of a single HM in a specific environment, but also the composite influence of multiple HMs. According to the I geo , the pollution degree rating is divided into seven levels (Table S2) [24]. The specific calculation equation used, Equation (1), is shown below: (1) where, C n is the actual concentration of an element n; B n represents the geochemical background value of an element n, taking the geochemical background values of element n in the Guangdong province as the fiducial values [25]. k is the coefficient set for changes in background values that may be caused by diagenesis (generally 1.5).

PCA-MLR Model
The PCA-MLR model is based on feature analysis and is used to extract the factor load matrix F and factor score matrix G of the receptor samples. The factor load is used to identify pollution sources, and the factor score is used to calculate source contribution rates [26,27]. The basic principle of PCA analysis is based on the least square method, that is, to seek F and G satisfying the minimum variance, which is described by the matrix as Equations (2) and (3): The multiple linear regression equation is described as Equation (4): where, Y is the total concentration of the pollutants, P is the number of factors extracted by the PCA, a h is the standardized regression coefficient of the h factor, X h is the score of the h factor extracted by PCA, and b is the regression constant term. The independent and dependent variables were standardized and regression analysis was performed by the standardized Equation (5).
where, B i is the multiple linear regression coefficient. The average contribution rate of source i (%) was calculated by Equation (6):

Data Processing
The descriptive statistics and principal component analysis-multiple linear regression (PCA-MLR) of the soil HMs concentrations were conducted by IBM SPSS Statistics 25 software (IBM, Chicago, IL, USA). Spatial distribution maps were accomplished by ArcGIS 10.2, and the statistical analysis by Origin 2021 software (OriginLab, Northampton, MA, USA).

Descriptive Statistics of HMs in PRD
The average contents of Cu, Pb, Zn, Cd, Ni, Cr, Hg, As, Co, and Mn were 16.45, 40.20, 45.10, 0.09, 12.93, 47.93, 0.13, 14.44, 5.68, and 199.66 mg/kg, respectively (Table 1). Compared with the background soil HM values of the Guangdong province, soil Pb, Cd, Hg, and As contents were 1.17, 1.61, 1.67, and 1.62 times higher, and the other soil HM contents were slightly lower in the PRD. Additionally, ten HMs were present in content values that were less than the control values stipulated by soil environmental quality standards for agricultural land in China [28]. This is due to that fact that this region has experienced the development process of pollution, reclamation, and improvement, with the environmental degradation index showing a downward trend within the PRD region during 2011-2015, which is a positive outcome of the efforts taken by the government to manage pollution to ensure that the environmental situation has been continuously improved [29]. The variation coefficients of soil Cd, Hg, and As exceeded 100%, indicating an obvious spatial heterogeneity which may be as a result of human activities [30]. Soil Pb, Cd, Hg, and As showed positive skewness, reaching values of 3.23, 4.87, 2.98 and 3.86, respectively, which may have been caused by the input of external sources [31]. The mean value of soil pH value was 5.10, indicating that the soil is acidic, as is consistent with previous studies [23].

Evaluation of Soil HM Pollution
Taking the soil background values of the study area as reference values, the geoaccumulation index (I geo ) was used to evaluate the soil HM pollution situation in the PRD. The spatial distribution of soil HMs showed that 36% sample sites were evaluated as sites without soil pollution, 32% as sites with slight pollution, 20% as sites with nearly moderate pollution, 9% as sites with moderate pollution, and 3% as site with serious pollution (Figure 3). Serious HM pollution was found in Foshan and south of Guangzhou city, with Cd and Hg being the main contributors, which was consistent with previous studies [32]. As a result, we evaluated the soil HM pollution in the PRD to be slight or non-polluted, with low environmental risk.

Spatial Distribution Characteristics of Soil HMs
The spatial distribution of soil HMs in the PRD are shown in Figure 4. It is shown that soil Cu, Pb, Zn, and Ni contents were much higher in Guangzhou and southern Foshan city compared with other regions, while soil Cd, Cr, Hg, and As contents were distributed evenly in the whole of the PRD except for the relatively higher values that appeared in Guangzhou, Dongguan, Foshan, Zhaoqing, and Huizhou city, which may have resulted from intensive industrial activities in these regions. Another reason may be the influence of natural geological factors, including the low mountains and hills that dominate in the PRD, causing HMs to be easily deposited in the valleys and basins of these areas.

Spatial Distribution Characteristics of Soil HMs
The spatial distribution of soil HMs in the PRD are shown in Figure 4. It is shown that soil Cu, Pb, Zn, and Ni contents were much higher in Guangzhou and southern Foshan city compared with other regions, while soil Cd, Cr, Hg, and As contents were distributed evenly in the whole of the PRD except for the relatively higher values that appeared in Guangzhou, Dongguan, Foshan, Zhaoqing, and Huizhou city, which may have resulted from intensive industrial activities in these regions. Another reason may be the influence of natural geological factors, including the low mountains and hills that dominate in the PRD, causing HMs to be easily deposited in the valleys and basins of these areas.

Pearson Correlation Analysis
Pearson correlation analysis showed that the soil element contents of the ten HMs could be explained by natural and anthropogenic factors [33,34]. We could see that the correlation between different HMs in the soil were strong, but some were not obvious, indicating that the ten HMs are partly caused by natural conditions, and partly, as a preliminary speculation, by long-term human activities ( Figure S1). The soil Cu content was found to be significantly positively correlated with Zn, Ni, and Cr (p < 0.01); Pb significantly positively correlated with Zn (p < 0.01); and Zn significantly positively correlated with Cu and Ni (p < 0.01). The soil As content had little relevance to the other nine types of HMs, indicating that it may have come from a different source.

Spatial Distribution Characteristics of Soil HMs
The spatial distribution of soil HMs in the PRD are shown in Figure 4. It is shown that soil Cu, Pb, Zn, and Ni contents were much higher in Guangzhou and southern Foshan city compared with other regions, while soil Cd, Cr, Hg, and As contents were distributed evenly in the whole of the PRD except for the relatively higher values that appeared in Guangzhou, Dongguan, Foshan, Zhaoqing, and Huizhou city, which may have resulted from intensive industrial activities in these regions. Another reason may be the influence of natural geological factors, including the low mountains and hills that dominate in the PRD, causing HMs to be easily deposited in the valleys and basins of these areas.

Source Apportionment by PCA-MLR Model
In order to further explore the sources of soil HMs in the PRD, principal component analysis (PCA) was used to reveal the potential relationships among different soil element contents. The Kaiser normalized orthogonal rotation method was used during the PCA analysis process, with Kaiser-Meyer-Olkin (KMO) (0.811 > 0.5) and Bartlett spherical tests (0.000 < 0.05) used. Three principal components whose eigenvalues were greater than 1 after rotation were obtained ( Figure 5, Table S3). The results indicated that the PCA led to a reduction of the initial dimension of the dataset to three components, explaining 72.32% of the variation in the data, the largest apportion of these components. It could be seen from the factor load after rotation (Table S4)  indicating that the ten HMs are partly caused by natural conditions, and partly, as a preliminary speculation, by long-term human activities ( Figure S1). The soil Cu content was found to be significantly positively correlated with Zn, Ni, and Cr (p < 0.01); Pb significantly positively correlated with Zn (p < 0.01); and Zn significantly positively correlated with Cu and Ni (p < 0.01). The soil As content had little relevance to the other nine types of HMs, indicating that it may have come from a different source.

Source Apportionment by PCA-MLR Model
In order to further explore the sources of soil HMs in the PRD, principal component analysis (PCA) was used to reveal the potential relationships among different soil element contents. The Kaiser normalized orthogonal rotation method was used during the PCA analysis process, with Kaiser-Meyer-Olkin (KMO) (0.811 > 0.5) and Bartlett spherical tests (0.000 < 0.05) used. Three principal components whose eigenvalues were greater than 1 after rotation were obtained ( Figure 5, Table S3). The results indicated that the PCA led to a reduction of the initial dimension of the dataset to three components, explaining 72.32% of the variation in the data, the largest apportion of these components. It could be seen from the factor load after rotation (Table S4)   The source apportionment of HMs could be inferred by combining PCA analysis with the spatial distribution characteristics. PC1 reflected an enrichment of Cu, Zn, Ni, Cr, The source apportionment of HMs could be inferred by combining PCA analysis with the spatial distribution characteristics. PC1 reflected an enrichment of Cu, Zn, Ni, Cr, Co, Mn, and a small amount of Cd. Compared with the soil background values of the Guangdong province, the soil Cu, Zn, Ni, Cr, Co, and Mn contents were all slightly lower except for the soil Cd, indicating that these HMs mainly from natural sources, potentially through the slow mineral soil formation processes of the parent rock [35,36]. Granite was widely distributed in the PRD, leading to quaternary soil with loose accumulation [23]. Previous studies have shown that limestone and sand shale are Mn-rich parent rocks [37,38]. The specific humid climate conditions present in the PRD accelerate the weathering process of parent rocks, and thus HMs are easily released from this natural source. On the other hand, as industrial sources may also account for the presence of soil HMs, we collected the spatial distribution information of 1904 enterprises associated with the surface treatment of metals and manufacturing in the PRD (Figure S2), based on which the kernel density of these enterprises was expressed ( Figure 6). This showed that the spatial distribution of enterprises associated with the surface treatment of metals and manufacturing, as well as their kernel density, was close to the Pearl River, leading to a potential accumulation site of soil HMs. Further, the spatial distribution patterns of soil Cu, Zn, Ni, Co, and Mn contents were highly correlated to the enterprises associated with the surface treatment of metals and manufacturing. It us suggested that these enterprises might be another important source apportionment of soil Cu, Zn, Ni, and Cr contents caused by atmospheric deposition [39,40]. Consequently, it is reasonable to conclude that PC1 was mainly a mixture of natural and industrial sources.
atmospheric deposition [39,40]. Consequently, it is reasonable to conclude that PC1 was mainly a mixture of natural and industrial sources.
PC2 was predominated by Pb, Hg, Cd, and Zn. The soil Pb and Hg may be correlated to the developed transportation industry in the PRD, which is an important automobile exhaust and production center [41][42][43]. A certain amount of soil Pb and Cd may be the result of waste gas emissions from vehicles, as well as the combustion of Pb-containing gasoline [4]. However, Cd occupied a large load proportion in both PC1 and PC2, reaching values of 0.405 and 0.400, respectively, indicating that there exists other Cd sources aside from natural and industrial sources. The geo-accumulation index results showed that Cd, Hg, and Pb were the main contributors to soil pollution in the PRD. As mentioned above, soil Cd and Pb might mainly come from vehicle emissions, oil leakages, cement pavement wear, rubber tires, and brake pad wear [44,45]. The larger contribution of PC2 to Pb occurred in Guangzhou city, which has the most intensive population and traffic [46]. As a result, we deduced that PC2 was traffic sources.  PC2 was predominated by Pb, Hg, Cd, and Zn. The soil Pb and Hg may be correlated to the developed transportation industry in the PRD, which is an important automobile exhaust and production center [41][42][43]. A certain amount of soil Pb and Cd may be the result of waste gas emissions from vehicles, as well as the combustion of Pb-containing gasoline [4]. However, Cd occupied a large load proportion in both PC1 and PC2, reaching values of 0.405 and 0.400, respectively, indicating that there exists other Cd sources aside from natural and industrial sources. The geo-accumulation index results showed that Cd, Hg, and Pb were the main contributors to soil pollution in the PRD. As mentioned above, soil Cd and Pb might mainly come from vehicle emissions, oil leakages, cement pavement wear, rubber tires, and brake pad wear [44,45]. The larger contribution of PC2 to Pb occurred in Guangzhou city, which has the most intensive population and traffic [46]. As a result, we deduced that PC2 was traffic sources. PC3 was occupied primarily by As. The average soil As content was significantly higher than that of the background value. Previous studies have considered livestock manure and some phosphate fertilizers as potential soil As sources [47,48]. Inorganic As may originated from agricultural activities, such as fertilizers and pesticides [49,50]. Phosphorus fertilizer is a commonly used fertilizer for agricultural activity due to soil P deficiency under acidic conditions in the PRD, which leads to the enrichment of As in the soil [51]. As-containing pesticides and feed additives had been added to the farmland historically in the PRD [52][53][54][55], which has potentially created excess As accumulation [56]. Therefore, we deduced that PC3 was significantly correlated to agricultural activities (fertilizers, pesticides, and livestock).
Overall, the results of this study showed that soil Cu, Zn, Ni, Cr, Co, Mn, and some amounts of Cd could be derived from natural and industrial sources, while soil Pb, Hg, Cd, and Zn originated from traffic sources, with As primarily originating from agricultural sources. This is in agreement with the previous literature carried out at the county-level of the Guangdong province, where source appointment results showed that Pb, Zn, and Cu mainly originated from vehicle emissions and atmosphere deposition; Hg and Cd originated from industrial activities; Cr and Ni mainly came from soil parent materials; and As mainly originated from agricultural inputs [13]. Similarly, Soffianian et al. (2015) applied geostatistical methods and also determined the effect of agricultural activities on HM concentrations and the spatial distribution using the Geographical Information System (GIS). The results showed that As, Cd, Zn, and Pb had a geological or agricultural origin and that Cr, Co, Ni, and V originated from bedrocks [57].
In order to quantify the contributions of the various sources obtained by PCA in the PRD, mixed logistic regression (MLR) of the sources (natural, industrial, transportation, and agricultural activities) was performed using Equation (6) [58]. The calculation results are showed (Figure 7, Table S5). These indicated that the contribution ratios of four sources occurred in the following order: natural and industrial sources (64.33%) > transportation sources (24.80%) > agricultural sources (10.87%).
Overall, the results of this study showed that soil Cu, Zn, Ni, Cr, Co, Mn, and some amounts of Cd could be derived from natural and industrial sources, while soil Pb, Hg, Cd, and Zn originated from traffic sources, with As primarily originating from agricultural sources. This is in agreement with the previous literature carried out at the countylevel of the Guangdong province, where source appointment results showed that Pb, Zn, and Cu mainly originated from vehicle emissions and atmosphere deposition; Hg and Cd originated from industrial activities; Cr and Ni mainly came from soil parent materials; and As mainly originated from agricultural inputs [13]. Similarly, Soffianian et al. (2015) applied geostatistical methods and also determined the effect of agricultural activities on HM concentrations and the spatial distribution using the Geographical Information System (GIS). The results showed that As, Cd, Zn, and Pb had a geological or agricultural origin and that Cr, Co, Ni, and V originated from bedrocks [57].
In order to quantify the contributions of the various sources obtained by PCA in the PRD, mixed logistic regression (MLR) of the sources (natural, industrial, transportation, and agricultural activities) was performed using Equation (6) [58]. The calculation results are showed (Figure 7, Table S5). These indicated that the contribution ratios of four sources occurred in the following order: natural and industrial sources (64.33%) > transportation sources (24.80%) > agricultural sources (10.87%). After, with the factor score variable as an independent variable and the standardized total HM contents as the dependent variable, Equation (7) was used (Table S6) where, Z is the total soil HM contents after data standardization and is the factor score variables by source apportionment .
To further obtain the contribution of each sampling site to the soil HM contents, Equation (8) was used to quantitively evaluate the source apportionment contribution from sampling site [59]. After, with the factor score variable as an independent variable and the standardized total HM contents as the dependent variable, Equation (7) was used (Table S6): where, Z is the total soil HM contents after data standardization and FS i is the factor score variables by source apportionment i.
To further obtain the contribution of each sampling site to the soil HM contents, Equation (8) was used to quantitively evaluate the source apportionment contribution from sampling site i [59].
where, Mean ∑ THM is the total mean HM contents and σ THM is standard deviation. The contributions of 80 sampling sites from the four sources are shown ( Figure 8, Table S7). The results suggest that the contributions of the different sampling sites to PC1 varied greatly, consequently suggesting that PC1 was attributed to by more than one source. This had again proved that PC1 has two different sources. However, no obvious fluctuations among different sampling sites occurred, suggesting that the sampling sites had contributed little to PC2 and PC3, indicating that just one source exists for both PC2 and PC3. This was consistent with the results of the source apportionment by the PCA-MLR model. It can be seen that the contribution of some individual sampling site appears to have a negative value, which is caused by nonnegativity restrictions in the solution process of the PCA-MLR model. fluctuations among different sampling sites occurred, suggesting that the sampli had contributed little to PC2 and PC3, indicating that just one source exists for bo and PC3. This was consistent with the results of the source apportionment by th MLR model. It can be seen that the contribution of some individual sampling site a to have a negative value, which is caused by nonnegativity restrictions in the s process of the PCA-MLR model.

Conclusions
In this study, 80 surface soil samples were collected in the PRD to analyze H tents as well as to trace possible HM sources through the use of the PCA-MLR mod soil environmental quality in the PRD was found to be generally good, though higher Pb, Cd, Hg, and As contents were found in the surface soil compared w background values in the Guangdong province. For the spatial distribution of so in the PRD, soil Cu, Pb, Zn, and Ni contents were much higher in Guangzhou and ern Foshan city, while soil Cd, Cr, Hg, and As contents were distributed evenly whole of the PRD except for the relatively higher values that appeared in Guan Dongguan, Foshan, Zhaoqing, and Huizhou city. The PCA-MLR model indicated soil HMs in the PRD mainly come from four sources in the following order: natu industrial sources (64.33%) > transportation sources (24.80%) > agricultural (10.87%). Soil contents of Cu, Zn, Ni, Cr, Co, Mn, and some Cd could be mainly ex by natural and industrial sources, whereas Pb, Hg, Cd, and Zn may originate from cles emissions, and As mainly originates from agricultural activities. This study p basic information to environmental managers and policy-makers on how to con gional soil HM pollution risks. Further investigation is still needed to build an a relationship between soil HM pollution and the rapid economic development leve PRD.

Conclusions
In this study, 80 surface soil samples were collected in the PRD to analyze HM contents as well as to trace possible HM sources through the use of the PCA-MLR model. The soil environmental quality in the PRD was found to be generally good, though slightly higher Pb, Cd, Hg, and As contents were found in the surface soil compared with the background values in the Guangdong province. For the spatial distribution of soil HMs in the PRD, soil Cu, Pb, Zn, and Ni contents were much higher in Guangzhou and southern Foshan city, while soil Cd, Cr, Hg, and As contents were distributed evenly in the whole of the PRD except for the relatively higher values that appeared in Guangzhou, Dongguan, Foshan, Zhaoqing, and Huizhou city. The PCA-MLR model indicated that the soil HMs in the PRD mainly come from four sources in the following order: natural and industrial sources (64.33%) > transportation sources (24.80%) > agricultural sources (10.87%). Soil contents of Cu, Zn, Ni, Cr, Co, Mn, and some Cd could be mainly explained by natural and industrial sources, whereas Pb, Hg, Cd, and Zn may originate from vehicles emissions, and As mainly originates from agricultural activities. This study provides basic information to environmental managers and policy-makers on how to control regional soil HM pollution risks. Further investigation is still needed to build an accurate relationship between soil HM pollution and the rapid economic development level in the PRD.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/su13179651/s1, Figure S1: Pearson correlation analysis of HMs in soil; Figure S2: Spatial distribution of 1904 enterprises associated with surface treatment of metals and manufacturing in the PRD; Table S1: Statistic description information of the quality control of soil HMs; Table S2: Geo-accumulation index and pollution levels; Table S3: Total variance of HMs explained by PCA; Table S4: Matrix of principal component analysis; Table S5: Contribution composition of soil HMs in the PRD; Table S6: The factor score variables of source i at 80 sampling sites; Table S7  Data Availability Statement: The data supporting the reported results in the present study will be available on request from the corresponding author or the first author.

Conflicts of Interest:
The authors declare no conflict of interest.