## Materials and Methods
### Data Evaluation  
#### 2.3.1 Evaluation of Soil Heavy Metal Pollution  
In order to compare the pollution degree of heavy metals in different soil layers of each sampling site, we used single-factor pollution index and Nemerow comprehensive index to calculate the heavy metal pollution degree at each sampling site with the following formulas. \[P_{i}=\frac{C_{i}}{B_{i}} \tag{1}\]  
\[P_{N}=\sqrt{\frac{\overline{P_{i}^{2}+P_{imax}^{2}}}{2}} \tag{2}\]  
In the Formulas (1) and (2), \(P_{i}\) is the single-factor pollution index, \(C_{i}\) is the measured concentration of pollutant \(i\), and \(B_{i}\) is the reference value of pollutant \(i\) in the \"GB15168-2018 Soil Environmental Quality Agricultural Land Pollution Risk Control Standard (Trial)\". \(P_{N}\) is the Nemerow comprehensive index at a particular sampling site, \(P_{i}\) is the average value of the single-factor pollution index, and \(P_{imax}\) is the maximum value of the single-factor pollution index. Figure 1: The distribution of sampling sites. 2.3.2. Soil heavy Metal Inverse Distance Weighting (IDW) and Comprehensive Evaluation  
Previous degrees of heavy metal pollution were based on a single sampling site, which cannot reflect the spatial proximity effect and spatial distribution characteristics of heavy metal indicators at each sampling site. Therefore, this study adopted IDW method to reveal the degree of heavy metal pollution from scattered sampling sites to continuous description and provides a theoretical basis for the pollution zoning of heavy metals. According to the soil quality standards of each functional area, a comprehensive evaluation was carried out. In order to eliminate the influence of the dimension of each heavy metal index, the minimum-maximum normalization method was used to normalize each heavy metal to \([0,1]\) before the comprehensive evaluation. The Raster Calculator tool in ArcGIS was used to comprehensively evaluate the pollution of heavy metals in each layer. In Formula (3), \(x_{max}\) represents the maximum value of the data for the samples, \(x_{min}\) represents the minimum value of data for the samples, \(x\) is the sample data, \(x^{\prime}\) is the normalized sample data. \[x^{\prime}=\frac{x-x_{min}}{x_{max}-x_{min}} \tag{3}\]  
2.3.3. Correlation Analysis  
Pearson correlation analysis was used to determine the relationship among the contents of different heavy metals as well as their relationship with soil physicochemical properties. 2.3.4. Positive Matrix Factorization Model (PMF)  
PMF modeling was performed using the US EPA PMF 5.0 program. As first proposed by [PERSON] et al. [13] in 1994, the least iterative squares algorithm decomposes the original matrix \(E_{ik}\) into two factor matrices \(A_{ij}\) and \(B_{jk}\) and a residual matrix \(C_{ik}\) with the following formula. \[E_{ik}=\sum_{j=1}^{p}A_{ij}B_{jk}+C_{ik}\quad(i=1,2\ldots,m;k=1,2\ldots,n) \tag{4}\]  
In Formula (4), \(E_{ik}\) represents the concentration of the \(k\)th heavy metal in the \(i\)th sample, \(A_{ij}\) represents the concentration (ng mg\({}^{-1}\)) of the \(i\)th sample from the \(k\)th source; \(B_{jk}\) represents the concentration of material from the \(j\)th source contributing to the \(k\)th sample, and \(C_{ik}\) represents the random error. PMF decomposes the original matrix through multiple calculations to obtain the optimal matrices A and B, so that the objective function Q reaches the minimum value. \[Q=\sum_{i=1}^{m}\sum_{k=1}^{n}\left(\frac{C_{ik}}{\sigma_{ik}}\right)^{2} \tag{5}\]  
In the Formula (5): \(\sigma_{ik}\) represents the uncertainty of \(E_{ik}\). During the calculation, both concentration data and the uncertainty data need to be loaded. When the heavy metal concentration is less than or equal to the corresponding method detection limit (MDL), the uncertainty is calculated with the following Formula (6). \[U_{nc}=5/6\times MDL \tag{6}\]  
When the heavy metal concentration is larger than the corresponding MDL, the uncertainty is calculated with the following Formula (7).