基于CA-Markov模型与ANUDEM内插法的崩岗侵蚀量预估

季 翔，黄炎和，林金石，蒋芳市，余明明，李绍鑫



基于CA-Markov模型与ANUDEM内插法的崩岗侵蚀量预估

季 翔1,2，黄炎和1※，林金石1，蒋芳市1，余明明1，李绍鑫1

（1. 福建农林大学资源与环境学院，福州 350002；2. 福建省土壤环境健康与调控重点实验室，福州 350002）

崩岗是中国南方最为严重的土壤侵蚀类型之一，产生的大量泥沙危害农业生产和生态环境，因此对其侵蚀量的预估是防治该现象的重要途径。崩岗面积较小且侵蚀剧烈，难以应用现有方法预估侵蚀量。该文应用CA-Markov模型和ANUDEM内插法对其高程级别模拟和空间内插，从而实现对崩岗侵蚀量的预估，并以福建省安溪县龙门镇的一处崩岗为例进行实证研究。结果表明：CA-Markov模型适用于对崩岗高程级别的模拟；ANUDEM内插法对崩岗地形的整体还原度较好，但对细节的刻画不够；以经过级别划分和内插处理的高程数据为基期底图计算得的崩岗侵蚀量较符合实际值，且实际侵蚀量越大，模拟精度越高；案例崩岗在一般年景、干旱年景和多雨年景中的年侵蚀量分别为：824.69、731.03和 924.57 m3，不同年景之间侵蚀量的最大差值为193.54 m3，因此在修建崩岗拦沙坝时需考虑不同降雨年景中侵蚀量的差异。研究结果不仅提供了预估崩岗侵蚀量的新思路，还可为崩岗侵蚀的防治工作提供参考依据。

侵蚀；模型；崩岗；侵蚀量；CA-Markov；ANUDEM

0 引 言

崩岗是指在水力和重力的作用下，山坡土石体受破坏而崩塌和受冲刷的侵蚀现象，主要集中在中国南方花岗岩发育的红壤丘陵区。根据相关研究[1-3]，南方红壤丘陵区共有各类崩岗239 125个，总侵蚀面积为1 220.05 km2，主要集中在江西、广东、湖南、福建、湖北、安徽、广西七省（自治区）。崩岗侵蚀产生的泥沙量巨大，不仅造成表土流失无法利用，而且泥沙被带到下游埋没农田、淤积河道水库，给农业生产、生态环境以及人民生活带来极大危害[4-5]。因此，预测和估算崩岗侵蚀量对减缓危害以及防治水土流失有重要作用。

土壤侵蚀量的估算方法众多，主要是根据通用土壤侵蚀方程（universal soil loss equation, USLE）计算土壤侵蚀模数[6-7]或通过流域模型估算土壤侵蚀量[8-9]。然而由于崩岗相对其他类型的土壤侵蚀空间尺度较小，且侵蚀剧烈地表变化迅猛[10-11]，因此较难应用该方法预估崩岗的侵蚀量。崩岗的侵蚀量可由其地形体积的变化表示[12-13]，如刘希林等曾利用3D激光扫描仪对广东五华县莲塘岗崩岗的演化过程进行空间监测，获得了精确的崩岗侵蚀 量[13]。因此，可通过对崩岗空间的动态模拟计算其地形变化体积，从而实现对崩岗侵蚀量的预估。由于崩岗侵蚀过程的复杂性限制了机理型模型的应用[14-17]，经验型模型更为适合。目前应用普遍的空间模拟模型是由马尔科夫链（Markov）和元胞自动机（Cellular automatic，CA）耦合的CA-Markov模型[18-20]。该模型适用于离散的类型数据，对于连续变量需要进行级别划分[21-22]。因此，应用该模型对崩岗地形变化进行模拟得到的结果为离散的高程级别，难以实现对崩岗侵蚀量的精确预估。

高程级别是将一定区间的高程值划分为一个级别，相邻级别边界线上的高程取值相同，这与等高线类似。因此，可通过空间内插将离散的高程级别转化成连续的高程数据，由此实现对崩岗侵蚀量的估算。空间内插法是一种将空间离散数据连续化的工具，在空间数据方面的应用较为成熟[23-25]，其中ANUDEM法在高程内插中的精度较高，应用成熟[24]。崩岗地表沟壑横生，若单纯根据等高线进行空间内插，难以体现地表细节。ANUDEM不仅可将等高线内插为精度较高的DEM数据，并且该方法可以根据地表径流、湖泊等辅助信息强化地表细节，较为适用于对崩岗高程的空间内插。

本研究按照“连续-离散-连续”的思路，将连续的崩岗高程进行级别划分使之适用于CA-Markov模型，并将预测得的崩岗高程级别进行ANUDEM空间内插，得到连续的崩岗高程预测值。通过与基期高程叠加得到土体体积变化量，从而实现对崩岗侵蚀量的估算。研究结果不仅提供了预估崩岗侵蚀量的新思路，还可为有针对性地防治崩岗侵蚀、减少河流泥沙提供参考依据。

1 研究方法

1.1 CA-Markov模型

1.1.1 Markov链

Markov链是基于Markov过程理论而形成的预测事件发生概率的一种方法，常用于具有无后效性特征地理事件的预测。一定区域内，不同空间类型具有相互转化的可能，各类型之间的转化过程有一些难以用函数关系准确描述的事件。在崩岗侵蚀过程中不同高程级别对应Markov过程中的“可能状态”，而不同高程级别间相互转换的面积数量或比例即为状态转移概率，可以利用如下公式对各高程级别的变化进行预测：

式中S1和S分别表示前一时刻和后一时刻的高程级别；m表示级高程转变为级高程的概率，%，且各行之和为1；为高程级别的数目。

m是由前后两期数据叠加而得，由于两期数据时间间隔难以控制，因此需要将转移概率进行转换，使其符合预测的时间长度，这在模拟相对匀速的空间变化时效果较好[20]。然而，崩岗侵蚀剧烈且在时间尺度上差异较大[26-27]，若采用时间长度对其转移概率进行转换，难以反映崩岗侵蚀特征，模拟结果精度不高。崩岗是受水力和重力共同作用的复合侵蚀，水力作用主要来自于降雨后的地表径流冲刷，重力作用同样由于降雨增加土体重力导致，由此可知降雨量对崩岗侵蚀起着决定性的作用[28-30]。因此，本研究采用降雨量对其转移概率进行转换，公式如下：

根据转移概率矩阵可以计算得下一时刻每种类型的面积，从而对CA的空间分配进行数量约束。

1.1.2 CA-Markov耦合模型

CA模型是一个空间分配的自动机制，其对变量没有要求，可以是离散数据也可以是连续数据，在土壤侵蚀方面得到过较多应用。该方法是通过将栅格图中每栅格视为一个元胞，每个栅格值即为元胞状态，根据转换规则计算元胞的下一个状态。CA-Markov模型是将Markov链的经验性与CA的空间分配结合起来，以崩岗的地形级别为基础数据，具体参数设置如下：

1）元胞：元胞大小即为栅格图像大小，本研究设定栅格大小为0.1 m×0.1 m。

2）元胞状态：在标准CA模型中，元胞状态为一个有限的离散集合。因此在对崩岗地形的模拟中，需要将连续的高程数据进行级别划分，高程级别即为元胞状态。

3）邻域：元胞的下一个状态由元胞自身及其邻域决定，需要计算邻域对该元胞状态转变的影响。由于崩岗高程为连续型变化，因此本研究采用元胞到目标状态最近距离表征邻域影响，距离越小则邻域影响越大，计算公式为：

式中N为某个级元胞转变为级的邻域影响；D为该级元胞到级的最近直线距离，m；将所有级元胞的到级最近直线距离组成集合，集合的最大值为max(D)，由此可将N的取值转化为0~1区间内。由于崩岗侵蚀的趋势为高程下降，越接近下一级别，越难保持不变，因此保持高程级别不变的邻域影响为到下一级别高程的 逆值。

4）适宜性图集：适宜性和邻域影响共同控制元胞状态的转变，由影响元胞转变的驱动因素计算而得，计算公式如下：

式中为适宜性，d和ω分别为第个驱动因素取值与该驱动因素的权重。

如前所述，崩岗受降雨后的地表径流冲刷和重力的双重作用。其中，崩岗的地表径流即为沟道，重力则主要作用在坡度陡峭的土体上[29-32]。因此，本研究选择沟道距离和坡度为表征崩岗侵蚀力的驱动因素。其中，降雨为崩岗侵蚀的关键驱动力，其作用远大与土体自身的重力[30-32]，结合前期研究将这两个驱动因素的权重设为0.75和0.25。

5）转换规则：转换规则是元胞自动机模型的核心，决定了元胞分配的动态过程。按照每个元胞的邻域影响与适宜性之和，从最高值开始赋值，然后次高，直到满足预测的数量。

1.2 ANUDEM内插

ANUDEM在对高程点或等高线的内插过程中，加入了如悬崖、河流等其他地形细节的数据作为约束，使内插后地形结构连续且地形细节准确，因此通过此方法能够构建出体现水文地貌关系的DEM[20]。该方法的输入数据包括基础数据和约束数据两部分。其中，基础数据为：等高线、高程点；约束数据为：洼地点、地表径流、湖泊、悬崖、海岸、边界以及不需要计算的区域。

崩岗的高程级别为面状数据，因此在进行ANUDEM时需要将其转换为线文件，然后根据高程级别的分级值对各边界线赋值，从而得到等高线数据，即ANUDEM内插的基础数据。例如，级高程同时与-1级和+1级相邻，将和-1级共用的边界线赋值为级和-1的分级值，将和+1级共用的边界线赋值为级和+1的分级值。崩岗侵蚀区的沟道对其土壤侵蚀以及崩岗地貌改变起着关键性作用，并且其空间位置较为稳定。因此，选择沟道（地表径流）作为应用ANUDEM对崩岗高程空间内插的约束数据。

1.3 崩岗侵蚀量

土壤侵蚀量可以通过质量或体积表示，考虑到在计算质量时，对土壤容重的选择会产生误差。因而，本研究采用崩岗的体积变化对崩岗侵蚀量进行表征，通过高程的栅格数据计算，公式如下：

1.4 降雨情景设置

降雨量对崩岗侵蚀起着决定性的作用，因此本研究在构建模拟崩岗高程级别变化的CA-Markov模型时，根据降雨量对转移概率矩阵进行转换。然而，降雨是一种具有不确定性的自然现象，以现有数据难以准确预测未来的降雨量。情景分析法是对未来可能出现情况的分析及预测。因此，本研究采用情景分析法，根据研究区域的历史降雨数据设置降雨情景，预测不同降雨情景下 崩岗地形的变化，从而估算出不同降雨情景时崩岗的侵蚀量。

2 结果与分析

2.1 案例点及数据

2.1.1 案例点选取

本研究的案例点位于福建省泉州市安溪县龙门镇（24°57¢N, 118°3¢E），该区域为亚热带季风气候，年均气温为19 ℃，降雨量高且集中，年降雨量1 546.2～2 023.8 mm，且多在5—9月。该崩岗面积较大，约5 300 m2，虽然曾被横断修路治理过，但由于该崩岗侵蚀剧烈，道路已无法辨别。因此本研究以此崩岗为案例点，如图1。

图1 案例崩岗示意图

2.1.2 数据来源与处理

1）高程及高程级别

分别于2017年3月11日、2017年7月21日和2017年12月2日，由RTK（Trimble R4 GNSS）设置控制点，无人机飞行器（DJI GO 4 APP）对案例崩岗进行低空摄影测量。经过拼接校正等处理，得到3个时期案例崩岗的DEM，其空间分布率为2.8 cm，投影坐标系为WGS-1984- UTM-zone-50N。过高的空间分布率将增加模型运行负担，因此本研究对3期DEM进行了重采样处理将其空间分辨率降为10 cm。

本研究中需要将连续的崩岗高程数据划分级别。一般来说，级别划分的越详细，模拟结果越准确，但数目过多会造成模型运行负担。因此，根据数据特征并结合现有研究成果[32]，将案例崩岗进行了高程级别的划分，如图2。

注：以2.25 m为间隔将案例崩岗划分为20个高程级别，其中1为最低级别高程250~252.25 m，20为最高级别高程292.75~295 m。下同。

2）降雨量

降雨量作为修正转移概率的指标在本研究中非常重要，与崩岗高程数据相应时间（2017年3月11日－2017年12月2日）的降雨量由放置在距离案例点约740 m的翻斗雨量计（Onset rg3-m）获取，其精度为2 mm。

降雨情景需要根据研究区域的历年降雨量记录进行设置，从而使其分析结果具有代表性。而翻斗雨量计测定雨量的时间较短，不符合情景设置的要求。因此，采用福建省安溪县第二次土壤普查成果《安溪土壤》[33]的降雨量记录进行降雨情景的设置。

2.2 精度验证

根据2017年3月11日和2017年7月21日的高程数据以及2017年3月11日－2017年7月21日和2017年7月21日－2017年12月2日的降雨量数据，对2017年12月2日的崩岗高程进行模拟并估算崩岗侵蚀量，然后与基于实测数据计算的崩岗侵蚀量进行对比，完成模拟精度和估算精度的验证。

2.2.1 CA-Markov模型

1）面积检验

首先计算2017年3月11日和2017年7月21日的转移概率矩阵，并根据该时段降雨量和2017年7月21日－2017年12月2日的降雨量进行转换，将2017年3月11日和2017年7月21日的转移概率矩阵转换为2017年7月21日－2017年12月2日的转移概率矩阵，如 表1所示。

表1 2017年7月21日－2017年12月2日的转移概率矩阵

然后根据转移概率矩阵计算2017年12月2日时各高程级别面积的模拟值，并计算与实测高程的各级别面积（图2）之间的误差，面积误差越小说明模拟精度越高，计算公式如下：

式中δ为面积误差，%；as和at分别为高程级别面积的模拟值和实测值，m2。各高程级别的模拟值与实测值的面积误差如图3所示。

图3中，各高程级别面积误差平均为1.47%，除第1级高程的面积误差较大（7.5%）外，其他级别的面积误差均小于5%，说明高程级别的模拟面积精度较好，可以用于下一步的空间模拟。

2）空间检验

借助ArcGIS和Matlab软件获得CA-Markov模型所需参数，并完成2017年12月2日时案例崩岗高程级别的空间模拟，然后计算其Kappa系数。Kappa系数是通过模拟值与实测值之间的误差矩阵计算得的，矩阵列表示实测类，矩阵行表示模拟类，矩阵中的数值表示像元数目，其计算公式如下：

式中KC为Kappa系数；为样本数目；P为第类被正确分类的样本数目；P为第类所在列的像元数目；P为第类所在行的像元数目。当KC≥0.75时，两图的一致性较高，精度较高；当0.4≤KC<0.75时，精度一般；当KC<0.4时，精度较差。

根据公式（8）计算得总Kappa系数为0.958，说明该模型的模拟结果可靠。另外，各高程级别的Kappa系数如图4。

图4 各高程级别的Kappa系数

由图4可以看出各高程级别的Kappa系数均大于0.85，说明各高程级别的模拟精度均可接受。

2.2.2 ANUDEM内插

1）崩岗高程

对岗高程级别的CA-Markov模拟结果进行ANUDEM空间内插，将离散的高程级别转化为连续的高程数据，即崩岗的模拟高程。然后与同期的实测高程进行对比，计算拟合度2和平均误差ME用以检验精度，拟合度2用以验证地形的整体还原度，平均误差ME则用以反映地形细节的一致性，计算公式如下：

根据式（9）和式（10）计算出拟合度2为0.998 7，说明ANUDEM内插后的崩岗地形整体上还原较好，精度可以接受；但平均误差ME为0.27 m，该值表示模拟高程与实测高程之间的平均差值，虽然数值较小，但相对于模拟期的崩岗地形变化较大。由此可以看出ANUDEM可以根据崩岗高程级别很好的还原崩岗的整体地貌形态，但对于崩岗地表细节刻画的准确度一般。

2）崩岗侵蚀量

将模拟高程与基期高程叠加，根据式（6）计算出崩岗侵蚀量的模拟值；同时，将相同时段的无人机实测DEM叠加，根据式（6）计算出崩岗侵蚀量的实际值。由于ANUDEM对基于高程级别还原崩岗整体地貌的精度较高，但对其细节刻画不够理想。因此，为了使基期高程与模拟高程对崩岗地形的还原程度一致，对原始DEM进行级别划分和ANUDEM内插处理。并将原始DEM和处理后的高程数据分别与模拟高程叠加，得崩岗侵蚀量的模拟值Ⅰ和模拟值Ⅱ。然后计算与实际值之间的误差，记为误差Ⅰ和误差Ⅱ。计算结果如表2。

由表2可以看出，2个时段的误差Ⅰ均较大，而2个时段的误差Ⅱ均较小，仅为17.55%和8.66%，这一方面确认了模拟值Ⅱ的精度可以接受，另一方面进一步确认了ANUDEM对还原崩岗整体地形的可靠性，但对其细节的刻画不够。另外，可以看出2017年7月21日－2017年12月2日时段的2个误差均大于2017年3月11日－2017年12月2日时段，这主要是由于后一时段的侵蚀量远大于前一时段，而大量的侵蚀导致崩岗地形的变化较为显著，减小了地形细节对计算侵蚀量的影响，从而提高的模拟精度。综合而言，模拟值Ⅱ的精度较高，因此在估算崩岗侵蚀量时，选择经过级别划分和ANUDEM内插处理的高程为基期数据；并且，崩岗的实际侵蚀量越大，估算精度越高。

表2 崩岗侵蚀量的模拟精度

注：模拟值Ⅰ和误差Ⅰ是以原始DEM为基期数据；模拟值Ⅱ和误差Ⅱ是以级别划分和内插处理后的高程为基期数据。

Note: Simulation Ⅰ and error Ⅰ are based on original DEM; Simulation Ⅱ and error Ⅱ are based on the elevation data processed by classification and interpolation.

2.3 崩岗侵蚀量预测

由于降雨是崩岗侵蚀的关键影响因素，因此本研究根据案例区历年降雨数据设置了一般年景、干旱年景和多雨年景3个情景，分别对应年降雨量的平均值1 785 mm、最小值量1 546.2 mm和最大值2 023.8 mm。然后根据CA-Markov模型和ANUDEM空间内插，估算出这3个情景中案例崩岗的2017年12月2日－2018年12月2日的侵蚀量，即该崩岗的年侵蚀量。

2.3.1 崩岗高程级别预测

首先根据2017年3月11日－2017年12月2日计算出转移概率矩阵，然后根据此时段的降雨量以及3个情景中2017年12月2日－2018年12月2日的降雨量对其进行转换，得到3个情景中2017年12月2日－2018年12月2日的转移概率矩阵。然后，应用CA-Markov模型对3个情景中2018年12月2日时崩岗高程级别进行空间模拟，如图5。

2.3.2 崩岗侵蚀量预测

应用ANUDEM对2018年12月2日时崩岗高程级别预测图进行内插，得到2018年12月2日时崩岗高程预测图，如图4。同时，对2017年12月2日时实测高程进行级别划分与内插处理。然后，分别与3个情景中2018年12月2日时崩岗高程预测图叠加，求得3个情景中，崩岗在2017年12月2日－2018年12月2日这一时段里的侵蚀量（年侵蚀量），如表3。

由表3可以看出，该崩岗在一般年景中的年侵蚀量为824.69 m3，即使在干旱年景时，该崩岗的年侵蚀量也有731.03 m3，说明该崩岗的泥沙产量巨大；该崩岗在2个极端情景之间的侵蚀量差值为193.54 m3，说明了不同降雨年景对崩岗泥沙产量的影响较大，在防治崩岗侵蚀修建拦沙坝时应考虑到不同降雨年景时侵蚀量的差异。

图5 不同情景下崩岗高程级别预测图和崩岗高程预测图（2018-12-02）

表3 3个情景中崩岗2017年12月2日－ 2018年12月2日的侵蚀量

3 讨 论

崩岗侵蚀剧烈，泥沙产量巨大，给当地的农业生产以及生态环境带来威胁。但目前对崩岗侵蚀的研究侧重于事后的治理工作，对崩岗侵蚀的预测、预警方面的研究相对较少。在实际的崩岗侵蚀治理工作中，多以修建堤坝拦截泥沙为主要途径，而这更加大了对预测崩岗侵蚀量的需求。由于崩岗面积较小，且侵蚀剧烈地表变化迅速[10-11]，难以应用现有方法对其侵蚀量进行预估。因此，本研究在应用CA-Markov对崩岗高程级别的进行预测的基础上，应用ANUDEM对高程级别进行空间内插实现对崩岗地形变化的预测，从而实现对崩岗侵蚀量的预估。不仅为预测和估算土壤侵蚀量提供新的思路和方法，并且为崩岗侵蚀的防治工作提供参考依据。

关于崩岗侵蚀量预测方面的研究相对较少，一般通过实测方式获取崩岗侵蚀量，如刘希林等[32]借助三维激光扫描仪测得莲塘岗崩岗（与案例崩岗相似）的3 a平均侵蚀量为833 m3，这与本所得结果非常接近。但由于不同崩岗的岩土及降雨等条件的影响，其侵蚀量也存在一定差异，且同一崩岗在不同时期其侵蚀量也不同。

CA-Markov为经验型模型，根据上一时段变化预测下一时段变化，即将模拟的生态过程视为匀速过程。而崩岗侵蚀过程的驱动力主要来自于在时间尺度上较为多变的降雨，为此本研究选择降雨量取代时间步长修正转移概率矩阵，且模拟精度较高。同时由于该模型的经验性，使经验期和模拟期中降雨量对崩岗侵蚀的作用机制一致，因而难以探讨不同降雨事件对崩岗侵蚀的作用。

有学者认为降雨强度对崩岗侵蚀也存在较为显著[34-35]的作用，如Jiang等[34]通过室内试验发现降雨强度与侵蚀泥沙颗粒的平均重量直径成正比。但降雨量是降雨强度与降雨时长的综合结果，与降雨强度相比，降雨量更适用于估算土壤侵蚀量[36]。但将降雨强度纳入该模型中将会提高模拟精度，因而这方面还有待进一步研究。

另外，崩岗地表土壤的理化性状以及如结皮等土壤方面的因素[37-38]，对崩岗侵蚀的空间差异性产生影响，但由于数据获取困难，本研究在进行崩岗侵蚀量预估时未能考虑该因素。

4 结 论

1）根据崩岗侵蚀特征，构建了适用于模拟崩岗地形变化的CA-Markov模型。并通过对案例研究对其进行精度验证：各高程级别面积误差的平均值为1.47%，空间分布的总Kappa系数为0.958，各高程级别的Kappa系数均大于0.85，由此验证了该模型对模拟崩岗高程级别变化的可靠性。

2）应用ANUDEM对崩岗高程级别进行内插，得到连续的崩岗模拟高程。通过与实测数据对比得到拟合度和平均误差，分别为0.998 7和0.27 m，说明该方法虽然对崩岗地表细节的还原度一般，但对崩岗整体地形的还原度很高。

3）选择了2种基期高程计算崩岗侵蚀量：原始DEM和经过级别划分和内插处理的高程数据。通过与基于实测高程计算的侵蚀量对比发现，采用后者为基期数据计算得的崩岗侵蚀量更接近实际值，误差为17.55%和8.66%；且实际侵蚀量越大，崩岗侵蚀量的估算精度更高。这主要是由于侵蚀越剧烈，地形改变越显著，对还原地表细节的需要越低。

4）由于降雨量较难以预测，本研究根据案例点所在区域的历史降雨数据设置了3种情景：一般年景、干旱年景和多雨年景。并对这3个情景中案例崩岗在2017年12月2日－2018年12月2日的侵蚀量（年侵蚀量）进行了预估，分别为824.69、731.03和924.57 m3。其中2个极端情景的侵蚀量之差较大，为193.54 m3，由此说明在修建崩岗拦沙坝时需考虑不同降雨年景时侵蚀量的差异。

[1] 水利部，中国科学院，中国工程院. 中国水土流失防治与生态安全[M]. 北京：科学出版社，2011：924－925.

[2] 梁音，宁堆虎，潘贤章，等. 南方红壤区崩岗侵蚀的特点与治理[J]. 中国水土保持，2009 (1)：31－34. Liang Yin, Ning Duihu, Pan Xianzhang, et al. Characters and control of collapsing hill erosion in South China[J]. Soil and Water Conservation in China, 2009 (1): 31－34. (in Chinese with English abstract)

[3] 牛德奎. 华南红壤丘陵区崩岗发育的环境背景与侵蚀机理研究[D]. 南京：南京林业大学，2009. Niu Dekui. Research on the Environmental Factors and Erosive Mechanism of Collapsing Hill in South China[D]. Nanjing : Nanjing Forestry University, 2009. (in Chinese with English abstract)

[4] 王礼先，孙保平，余新晓，等. 中国水利百科全书——水土保持分册[M]. 北京：中国水利水电出版社，2004：48－49.

[5] 王学强，蔡强国. 崩岗及其治理措施的系统分析[J]. 中国水土保持，2007(7)：29－31. Wang Xueqiang, Cai Qiangguo. Systematic analysis on dilapidated granite and measures for management[J]. Soil and Water Conservation in China, 2007(7): 29－31. (in Chinese with English abstract)

[6] 赵明松，李德成，张甘霖，等. 基于RUSLE模型的安徽省土壤侵蚀及其养分流失评估[J]. 土壤学报，2016，53(1)：28－38. Zhao Mingsong, Li Decheng, Zhang Ganlin, et al. Evaluation of soil erosion and soil nutrient loss in Anhui Province based on RUSLE Model[J]. Acta Pedologica Sinica, 2016, 53(1): 28－38. (in Chinese with English abstract)

[7] Kazuya Yoshimura, Yuichi Onda, Hiroaki Kato. Evaluation of radiocaesium wash-off by soil erosion from various land uses using USLE plots[J]. Journal of Environmental Radioactivity, 2015, 139: 362－369.

[8] 贾宁凤，段建南，李保国，等. 基于 AnnAGNPS模型的黄土高原小流域土壤侵蚀定量评价[J]. 农业工程学报，2006，22(12)：23－27. Jia Ningfeng, Duan Jiannan, Li Baoguo, et al. Soil erosion quantitative evaluation of small watershed in Loess Plateau based on AnnAGNPS model[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2006, 22(12): 23－27. (in Chinese with English abstract)

[9] 于国强，李占斌，李鹏，等. 黄土高原小流域重力侵蚀数值模拟[J]. 农业工程学报，2009，25(12)：74－79. Yu Guoqiang, Li Zhanbin, Li Peng, et al. Numerical simulation on gravitational erosion of small watershed system in Loess Plateau[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2009, 25(12): 74－79. (in Chinese with English abstract)

[10] 阮伏水. 福建省崩岗侵蚀与治理模式探讨[J]. 山地学报，2003，21(6)：675－680. Ruan Fushui. Study on slump gully erosion and its control in Fujian Province[J]. Journal of Mountain Science, 2003, 21(6): 675－680. (in Chinese with English abstract)

[11] 林敬兰，陈志明，黄炎和，等. 安溪县崩岗侵蚀空间分布特征探讨[J]. 水土保持研究，2009，24(6)：63－68. Lin Jinglan, Chen Zhiming, Huang Yanhe, et al. Study on the charactistics of spatial distribution of slope disintegration erosion in Anxi County[J]. Research of Soil and Water Conservation, 2009, 24(6): 63－68. (in Chinese with English abstract)

[12] Marco Cavalli, Barbara Goldin, Francesco Comiti, et al. Assessment of erosion and deposition in steep mountain basins by differencing sequential digital terrain models[J]. Geomorphology, 2017, 291：4－16.

[13] 刘希林，张大林. 基于三维激光扫描的崩岗侵蚀的时空分析[J]. 农业工程学报，2015，31(4): 204－211. Liu Xilin, Zhang Dalin. Temporal－spatial analyses of collapsed gully erosion based on three－dimensional laser scanning[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2015, 31(4)：204－211. (in Chinese with English abstract)

[14] Aliakbar Nazari Samani, Hassan Ahmadi, Mohammad Jafari, et al. Geomorphic threshold conditions for gully erosion in Southwestern Iran (Boushehr－Samal watershed)[J]. Journal of Asian Earth Sciences, 2009, 35(2): 180－189.

[15] Kakembo Vincent, Xanga WW, Rowntree Katharine. Topographic thresholds in gully development on the hill slopes of communal areas in Ngqushwa Local Municipality, Eastern Cape, South Africa[J]. Geomorphology, 2009, 110 (3/4): 188－194.

[16] Antonina Capra, Paolo Porto, Baldassare Scicolone, et al. Relationships between rainfall characteristics and ephemeral gully erosion in a cultivated catchment in Sicily (Italy)[J]. Soil and Tillage Research, 2009, 105(1): 77－87.

[17] 蒋芳市. 花岗岩崩岗崩积体侵蚀机理研究[D]. 福州：福建农林大学，2013. Jiang Fangshi. Erosion Mechanism of Colluvial Deposits in Granite Benggang[D]. Fuzhou: Fujian Agriculture and Forestry University, 2013. (in Chinese with English abstract)

[18] Kuswantoro Marko, Faris Zulkarnain, Eko Kusratmoko. Coupling of Markov chains and cellular automata spatial models to predict land cover changes (case study: upper Ci Leungsi catchment area)[C]// IOP Conference Series: Earth and Environmental Science. Bristol, England: IOP Publishing Ltd, 2017: 012032.

[19] Zhou De, Lin Zhulu, Liu Liming. Regional land salinization assessment and simulation through cellular automaton－markov modeling and spatial pattern analysis[J]. Science of the Total Environment, 2012, 439: 260－274.

[20] 季翔，刘黎明，李洪庆. 基于生命周期的乡村景观格局演变的预测方法——以湖南省金井镇为例[J]. 应用生态学报，2014，25(11)：3270－3278. Ji Xiang, Liu Liming, Li Hongqing. Prediction method of rural landscape pattern evolution based on life cycle: A case study of Jinjing Town, Hunan Province, China[J]. Chinese Journal of Applied Ecology, 2014, 25(11): 3270－3278. (in Chinese with English abstract)

[21] Liu Ruimin, Men Cong, Wang Xiujuan, et al. Application of spatial Markov chains to the analysis of the temporal-spatial evolution of soil erosion[J]. Water Sci Technol, 2016, 74(5): 1051－1059.

[22] Rhonda Skaggs, Soumen Ghosh. Assessing changes in soil erosion rates: A Markov chain analysis[J]. Journal of Agricultural and Applied Economics, 2009, 31(3): 611－622.

[23] 王俊，江帆，朱长青，等. 利用等高线栅格化的DEM精度评估[J]. 武汉大学学报·信息科学版，2009，34(8)： 952－955. Wang Jun, Jiang Fan, Zhu Changqing, et al. Estimate the accuracy of DEM based on reconstructed contours[J]. Geomatics and Information Science of Wuhan University, 2009, 34(8): 952－955. (in Chinese with English abstract)

[24] 张笑，陈刚，段淼然，等. 南方复合地貌区DEM内插算法的对比研究[J]. 水土保持研究，2015，22(5)：106－112. Zhang Xiao, Chen Gang, Duan Miaoran, et al. Comparative study of different DEM interpolation algorithms in the southern complex geomorphologic zone[J]. Research of Soil and Water Conservation, 2015, 22(5): 106－112. (in Chinese with English abstract)

[25] 徐亚莉，罗明良，梁倍瑜，等. DEM 空间插值方法对土壤侵蚀模拟的影响研究：以USPED分析干热河谷典型冲沟为例[J]. 地理科学进展，2016，35(7)：870－877. Xu Yali, Luo Mingliang, Liang Beiyu, et al. Effects of different DEM spatial interpolation methods on soil erosion simulation: A case study of a typical gully of dry－hot valley based on USPED[J]. Progress in Geography, 2016, 35(7): 870－877. (in Chinese with English abstract)

[26] 刘希林，张大林. 崩岗土体剖面水分分布特征及变化规 律：以广东五华县莲塘岗崩岗为例[J]. 热带地理，2015，35(3)：291－297. Liu Xilin, Zhang Dalin. Distribution characteristics and spatial variation of Benggang soil moistures: A case study of Liantanggang in Wuhua County, Guangdong[J]. Tropical Geography, 2015, 35(3): 291－297. (in Chinese with English abstract)

[27] 史德明. 我国热带、亚热带地区崩岗侵蚀的剖析[J]. 水土保持学报，1992，6(4)：32－37. Shi Deming. Analysis of collapsing erosion in tropical and subtropical regions of China[J]. Journal of Soil and Water Conservation, 1992, 6(4): 32－37. (in Chinese with English abstract)

[28] 林敬兰，黄炎和，张德斌，等. 水分对崩岗土体抗剪切特性的影响[J]. 水土保持学报，2013，27(3)：55－58. Lin Jinglan, Huang Yanhe, Zhang Debin, et al. Influence of soil moisture content on shear characteristics of Benggang[J]. Journal of Soil and Water Conservation, 2013, 27(3): 55－58. (in Chinese with English abstract)

[29] Vincent Chaplot, James Brown, Phumlani Dlamini, et al. Rainfall simulation to identify the storm-scale mechanisms of gully bank retreat[J]. Agricultural Water Management, 2011, 98: 1704－1710.

[30] 蒋芳市，黄炎和，林金石，等. 坡度和雨强对崩岗崩积体侵蚀泥沙颗粒特征的影响[J]. 土壤学报，2014，51(5)： 974－982. Jiang Fangshi, Huang Yanhe, Lin Jinshi, et al. Effects of slope gradient and rainfall intensity on particle size composition of erosion sediment form conlluvial deposits of Benggang[J]. Acta Pedologica Sinca. 2014, 51(5): 974－982. (in Chinese with English abstract)

[31] 刘希林，张大林，贾瑶瑶. 崩岗地貌发育的土体物理性质及其土壤侵蚀意义：以广东五华县莲塘岗崩岗为例[J]. 地球科学进展，2013，28(7)：802－811. Liu Xilin, Zhang Dalin, Jia Yaoyao. Soil physical properties of collapsing hill and gully and their indications for soil erosion: An example of Liantanggang collapsing hill and gully in Wuhua county of Guangdong[J]. Advances in Earth Science, 2013, 28(7): 802－811. (in Chinese with English abstract)

[32] 刘希林，张大林. 崩岗地貌侵蚀过程三维立体监测研 究：以广东五华县莲塘岗崩岗为例[J]. 水土保持学报，2015，29(1)：26－31. Liu Xilin, Zhang Dalin. Study on erosion process of collapsing hill and gully by Three-dimensional Monitoring: An example of Liantanggang in Wuhua County of Guangdong Province[J]. Journal of Soil and Water Conservation, 2015, 29(1): 26－31. (in Chinese with English abstract)

[33] 福建省安溪县土壤普查办公室. 安溪土壤[R]. 1984.

[34] Jiang Fangshi, Huang Yanhe, Wang Mingkuang, et al. Effects of rainfall intensity and slope gradient on steep colluvial deposit erosion in Southeast China[J]. Soil Science Society of America Journal, 2014, 78(5): 1741－1752.

[35] 熊传祥，王涛，鲁晓兵. 降雨作用下崩岗形成细观机理模拟[J]. 山地学报，2013，31(6)：710－715. Xiong Chuanxiang, Wang Tao, Lu Xiaobing. Meso-mechanical simulation of slope disintegration erosion under rainfall[J]. Journal of Mountain Science, 2013, 31(6): 710－715. (in Chinese with English abstract)

[36] 彭双云，杨昆，洪亮，等. 基于 USLE 模型的滇池流域土壤侵蚀时空演变分析[J]. 农业工程学报，2018，34(10)：138－146. Peng Shuangyun, Yang Kun, Hong Liang, et al. Spatio- temporal evolution analysis of soil erosion based on USLE model in Dianchi Basin[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2018, 34(10): 138－146. (in Chinese with English abstract)

[37] 程琴娟，蔡强国，李家永. 表土结皮发育过程及其侵蚀响应研究进展[J]. 地理科学进展，2005，24(4)：114－122. Cheng Qinjuan, Cai Qiangguo, Li Jiayong. Summarization on study of soil surface crust or sealing and its effects on erosion[J]. Process in Geography, 2005, 24(4): 114－122. (in Chinese with English abstract)

[38] 蒋芳市，黄炎和，林金石，等. 花岗岩崩岗崩积体颗粒组成及分形特征[J]. 水土保持研究，2014，21(6)：175－180. Jiang Fangshi, Huang Yanhe, Lin Jinshi, et al. Soil particle size distribution and fractal dimensions of colluvial deposits in granite Benggang[J]. Research of Soil and Water Conservation, 2014, 21(6):175－180. (in Chinese with English abstract)

Estimation of erosion amount in collapsed gully based on CA-Markov model and ANUDEM interpolation

Ji Xiang1,2, Huang Yanhe1※, Lin Jinshi1, Jiang Fangshi1, Yu Mingming1, Li Shaoxin1

(1.350002,; 2.350002,)

Collapsed gully is a kind of soil erosion in South China threating agriculture production and ecological environment due to its considerable sediment yield. Therefore, the prediction and estimation of erosion amount play a significant role in its prevention and control. However, collapsed gullies normally cover a relatively smaller area and eroded more severely than other kinds of soil erosion, which made existing methods and approaches hardly to apply on it. In this research, cellular automation-Markov (CA-Markov) model and Australian National University digital elevation model (ANUDEM) are employed to estimate and predict the erosion amount of collapsed gully. The former is taken to simulate the transition of collapsed gully elevation classes, and the latter is taken to spatial interpolate on the simulation results of collapsed gully elevation classes. Based on the two processes, the continuous elevation data of collapsed gully can be simulated and predicted. Then the erosion amount is computed according to the prediction result and base time data. To test this approach, a collapsed gully located in the Longmen Town, Anxi County, Fujian Province is chosen as the case study. And it contains two parts, which are the confirmation of simulating accuracy and the prediction under 3 types of scenes for annual erosion amount. The results show that the modified CA-Markov model performances well on the simulation of collapsed gully elevation classes, which resulted from a lower area error between simulation results and actual value (1.47%), and a higher Kappa coefficient of 0.958. The goodness of fit (2) of the interpolation result using ANUDEM is 0.998 7 and its mean error (ME) is 0.27 m, which means ANUDEM is good at the integral restoration of collapsed gullies terrain but not sufficient at details. As a consequence, two types elevation at the base time are chosen to estimate the erosion amount, one is original DEM, and the other one is the elevation data through classes division and interpolated processing. The erosion amount calculated from the second type of base data performances difference rates of 17.55% and 8.66%, which are much lower than the one derived from the first type of base data. Thus, the second type of base data is accessible for computing collapsed gully erosion amount, and the huger the real erosion amount, the higher the accuracy. After the confirmation of modeling accuracy, 3 types of rainfall scenes were set based on the historical records of rainfall, which are normal years, drought years and rainy years. The annual erosion amounts of case collapsed gully under the 3 types of rainfall scenes are 824.69, 731.03 and 924.57 m3, and the difference of the two extreme rainfall scenes is 193.54 m3, which means that the difference of erosion amount between different rainfall scenes should be taken into account when building a sediment storage dam to control the damage of a collapsed gully. The methods proposed in this research can not only provide a novel way to study the erosion process and erosion amount of collapsed gully, but also supply a reasonable reference for the prevention and control of collapsed gully erosion.

erosions; models; collapsed gully; erosion amount; CA-Markov; ANUDEM

10.11975/j.issn.1002-6819.2018.21.016

S157

A

1002-6819(2018)-21-0128-09

2018-05-18

2018-09-17

国家自然科学基金项目（41601557; 41571272）；国家科技支撑计划项目（2014BAD15B0303）

季 翔，讲师，博士，主要从事景观生态及水土保持研究。 Email：jixiangss@126.com

黄炎和，教授，博士，主要从事水土保持和崩岗侵蚀研究。 Email：yanhehuang@163.com

季 翔，黄炎和，林金石，蒋芳市，余明明，李绍鑫. 基于CA-Markov模型与ANUDEM内插法的崩岗侵蚀量预估[J]. 农业工程学报，2018，34(21)：128－136. doi：10.11975/j.issn.1002-6819.2018.21.016 http://www.tcsae.org

Ji Xiang, Huang Yanhe, Lin Jinshi, Jiang Fangshi, Yu Mingming, Li Shaoxin. Estimation of erosion amount in collapsed gully based on CA-Markov model and ANUDEM interpolation[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2018, 34(21): 128－136. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2018.21.016 http://www.tcsae.org