灌区土壤盐分空间变异及多因素响应关系

2020-06-15谷晓伟韩金旭张凤燃刘姝芳

水土保持通报 2020年2期

谷晓伟, 韩金旭, 张凤燃, 张丹, 刘姝芳

(黄河水利委员会黄河水利科学研究院, 河南郑州 450003)

土壤盐碱化是灌溉农业的伴生问题，据FAO统计全球约1/3的灌溉土地因盐碱化而生产力降低。中国盐渍化土壤面积约3.60×107hm2，占全国可利用土地的4.88%，其中耕地盐渍化面积9.21×106hm2，占全国耕地面积6.62%[1]。土壤受母质、地形等自然及人为因素共同作用，具有高度空间异质性[2]。了解土壤空间特性是盐渍化防治的重要前提。刘鑫等[3]基于ArcGIS对河套灌区土壤盐碱化指标进行了等级划分；徐英等[4]研究发现地下水对地表返盐的影响具有滞后性，河套灌区4月底地表发生中、轻度盐渍化的地下水临界埋深分别为2.0 m，2.5 m。王卓然[5]、吴亚坤[6]、徐存东[7]等基于多源数据，分析了不同区域土壤盐分空间变异特征。从目前成果来看，土壤盐分空间变异研究较多关注盐分单一要素的空间异质性分析，对多参数空间协同关系研究尚待进一步加深。黄河南岸灌区是鄂尔多斯节水改造的重点区域，也是受盐渍化影响的典型灌溉绿洲灌区。本文以吉格斯太灌域为例，根据实测数据，利用地统计学方法解析灌域土壤盐分空间分布特征，研究土壤盐分空间格局与地下水、土壤物理特性之间空间响应关系，成果可为灌区土壤盐渍化防控提供理论依据。

1 试验设计、材料与方法

1.1 试验区概况

试验区位于鄂尔多斯市北部黄河南岸灌区吉格斯太灌域(110°33′29″—110°42′30″E，40°18′25″—40°19′8″N)扬黄节水改造区，面积约2 500 hm2。灌域东以呼斯太河为界，西以东柳沟为界，北至黄河，南至毛乌素沙漠，为井渠结合灌域，总灌溉面积6 800 hm2。区域多年平均气温6.8 ℃，平均无霜期160 d，平均日照时数3 159 h，平均降水量311 mm，平均蒸发量2 115 mm。

1.2 样点布设及样品采集

样品采集于2016年5月中旬进行，试验区按照网格化结合土地利用类型布置采样点29个，地下水观测井11眼(如图1所示)。按照0—20 cm(以下称表层)、20—40 cm和40—60 cm 3个深度人工取土，采用环刀法测定土壤干容重及含水量，同步观测地下水埋深和含盐量。风干土样过2 mm筛，表层土壤测定颗粒组成，其余按5∶1配置土壤溶液测定电导率，并换算含盐量。

图1 研究区样点分布

根据监测情况，试验区地下水埋深在1.13～3.27 m，含盐量为1.12～2.86 g/L，地下水埋深自西北向东南逐渐加大，含盐量逐渐降低。试验区表层土壤容重1.40～1.79 g/cm3，砂粒含量占比86.58%～96.68%，粉粒含量占比2.65%～13.17%，黏粒占比0.16%～3.03%。按照国际土壤质地分类标准，29个采样点全部为砂土(表1)。

表1 表层(0-20 cm)土壤物理特性参数统计表

注：黏粒<0.002 mm,粉粒0.002～0.02 mm,砂粒0.02～2 mm。

1.3 数据分析方法

本次采用地统计学方法分析土壤含盐量空间分布特征，采用互相关函数法研究土壤盐分与地下水及表层土壤物理特性间的空间响应关系。其中，土壤物理特性参数选择黏粒、粉粒、砂粒含量、干容重、体积含水量、导热率与热容量等指标，热容量、导热率采用相关文献[8-9]推荐公式计算。地下水参数包括埋深和含盐量。主要统计分析方法如下：

(1) 经典统计学方法。

(1)

(2) 地统计学方法。

(2)

式中：γ(h)为空间变量半方差函数值，用以表示变量空间变异结构；h为两样本点空间距离；Z(xi)，Z(xi+h)分别为变量在点xi和xi+h处观测值；N(h)表示滞后距离为h时的样本对数。本文中半方差函数主要采用如下模型：

高斯模型：γ(h)=C0+C〔1-exp(-h/a)2〕

(3)

指数模型：γ(h)=C0+C(1-e-h/a)

(4)

式中：C0为块金值；a为变程(km)；C为拱高；C0+C为基台值。

区域化变量空间变异由结构性因素和随机性因素共同引起，气候、地形、母质等结构性因素使变量具有空间连续性，耕作、灌溉等随机因素则会弱化其连续性[4]。一般用块基比〔C0/(C0+C)〕表示随机变异占总变异的大小，块基比0～25%属较强空间相关性，25%～75%属中等空间相关性，75%～100%为弱空间相关性。

(3) 空间自相关性。

(5)

(6)

(7)

式中：Z(I)表示I的标准正态分布值，当Z>1.96或Z<-1.96(α=0.05)时，表明变量在空间上具有显著自相关性；Z位于[-1.96,1.96]则观测值呈独立随机分布；E(I)为理论上的数学期望； var(i)为方差。

(4) 互相关函数法。

互相关函数可以较好地表达两变量间空间相关关系[9-12]，计算公式为：

(8)

2 结果与分析

2.1 土壤盐分统计特征

计算土壤盐分统计参数如表2所示。参考相关成果[3-4,13]，含盐量<2 g/kg为非盐化土；2～4 g/kg为轻度盐化土；4～6 g/kg为中度盐化土；6～10 g/kg为重度盐化土。由表2可知，各土层平均含盐量1.96～2.46 g/kg，属非—轻度盐化土。表层土壤平均含盐量最高，变幅最大，极差达到8.29 g/kg。各层含盐量Cv值为0.52～0.61，属中等变异强度。

表2 土壤含盐量统计特征

2.2 土壤盐分空间分布与变异特征

2.2.1 土壤盐分空间变异性根据表2偏度、峰度及K-S检验结果，土壤含盐量近似服从对数正态分布，满足地统计分析要求，利用GS+9.0对土壤含盐量进行半方差分析，采用GeoDa 1.14计算Moran’sI指数[14-15]，结果详见表3。

根据分析结果，土壤含盐量可采用高斯模型和指数模型较好拟合，R2在0.74～0.82。根据表3各土层含盐量块基比在25%～75%，结合公式(2)—(4)可知土壤含盐量空间分布受到的人为、自然因素的共同作用，属于中等空间变异程度。不同土层Moran’sI值为-0.088～0.044，Z得分为-0.488～0.744，位于[-1.96,1.96]，土壤盐分在空间上呈相对独立的随机分布。

表3 土壤盐分地统计分析结果

2.2.2 土壤盐分空间分布特征利用Surfer 12.0基于Kringing插值计算土壤含盐量空间分布如图2所示。

由图2可知，土壤含盐量呈斑块变化，总体自西北向东南逐渐减小。除中东部存在非盐渍化区域外，其余均达到或超过轻度盐化阈值，西北局部达到重度盐化程度。随土层加深，盐分峰值总体有所降低，东南部非盐化土范围不断扩大，西部重度盐化土面积减小。按表层土壤划分，试验区非盐渍化面积(<2 g/kg)占25.6%，轻度盐渍化面积(2～4 g/kg)占64.4%，中度以上盐渍化面积占10%(>4 g/kg)，总体处于非—轻度盐渍化状态。

2.3 土壤含盐量多因素响应关系

2.3.1 土壤含盐量多因素相关性分析根据半方差分析结果，试验区土壤盐分受结构性和随机性因素共同作用。地下水是耕层土壤的下边界，表层土壤物理特性除影响本层外，也作为上边界影响下层土壤，两者均可视为影响土壤盐分的结构性因素。为便于分析，将土壤划分为表层(0—20 cm)和中下层(20—60 cm)，中下层土壤含盐量采用20 cm以下含盐量均值。表层土壤物理特性及地下水参数与含盐量相关性分析结果如表4所示。由于地下水井与采样点无法一一对应，本次地下水数据采用空间差值结果。由表4可知，表层土壤含盐量与其自身物理特性之间相关性不显著，中下层土壤含盐量与表层土壤含水量显著正相关。地下水对土壤含盐量影响明显，土壤含盐量与地下水埋深显著负相关，与地下水含盐量显著正相关。

表4 不同土层土壤含盐量与地下水、表层土壤物理特性相关性

注：*表示在p>95%水平上相关性显著；地下水埋深及含盐量数据为空间差值计算结果。

2.3.2 土壤含盐量多因素空间相关性分析双因子相关性分析(表4)表明土壤含盐量受地下水影响显著，但除含水量外，土壤含盐量与表层土壤物理特性无显著相关性。作为空间连续体，土壤盐分与其相关因素之间还可能存在空间依赖性。采用互相关函数[9-11]分析土壤含盐量与地下水、表层土壤物理特性间的空间相关性如图3—4所示。表层土壤含盐量与地下水及其自身物理参数空间相关性分析结果如图3所示。由图3可知，表层土壤盐分与其自身黏粒含量、容重、水热参数在2～6 km范围内显著正相关，与砂粒含量在2～4 km范围内显著负相关，表明在相关范围内土壤含盐量随黏粒含量、容重、含水量、热容量及导热率的增加而增大，随砂粒含量增加而减小。表层土壤盐分与地下水含盐量显著正相关，与地下水埋深显著负相关，地下水对土壤盐分影响较大。

对比表4和图3可知，土壤盐分与物理特性和地下水之间存在空间响应关系，物理特性对盐分的影响具有空间滞后性。

注：图中虚线为95%置信线，相关函数值超过置信线即为相关性达到显著。下同。

由图4可见，中下层土壤盐分与地下水及表层土壤物理特性存在较显著的空间相关关系。其中，中下层土壤盐分与地下水参数、表层土壤黏粒、砂粒含量、热容量及导热率的空间相关范围和属性均与表层土壤基本相当，与表层土壤含水量正相关范围扩大至0～6 km，容重和粉粒含量则未达显著级别。由此可见，表层土壤物理性质与其下层土壤盐分存在空间协同性，表现出一定的依赖关系。

图4 中下层土壤含盐量与地下水、表层土壤物理参数空间相关性

2.3.3 土壤含盐量与土壤物理特性及地下水空间分布关系绘制地下水埋深、表层土壤物理特性及不同土层盐渍化风险空间分布如图5所示。由图5可见，土壤砂粒含量、容重及水热特性呈较一致的斑块分布。在试验区西北—北、西南和东南部表层土壤各物理特性存在明显的极值斑块，其中北部轻度盐渍化高风险区与土壤物理参数空间分布一致性较明显，但在东南部上述因素与盐渍化风险分布一致性减弱，这可能与东南部地下水埋深明显加大有关。总体而言，从斑块分布来看土壤物理参数和含盐量具有较明显的响应。研究区地下水埋深h由北向南逐渐加大，北、西北部受黄河侧向补给h<1.6 m，南部h逐渐加大至3 m以上。地下水含盐量呈与埋深相类似的带状分布，自西北向东南逐渐降低。结合土壤盐渍化概率分布，地下水埋深较浅(<1.6 m)、含盐量较高(>2.4 g/L)的灌域西、北部轻度盐渍化风险较高，地下水与土壤盐分间存在较明显的空间响应关系，与前述分析结果一致。

3 讨论与结论

3.1 讨论

土壤盐分变化受到立地属性、耕作条件等自然和人为因素的共同影响。在中等尺度上研究土壤盐分空间变化，分析其与土壤物理组成、水热效应和地下水间的关系，有助于在较大尺度范围内认识盐渍化发生的空间规律，为区域控盐提供科学依据。

图5 表层土壤物理特性、地下水埋深及轻度盐渍化风险空间分布

(1) 土壤特性的空间属性可以采用地统计学方法描述。王全九等[9,14,16-17]在不同地区开展的研究均发现土壤含盐量在一定范围内具有空间结构特征，半方差函数可以用高斯模型和指数模型拟合，与本次研究结论基本一致。结合前人成果可见，不同区域变量空间自相关性存在差异，王维维[14]等发现焉耆盆耕作层土壤盐分Moran’s I接近于0，即呈相对独立的随机分布，与本次研究结论一致；张飞[15]等则发现精河绿洲0—20 cm土壤电导率表现出较强的空间正相关性。

(2) 潜水蒸发是耕作层土壤盐分的主要来源。根据河套灌区开展的相关研究，窦旭等[2]认为埋深大于1.6 m时土壤含盐量随埋深增加变幅逐渐减小；徐英[4]等研究发现灌区4月底发生盐渍化的临界地下水埋深为2.0～2.5 m。本研究表明，地下水参数与土壤含盐量存在显著的相关关系，地下水埋深<1.6 m的区域发生轻度以上盐渍化概率较高，与前述相关研究结论基本一致。

(3) 土壤物理性质与盐分关系尚无统一认识。吕廷波[18]、毛海涛[19]发现细颗粒含量与土壤最终积盐量成正比；刘继龙[16]认为0—20 cm土壤电导率与基本物理特性相关性不显著，但20—40 cm土层两者显著相关。本次通过相关性分析(表4)认为，地下水与土壤盐分显著相关，与土壤物理特性相关性不显著。但相关性分析只能代表各因素在样本点上的关系，对空间响应关系无法表达，本次采用互相关函数分析发现土壤盐分与土壤物理特性和地下水等多因素之间存在空间相关性，更能体现区域变量在空间上的协同性。目前对土壤盐分与其他因素空间相关性研究尚不多，王全九[9]研究粉壤土灌区得出土壤盐分空间分布与含水量、热容量及导热率等水热条件的空间关系不显著，与本次研究结论不一致，这一差异是否与土壤质地条件相关尚需开展进一步的比较研究。本次研究发现表层土壤物理特性对中下层土壤盐分空间分布存在影响，这一影响在纵向上可以理解为表层土壤作为上边界条件对下层土壤的作用，樊会敏[20]在研究中也发现土壤物理性质和盐分在上下层之间存在异位影响，但没有从空间分布角度对异位影响进行探讨。