不同三维模型对铅酸蓄电池污染场地土壤Pb空间分布预测的影响

刘 庚,王尚义,牛俊杰,张 朝,赵 鑫,郭观林*(.太原师范学院,汾河流域科学发展研究中心,山西太原 0300；.中国环境科学研究院,环境基准与风险评估国家重点实验室,北京 000)

不同三维模型对铅酸蓄电池污染场地土壤Pb空间分布预测的影响

刘 庚1,王尚义1,牛俊杰1,张 朝2,赵 鑫2,郭观林2*(1.太原师范学院,汾河流域科学发展研究中心,山西太原 030012；2.中国环境科学研究院,环境基准与风险评估国家重点实验室,北京 100012)

为准确界定污染场地土壤中污染物三维分布范围和受污染土壤土方量并提供相关技术方法和思路,以某典型铅酸蓄电池污染场地为例,对比研究了不同三维插值模型(三维克里格、反距离加权、最邻近点)在不同水平垂直向异性系数设置,对土壤铅污染评价和污染边界分布的影响.研究结果显示,不同模型预测的精度不同,随着水平垂直比值增大,平均误差和均方根误差有增大趋势;通过对比修复目标表明,反距离加权统计的受污染土方量最大,其次是三维克里格和最邻近点,统计值相差约为 13%;从污染评价结果和预测精度来看,三维克里格模型和反距离加权要优于最邻近点,选择适宜的模型和合理的参数设置对准确界定污染边界和降低不确定性具有重要影响.

污染场地；三维插值；污染评价；不确定性

基于有限的土壤钻孔样点数据,采用不同的空间预测模型插值计算是描述污染场地污染物空间分布制图的主要手段[1-2],是确定修复边界以及修复治理决策制定的重要基础[3-4],同时在认识污染物环境行为过程上也具有重要作用和意义.目前基于机理模型和空间统计模型的二维空间插值技术应用相对已比较成熟[5-6],但其预测结果割裂了变量在三维空间上的连续性分布特征,因此,很多学者在二维插值模型基础上,探索适合重建三维空间场的预测方法.常用的三维空间插值模型有三维克里格模型、反距离加权模型、邻近点取样模型、高次函数逼近模型以及分形插值模型等[7-9].这些三维插值模型已应用到地层建模、土壤盐分、渗透系数、土壤质地等不同领域[10-13].据现有研究报道,三维随机模拟方法和三维反距离加权模型结合土壤电导率能够真实反映土壤盐分的空间分布以及对土壤盐分进行精确解译和评估[14-15];对土壤质地层次三维空间分布模拟时,多点地统计学三维模拟方法在精度上要优于传统顺序指示模拟[16];在对土壤中元素和土壤污染物三维分布模拟时,三维克里格模型预测精度更高,能够更好表征局部区域变量的分布规律[17-18].不同模型的适用范围和精度优势并不相同.受模型参数以及对预测数据要求等因素影响,其中三维克里格、反距离加权和邻近点是这些三维插值模型中使用频率较高、应用较为广泛的几种方法.

在污染场地土壤污染评价以及污染物三维空间分布表征时,重点关注的是污染物空间分布范围和受污染土方量的统计.污染边界的准确界定能够合理节约场地修复治理成本,提高场地环境管理效率.受污染物累积释放和迁移转化因素的影响,污染物含量在不同土层中具有各向异性特征,缺乏对变量的空间各向异性空间特征分析,会影响三维插值计算的精度,插值结果不能真实反映变量的三维空间分布特征.三维预测模型虽已应用于不同变量的插值计算,但在污染场地土壤中污染物的三维分布表征方面的应用较少,同时现有三维模型在空间分布预测过程中,也缺少对变量空间向异性特征的分析.因此,本文选择某典型铅酸蓄电池污染场地土壤Pb污染为研究案例,在分析污染物各向异性空间结构特征基础上,对比研究常用的三维克里格、反距离加权和邻近点模型在不同各向异性特征参数设置下的预测精度,分析对受污染土方量计算和污染评价制图的不确定性,为指导该场地的修复治理提供依据.

1 材料与方法

1.1 目标污染场地概况

目标污染场地曾是我国西南地区较大的铅酸蓄电池生产厂家,始建于 1943年,场地面积为13.6hm2,主要生产军用通信干电池、铅酸蓄电池和制铅粉等产品,生产主要原材料包括铅、锌、硫酸等.生产过程中涉及到重金属污染.根据工程地质测绘及钻探揭露,该场地地层有素填土(Q4ml)、粉质粘土(Q4el+dl)、砂岩和泥岩.该场地所在区域属亚热带湿暖气候,无霜期长,日照时间短,空气湿度大,年均气温约为 18.5℃,最高气温近 40℃,年均降水量为 1187mm.通过对场地的生产历史、工艺特点和布局、现阶段场地现状等资料进行收集和分析,对该场地展开了初步污染调查并设置了土壤监测点,明确了该场地存在明显的铅污染问题.

1.2 样品采集与分析

在现场采样时,该厂区部分厂房已经拆除,表层土中含有建筑垃圾,场地实际情况较为复杂,部分区域难以采集到有效样品,根据前期初步环境调查结果,结合该厂区实际现状,本次采样采用随机布点和判断布点相结合的原则,重点关注已经暴露的土壤污染地点和污染物特征并确保取样点对整个场地有合理的覆盖,对前期已经察明的重污染区域加密布点,对已知未污染或污染较轻区域样点间距适当加大.结合场地水文地质背景资料以及实际现状,将该场地从表层到下依次划分为建筑垃圾层(0～0.5m)、杂填土层(0.5～1m)和原土层(1m以下),土壤样品按照上述分层进行3层采样,采样深度为0～9.4m,共布设有效土壤钻孔样点79个,每个钻孔基本按照分层的深度取样,每个孔位取样品 3～5个,经分析化验后,取本层中污染物含量的最大值来进行插值计算.所采集样品的数量和代表性能够满足本文所使用三维插值模型的要求.场地的三维地形和钻孔样点分布如图 1所示.因场地搬迁与拆除造成部分建筑垃圾堆砌在地面,钻孔采样采用挖掘机和钻机相结合的方式,在预先设定的采样点位钻探,获取该点位的土柱样品,将采样管移出地面后,将筛选出的土样混合均匀后立即装入250mL直口玻璃瓶中,并正确密封、填写标签、放入装有适量低温蓝冰的保存箱中,样品采集完毕后送回实验室并进行检测.

土壤中Pb含量采用USEPA 6010C-2007方法进行测定[19],具体分析检测流程:将采集到的土壤样品经自然风干,粗磨除去土壤中的碎石等异物,过 20目尼龙筛,混匀后用四分法缩分至约100g,再用玛瑙研磨,过100目尼龙筛,混匀后备用,准确称取0.100g过100目筛的样品(100℃烘干4h),放入聚四氟乙烯高压釜内,加入 3mL优级纯硝酸、1mL次氯酸、1mL氢氟酸,放入不锈钢套内,拧紧盖,烘箱内160℃加热4h,冷却后取出样品,定容到10mL在日本日立Z-2000原子吸收光谱仪上测定土壤中Pb含量.

图1 场地土壤钻孔样点分布Fig.1 The distribution of drilling soil samples in sites

1.3 插值方法及参数设定

选择三维插值模型中应用较为普遍的三维克里格模型(3D-OK)、反距离加权模型(IDW)和邻近点模型(NN)进行该场地土壤铅的三维空间插值预测,每种模型的算法和原理见参考文献[20-22].污染物含量在不同土层中具有各向异性空间结构特征,在插值计算时,邻域点集内水平方向的点与垂直方向的点对预测点含量值的影响权重不同,具体表现在水平垂直各向异性比系数上,水平垂直各向异性比系数一般默认值为10,即垂向上 1个长度单位的采样点影响权重,与水平方向上10个长度单位的采样点权重相同.为揭示不同水平垂直比值对不同模型预测精度和污染评价的影响,不同模型设定的水平垂直比值、搜索半径等参数情况如表 1所示.在应用三维克里格模型前,需要对原始数据集进行对数正态变换处理,根据不同水平垂直比值设置后拟合的最优三维半变异函数理论模型进行插值计算,最后将预测结果根据正态变换公式进行逆变换回推到原始单位.

表1 三维插值模型及参数设置Table 1 Three dimensional interpolation models and parameters

1.4 数据处理方法

本研究中Pb在不同土层含量的描述性统计分析采用软件SPSS16.0,三维半变异函数拟合采用软件S-GeMS2.1,三维空间插值计算及可视化表征采用软件 EVS. Pro9.42.不同插值方法预测结果精度验证采用交叉验证法中常用的评价指标平均误差(ME)和均方根误差(RMSE)进行评价[23],ME越接近于0,RMSE越小表明插值精度就越高.

式中:u( χi)为预测值,u*(χi)原始样点值.

2 结果与讨论

2.1 场地不同地层土壤中 Pb含量的描述性统计分析

从图2可以看出,在不同土层中钻孔样点含量数据的极差较大,从第1层到第3层的最小值分别为 23.6,21.8,19.2mg/kg,最大值分别为166000,35200,22900mg/kg,最小值与最大值相差数千倍,表明在场地局部区域污染严重,样点含量值极高.从数据频率分布直方图发现,不同层样点含量数据均存在右偏尾现象,不符合正态分布特征,且不同层样点含量数据的偏度和峰度值都较大,表明样点含量数据集在空间上具有一定的离散特性.不同层样点的变异系数均超过了150%,具有很强的空间变异特征.对同一钻孔不同层的含量值进行比较分析,从表层到底层的含量值呈减小趋势,表明场地土壤中污染物的累积释放受迁移转化和人为干扰因素的影响.不同土层的样点含量数据集中均含有异常真实高值点,远超过北京住宅用地标准中规定的 400mg/kg,目标场地已经存在较大的环境风险.将钻孔样点叠加到厂区平面布置图上可以看出,含量超标严重的样点主要分布在废铅堆放区,五车间和配件厂、二车间、一车间和四车间等区域内,赋存在土壤中的污染物与污染源分布以及不同生产工艺产生污染物的情况较为一致,上述车间在生产过程中对局部区域土壤产生较为严重的Pb污染.

图2 场地土壤中Pb污染数据的基本统计特征Fig.2 Statistical characteristics of Pb concentration in soils (from sites)

2.2 不同模型对场地土壤中 Pb含量的三维空间分布预测精度比较

不同的插值模型可以应用于该场地土壤 Pb的三维污染分布表征,但不同模型的精度优势和插值效果不同.本文选用交叉验证法中常用的均方根误差(RMSE)和平均误差(ME)2个评价指标来分析不同模型的预测精度,RMSE值越小、ME值越接近 0,表明该模型插值结果精度越高.交叉验证结果如表 2所示,不同模型随着水平垂直向异性比值增大,其 RMSE和 ME都呈增大趋势. 3D-OK5的 ME和 RMSE最小,分别为 42.7, 221.5mg/kg,NN20的ME和RMSE最大,分别为89.7,892.7mg/kg.总体上来看,3D-OK模型在不同参数设置下的预测精度优势较为明显,其次是IDW模型,NN模型预测精度最低.3D-OK在预测过程中考虑了污染物含量的空间结构信息,虽具有一定的平滑效应,但基本上还能较好地反映出污染物三维空间分布信息,IDW模型和NN模型在插值过程中只依据数据的几何结构特征进行计算,不能反映出变量在空间上的变异特征,因此预测精度相对较低.

表2 不同插值模型的预测误差Table 2 ME and RMSE for different interpolation models

2.3 不同模型估算场地不同地层受污染土方量分析

参照北京市住宅标准所规定的 Pb含量值400mg/kg,将高于 400mg/kg的区域界定为污染范围,对不同模型的三维分布预测结果受污染土方量进行统计,不同模型对不同土层以及总的受污染土方量如表3所示.从统计结果可以看出,不同模型和同一模型不同参数计算出的受污染土方量都不一致.受表层真实高值样点的影响,3种模型随着水平垂直比值的增大,第 1层计算出受污染土方量基本呈减小趋势,第2层和第3层受污染土方量呈增大趋势.3D-OK模型和NN模型统计的总污染土方量随水平垂直比值增大而增大,IDW模型计算出总的受污染土方量随水平垂直比值的增大而减小.所有模型中IDW5方法统计出的受污染土方量最大,为 263402m3,NN5方法统计的受污染土方量最小,为 215379m3.总体来看,3D-OK模型和IDW模型在不同参数下计算出的受污染土方量明显大于 NN模型.从分层统计结果来看,表层受污染土方量最大,其次是第2层,第3层受污染土方量最小,表明该场地表层土壤污染最为严重.3D-OK模型和IDW模型在不同参数下计算出的表层受污染土方量较为接近,NN模型在不同参数下对表层的计算结果低于另外 2种模型.从分层和总的受污染土方量统计结果来看,第1层中受污染土方量总体呈现为NN＜IDW＜3D-OK,第2层为NN＜3D-OK＜IDW,第3层为3D-OK＜NN＜IDW,不同模型对受污染总土方量统计结果为 NN＜3D-OK＜IDW,同一模型在不同参数设置下的统计结果也略有差异.从上述分析结果可以看出,不同模型受其适用原理以及样点含量数据集特征的影响,预测结果中统计的受污染土方量并不一致.在进行钻孔采样时,该场地已进行了部分拆除,建筑垃圾和表层杂填土已经混为一体,因此,本文基于所有钻孔样点含量数据采用不同模型进行受污染土方量统计计算时,包含表层杂填土中的建筑垃圾部分,在对建筑垃圾中铅浸出含量分析时发现未超标,在估算修复成本时,应根据建筑垃圾所占比例将其方量剔除.

2.4 不同模型界定的 Pb含量三维污染范围评价比较

在对不同土层钻孔样点含量数据规范化以及预处理后,采用表1中的不同模型和设定的模型参数对数据集三维插值计算并进行污染制图和可视化表征.使用 3D-OK模型前需要三维半变异函数的拟合,根据拟合后的最优理论模型进行三维空间分布预测.参照北京市住宅标准所规定的Pb含量值400mg/kg,将不同模型的三维分布预测结果转化为栅格数据处理,将高于400mg/kg的区域界定为污染范围,最终形成的污染超标范围如图3所示.

表3 不同模型统计的受污染土方量(m3)Table 3 The contaminated earth volume analyze by different interpolation models (m3)

从图 3可以看出,不同三维模型以及同一模型不同参数预测的污染范围总体趋势相似,但局部区域有些差异.从表层预测的污染范围来看,3D-OK模型预测局部高污染区域范围最大,其次是NN模型,IDW模型预测的局部高污染区域范围较小.对预测结果进行不同方向上的切片处理可以看出,不同模型预测的表层土壤污染最为严重,其次是第2和第3土层,这与钻孔样点含量值的分布特征较为一致.IDW 模型随着水平垂直向异性比值增大,局部污染超标严重区域有增大趋势,表明污染空间分布预测结果受不同模型和模型参数的影响较大.污染超标严重区域是本场地风险评估和修复治理需要重点关注的位置,不同模型预测的污染超标范围分布较为相似,但污染超标严重区域主要分布在场地的中下部位置,结合该场地原厂区平面车间分布图以及生产工艺和历史生产活动可知,在厂区的中下部主要有铅堆放区,五车间和配件厂、二车间、一车间和四车间等生产或存储车间,在生产过程中涉及到 Pb的存储和使用,是产生土壤Pb污染的最主要原因,个别车间在生产、存储过程中的泄露、遗洒等原因使得局部地区污染超标严重.场地土壤Pb污染程度从表层到底层基本为减小趋势,表明存留在土壤中的污染物分布规律受污染物在不同土层的迁移转化特征和污染源分布等因素影响.

图3 场地土壤Pb含量空间分布Fig.3 The spatial distribution of Pb in soil

3 结论

3.1 受不同插值模型适用原理和原始数据集数据特征影响,该场地土壤Pb含量的三维插值误差相对较大,不同模型预测的ME均大于40mg/kg,插值的RMSE均大于200mg/kg.不同模型的插值误差随水平垂直比值的增大呈增大趋势.水平垂直向异性比值越小,预测精度越高,但不能真实反映污染物含量的空间向异性特征,从预测结果的切片和揭层显示来看,比值设为 15时,预测结果比较符合场地的实际污染状况.

3.2 从不同模型三维空间分布预测精度和污染评价结果看,3D-OK模型和IDW模型要优于NN模型,后者界定的污染范围与实际污染情况差异较大,不适合该场地土壤Pb的三维污染评价.

3.3 局部污染严重区域的样点含量具有很强的空间离散特征,对样点高值区域和高值向低值过度区域加密采样,则能够更好地反映污染物空间向异性特征,以提高污染物空间分布预测精度和更准确界定污染分布范围.

[1] Bishop T F A, McBratney A B. A comparison of prediction methods for the creation of field extent soil property maps [J]. Geoderma, 2001,103(1/2):149-160.

[2] Zhao F J, McGrath S P, Merrington G. Estimates of ambient background concentrations of trace metals in soils for risk assessment [J]. Environmental Pollution, 2007,148(1):221-229.

[3] 郭观林,王 翔,关 亮,等.基于特定场地的挥发/半挥发有机化合物(VOC/SVOC)空间分布与修复边界确定[J]. 环境科学学报, 2009,29(12):2597-2605.

[4] 刘 庚,毕如田,张 朝,等.某焦化场地苯并(a)芘污染空间分布范围预测的不确定性分析 [J]. 环境科学学报, 2013,33(2):587-590.

[5] Kamaruddin S A, Sulaiman W N A, Rahman N A, et al. A review of laboratory and numerical simulations of hydrocarbons migration in subsurface environments [J]. Journal of Environmental Science and Technology, 2011,4(3):191-214.

[6] Robinson T P, Metternicht G. Testing the performance of spatial interpolation techniques for mapping soil properties [J]. Computers and Electronics in Agriculture, 2006,50(2):97-108.

[7] Franke R, Nielson G. Smooth interpolation of large sets of scattered data [J]. International Journal for Numerical Methods in Engineering, 1980,15(11):1691-1704.

[8] Honda H, Haseyam M, Kitajlma H. Fractal interpolation for natural images [J]. Image Processing, 1999,3:657-661.

[9] Li L, Revesz P. Interpolation methods for spatio-temporal geographic data [J]. Computers, Environment and Urban Systems, 2004,28(3):201-227.

[10] 郭艳军,潘 懋,王 喆,等.基于钻孔数据和交叉折剖面约束的三维地层建模方法研究 [J]. 地理与地理信息科学, 2009,25(2): 23-26.

[11] 刘广明,吴亚坤,杨劲松,等.基于光谱指数的区域土壤盐分三维空间变异研究 [J]. 光谱学与光谱分析, 2013,33(10):2758-2761.

[12] Liu C W, Jang C S, Chen S C. Three-dimensional spatial variability of hydraulic conductivity in the Choushui River alluvial fan, Taiwan [J]. Environmental Geology, 2002,43(8):48-56.

[13] 陈 冲,胡克林,贺 勇.农田土壤质地空间分布的三维随机模拟及其不确定性评价 [J]. 土壤, 2013,45(2):319-325.

[14] 姚荣江,杨劲松,赵秀芳,等.沿海滩涂土壤盐分空间分布的三维随机模拟与不确定性评价 [J]. 农业工程学报, 2010,26(11): 91-97.

[15] 吴亚坤,刘广明,杨劲松,等.基于反距离权重插值的土壤盐分三维分布解析方法 [J]. 农业工程学报, 2013,29(3):100-106.

[16] 贺 勇,胡克林,李保国,等.区域土壤质地层次三维空间分布的地统计模拟方法比较 [J]. 土壤, 2010,42(3):429-437.

[17] Meirvenne M V, Maes K, Hofman G. Three-dimensional variability of soil nitrate-nitrogen in an agricultural field [J]. Biology and Fertility of Soils, 2003,37(3):147-153.

[18] Liu C W, Jang C S, Liao C M. Evaluation of arsenic contamination potential using indicator kriging in the Yun-Lin aquifer (Taiwan) [J]. Science of The Total Environment, 2004,321(1-3):173-188.

[19] USEPA (United States Environmental Protection Agency). Method 3050B: Acid Digestion of Sediments Sludges and Soils [Z]. 1996.

[20] Li L, Revesz P. Interpolation methods for spatio-temporal geographic data [J]. Computers, Environment and Urban Systems, 2004,28(3):201-227.

[21] Lu G Y, Wong D W. An adaptive inverse-distance weighting spatial interpolation technique [J]. Computers and Geosciences, 2008,34(7):1044-1055.

[22] Pan G C, Moss K, Heiner T, et al. A Fortran program for three dimensional Cokriging with case demonstration [J]. Computer and Geosciences, 1992,18(5):557-578.

[23] Robinson T P, Metternicht G. Testing the performance of spatial interpolation techniques for mapping soil properties [J]. Computers and Electronics in Agriculture, 2006,50(2):97-108.

The influence of different three dimensional models on Pb distributing prediction in lead-acid battery

contaminated sites.

LIU Geng1, WANG Shang-yi1, NIU Jun-jie1, ZHANG Chao2, ZHAO Xin2, GUO Guan-lin2*

(1.Research Center for Scientific Development in Fenhe River Valley, Taiyuan Normal University, Taiyuan 030012, China；2.State Key Laboratory of Environmental Criteria and Risk Assessment, Chinese Research Academy of Environmental Sciences,Beijing 100012, China). China Environmental Science, 2014,34(12)：3157～3163

In order to define the three dimensional distribution of pollutant and contaminated earth volume in contaminated sites, and provide the relative technical methods and ideas, a typical lead-acid battery contaminated site was studied as an example. The influence of different three dimensional interpolation models (Krig-3D, Inverse Distance Weighted (IDW), Nearest Neighbor (NN)) in different horizontal and vertical anisotropy coefficient setting conditions on Pb pollution assessment and boundary distribution was examined. Different models had different predictive accuracy, the mean error and root mean square error raised, with the increasing of vertical and horizontal ratio. Furthermore, comparison between the results and the remediation target indicated that the contaminated earth volume was largest by IDW statistical model, followed by Krig-3D and NN models, the statistical values differed by 13%. In terms of the contaminated evaluation results and prediction accuracy, the Krig-3D and IDW were better than NN models. Thus, the selections of the appropriate models and rational parameters had very important effects on the accuracy definition of pollution boundary and reducing uncertainty. This result provided important basis for guiding the remediation and management of the contaminated sites.

contaminated site；3-D interpolation；pollution evaluation；uncertainty

X53

A

1000-6923(2014)12-3157-07

刘 庚(1981-),男,山东微山人,讲师,博士,主要从事区域土壤行为过程及其环境效应、土壤有机污染物分布表征和不确定性研究.发表论文近20篇.

2014-06-18

国家自然科学基金青年基金项目(41401236,40901249);国家“863”计划(2013AA06A206)

* 责任作者, 研究员, guogl@craes.org.cn