典型黑土区村级尺度侵蚀沟演变

李浩，张兴义，刘爽，Yury Kravchenko，Kateryna Ivanova，李续峰，陈强

(1.中国科学院东北地理与农业生态研究所，150081，哈尔滨;2.中国科学院研究生院，100049，北京;

3.乌克兰国立生命与环境科学大学，03041，基辅，乌克兰)

东北黑土区是我国最大的商品粮基地，典型地貌类型为漫川漫岗，坡缓坡长，60%为坡耕地。大部分的黑土区已有几十年甚至上百年的耕种历史［1］，尤其是新中国成立以来，大规模的过度垦殖和不合理耕作，昔日有“北大仓”之称的千里沃野已存在不同程度的水土流失，沟蚀严重，且逐年向严重和剧烈侵蚀发展，密度增加，强度加强。沟蚀的危害表现为吞噬耕地，毁坏道路桥梁，破坏原本整齐的坡耕地，造成土地支离破碎，影响机械耕作，损害黑土区农业综合生产能力，对我国的粮食安全问题构成了严重威胁［2］。

利用遥感、地理信息系统技术，结合监测手段，使研究侵蚀沟分布及动态变化已成为现实。目前，国内外对沟蚀演变的相关报道，主要集中于采用卫片判读监测数十年跨度中侵蚀沟的演变及GPS 野外测量监测侵蚀沟短期变化等方面，研究的空间尺度一般为几km2到几百km2，分辨率多在5 m 以上［3－6］;然而，遥感手段在沟蚀演变的应用中也有很多的不足。一方面影像分辨率与研究对象尺度未能完全匹配，受影像分辨率的限制，尺寸低于影像分辨率的实际地物很难在影像中被解译。前人所作研究采用的遥感影像分辨率多在5 m 以上，而野外调查发现，黑土区沟蚀中相当比例的浅沟宽度在5 m 以下;因此，在计算侵蚀沟密度时未能将浅沟充分考虑，在一定程度上减小了侵蚀沟密度，降低了对沟蚀潜在危害的估计，也增大了解译误差。另一方面，21世纪初以来，高分辨率甚至亚米级卫星遥感影像的广泛应用，相近空间分辨率历史影像的缺乏也限制了遥感影像在沟蚀演变方面的应用。笔者以1968年9 月分辨率为1∶1.5 万的黑白航摄影像及2009年6 月分辨率为0.61 m 的Quickbird 模拟假彩色合成影像为主要数据源提取沟蚀信息(根据航摄规范［7］，这2 种遥感影像具有相似的空间分辨率)，并在遥感及GIS 技术支持下，应用亚米级影像，对41 a来典型黑土区村级小尺度研究范围内侵蚀沟的动态变化进行调查研究，旨在明确村级尺度沟蚀现状、演变及其驱动因素，为防治沟蚀和促进农业可持续发展提供科学依据。

1 研究区概况

研究区为村级尺度，总面积24.2 km2，位于黑龙江省海伦市前进乡光荣村(E 126°50'1.42″，N 47°21'12.61″)，属寒温带大陆性季风气候区，一年一熟制，冬季寒冷干燥，夏季温热多雨，年均气温1.5 ℃，极端最高温度为37 ℃，极端最低温度为－39.5 ℃，多年平均降水量530 mm。土壤类型为典型黑土，地形为漫川漫岗，平均坡度2.55°。土地利用类型结构较为单一，70%以上开垦为农田，多为缓坡耕地，林草覆盖率低。

2 数据来源及研究方法

2.1 数据来源

常用的TM 及SPOT 影像受分辨率的限制，很难准确地反映侵蚀沟尤其是宽度低于影像分辨率的浅沟的分布及演化过程。笔者采用的跨度41 a 的影像数据均为亚米级或米级影像，包括1968 年9 月16 日航摄黑白影像，实际航摄比例尺1∶1.5 万，绝对航高1 290 m，单幅影像覆盖地400 m×400 m，航向重叠55%～65%，旁向重叠18%～50%。将14景航摄黑白影像融合生成分辨率为0.6 m 的光荣村1968 年遥感数据。2009 年6 月Quickbird 卫星影像包括多光谱(空间分辨率1.44 m)和全色波段(空间分辨率0.61 m)，为了突出区域的侵蚀沟信息，对Quickbird 影像进行波段合成、图像融合处理，生成空间分辨率为0.61 m 的模拟真彩色合成影像。利用光荣村1∶1万地形图对上述影像进行几何精纠正，误差控制在2 个像元以内，最后统一到西安1980坐标系下。

2.2 数据获取

1)首先进行遥感影像室内初步预判。由于该典型黑土区以坡耕地为主，土地利用结构较为单一，因此，将土地利用类型划分为耕地、草地、林地、柳条通、侵蚀沟、居民地、废弃地。

2)在预判的基础上建立先验知识，发现存在的问题及疑点，然后于2010 年7 月对光荣村进行实地踏勘，详细调查了土地利用状况，侵蚀沟的沟宽、沟长、沟内沉积物及沟内外植被生长状况。

3)对于在1968 年的航摄黑白影像上难以判断的斑块，采取咨询当地人的方式解决。

4)量测部分耕地斑块的垄向坡度、坡长。

5)在建立了侵蚀沟及植被类型的解译标志后，回到室内详细判读，对41 a 的侵蚀沟进行线文件及面文件的矢量化，解译土地利用类型。

6)判读结束后，进行侵蚀沟及土地利用类型的验证及补判工作，最终获得光荣村41 a 侵蚀沟数据、土地利用类型分布数据及部分耕地斑块的坡长垄向坡度数据。

2.3 研究方法

侵蚀沟密度以单位面积上的侵蚀沟总长度来度量，用以反映地表的破碎程度和土壤侵蚀的剧烈程度［8］。在村级尺度内，气象、植被、地质等因素具有一致性，对土壤侵蚀影响最大的自然因素主要为地形地貌，人为因素为耕作措施。采用侵蚀沟密度作为侵蚀动态分析指标，从坡度、坡长等自然因素，以及土地利用变化、垄向坡度等人为因素来分析41 a 来光荣村的沟蚀演变，并探讨自然因素及人为因素对沟蚀变化影响的特点。将其中的垄向坡度定义为在沿着垄向上坡面的垂直高度和水平宽度的比。

3 结果与分析

3.1 沟蚀现状及其演变

在土壤侵蚀分类中，朱显谟［9］根据沟蚀的发展阶段，分为细沟侵蚀、浅沟侵蚀和切沟侵蚀，认为浅沟是细沟发展为切沟的过渡形态，能横向耕作，但不能消除它的痕迹，而切沟侵蚀是沟蚀的最后阶段，以不能横过耕作为其主要特征，沟床的下切至少在1 m 以上，深度也可达20 m 至100 m，横剖面为V型或长U 型，纵断面可以和坡面不相一致。笔者所计算的侵蚀沟密度包含切沟、浅沟，采用窗口移动法，利用ArcGIS 9.2 desktop 中的Focal 函数，对1968 及2009 年的侵蚀沟线状分布进行处理。首先将研究区域均分为300×300 个方格，方格中心间距为25 m，计算单位方格区域内侵蚀沟密度后，将值赋予该方格。栅格化后以像素为计算单元，搜索边长为1 000 m 的正方形内的侵蚀沟密度均值，将该值赋予中心像元，并移动窗口，得到1968 和2009 年侵蚀沟密度分布图(图1)，由二者的差值获得41 a 侵蚀沟密度变化图(图2)。

图1 光荣村1968 与2009 年侵蚀沟密度分布图Fig.1 Distribution of gully erosion density in Guangrong village in 1968 and 2009

3.1.1 侵蚀沟密度分布 计算结果显示，1968 与2009 年侵蚀沟平均密度分别为3.33 和4.96 km/km2，远大于以往在黑土区前人的研究结果［10］。根据黑土区水土流失综合防治技术标准［8］，在1968年，光荣村已有40%的面积侵蚀沟密度大于4 km/km2，即已处于剧烈程度。侵蚀沟宽度大于5 个像元，即宽度大于3 m 的侵蚀沟密度，1968 和2009 年分别为1.68 和2.27 km/km2，沟蚀强度分别达到极强烈和剧烈程度，二者相差2 倍左右。宽度在3 m以下的浅沟、切沟侵蚀沟密度占研究区总侵蚀沟密度的50%，表明亚米级高分辨率遥感影像在沟蚀调查上具有较好的应用潜力。

图2 光荣村41 a 浅沟切沟变化图Fig.2 Development of gully in Gurangrong village during 41 years

根据侵蚀沟密度现状及其变化程度，将研究区侵蚀沟密度分为0 ～2、2 ～4、4 ～6、6 ～8、8 ～10 和＞10 km/km26 个等级。从图1 可以看出:1) 1968和2009 年2 年的高值区均为连片分布，且分布在坡度较大的区域;2)1968 和2009 年平均侵蚀沟密度均大于4 km/km2，沟蚀强度处于剧烈程度［11］的面积比例分别为37%和60%;3)41 a 来侵蚀沟密度高值区面积在扩大，大于10 km/km2的区域从1 片增加到3 片。

3.1.2 侵蚀沟演变 为更好地了解光荣村41 a 沟蚀的演变以及侵蚀沟密度变化的区域差异，将1968和2009 年研究区浅沟切沟的矢量面文件叠加，生成了浅沟切沟变化图(图2)，在侵蚀沟密度空间分布图的基础上，生成41 a 侵蚀沟密度变化图(图3)和侵蚀沟面积与侵蚀沟密度关系图(图4)。将侵蚀沟密度变化值分成－3 ～－1、－1 ～0、0 ～1、1 ～2、2 ～3、3 ～4 和4 ～5.5 km/km27 个等级，以分析侵蚀沟密度的变化趋势。

可知，41 a 间，光荣村侵蚀沟面积从0.86 km2增长到了2.22 km2，所占研究区域总面积的比例由3.5%增长到9.1%，总侵蚀沟密度值增加了1.6 km/km2;但从整个区域来看，侵蚀沟密度变化程度很不均衡，有20%区域沟蚀程度下降，20%区域侵蚀沟密度增加值超过了3 km/km2。在侵蚀沟密度值增加的区域中，随着侵蚀沟密度变化值的增加，相应面积也在减少，且呈现出极强的线性关系，相关性达－0.983。

图3 光荣村41 a 侵蚀沟密度变化图Fig.3 Dynamic of gully erosion density in Guangrong village during 41 years

图4 侵蚀沟面积与侵蚀沟密度关系图Fig.4 Relationship between gully erosion area and gully erosion density in Guangrong village during 41 years

从研究区侵蚀沟密度的剧烈变化情况可知，41 a来，光荣村加大了土地的开垦力度，破坏了原本整齐的坡耕地，造成土地支离破碎。从土地利用转换的解译结果可知，41 a 间，1.41 km2的草地转成了坡耕地，同时却有1.16 km2的坡耕地被侵蚀沟吞噬，耕地面积从19.73 km2降低到了19.46 km2，耕地地块数也从27 块增加到68 块，说明侵蚀沟的发展造成了坡耕地斑块的破碎化。41 a 来，侵蚀沟密度高值区面积在扩大，1968 年，2 块侵蚀沟密度大于8 km/km2的高值区在2009 年已连成一体。

3.2 侵蚀沟动态变化影响因素分析

侵蚀沟是在自然和人为因素共同作用下形成的，另一方面，人为因素也起到了加速或减缓侵蚀沟发展的作用。为了研究典型黑土区自然和人为因素对侵蚀沟发展演化的影响，在自然因素中选取了坡度、坡长因子，在人文因素中选取了垄向坡度、土地利用变换因子，分析各个因子与侵蚀沟演变之间的关系。利用ArcGIS 9.2 desktop 的空间分析模块，结合实地测量，提取了光荣村的坡度、坡向、坡长及垄向坡度等数据。

3.2.1 侵蚀沟密度与坡度关系分析 东北黑土区的地形特点是漫川漫岗，坡度较缓［8］。在Arcgis9.2 desktop 的支持下，利用1∶1万地形图数据建立数字高程模型。汤国安等［10］研究发现，对于1∶1万比例尺的DEM，5 m 是保证地形描述精度的理想分辨率尺度;因此，利用Slope 函数生成了研究区分辨率为5 m 的坡度图(图5)，平均坡度为2.55°。根据黑土区水土流失综合防治技术标准［11］，将研究区坡度分成≤0.25°、0.25°～1.5°、1.5°～3°、3°～4°、4°～5°、5°～8°和＞8°等7 个坡度带，分别计算不同坡度带上2 个时期的侵蚀沟密度以及41 a 侵蚀沟密度变化值。2 个时期的侵蚀沟密度最大值均出现在坡度值为8°的区域(图6)，但在2009 年，在坡度大于8°的区域侵蚀沟密度减少，反映出20 世纪80 年代以来荒岗荒坡地退耕还林的水土保持效益显著。侵蚀沟密度增加比例最大的区域出现在0.25°～1.5°及4°～5°之间，且密度增加绝对值最大区域出现在4°～5°之间。由光荣村土地利用类型分布现状可知，上述3 个区域均是坡耕地。鉴于该区域坡耕地面积占总面积的4/5 强，参照水利部标准，对2 个时期的坡耕地进行坡度分级，分为≤0.25°、＞0.25°～1.5°、＞1.5°～3°、＞3°～4°、＞4°～5°和＞5°等6个等级(表1)。

图5 光荣村坡度分布图Fig.5 Slope distribution in Gurangrong village

图6 光荣村1968 及2009 年侵蚀沟密度与坡度关系Fig.6 Relationship between gully erosion density and slope steepness in Guangrong village during 41 years

表1 光荣村1968 及2009 年坡度与坡耕地侵蚀沟密度、面积关系Tab.1 Relationship between slope steepness and gully erosion density of cultivated land in Guangrong village during 41 years

由表1 可以看出:41 a 间研究区的耕地总面积减少了0.3 km2，但在不同坡度分级范围内坡耕地的面积增减是不一致的;在0.25°～1.5°范围内坡耕地面积增加了0.3 km2，增加比例为15%，在0.25°～4°范围内坡耕地的面积几乎没有变化，而在坡度大于5°的区域内耕地面积减少了0.3 km2，减少比例为12.5%。不同坡度分级范围内，坡耕地面积的增减反映出20 世纪80 年代以来对河套草甸开荒及荒山荒坡地退耕还林的结果。研究区任何坡度等级的侵蚀沟密度变化量都大于1 km/km2，说明41 a 间该区域加大了对所有坡耕地的开垦力度。耕地多分布在0 ～4°范围内，在该区域内侵蚀沟密度变化最大的坡度地带是0.25°～1.5°(图6)，反映出41 a 来该区域尤其加剧了对缓坡土地的开垦，导致侵蚀情况恶化。

3.2.2 侵蚀沟密度与坡长关系分析 东北典型黑土区地貌特征为漫川漫岗，使得该区坡度平缓，坡长较长，因此，以往部分研究指出，坡长是该区侵蚀产生的主要原因［12］。将研究区坡长分成0 ～99，100 ～199，200 ～299，300 ～399，400 ～499 和＞500 m 等6级，分析坡长因素对侵蚀沟演化的影响。

图7 为光荣村41 a 间坡长与侵蚀沟密度的关系图。可以看出，2 个时期侵蚀沟密度在不同坡长范围下的分布具有一定的相似性，侵蚀沟密度第1个峰值均出现在坡长300 m 附近，第2 个峰值出现在500 m 处。在0 ～500 m 内，随着坡长的增加，1968 年侵蚀沟密度呈线性增加，超过500 m 后侵蚀强度降低。2009 年也有类似相关性，但坡长超过400 m 后侵蚀沟密度开始降低，这可能是受到耕作措施的影响。侵蚀沟密度最大值及密度最大增加值均出现在300 m 附近，与前人的研究结果［13］相似，也说明了并不是坡长越长侵蚀就越强。

图7 光荣村41 a 间坡长与侵蚀沟密度的关系Fig.7 Relationship between slope length and gully erosion density in Guangrong village during 41 years

3.2.3 侵蚀沟密度与垄向坡度关系分析 黑土区的坡耕地具有坡缓且长的特点，合理的耕作措施可有效防治坡耕地水土流失。实践证明，等高耕作可以显著减少水土流失。实际操作中，为了防治坡耕地水土流失，人们通常将顺坡垄作改为横坡垄作来降低地表径流对耕层土壤的冲刷力［14］;但由于微地形的起伏变化，横坡垄很难做到绝对等高，垄向仍然存在不同程度的坡度［14－15］。野外调查也发现，本研究区域顺坡垄较少，垄向多为介于顺坡垄与横坡垄之间。对比1968 与2009 年影像，41 a 间部分区域耕地的垄向也发生了一定的变化;因此，利用2 个时期的坡向、坡度及垄向分布图层，绘制1968 和2009年坡耕地的垄向坡度分布图，将垄向坡度进行分级，并以此计算侵蚀沟密度及变化与坡耕地垄向坡度之间的关系，结果见图8 和表2。

图8 光荣村41 a 间侵蚀沟密度与坡耕地垄向坡度的关系Fig.8 Relationship between gully erosion density and gradient along the ridge of cultivated land in Guangrong village during 41 years

表2 光荣村41 a 间侵蚀沟密度与坡耕地垄向坡度分布Tab.2 Distribution of gully erosion density and slope steepness along the ridge of cultivated land in Guangrong village during 41 years

由图8 可知，2 个时期侵蚀沟密度在不同垄向坡度范围内的分布趋势以及41 a 间的侵蚀沟密度变化具有一定的相似性。由表2 可知，2 个时期超过60%的坡耕地的垄向坡度介于0 ～0.5°之间，在此坡度区间内，侵蚀沟密度随着垄向坡度的增加而升高，且侵蚀沟密度最大值均出现在垄向坡度0.5°处，在0.5°～1.25°区间内，侵蚀沟密度随着垄向坡度的增加而降低，这可能与坡长的增加有关。研究区域1968 年垄向坡度平均值为0.44°，2009 年垄向坡度平均值为0.40°，其中垄向坡度大于0.25°的坡 耕地面积均在减少，而垄向坡度在0 ～0.25°之间的 面积大幅增加，反映了近年来人们为了减少水土流 失而采取的改变垄向来降低垄向坡度的效果。

3.2.4 光荣村41 a 土地利用变化 1968—2009年，随着人口的增加及对研究区开发力度的增大，光 荣村的土地利用结构发生了一定的变化。利用光荣 村土地利用类型转移矩阵(表3)，对侵蚀沟密度变 化的原因做出解释。

表3 光荣村41 a 间土地利用类型面积转移矩阵Tab.3 Land use area transformation matrix in Guangrong village during 41 years km2

由表3 可知，41 a 间光荣村草地面积从1.77 km2下降到0.05 km2，其中1.41 km2，即80%的草地被开垦成耕地，且废弃地的消失说明该区域近年来加大了对土地的开垦。光荣村从19 世纪末开始垦殖，到1968 年林地只剩0.11 km2。为了恢复当地的植被，自20 世纪80 年代实施退耕还林措施以来，林地面积增加到了0.97 km2，其中0.69 km2来自于耕地。41 a 间侵蚀沟面积从0.86 km2增加到2.22 km2，其中有1.16 km2来自于耕地，0.10 km2来自于草地。

41 a 来，光荣村加大了土地的利用强度，草地被大量开垦，从而导致侵蚀沟面积的快速发展。另一方面，林地面积的大大增加得益于退耕还林政策。

4 结论与讨论

1)应用亚米级遥感影像解译的侵蚀沟密度远高于以往研究结果。2009 年光荣村侵蚀沟密度达到了4.96 km/km2，是以往相关文献中数值的2 倍左右。实地调查发现，该区域宽度在3 m 左右的浅沟占沟蚀的1/2 以上。在以往侵蚀沟研究中，沟宽低于遥感影像分辨率的部分在本文中被解译，是侵蚀沟密度显著提高的原因，说明黑土区沟蚀强度需要被谨慎对待，而未被评估的浅沟可反映潜在沟蚀危害。

2)典型黑土区沟蚀仍呈加剧的态势。41 a 来，光荣村侵蚀沟密度从3.33 km/km2增加到了4.96 km/km2，但侵蚀沟密度变化并不均匀。与1968 年相比，2009 年高值区扩展、集中连片，同时也有20%的区域侵蚀沟密度下降。侵蚀沟密度的增加造成了耕地斑块的破碎化，在耕地总面积几乎不变的情况下，耕地块数增加了2 倍。

3)坡度是影响侵蚀沟演变的重要自然因素。随着坡度的增加，侵蚀沟密度增加值的最大值出现在坡度0.25°处，说明近年来光荣村加大了对草地及坡耕地的开垦力度，加速了缓坡地的水土流失。在坡长方面，侵蚀沟密度最大值及密度增加值的最大值均出现在300 m 处。

4) 农耕是加速沟蚀的重要人为因素。41 a 间，坡耕地不同垄向坡度的侵蚀沟密度变化比较一致，光荣村坡耕地的平均垄向坡度在降低。土地利用类型发生了一定的变化，主要表现为大部分草地向耕地转化，坡度较大区域的坡耕地被退耕还林，新增加的侵蚀沟面积中，大部分来源于耕地及草地。

本文的不足之处在于仅考虑了侵蚀沟的长度及面积，没有加入侵蚀沟的深度及体积数据。在接下来的工作中，需要将研究重点放在如何通过遥感与GIS 手段获得这2 种数据，这对于估测侵蚀沟侵蚀量具有重要的意义。

［1］ 汪景宽，王铁宇，张旭东，等.黑土土壤质量演变初探［J］.沈阳农业大学学报，2002，33(1):43－47

［2］ 水利部，中国科学院，中国工程院.中国水土流失防治与生态安全:东北黑土区卷［M］.北京: 科学出版社，2010:67－72

［3］ Vandekerckhove L，Poesen J，Govers G.Medium－term gully headcut retreat rates in Southeast Spain determined from aerial photographs and ground measurements［J］.Catena，2003，50:329－352

［4］ Martinez－Casasnovas J A.A spatial information technology approach for the mapping and quantification of gully erosion［J］.Catena，2003，50:293－308

［5］ 胡刚，伍永秋，刘宝元，等.GPS 和GIS 进行短期沟蚀研究初探［J］.水土保持学报，2004，18(4):16－41

［6］ 闫业超，张树文，岳书平.基于Corona 和Spot 影像的近40 年黑土典型区侵蚀沟动态变化［J］.资源科学，2006，28(6)，154－160

［7］ 国家测绘局.GB/T 15661—1995 1∶5 000，1∶10 000，1∶25 000，1∶50 000，1∶100 000 地形图航空摄影规范［S］.北京: 测绘出版社，1996

［8］ 唐克丽，史立人，史德明，等.中国水土保持［M］.北京:科学出版社，2004:146

［9］ 朱显谟.黄土区土壤侵蚀的分类［J］.土壤学报，1956，4(2):99－114

［10］汤国安，赵牡丹，李天文，等.DEM 提取黄土高原地面坡度的不确定性［J］.地理学报，2003，58(6): 824－830

［11］中华人民共和国水利部.SL 446—2009 黑土区水土流失综合防治技术标准［S］.北京:中国水利水电出版社，2009:5－6

［12］张晓平，梁爱珍，申艳，等.东北黑土水土流失特点［J］.地理科学，2006，26(6):687－692

［13］刘宝元，阎百兴，沈波，等.东北黑土区农地水土流失现状与综合治理对策［J］.中国水土保持科学，2008，6(1):1－8

［14］张兴义，张少良，刘爽，等.严重侵蚀退化黑土农田地力快速提升技术研究［J］.水土保持研究，2010，17(4):1－5

［15］杨爱民，沈昌蒲，刘福，等.坡耕地垄作区田水土保持效益的研究［J］.水土保持学报，1994，8(3):52－58