黄土高原切沟空间分布特征

2023-01-09苏建华王春梅庞国伟杨勤科杨丽娟杨安南刘宝元

水土保持研究 2023年1期

苏建华, 王春梅, 庞国伟, 杨勤科,仲原, 杨丽娟, 杨安南, 刘宝元

(1.陕西省地表系统与环境承载力重点实验室, 西安 710127; 2.西北大学城市与环境学院, 西安710127; 3.北京师范大学自然科学高等研究院, 广东珠海 519087; 4.北京师范大学地理科学学部地表过程与资源生态国家重点实验室, 北京 100875)

黄土高原是全球土壤侵蚀最为剧烈的地区之一，切沟是土壤侵蚀发育到最严重阶段的表征，特别是处于发育活跃期的切沟，是黄土高原重要的产沙来源，其产沙量占整个流域的一半及以上[1-2]，对黄河中游流域生态治理带来了极大的风险与挑战[3]。在区域尺度明确切沟空间分布特征对黄土高原切沟侵蚀治理具有重要意义。

切沟是不能被普通耕作工具横跨的侵蚀沟[4]，沟宽和沟深一般都超过50 cm，其纵剖面与所在坡面基本一致[5]，刘宝元将宽度大于10 m的切沟定义为大切沟，小于10 m的为小切沟[6]，大、小切沟形态、发育速率、空间分布差异均较大，小切沟经过沟头的溯源侵蚀，沟底下切，沟坡冲淘等侵蚀形式[7]可逐渐发育成为大切沟。近年来国内外学者就切沟区域尺度空间分布的研究，受限于基础数据与方法鲜有报道，仅有为数不多的国内外学者对其进行了研究。国内方面，杜国明等以黑龙江省宾县漫岗区为研究区，结合实地调查，以SPOT5遥感影像为基础数据[8]，基于人机交互方式提取切沟，分析了东北典型黑土漫岗区切沟的空间分布特征。董一帆等使用Google Earth高清影像结合ArcGIS的方法[9]，基于可见性的基本原则，对横断山区开展了侵蚀沟抽样调查，发现横断山区切沟和冲沟分布广泛，发育活跃。李镇等基于QuickBird影像选取了陕北黄土区吴起县合沟与绥德县桥沟流域作为研究区[10]，发现QuickBird影像目视解译合沟和桥沟小流域切沟面积、周长的平均相对误差在5%左右，对于小流域切沟监测具有准确、便捷的特点。国外方面，希腊panagos团队基于Google Earth影像[11]，对希腊全国耕地进行随机抽样调查，利用GIS技术进行目视解译，调查发现希腊北部和东部出现了较严重的沟蚀空间集群，对希腊耕地上发育的侵蚀沟在国家尺度有了新的认识。综上，基于高分辨率遥感影像调查的方法，为侵蚀沟空间分布特征研究提供了新的机遇，目前黄土高原切沟区域尺度空间分布认识不足，尤其是对大、小切沟空间分布差异性研究缺乏，严重阻碍了黄土高原侵蚀沟的系统治理。

为了在区域尺度回答黄土高原切沟的空间分布格局问题，本研究基于Google Earth高清影像，运用GIS空间分析与统计方法，以小流域为单元进行系统抽样，通过目视解译的方式，提取切沟密度、长度、宽度、距分水岭距离、切沟发育所在坡面土地利用类型等。本研究将有助于进一步明确黄土高原的切沟空间分布特征，以支持黄土高原的切沟治理，促进黄土高原生态保护和高质量发展。

1 材料与方法

1.1 研究区概况

黄土高原位于中国中部偏北部，区域范围：33°41′—41°16′N，100°52′—114°33′E，面积约64.6万km2，属大陆性季风气候区，是黄河流域的主要产沙区，黄河中游的主要侵蚀区，区域内地形破碎，沟壑纵横[12]，极易发生切沟侵蚀。本研究通过系统抽样方法，在黄土高原均匀布设256个抽样单元，每个抽样单元面积约为0.2 km2，对于地形起伏较为平坦的抽样单元(如居民用地、生活用地、未发育切沟的农地)，选择在抽样单元的中心位置绘制500 m×500 m的矩形区域，见图1。

1.2 基础数据

本研究基础数据主要分为影像数据、高程数据、标准地图数据。影像数据为Google Earth影像，亚米级分辨率，影像选取过程遵循两个原则：(1) 选择近期拍摄且切沟辨识度较高的影像(多数影像年份在2019—2021年期间，极少数区域近年影像因天气原因导致切沟辨识度较低的，选择2019年以前切沟辨识度较高的影像)。(2) 为避免植被影响切沟解译，主要选择冬春两季影像；高程数据为SRTM数据集，空间分辨率为30 m，来源于网站(http:∥srtm.csi.cgiar.org)，用于提取流域边界与坡度；标准地图数据包括黄土高原边界及其边界内省界、省会城市、地市、主要河流等数据，边界数据来源于黄土高原科学数据中心网站(http:∥loess.geodata.cn)，省界、省会城市、地市、主要河流等数据来源于中国自然资源部标准地图服务中心网站(http:∥bzdt.ch.mnr.gov.cn)。

1.3 切沟解译方法与流程

1.3.1 解译方法本研究基于人工目视判读方式解译切沟，主要包括确定流域单元、切沟形态特征、土地利用类型3个方面的解译。确定流域单元包括分水岭、沟沿线解译。切沟形态特征解译包括切沟长度、宽度、距分水岭距离等解译。土地利用类型解译主要是指切沟发育所在坡面的土地利用方式。

图1 研究区抽样单元分布及抽样单元示例

使用SRTM数据集生成分水岭，同时查看Google Earth遥感影像，沿着流域山脊线进行手动勾画可获得抽样单元流域边界，即流域分水岭。通过SRTM数据集生成的流域坡度，同时结合Google Earth影像，沿着沟边坡度发生剧烈变化处进行手动勾画可获得流域沟沿线。从切沟沟头位置出发，沿切沟汇流路径进行手动勾画，到切沟沟尾处结束，可获得切沟的长度。分别在切沟上部、中部、下部各测量一次宽度，计算三次测量宽度均值作为切沟宽度。从切沟沟头上方分水岭出发，沿垂直等高线方向，从上往下手动勾画至切沟沟头位置处结束，可获取切沟距分水岭的距离，对于矩形抽样单元中的切沟，则按照矩形单元所在的完整流域确定流域的分水岭进而获取切沟距分水岭距离。切沟所在坡面的土地利用类型解译按照GLC-30土地利用解译方法及标准[13]，以人工目视解译方式判读获取。切沟解译示意见图2。

1.3.2 解译流程与质量控制为保证切沟解译的精度，本研究采用了严格的解译流程和质量控制。基本流程包括：第一步，基于切沟定义分析[4]及专家论证，制定完整解译标准。第二步，随机筛选30个流域单元进行预解译，逐一对比解译结果，进一步统一切沟的提取标准，确保不同人员对切沟抽样单元的解译结果相似度达到95%以上，同时依据解译过程中出现的问题进行解译标准优化。第三步，完成所有抽样单元的切沟解译。第四步，对所有解译数据进行质检，解译结果进行三层质检后，合格数据保存至数据库。三层质检具体指：(1) 解译人员组内互查，发现有错误的解译数据，及时进行修正。(2) 解译工作小组质检员对切沟解译数据进行100%检查，发现问题反馈给解译人员进行修正，并对修正后数据重新进行检查，直至通过质检员检查。(3) 组织专家对解译数据进行20%抽查，如错误率超过5%，则对数据进行逐一修订并再次进行三层质检(图3)。

注：图2所示为陕西省榆林市子洲县王武沟村附近典型切沟小流域(37°43′15.48″N，109°57′44.33″E，海拔1054.32 m)。

1.4 统计与分析方法

本研究获取了计算切沟各空间分布参数的原始数据，聚类和异常值分析(Anselin Local Moran I)方法[14]，可对给定的一组加权要素，识别具有统计显著性的热点、冷点和空间异常值，并划分为4种空间分布类型，具体如下：

(1)

(2)

式中：n等于要素的总数目。统计数据的zIi得分的计算方法如下:

(3)

其中：

(4)

(5)

若zIi得分是较高的正值，则表示周围要素具有相似的高值或者低值，会分别将具有统计显著性高值、低值聚类表示为HH，LL空间聚类，如果zIi得分是较低的负值，则表示存在一个具有统计显著性的空间数据异常值，当低值被高值包围时，表示为LH空间异常，相反，则表示为HL空间异常。

图3 切沟解译流程与质量控制

2 结果与分析

2.1 切沟参数统计特征

本研究在黄土高原共计布设256个抽样单元，统计得出：存在切沟、大切沟、小切沟抽样单元占比分别为35.94%，9.77%，35.16%。切沟密度主要介于0.01～4 km/km2，整体密度均值为1.47 km/km2，其中存在切沟的92个流域平均切沟密度为4.08 km/km2，主要为小切沟(表1—2，图4)。

切沟长度一般小于60 m，平均长度43.53 m，小切沟平均长度与切沟接近，相比小切沟，大切沟长度更长，平均值达70.90 m。超过一半的切沟，宽度在8 m以内，以小切沟为主，绝大部分大切沟宽度介于10～14 m，平均宽度13.49 m。与小切沟相比，大切沟距流域分水岭更远，平均距离达89.27 m，小切沟则为70.31 m。

表1 切沟密度特征统计

2.2 切沟空间分布格局

2.2.1 黄土高原土壤侵蚀二级分区切沟空间分异分析黄土高原以水力侵蚀为主，是我国土壤侵蚀发生最为严重的地区，绝大部分属黄河中游。黄秉维依据土壤侵蚀营力、类型、发展途径及治理趋势，将黄河中游地区划分为9个二级侵蚀区[16]，其中黄土丘陵沟壑区面积最大，下设5个副区(黄土丘陵I副区、Ⅱ副区、Ⅲ副区、Ⅳ副区、V副区)，本研究根据划分的二级侵蚀分区统计得出：绝大多数切沟抽样单元分布在黄土丘陵沟壑区，该地区切沟密度均值最大，达2.67 km/km2。黄土高塬区切沟长度最大，均值达114.23 m。干旱草原区虽然有较多抽样单元分布，但是存在切沟抽样单元较少，主要分布区靠近黄土丘陵I副区(榆林—东胜一带)。具体数据参看图5及表3。

表2 切沟形态特征统计

图4 切沟参数频率分布

注：林区—石质山地与林区—黄土丘陵统一归为林区。

2.2.2 密度空间分布热点分析结合反距离插值方法可以快速判断出黄土高原切沟密度变化较大地区的空间分布情况。黄土高原切沟主要沿400 mm等降雨量线分布，特别是延安及其以西至固原一带，榆林及其以北至东胜一带(图6)。在95%～99%置信区间内，黄土高原西部(西宁—同仁一带)、中部(西峰及其以北)、东部(晋陕黄河两岸延安—东胜一带)，切沟密度呈现显著的空间热点分布，发生高值聚集现象，且这些地区切沟密度值大于7 km/km2。在90%～95%置信区间内，西宁及其以南、固原及其以北、延安及其以西，切沟抽样单元密度也呈现出次热点空间分布，为切沟密度第二梯度地区，这些地区切沟密度介于3～7 km/km2。大切沟分布较为零散，主要的集中区为天水—定西一带，在黄土高原东北部及西峰至榆林沿线一带及西宁、呼和浩特附近也有分布。小切沟主体分布在400 mm等降雨量线两侧，在95%～99%置信区间内，黄土高原中部(西峰及其以北)、东部(东胜及其以西)，小切沟抽样单元密度呈现显著空间热点分布，且具有中间高四周低的特点，这些地区小切沟密度值均大于7 km/km2。在90%～95%置信区间内，延安以西的小切沟密度呈现出次热点空间分布特征，该地区小切沟密度介于3～7 km/km2。

表3 黄土高原土壤侵蚀二级分区切沟参数统计值

2.2.3 切沟长度、宽度、距分水岭的距离空间格局聚类与异常值分析结合反距离插值方法可以快速判断出黄土高原切沟长度较大地区的空间分布。黄土高原中部(西峰及其以南)、东北部(朔州—忻州一带)、西南部(定西—天水一带)，切沟长度在空间上呈现明显的HH型分布(图7)，发生高值聚类现象，且这些地区的切沟平均长度值介于50～150 m。黄土高原西南部(天水一带)，大切沟长度也呈现明显的HH型空间分布，该地区大切沟长度值介于80～150 m。小切沟长度较大地区主要位于黄土高原南部(铜川—宝鸡一带)、东北部(离石—忻州一带)。

黄土高原西南部(天水及其以南)，切沟宽度最大，其值大于12 m，其次为中部(榆林一带)、东北部(朔州及其以西)，切沟宽度介于7～10 m，东部(晋陕黄河两岸延安—东胜一带)，切沟宽度值也较大(图8)。大切沟在天水一带宽度最大。小切沟在东部(晋陕黄河两岸延安—东胜一带)、南部(铜川及其以西)宽度值较大且呈现明显的空间聚集特征，其值介于5～7 m。

黄土高原西南部(定西—天水一带)，切沟距分水岭距离呈现明显的空间聚集现象，且数值较大，均大于150 m，其次为南部(铜川及其以西)、中部(西峰及其以东)、东南部(朔州西北地区)，其值介于70～150 m，局部地区大于150 m(图9)。大切沟距分水岭距离较大地区主要集中在黄土高原西南部(定西—天水一带)，这些地区值介于100～150 m，局部地区大于150 m。小切沟距分水岭距离较大地区主要位于黄土高原南部(铜川及其以西)、中部(西峰以东与固原以西)、东北部(朔州以西)，其值介于70～100 m。

2.2.4 切沟所在坡面土地利用类型目前黄土高原土地利用类型主要以草地(41.33%)和耕地(32.60%)为主[17](图10)。本研究主要探讨切沟发育所在坡面的土地利用类型，按照土地利用一级分类[18]对2 154条切沟分析发现：目前切沟所在坡面主要土地利用类型为草地，对应存在切沟抽样单元占比接近一半，其次为耕地(29.76%)和林地(17.27%)。

分别对大、小切沟目前发育所在坡面土地利用类型按照单沟统计发现：大、小切沟发育所在坡面的主要土地利用类型均为草地，其次为耕地与林地。具体统计见表4。

3 讨论与结论

3.1 讨论

本研究基于Google Earth亚米级影像，采用系统抽样方法，以256个小流域为调查单元对切沟进行目视解译，得到了区域尺度黄土高原切沟空间分布，弥补了以往研究认识的不足。在黄土高原区域尺度，已有研究者对沟谷密度进行过探讨，其主要研究对象为冲沟、干沟甚至河沟等大型沟谷，沟长一般大于100 m，与本研究所研究对象切沟有较大差异。如景可等和田剑等分别得出了黄土高原区域尺度沟谷密度的空间分布，且整体以晋陕黄河两岸多沙粗沙区密度最大，并以此为中心向外围密度递减[19-20]。本研究表明切沟主要沿400 mm等降雨量线分布，在晋陕黄河两岸也是密度高值区，但不是唯一的高值区，与沟谷密度空间分布存在一定差异。沟谷密度与土壤侵蚀有一定关联，因为沟谷为主要的泥沙输运路径[20]，但切沟是重要的产沙来源，从侵蚀治理角度，切沟是黄土高原沟道侵蚀防治的重点对象，是土壤侵蚀达到严重程度的重要表征[4]，本研究得出的切沟空间分布特征可为后续黄土高原切沟治理工作提供一定依据。

图6 切沟、大切沟、小切沟密度空间分布

400 mm等降雨量线是中国半湿润与半干旱地区的分界线[21]，本研究发现黄土高原切沟主要沿这一等降雨量线分布，可能的原因是在更少的降雨量地区(如黄土高原西北部)虽然植被覆盖度更低，但降雨侵蚀力较弱，侵蚀动能不足以使这一地区切沟发育更广泛。而在降雨量更大的地区(如黄土高原东南部)植被覆盖度增加，很大程度上抑制切沟发育[22]。因此，切沟的分布是降雨侵蚀动能与地表抵御侵蚀双重影响下的结果。

极端降雨事件增加[23]可直接影响切沟的产生及发育空间分布格局，其产生的大量地表径流迅速汇集使得侵蚀动能快速增大[24]，会引起已发育切沟的沟头溯源、沟底下切等侵蚀形式的快速扩张。本研究中位于黄土高原南部(定西—天水一带)、中部(西峰以东)的切沟，其沟头距分水岭距离较远，沟头上方汇水面积较大，地表径流到达切沟沟头的侵蚀动能较强，若后续该地区降雨增大或遭遇暴雨等极端天气时，势必会引起切沟发育的强烈响应[25]，为防止极端天气使得切沟发育增强，黄土高原切沟防治监测过程应重点关注该地区。

图7 切沟、大切沟、小切沟长度空间分布

图8 切沟、大切沟、小切沟宽度空间分布

本研究发现，在退耕还林(草)工程实施背景下，目前黄土高原切沟所在坡面土地利用类型主要为草地，但仍有近30%的切沟发育在耕地上。前人研究发现，当植被覆盖度大于60%时[26-27]。发育在林地与草地坡面上的切沟受到植被抑制作用明显，大部分较为稳定，而耕地上的切沟因翻耕、缺少有效植被保护等，土壤侵蚀更为剧烈，切沟发育危险性更大[28-29]。因此后续治理过程应当重点关注发育在耕地上的切沟。

图9 切沟、大切沟、小切沟距分水岭距离空间分布

黄土高原区域面积广大，基于无人机、卫星、航空获取的遥感影像，由于数据源成本和调查方法的限制，黄土高原区域尺度的切沟空间分布研究较少，但近十几年来随着Google Earth亚米级遥感影像基本覆盖全球，为区域尺度切沟空间分布研究提供了重要数据基础。很多学者也基于此数据基础展开了大尺度切沟研究，如董一帆等调查了横断山区侵蚀沟空间分布[9], Karydas C等对希腊全国侵蚀沟进行了调查[11]。从方法上来看，人工目视解译是目前较为可信的调查方法，但也有其他方法需要在进一步研究中深入探讨与应用，如近年来机器学习算法在切沟空间分布研究中取得了较大进展[30]，外国学者的研究区域虽不是黄土高原，但是其思路与方法可以较好的推广至黄土高原地区，今后通过充分运用切沟自动化提取算法可以快捷、便利的获取黄土高原的切沟空间分布特征，而本研究的相关结果也将为下一步切沟自动化算法研究模型的构建和验证提供数据支撑。