全 文 ：应用转向带法进行土壤铜和铅含量的条件模拟*
王学军 (北京大学城市与环境学系, 地表过程分析与模拟教育部重点实验室, 北京 100871)
摘要 以北京东郊污灌土壤中铜和铅为例, 探讨了条件模拟中的转向带方法在研究土壤中重金属空间
分布特征的应用, 从统计特征、空间结构等方面考察了模拟结果 ,分析了误差产生的原因. 结果表明, 空间
模拟方法可以很好地再现土壤重金属含量的空间结构特征, 为揭示土壤重金属含量的空间结构特征提供
了一个有力的工具. 应用该方法有助于了解土壤、水体和大气环境中各种物质的来源、迁移转化路径以及
最终的分布特征, 从而为环境评价、规划及解决有关问题提供技术支持.
关键词 土壤 转向带法 条件模拟 重金属
文章编号 1001- 9332( 2002) 12- 1667- 04 中图分类号 X524 文献标识码 A
Application of turning band method in conditional simulation of soil heavy metals. WANG Xuejun ( Depar t
ment of Urban and Environmental Science , Peking University , Beij ing 100871) . Chin. J . A pp l . Ecol . , 2002,
13( 12) : 1667~ 1670.
Geostatistical simulation techniques are widely used in studying spatial phenomena. The purpose of this study is
to make an experimental study in apply ing this technique to understand the spatial distribution of copper and lead
in w astew ater irr igated ag ricultural so il in the eastern suburb of Beijing . The results showed that the simulation is
efficient, because it fits to the main char acteristics of the r evealed reality of the spatial distribution of these ele
ments. Geostatistical simulat ion could play an important role in regional environmental and ecolog ical studies.
Key words Soil, Turning band method, Conditional simulation, Heavy metal.
* 国家自然科学基金资助项目( 40041002, 40024101) .
* * 通讯联系人.
2001- 06- 26收稿, 2002- 03- 01接受.
1 引 言
地统计模拟和 Krig ing 插值是区域空间特征研
究中常用的方法[ 1, 5] . Kriging 插值可以给出有限区
域内区域化变量的最佳线性无偏估计量, 因此在地
质学、生态学等领域得到广泛应用. 但是, 该方法也
有其缺点, 其中之一就是不能很好地再现区域化变
量的空间结构, 因为 Krig ing 法平滑了真实数据. 应
用条件模拟的方法,则能正确反映区域化变量的空
间波动特点[ 3, 8] .
地统计模拟的模拟现实具有与实测数据相同的
频率直方图.更重要的是,模拟现实与真实现实具有
相同的空间自相关结构. 这实际上是构造地统计模
拟各种方法的基本要求[ 6] . 因为区域化变量在空间
各点的值既有差异又相互影响, 形成了特定的空间
结构,模拟现实只有再现了真实现实的空间结构, 才
能在一定程度上代表真实现实,从而具有预测功能.
地统计模拟自提出以来发展很快, 产生了多种
方法,应用领域也不断拓宽,在石油勘探、采矿、水文
地质、工程地质等方面均有成功应用的例子. 例如,
Journel 和 Issaks[ 4] 用 指 示 变 量 法 对 加 拿 大
Saskatchew an 省一铀、砷交错矿进行条件模拟,满足
了生产提出的各项要求, 节约了开支, 效果很好. 又
如,在石油勘探开发中,孔隙度和渗透率具有重要意
义.周叶、王家华[ 10]讨论了序贯指示模拟在估计渗
透率空间分布中的应用原理和特点, 进行了模拟重
构实验,将处理结果和普通 Kriging 法、指示 Krig ing
法的处理结果进行了对比分析. 该方法在含油气储
集层的评估方面应用前景甚广. 但是,地统计模拟方
法在环境学和生态学研究中的应用尚不多见.
地统计模拟研究的开端是以转向带法的提出为
标志的[ 7] . 在该法的基础上, 出现了谱方法、指示模
拟方法、序贯方法等多种方法,使转向带法应用领域
不断拓宽,但转向带法仍是最主要的模拟方法.
2 材料与方法
2 1 样品采集与分析
土壤样品采自北京东郊通惠河畔的污灌区, 采样时间为
1995 年 7 月. 该区土壤为潮土, 中壤土, 呈弱碱性反映, pH
值在 7. 1~ 8. 5, 有机质含量在 0. 8% ~ 3. 7%之间. 种植水稻
和蔬菜, 采样区总面积约 0. 5km2 . 自 20 世纪 60 年代以来,
一直利用通惠河和南北灌溉渠进行污水灌溉, 地形平坦, 土
壤质地也较为均一. 采样按均匀布点方式进行,共采集 10 行
19 列共计 190 个表土样品, 采样深度为 20cm, 南北向采样点
应 用 生 态 学 报 2002年 12 月 第 13 卷 第 12 期
CHINESE JOURNAL OF APPLIED ECOLOGY , Dec. 2002, 13( 12)!1667~ 1670
间距为 50m,东西向采样点间距为 40m. 土样在室内风干, 过
1mm 尼龙网筛, 用重铬酸钾法测定有机质含量, 用比重计法
测定< 0. 01mm 颗粒组成.少量样品用玛瑙研钵碾至 100 目
以下,用硝酸高氯酸氢氟酸对土壤样品进行消化, 提取液
中铜和铅含量用岛津( 18080)原子吸收分光光度计测定[ 9] .
各元素含量和土壤理化参数数据多呈对数正态分布. 用
Grubbs 法进行异常值检验. 为减少信息损失, 显著性水平取
5% , 仅有少量异常值被剔除.
22 研究方法
地统计模拟一般是指条件模拟,即模拟值以临近点的实
测数据为条件.如果模拟点与某一实测点临近,则模拟值与
该点的实测值相近;如果模拟点与某一实测点重合, 则模拟
值与该点的实测值相等. 这样, 模拟现实不仅空间变化规律
上与真实现实相符,而且量值与实际相近, 因而对特定问题
更具实际意义.
条件模拟中的转向带法不直接模拟二维(或三维)过程,
而是先沿若干条直线独立地模拟一维过程.然后,对二维(或
三维)空间的每一点, 用已模拟的若干个一维过程的相应值
的加权和作为其模拟值. 使用该方法有 3 个前提假设, 一是
待模拟的随机场是二阶平稳的;二是各点的值服从均值为零
的正态分布.若此假设不真, 可先作标准正态变换;三是待模
拟的随机场协方差函数已知.
本文应用条件模拟方法中的转向带法研究了铜和铅的
含量.对 10m ∀ 10m 规则网格进行了 10 次条件模拟, 模拟点
为 46 行 73列共计 3358 个.
模拟过程中半变异函数拟合利用了美国环保局开发的
GEOEAS 1. 2. 1 软件[ 2] , 其余计算全部用 Visual Basic 3. 0
语言编制程序, 在 PC 机上运行通过. 结果图件分别用 MS
Excel 5. 0 和 Surfer 5. 01 绘制.
对模拟过程中几个问题的说明: 1)对漂移的处理 漂移
是指非平稳点随机函数的非平稳数学期望, 是位置的函数.
本文研究的变量存在较明显的漂移,因此用趋势面拟合其漂
移部分,然后对残差进行模拟. 模拟完成后,将残差模拟值与
趋势值相加 ,得到原始数据的条件模拟.趋势面采用多项式
曲面. 多项式的次数考虑两个因素之后确定, 其一是实际需
要,因为现在的问题是要把局部构造从总体趋势中分离出
来,得到较平稳的残差, 并不是要探讨区域性的趋势变化, 故
而可根据对拟合精度的要求取较高次数; 其二是拟合精度,
通过计算 3~ 7 次趋势面, 比较残差平方和,决定采用 5 次曲
面进行拟合.
2)正态变换 因为转向带法是将确定在若干直线上的
很大数目(理论上为无穷多个, 本文中用 6 个)的独立的现实
相加, 因此根据中心极限定理, 得到的随机函数的模拟现实
将具有高斯分布特点.文中需要模拟的随机函数的一元分布
近似对数正态分布,因此先把数据进行高斯变换, 根据由变
换后的数据计算的协方差函数进行模拟,并用变换后的数据
对非条件模拟结果进行条件化,最后再对条件模拟结果进行
逆高斯变换.对趋势面残差进行 BoxCox 变换.进行 BoxCox
正态变换时要求数据为正, 因此根据其极值分布情况, 使残
差统一增大一定值(不会改变数据的自相关结构) , 相应地,
把趋势面向下平移一定值.
3)结构分析 分析变量空间结构的主要工具是半变异
函数. 半变异函数最常见的模型有球状模型、高斯模型、指数
模型、线性模型等.本文对两种元素半变异函数的拟合采用
球状模型. 铜的基台值和变程分别为 3. 3 和 73,铅的基台值
和变程为 2. 5 和 60. 所研究的变量可近似为各向同性.
4)计算滑动平均的权重因子 根据协方差函数计算滑
动平均的权重因子, 采用 NewtonRaphson 迭代方法. 求得权
重因子之后, 进行一维模拟.
5)转向带法迭加 首先采用滑动平均在 0#, 30#, 60#,
90#, 120#和 150#6 条直线上分别独立地产生一维模拟, 然后
通过矢量内积求得模拟点在各条直线上的投影位置,并求出
该直线在该位置上的一维模拟值, 最后进行转向带迭加.
3 结果与分析
图 1显示了模拟结果及原始数据的描述数据离
散特征的基本统计量, 从中可以看出,模拟结果极端
值范围比原始数据和 Krig ing 插值结果的相应范围
明显宽些,但大部分数据(图中反映了一半的数据)
所处的范围二者相近; 模拟值和原始值的均值及中
位数相近;另据检验,在控制点上模拟值与实测值严
格相等. 以上特点显示出条件模拟在实测数据控制
下的波动性, 并反映出模拟有效地再现了原始数据
的基本统计特征.
限于篇幅,图 2给出了土壤中铜和铅的部分模拟
图 1 铜和铅的基本统计量
Fig. 1 Basic stat istics of copper and lead.
1)最小值M inimum; 2 )第一四分位数 First quart ile; 3 )中位数 Medi
an; 4)均值 M ean; 5) 第三四分位数 Third quart ile; 6) 最大值 Maxi
mum .
1668 应 用 生 态 学 报 13卷
图 2 铜、铅的原始数据和几次模拟结果的对比
Fig. 2 Comparison of original data and simulated results of copper and lead.
a)原始数据Original data; b)模拟 2 Simulat ion No. 2; c)模拟 3 Simulat ion No. 3; d)模拟 5 Simulat ion No. 5; e)模拟 7 S imulation No. 7; f )模拟 8 S im
ulat ion No. 8; g)模拟 9 Simulat ion No. 9.
结果.它反映了原始数据的空间结构特征,再现了两
种重金属元素含量的空间变化特点.
图 3给出了模拟结果的半变异函数, 由此可以
看出, 模拟结果的半变异函数与原始数据符合的很
好,而 Kriging 插值的半变异函数则与原始数据有
较大差异.
是否符合原始数据的半变异函数, 正是模拟与
插值的重要不同之处. 模拟与估计(插值)是研究区
域化变量时常用的两种方法, 二者在理论意义和应
用领域上均有不同.估计的主要目的是在每一点 x
提出一个尽可能接近未知的真实值的估计量, 估计
的质量用无偏性和最小均方误差或者估计方差来衡
量;估计有明显的∃平滑∃效应,改变了变量的空间结
构.模拟的目的主要是在模拟值接近真实值的前提
下,再现变量的空间变化特点,而模拟值并不是最优
估计值. 因此当空间变化对所研究的问题有重要意
义时,就要用模拟值,而不是用估计值.
模拟与估计的另一点不同在于, 对研究空间的
一点(不是控制点) ,估计值是一定的,模拟值则可以
有多个,事实上, 模拟现实的个数没有任何限制. 因
此,经常把随机场的模拟现实作为Monte Carlo 模拟
或其它模型的输入数据(或参量) , 由同一数据系列
的多个条件模拟现实得到多个模型输出, 可计算输
出的统计规律, 也可分析输入(或参量)与输出的统
166912 期 王学军: 应用转向带法进行土壤铜和铅含量的条件模拟
图 3 铜和铅在不同方向上的半变异函数
Fig. 3 Variograms of copper and lead in different directions.
a)东西方向EW direct ion; b)南北方向 SN direct ion. % .模拟现实 S imulated real it y; & . Kriging 插值Kriging interpolat ion ; ∋ .实测数据 Measured
data.
计关系,空间不确定性对输出的影响,以及模型的灵
敏度等.
因此,相对而言,应用条件模拟的方法更能正确
反映区域化变量的空间波动特点.这对于研究土壤
中重金属的来源、迁移转化以及最终分布具有重要
意义.相比之下, Krig ing 插值则更适于进行物质含
量水平的空间估计.
就本研究而言, 采样区虽然面积很小,但污灌口
有多个,由于土壤中大量存在的有机质以及粘土矿
物对铜、铅等污染物质的吸附作用.这些污染物在土
壤中的迁移路径很短, 造成了其含量在空间上的明
显差异.例如, 采样区东北角的污灌口造成了该处
铜、铅含量的偏高.
4 结 论
利用转向带法对北京东郊污灌区污灌土壤中重
金属的空间分布进行了模拟研究, 从基本统计特征、
空间结构等方面考察了模拟结果,分析了产生误差
的原因.研究结果表明, 与 Kriging 插值方法相比,
空间模拟方法可以很好地再现土壤重金属含量的空
间结构特征,为揭示土壤重金属含量的空间结构特
征提供了一个有力的工具.但是,如果研究目标是估
计土壤物质含量分布,则 Kriging 方法更为可取. 应
用空间模拟手段有助于了解土壤、水体和大气环境
中污染物的来源、分布, 从而为环境评价、规划及解
决有关问题提供技术支持.
参考文献
1 Black T C, Freyburg DL. 1990. Simulat ion of onedimensional corre
lated f ields using a mat rixf actorisat ion moving average approach.
Math Geol , 22( 1) : 39~ 62
2 Englund E, Sparks A. 1991. GEOEAS 1. 2. 1 User( s Guide. Las
Vegas: U. S. EPA.
3 Journel AG, Alabert F. 1990. New method for reservoir mapping. J
Pet Tech nol , 42( 2) : 212~ 218
4 Journel AG, Isaaks EH. 1984. Condit ional indicator simulation: Ap
plicat ion to a Saskatchew an uranium deposit . Math Geol , 16 ( 7 ) :
685~ 718
5 Kentw ell DJ, Bloom LM , Comber GA. 1999. Im provements in
grade tonnage curve predict ion via sequent ial Gaussian f ractal simu
lat ion. Math Geol , 31( 2) : 311~ 325
6 Kentw ell DJ , Bloom LM , Comber GA. 1999. Geostat ist ical condi
tional fractal simulat ion w ith irregularly spaced data. Math Comp u
S im , 48: 447~ 456
7 Matheron G. 1973. The int rinsic random fun ctions and their appl i
cat ions. Ad v Ap pl Prob, 5: 439~ 468
8 Sen Z. 1990. Spacial simulat ion of geological variables. Math Geol ,
22( 2) : 175~ 188
9 Wang XJ (王学军 ) , Deng BS (邓宝山) , Zhang ZP (张泽浦 ) .
1997. The spat ial st ructure of t race element in surface soils in east
ern suburb of Beijing. A cta Sci Cir c( 环境科学学报) , 17 ( 4) : 412
~ 416( in Chinese)
10 Zhou Y(周 叶) , Wang JH (王家华) . 1991. Method for Est imat
ing the Spat ial Dist ribution of Permeabilit y. Beijing: Geology Press.
145~ 156( in Chinese)
作者简介 王学军, 男, 1964 年生, 教授, 博士生导师, 主要
从事环境地学研究, 发表论文 100 余篇, 专著 5 部. T el: 010
62759190, Email: x jw ang@ urban. pku. edu. cn
1670 应 用 生 态 学 报 13卷