东北黑土漫岗区春耕期土壤水分空间变异及地形影响

靖亭亭, 胡海棠, 邱春霞, 李存军, 周静平, 张巧玲, 白 翠

(1.西安科技大学, 西安 710054; 2.北京农业信息技术研究中心,北京 100097; 3.国家农业智能装备工程技术研究中心, 北京 100097)

东北黑土区是我国重要的粮食主产地，大部分处于中深度季节性冻土区，春季水热条件是东北地区粮食生产主要限制因子[1]。漫川漫岗坡耕地占黑土区耕地面积的60%，春季黑土冻融和积雪融化及降雨造成了水土流失，在春耕期容易导致农田土壤水分饱和甚至积水，农机无法进地作业，影响春耕整地和春播进程。2013年黑龙江地区由于低温春涝，作物播种较常年推迟5～15 d[1]，是当年大豆和马铃薯减产8.5%和26%的重要影响因素。由于受漫川漫岗地形的影响，土壤水分呈现显著的空间异质性。因此，研究东北黑土漫岗区土壤水分分布规律，对于黑土区土壤水分管理、防治水土流失、合理安排春耕春播，保障粮食产量具有重要的理论和现实意义。

土壤水分对于农业、水文、土壤侵蚀、土壤养分变化以及局部气候等都有重要影响[2-3]，在作物生长和农业生产上起到决定性作用[4]。土壤水分在不同时间、地点和土层深度上表现出明显的时空异质性，这是多重尺度上地形、土壤、土地利用及植被、气象(降雨)、人为活动等多因子综合作用的结果[2]。由于土壤水分在区域和全球尺度的重要性，过去20 a来，国内外学者对土壤水分进行了大量不同尺度的研究[4]。邱扬等[2]对不同研究中不同区域土壤水分主要影响因子进行了总结，主要包括相对海拔、坡向、坡度、坡形、汇流面积、湿度指数、植被类型、导水率、黏粒含量等，但是其主控因子因地因时因尺度而异。Xuemei Mei[5]对黄土高原不同时期不同土层土壤水分的主控因子分析也证实了这一点。Entin等[6]认为土壤水分的时空变异尺度可分成大尺度和小尺度2个组分，大尺度由大气控制主要决定于降雨和蒸发格局，小尺度主要决定于土壤、地形、植被等。Western A W[7]使用半方差函数研究土壤水分的空间变异，并定义变异尺度类型是由间距、范围、点大小3种尺度构成，证明半方差图是可以应用于土壤水分空间变异，且其变异尺度是一个重要性问题。

东北典型黑土区以漫岗丘陵为主，坡长坡缓，气候低温湿润，其独特的冻融侵蚀是水土流失的重要形式。近年来随着对黑土区水土流失和退化问题的重视，东北黑土区土壤水分的研究逐渐兴起，坡面尺度土壤水分异质性及主控因素的研究成为关注的热点[8-12]。葛翠萍[11]、苏子龙[12]、胡向辉[13]等主要探讨了坡位、坡向、坡度、海拔等地形因子和植被景观对土壤水分的影响。黑土区土壤水分研究主要集中在生长季或收获后，但是春耕期土壤水分状况对春耕春播农机科学调度加快春耕春播有重要指导作用，因受冻融、积雪融化及雨雪、大风的综合影响，土壤水分异质性更加复杂，且土壤水分变化剧烈，是研究的难点。一方面对于春耕期土壤水分分布特征及影响因素还尚不明确，另一方面已有研究通常基于田间离散采样获取土壤水分数据进行地形因子相关性分析，样点空间代表性和连续性不足。对于山地等土壤水分空间异质性显著的地区，基于遥感尤其是Sentinel-1等雷达数据的土壤水分反演方法[14-15]比传统空间插值方法，将能更准确反映土壤水分的空间分布特征[4]。

本研究以东北黑土区黑龙江赵光农场为研究区，基于Sentinel-1雷达数据反演的土壤水分以及利用DEM提取的坡度、坡向、坡位、高程、地形湿度指数等地形因子数据，采用半方差函数、相关性分析结合集成推进树(Aggregated boosted tree, ABT)方法，揭示东北黑土漫岗区春耕期土壤水分空间变异尺度特征及关键地形因子，旨在为黑土区坡耕地水分管理、水土流失调控和春耕春播农机科学调度提供依据。

1 材料与方法

1.1 研究区概况

赵光农场位于黑龙江省北安市中部，东经126°26′—127°6′，北纬47°54′—48°12′。赵光农场是我国第一个国营机械化农场，地处小兴安岭南麓边缘，总的地势是西南高，东北低。地处中高纬度，属寒温带大陆性季风气候，冬季寒冷干燥，夏季湿热多雨，年降雨量570 mm，全年无霜期为120 d左右。耕地面积2.44万hm2，种植作物主要为大豆和玉米，垄向以顺垄为主。地形具有典型黑土区漫川漫岗的特征，坡长坡缓，坡度多分布在3°～5°。东北地区冬季土壤上冻，春季耕层冻土开始融化适宜农机进地耕作。

1.2 数据源及预处理

1.2.1 地面采样及Sentinel-1数据处理 研究区于2018年4月19日开始春耕，地面水分样点于2018年4月24日在研究区均匀布设，每个样点使用GPS记录坐标、高程等信息，采用五点混合取样法采集土壤表层0—20 cm土壤样品，取部分鲜土放入铝盒内称重，之后烘干得到125个土壤重量含水量采样点，剔除异常值后共计120个样点。土壤水分遥感反演利用Sentinel-1卫星数据，该卫星4月24日过境，保证了地面观测数据与遥感影像的同步性。

Sentinel-1卫星搭载C波段SAR传感器，具有双极化、重访周期短等特点[16]，本文选取2018年4月24日Sentinel-1地距(ground range detected, GRD)产品，其包括垂直发射水平接收(VH)和垂直发射垂直接收(VV)两种极化方式，本试验数据在欧空局官网(https:∥scihub.copernicus)下载。利用欧空局提供的SNAP软件对Sentinel-1 SRA数据进行斑点滤波、辐射定标、几何校正和地理编码，获取影像的后向散射系数。

1.2.2 DEM数据及预处理 本研究在NASA中下载得到研究区分辨率为12.5 m的DEM数据，借助GIS软件，使用DEM数据提取高程、坡度、坡向、坡位、地形湿度指数(Topographic Wetness Index, TWI)等地形因子。其中，坡度、坡向均在Arc GIS空间分析模块中从DEM模型中直接提取，其中坡向按平地(-1°～0°)，阴坡(0°～45°，315°～360°)，半阴坡(275°～315°)，半阳坡(45°～135°)，阳坡(135°～275°)分类[17]。地形湿度指数计算见文献[18]，在SAGA GIS中提取。

坡位提取见1.3研究方法，各地形因子描述见表1。根据东北典型黑土区长坡长、缓坡度的地形特点，水土流失严重区的地表坡度一般集中于3°～5°农耕地上，个别地方可达8°以上[19]。3°是研究区农耕地的代表性坡度，5°是该区侵蚀严重农耕地的代表性坡度[20]。考虑到该区特殊的地形特征，故选择了 3°和 5°作为试验研究农耕地的坡度，坡度进行划分：缓坡(<1°)、缓中坡(1°～3°)、中坡(3°～5°)、陡坡(>5°)。

表1 地形因子及描述

1.3 研究方法

1.3.1 土壤水分遥感反演 黑土春耕期均为裸土，雷达后向散射系数不受植被影响，用其建立土壤水分模型。建模点与验证点比例为6∶1，以决定系数(R2)和均方根误差(RMSE)作为拟合模型的精度评价标准。基于120个土壤采样点水分值与VH、VV后向散射系数的相关性分析可得，土壤水分与VV相关性为-0.427，与VH相关性为-0.008，VV极化后向散射系数与土壤水分相关性较高且通过ρ<0.01的显著水平，VV极化后向散射系数对于土壤水分变化更敏感。使用VV极化建立反演模型(表2)，经试验表明指数模型R2最大RMSE最小，拟合效果相对较好，最终采用4折交叉验证法[21]验证所建立的模型。经交叉验证后的反演模型平均R2为0.576 8，比拟合模型增大了0.002 6；RMSE为1.752 4，比拟合模型减小了0.020 9，选择的模型可以使用此研究。

1.3.2 土壤水分空间变异分析 半变异函数γ(h)(公式1)是通过计算两点空间属性的方差与两点间距的关系来描述要素的分布特征[22]。半变异函数主要有以下参数：块金值(C0)、基台值(C0+C)和变程(A)。块金值主要反映随机因素对空间变异的影响；基台值主要反映结构性因素对空间变异的影响；变程主要反映影响空间变异的尺度大小，变程范围内存在空间相关性，以外的距离不再有相关性[23]。土壤水分半方差函数计算利用ENVI中统计功能完成，并用Matlab实现半方差图构建。本研究基于土壤水分遥感反演结果，分别以整个赵光农场研究区以及典型地块为对象，主要利用半变异函数(semi-variogram)分析土壤水分农场尺度和地块尺度的空间变异特征。

(1)

式中:γ(h)为半方差函数;Z(Xi)为像元所对应的土壤属性值；Z(Xi+h)则为步长为h的另一像元所对应的土壤属性值;h为步长；N(h)为步长为h的像元对数。

1.3.3 组合尺度的坡位因子计算 坡位采用Weiss A[24]提出的方法，利用地形位置指数(Topographic Position Index,TPI)和坡度结合划分某点不同坡面的位置。TPI反映的是某点高程与其一定邻域内高程的差值(公式2)，其中邻域阈值不同反映了不同尺度上的坡位特征，含有尺度效应。TPI利用Arc GIS扩展工具分析套件Land Facet Corridor Designer中的Topographic TPI工具进行计算，其邻域阈值可以根据研究区土壤水分的空间变异特征来确定。大尺度上的谷地，也可能是局部的漫岗。将大、小尺度上的坡位进行组合生成组合尺度坡位，能够反映在全局和局域上的坡位特征。

(2)

式中:Z0表示中心点高程；R表示预设的邻域；Zi表示邻域内高程；n表示邻域内高程点数。

根据土壤水分的空间变异(变程)特征分析，分别设定不同邻域阈值计算不同尺度地形位置指数(TPI)，并将两种尺度坡位进行组合，从而表征每个像元在区域(农场)和局域(地块)的坡位，邻域大小根据土壤水分空间变异特征设定。如组合坡位值为1f，代表在农场尺度处于谷地，而在地块内部位于漫岗区，具体阈值情况见表3—4。

表3 大尺度(农场)坡位分类阈值

1.3.4 影响因子分析方法 为探明土壤水分在农场尺度和地块尺度的空间异质性形成机制，利用Arc GIS在研究区740个地块上随机生成6 000点，提取各点对应的土壤水分和地形因子(坡度、坡向、高程、地形湿度指数、农场尺度坡位、地块尺度坡位、组合尺度坡位)的值。利用SPSS软件对数据结果进行双变量Spearman相关性分析，基于R语言gbm软件包进行集成推进树(Aggregated boosted tree, ABT)分析，根据重要性排序评价坡度、坡向、高程、坡位、地形湿度指数等地形因子对土壤水分的相对影响[25]，方差分析结果采用 Origin绘图表示。

表4 小尺度(地块)坡位分类阈值

集成推进树算法(Aggregated boosted tree, ABT)是基于 PAC模型(Probably Boosting Stacking, PAC)，通过相对重要性图来定量评价和描述自变量对因变量的相对影响大小，重要性排序越靠前，其对土壤属性的推测能力越强。相应地，重要性排序越靠后的变量对土壤属性的推测能力越弱。

2 结果与分析

2.1 土壤水分分布及空间变异特征

研究区土壤水分空间分布见图1，为使该数据与DEM像元大小保持一致将土壤水分反演数据重采样到12.5 m。图中可以看出研究区土壤含水量在25%～37%，均值为30.65%，属于湿润级别[26]。土壤水分空间分布呈现东南部高，西部相对较低的格局。基于土壤水分反演数据，计算在不同尺度变异系数(表5)，地块尺度上以地块为单位计算了农场上每个地块的变异程度，即地块内部变异为5.81%(取所有地块变异系数的均值)；农场尺度上每个地块的均值为基础值统计地块之间的变异为4.16%，表明黑土漫岗区土壤水分地块内部的变异大于地块之间的变异。单样本K-S检验(α=0.05)结果显示，土壤水分在各个尺度上均服从正态分布。

表5 不同尺度土壤水分统计特征

图1 土壤水分空间分布

通过式(1)计算土壤水分半方差函数，半方差图采用y轴表示半方差值，x轴表示不同的步长。土壤水分在农场尺度及地块尺度上变异特征曲线相似，但经验半方差值达到了稳定的常数(基台值)的像元个数不同。在农场尺度上，在250个像元前达到稳定值，即有效变程为3 000 m；地块尺度上，在25个像元前基本达到稳定值，即有效变程为300 m。在本研究中，农场尺度上主要呈现了3 000 m左右范围内变异特征，地块尺度上主要反映了300 m左右的变异特征。该变程范围为上一节中坡位的划分邻域提供参考依据。

2.2 土壤水分与地形因子的相关性及主控因素识别

Spearman相关性分析显示，土壤水分与地形因子之间均呈显著相关(p<0.01)，其中与地形湿度指数呈正相关，与高程、坡度、坡向、坡位呈负相关。表明地势越平缓，地表径流量越大，地块的土壤水分含量越高，其中在不同尺度上的坡位均影响水分的分布(表6)。

表6 土壤水分与地形因子Spearman相关性分析

集成推进树算法(Aggregated boosted tree, ABT)显示了各个地形因子对土壤水分分异性的相对重要性(图2)。坡位是影响土壤水分分的首要因素，其相对影响大小为36.28%，其次是坡度，其相对影响大小为26.54%，坡向15.57%，地形湿度指数为12.72%，高程为8.89%。

图2 地形因子对土壤水分分异的相对影响

2.3 土壤水分主要影响因子

根据ABT分析的结果，坡位、坡度和坡向是影响土壤水分分布的关键地形因子，根据研究区地形的划分对在不同坡位、不同坡度以及不同坡向条件下的土壤水分分布情况做定量研究，进一步探索土壤水分的空间分布特征。

2.3.1 坡位 图4展示了在不同尺度上土壤水分含量随坡位的变化情况。五组柱状图分别以不同填充代表农场尺度和地块尺度坡位的影响，横轴为农场尺度的坡位代码对应表3，虚线是农场尺度各个坡位上土壤水分的均值，作为组间的参照值。农场尺度上，谷地、下坡位、中坡位、上坡位和漫岗平均含水量分别为33.23%，33.07%，32.76%，32.47%和32.71%。可以看出，除第五组即漫岗地块外，土壤水分随农场尺度和地块尺度坡位的提升都呈现梯度下降趋势，地块尺度坡位的组内差异显著大于农场尺度坡位之间的差异。

2.3.2 坡度 在不同坡度下(图4A)，土壤水分含量随坡度的升高呈梯度下降趋势，缓坡、缓中坡、中坡、陡坡平均含水量为32.99%，32.71%，32.05%和31.53%。缓坡>缓中坡>中坡>陡坡，即坡度越大，土壤含水量越低，陡坡和缓坡平均含水量相差1.46%。

注：横轴代表农场尺度坡位代码用不同填充显示，类别与表3对应；每组相同填充的6个柱行表示地块尺度坡位，类别与表4对应，依次对应编号类型a—f。

2.3.3 坡向 不同坡向上，土壤水分含量有较大差异，平地、阴坡、半阴坡、半阳坡和阳坡平均含水量分别为32.92%，32.70%，32.48%，32.92%和32.37%，半阳坡>平地>阴坡>半阴坡>阳坡(图4B)。半阳坡与平地相差不大，属于水分最大的区域，阳坡处最小两者相差0.56%。在不同坡向上，土壤平均含水量差异较小，表明了这一时期坡向相比坡度来说对于土壤水分分布影响较小。

注：不同字母表示差异显著，相同字母表示差异不显著。

3 讨 论

土壤水分空间分布表明，东南部区域水分含量相对较高，西部区域较低，且空间变异在不同尺度上存在差异。通过每个尺度上土壤水分变异系数发现，每个地块内部以及相邻地块之间土壤水分都存在空间分异。因此推断，土壤水分可能具有地块尺度和农场尺度两个尺度上的空间变异特征，可能原因是受到地形、田间田埂沟渠结构影响。赵光农场是大机械化生产，每个地块面积通常约30万m2，地块之间常有田埂、机耕道、沟渠、防风林带分布。因此，一方面大面积地块，受地形影响地块内部存在明显空间差异；另一方面地表径流等土壤水分运移过程侧重于以地块为相对封闭单元发生。对整个区域而言，又存在着集水区尺度的水分运移。整个农场的土壤水分横向迁移表现为以顺垄地块的地块尺度水分运移为主，以地块间水分迁移为辅。通过半方差函数计算发现，土壤水分在整个农场上空间变异有效变程为3 000 m，地块上空间变异有效变程为300 m。

相关性分析结果表明地形因子(坡度、坡位、坡向、高程、地形湿度指数)和土壤水分含量呈显著相关(p<0.01)。土壤水分与坡度、坡向、坡位、高程呈负相关关系，与地形湿度指数呈正相关关系。本研究区春耕期处于土壤水分冻融期，土壤水分主要来源于大气降水、冻融和积雪融化。高程影响降雨的再分配以及径流量大小，随海拔的升高土壤水分含量降低[27]。曹伟等[28]通过研究坡面冻土表明，受重力作用坡度越大，土壤径流越大，壤中流减少，土壤水分含量减少，与本研究结果相似。坡向影响接收太阳辐射大小，从而造成不同程度的水分蒸发[29]。不同坡位体现了地表径流及汇流累积过程造成土壤水分的空间差异[11]。邓慧平等[30]也发现土壤相对含水量(干湿情况)与地形湿度指数呈线性关系，且两者的空间变化呈正比。综上，土壤含水量随坡度的增大、太阳光辐射的增强、坡位的升高而减小，随地形湿度指数的增大而增大。

基于ABT算法的因子相对影响程度分析表明，坡位的相对影响最大，坡度、坡向次之，地形湿度指数和高程最小，与朱华德等[31]研究土壤水分湿润时期结果一致，但与张继光等[32]、葛翠萍等[11]、郭欣欣等[9]研究结论有所差异。土壤水分的主控因素是随时间和空间变化的[5，33]。例如：葛翠萍研究作物生长期表明坡度是主要影响因子，本研究中坡位是该时期土壤水分异质性的主要影响因子。可能原因一是春耕时期为裸土，地块上土壤水分的空间分布没有植被的干扰，坡位成为主导因素[34]；二是春耕时期研究区三天前发生过一次雨雪过程，叠加冻融过程，共同形成的地表径流造成土壤水分空间上的重分配且坡位划分时考虑了坡度的变化且分类较细致。研究中使用两种尺度的坡位进行组合分析对土壤水分的影响，且在各个尺度表现出明显规律性。农场尺度上，土壤水分含量沿坡长方向呈先降低后升高的趋势，谷地到达水分最大值。这种趋势规律与黑土区前人研究成果相似[11-12，35]。研究中发现漫岗(坡顶)土壤水分含量较大，一方面因为坡顶漫岗处坡度较缓、径流作用不强；另一方面因为研究区漫岗处的地块主要以横垄垄作方式。相关研究表明横坡垄比顺坡垄更能减缓坡面径流的下切作用，有效拦截水分进而增加土壤含水量[9，36]，表6也反映出垄向对土壤水分分布有显著影响。而上坡位、中坡位坡度较大，排水通畅水分在此处留滞较少，谷地(坡脚)处则有更多的水分累积。地块尺度上，土壤水分呈现随坡位下降而升高的显著规律性。比较农场尺度和地块尺度不同坡位的土壤水分差异发现，地块尺度坡位引起的土壤水分差异大于农场尺度坡位，说明农田土壤水分的空间重分配过程主要发生在地块内部。研究区地块边界有田埂，机耕道、沟渠、防风林带等分布，形成了地块之间土壤水分横向迁移的屏障，研究土壤水分的空间分布不能将整个区域作为一个连续变化的地表来处理。虽然有田埂等屏障，地块也不是完全封闭的单元，除了壤中流，少数地块地表径流量较大，可能会溢出田埂或沿着侵蚀冲沟流动，形成地块间的水分交换。研究区在4月19～20日发生的雨雪过程，以降雪为主。冻土化冻和融雪过程都是缓慢而持续的过程，大多数地块形成的地表径流很小，未能形成跨越田埂的地表径流。

不同坡度上土壤水分含量大小为，缓坡>缓中坡>中坡>陡坡，随坡度增加水分减小。坡度主要影响降雨的入渗速率[37]，研究区耕地坡度较缓，不同于喀斯特[32]、黄土地区[38]坡度较陡在地形要素中对径流影响最强。坡向在各个方向上接收太阳辐射强度和日照时间不同，表现为阳坡>半阳坡>半阴坡>阴坡[39]。不同坡向土壤水分的蒸发量不同，本研究发现土壤水分含量在坡向上的趋势为半阳坡>阴坡>半阴坡>阳坡，与以往经验不同，可能原因是不同坡向在各个方向上接收太阳辐射强度和日照时间不同，造成土壤水分的蒸发量不同。而研究区当时气温低，春季主要风向为西北风，当日及前两天风向为西北风3～5级，阴坡和半阴坡属于迎风坡，水分蒸发较快，而半阳坡属于背风坡，风干较慢，土壤含水量相对较高。地形湿度指数对研究区土壤水分的贡献度偏低，它的提出是基于一个比较理想的连续变化的流域表面，根据汇流面积与坡度计算，且计算前将地形数据进行了填洼处理忽视了微小地形的影响[40]。本文研究区地形起伏较大且地块中有冲沟分布，壤中流的速度不稳定，造成TWI对土壤水分的空间分布解释较低，这与Ladson A[41]的研究结果相似。王洪明[42]研究表明地形湿度指数对表层土壤水分分布的影响较小。因此，尽管地形湿度指数是土壤水分空间分布研究中被广泛应用的指标，但是它对本文中土壤水分的解释并不太有效。本研究提出的结合农场尺度和地块尺度的坡位，更好地反映了东北漫岗区春耕期土壤水分的运移和分布规律。

4 结 论

(1)赵光农场于2018年4月24日春耕期土壤水分含量分布在25%～37%，均值为30.65%，整体属于偏湿润。东南部区域含水量相对较高，西部区域较低。农场地块内部土壤水分变异系数(5.81%)大于地块之间变异系数(4.16%)。土壤水分表现出在地块尺度和农场尺度两个尺度上的空间变异特征，土壤水分在农场上空间变异的变程为3 000 m，地块上空间变异的变程为300 m。

(2)土壤水分含量与地形湿度指数呈显著正相关关系，与坡度、坡向、坡位、高程呈显著负相关关系。

(3)利用集成推进树方法探测地形因子对土壤水分的相对影响表明，坡位、坡度、坡向是影响春耕时期土壤水分分异的主要地形因子，累计相对解释度超过了70%，其中坡位的相对解释度达到36.28%，是主要的影响因子。