基于小流域尺度的矿区农田土壤重金属源解析

杨宇，周佳俊，郭婷婷，丁洪伟，刘孝利*

（1.湖南农业大学资源环境学院，长沙 410128；2.湖南省灌溉水源水质污染净化工程技术研究中心，长沙 410128）

当前，我国耕地土壤重金属污染情况不容乐观，全国首次土壤污染状况调查显示，耕地土壤点位超标率为19.4%[1]。当种植在受污染土壤中的作物对重金属的吸收和富集达到一定限值时，则会危害人体健康[2]。近年来，国内外学者就耕地土壤重金属污染问题展开了包括但不限于土壤-作物系统累积特征[3-4]、分布特征[5-7]、生态与健康风险评价[8-9]、来源分析[10-11]等方面的研究。其中，明确土壤重金属污染来源是从源头上防治土壤重金属污染的根本措施，对于提高土壤环境质量具有重要意义。

目前国内外在土壤领域应用较多的污染源解析方法主要有同位素示踪法、多元统计法、UMIX 模型法、化学质量平衡法（CMB）、正定矩阵因子分析法（PMF）和绝对因子得分-多元线性回归法（APCSMLR）等，不同源解析方法根据模型原理或假设不同，所需样本量和数据处理量也有所差异，最终解析出的结果可能不尽相同[12]，因此常需要根据研究区域的实际情况联合应用多种方法，取长补短。如马杰等[13]将APCS-MLR 和PMF 模型相结合使用，发现两种方法解析出的结果基本一致。郑宇娜等[14]结合地统计、多元统计和PMF 模型3 种方法，彼此佐证，准确识别了台州市电子垃圾拆解场地周边农田土壤中重金属污染来源。Fei 等[15]将PMF 模型与地理探测器模型（GeogDetector）联合应用，结果证实了这种新的合成模型可以很好地用于估计土壤重金属的来源分布。Lv[16]使用地统计学描述比较了APCS-MLR 和PMF 两种方法的解析结果。

浏阳七宝山矿区矿产资源丰富，采矿历史悠久，长期的采矿活动以及数百家大小矿山在历史上采取的粗放型管理方式导致了严重的土壤重金属污染问题，极大影响了农产品的安全生产，针对这样的典型涉矿区域重金属污染进行研究极有必要。目前关于土壤重金属污染研究的尺度多集中在较大的省、市等，或田块、采矿场地及周边等较小尺度，以农业小流域特别是典型矿区所在流域作为尺度的研究相对较少。流域作为陆地生态系统集水与水循环的基本单元，是污染综合治理的最佳研究与管理尺度。从小流域尺度开展耕地土壤重金属污染来源及其空间异质性特征研究，对了解耕地土壤重金属污染来源、扩散和空间迁移特征，以及在流域尺度实现耕地的分级、分区污染综合防治具有重要的科学意义和实践价值。因此，本研究选取浏阳七宝山矿区小流域为研究对象，利用ArcGIS 进行网格化布点，采集水田耕层土壤（0～20 cm）样品，运用地统计学、多元统计和PMF 模型相结合的方法，解析小流域内7 种重金属（Cd、As、Fe、Mn、Cu、Zn 和Pb）的污染特征，有效识别其污染来源，以期为涉矿区域农业小流域耕地重金属污染综合防控与耕地安全利用提供科学参考。

1 材料与方法

1.1 研究区概况

研究选取的小流域位于湖南省浏阳市，介于北纬28°14'～28°20'、东经113°50'～113°53'之间，总面积约为72 km2。流域总体为一片向西开阔的狭长谷地，南北两侧为林地，地貌以剥蚀丘陵和溶蚀堆积谷地为主[17]，成土母质主要有砂页岩、石灰岩风化物和中性冲积物。流域内主要河流为浏阳河的二级支流，其也是该地农业生产的主要灌溉水源。区内矿冶开采活动历史悠久，目前仍有3 处正在营业的厂矿企业，分别为铜锌矿厂、硫铁矿厂和磷肥厂（图1）。

图1 研究流域采样点位、主要道路、工矿及水系分布情况Figure 1 Sampling points，main roads，and distribution of industrial mines and water systems of study area

1.2 样品采集与分析

首先利用ArcGIS 10.6 软件进行网格化布点，布点原则遵循《土壤质量 土壤采样技术指南》（GB/T 36197—2018），布点精度为800 m×800 m，根据前期收集到的有关资料，结合该小流域的地形地貌和工业分布情况，在城镇及矿冶活动区进行了点位加密处理[18]。实地调查采样时使用“两步路”软件按照预设样点的经纬度进行定位搜索并记录实际采样点的经纬度，最终共采集60个土壤样品（图1）。采样工具选用直径为80 mm 的手动螺纹钻，垂直采集0～20 cm 的表层土壤，每个样点选取同块田地内3～5个彼此接近的较小点的样本进行充分混合制成混合样，共计土壤质量约1.5 kg。采样选在4 月进行（早稻种下前），采样同时现场记录样点的土地利用类型、历史耕作情况、灌溉水来源等有关信息，采集的样品装入聚乙烯密封袋中带回实验室保存。

样品在室内避光处自然风干，挑出石砾、植物根系等杂物，研磨后依次过10目和100目分样筛。首先用HCl-HNO3-HClO4三酸消解法消煮，经0.45 µm 滤膜过滤后使用电感耦合等离子发射光谱仪（ICP-OES Optima 8300）测定Cd、As、Fe、Mn、Cu、Zn 和Pb 的全量，每批次样品消解全程均设置空白对照和3 组平行样，并使用国家标准土样（GSS-3）进行质量控制。土壤pH 值按《土壤pH 值的测定电位法》（HJ 962—2018）使用pH计（雷磁PHS-3C）测定。

1.3 地统计学和半方差函数

地统计学是一门以区域化变量理论为基础，将克里金法作为基本分析方法，研究在空间分布上具有变异性或结构性的现象及其规律的学科[19]。半方差函数又称半变异函数，可以描述区域化变量的结构性和随机性变化[20]。该函数常被用来进行各种空间变异分析，为空间预测提供基础，其计算公式如下：

式中：N（h）代表空间上步长为h时离散点的数量；Z（xi）和Z（xi+h）分别代表区域化变量在位置xi和xi+h时的值。

半变异函数在计算过程中会产生4 个重要的参数：块金值（C0）、基台值（C）、变程（A）、结构方差（C1）。块金值用来衡量在实际工作中客观存在的测量误差和空间变异导致的块金效应，基台值（C=C0+C1）反映了某区域化变量在研究范围内变异的强度，其中C1代表由于样本数据中存在空间相关性而引起的方差变化范围。块金值与基台值的比值可表示空间变异性程度，比值高说明由随机部分引起的空间变异性程度较大，反之说明由空间自相关部分引起的空间变异性程度大[21]。

1.4 PMF模型

PMF 模型最早由Paatero 等[22]在1994 年提出用于大气污染源解析，是一种以受体模型为基础并根据源成分来量化每种源对样品数据贡献量的数学方法，近年来也被广泛应用于土壤重金属污染源解析。PMF模型将实测数据矩阵（Xij）分解为因子贡献矩阵（gij）和因子成分矩阵（fij），目标是通过最小二乘法经过多次迭代运算得到目标函数Q，函数Q的定义如下：

式中：uij表示第i个样品的第j个元素的不确定值。

uij的计算方法如下：

式中：RSD为相对标准偏差；Cij为第i个样品的第j个元素的含量；MDL为方法检出限。

1.5 数据处理与分析

数据使用Excel 和SPSS 23 进行预处理，使用Origin 2018 进行相关性和主成分分析，地统计分析、半方差函数模型拟合利用GS+ 9.0 软件进行，使用Arc-GIS 10.2 进行克里金插值分析，绘制空间分布图，使用EPA PMF 5.0软件进行定量源解析。

2 结果与讨论

2.1 重金属含量的描述性统计

研究流域内土壤pH 为4.32～6.41，整体呈酸性，与全国表层土壤pH 中位值（7.88）相比，有明显的酸化特征。表1 为重金属含量描述统计分析结果，7 种重金属含量数据均呈现高峰和正偏形态，说明含量高的点偏多且数据分布比较集中，具有整体性。Cd、As、Cu、Zn和Pb的均值分别为湖南省土壤背景值[23]的10.0、2.1、4.5、2.0 倍和4.8 倍，且除As 外都超过了《土壤环境质量农用地土壤污染风险管控标准（试行）》（GB 15618—2018）规定的风险筛选值，说明5 种元素均有不同程度的富集，且以Cd、Cu 和Pb 较为严重。就变异系数而言，Fe、Mn和Zn属于中等变异，Cd、As、Cu 和Pb 达到强变异程度，说明研究流域内这4 种重金属可能受人类活动影响较为显著。

表1 土壤重金属含量统计结果Table 1 Statistical results of heavy metals in soil

2.2 半方差分析

半方差函数拟合要求数据服从正态分布或近似于正态分布[24]，因此先对原始数据进行剔除异常值、对数转换等处理，经Kolmogorov-Smirnov 检验P＞0.5，说明处理后的数据符合正态分布，在此基础上再进行拟合。计算7 种元素的函数值，根据残差平方和（RSS）最小和决定系数（R2）最大的原则获取最佳拟合模型的参数[25-26]（表2）。结果显示7 种元素的决定系数均在0.7 以上，较好地反映了各自的空间结构特征，符合最佳拟合模型，其中Cd、Mn、Cu 和Pb 符合指数模型，As 和Fe 符合球形模型，Zn 符合高斯模型。变程描述该间隔内样点的空间相关特征，若两点间的距离超过了变程，那么其不具有空间相关性。7 种元素拟合得到的变程均大于采样间距（800 m），说明采样数据间存在空间相关性，可以用于内插或外推。根据张红桔等[27]的研究，结构性因素（如成土母质、母岩、地形地貌和土壤矿物等）和随机因素（如工业生产、农业活动和土地利用强度等人类活动）是影响土壤性质发生空间变异的两大决定性因素。而块基比能够反映由随机因素引起的空间变异占总体变异的比例，块基比＜25%表示变量具有显著的空间自相关性，其空间变异主要受到结构性因素的影响；块基比＞75%表示变量间的自相关性不显著，受随机因素的影响大于结构性因素；块基比介于25%～75%之间时，变量具有中等自相关性，其空间变异特征由两种因素共同决定[24]。7 种元素的块基比结果表明（表2），Fe、Mn 和Zn 主要受结构性因素影响，其空间变异主要是由自然活动引起的，Cd、As、Cu 和Pb 受两种因素共同影响，其中Cd、Cu 和Pb 的块基比约为50%，且样点间表现出很强的随机性，说明其受到较高程度的人为活动的影响。

表2 半方差函数拟合结果Table 2 Fitting results of semi-variogram

2.3 耕地土壤重金属来源分析

2.3.1 多元统计分析

研究显示相关系数越大，元素间存在同源、伴生或复合污染的可能性越大[28-29]，聚类分析则可以表明元素间的分类情况。7 种元素的相关性分析结果显示（图2），Mn与Fe极显著相关（P＜0.01），与As显著相关（P＜0.05），与其他元素不相关；Fe 与除Cd 外的5 种元素都极显著相关；Cd、As、Cu、Zn 和Pb 两两之间极显著相关。而聚类分析将7 种元素分成了三大类（图2）：第一类为Fe 和Mn，第二类为As 和Pb，第三类为Cd、Cu 和Zn。由结果可初步推断，Fe 和Mn 可能具有同源性，其中Fe 的来源更为复杂；As 和Pb 的主要来源可能相同，Cd、Cu、Zn 和Pb 可能存在复合污染，其中Cd、Cu和Zn的来源可能相似。

图2 相关性及聚类分析结果Figure 2 Results of correlation and cluster analysis

为进一步确定研究流域耕地表层土壤7 种重金属的来源，利用SPSS进行了主成分分析（表3）。经过正交变换后，提取出的两个主成分累积解释了75.94%的总方差，能够说明原始数据的大部分信息，符合预期要求。Cd、As、Cu、Zn 和Pb 在第一个主成分（PC1）上的载荷，与Fe 和Mn 在第二个主成分（PC2）上的载荷均大于0.8，说明PC1 和PC2 分别能够反映上述两组元素的绝大部分信息。该结果与相关性和聚类分析的结果一致，进一步说明了Fe 和Mn 的主要来源明显区别于其余5种元素。

表3 提取成分矩阵Table 3 Extract component matrix

2.3.2 PMF模型解析

将几种多元统计方法结合在一起能够相互对比验证，定性分析出不同重金属的主要污染源，但不能定量解析出每类源的贡献率和各个元素对源的贡献度。因此，在多元统计分析结果的基础上，使用PMF模型设置因子数为3～5个，调整不确定性参数进行迭代运算20次，发现当因子数设置为3时，包含所有点的拟合参数（Qture）和剔除高残差样点后的拟合参数（Qrobust）非常接近，除Fe外的6种元素拟合度R2均大于0.8（Fe元素可能来源较为复杂），且样本残差大部分处于-3～3之间。同时对3个因子在各样点上的贡献率进行插值，以便揭示每种污染源在空间上的分布特征。

结果显示，As、Pb 和Cu 在源1 上有很高的贡献率，分别达到了78.2%、64.7%和59.1%，其中As和Pb的贡献率远高于其他重金属（图3）。由图4可知，源1自流域中游起向下至出口处沿干流（主干路）左右两侧有较高的贡献率。有研究表明，As 和Cu 是农药的主要组成成分[30-31]，农业投入品的不合理使用会造成土壤中As、Cu 和Zn 等的积累[32]。根据实地调查资料，研究流域内有一磷肥厂（图1），Zhang 等[33]研究发现，施用磷肥会促进As 的协同作用，致使As 以砷酸盐和亚砷酸盐的形式在土壤中富集。Pb 经常被认为是机动车燃油和汽车尾气的标志元素[34-36]，张恬雨等[37]研究发现Cu 和Zn 也常与交通运输有关，汽车轮胎和路面磨损、运输物泄露和尾气排放等产生的含Cu 和Zn 的粉尘会在沉降作用下在土壤中累积[38-40]。除此之外，Cd 和Zn 在源1 上也分别有18.2%和28.9%的贡献率，Fe和Mn则较少。据调查，该流域内的干流沿一条贯穿流域的主要道路由东向西而下，为其两侧耕地主要的灌溉水源（图1）。戴塔根等[41]早期对该流域进行调查时发现，河流两侧农田重金属的积累与长期使用受污染的河水灌溉有很大关系。杨宇等[6]的研究结果表明，该流域内距离河流近的比距离较远不受河流灌溉的农田的土壤生态风险等级高。到目前为止，该地已开展过多次矿区地质环境改善防治工作[42]，河道水质得到提高，污水灌溉已成为流域土壤重金属污染的历史源。因此，可以认为源1主要为施用农药化肥、交通运输和污水灌溉（历史遗留）的混合源。

图3 土壤重金属PMF因子贡献率Figure 3 Factor contribution from PMF of soil heavy metals

图4 污染源空间分布特征Figure 4 Spatial distribution of pollution sources

源2 的主要贡献元素为Cd（80.7%），其次为Zn（52.8%），Fe、Cu、Pb、Mn 和As 分别有32.2%、28.2%、22.8%、14.2%和5.0%的贡献率。根据每种元素在源2 上的贡献率，初步推断Cd 为源2 的标识性元素。从空间分布特征来看，源2 在两处采矿区尤其是金属矿厂附近有较高的贡献率（图4）。有研究表明，在无工业活动干扰时，Cd 常被看作农耕活动的标志性元素[43-44]，但在人为活动干扰下，重金属极易在土壤富集，造成空间分布上的不均匀性[45]。研究发现，Cd 常在矿石冶炼、精炼的过程中伴生出现[46-49]，而流域内金、银、铜、铅、锌和铁锰矿产资源丰富，矿冶开采活动历史悠久。采矿常伴随着“三废”的产生，在绿色矿山开发模式出现前，废气可以通过干湿沉降向四周扩散，废水未经处理直接排入铁山河，发生降雨时堆积在地表的废渣则可能随降雨径流或流失的水土向周围的环境迁移，由此造成“三废”中的重金属在空间上产生富集。因此判断源2主要为工业源。

源3 的主要贡献元素为Mn（85.5%）和Fe（65.6%），主成分分析中将Fe 和Mn 提取为PC2（表3），聚类分析将Fe和Mn单独分为一类（图2），说明二者的主要来源区别于源1和源2。图4显示，源3主要在流域上游贡献率较高，其高贡献率的位置特征在空间上与河流和矿厂的分布不相关。Mn作为地壳的组成元素之一，在自然土壤中有较高的丰度。Zhou等[50]发现在以石灰岩为基础风化发育的土壤中Cu 和Mn的含量较高，而研究流域的主要成土母质之一就是石灰岩。同时，丰富的铁矿资源也会使含Fe 的硫化物和氧化物在岩石风化和雨水冲蚀的作用下迁移到土壤中。因此可以推断源3主要为成土母质。

3 结论

（1）研究流域内耕地土壤有明显的酸化特征，Cd、As、Cu、Zn 和Pb 的均值高于湖南省土壤背景值，受人为活动影响显著，Cd、Cu、Zn 和Pb 的均值超过农用地土壤污染风险筛选值（GB 15618—2018）。

（2）本研究的采样间距合适，半方差函数拟合结果较好，所采样点能够用于内插或外推进而反映研究流域耕地土壤重金属污染的整体特征与空间差异。

（3）表层土壤中As 和Pb 的主要来源为农药化肥等农业活动和交通运输；Cd 主要受工业活动如矿石冶炼的影响；Fe 和Mn 主要受成土母质风化影响；7 种重金属在耕地表层土壤中的积累均在不同程度上受到历史污水灌溉、干湿沉降的影响。