Distribution, Origins and Hazardous Effects of Polycyclic Aromatic Hydrocarbons in Topsoil Surrounding Oil Fields: A Case Study on the Loess Plateau, China

The Loess Plateau has one of the most vulnerable ecological environments in the world, but it also contains abundant oil and gas resources that are regularly exploited, which has resulted in serious environmental problems. Therefore, it is important to analyze the polycyclic aromatic hydrocarbons (PAHs) present in the topsoil of this region. The ∑16PAHs concentrations between 1980–1999 and 2000–2019 ranged from 1134.20–15871.04 and 1010.67–18,068.80 µg kg−1, with average values of 5021.30 and 5662.82 µg kg−1. All samples displayed heavy pollution levels according to European soil quality standards. In addition, among the measured physicochemical properties, the soil organic carbon (SOC) had the greatest influence on PAHs, while soil particle size distribution had the smallest effect. Source apportionment indicated that the two main sources were petroleum source (37.57%) and vehicular traffic source (25.88%). Lastly, an assessment of the carcinogenic risks illustrated that more focus should be placed on the dermal pathway in which the human body is exposed to soil PAHs. Overall, the carcinogenic risks in different populations did not exceed 10−4, but there was still a potential carcinogenic risk in some age groups, especially in adult women.


Introduction
Polycyclic aromatic hydrocarbons (PAHs) are a major class of persistent organic pollutants (POPs) that are primarily generated by human activities, including the incomplete combustion of fossil fuels and biomass fuels, automobile exhaust, and oil spills [1][2][3][4][5]. They are widely distributed throughout the atmosphere, soil, water, sediment, and other environmental media. Among these, soil is the most important reservoir of PAHs, and is also a stable index used to reflect the status of environmental pollution [6]. Among the hundreds of known PAHs, the United States Environmental Protection Agency (USEPA) has identified sixteen as priority pollutants (16PAHs), seven of which are classified as carcinogenic PAHs (7PAHs) [7,8]. Due to their toxicity, mutagenicity, and carcinogenicity, the origin, distribution, and hazards of 16PAHs in soil must be analyzed [9][10][11].

Soil Sampling and Preparation
The sample collection was undertaken in July 2017 during various trips. The places for sampling belonged to the regions of Yan'an, Yulin, and Qingyang producing oils ( Figure 1). All soil samples were taken near oil wells, but away from the apparent oil flows into the soil. A stainless-steel soil auger was used to collect 5 different samples at a depth of 0-10 cm soil layer and a composite sample was achieved by mixing the 5 samples. The above sampling process was repeated three times, that is, there were three composite samples at each sample point. After being taken to the laboratory, grounding of the samples was carried out followed by sieving with a sieve of stainless-steel with 60 mesh sieve and storage at 4 • C prior to analysis.

Analysis of Soil Physicochemical Properties
A dichromate-based oxidation method was used for measuring the soil organic carbon (SOC) content of the samples. The pH was measured in a suspension (1:2.5, soil/water) using a pH meter (Lei-ci PXSJ-216F, lei-ci, Shanghai, China). In the wet measurement mode, the soil texture was analyzed by laser diffraction using a Mastersizer 2000 (MS-2000, Malvern Panalytical, Malvern, UK).

Analysis of Soil Physicochemical Properties
A dichromate-based oxidation method was used for measuring the soil organic carbon (SOC) content of the samples. The pH was measured in a suspension (1:2.5, soil/water) using a pH meter (Lei-ci PXSJ-216F, lei-ci, Shanghai, China). In the wet measurement mode, the soil texture was analyzed by laser diffraction using a Mastersizer 2000 (MS-2000, Malvern Panalytical, Malvern, UK).

Extraction and Analysis of PAHs
The steps of analysis of PAHs mainly included dichloromethane extraction, QuEChERS reagent purification, and organic filter membrane filtration. The detailed experimental methods were described previously by Wang et al. [12]. Gas chromatography-tandem mass spectrometry (GCMS-TQ8040, Shimadzu, Kyoto, Japan) in combination with the QuEChERS method was used to determine 16PAHs. Chromatographic separation was actualized with a capillary column Rxi-5Sil Ms (30 m × 0.25 mm × 0.25 µm). He (carrier gas) flow rate was 1.0 mL min −1 with temperature programming as follows: increase to 50 • C with a 2 min waiting time followed by an increase at a rate of 20 • C min −1 up to 250 • C with a 3 min waiting time, and a final increase to 300 • C at a heating rate of 5 • C min −1 with a 5 min waiting time. Identification of PAHs was implemented by comparison of retention time and ion characteristics with standards (Supplementary Table S1).

Detection of PAHs
The external standard method having seven points was used for 16PAHs described by a linear equation having R 2 > 0.996 (Supplementary Table S2). The limit of detection (LOD) was found to be thrice to that of the standard deviation of the blank. The LOD of PAHs was in the range of 0.02-0.80 µg kg −1 . Any concentrations below the LOD were defined as non-detected (N.D.). The 16PAHs had an average recovery of 65-119% with 0.5-9.5% being the range of relative standard deviation (RSD). The recoveries and RSDs were found to be in accordance with the experimental standard (Supplementary Table S3).

Data Analysis
2.4.1. Source Apportionment PMF (Positive Matrix Factorization) is an analysis tool of the source apportionment of pollutants recommended by USEPA [14]. The tool is applicable for different environmental media such as soil/sediment, water, and air samples. The PMF model is aimed at calculating the lowest objective function Q for 'm' number of sources containing 'n' samples in terms of the residuals (e ij ) and uncertainty (u ij ): where, for measurement of a source 'j' in sample 'i' in the presence of 'p' number of factors, u ij is the uncertainty estimate, the contribution of each factor towards the individual sample is g ik , and the profile of species in each source is given by f kj .
For each sample, the uncertainties were a combination of uncertainty in measurement and detection limits of the method (MDL). For samples having a concentration equal to or below MDL, the uncertainty U nc was taken as a fixed fraction of the MDL: but for higher sample concentrations (>MDL), U nc was based on our provided fraction of the concentration and MDL: Starting from different points and choosing 3-5 variables, the PMF model was run after selecting a random seed mode having 33 starting points chosen randomly. The number of factors corresponding to PMF was determined after comparison of the values Q true and Q robust, where Q true is a measure of closeness of the input data to the fit parameter and Q robust is achieved after excluding the outliers [15]. The Q robust value (2165.9) calculated by PMF was close to the Q true value (2167.2). Therefore, the analytical results of PMF were reasonable. Hence, optimal solution based on four factors indicated a strong relation between observed and predicted PAH concentrations (r 2 > 0.90, p < 0.01). Three different types of PAHs were identified by the PMF model: strong (FLU, ANT, FLA, PYR, and BaA), weak (NAP, ACE, PHE, CHR, BbF, BkF, BaP, InP, DBA, and BgP), and bad species (ACY).

Health Risk Assessment
According to the appraisal of health risk, different PAHs have different toxicities [2]. In this work, the direct exposure risk of PAHs towards cancer was investigated using the incremental lifetime cancer risk (ILCR). This model based on ILCR was used in combination with a toxic equivalent method based on toxicity equivalency factors (TEFs). High carcinogenicity of BaP influenced its choice as a reference compound in the TEF method [16]. The ILCRs included ingestion, dermal contact, and inhalation [17]. They were calculated through the following formula: where PAH i denotes the average soil concentration of a single PAH (mg kg −1 ), TEF i is the toxic equivalency factor for an individual PAH in relation to BaP, and CS denotes the summation of converted PAH concentrations based on TEFs. CSF denotes the carcinogenic slope factor (mg kg −1 day −1 ) −1 , where CSF ingestion , CSF dermal , and CSF inhalation were assigned the values of 7.3, 25, and 3.85 (mg kg −1 day −1 ) −1 , respectively [18]. BW refers to the body weight (kg), EF refers to the exposure frequency (d year −1 ), ED is the duration of exposure (year), and InhR and IngR denote the rate of inhalation (m 3 d −1 ) and the rate of soil ingestion (mg d −1 ), respectively. SA corresponds to the surface area of dermal exposure (cm 2 day −1 ), the dermal adherence factor is referred to AF (mg cm −2 ), AT denotes the average life span (day), ABS is the dermal adsorption factor, and PEF corresponds to the particle emission factor (m 3 kg −1 ). People were divided in terms of age and gender in six groups, such as: (1) adult male (18-70 years age), (2) adult female (18-70 years age), (3) adolescent male (11-17 years age), (4) adolescent female (11-17 years age), (5) male child (2-10 years age), and (6) female child (2-10 years age) [19]. The corresponding parameter values are given in Supplementary Table S4.

Concentrations and Compositions of PAHs
The concentration profiles of PAHs around oil wells of different initial production time are shown in Table 1 , respectively. In general, the concentration of HMWPAHs is generally higher than those in LMWPAHs because of the more significant degradation profile of LMWPAHs [20]. Moreover, the relatively high resistance of HMWPAHs to degradation contributes to their accumulation in soil [21]. It is evident from Table 1 that the pollution of soil PAHs in 2000-2019 is more serious than that in 1980-1999. One reason is that the utilization rate of oil wells is low in the first two decades (some oil wells have been abandoned), and the pollution degree of soil PAHs has been alleviated after a long period of natural degradation and remediation. The other reason is that the utilization rate and oil recovery rate of newer oil wells are relatively high, so the leakage of oil and the emission of exhaust gas from oil production vehicles will aggravate the pollution of soil PAHs in the process of oil exploitation.
The soil PAHs pollution was separated into four grades with reference to the soil quality standard proposed by Maliszewska-Kordybach [22]. The total concentration of 16PAHs below 200 µg kg −1 indicates non-polluted, with lightly polluted being indicated by a concentration of 200-600 µg kg −1 , and a concentration of 600-1000 µg kg −1 representing moderately polluted. When the concentration exceeds 1000 µg kg −1 , the area is heavily polluted. Compared to the above classification standards, all sites were heavily contaminated ( Figure 2).
According to the number of aromatic rings, the 16PAHs were classified into three different groups, one with 2-3 rings, another with 4 rings, and last one with 5-6 rings [23]. The compositions of soil PAHs around oil wells with different production times are shown in Figure 3. In 1980-1999, the proportion of PAHs having 2-3 rings, 4 rings, and 5-6 rings in this region was 28.71−63.06% (mean 47.80%), 30.67-49.48% (mean 40.35%), and 4.24-21.81% (mean 11.86%), respectively. Similar results were found in 2000-2019, the highest ratio value of soil PAHs was found in 2-3 rings (26.20-65.33%, mean 43.97%), followed by 4 rings (17.67-54.38%, mean 39.57%), and 5-6 rings (4.72-51.95%, mean 16.47%). It can be found that the proportion of LMWPAHs is large, which are mainly related to combustion products of pyrogenic processes at low temperatures [24]. However, HMWPAHs are more likely to accumulate in soil because of their low volatility, and their sources are different from those of LMWPAHs. Generally, LMWPAHs are produced by substances like petroleum and its derivatives, while HMWPAHs mainly come from the combustion of coal, organic compounds, and higher plants [25]. The specific sources of PAHs will be discussed below. ∑16PAHs-summation of concentrations of sixteen PAHs. ∑7PAHs-summation of concentrations of seven carcinogenic PAHs including BaA, CHR, BbF, BkF, BaP, DBA, and InP. ∑LMWPAHs-the concentration sum of low molecular weight PAHs (2-3 rings). ∑HMWPAHs-the total concentrations of high molecular weight PAHs (≥4 rings). All the other abbreviations used are explained in Supplementary Table S1. N.D.-Not Detected.

Relationships between PAHs and Environmental Factors
As shown in Figure 4, PCA provided a single two-dimensional model, which would explain 72.53% of the variance in the data. It can be found that there is a high correlation between individual PAHs. Reflection of calculated Pearson correlation coefficients for the individual PAHs showed interrelation of most variables at level 0.01 (Supplementary Table S5).
Commonly, soil properties such as pH, SOC, and texture, are important factors affecting the concentrations of PAHs [26]. In our study, the vector magnitude of environmental factors indicated that SOC was the most significant factor. The soils with high SOC had a strong adsorption ability to PAHs, which affected the distribution and migration of PAHs in the environment. Poor correlation existed between pH and PAHs, indicating that pH was a non-critical parameter affecting PAHs. For soil particle size distribution, we found that the influence of sand, silt, and clay on the distribution of PAHs was basically the same. The difference was that most individual PAHs showed positive correlations with soil clay, whereas negative correlations with soil sand, suggesting that PAHs were more likely to accumulate in fine soil particles. The reason is that small particles revealed stronger affinities of PAHs than larger particles due to the higher surface area of smaller particles.

Relationships between PAHs and Environmental Factors
As shown in Figure 4, PCA provided a single two-dimensional model, which would explain 72.53% of the variance in the data. It can be found that there is a high correlation between individual PAHs. Reflection of calculated Pearson correlation coefficients for the individual PAHs showed interrelation of most variables at level 0.01 (Supplementary Table S5).
Commonly, soil properties such as pH, SOC, and texture, are important factors affecting the concentrations of PAHs [26]. In our study, the vector magnitude of environmental factors indicated that SOC was the most significant factor. The soils with high SOC had a strong adsorption ability to PAHs, which affected the distribution and migration of PAHs in the environment. Poor correlation existed between pH and PAHs, indicating that pH was a non-critical parameter affecting PAHs. For soil particle size distribution, we found that the influence of sand, silt, and clay on the distribution of PAHs was basically the same. The difference was that most individual PAHs showed positive correlations with soil clay, whereas negative correlations with soil sand, suggesting that PAHs were more likely to accumulate in fine soil particles. The reason is that small particles revealed stronger affinities of PAHs than larger particles due to the higher surface area of smaller particles.

Source Analysis of PAHs
As this study was carried out using samples from locations around oil wells, we can relate the main sources of 16PAHs with oil exploitation and human-related activities. But, it is also required to take into account other sources. For further identification of PAHs sources, modeling of the data was carried out by PMF [27]. Figure 5a presents the average contributions of PAHs species towards four PMF factors.
A major portion (37.57%) of the total measured PAHs was due to the factor 1, which was loaded predominately by NAP and ACE. Among all the different PAHs, NAP was the heaviest and was usually used as an indicator for the leakage of petroleum substances [28][29][30], whereas ACE was mainly from the combustion of petroleum-based fuel at lower temperatures. Earlier investigations have indicated that origins of LMWPAHs are petroleum products related to petrogenic sources [6,31,32]. Therefore, factor 1 was marked as related to petroleum sources. ANT, PHE, and FLU, which are believed to be the tracers of wood/grass burning, dominated factor 2, which accounted for 18.75% of the total PAHs [33][34][35][36][37]. Hence, factor 2 was identified as originating from burning of biomass. Factor 3 explained 17.80% of the total PAHs and heavily loaded on BaA, followed by ACE, ANT, PYR, NAP, InP, and BaP. Among them, BaA was not only the characteristic indicator of natural gas, but also the product of coal combustion [38,39]. Other PAHs were related to the release from burning of coal and coke combustion [40,41]. Thus, factor 3 was related to the combustion of coal/coke. Factor 4 explained 25.88% variance of the data and was high loaded predominately by molecular weight PAHs with a high loading on InP and BaP and moderately by BgP and DBA. InP and BgP were characteristic indicators of diesel and gasoline combustion, respectively [42]. Meanwhile, particulate matters of gasoline vehicle exhaust showed high emissions of InP and DBA [43,44]. Consequently, factor 4 was related to the source of vehicular traffic.
The PMF results showed the following order for relative contributions to the total soil PAHs burden: petroleum source (37.57%) > vehicular traffic source (25.88%) > biomass combustion source (18.75%) > coal/coke combustion source (17.80%) (Figure 5b). The rationalization proposal is to pay attention to the leakage and sprinkling of petroleum in the process of petroleum exploitation.

Source Analysis of PAHs
As this study was carried out using samples from locations around oil wells, we can relate the main sources of 16PAHs with oil exploitation and human-related activities. But, it is also required to take into account other sources. For further identification of PAHs sources, modeling of the data was carried out by PMF [27]. Figure 5a presents the average contributions of PAHs species towards four PMF factors.
A major portion (37.57%) of the total measured PAHs was due to the factor 1, which was loaded predominately by NAP and ACE. Among all the different PAHs, NAP was the heaviest and was usually used as an indicator for the leakage of petroleum substances [28][29][30], whereas ACE was mainly from the combustion of petroleum-based fuel at lower temperatures. Earlier investigations have indicated that origins of LMWPAHs are petroleum products related to petrogenic sources [6,31,32]. Therefore, factor 1 was marked as related to petroleum sources. ANT, PHE, and FLU, which are believed to be the tracers of wood/grass burning, dominated factor 2, which accounted for 18.75% of the total PAHs [33][34][35][36][37]. Hence, factor 2 was identified as originating from burning of biomass. Factor 3 explained 17.80% of the total PAHs and heavily loaded on BaA, followed by ACE, ANT, PYR, NAP, InP, and BaP. Among them, BaA was not only the characteristic indicator of natural gas, but also the product of coal combustion [38,39]. Other PAHs were related to the release from burning of coal and coke combustion [40,41]. Thus, factor 3 was related to the combustion of coal/coke. Factor 4 explained 25.88% variance of the data and was high loaded predominately by molecular weight PAHs with a high loading on InP and BaP and moderately by BgP and DBA. InP and BgP were characteristic indicators of diesel and gasoline combustion, respectively [42]. Meanwhile, particulate matters of gasoline vehicle exhaust showed high emissions of InP and DBA [43,44]. Consequently, factor 4 was related to the source of vehicular traffic.
The PMF results showed the following order for relative contributions to the total soil PAHs burden: petroleum source (37.57%) > vehicular traffic source (25.88%) > biomass combustion source (18.75%) > coal/coke combustion source (17.80%) (Figure 5b). The rationalization proposal is to pay attention to the leakage and sprinkling of petroleum in the process of petroleum exploitation. Moreover, car drivers are advised to replace gasoline and diesel with cleaner fuels, or to use electric vehicles, which can reduce the burden of PAHs in the soil environment. Moreover, car drivers are advised to replace gasoline and diesel with cleaner fuels, or to use electric vehicles, which can reduce the burden of PAHs in the soil environment.  Supplementary Table S1.

Assessment of Health Risk
The results of lifetime cancer risk levels are given in Table 2. Normally, values of ILCRs indicate the magnitude of risk with negligible risks, potential health risks, and higher risks being denoted by the values of 10 −6 , between 10 −6 to 10 −4 , and greater than 10 −4 , respectively [45]. In the current study, it is clear that all soil samples had the same level of cancer risk (10 −12 -10 −10 ) through inhalation in all soil

Assessment of Health Risk
The results of lifetime cancer risk levels are given in Table 2. Normally, values of ILCRs indicate the magnitude of risk with negligible risks, potential health risks, and higher risks being denoted by the values of 10 −6 , between 10 −6 to 10 −4 , and greater than 10 −4 , respectively [45]. In the current study, it is clear that all soil samples had the same level of cancer risk (10 −12 -10 −10 ) through inhalation in all soil samples, indicating that there is no cancer risk during different periods. By contrast, the levels of carcinogenic risk to the human body are large through ingestion and dermal contact pathways. For the ingestion pathway, the average ILCR ingestion (1980)(1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999) for male and female were estimated to be 8.54 × 10 −7 , 5.09 × 10 −7 , and 9.26 × 10 −7 , and 8.78 × 10 −7 , 5.26 × 10 −7 , and 1.01 × 10 −6 , for children, adolescents, adults respectively, suggesting that adults were most sensitive to PAHs contamination. In addition, adults have the highest risk (mean value: 1.65 × 10 −6 -1.79 × 10 −6 ) of exposure to soil PAHs through the dermal pathway, followed by adolescents (mean value: 1.27 × 10 −6 -1.31 × 10 −6 ) and children (mean value: 1.06 × 10 −6 -1.09 × 10 −6 ). Accurately, the carcinogenic risk of male and female in children and adolescents were basically the same, while in adults, the carcinogenic risk in women was higher than that in men. The higher inhalation rate with larger dermal exposure area along with exposure duration for adults increases carcinogenic risk. Moreover, the weight difference between men and women is responsible for higher carcinogenic risk in women. These results indicated that the dermal exposure of the human body towards soil PAHs is the most important in this region, causing great concern to us. In 2000-2019, the toxin effect of PAHs from soil to all the groups through direct ingestion, dermal contact, and inhalation was basically the same as that in 1980-1999. However, the mean value of ILCRs in 2000-2019 is higher than that in 1980-1999 owing to the rapid demand for petroleum during the recent two decades, resulting in the soil being seriously polluted. So, administrators in petroleum production areas on the Loess Plateau need to take some steps to reduce the risks posed by PAHs.

Conclusions
Understanding the distribution, origins, and potential health hazards of 16PAHs in soil around oil wells on the Loess Plateau is critical to appropriately manage PAHs in the soil matrice. The pollution of soil PAHs in 2000-2019 is more serious than that in 1980-1999. Among the measured PAHs, the ratio of soil PAHs decreased in the following order with respect to the number of rings: 2-3 > 4 > 5-6. SOC was the most significant factor, followed by pH and soil particle size distribution. Four different sources were identified by PMF. Source apportionment revealed that origins of PAHs were mainly from petroleum and vehicular traffic. Biomass and coal/coke combustion also contribute to the PAHs recorded in this area. The ILCRs results indicated that the exposure pathway of dermal contact was more noteworthy than that of ingestion and inhalation. Moreover, the health risk of adults is higher than that of adolescents and children.
Supplementary Materials: The following are available online at http://www.mdpi.com/1660-4601/17/4/1390/s1, Table S1: Retention time and ion characteristics of the selected PAHs in the MRM mode; Table S2: Equations of calibration curve of 16PAHs; Table S3: Performance and validation of the analytical method (n = 6); Table S4: Parameters used in the incremental lifetime cancer risk assessment; Table S5: Correlation analysis between PAHs and environmental factors.