不同方法下土壤重金属元素空间最大变异距离研究——以浙江上虞为例_刘永生.pdf

第 26 卷第 4 期 2007 年 8 月 岩矿测试 ROCK AND MINERAL ANALYSIS Vol． 26，No． 4 August， 2007 文章编号 02545357 2007 04033105 不同方法下土壤重金属元素空间最大变异距离研究 以浙江上虞为例 刘永生，李瑞敏 中国地质环境监测院，北京100081 摘要 从地质统计学和全域型空间自相关分析两种途径探讨浙江上虞市土壤重金属的空间变异情况。研究 结果表明， 空间分析法得到的上虞市重金属空间变程比地质统计学法得到的空间变程更为直观和稳定; 但空间分 析法的准确程度依赖于空间间隔的大小， 不同的空间间隔计算量有很大的不同， 随着间隔的逐渐减小， 计算量成 倍增加。地质统计学法能快速得到空间变程， 还能推测空间变异程度， 用优化的拟合模型得到的结果更为准确。 关键词 重金属元素; 地质统计学; 空间变异; 空间自相关 中图分类号 O213 P59; S151． 9文献标识码 A 收稿日期 2007- 02- 26; 修订日期 2007- 05- 15 基金项目 国土资源地质大调查 浙江省农业地质环境调查项目资助 20021010; ZD －17 作者简介 刘永生 1976 ， 男， 四川什邡人， 硕士， 工程师， 从事生态地球化学与地下水资源环境方面的研究。 E- mail liuys mail． cigem． gov． cn。 Investigation on Maximum Spatial Variation Intervals of Heavy Metal Elements in Soils with Different s Example from Shangyu City of Zhejiang Province LIU Yong- sheng，LI Rui- min China Institute of Geo- Environment Monitoring，Beijing100081，China Abstract Heavy metal spatial variation in soils was studied in Shangyu City of Zhejiang Province with s of geostatistics and global spatial autocorrelation． The results show that the of spatial autocorrelation is more directly- viewing and stable than geostatistic ． But the precision of the spatial autocorrelation depends on spatial intervals and the computational com- plexity is remarkably increased with decrease of the spatial intervals． Whereas，geostatistic can not only get spatial variation intervals very quickly， but also can speculate the spatial variability degree． And more accurate results can be obtained by geostatistic when fitting model is further optimized． Key words heavy metal element; geostatistics; spatial variability; spatial autocorrelation 在涉及地学空间信息的处理上， 利用传统的数理统计 方法对其进行研究， 由于空间变量特有的空间特性， 不能满 足传统数理统计所要求的随机条件而可能产生误差， 许多 关于时间序列方面的数据分析方法也不适用于空间数据的 处理。因而空间统计学作为一个研究空间变量的学科分支 而产生， 许多学者在这方面进行了深入研究 ［1 －6 ］。 元素作为地学上一个反映土壤的变量， 它与数学上随机 变量的显著区别在于它是一个空间变量， 因而空间统计学在 土壤调查中以及土壤环境质量监测方面有着广泛的应用 ［7 ］。 探讨一地区某一元素与另一元素的相关关系， 假设它们理论 上是线性关系， 然而不同的采样方法得出的线性关系有很大 的不同， 或者线性拟合效果很差。导致这种结果的一个重要 原因是调查区中所获得的元素变量并不满足数理统计学上 要求随机变量间相互独立的条件， 而是相互之间存在着空间 自相关。在土壤环境监测过程中， 面临一个重要的问题是监 测网点的布置问题， 也就是研究区样本容量的确定。人们究 竟在该区布置多少监测点， 这些监测点所获得空间信息在多 大程度上能反映该区的实际情况。这些问题都涉及到对空 间信息的最大相关尺度的探讨， 这个最大空间相关尺度即为 本文要探讨的变程。在变程之内的采样会给分析的变量之 间带来自相关干扰。本文以浙江省上虞市为例， 利用地质统 计学和空间统计学来探讨土壤重金属的空间分布模式， 得到 土壤重金属的最大空间变程。 1研究区域概况与样品采集 1． 1研究区域概况 上虞市地处浙江省东北部， 宁绍平原中部， 曹娥江 下游， 杭甬铁路中段， 东邻余姚市， 南接嵊州市， 西连绍兴 133 ChaoXing 县， 北濒杭州湾， 隔水与海盐县相望。界于东经 1203623″ ～12169″、 北纬 294338″ ～ 301617″。市域南北最长 60 km， 东西最宽 46 km， 总面积 1427． 5 km2。该市地处宁 绍平原中部， 位于上海、 杭州、 宁波三城市形成的三角地带 之中心， 直接受到这三个大城市经济发展辐射的影响， 水、 陆交通十分便利 图 1 ， 杭甬铁路、 杭甬运河、 杭甬高速公 路、 上三高速公路、 104 和 329 国道纵横贯穿全境。上虞市 地形南高北低， 属亚热带季风气候， 气候温和， 四季分明， 湿 润多雨。年均气温16． 4℃， 全年无霜期 250 d 左右。一般 年降雨量约 1 400 mm， 年平均降水天数 160 d， 占全年的 42． 7， 最大年降水量 2 116． 6 mm， 最小年降水量 940． 6 mm， 常年主导风向为南风。 图 1上虞市区位交通图 Fig． 1The traffic map in Shangyu City 1． 2样品采集与分析 以 1 ∶ 5 万地形图为依据， 按1 个样/km2的采样密度在 上虞全区进行土壤样品采集， 共获得1040 个。土壤样品经 过风干， 剔除杂质， 研磨， 过筛进行分析。土壤中重金属元 素测定采用的分析方法及其检出限见表 1。 表 1重金属分析方法 Table 1Analysis s for heavy metals 重金属元素分析方法检出限 Cd等离子体质谱法 ICP － MS20 ng/g As原子荧光光谱法 AFS1 μg/g Hg原子荧光光谱法 AFS2 ng/g Cr压片 － X 射线荧光光谱法 XRF5 μg/g Cu压片 － X 射线荧光光谱法 XRF1 μg/g Ni压片 － X 射线荧光光谱法 XRF2 μg/g Pb压片 － X 射线荧光光谱法 XRF2 μg/g Zn压片 － X 射线荧光光谱法 XRF2 μg/g 2数据分析方法 2． 1地质统计学分析 土壤重金属在空间上的变异也仅仅表现在一定的范围 内， 超出某范围后， 就表现出其相关性。该范围就是本文需 要计算的变程 图 2 。 假设区域化变量满足二价平稳和本征假设 ［8 ］， 则半变 异函数［ γ h ］ 可用式 1 表示 γ h 1 2N d ∑ N d i 1 Z xi － Z xi h[] 2 1 式 1 中， h 表示空间欧式距离; N d 表示被间隔为 h 的样 点数; xi表示第 i 个数据; Z xi 和 Z xi h 分别表示 Z x 在 xi和 xi h 处的值。h 0 时的变异称为块金值， 通常表 示由实验误差或小于实验取样尺度引起的变异，较大的块 金值说明了在较小尺度下存在着重要的地球化学过程。 基台值表示系统内总的变异。 图 2半变异函数图 Fig． 2Diagram of semi- variation functions 2． 2空间自相关分析 在对空间数据的处理研究中， 必须面对的一个很重要 的问题就是数据之间广泛存在空间自相关， 而经典线性回 归模型 CLRM 有 10 项假设条件 ［9 ］， 其中之一就是要求变 量与变量之间是相互独立的， 由于空间信息的自相关特性， 因而经典线性回归模型对此无能为力。对于具有空间自相 关关系的空间数据， 有很多学者对此进行了研究， 提出了一 些空间分析模型 ［10 －12 ］。 为了研究变量的空间模式， 即空间关系和空间自相关， 一些标准的全域和新的局域空间统计方法被开发出来。这 些方法包括 Moran I 法 ［13 ］、 Geary C 法［14 －15 ］、 G statistics 法 ［16 ］、 LISA 法［10 ］和 GLISA 法［17 ］。所有这些空间分析技术 通常包括两方面。第一， 空间模式都是以随机分布为前提; 第二， 所得到的空间模式、 空间结构或者是空间依赖关系的 形式都只是从空间数据中获取的， 而并没有先前预想的 理论概念。 计算空间自相关的方法有很多种。用得最广的主要有 Morans I 法、 GearyC 法、 Getis 法 ［16， 18 －19 ］、 Join count 法［15 ］等 等， 但这些方法作用不是很相同。一般来说， 这些方法从功 用上来说大致可以分为两大类 一为全域型; 另一为局域型。 全域型方法是探索整个空间上属性数据的分布模式以及空 间结构情况; 而局域型空间分析方法可以进一步探索属性数 据在空间上的聚集行为、 聚集位置等空间信息。 在许多研究中， MoranI 空间自相关分析法和 Getis 空 间自相关分析法是最被经常使用的方法。本文主要利用 MoranI 空间自相关分析法来探讨上虞市土壤重金属的空 间变异特征。MoranI 法计算公式见式 2 I d n ∑ n i 1∑ n j≠1Wij ∑ n i 1∑ n j≠1Wij xi － 珋x xj－ 珋 x ∑ n i 1 xi － 珋x 2 2 式 2 中， I d 表示随空间权重间距变化的MoranI值， d 为 空间权重间距; xi和 xj分别为位置 i 和位置 j 的观测值， 珋x 为变量{ xi} 在 n 个位置的平均值; Wij为空间中位置 i 和 位置 j的空间权重矩阵。若 i 与 j 相邻， 则 Wij 1， 否则 Wij0。 233 第 4 期 岩矿测试 http ∥ykcs． i3t． com． cn/ 2007 年 ChaoXing 样本的 MoranI 值的数学期望和方差 s 的计算主要是 依据假设的空间数据分布模式 ［11， 15 ］。 全域型 MoranI 法的统计性质如下。 由于元素在土壤中的分布近似于正态分布， 因而 MoranI利用正态分布模式进行计算 En I － 1 n －1 VarN I n2s1－ ns23s2 0 s2 0 n 2 －1 － E2n I 若要不考虑单位， 则 MoranI 标准化的公式如式 3 Z d I d － E I Var I 槡 3 式 3 中， Z d 为 MoranI 的标准化值， 即将 MoranI 值转化 标准正态分布; I d 表示不同距离下 MoranI 值; E I 表示 MoranI值的数学期望值; Var I 表示 MoranI值的方差。 如果某些位置的观测值在距离 d 的范围内倾向于相 似， 那么 MoranI 是显著的正相关; 如果不相似则为负相关。 如果观测值在空间上完全独立并呈随机分布， MoranI 的值 接近为 0。 2． 3数据处理软件 本文数据处理和图件制作均是在开源自由统计软件 R 2．1．1 环境下进行编程处理。R 是一种为统计计算和图形显 示而设计的语言及环境， 提供了一系列统计和图形显示工 具， 它和贝尔实验室 John Chambers 等开发的 S 系统相 似 ［20 －21 ］。地统计计算利用基于 R 的软件包 gstat －0． 9 －18， 而空间自相关利用软件包 spdep －0．3 －10 进行计算 ［12， 22 ］。 3结果与分析 3． 1土壤重金属地质统计学特征 块金系数为块金值与基台值之比， 表示由随机性因素 引起的空间异质性占系统总变异的比例。块金系数 ＜25、 25 ～75、 ＞75分别表明变量的空间相关性强 烈、 中等、 较弱 ［7 ］。 与实验半变异函数拟合程度较好的理论模型为指数模 型 图 3 ， 拟合度都在 0． 76 以上。从表 2 土壤重金属元素 的半变异函数参数中可以看出， 重金属 Pb 的块金系数值最 小 0． 093 ， 表明 Pb 具有强烈的空间自相关性。其他 7 种 重金属元素的块金系数值均在 25 ～75， 表明其空间自 相关性中等。 表 2土壤重金属元素半变异函数参数 Table 2Parameters of semi- variation function of heavy metals in soils 重金属 元素 理论模型块金值基台值 变程①/km 块金系数拟合度 As指数模型0． 0800． 28015． 1880． 2850． 896 Cd指数模型0． 0720． 1267． 7140． 5690． 763 Cr指数模型0． 0500． 13137． 3740． 3790． 820 Cu指数模型0． 0860． 18213． 0210． 4710． 894 Hg指数模型0． 2200． 90724． 7600． 2420． 910 Pb指数模型0． 0230． 24726． 6790． 0930． 973 Ni指数模型0． 0520． 15523． 6420． 3370． 876 Zn指数模型0． 0390． 08317． 4240． 4650． 871 ① 指数模型计算公式 γ h c 1 － e－ h a ， 变程 3a。 图 3重金属元素的半变异函数图 Fig． 3Semi- variogram of heavy metal elements in soils 3． 2土壤重金属空间自相关特征 图 4 和图 5 是以空间距离间隔为 1 km 计算的上虞土 壤中 8 种重金属元素在 60 km 内的空间自相关特征曲线 图。其中图 4 以系数 MoranI 值I d 为纵坐标。在公式 2 中可以看出， MoranI 计算公式的分子为空间中相邻单 元属性数据的离差乘积。若空间中相邻单元属性数据均大 于平均值或均小于平均值， 则 I d 值均大于 0， 这代表空间 中相邻单元有相似的空间属性， 空间中有聚集现象。若 I d 值小于 0， 表示空间中相邻单元空间属性差异大， 空间 上数据的数据分布为高低间隔分布的状态。I d 值越趋 近于 0， 表示相邻空间单元相关性越低， 在空间上越呈随机 分布。从图 4 可以看出， 随着空间距离的逐渐增大， 8 种重 金属元素的空间自相关强度都逐渐降低， 到接近 0 时变为 随机状态。随着空间距离的进一步扩大， 逐渐过渡为空间 负相关。重金属元素 Cu 和 Cr 在空间间隔为 20 km 左右出 现负相关; Hg 在 27 km 出现负相关; As 在32 km出现负相 关。Cr、 Pb、 Ni、 Zn 等 4 种重金属元素随着空间间隔的增大 逐渐趋于随机分布状态， 因而在整个上虞境内这 4 种元素 可以认为在空间上的分布具有单核心的分布状态; 而 Cu、 Cr、 Hg、 As 等 4 种元素的空间分布特征是多核心的。 333 第 4 期刘永生等 不同方法下土壤重金属元素空间最大变异距离研究 以浙江上虞为例第 26 卷 ChaoXing 图 5 是利用正态模型将 MoranI 值进行标准化［见式 3 ］ 后作为纵坐标的自相关变化图， 是对MoranI值进行了 显著性检验的值， 图中在 Z I 等于 0 的两边两条虚线是在 显著性水平 0． 05 下的区间线。在这两条线构成的区间外 的空间距离上空间相关性可以说具有显著性， 即 0． 05 的显 著性水平下， |Z I |≥1． 96。 图 4重金属元素的空间自相关变化图 Morans I Fig．4Spatial autocorrelation variogram of heavy metal elements Morans I 图 5重金属元素的空间自相关标准化图 1 km 间隔 Fig． 5Spatial autocorrelation standard diagram for heavy metal elements 从图 5 中可以看出 8 种元素的空间变程。Cr 和 Ni 的 空间自相关图很相近， 开始空间关系表现为正相关， 在间隔 24 km 处空间自相关达到最强 曲线驻点， 数据变化缓慢的 地方 ， 随后自相关性逐渐降低， 数据分布呈现随机分布状 态。重金属元素 Zn 在空间间隔 13 km 左右空间自相关达 到最强， 随后空间自相关性降低， 到空间间隔为 31 km 处达 到最低， 随后空间自相关性又随空间间隔的增大而增强; Cd 和 Hg 空间自相关分别在8 km 和16 km 处达到最强; Pb 的空间自相关曲线的驻点出现在 26 km 处， 也就是说， 在空 间显著相关的状态 在 95 区间线之外 下， Pb 在 26 km 处空间自相关开始变得比较缓慢， 也就是达到了空间自相 关的最大变程， 因而 26 km 为 Pb 的空间最大变程。同理 As 和 Cu 分别在 15 km 处和 12 km 处曲线变化， 相应的变 程也为 15 km 和 12 km。 表3 是从地质统计学法和空间分析法两种途径得到的 上虞市土壤重金属的空间最大变程， 其中 As、 Cd、 Cu、 Pb 和 Ni 等5 种重金属元素两种途径得到的空间最大变程都很接 近， 而 Cr、 Hg 和 Zn 等3 种重金属元素两种方法得到的空间 最大变程相差较大， 其中 Cr 相差最大， 然而从 Cr 的空间自 相关图上， 很容易得出 Cr 的最大变程是 24 km 的结论。 对于地质统计学法， 虽然其拟合曲线的决定系数能达到 0．82， 却得到的是37．4 km， 两者差距比较大。这应该是地质 统计学法的拟合模型不是很好而产生较大的误差所致。Zn 元素的空间自相关图显示有两个峰， 曲线在 0． 05 的显著性 水平线上有3 个驻点 平缓点 ， 即13 km、 30 km 和49 km处。 这3 个驻点可以在不同尺度上表征 Zn 的空间变程。 表 3不同途径求得的重金属元素空间变程 Table 3Spatial variability distance of heavy metal elements by different ways 重金属 元素 变程/km 地质统计学法 空间分析法 重金属 元素 变程/km 地质统计学法 空间分析法 As15． 215Hg24． 816 Cd7． 78Pb26． 726 Cr37． 424Ni23． 624 Cu13． 012Zn17． 413 433 第 4 期 岩矿测试 http ∥ykcs． i3t． com． cn/ 2007 年 ChaoXing 4结果与讨论 在研究空间数据的过程中， 地质统计学法和空间自相 关分析法是两种重要的面向数据的处理方法， 两种方法从 不同的途径去探索空间信息。本文从两种途径探讨了上虞 市土壤重金属空间变异特征。从两种途径得到的空间变程 来看， 空间自相关分析法得到的空间变程最为稳定， 也最为 直观， 而地质统计学法得到的空间变程随着拟合模型的不 同其结果变化较大。两种方法可以互相比较， 都有各自的 优缺点。地质统计学法是通过对模型的拟合来得到重金属 元素的空间最大变程， 其优点是若模型选择得好， 拟合效果 就比较好， 得到的最大变程也就越准确; 但由于空间数据的 复杂性， 往往很难找到能够拟合的非常好的模型， 有时常用 的经验模型都拟合得很差， 而导致最后结果不准确。地质 统计学法除了能得到空间变程信息外， 还能根据模型得到 的基台值和块金值推测重金属元素空间变异程度; 全域型 空间自相关分析法是通过空间单元间的邻接关系来分析空 间数据的特征， 它不但能得到空间属性数据的正相关尺度， 还能得到负相关尺度; 但由于其计算方法复杂， 得到的属性 数据空间变程的精度依赖于计算的空间距离， 计算的空间 距离越短， 数据越多， 计算时间也就越长。 5参考文献 ［ 1］Anselin Luc，Getis Arthur． Spatial Statistical Analysis and Geo- graphic Ination Systems［J］ ． Annals of Regional Science， 1992 26 19 －33． ［ 2］Goodchild M F， Haining R P， Wise S． Integrating GIS and Spatial Data Analysis Problems and Possibilities［ J］ ． International Journal of Geographical Ination Systems， 1992， 6 5 407 －423． ［ 3］Haining R P． Spatial Data Analysis in the Social and Environ- mental Science［M］． Cambridge Cambridge University Press， 1990 291 －312． ［ 4］Odland J． Spatial Autocorrelation［M］ ． California Sage Publica- tions Inc， 1988 25 －55． ［ 5］Paelinck，Jean H P，Klaassen，et al． Spatial Econometric［ M］ ． England Saxon House， 1979 211． ［ 6］Yue Hone Chou． Map Resolution and Spatial Autocorrelation ［ J］． Geographical Analysis， 1991， 23 3 228 －246． ［ 7］Cambardella C A，Oorman A T，Novak J M，et al． Field- scale Variability of Soil Properties in CentralIowa Soils［J］ ． Soil Sci Soc Am J， 1994 58 1501 －1511． ［ 8］赵鹏大． 定量地学方法及应用［M］． 北京 高等教育出版社， 2004 78 －82． ［ 9］Gujarati D． 计量经济学 Basic Econometrics ［M］． 北京 中国 人民大学出版社， 2000 47 －55． ［ 10］Anselin Luc． Local Indicators of Spatial Association- LISA［J］ ． Geographical Analysis， 1995， 27 2 93 －115． ［ 11］Goodchild M F． Spatial Autocorrelation［M］． CATMOG 47， Geobooks Norwich UK， 1986 6 －25． ［ 12］Pebesma E J， Wesseling C G Gstat． A Program for Geostatistical Modelling，Prediction and Simulation［J］ ． Computers ＆ Geosci- ences， 1998， 24 1 17 －31． ［ 13］Moran P A P． Notes on Continuous Stochastic Phenomena［J］ ． Biometrika， 1950 1/2 17 －23． ［ 14］Cliff A，Ord J Keith． Spatial Autocorrelation［M］． London Pion， 1973 7 －17． ［ 15］Cliff A，Ord J Keith． Spatial Processes Models and Applica- tions［M］． London Pion， 1981 6 －22． ［ 16］Getis Arthur，Ord J Keith． The Analysis of Spatial Association by Use of Distance Statistics［J］． Geographical Analysis， 1992， 24 3 189 －206． ［ 17］Bao Shuming，Mark S Henry． Heterogeneity Issues in Local Measurements of Spatial Association［ J］ ． Geographical Systems， 1996 3 1 －13． ［ 18］Ord J Keith，Getis Arthur． Testing for Local Spatial Autocorre- lation in the Presence of Global Autocorrelation［J］ ． Journal of Regional Science， 2001， 41 3 411 －432． ［ 19］Ord J Keith，Getis Arthur． Local Spatial Autocorrelation Statis- tics Distributional Issues and Application［J］． Geographical Analysis， 1995， 27 4 286 －306． ［ 20］Stephen P Kaluzny，Silvia C Vega，Tamre P Cardoso，et al． S SpatialStats Users Manual for Windows and UNIX［M］ ． New York Springer， 1998．327． ［ 21］Venables W N，Ripley B D． Modern Applied Statistics with S ［ M］ ．4th ed． New Jersey Springer， 2002 41 －66． ［ 22］Ripley B D． Spatial Statistics in R［M］ ． New York Wiley- Interscience， 2001 Vol．1/2  14 －15． 2007 全国地质分析学术报告会通知 为推动全国地质分析技术发展， 促进学术交流， 中国地 质学会岩矿测试技术专业委员会定于 2007 年 10 月举办全 国地质分析学术报告会。 本届学术报告会将以形态分析及环境地球化学分析为主 题， 关注环境和人类健康， 并涵盖地质、 环境、 材料等诸多学科 1 形态分析技术与应用; 2 地质调查及生态环境地球化学评价中的分析技术; 3 有机地球化学分析; 4 岩石矿物、 材料分析技术; 5 环境与生物分析技术。 大会将邀请国内外著名学者作特邀报告， 欢迎踊跃 参加。 联 系 人吴晓军 电 话 010 －68999770传真 010 －68998605 电子邮件wuxiaojun cags． net． cn 邮 编 100037 地 址 北京百万庄大街 26 号 国家地质实验测试中心 中国地质学会岩矿测试技术专业委员会 533 第 4 期刘永生等 不同方法下土壤重金属元素空间最大变异距离研究 以浙江上虞为例第 26 卷 ChaoXing