黄土高原是世界上水土流失最严重的地区之一^[1], 是黄河泥沙的主要来源区^[2—3], 其中浅沟侵蚀是土壤侵蚀的重要类型^[4—5], 其侵蚀量通常占坡面总侵蚀量的26.6%—59.2%^[4]。退耕还林(草)工程实施二十多年来黄土高原土壤侵蚀得到有效控制^[6—7], 但部分地区浅沟侵蚀问题仍需重视。从区域尺度阐明目前浅沟空间分布特点和防治重点是进一步优化黄土高原水土保持措施布设的重要基础。

浅沟是指能被普通耕作工具跨过和填充但不能被其完全消除的侵蚀沟^[8—9], 也被称作临时性切沟(ephemeral gully), 其宽度和深度介于细沟与切沟之间^[10], 一般深入犁底层(20 cm), 宽约30—50cm^[9]。已有研究表明, 浅沟的发生与耕作有关^{[8—9, 11—13]}。在黄土高原, 浅沟发生一般伴随着瓦背状地形的形成, 浅沟槽形部位的横断面因再侵蚀和再耕作呈弧形扩展, 形成瓦背状地形^[9], 耕作过后, 浅沟地势看似与周围齐平, 但是在下一次侵蚀事件中, 浅沟还会在瓦背状的槽底部相同位置再次发育 ^{[8—9, 12]}。因浅沟发生的临时性, 一般浅沟调查多针对瓦背状地形底部开展^[14—15]。已有报道对浅沟分布特征^{[14, 16—18]}、浅沟发育过程^{[4, 15, 19—20]}、浅沟侵蚀发生的临界地形条件^[21—26]等问题进行了研究, 其中对黄土高原浅沟侵蚀分布特征的研究多集中于流域尺度, 区域尺度的研究仍较为缺乏。同时, 对于退耕二十年后浅沟目前分布的土地利用, 尤其是仍分布在耕地上的浅沟空间分布认识也亟需加强。

近十几年来高分辨率遥感影像的发展为浅沟的区域尺度调查提供了新的思路和数据基础^[27]。本研究在黄土高原选取137个小流域调查单元, 基于亚米级Google Earth遥感影像并综合运用GIS方法对浅沟以及各调查单元土地利用进行人工目视解译和空间分析, 在区域尺度调查和分析黄土高原浅沟分布特征及浅沟所在坡面土地利用分布, 并明确目前仍处在耕地上的浅沟空间分布格局, 以期为黄土高原侵蚀治理和黄河流域生态保护与高质量发展提供科学依据。

1 材料与方法 1.1 研究区概况

本研究在黄土高原开展, 黄土高原范围参考陈永宗1988年划定的范围^[28](图 1), 即从小浪底沿沁河西分水岭, 接汾河东分水岭, 向北连清水河分水岭, 向西经托克托过黄河, 由内蒙古准格尔旗东部转向西南沿长城至中宁, 再沿黄河至兰州北, 经乌鞘岭至日月山东坡, 由贵德转洮河分水岭, 然后向东沿秦岭北坡直至伊、洛河的北分水岭于小浪底闭合, 总面积约38万km², 属于季风气候区, 多年平均气温在3.6—14.3℃, 年降雨量在700 mm左右, 降雨年际变化较大, 降雨的区域分布不均匀。

采用系统抽样方法, 纬度和经度均取0.5°为间隔, 选取了137个调查单元。基于1弧秒分辨率SRTM (Shuttle Radar Topography Mission) 数据计算调查单元平均坡度, 如大于2°, 地形起伏相对较大, 调查单元为面积0.3 km²左右的小流域；否则, 地形相对平坦, 调查单元为0.5 km×0.5 km左右的矩形范围(图 1)。其中63个调查单元分布在黄土高塬沟壑区, 41个分布在黄土丘陵沟壑区, 24个分布在河谷平原区, 6个分布在土石山区, 3个位于其它区域。同时选取无定河中游典型小流域(图 1, 面积0.44 km²)作为野外实测样区, 对45条浅沟进行野外实测, 进行浅沟解译精度验证。

1.2 基础数据

基础数据名称、来源、分辨率、获取时间及主要用途见表 1。

1.3 浅沟解译方法与质量控制

浅沟调查于2020年5月开始, 至9月底完成。基于浅沟的定义、发育过程、影像特征等进行分析及专家论证, 制定浅沟解译标准, 由解译人员分组进行解译；解译完成后进行两层质检: 由解译小组负责人进行100%复检, 发现问题后由解译人员进行全面修改；由侵蚀沟研究领域专家, 对通过复检的解译数据进行20%抽检, 正确率在95%以上为合格, 抽检合格后入库, 如不合格则重复进行第一步质检。

浅沟解译方法: 在ArcGIS 10.2软件环境下对每个调查单元浅沟长度、沟头汇水路径、条数等参数进行调查。影像时间主要选在春秋季, 以尽量减少植被、冰雪覆盖等的影响, 对于每个调查单元均利用多期Google Earth影像数据, 判断影像的可辨识程度, 保留有最大浅沟数量的影像。解译标准为: 浅沟宽约30—50 cm^[9], 一般呈现为按一定规律排列的浅洼地^[11]；在影像上, 浅沟较明显的亮或暗于周围地物, 形状较直；横断面呈弧形扩展, 呈现瓦背状地形, 无明显沟缘；以瓦背顶部分水岭为分界, 则其分水岭之间间隔一般可达十几米到几十米^[9]; 沿着同一汇水路径但中间断开几米的浅沟原则上被识别为同一条浅沟; 每条浅沟均从其分水岭位置数字化至浅沟末尾, 沟头以上部分为沟头汇水路径(图 2), 其长度为浅沟发育的临界坡长^[11]。

流域边界获取方法: 基于1弧秒分辨率SRTM数据, 利用ArcGIS水文分析模块进行大于0.1 km²的小流域边界粗提取, 以此为参考, 基于亚米级影像进行人工修订。

浅沟所在坡面土地利用解译: 基于Google Earth影像, 对每条浅沟发育所在坡面土地利用进行人工目视解译, 按1:10000解译精度要求进行, 并基于野外考察对解译结果进行了修订^[30]。

1.4 解译精度验证

采用GNSS RTK对所选典型小流域内的所有浅沟进行实地测量, 按平均0.5 m的间距测量浅沟各点的三维坐标, 同时在地形突变处适当增加测量点, 得到浅沟实测长度。基于Google Earth影像解译浅沟, 并与实测数据作对比。利用SPSS软件进行配对样本T检验, 判断实测浅沟长度、解译浅沟长度是否有显著性差异。为评估Google Earth影像解译浅沟的精度, 计算单条浅沟长度相对误差(式1)、浅沟长度相对误差的平均值(式2)、浅沟密度相对误差(式3)。

(1)

(2)

(3)

式中, R_Li为单条浅沟长度相对误差, R_L为浅沟长度相对误差的平均值, R_ρ为浅沟密度相对误差, L_GEi、L_RTKi分别为基于Google Earth影像解译和基于实测得到的单条浅沟长度(m)。ρ_GE、ρ_RTK分别为解译和实测得到的浅沟密度(km/km²)。

1.5 空间统计与分析方法

调查单元基本特征统计: 计算每个调查单元的浅沟密度(单位面积沟长度, km/km²)、浅沟沟数密度(单位面积沟条数, 条/km²)^[31]、浅沟平均长度(m)。参考土壤侵蚀分类分级标准^[32], 进行浅沟侵蚀强度分级。

空间格局统计方法: 利用ArcGIS 10.2对浅沟密度、浅沟平均长度、临界坡长统计结果进行制图, 得到黄土高原整体浅沟密度、平均长度、临界坡长的空间分布。使用莫兰指数(Moran's I)、热点分析(Getis-Ord Gi*)分析浅沟密度的空间分布特征。莫兰指数(Moran's I)用于分析空间分布的相关性, 其值介于-1和1之间, Moran's I > 0表示数据的空间整体分布呈现正相关, 其值越大代表空间上的相似性越明显；Moran's I < 0表示数据呈现高低间隔分布的状态, 其绝对值越大代表空间差异性越大；Moran's I=0, 表示数据的空间分布呈随机性^[33]；热点分析中热点是指以特定指标的高密度簇为特征的区域, 冷点是周围环绕着同一指标的低密度簇^[34]。通过热点分析可以呈现出浅沟密度高值的集中区域。

2 结果与分析 2.1 浅沟解译精度验证

对验证小流域的45条浅沟通过野外实测, 得到浅沟总长度为1067.87 m, 平均沟长23.73 m, 密度为2.43 km/km²。将实测浅沟长度与基于Google Earth影像解译的结果进行对比(图 3), 两种方法得到的浅沟长度基本分布在1:1线附近, 单条浅沟长度相对误差R_Li介于0.07%至11.66%之间, 占总数量约58%的浅沟R_Li在5%以下, 仅有不到7%的浅沟R_Li大于10%；其均值R_L为4.99%；浅沟密度相对误差R_ρ为1.23%；在0.05的置信水平上, 两种方法所得到的浅沟长度无显著性差异(P > 0.05)。

2.2 黄土高原浅沟特征统计

在黄土高原137个调查单元中, 有46个调查单元出现浅沟, 占总调查单元的33.6%。出现浅沟的调查单元中, 解译得到浅沟587条, 表 2统计结果表明, 平均浅沟密度为3.41 km/km², 最大为21.92 km/km², 最小为0.07 km/km²；沟数密度介于2—311条/km²之间, 平均为65条/km²；变异系数均大于100%。黄土高原浅沟的长度、临界坡长的范围分别为7.71—237.15 m和0—285.29 m。浅沟长度平均值为63.31 m, 中值为51.64 m；临界坡长平均值为56.20 m, 中值为47.59 m；变异系数均大于60%。

图 4显示了浅沟侵蚀达到中度及以上(浅沟密度 > 2 km/km²)的调查单元占45.6%；约87%的调查单元沟数密度小于150条/km²；约80%的浅沟长度和临界坡长介于0—80 m；浅沟长度主要分布在20—60 m之间, 约占55%；临界坡长分布在40—60 m的最多, 约占29%。

2.3 浅沟分布的空间格局

从浅沟密度的空间分布上看(图 5), 出现浅沟的调查单元主要位于黄土高原的中北部(陕西延安市、榆林市, 山西吕梁市等)和西部(甘肃庆阳市、定西市、白银市会宁等), 由东北向西南呈现条带状分布。其中浅沟密度较大的区域主要有3个, 陕西延安市西北部吴起、志丹、榆林市靖边交界处, 延安市甘泉、安塞一带, 其中吴起、志丹及靖边一带的浅沟密度高达21.92 km/km², 甘泉附近的浅沟密度为15.82 km/km², 安塞附近的浅沟密度为15.54 km/km²；山西吕梁、陕西榆林市绥德县一带黄河两岸浅沟密度也较大, 为4—8 km/km²；甘肃定西市、白银市会宁一带以及庆阳市, 浅沟密度为5—7 km/km²。浅沟密度莫兰指数为0.37, 由热点分析可以看出, 浅沟密度的热点区主要为延安市西北部, 呈现聚集分布(图 6)。

浅沟平均长度较大的区域在延安市富县、甘泉县、安塞区一带, 浅沟平均长度可达85—105 m；铜川市宜君县、白银市会宁一带也较长, 浅沟平均长度在80—85 m之间, 其余地区浅沟一般在80 m以内, 其中无定河及晋陕黄河两岸主要在40—60 m, 宁夏东南部和甘肃兰州一带主要在20—40 m(图 7)。浅沟临界坡长反映了浅沟沟头距分水岭沿流水线的距离, 在延安西北部、晋陕黄河两岸呈聚焦性低值, 主要在20—60 m, 定西到天水沿渭河一带出现浅沟的调查单元临界坡长较长, 一般在100 m以上, 宁夏固原、陕西铜川及延安南部一带临界坡长一般也在100 m以上(图 8)。

各治理分区浅沟参数有较大差异(表 3)。河谷平原区未发现浅沟, 黄土高塬沟壑区24个调查单元出现浅沟, 约占该区总调查单元的38%, 平均浅沟密度为2.07 km/km², 属于浅沟侵蚀中度(2—3 km/km²)范畴。黄土丘陵沟壑区20个调查单元出现浅沟, 约占该区总调查单元的49%, 平均浅沟密度为5.32 km/km², 属于浅沟侵蚀极强烈(5—7 km/km²)范畴。黄土丘陵沟壑区平均浅沟密度最大, 其次是黄土高塬沟壑区；同时黄土丘陵沟壑区浅沟长度平均值最大, 为65.74 m, 其次是黄土高塬沟壑区, 为60.11 m；土石山区(此区仅发现一条浅沟)临界坡长最大, 为69.86 m, 其次是黄土高塬沟壑区, 为67.24 m。黄土丘陵沟壑区临界坡长平均值最小, 为46.58 m, 这在一定程度上反映了该区发生浅沟侵蚀所需坡长较短。

黄土高塬沟壑区、黄土丘陵沟壑区浅沟密度较大, 对其出现浅沟的调查单元进一步统计得到(图 9): 黄土高塬沟壑区浅沟密度小于2 km/km²约占70%, 未发现浅沟侵蚀剧烈(> 7 km/km²)的单元。黄土丘陵沟壑区浅沟侵蚀达到强烈及以上(> 3 km/km²)的约占55%, 浅沟侵蚀极强烈(5—7 km/km²)的调查单元最多, 约占30%。黄土高塬沟壑区沟数密度小于20条/km²的约占63%, 黄土丘陵沟壑区沟数密度主要集中在50—150条/km², 约占55%。

黄土高塬沟壑区约70%的浅沟长度介于0—80 m, 浅沟长度为20—40 m的最多, 约占37%；黄土丘陵沟壑区约80%的浅沟长度介于0—80 m, 浅沟长度为40—60 m的最多, 约占28%。黄土丘陵沟壑区长度大于140 m的浅沟均多于黄土高塬沟壑区。黄土高塬沟壑区约70%的浅沟临界坡长介于0—80 m, 临界坡长为40—60 m的最多, 约占27%；黄土丘陵沟壑区90%以上的浅沟临界坡长介于0—80 m, 临界坡长为40—60 m的最多, 约占32%, 此区未发现临界坡长为140 m以上的浅沟。

2.4 浅沟所在坡面土地利用分析

通过统计每条浅沟所在坡面土地利用可以看出, 目前黄土高原38.3%的浅沟分布在耕地, 35.3%分布在草地, 22.8%分布在林地(图 10)。其中耕地上浅沟平均密度最大, 为5.46 km/km²；其次是草地, 为5.32 km/km²；林地为4.99 km/km²。

在黄土高原所有调查单元耕地上共发现浅沟225条, 图 11为耕地上浅沟密度的分布, 其中密度较大的区域为北洛河中上游延安市中部以及西北部, 其中延安市中部耕地浅沟密度达到20 km/km²以上, 延安市西北部耕地浅沟密度约为7 km/km²；定西、白银市会宁交界处耕地浅沟密度也较大, 达14 km/km²；无定河下游黄河两岸榆林市绥德、吕梁市一带也是重要的耕地浅沟分布区域, 密度为10 km/km²左右。以上三个区域为目前黄土高原耕地浅沟需重点关注的区域。

3 讨论

本研究基于遥感、GIS方法, 以小流域为抽样调查单元, 结合野外实测方法, 回答了黄土高原区域尺度浅沟统计特征和空间格局特征, 指出浅沟目前所在坡面的土地利用, 并重点分析了仍处在耕地上的浅沟空间分布格局, 对于进一步浅沟的治理决策有重要意义。

近年来高分辨率Google Earth影像在沟道侵蚀调查中应用越来越多, 为大尺度侵蚀沟制图提供了重要数据支撑。本研究结果也说明基于Google Earth影像的浅沟调查结果与野外高精度实测结果无显著性差异, 前者对于浅沟调查应用是适用的。在一些研究中, 如在Karydas和Panagos^[27]对希腊全国尺度浅沟调查中应用了这一数据源。但这一数据源目前常用于对浅沟位置的调查, 本研究中浅沟也特指了瓦背状地形底部这一地貌特征, 而对于浅沟发育动态的调查研究工作中满足需求的影像分辨率水平和大尺度调查方法尚需进一步探索。

本研究对黄土高原浅沟的调查结果与已有研究在某些小流域上的研究结论基本一致。20世纪80年代, 刘元保^[35]基于野外调查得到志丹杏河庙坬浅沟密度19.80 km/km²、志丹杏河米老家庄20.00 km/km²、安塞茶坊13.98 km/km², 其密度较大的可能原因为志丹、安塞地处黄土丘陵沟壑区, 退耕之前的耕地分布也较广^[36—37], 加上降水、植被等条件共同作用导致浅沟密度较大^[38]。张岩等^[15]基于QuickBird影像得到吴起县合沟小流域浅沟密度为36.21 km/km², 而本文调查出吴起附近浅沟密度最大值为21.92 km/km², 可能的原因为合沟小流域为当地浅沟典型分布的小流域, 而本研究采用系统抽样调查, 但这一区域整体浅沟密度较大的结论是一致的。

在各治理分区中, 黄土丘陵沟壑区的浅沟长度平均值最长, 本文得到黄土丘陵沟壑区的浅沟长度平均值为65.74 m, 与姜永清等^[17]利用航空影像分析黄土丘陵沟壑区瓦背状浅沟分布特征得到浅沟平均长度为64.7 m的结果较接近。黄土丘陵沟壑区以峁状、梁状丘陵为主, 地形起伏大, 植被盖度低, 耕作方式落后, 为浅沟的发育提供了条件, 导致该区的浅沟长度较大。

本研究表明黄土高原浅沟临界坡长多介于40—60 m, 与刘元保等^[5]基于野外调查得到的黄土高原浅沟侵蚀凸形坡上临界坡长为58.70 m、直线坡为44.30 m、凹形坡为32.24 m, 张科利等^[11]基于野外调查得到的黄土高原浅沟侵蚀临界坡长平均为40 m左右均较接近。在不同治理分区中, 黄土丘陵沟壑区浅沟密度最大且临界坡长平均值最小, 为46.58 m, 此结果在一定程度上反映了该区发生浅沟侵蚀所需坡长较短, 浅沟发育所需汇水条件较其它治理区低, 这与其地形较陡、土壤性质易蚀、植被覆盖条件较差有关。

浅沟是长期耕作-侵蚀的结果^{[8—9, 11—13]}, 理论上发育于耕地坡面。黄土高原自上世纪末开始的大规模退耕还林(草)工程, 对浅沟的治理起到了重要作用, 目前有60%以上的浅沟已不在耕地坡面上, 但仍有38.3%的浅沟仍在耕地坡面上, 这部分浅沟在暴雨条件下仍可能发生侵蚀-耕作填充-再侵蚀的过程, 极易继续发展从而造成较大侵蚀量, 并有可能进一步成为切沟^[8], 造成更严重的土壤侵蚀, 因此目前仍分布在耕地上的浅沟具有进一步发展的潜力, 对于耕地上浅沟密度较大的区域需重点关注。

在进一步研究中, 亟需对黄土高原区域尺度浅沟分布的差异性原因进行更为深入的分析, 对浅沟主要分布区及耕地上浅沟密度较大的区域降雨、地形、土壤、土地利用、水土保持措施等影响因素的进一步研究亟待开展。同时综合分析多空间尺度浅沟相关研究, 从机理、模拟、格局多角度综合开展黄土高原浅沟研究也将促进浅沟科学治理。

4 结论

本文以亚米级遥感影像为数据基础, 综合运用遥感和GIS技术, 基于137个调查单元对退耕背景下黄土高原的浅沟分布及其特征以及防治重点问题进行研究。主要结论如下:

(1) 亚米级Google Earth影像可以用于浅沟(本文特指瓦背状地形底部)解译, 解译结果与GNSS RTK野外实测结果无显著性差异。

(2) 黄土高原33.6%的调查单元存在浅沟, 出现浅沟的调查单元中, 平均浅沟密度为3.41 km/km², 最大为21.92 km/km², 平均沟数密度为65条/km²。浅沟长度多介于20—60 m, 临界坡长多介于40—60 m。浅沟密度较大的区域为延安市西北部吴起、志丹、榆林市靖边一带以及延安市甘泉、安塞, 榆林市绥德县、吕梁市一带, 定西市、白银市会宁一带以及庆阳市。不同治理分区中, 黄土丘陵沟壑区平均浅沟密度最大, 其次是黄土高塬沟壑区。

(3) 目前黄土高原38.3%的浅沟分布在耕地, 35.3%分布在草地, 22.8%分布在林地。其中耕地上浅沟平均密度最大, 其次是草地、林地。耕地浅沟密度较大的区域有北洛河中上游延安市中部与西北部, 定西市、白银市会宁交界处, 无定河下游黄河两岸榆林市绥德、吕梁市一带。