PMF模型的影响因素考察——以某铅锌矿周边农田土壤重金属源解析为例

魏迎辉，李国琛，王颜红，张 琪，李 波，王世成，3，崔杰华，张 红，周 强

（1.中国科学院沈阳应用生态研究所，沈阳 110016；2.中国科学院大学，北京 100049；3.辽宁省农产品环境及质量安全控制工程技术研究中心，沈阳 110016）

土壤重金属污染不仅直接改变土壤理化性状、降低土壤生物活性、阻碍养分有效供应，从而降低农产品品质，更重要的是通过食物链数十倍富集[1-2]，严重威胁人类健康[3-4]。2016年国家印发《土壤污染防治行动计划》，加大土壤防治力度。寻找并阻断污染源能从根本上治理土壤污染防治[5]。土壤重金属可分为自然源和人为源：自然源主要为成土母质，人为源主要有大气沉降、农业活动、交通污染以及工矿业活动等[6-7]。土壤重金属往往受多个污染源共同作用，且来源不明，各污染源贡献高低无法量化表征。基于同位素比值的溯源，要求必须采集端源物质并测定其同位素比值[8]，因而并不适用于污染源未知的情况。未知污染源解析主要有两个步骤：源分类和贡献率计算[9-10]。主成分分析（Principal component analysis，PCA）和因子分析（Factor analysis，FA）只能对污染源进行分类，无法直接给出完整的源贡献率[11]。正定矩阵因子分析法（Positive matrix factorization，PMF）和UNMIX受体模型不仅可以给出污染源个数、每个污染源贡献的元素，还可以给出各污染源对某种元素的贡献率[11-12]。

PMF是由Paatero在1994年首次提出，被美国环保署认可用于大气污染源解析[13]。PMF运算原理是基于最小二乘法进行迭代运算，目标是解决所测量的物质浓度和来源之间的化学质量平衡（CMB）[14-15]。近年来，PMF模型被国内学者应用于土壤重金属污染源解析[16]，与PMF用于大气颗粒物源解析经过广泛论证不同[17]，PMF模型用于土壤重金属源解析尚未进行系统评价。与大气样本往往采自同一观测点的不同时间不同，土壤样本通常采自同一时间的不同点位。田莉[18]在利用PMF模型进行土壤重金属源解析时，提出应考虑样本量个数对PMF模型源解析结果的影响。支裕优[19]在应用PMF模型对浙江长兴煤山盆地土壤重金属进行污染源解析时，首次提出源解析结果受异常值影响。此外成分谱的元素和数量也能影响源解析结果[20-21]。总之，数据处理（成分谱、异常值检验、数据替换、不确定度加权）和样本特征（样本量大小、含量高低、是否存在零含量）均可对模型结果产生扰动。本研究拟考察成分谱和异常值两个因素的影响。

水口山铅锌矿是我国重要的铅锌生产基地之一，大规模矿产开采及冶炼活动已有100多年的历史，导致其周边土壤被严重污染，现已形成以开采、冶炼为主，其他污染为辅的综合污染区[22]。本研究以水口山铅锌矿周边农田土壤为例，探究异常值剔除和地壳元素引入对PMF模型源解析结果的影响，结合元素浓度空间分布图，初步确定重金属的污染来源，为今后应用PMF模型对土壤重金属污染源解析的数据处理和成分选择提供参考。

1 材料与方法

1.1 研究区概况

水口山铅锌矿研究区位于湖南省常宁县水口山区，沿湘江南岸东西向分布，地处东经112.606 421°，北纬26.588 913°。气候属于亚热带季风性湿润气候，四季分明，雨量充沛，春秋天气多变，各年度降水量、气温、日照等气象要素变化无常。以黄壤土为主，主要的作物类型有水稻、蔬菜等。

1.2 样品采集与分析方法

2016年9月，从水口山铅锌矿矿区核心位置出发，由近及远采集土壤，在矿区周边农田分布较多的西南方向（东经 112.576 508°～112.606 421°，北纬26.553 302°～26.597 315°）呈放射状布设采样点，整个研究区采样点尽量覆盖周边的农田区。共布设44个采样点，采集0～20 cm的表层土壤，并使用GPS记录下各采样点的地理位置。采样点分布如图1所示。

图1 研究区农田土壤采样点分布图Figure1 Distribution of sampling sites in study area farmland soil

土壤样品经风干后，四分法取样，过100目筛，放入聚乙烯自封袋中备用。使用电感耦合等离子体质谱仪（ICP-MS）测定样品中 Al、Cd、Cr、Co、Cu、Fe、K、Mg、Mn、Mo、Na、Ni、Pb、Sb、Tl、V、Zn 17种元素含量，Ca含量用原子吸收分光光度计（AAS）测定，As和Hg含量用原子荧光光度计（AFS）测定。每个样品设定3个平行，做试剂空白，采用国家标准土壤样品GSS-14和GSS-16进行质量控制。

表1 土壤组分相关性Table1 Correlation of soil components

1.3 数据分析方法

元素浓度数据分析和空间图采用Excel 2010和Origin 8.0，相关分析及双样本T检验采用SPSS 18.0，PMF模型采用EPA PMF 5.0。

2 结果与讨论

2.1 元素相关分析

相关性分析是土壤重金属污染来源分析的简单方法之一，相关系数越高，说明同源可能性越大[23-27]。相关分析显示（表1），Cu、Zn、Pb、Cd、As、Hg和Sb两两之间，Cu、Cr、Mo、Co、Ni、Tl和V两两之间，Al、Mg、K、Na和Ca两两之间，Mn和Fe之间均呈显著正相关（P＜0.05），Pearson相关系数大于0.30。其中，Zn、Pb和Cd两两之间，Co、Ni和Mo两两之间，Mg和Ca之间呈极显著相关（P＜0.01），且Pearson相关系数在0.80以上，表明元素间同源性很强。

图2 主要元素浓度数据散点图Figure2 Scatter plot of the concentration data of main elements

2.2 不同因素处理

2.2.1 散点图剔除异常值

散点图是一种常用的离群值检验方法，以一种元素（通常选择含量高的）为x值，其余元素依次为y值，绘制散点图。在本研究中，离群值即被认定为疑似异常值。结合所有散点图异常点的样本编号，确定异常样本的数量并予以剔除。本研究土壤采自铅锌矿周边农田，因此以含量高的Zn作横坐标，其余19种元素分别为纵坐标，散点图有3种不同分布，离群值有2个：SK016-T1-002和SK016-T2-014，选择有代表性的元素散点图见图2，方框圈出的点为异常值。图2a、图2c和图2e中除SK016-T1-002和SK016-T2-014号点外，其余样本点均匀分布在一条线周围，说明As、Cd和Hg 3种元素与Zn可能来自于同一污染源；图2b中除SK016-T1-002号点外，其余样本点分别位于两条线周围，说明Ca与Zn可能来自不同的污染源；图2d和图2f中除SK016-T1-002号点外，少部分样本点单独在一起，说明它们可能有其他来源。对离散严重的SK016-T1-002和SK016-T2-014号点进行考察，确认其距离矿区较近，推测这两个土壤样本混入矿石粉成分较高，因此将其剔除。

2.2.2 地壳元素引入

土壤没有被外来污染源污染时，其元素含量由成土母质中元素含量以及成土过程中元素的一系列物理化学变化（风化、淋溶等）共同决定[28]。外来源的输入使土壤中的成分含量发生不同程度的扰动。因此，在成分谱中引入相对稳定且含量高的地壳元素，可以增加对污染源的判断精度，一定程度上能够方便对自然源与人为源的区分。本研究选择Al、Ca、Fe、K、Mg和Na共6种金属类地壳元素。

表2 PMF模型解析出A样本集各元素源贡献率Table2 Source contribution rate of elements of the A sample set by PMF

表3 PMF模型解析出B样本集各元素源贡献率Table3 Source contribution rate of elements of the B sample set by PMF

2.2.3 异常值剔除和地壳元素引入对PMF模型源解析结果的影响

PMF模型根据元素浓度和元素不确定度计算信噪比（S/N）[29]。初步将元素数据载入后信噪比（S/N）都大于2，定义为“strong”。运行模型后，多次调试元素的“strong”、“weak”以及因子个数，以实现元素实测/模拟浓度拟合系数大于0.5，Q值与理论Q值的差值小于10%。

本研究主要考察异常值和地壳元素两个因素，共建立4个样本集：样本集A，未剔除异常值且未引入地壳元素（表2）；样本集B，剔除异常值但并未引入地壳元素（表3）；样本集C，未剔除异常值但引入地壳元素（表4）；样本集D，剔除异常值并引入地壳元素（表5）。

结果显示，异常值剔除和地壳元素引入使源解析结果发生了改变。样本集A解析出4个污染源，源1对Zn、Pb、Cd和Hg有重要贡献，源2对Hg有重要贡献，源3对Cr、Co和V有重要贡献，源4对Mn和Mo有重要贡献。样本集B解析出4个污染源，源1对Zn、Pb、Cd和Hg有重要贡献，源2对Hg有重要贡献，源3对Cr、Co和V有重要贡献，源4对Mn、Ni和Mo有重要贡献。样本集C解析出5个污染源，源1对Zn、Pb、Cd、Sb和Hg有重要贡献，源2对Hg和As有重要贡献，源3对Cr、Co、Mo和Ni有重要贡献，源4对Mn和Fe有重要贡献，源5对Al、Ca和Mg有重要贡献。样本集D解析出5个污染源，源1对Zn、Pb、Cd和Sb有重要贡献，源2对Hg和As有重要贡献，源3对Cr、Co、Mo和Ni有重要贡献，源4对Mn和Fe有重要贡献，源5对Al、Ca和Mg有重要贡献。

A与B相比、C与D相比，异常值剔除前后，污染源成分谱相同，解析出污染源的数量一致，但各污染源的总贡献率发生变化。Norris等[30]运用PMF模型进行大气颗粒物污染源解析时发现，源解析结果容易受异常值的影响，为了使源解析结果合理反映大气颗粒物的普遍污染来源，必须对异常值进行剔除。在运用PMF模型进行土壤重金属源解析时，支裕优[19]首次针对异常值进行研究，他将部分异常高值剔除，结果表明，剔除异常值后，PMF模型能够给出较为合理的源解析结果。分析其原因主要包括两个方面：其一，异常高值中包含比正常值更多的污染源信息；其二，基于PMF模型的求解原理，会优先拟合异常值使目标函数最小，导致污染源的总贡献率趋向异常值。本研究剔除了两个异常样本，其Ca、Fe、Cd和Zn等元素含量均偏高。异常值剔除后，这几种元素所代表的污染源总贡献率变小，进一步验证了异常值对源解析结果的影响，使整体的源解析结果能够更合理地反映当地的污染源情况。

表4 PMF模型解析出C样本集各元素源贡献率Table4 Source contribution rate of elements of the C sample set by PMF

表5 PMF模型解析出D样本集各元素源贡献率Table5 Source contribution rate of elements of the D sample set by PMF

A与C相比、B与D相比，引入地壳元素前后，污染源成分谱不同，解析出的污染源类型和个数均发生变化。土壤中各个元素都是来源于一个或多个污染源，在进行污染源解析时，需选择不易造成污染的元素来反映当地的土壤本底情况。艾建超等[31]利用UNMIX模型对夹皮沟金矿区土壤重金属进行源解析的研究表明，引入Al、Mg和Ca等元素，源解析结果有利于后期污染源类型的辨识。本研究在此基础上，将土壤中含量极高且相对稳定的地壳元素（Al、Ca、Fe、K、Mg、Na）引入PMF模型，由于其几乎不受外来污染的影响，使受土壤本底影响大的元素可以更好地与其他元素进行分离，从而更容易确认各污染源性质。

综上，B与A比较，仅剔除2个异常值，未改变分类结果（污染源仍保持4个），只改变贡献率；C与A比较，仅引入地壳元素，改变分类结果（污染源由4个变成5个）。C与D比较，仅剔除2个异常值，未改变分类结果（污染源仍保持5个），只改变贡献率；B与D比较，仅引入地壳元素，改变分类结果（污染源由4个变成5个）。证明引入地壳元素改变分类结果，异常值剔除改变贡献率。D与A相比，既剔除了异常值，又引入地壳元素，贡献率改变，且分类由4个变成5个，并未进一步增加，表明两者的共同作用只是普通的加和，并不是协同作用。

由样本T检验可知（表6）：A与B解析出的源1显著性（Sig.）小于0.05，说明剔除异常值对源1有显著影响；而源2、源3和源4的Sig.大于0.05，说明剔除异常值后对解析结果的源2、源3和源4无显著影响。C与D相比，5个污染源的Sig.均大于0.05，说明在剔除异常值后，对解析出的5个污染源均无显著影响。因此，引入地壳元素后（C和D），源解析结果受异常值的影响较小，稳定性更强。由于PMF模型对数据个数要求比较高，需要大量的数据来确保化学成分间的不确定度范围[17，32]。因此，引入地壳元素后，异常值占总数据量的比例低于引入前所占比例可能是其结果相对稳定的原因之一。

表6 成对T检验Table6 The factors′paired samples test

2.3 水口山铅锌矿周边农田土壤重金属的污染来源分析

通过对2.2部分的4种不同样本集的分析可知，D样本集的源解析结果最优。结合浓度空间分布图（图3），对D样本集进行污染源分析：因子1对Pb、Zn、Cd和Sb的贡献率高达50%以上，对As和Hg也有较高的贡献率。有研究显示，水口山铅锌矿区的冶炼和开采活动产生的废水、废气、废渣中含有大量的Zn、Pb、Cd、As和Hg等元素，排放后通过各种途径进入土壤，导致土壤污染[33]。Sb与这几种元素呈显著正相关，具有同源性。并且Pb和Hg的浓度越靠近矿山和冶炼厂越高。因此，推断源1为采选矿及冶炼等工业活动源；因子2对Hg和As的贡献率较高。浓度空间分布图显示，Hg除矿山附近有高值外，在靠近湘江支流的位置也有高值区。根据前期的实际调查，研究区主要种植水稻和蔬菜，用于灌溉的水大多来自湘江支流，在种植过程中喷洒以杀虫剂为主的农药以及施用各种类型的肥料。农业投入品（化肥、农药等）中含有的Hg、As、Cu、Zn、K等元素会在土壤中残留，仅就磷肥而言，As的含量一般在20～50 mg·kg-1，其不合理施用会使土壤中As含量升高[14，34-35]。因此，推断源2为污水灌溉和农药化肥施用等农业活动源；因子3对Mo、Co和Ni的贡献率均超过50%。Ni的浓度空间分布图显示，除在公路附近有高值区外整体平缓，且从表7中可以看出，这3种元素的变异系数相比其他元素均偏小，说明3种元素的空间变异程度小，表明其受人为活动影响不明显。同时3种元素的平均值略高于湖南省背景值。因此，推断源3为土壤母质源；因子4对Fe和Mn的贡献率比较高，分别为63.57%和51.21%。从浓度空间分布图可以看出，Fe有两个高值区，主要是因为水口山铅锌矿选矿副产的硫铁矿精矿，有一部分会分别销售于松柏化肥厂和水口山一厂，在运输中铁矿石颠落会导致Fe和伴生组分Mn在道路两旁土壤中含量增高[31，36]。同时因子4对Pb、Ni和Mo也有一定的贡献。在交通运输方面，Pb被作为汽车尾气的指示性元素[18]，机动车所用的润滑油在高温下与空气中的氧气发生氧化反应生成有机酸、醛、酮、醇及其他有机化合物，这些化合物能够腐蚀与其接触的含Ni和Mo油泵[37]，导致其向环境中释放金属元素。并且研究区内交通网较发达，有两条高速公路贯穿。因此，推断源4为铁矿石开采和交通运输源；因子5对Ca和Al的贡献明显，Ca的浓度分布和Pb的浓度分布比较相似，在矿山附近偏高，但其浓度跨度比Pb的跨度要小。矿山周边岩石经过风化与生物作用，发生一系列物理作用变为细小的颗粒，经过雨水冲刷和搬运后导致各元素在土壤表层迁移并积累[38-39]。表7中，Ca的含量高于湖南省背景值，Al的含量接近湖南省背景值。路睿[40]对水口山铅锌矿的岩石进行研究发现，研究区主要为花岗岩，Al在蚀变过程中含量相对稳定，Ca由于后期大量碳酸盐化蚀变而含量升高。因此，推断源5为岩石风化源。

图3 土壤主要污染源代表组分空间分布Figure3 Spatial distribution of representative components of main pollutant sources of soils

表7 研究区农田土壤金属元素分析统计表Table7 Statistical analysis of elements in the study area farmland soils

3 结论

（1）异常值剔除前后，PMF模型解析出的源数量和源成分谱相同，但各污染源对元素的贡献率和总贡献率均发生变化。

（2）地壳元素的引入对PMF模型解析的污染源类型和个数产生影响，使各污染源的性质更容易确认；并且引入地壳元素后，源解析结果受异常值的影响较小，稳定性更强。

（3）选择最优的D样本集源解析结果，结合浓度空间分布图得出5个污染来源：采选矿及冶炼等工业活动源（26.81%）、污水灌溉和农药化肥施用等农业活动源（14.68%）、土壤母质源（24.41%）、铁矿石开采和交通运输源（16.39%）、矿石风化源（17.72%）。