不同DEM分辨率下的坡面土壤侵蚀模型的坡长因子提取对比研究

孔锋　王一飞　吕丽莉　闫绪娴

摘要 采用基于GIS的非累计流量的坡长直接提取算法（NCSL）和空间分析提取法（SAC）对陕北省安塞县10个样区的坡长进行提取，并计算相应的坡长因子。对不同DEM分辨率下的坡长提取结果对比，结果表明：NCSL对坡长的提取精度明显好于SAC的提取结果，其中DEM在5 、10 m分辨率下的提取精度最好，且二者计算的坡长因子值差异不大，因此，可采用10 m分辨率NCSL方法提取坡长。

关键词 DEM；分辨率；坡长；精度差异；土壤侵蚀模型

中图分类号 S157 文献标识码 A 文章编号 0517-6611（2017）29-0104-05

Comparison of Slope Length Factor Extraction in Hillslope Soil Erosion Model with Different DEM Resolution

KONG Feng1，2，3，4，5， WANG Yifei1， L Lili1，2 et al

（1.China Meteorological Administration Training Center， Beijing 100081；2. China Meteorological Administration Development Research Center， Beijing 100081；3.State Key Laboratory of Earth Surface Processes and Resource Ecology， Beijing Normal University， Beijing 100875；4.Academy of Disaster Reduction and Emergency Management， Ministry of Civil Affairs & Ministry of Education， Beijing 100875；5. College of Management Science and Engineering， Shanxi University of Finance&Economics，Taiyuan， Shanxi 030006 ）

Abstract In this research，noncumulative slope length （NCSL） calculation method and spatial analytical calculation （SAC） method were respectively applied to extract slope length and slopelength factor from 10 sample areas， which were located in Ansai County， north Shaanxi Province. The comparison of computation precision among variable DEM resolutions showed that NCSL was superior to SAC entirely. And the results were best when DEM resolutions was 5m and 10 m；besides，the result of slopelength factor was nearly the same under the same condition. So DEM of 10 m resolution can be used to extract slopelength.

Key words DEM；Resolution；Slope length；Precision differentiation；Soil erosion model

基金項目 中国气象局气象软科学重点项目“基于综合风险防范视角的中国及周边国家安全和全球战略研究”（2017〔21〕）；中國气象局气象软科学自主项目“新常态下中国自然灾害风险时空格局和综合防灾减灾工作的现状、趋势、挑战及战略对策范式研究”（2017〔36〕）；中国气象局气象软科学自主项目“中国气象灾害防御能力评估及政策建议”（2017〔35〕）。

作者简介 孔锋（1986—），男，山西临汾人，助理研究员，博士，从事自然灾害与环境演变研究。

收稿日期 2017-08-02

快速城镇化背景下人类的不合理生产活动使土壤侵蚀的强度日渐加剧，由此引发严重的水土流失，不仅造成全球范围内的土壤肥力下降，土地生产力丧失，而且导致诸多江河的河床抬高，水库库容减少以至报废，进而加剧洪涝等灾害的发生。土壤侵蚀已经成为人类生存与可持续发展的世界性重大环境问题之一。

1 土壤侵蚀模型与坡长提取

1.1 土壤侵蚀模型进展

为了更好地掌握土壤侵蚀发生的规律和机理，减小因严重的水土流失而带来的灾难和损失，世界各国的科学家相继在该领域开展了大量的研究工作。建立定量的模型是研究土壤侵蚀的重要手段之一。 Wischmeier等[1]于1959年提出了以下通用土壤流失方程USLE（universal soil loss equation）：

A=R×K×LS×C×P （1）

式中，各因子分别代表降雨侵蚀力、土壤可蚀性、坡长、坡度、覆盖和管理以及水土保持措施。该公式发表后，世界各国学者又不断地修订并改进上述方程中的各因子，这些改进成果被统称为RUSLE（revised universal soil loss equation）[2-3]。基于上述思路，世界许多国家相继建立了更为适合本国实际的土壤侵蚀模型。例如我国学者刘宝元等[4]曾提出中国土壤流失预报方程，将上述公式中的覆盖和管理因子以及水土保持措施因子变为水土保持生物措施（B）、工程措施（E）与耕作措施（T）3个因子：

A=R×K×LS×B×E×T （2）

RUSLE是目前应用最多的水蚀预报经验模型。近年来，该领域又有新的进展，特别值得一提的是基于土壤侵蚀机理过程（水蚀和风蚀等）的物理模型的出现，开辟了土壤侵蚀定量研究的新纪元。

坡长指从地表径流原点到坡度减小直至有沉积出现的地方之间的距离，或到一个明显渠道之间的水平距离，这是关于完整坡的定义[1-3]。然而在绝大多数情况下，实际的坡面形态是非常复杂的，故上述完整坡长的定义不再适用。为此，Foster等[3]提出分段处理不规则坡面的方法，该方法视每一段坡长为上游各个分段坡长的累加，我国学者汤国安等[5-7]据此将坡长定义为地面上一点沿水流方向到其源点的最大距离在水平面上的投影。该定义是基于DEM坡长自动提取的理论基础。在上述通用土壤流失方程USLE和改进后的通用土壤流失方程RUSLE中，坡长是以坡长因子的形式出现的，用来描述坡长与土壤侵蚀量的关系。在RUSLE中，坡长因子的计算公式为：

L=（λ/22.13）m （3）

式中，λ为实际坡面的水平投影距离；22.13是标准化坡长的水平投影长度，其单位是m；m为可变的坡长因子指数，其取值与细沟侵蚀与细沟间侵蚀的比例有关。对于不规则坡和分段坡，人们作分段处理，将每个分段内的地形特征看作是均匀的，则从坡顶开始第i分段对应的坡长因子计算公式为：

Li= λm+1i-λm+1i-1 （λi-λi-1）22.13m （4）

式中，λi为第i段坡底到上游边界的长度。这种方法考虑坡面的复杂性，为坡长因子量化提供了依据。如果考虑到水流的二维性，坡面上某点的径流量和土壤侵蚀量并非完全取决于坡长值，而取决于向该点汇水的单位汇水面积，这一观点更能体现径流的分散和汇集，其水文学意义更明显。接着，Desmet等[8]提出了在二维情况下土壤侵蚀模型中的坡长因子应由坡长长度转为单位汇水面积，例如单位等高线长度的上游汇水面积。基于上述思路，Kinnell[9]对USLE中的地形因子的计算方法进行了改进，得出了一种基于DEM自动提取的新算法。该算法用单位汇水面积A代替坡长λ：

Li= Am+1i，j-out-Am+1i，j-in （Ai，j-out-Ai，j-in）22.13m （5）

式中，Ai，j-in表示流入i，j格网中的上游汇水面积，Ai，j-out代表从该网格流出的汇水面积，后者由于包含网格本身的面积而大于前者，上式又可改为：

Li，j= （Ai，j-in+g2）m+1-Am+1i，j-in g2xi，j（22.13）m （6）

式中，g为格网分辨率，xi，j为该格网相对应的有效等高线长度系数。对于RUSLE[1-3]，坡长因子计算公式为：L=（m+1）（λ/22.13）m。转化为单位汇水面积形式为：

Li，j=（m+1） 2Ai，j-in+g2 2gxi，j（22.13） m （7）

在实际应用中，人们采用GRASS-GIS进行坡长因子的计算[10-12]。国内有人曾经使用该算法求得坡长因子，然后带入公式（1）中反演实际坡长，结果表明这种方法只能整体上反映区域内的坡长分异规律，而坡长值的提取结果精度较低[13-16]。

1.2 4种坡長提取算法对比

基于GIS的地形因子提取多研究坡度、地形等因子，而对其在坡长提取中的应用研究尚不很成熟。目前，非累计流量的直接计算方法（noncumulative slope length，NCSL）、基于累计流量的单位汇水面积计算法（specific contribution area，SCA）、

基于水流强度指数的间接计算法（stream power index）、空间提取算法是比较常用的、也是相对较成熟的基于GIS的坡长自动提取算法，被广泛应用于水利、国土、交通、环保和农业等部门。

1.2.1 非累计流量的直接计算方法（noncumulative slope length，NCSL）。该算法主要基于“水流方向与最大坡降的方向一致”這一假设。其主要原理是利用“八联通”法则，水流由源头（如山脊线，山峰等）像素出发，朝着最大坡降方向像素流动，依次类推， 得出水流路径，进而求得坡长。该算法的特点是算法简单，效率较高，尤其是对洼地和平坦区域有较强的处理能力。但是该算法认为水流方向与最大坡降一致，实地试验表明这与现实的地形和水流情况并非完全吻合[15]。

1.2.2

基于累计流量的单位汇水面积计算法（specific contribution area，SCA）。该算法认为水流分布具有分散性质，其基本思想是用多流向算法得到单位汇水面积，进而取代方程中的坡长因子。由于汇水面积的计算是通过流量累积的方式得到，故称之为基于流量累积的坡长计算方法。该算法解决了方法1中的水流流向问题，但是算法设计极为复杂，并且要考虑诸多特殊情况，因此虽然该算法能够精确地算出汇水面积，但是程序复杂且运行可靠性不高而应用受限[14]。

1.2.3

基于水流强度指数的间接计算法（stream power index）。该算法是基于流量指数的坡长，非直接物理意义上的坡长，而是计算坡长坡度的合成因子（USLE中的LS因子），并认为LS因子是地表径流输沙能力的度量，从而将代表地面曲面形态的LS因子的计算解译为流量和坡度呈非线性函数关系的无量纲输沙能力的指数的计算。值得注意的是，基于水流强度指数的间接计算法得到的坡长其意义已经不是地理学中的物理坡长。该算法较少应用于地理学领域，而在水土保持和水文学方面利用较多[16]。

1.2.4

空间提取算法。除上述需编程实现的算法之外，GIS软件的水文分析功能有直接的坡长提取算法，我们称之为空间分析提取法，其大致原理是如果坡地近似满足坡面水流方向与山脊线垂直，则只要计算出每个点沿垂直方向到山脊线的水平距离可以作为该点坡长的近似值，相对以上3种提取方法而言，该算法较为简单，其大致原理是先采用负地形法提取出研究区域的山脊线，提取出山脊线以后，求取格网单元到最近邻的山脊线的垂直距离作为近似坡长。由于该方法操作简便，计算快捷，又称为快速计算方法。快速计算方法的优点是计算快捷，操作简单。缺点是适用范围比较窄——仅仅局限在等高线与山脊大致平行的地貌（例如川和梁等）或者坡面较均匀、坡度较缓的坡地地貌，计算精度也相对不是很高[6，15]。

综上可知，坡度、坡长和坡向是坡面的3个最基本的地形要素。然而相对于诸如坡度、坡向、沟壑密度、坡面曲率等其他地形要素，人们对坡长的定义仍然有很大争议；有关坡长及坡长因子的研究还不够深入。已有的不同提取算法也有各自的优劣，而对于这些优缺点和方法的适应性并未给出明确的限定和详细的讨论，生成的坡长往往误差较大，精度也不好控制。此外，对于基于GIS的提取算法，在DEM（数字高程模型）高分辨率（<20 m）条件下尚无研究，因DEM分辨率不同而造成的提取结果精度的差异。以足够大比例尺地形图为源数据并数字化，利用功能强大的GIS软件可以很方便地提取到足够精确的DEM。目前已知从1 ∶10 000比例尺地形图中提取得到的DEM的理想分辨率为5 m[5-7]。在源数据精度得以保证的情况下，DEM的分辨率越高，与真实地形的吻合程度也越高，相应地形要素的提取效果也越好。但是利用程序从高分辨率DEM中提取坡长意味着高耗时、低效率，这一缺陷在处理海量数据的时候体现得尤为明显。在DEM高分辨率条件下，其高低对坡长的提取和坡长因子的计算精度有无显著性影响目前尚无定论。

2 研究区概况与方法

2.1 研究区概况

选择的试验区位于陕西省延安市安塞县境内（图1），该县地处陕北黄土高原核心区域，鄂尔多斯盆地南缘，黄河一级支流的延河流域，县域面积2 950 km2。区域地貌類型复杂多样，沟壑密度达4.7 km/km2，平均海拔1 371.9 m，整个区域均被黄土覆盖，厚度大，基本为全新世黄土。气候属中温带大陆性半干旱季风气候，年均气温9.0 ℃，年均降水量约500 mm。水土流失严重，长期观测数据表明，该区域年平均侵蚀输沙模数高达9 370 t/（km2·a），自然条件恶劣，生态环境恶化。

2.2 数据来源

数据主要涉及从1 ∶10 000大比例尺地形图上进行试验小区的选择和DEM的生成。试验小区选择的好坏直接关系到研究结果的可靠性。基于科学性原则、代表性原则和信息全面原则，同时考虑该研究的实际情况，选取典型流域小区共10块，这些地块的覆盖面积介于0.5～1.0 km2，以这些小区的1 ∶10 000比例尺的矢量地形图作为源数据（表1）。以1 ∶10 000地形图为源数据所生成的DEM的理想分辨率为5 m，是分辨率的上限。杨勤科等[17]在榆林、延安和秦岭地区的3个样区基于10、50 m分辨率DEM提取的坡长表明，低分辨率下的提取结果是高分辨率下提取结果的近4倍，由此可见，对于坡长提取尤其是小流域的坡长提取来说，50 m分辨率设置实在太低，提取结果必然误差很大。为探讨DEM分辨率差异带来的计算结果的差异，该研究以5 m分辨率为上限，20 m分辨率为下限，每隔5 m选取一个分辨率值，共选择4组。对每一块试验小区均生成4组由高到低不同分辨率的DEM，然后利用算法提取坡长，分析结果。

2.3 研究方法

利用2种不同的既有坡长提取算法，以若干具有代表性的小流域试验区的不同分辨率DEM为源数据，计算相应小区的相应分辨率的坡长和坡长因子，然后采用数理统计分析法进行结果数据处理，得出相应指标，确定DEM分辨率与提取精度的关系，并分析产生误差的原因。考虑到算法实现的复杂性和论文篇幅限制，采用非累计流量的直接计算方法（noncumulative slope length，NCSL）和空间分析提取法（快速计算方法）提取坡长。如果DEM分辨率高低对提取结果的精度没有显著性影响，则在保证精度的情况下可以通过降低DEM分辨率来提高计算效率，这在处理海量数据时是有一定意义的。采取的统计方法主要是常规数理统计分析法。主要包括以下4种：①绝对误差和相对误差：提取值和标准坡长值的差值的绝对值为绝对误差；绝对误差与标准坡长值的比值乘以100%为相对误差。②合格率：设30%为相对误差上限，当提取得到的坡长值与标准值的误差小于30%为合格点，反之为不合格。以合格点数目与总数的比值乘以100%即得到合格率。③相关系数。④确定性系数：该系数由Nash和Sutcliffe于1970年提出，用以分析提取值与实际值的相关程度。该系数反映的是计算值与实测值1 ∶1线的接近程度，其取值范围是（-∞，1）。当确定性系数大于0时，表明提取值与标准值相当吻合；当确定性系数越接近1，提取值越精确；当确定性系数为0时，表示计算结果与实测值之差的平方和恰好等于实测值与其均值的平方和，表明用实测值的平均值来估算与计算值的效果一样。此外，采用单因素方差分析法以确定DEM分辨率是否对坡长因子计算结果有显著性影响。

3 结果与分析

3.1 2种提取方法比较

由于在自动选择100 m见方的中心样点时有很大的随机性，故可能会出现样点落在山脊线或山顶上，以及有些点落在沟谷的底部的情况。首先应该对参与分析的样点进行筛选，剔除那些不宜与坡长真值作比较的样点，最终确定参与比较的样点的数目为218个（原始的样点数目为292个）。统计结果（表2、3）显示，2种算法的提取精度均随着DEM分辨率的降低而降低，说明DEM分辨率越高，坡长的提取效果越好；非累计流量的直接计算法的坡长提取效果远好于空间分析提取法。

3.2 DEM分辨率对坡长因子的影响

通过单因素方差分析确定DEM分辨率对坡长因子的计算结果是否有显著影响。当DEM分辨率为5、10 m时，非累计流量的直接计算法提取的坡长精度较为理想，均为80%以上，而且二者精度相差不明显；但当DEM分辨率为15、20 m时，提取到的精度已经不够理想。坡长决定坡长因子，坡长因子的计算精度依赖于坡长的提取精度，故当DEM的分辨率大于10 m时，坡长因子的提取精度必然也不高。因此，只需要确定在DEM分辨率为5、10 m这2种条件下对坡长因子的提取有无显著影响。由于在选择样区和样点时均考虑随机性和典型性，可以将参与统计分析的218个样点的坡长因子值看作是来自正态总体的一个样本。考虑DEM分辨率分别为5、10 m的2个水平，取显著性水平α=0.05，检验在这2种水平下坡长因子的均值有无显著性差异，进而可以確定当DEM的分辨率

高于10 m时分辨率对坡长因子的计算效果有无显著影响。

建立给定问题的原假设和备择假设：H0：μ1=μ2，H1：μ1≠μ2（μ1、μ2分别为2组数据的期望值）。F（0.02）α（0.05）（表4）。故接受H0，认为2组坡长计算结果没有显著性差异。经过上述分析可以得出，从5 m分辨率的DEM和10 m分辨率的DEM中提取到的坡长因子没有明显差别，也就是说当DEM分辨率高于10 m時分辨率对坡长因子的计算效果没有显著性影响。

4 结论与讨论

4.1 结论

（1）非累计流量的直接计算法提取出的坡长的效果整体上远好于空间分析法的提取效果。之所以如此，是因为非累计流量算法有明确的物理意义，其设计思路有严密的水文地貌学依据，并且依照分段坡的计算方法计算坡长，与真实的地形特征较为接近。此外，当DEM的分辨率足够大（高于10 m）时，非累计流量直接计算法的提取精度已经相当理想。所以，非累计流量直接计算法已经足够成熟，能够在坡长提取中推广应用。

（2）当DEM的分辨率分别为5、10 m时，利用非累计流量的直接计算法提取的坡长和坡长因子的精度均没有显著性差异。可见，在提取坡长和坡长因子时，一般情况下将DEM的分辨率设置为10 m即可，没有必要设置为5 m。降低DEM的分辨率很大程度上意味着提高程序的计算效率。

（3）空间分析计算法的整体精度较差。空间分析法提取精确的前提条件是山脊线与等高线基本平行，这显然与一般的实际地形不符甚至差别很大；还有该算法在提取山脊线时阈值的设置没有公认的可靠标准，存在盲目性和随意性，极易导致山脊线的误判和漏提， 故造成计算精度很低，尤其是在以低分辨率DEM为源数据时计算结果误差极大，无实际价值，表明该算法在坡长提取方面的应用尚不成熟，亟待改进。

4.2 讨论

首先，受时间等客观条件的制约，试验小区的选择数量偏少。要想得到更为可靠和令人信服的结论，按照以往研究的惯例，仅安塞县境内就需要选择试验小区近100个。在试验小区数量少的情况下，很难反映出该区域的地域分异规律。其次，由于该研究水平十分有限，仅从已有的4种常用算法中选择了2种较为简单的进行提取效果的对比，未能开阔思路，发掘出精度和效率更高的、更适合地理学坡长提取的算法。此外，没有系统地分析误差，未能给出造成计算误差的原因，更没有提出减少和校正误差的可行的办法。这些都是很值得进一步研究的内容。最后，算法的精度评价是一个精细而又复杂的过程，尤其对坡长这一原本就有争议的地形要素来说，其精度评价必然存在着很多尚未解决的问题，涉及这些问题的难度远远超出了该研究的实际水平。此外，整个研究中从地形图的数字化到坡长真值的量取，均不可避免地存在误差；设置的评价标准也较为单一。在这样的条件下得出的结论是否真正可靠有待进一步研究。

参考文献

[1] WISCHMEIER W H.A rainfall erosion index for a universal soil-loss equation[J].Soil science society of America journal，1959，23（3）：246-249.

[2] ONORI F，BONIS P D，GRAUSO S.Soil erosion prediction at the basin scale using the revised universal soil loss equation （RUSLE） in a catchment of Sicily （southern Italy）[J].Environmental geology，2006，50（8）：1129-1140.

[3] FOSTER G R，MCCOOL D K，RENARD K G，et al.Conversion of the universal soil loss equation to SI metric units[J].Journal of soil & water conservation，1981，36（6）：355-359.

[4] 刘宝元，谢云，张科利，等.土壤侵蚀预报模型[M].北京：中国科学技术出版社，2001.

[5] 汤国安，陈正江，赵牡丹，等.Arc View地理信息系统空间分析方法[M].北京：科学出版社，2002.

[6] 汤国安，刘学军，房亮，等.DEM及数字地形分析中尺度问题研究综述[J].武汉大学学报（信息科学版），2006，31（12）：1059-1066.

[7] 汤国安，杨昕.ArcGIS地理信息系统空间分析实验教程[M].北京：科学出版社，2006.

[8] DESMET P J J，GOVERS G.GISbased simulation of erosion and deposition patterns in an agricultural landscape：A comparison of model results with soil map information[J].Catena，1995，25（1/2/3/4）：389-401.

[9] KINNELL P I A.The effect of slope length on sediment concentrations associated with side-slope erosion[J].Soil science society of America journal，2000，64（3）：1004-1008.

[10] 李俊.基于DEM的黄土高原坡长的自动提取和分析[D].西安：西北大学，2007.

[11] HOYOS N.Spatial modeling of soil erosion potential in a tropical watershed of the Colombian Andes [J].Catena，2005，63（1）：85-108.

[12] DE VENTE J，POESEN J，VERSTRAETEN G，et al.Spatially distributed modelling of soil erosion and sediment yield at regional scales in Spain [J].Global and planetary change，2008，60（3/4）：393-415.

[13] 周貴云，刘瑜，邬伦.基于数字高程模型的水系提取算法[J].地理学与国土研究，2000，16（4）：77-81.

[14] 蔡强国，刘纪根.关于我国土壤侵蚀模型研究进展[J].地理科学进展，2003，22（3）：242-250.

[15] 曹龙熹，符素华.基于DEM的坡长计算方法比较分析[J].水土保持通报，2007，27（5）：58-62.

[16] 张启旺，安俊珍，王霞，等.中国土壤侵蚀相关模型研究进展[J].中国水土保持，2014（1）：43-46.

[17] 杨勤科，赵牡丹，刘咏梅，等.DEM与区域土壤侵蚀地形因子研究[J].地理信息世界，2009，7（1）：25-31.