磷化工区农田土壤重金属的空间变化与源解析

杨光，唐翔宇，关卓，崔俊芳，程建华

（1.中国科学院、水利部成都山地灾害与环境研究所，成都 610299；2.中国科学院大学，北京 100049；3.浙江农林大学亚热带森林培育国家重点实验室，杭州 311300）

土壤重金属污染具有潜伏性、不可逆性和长期性，且其污染后果严重，不仅能够引起土壤组成、结构和功能的变化，还会抑制农作物生长，造成经济损失[1]。更重要的是，重金属可通过食物链富集后进入人体，严重危害人体健康。镉米、砷毒、血铅等环境公害事件常见诸报道，土壤重金属污染已经成为备受关注的环境问题之一[2]。2014 年公布的《全国土壤污染状况调查公报》显示[3]，全国土壤总调查点位超标率为16.1%，中南、西南地区土壤重金属污染较为突出。其中，重金属污染物Cd、Ni、Cu、Hg、Pb、Cr、Zn 的点位超标率分别为7.0%、4.8%、2.1%、1.6%、1.5%、1.1%、0.9%，以Cd的污染超标点位比例最高。

不同形态的重金属污染物会产生不同的环境效应，从而影响其生物毒性和迁移性[2，4]。辨识土壤中重金属污染物的来源，查明其赋存形态及其空间变化规律，可为土壤重金属污染的针对性治理提供参考。源解析技术是对污染来源进行定性（源识别）或定量（源解析）判断的方法。源解析模型主要分为受体模型和扩散模型两种类型。扩散模型以污染源为研究对象，根据污染物排放量、研究区与排放源距离、污染物的理化性质以及气象条件等参数，应用扩散方程计算源贡献量。但扩散模型无法应用于源强难以确定的开放源[5]。

受体模型主要有富集因子法、同位素标记法、UNMIX 法、主成分/绝对主成分-受体模型（PCA/APCS-MLR）和正定矩阵因子分解法（PMF）等[6]。其中，PMF是由美国环保署（USEPA）组织开发的源解析模型，可在源未知的情况下应用非负因子分析来弱化源解析结果的不确定性，从而提高其精确性[7-8]。与扩散模型相比，受体模型不依赖于污染源排放条件、气象、地形等因素，不用追踪污染物迁移过程，避开了扩散模型的缺点。近年来，已有许多学者应用PMF模型进行土壤重金属源解析研究，并获得了准确可信的结果[9-11]。

成都平原是我国西部工业集中区分布较多的区域，磷矿资源丰富，其开采及生产过程会带来环境重金属污染问题。李杰等[12]调查发现成都、德阳、蒲江彭山3 个典型城市土壤中的重金属Cd 均超过土壤二级质量标准，其中德阳地区超标最为严重；李佳宣等[13]对什邡市磷石膏堆及周围土壤重金属含量的研究发现，长期堆放磷石膏使周围土壤重金属含量升高，其中Cd、Pb 尤为明显。前人研究多集中在土壤、植物和水体污染情况及污染物质重金属含量特征方面，而有关工业集中区周边农田土壤重金属污染评价及其源解析的研究工作相对较少，特别是对土壤Cd赋存形态及其空间变化规律的认识尚有所欠缺。本研究运用地统计学方法分析农田土壤金属污染物（包括不同形态Cd）的空间分布特征，基于污染负荷指数法和潜在生态危害指数法评价其污染现状，采用PMF模型进行来源解析，以期为磷化工区周边农田土壤环境质量的全面评估、重金属污染源的综合治理以及土壤资源的合理利用提供科学依据。

1 材料与方法

1.1 研究区概况

研究区位于成都平原西北部（31°09′N，104°14′E），其概况如图1 所示。该地区属亚热带季风气候，年平均气温16 ℃，年平均降雨量940 mm，夏季盛行风向为东北风。什邡处于龙门山地槽边缘坳陷带中南段之什邡-绵竹复式褶皱带，山区约占全市总面积的60%，耕地约占总面积的30%，土壤以水稻土为主。什邡是我国六大磷矿基地之一，研究区的东北部为磷化工区，产品主要包括黄磷、磷酸以及磷酸盐（化肥、工业材料等）。农业灌溉水主要来自地表水（湔江），施肥为有机肥和化肥综合使用。

图1 研究区采样图Figure 1 Sampling map of study area

1.2 样品采集与测定

在磷化工区的下风向农田区，沿3 个方向（西南偏西、西南、南）在距磷化工区20、500、1 000、1 500、2 000、3 000 m 和4 000 m 处分别采集农田表层（0～20 cm）土壤（图1）。在每个点位，设立3个10 m×10 m 的样方，采样区相距约10 m。在每个样方中，随机采集3个土壤样品，将其制为一个混合样。2019年2月，在21 个点位共采集了63 份土壤样品并制成21 个混合样，用冷藏箱运回实验室。将土壤样品风干后研磨并过10目和200目筛以测定不同指标。

用pH 计（PHS-3C，INESA Instrument，中国）以1∶5（m/V）的土水比测量过10 目筛的土壤pH 值。在分析元素含量之前，用混合强酸（HF、HNO3和HClO4）消煮过200 目筛的土壤样品。样品中主要元素（Al、Ca、Fe、K、Mg、Mn 和P）的含量用电感耦合等离子体原子发射光谱仪（ICP-AES）（Profile DV，利曼-徕伯斯，美国）测定，其他（重）金属元素（Cd、Cr、Cu、Zn、Ni、Pb、Mo、V 和Ti）的含量采用电感耦合等离子体质谱仪（ICP-MS）（7700x，安捷伦，美国）测定。

1.3 研究方法

1.3.1 污染指数法

单因子污染指数法是对单种污染物的污染程度进行评价的方法，计算公式如下[14-15]：

式中：Pi为单因子污染指数；Ci为实测值；Bi为背景值，选择四川省成都经济区土壤背景值作为标准[16]。每种污染物的污染指数分为无污染（Pi＜1）、轻度污染（1≤Pi＜3）、中度污染（3≤Pi＜6）和重度污染（Pi≥6）。

污染负荷指数法是反映综合污染程度的一种方法，用于描述单一位点多种污染物的综合污染程度，计算公式如下[17-18]：

式中：PLI为污染负荷指数；n为污染物个数。污染负荷指数评价标准：PLI≤1，无污染；1＜PLI≤2，轻度污染；2＜PLI≤3，中度污染；PLI＞3，重度污染。

1.3.2 潜在生态风险指数法

Hakason 提出的潜在生态风险指数（Potential ecological risk index）是计算土壤污染物潜在风险的有效定量方法[15，19]。该方法综合考虑了土壤中污染物的含量和毒性水平，可评估一种或多种污染物的潜在生态危害程度，被广泛应用于沉积物及土壤重金属污染评价[20-21]，计算公式如下：

式中：RI为土壤污染物的潜在生态风险指数，根据RI值，可将土壤潜在风险程度划分为轻微（RI＜150）、中度（150≤RI＜300）、强度（300≤RI＜600）、很强（RI≥600）4 个等级；Eir为第i种污染物的潜在生态风险系数，根据Eir值将风险级别划分为轻微（Eir＜40）、中度（40≤Eir＜80）、强度（80≤Eir＜160）、很强（160≤Eir＜320）、极强（Eir≥320）5 个等级；Ti为第i种污染物的毒性系数，参考相关研究[15，22]，各种（重）金属元素的毒性系数分别取Mn=Zn=1、V=Cr=2、Ni=Pb=Cu=5、Mo=15、Cd=30；Ci为污染物i的实测含量值；Cin为污染物i的背景值，取值与1.3.1中背景值相同。

1.3.3 土壤重金属形态分析

通过初步分析确定该地区污染最严重的重金属元素，并对其进行形态分析，土壤重金属形态分析采用Tessier五步连续提取法[23]，该方法将土壤重金属赋存形态分为可交换态、碳酸盐结合态、铁锰氧化物结合态、有机结合态和残渣态。采用ICP-MS 测定各形态含量。

1.3.4 PMF模型

式中：xij表示样本i中第j个污染物的含量；gik为样本i中第k个污染源的贡献率；fkj为污染源k对第j个污染物的特征值。

eij由目标函数Q的最小值计算：

式中：uij是i样本中第j个污染物的不确定度。

PMF模型引入了不确定度参与建模，若实测含量大于方法检出限（MDL），则采用式（6），若实测含量小于MDL，则采用式（7）。

式中：σ是相对标准偏差。

本研究利用PMF 5.0模型分析土壤中重金属的来源，对采样数据进行20 次迭代处理，寻找最小最稳定的Q（robust），以确定因子（源成分），求得其贡献率。

2 结果与讨论

2.1 农田土壤重金属含量的空间变化规律

研究区所有农田采样点的土壤pH 范围为5.23～6.90，平均值为6.04，变异系数为8.57%，呈弱酸性。与表1 中四川省土壤背景值相比，研究区土壤Cr、Ni、Mo、Cd、Pb、V、Zn含量的平均值超过背景值，其中，Cd和Zn 的含量分别达到背景值的16 倍和1.8 倍，且所有采样点的含量均超过背景值，表现出明显的积累，其他金属元素的含量略超出背景值，对土壤环境质量影响较小。研究区所有采样点土壤Cd 含量均高于《土壤环境质量农用地土壤污染风险管控标准（试行）》（GB 1561—2018）的风险筛选值（0.3 mg·kg-1），最高超标13.5 倍，且部分点位超出风险管控值[3 mg·kg-1（6.5＜pH≤7.5）]；有4 个采样点Zn 超标（200 mg·kg-1），超标率为19%，最高超标1.62倍；其他种类重金属的含量均不超标。

表1 农田土壤（重）金属统计特征Table 1 Statistical characteristics of（heavy）metals in farmland soil

一般可用变异系数反映污染物的平均变异程度，若其大于75%，则表明污染物的空间分布不均，有外源输入[24]。如表1 所示，研究区土壤Cd 和Zn 的变异系数分别达93.51%和64.27%，表明其受人为干扰影响较大。其他金属元素含量的变异程度较小（均不超过30%），说明它们的空间分布较为均匀，受外源的影响较小。由图2可见，Cd含量随与磷化工区距离的增大而显著降低。综上，研究区主要受到重金属Cd 的污染，且污染程度较重，受到的Zn污染较轻。可以推断，农田土壤的Cd 和Zn 污染可能主要与当地磷化工生产有关。此外，该研究区土壤并未受到其他金属元素的污染。

目前防治病毒病尚无特效药物，应预防为主，根据现有的药物配合施用，尽量抑制病毒病发生。可用20%病毒A500倍液、抗病毒1号400倍液，40%植物病毒灵1000倍液。结合叶面喷施牛奶、葡萄糖、含P、K、Zn的叶面肥，增加抗病力，既能防病又丰产。

图2 与磷化工区不同距离采样点土壤的重金属含量（3个等距采样点的平均值±标准差）Figure 2 Contents of heavy metals in soil samples collected at different distances from the phosphorus chemical industrial zone（average value±standard deviation for three sampling locations at the same distance）

2.2 农田土壤Cd赋存形态的空间变化规律

如图3a 所示，主要污染物Cd 的不同赋存形态的占比随与磷化工区距离的增大而发生规律性变化。在距离≤500 m 的采样点土壤中，交换态Cd 和碳酸盐结合态Cd 之和占比超过80%，不同形态Cd 占比大小顺序为交换态Cd＞碳酸盐结合态Cd＞铁锰氧化物结合态Cd＞残渣态Cd＞有机结合态Cd。交换态Cd 和碳酸盐结合态Cd 与土壤结合能力较弱，在酸性环境下易于释放，有较高的生物活性，在土壤中迁移性较强，易被作物富集，进而危及人体健康[25]。在距离＞500 m的采样点土壤中，交换态Cd 和碳酸盐结合态Cd 之和的占比大幅低于≤500 m 范围内的采样点土壤，且其随距离增大而降低的趋势变缓并在距离4 000 m处达到最低值（33%）。铁锰氧化物结合态Cd 是中等活性态Cd，在还原条件下易于释放，在距离≤500 m 的采样点土壤中，其占比约为10%，在距离＞500 m 的采样点土壤中其占比明显高于≤500 m 范围内的采样点土壤，且其含量随距离的增大基本保持不变（介于23%～24%）。有机结合态Cd 和残渣态Cd 在土壤中稳定性较好，是作物难以利用的形态，对生态系统危害小[26]。有机结合态Cd 和残渣态Cd 的占比均随与磷化工区距离的增大而上升并在4 000 m 处达到最高值，两者占比之和达43%。在距离磷化工区最远的轻污染采样点土壤中，不同形态Cd占比大小顺序为残渣态Cd＞铁锰氧化物结合态Cd＞交换态Cd＞碳酸盐结合态Cd＞有机结合态Cd，与邻近磷化工区的重污染采样点土壤的Cd赋存特征有明显不同。

图3 与磷化工区不同距离采样点土壤的Cd赋存形态Figure 3 Speciations of Cd in soil samples collected at different distances from the phosphorus chemical industrial zone

图3 b 所示土壤中不同赋存形态Cd 的含量随与磷化工区距离的变化规律与图3a 所示其占比的变化规律相似。在距离为20 m 处，土壤中交换态Cd 含量最高，达2 mg·kg-1，且其随着距离增大而急剧降低，在4 000 m处达到最低值（0.10 mg·kg-1）。碳酸盐结合态Cd 和铁锰氧化物结合态Cd 的变化趋势与交换态Cd相似，其中碳酸盐结合态Cd 最大值为0.27 mg·kg-1，而最小值为0.09 mg·kg-1，铁锰氧化物结合态Cd 最大值为0.24 mg·kg-1，而最小值为0.12 mg·kg-1。有机结合态Cd 以及残渣态Cd 与前3 种Cd 形态的变化趋势相反，总体表现为随着距离增大而逐渐增加。有机结合态Cd 和残渣态Cd 的最大值分别为0.05 mg·kg-1和0.21 mg·kg-1。

综上，磷化工区的存在使研究区土壤中主要重金属污染物Cd 的分布呈现明显的空间变化趋势，这与已有文献的报道相一致[27-28]。外源输入是研究区农田土壤中Cd 污染的主要原因，可能的外源主要包括工业活动、交通、化肥施用以及污水灌溉等。研究区位于磷化工区的下风向，因此会受到污水灌溉和污染降尘输入的影响。此外，该地区磷矿资源丰富，土壤磷素的地质背景含量高，存在伴矿重金属（Zn、Pb、Cu、Co、Ni、Cd）[29]。邻近工业区采样点土壤的高Cd含量则可能主要与磷化工生产排污有密切关系。

2.3 土壤污染及生态风险评价

2.3.1 污染负荷指数法评价

研究区农田土壤单因子污染指数平均值的排序从大到小依次为Cd（9.86）、Zn（2.17）、Pb（1.21）、Mo（1.17）、V（1.12）、Cu（1.08）、Ni（1.05）、Cr（0.93）、Mn（0.62）。依据该评价方法，重金属Cd 属于重度污染，Zn、V、Ni、Mo、Pb、Cu 属于轻度污染，Cr、Mn 属于无污染。由表2 可知，所有采样点土壤中Cd 的单因子污染指数属于中度污染和重度污染，少数采样点土壤中Zn 属于中度污染和重度污染，其他金属属于轻度污染或无污染。依据污染负荷指数，所有采样点土壤均属于轻度污染。

表2 不同污染级别采样点的数量占比（%）Table 2 Proportions of the number of sampling point at different pollution levels（%）

由图4 可见，土壤中主要污染元素Cd 在与磷化工区距离20 m 和500 m 处的单因子污染指数的平均值分别为39.96 和32.59，其最大值分别高达51.13 和51.26，随着距离的增大逐渐降低至7.33，但依然属于中度污染级别。其他金属元素的单因子污染指数也有相似的空间变化规律。类似地，综合污染指数也呈现此变化规律，其最大值的最高值出现在500 m 处（1.87），其平均值随与工业区距离增大而逐渐降低，在4 000 m处降为1.23。

图4 与磷化工区不同距离处采样点土壤的单因子污染指数和综合污染指数Figure 4 Single factor pollution indexes and comprehensive pollution indexes of soil samples collected at different distances from the phosphorus chemical industrial area

2.3.2 潜在生态危害指数法评价

研究区农田土壤的单项潜在风险系数平均值的排序由大到小依次为Cd、Mo、Pb、Cu、Ni、V、Zn、Cr、Mn。其中，Cd 的单项潜在生态风险指数均大于80，有7 个采样点处于很强生态风险级别，有4 个点位处于极强生态风险级别；其他金属元素的单项潜在风险指数均小于40，处于轻微生态风险级别（表3）。研究区仅有2 个采样点的土壤重金属综合潜在生态风险指数小于150，其中71.43%的采样点存在中等生态风险，19.05%的采样点存在很强生态风险。由表4 可知，潜在生态风险指数也与采样点位置存在明显的关系。Cd 的单项潜在风险指数在与磷化工区距离500 m 范围内最高，指数超过800 的采样点均在这个范围内；在距离＞500 m 的采样点，指数均小于300，在距离4 000 m 处，指数低至150 以下。综合潜在生态风险指数也呈现相似的空间变化规律，指数超过800的点位均在距磷化工区500 m 范围内，随着距离增大该值逐渐降低，在4 000 m处时，降低至200以下。此外，s3和s8点位所在线上Cd的潜在生态风险指数变化不同于其他两条，s3 和s8 所在的点位线位于农田旁，该线上土壤为进行农业活动而进行过一定的修复与改良，所以与另两条线上土壤潜在生态风险指数变化有差别，但整体上也遵循随距离增加Cd潜在生态风险逐渐降低的规律。

表3 不同生态风险级别采样点的数量占比（%）Table 3 Proportions of the number of sampling point at different ecological risk levels（%）

表4 与磷化工区不同距离处采样点土壤的潜在风险指数Table 4 Potential risk indexes of soil samples collected at different distances from the phosphorus chemical industrial area

综上，应用单因子污染指数和单项潜在生态风险指数评价研究区农田土壤污染状况时，均发现Cd 污染严重，而对于其他金属元素，不同评价方法的结果有一定差异。依据单因子污染指数，个别采样点Zn污染较为严重；依据单项潜在生态风险指数，该研究区所有采样点的Zn污染处于轻微生态风险。依据污染负荷指数，整个研究区均处于轻度污染；然而，依据潜在生态风险指数，该研究区存在中等到很强级别的潜在生态风险。这是因为污染负荷指数法作为一种化学评价方法，重点突出采样点含量最高的污染物对环境的影响，而潜在生态危害指数法则凸显高毒性污染物的生态环境风险。前者考察的是污染状况和背景之间的差异，而后者通过引入生态毒理系数表征不同污染物对生态系统的潜在影响。将两种方法组合应用可获得更全面的评价结果。

2.4 土壤污染源解析

Pearson 相关性分析发现，研究区土壤中多种元素之间存在相关性，结果如图5 所示，Cd与Zn、Mg呈显著正相关（P＜0.01），相关系数分别为0.60和0.57；Ni与K、Cr、V 呈显著正相关（P＜0.01），相关系数分别为0.56、0.91、0.54；Cr 与Pb、Fe、V 呈显著正相关（P＜0.01），相关系数分别为0.63、0.58、0.62；Pb与Fe、Cu、V呈极显著正相关（P＜0.001），相关系数分别为0.67、0.68、0.68；Fe 与V、Al、Ti 呈极显著正相关（P＜0.001），相关系数分别为0.88、0.72、0.71；Ca与V、Al、Ti呈显著负相关（P＜0.01），相关系数分别为-0.62、-0.80、-0.74。具有显著相关性的元素可能来源相同，Cd 与Zn 是该区超标的重金属且相关性较强，说明土壤中两种元素可能来自相同的外部输入源；Al、Fe、Ti、V 等元素相关性较高，且它们都是地壳元素，说明它们可能都是土壤内源产生的。

图5 土壤污染物Pearson相关性分析结果Figure 5 Results of Pearson correlation analysis for soil pollutants

通过PMF模型分析，最终确定4因子，其Q（robust）值最小，元素信噪比（S/N）均大于5，大部分残差处于-3～3之间，结果较理想，源成分贡献率如图6所示。

因子3对地壳元素Al、Fe、Ti、V的贡献率较高，分别为37.8%、35.1%、38.0%、37.1%，对重金属元素Cu、Pb、Ni、Cr、Zn的贡献率也较高，分别为43.4%、43.0%、33.0%、34.8%、36.0%，且Al、Fe、Ti、V 的变异系数较低，分别为14.46%、15.02%、6.84%、8.58%，表明这些元素在研究区分布较为均匀，外源贡献较低。因此，因子3被识别为土壤母质源。

因子2 对Cd 和Zn 的贡献率最高，分别为56.2%和45.9%，对其他元素的贡献率都不超过20%，说明该因子可作为Cd 和Zn 的识别因子。该结果与前述发现相一致：变异系数结果显示土壤中的Cd 和Zn 主要由外来源引入（表1），且Cd含量随与磷化工区距离的增大而快速下降（图2）。该区域已有工业废弃物调查研究发现，磷渣和沟底淤泥中Cd 和Zn 的含量均超标，其中污泥的Cd 含量高达123 mg·kg-1[30]。该地区磷化工产业较多，产生大量磷石膏副产物，磷石膏中含有较高含量的Cd、Zn、Cu、Pb等[31]，其暴露在自然界中受溶出、雨水淋滤、风力搬运等作用影响会不断富集在土壤和水体中，且有研究发现磷石膏中Cd 的存在形态主要为交换态，生态风险较大[32]。因此，确定因子2为磷化工源。

因子4对肥料元素Ca、Mg、K、P的贡献率较高，分别达49.8%、41.8%、33.7%、48.2%，这些元素是重要的肥力因子，主要来源于化肥，长期施用会造成其在局部土壤中的含量较高。此外，该因子对重金属Cu 和Zn也有一定的贡献率（分别为25.6%和18.1%）。化肥种类繁多，有研究表明，部分肥料中重金属Cu和Zn的含量较高，土壤中的Cu 和Zn 与农业化学品投入关系密切[33]。因此，因子4确定为农业化学投入品源。

因子1 对Cd 和Pb 的贡献率较高，分别为43.8%和22.3%。采样区土壤位于工业区下风向，重金属污染物经过大气迁移后易沉降在附近土壤，有研究表明，该区域大气降尘的Cd含量超标，且部分点位超标严重[28]，说明大气沉降是Cd 的一个重要输入源。近年来，随着经济发展，私家车数量激增，带来的尾气排放也是逐年增加。Pb是交通运输的主要标志物[34]，农田采样点周边道路汽车尾气所含的Pb，会随着降尘飘落至表层土壤中。成都平原是我国西部工业集中区分布较多的区域之一，研究区周边工业生产活动产生的含Pb 灰尘会随着干湿沉降进入土壤表层[35]。因此，判定因子1为大气沉降和尾气排放源。

3 结论

（1）磷化工区下风向农田主要重金属污染物为Cd 和Zn，以Cd 的超标较为严重，距离磷化工区越近，Cd超标倍数越高，主要是外源输入的结果。

（2）土壤中Cd 的赋存形态随距离增大呈现空间变化规律，残渣态Cd 含量总体增加，而交换态Cd 含量大幅下降。

（3）污染负荷指数法和潜在生态危害指数法评估表明，靠近磷化工区（距离≤500 m）的采样点土壤的Cd污染级别及生态风险均相对较高。

（4）PMF 源解析结果表明，对于污染最重的Cd，工业污染源的外源转入贡献率为56.2%，大气沉降和尾气排放源的外源转入贡献率为43.8%；对于污染较轻的Zn，土壤母质源的内源贡献率为36.0%，工业污染源的外源转入贡献率为45.9%，农业化学投入品源的外源转入贡献率为18.1%。因此，有必要对当地涉重工业企业进行有效监管，加强农业环境的综合治理，以确保粮食生产安全。