莱州湾南岸土壤盐分空间变异性指示克立格阈值研究

（河南工程学院，河南郑州市，491191） 毕延凤

莱州湾地区有浓度最高、分布最广、储量最大的滨海相地下卤水，是我国重要的海盐生产基地和海水养殖基地，还是山东省粮食高产区和全国蔬菜市场的重要基地。这些经济农业活动扩大了咸水入侵范围，2007年7月开始，在国家908专项支持下，对莱州湾沿岸海水入侵严重地区进行了为期两年的咸水入侵灾害监测，结果显示莱州湾地区咸水入侵范围已达2 500km2。咸水入侵导致大量土地发生盐渍化，植物群落和生态系统相应发生改变。

土壤盐分高值区以特异值形式存在，数据往往不呈正态分布，导致变异函数不稳健[1]，传统的高斯克里格估值存在困难。Journel于1983年创立的指示克里格(IK)，不要求数据服从正态分布且不受特异值的影响，它对连续数据通过阈值变换转变成0和1的二进制数据，并以指示变异函数为基本结构分析工具，预测精度超过了普通克里格法[2]。IK法自创立以来广泛应用于地质[3]、土壤[4、5]、水文[6、7]等领域的风险预测研究中，但多为局部估计的初步探讨，对IK法阈值的选择及与函数关系的系统研究尚少见[8、9]。

以莱州湾南岸潍坊北部平原的土壤盐分空间变异为案例，以ArcGIS9.3、GS+7.0为操作平台，重点分析了不同分位数阈值下盐分空间结构、IK插值精度及指示阈值与条件概率分布函数的相互关系和取值规则，为IK法的应用提供参考。

1 研究区概况

研究区位于渤海莱州湾南岸，潍坊北部地区，面积1 800km2。地势南高北低，地形平坦广阔，属滨海堆积平原地貌。自西向东依次有弥河、白浪河、潍河和胶莱河由南向北注入渤海，海岸线为东南－西北走向，呈弧形曲线状。本区属干旱半干旱气候，年均降雨量652.8mm，年均蒸发量1 774mm，蒸发量远大于降水量。

自南至北分布着棕壤、褐土、潮土、砂姜黑土和盐土5大土类，下分15个亚类、34个土属、110个土种。其中盐土类面积约149×104亩，占可利用土壤面积的7.43%，集中分布在市域北部滨海地带，即寿光、寒亭、昌邑3县、区的北部。盐土土类下分滨海潮盐土亚类及2个土属、8个土种。

2 材料与方法

自西向东沿垂直海岸线方向共布设了13条监测剖面，每条剖面上设置监测点5～7个不等，均由海水入侵区过渡至非海水入侵区，共计80个，各相邻剖面及相邻点之间的直线距离约2～3km（图1）。取样点大部分选在农田，亦有部分采样点选在盐荒地、林地及草地。土样的采集利用手动土壤采样器，每个监测点取地表以下1m土样，分6层，分别为0～10cm、10～20cm、20～40cm、40～60cm、60～80cm及80～100cm，每份样品重量约1kg，编号，带回实验室，风干、研磨、过筛，进行土壤全盐含量的测定[10]。

指示克里格法的原理详见文献[5、7]，首先确定评价指标的阈值，在此选取20、30、40、50、60、70、80分位数的盐分阈值，并为评价指标确定指示函数，将采样点数据进行二态指示变换（指示变换值1或0），然后将指示变换值在GS+7.0中进行变异函数模拟计算，最后将得到的最佳变异函数模型参数输入ArcGIS9.3进行IK插值，得到土壤盐分满足一定阈值下的概率空间分布图。

经典统计采用SPSS，变异函数的模型拟合采用GS+7.0，概率空间分布图的绘制以及平均概率、预测概率均方根误差（RMSE）的计算采用ArcGIS9.3。

3 结果与分析

3.1 土壤盐分含量的描述性统计

从表1中K-S_P值可知，20～40cm和80～100cm土壤盐分经自然对数转换后K-S检验的p值＞0.05（双截尾），服从对数正态分布，其余土层以及0～100cm加权平均含盐量（表1最后一行）的分布均呈明显的偏态分布。土壤剖面中，上层0～10cm，10～20cm和底层80～100cm土层土壤盐分含量的最大值和最小值差别最大。上层0～10cm和10～20cm盐分均值为中度盐渍化，往下层盐分有增加的趋势，呈现重度盐渍化。标准差也呈现往下层增加的趋势。从各层变异系数来看，范围为2.767～3.699，均为强变异，20～40cm变异系数最大，其他各层总体呈现往下层递增的趋势。

在秋季持续干旱的地表蒸发作用下，表层0～10cm全盐量均值大于10～20cm，而总体剖面呈底聚类型，这与莱州湾南岸广大平原区除紧邻莱州湾的沿海狭窄条带，由于卤水开发逐年降低的地下水位，干旱少雨的气候和灌溉水的紧缺密切相关。同时也体现了盐渍化程度正趋于和缓。

3.2 不同阈值下土壤盐分空间结构的变化

用0～100cm加权平均含盐量为研究对象，由于数据不呈正态分布，拟对其进行指示变异函数的研究。

利用表2中的盐分阈值，分别对土壤盐分实测值进行二态指示变换，在GS+7.0中设有效滞后距67 667.22m，步长大小为6 766.72m，步长数为10，模拟获得不同阈值下的盐分指示变异函数最佳模型及其参数。

从表中可以看出，土壤盐分的空间相关性总体上随着盐分阈值的增大而增强，在第50百分位数之前均呈中等强度的空间相关性，盐分从第60百分位数开始大多呈现强烈的空间相关性。RSS都很小，尤其是在第50分位数之前，在第30分位数处达最小。模型的决定系数R2在第50百分位数之前在0.758～0.995之间，模型拟合精度高，盐分第60百分位数之后R2在0.220～0.386之间，拟合精度较低。

表1 1m以内土体各层及各层平均土壤盐分含量的经典统计

表2 不同盐分阈值下土壤盐分的指示变异函数模型

表3 土壤盐分指示阈值与指示变异函数的关系

3.3 指示阈值与指示变异函数的关系

由表3可见，当盐分阈值由小逐渐增大时，指示变异函数值有一种由小到大、再到小的趋势，转折点在第50分位数（盐分阈值0.13%）处，可见中位数阈值附近的指示变异函数结构性最差，其反应的指示概率的空间变化性最大，由此推测其指示概率的预测误差也将最大。

3.4 指示阈值与指示概率函数的关系

指示概率均值随着盐分阈值的增大逐渐减少，由83.4%减少到27.6%，从第60分位数（盐分阈值0.16）开始减小速度变慢；均方根预测误差（RMSE）有先增加后减小的趋势，由38.4%迅速增加至45.9%，后随着阈值增大缓慢减小至38.7%，在第50分位数（盐分阈值0.133）处概率预测误差最大，与上文提及的指示变异函数在此值处最大一致，采用此阈值会降低风险评价与IK法估计的可靠性，在第20分位数处（盐分阈值0.09）概率预测误差最小。可见，单从IK法插值精度考虑，中位数盐分阈值并不是最佳阈值。这点与杨奇勇等对山东禹城市[9]、刘全明等对内蒙古河套灌区[8]利用IK法进行的土壤盐分空间变异性研究结果一致。但在杨与刘等的研究中，由于其研究区中位数以下的盐分阈值小于0.1%，属非盐渍化，其概率分布图对盐渍化改良参考意义不大，因此建议研究区盐分阈值应尽量在精度容许情况下选取超过中位数的偏大值。而本文的研究中，中位数以下的30和40分位数处(盐分阈值分别是0.10%和0.11%)均属轻度盐渍化范畴，尤其是30分位数处，半变异函数模型拟合精度最高，土壤盐分呈中等强度的空间相关性；超过中位数阈值的RMSE相对较小的第70、80百分位数处，土壤盐分呈强烈的空间相关性，而其半变异函数的拟合精度较低。

3.5 指示概率分布图

在盐碱土壤改良过程中，更感兴趣的往往不是某一点处土壤盐分的具体值，而是盐分大于某一阈值在空间上的分布概率。用IK法做出不同盐分阈值下的指示概率分布图如图2。

在图2中，随着盐分阈值的增大，高概率面积逐渐减小，低概率面积逐渐增大，导致平均概率逐渐由盐分阈值0.09%的83.1%减小到盐分阈值0.35%的27.6%，概率值分布的偏度由盐分阈值0.09%的负偏（偏度为-0.9 775）过渡到盐分阈值0.35%的正偏（偏度为0.6354），整体呈逐渐增大趋势；概率值分布的峰度在盐分阈值0.09%、0.10%、0.11%和0.13%处普遍偏大，在中位数阈值0.13%处达到最大，而在盐分阈值0.16%和0.22%处急剧降低，在盐分阈值0.35%处又开始回升，峰度随着盐分阈值的变化规律与研究区土壤盐分的分布状况密切相关。从概率值分布的偏度和峰度看，在中位数阈值0.13%处最接近正态分布（偏度-0.929，峰度3.252）。随着盐分阈值的增大，主要概率分布区间也由盐分阈值0.09%的0.9～1.0过渡到盐分阈值0.35%的0.0～0.1。当盐分阈值0.10%时，高概率区主要分布在广大的西部以及北部地区，潍坊北部是主要盐场所在地，而西部则主要与广泛分布的大棚蔬菜种植区导致的次生盐渍化密切相关。在盐分阈值0.22%处，盐分高概率区主要分布在北部盐田区以及胶莱河和白浪河河谷地带。某一盐分的高概率分布区域一般包含在小于它的上一个盐分阈值的高概率分布范围内，相应的低概率分布区域一般包含在大于它的上一个盐分阈值的低概率分布范围内。盐分阈值可根据不同阈值目标下特定区域的概率与风险进行选择。

莱州湾南岸地区盐渍土的形成既有自然因素也有人为因素，较高的气温及较大的蒸发量使得受海水（或卤水）浸染的矿化度较高的地下水通过蒸发作用及毛细作用将水中盐分聚积在地表及地表下的土层中，形成了盐渍土。除自然因素外，人类活动是造成土壤次生盐渍化的重要原因。莱州湾是山东省粮食高产区和全国蔬菜市场和水产养殖的重要基地，自然条件优越。但过大的人口密度（是全国人口平均密度的4倍多）及非农业用地比例的增加使得人口环境压力十分突出。近年来随着沿海改革开放战略的实施，在各县市经济快速发展的同时，水资源的需求量猛增。由于过度开采地下水，使得地下水位以1～3m/yr的速度迅猛下降，形成大于2 500km2的地下水位负值区。此外莱州湾地区是我国重要的海盐生产基地和海水养殖基地，海水养殖和抽卤晒盐等经济活动把大量海水和卤水引入陆地，人为地促进海水或咸水向内陆运移5～15km2，扩大了咸水入侵范围。咸水入侵一方面直接浸渍土颗粒使其盐渍化；另一方面，为了农田保产，人们不得不使用一部分含盐量较高的地下半咸水进行灌溉，在没有良好的排灌系统的情况下长期使用高矿化度的水进行灌溉，盐分不断在土壤表层聚集，使耕层土壤盐分含量增加，必然导致了次生盐渍化的发生。据各县市区的统计资料，目前由海水入侵造成的次生盐渍化面积已达2.5×104hm2，受害土地面积达1.01×105hm2。土壤盐渍化的结果是土壤结构变差，理化性能降低，土地资源退化，严重的丧失生产能力，部分农田因地下水变咸而被弃耕成为荒地，栽培作物逐渐为盐生植被所替代，从而进一步加剧了耕地资源的危机，也增加了区域生态恢复和生态开发的难度。

4 结论

莱州湾南岸潍坊北部的广大平原地区土壤盐分呈偏态分布，利用IK法，选择不同分位数盐分阈值，对土壤盐分空间结构、盐分阈值与指示概率和指示概率均方根预测误差的关系以及IK指示概率插值图进行了研究。

随着土壤盐分阈值的增大，盐分空间相关性总体呈增强趋势，拟合精度降低，指示变异函数值有一种由小到大、再到小的趋势，转折点在第50分位数处，可见中位数阈值附近的指示变异函数结构性最差，其反应的指示概率的空间变化性最大，导致其指示概率均方根预测误差也最大，采用此阈值会降低风险评价与IK法估计的可靠性。中位数以下的30和40分位数处(盐分阈值分别是0.10%和0.11%）均属轻度盐渍化范畴，尤其是30分位数处，半变异函数模型拟合精度最高，土壤盐分呈中等强度的空间相关性；超过中位数阈值的RMSE相对较小的第70、80百分位数处，土壤盐分呈强烈的空间相关性，而其半变异函数的拟合精度较低。

从不同盐分阈值下的指示概率分布图可见，随着盐分阈值的增大，高概率面积逐渐减小，低概率面积逐渐增大，盐分预测平均条件概率随着阈值的增大逐渐减少。主要概率分布区间由盐分阈值0.09%的0.9～1.0过渡到盐分阈值0.35%的0.0～0.1。概率值分布的偏度由盐分阈值0.09%的负偏过渡到盐分阈值0.35%的正偏，峰度先增大在中位数阈值0.13%处达到最大，后急剧降低，在盐分阈值0.35%处又开始回升，在中位数阈值0.13%处最接近正态分布。某一盐分的高概率分布区域一般包含在小于它的上一个盐分阈值的高概率分布范围内，相应的低概率分布区域一般包含在大于它的上一个盐分阈值的低概率分布范围内。盐分阈值可根据不同阈值目标下特定区域的概率与风险进行选择。