杨 通，张智杰，王效科※，郭旭东，于 潇

（1.中国科学院生态环境研究中心，城市与区域生态国家重点实验室，北京 100085；2.中国科学院空天信息研究院，北京 100094；3.中国国土勘测规划院，自然资源部土地利用重点实验室，北京 100035）

0 引 言

土壤侵蚀是中国干旱半干旱地区面临的主要生态问题之一，严重的土壤侵蚀一方面迅速恶化生态环境，另一方面加剧当地人民群众的贫困程度[1]。现阶段土壤侵蚀监测方法包括建立观测站、样地调查、立体摄影、高精度GPS、三维激光测量、元素示踪、模型模拟以及遥感监测等[2-4]。遥感方法具有宏观性、综合性、可视性、实用性、经济性等特点，突破了传统的样地调查微观机理研究范畴，将研究层面拓展到面上研究和动态研究[5]，对土壤侵蚀研究发挥重要作用[6]。现行遥感监测土壤侵蚀方法属于光学遥感范畴，监测地表正射水平面的二维特征，以像元流失作为土壤侵蚀的主要依据[7]。然而，土壤侵蚀的二维变化特征并不显著，光学遥感传感器的分辨率也不足以捕捉以毫米计量的土层剥蚀特征，从而低估土壤侵蚀的实际情况，误以为土壤侵蚀很少或侵蚀程度较轻。

土壤侵蚀的主要特征是土层厚度的剥蚀，是高程向的累计形变。水利部发布的《土壤侵蚀分类分级标准SL190—2007》规定：水力侵蚀的土壤平均流失厚度介于5.9～11.1 mm/a，其侵蚀强度为极强烈；风力侵蚀的风蚀厚度介于10～25 mm/a，其侵蚀强度为中度。理论上，若高程向监测精度达到厘米级，即可监测到中度以上风蚀和极强烈水蚀。

雷达遥感对地物的高程变化信息敏感，监测尺度由二维空间转换为三维空间，高程向监测精度可达厘米甚至毫米级[8]，理论上具备识别土壤侵蚀的能力。其中，永久散射体合成孔径雷达干涉测量 PS-InSAR（Persistent Scatterer Synthetic Aperture Radar Interferometry）是雷达微波成像和电磁波干涉的集成技术[9]，是一种新的D-InSAR（Differential SAR Interferometry）技术。该技术通过合成孔径雷达在不同位置对观测区域的相位差等信息分析处理，获取观测区域的三维空间信息。该技术以其覆盖范围大、不受云雾影响、厘米/毫米级监测精度的特点，成为监测地表形变最有前景的空间大地测量技术之一，并广泛应用于火山活动[10]、地震[11]、滑坡[12]、不均匀沉降[13]、湿地变化[14]、地形测绘等研究领域[15-16]。PS-InSAR通过分析时间序列影像，提取时序内相位和幅度保持稳定的离散点，即PS点（Persistent Scatterer），反演地表缓慢形变[17]。与传统的 InSAR技术相比，PS-InSAR克服了空间去相干、时间去相干、大气延迟，显著提高了测量精度和信噪比，其高程向测量精度优于厘米级甚至达到毫米级[8,12]。

目前，国内外少见PS-InSAR监测土壤侵蚀的研究，有学者指出土壤侵蚀造成了上海市的地面沉降[2]，也有学者注意到浙江四明山自然地表在雨季发生沉降，推测为土壤侵蚀所致[18]；国外在伊朗西部[19]、意大利中部[20]有应用InSAR技术监测土壤侵蚀的研究案例。本研究拟利用哨兵一号卫星平台数据（Sentinel-1A），基于改进后的PS-InSAR方法，对内蒙古和林格尔县进行土壤侵蚀的面上观测，结合多源数据对监测结果进行评估和验证，检验PS-InSAR技术监测土壤侵蚀的可行性，为该方法在土壤侵蚀监测领域的推广提供助益。

1 研究方法与数据

1.1 研究区概况

研究区位于内蒙古自治区中南部（39°58′～40°41′N，111°26′～112°18′E），总面积 3 436 km2，海拔 1 400～2 028 m，年平均降雨量为 376 mm，年均蒸发量为2 073 mm。和林格尔县属内蒙古高原和黄土高原的过渡地带，东南部山区、黄土丘陵区的面积占县域面积的77.7%；西北部属土默川平原，面积占县域面积的22.3%。研究区超过一半区域属于黄土高原。黄土高原是世界上土壤侵蚀最严重区域，每年入黄泥沙量达 16亿 t，且每年由于人类活动新增水土流失面积900～1 100 km2。自大规模农耕活动开始，人类活动就已经与水力、风力并列，成为土壤侵蚀的主要营力[1]。根据《土壤侵蚀分类分级标准SL190—2007》分类体系，研究区既属于西北黄土高原水力侵蚀区，又属于内蒙古高原草原中度风蚀水蚀区，水力侵蚀与风力侵蚀并存。

1.2 数据概况

监测数据为2017至2018年的17期Sentinel-1A雷达数据，每期1景，用于监测土壤侵蚀。Sentinel-1A的轨道类型为近圆形的太阳同步轨道，其轨道高度为693 km，倾角为98.18°。卫星重访周期为12 d，内有175个轨道。为定向监测土层厚度变化，选用干涉宽模式（Interferometric Wide，IW）下的S1 TOPS-modes单视复数SLC（Single Look Complex）数据，幅宽250 km，入射角29.1°～46°。该数据为C波段的升轨数据，C波段属于长波，可以穿透稀疏植被直接作用于地表。SLC数据处理后的距离向分辨率为2.3 m，方位向分辨率为13.9 m。

辅助数据为同步观测的哨兵2号（Sentinel-2A）数据，用于主被动遥感数据协同观测；DEM数据用于控制地形影响；土地利用数据用于掩膜不透水面；降雨量和土壤数据用于信息复合。以雷达数据序列为基准数据集，数据源见表1。

表1 数据源一览表Table 1 Data source list

1.3 研究方法

本研究的总体技术方案如图1所示，主要技术流程为：1）预处理流程，包括数据导入、精轨数据更新、主影像选择、DEM主辅影像配准、ESD配准，burst拼接，子条带拼接等；2）差分干涉流程[21]，包括ROI截取、干涉相位计算、去除平地和地形相位、生成差分干涉图等；3）PS点提取流程，包括累计形变点（即PS点）候选点提取、构建Delaunay三角网、选择参考点求取PS点的参数信息、求解二次差分相位模型、残余相位解缠、大气相位估计、累计形变量计算、地理编码等[9]；4）多尺度滤波流程，包括空间域的不透水面掩膜、稳定地表分离，时间域的畸变滤波、线性滤波等；5）侵蚀特征分析流程，包括土壤侵蚀的空间特征、时间特征等。

需要指出的是，PS-InSAR方法监测到的地表累计形变，属于多种人类扰动和自然因素的混合效应，如何从混合效应中分离出土壤侵蚀的贡献，将监测目标锚定于土壤侵蚀，是方法有效性的决定因素。本研究对PS-InSAR方法进行了针对性的改进，基于土壤侵蚀的物候特征，设计了空间域和时间域的多尺度滤波环节，定向提取土壤侵蚀造成的累计形变。

1.3.1 数据预处理

首先，生成 SLC数据，同时利用各子条带及 burst的重叠信息生成mosaic SLC和mosaic MLI（Multi Looks Intensity），多视参数默认为 10和 2。其次，利用 ESA（European Space Agency）提供的精轨 POD数据，（https://qc.sentinel1.eo.esa.int/）修正轨道信息查询列表，减小定时及轨道误差对干涉相位的影响。再次，选择同时位于时间基线中部，及时序影像多普勒质心频率中心的时段的数据作为主影像（表2）。然后，将 DEM与主影像配准，得到每个像元的三维坐标信息。之后，利用ESD方法生成主影像与副影像的几何关系（查询列表），确保配准精度达到1/1 000像素的要求[22]。最后，将分块存储的Sentinel-1A数据根据研究区的位置提取出来，以便缩小处理数据量加快处理时间。此外，对光学遥感数据、DEM、土地利用、土壤、地面观测站数据等，执行校正和配准等处理步骤，生成多源可比数据。

表2 合成孔径雷达数据序列Table 2 Time series synthetic aperture radar data set

1.3.2 差分干涉处理

主影像分别与辅助影像生成17期差分干涉图，并结合卫星轨道参数、成像模型、参考DEM，去地形相位和平地相位。每一个干涉对都需要进行配准和重采样，形成干涉相位。利用轨道和外部的DEM去除平地效应和地形效应，生成时间序列的干涉相位。具体包括感兴趣区选择、干涉处理、地形相位去除、差分模型构建 4个步骤。主辅影像干涉相位可表示为公式（1）

式中φtopo为地形相位；φdef为形变相位；φatm为大气延迟相位；φflat为平地相位；φnoise为噪声相位。

1.3.3 提取永久散射体

PS候选点选取：选择后向散射较强，并且在时序上较稳定的地物目标，如建筑物与岩石等，将相干性差或者一些噪声点去除，留下高相干且稳定的点进行分析。研究区有150个行政村，并有蒙牛、蒙草等大型企业，建成区的人工建筑可以部分替代强散射点的作用。为了保证PS候选点密度，本文采用平均相干性阈值法选点，包含一些中等相干性点目标，以增加PS候选点的选点数量。

构建三角剖分网：利用Delauany三角化方法构建三角剖分网络，以解决PS离散点的图像重构问题。同时，利用三角网络的冗余测量边，以区域网络平差方法控制形变测量误差，削弱大气相位影响。Sentientl-1A时空基线较短，轨道误差较小，通过构建三角网及相邻PS点的相位双差模型等方法，可以削弱大气相位的干扰。

计算二次差分相位：基于PS点的三角剖分网，计算网络上每条边的相邻PS 点的二次差分相位。通过求解二次差分相位模型，获得PS点对的形变速率差和DEM修正差。

相位解缠：相位解缠后可得每个PS点的绝对形变速率和DEM修正量，然后可得差分相位的形变相位和残余地形相位。

大气滤波：大气延迟相位在时间域上表现为随机高频信号，在空间域上则表现为平滑低频信号；非线性形变在时间域上表现为低频信号，在频率域表现为低频信号；噪声相位则为时域高频，频域高频。根据各相位成分的特点，通过时空域高通、低通滤波分离各相位成分，从而得到形变相位，再由地理编码得出累计形变量。

1.3.4 多尺度滤波

地表累计形变是多种因素共同驱动的叠加效应，包括人类扰动和自然变化。自然变化又包含植被生长、土壤侵蚀、滑坡、泥石流、冻土、岩石风化搬运、枯落物堆积、地质运动等多种类型。本研究设计了空间域和时间域的多尺度滤波，剥离人类扰动和自然变化影响，将监测对象锚定为土壤侵蚀引发的地表形变。

空间域滤波又分为植被、稳定地表滤波、人类活动滤波。首先，通过计算主末影像的相位标准差，加入植被指数的限制因素，振幅小的信号通过，振幅大的剔除，滤除茂密植被，将感兴趣区限制于植被稀疏区域[23]。雷达C波可以穿透稀疏植被，直接作用于地表。其次，设置稳定地表值域，对永久散射体的累计形变信号进行高通滤波，滤除无变化或变化极小的区域。然后，根据土地利用数据掩膜城市村镇的不透水面及活跃耕地，并设置阈值分割上限，对永久散射体的形变信号进行低通滤波，滤除人类活动干扰。

时间域滤波分为线性滤波和畸变滤波。不同自然因素驱动的地表形变在时间序列上有其规律性，如土壤侵蚀表现为雨季的非线性形变；滑坡、泥石流表现为突发性形变；地质隆起或沉降表现为平滑线性形变；冻土表现为双向形变等。根据各种自然因素在时间序列上的形变方向、形变区间、形变幅度差异，可以初步区分主要形变类型；再结合高分光学遥感数据的光谱信息、几何信息及动态变化进行检验。选取侵蚀特征显著的样本点，统计其时序变化参数，迭代更新初始参数，修正分离区间和分离阈值，进而实现土壤侵蚀的定向识别。

2 结果与分析

2.1 形变速度场

基于时序差分干涉结果，提取后向散射系数较高、并且在时序上较稳定的永久散射体。对PS候选点进行迭代更新、解缠等处理后，在研究区共计提取到 1 007 958个累计形变点，样本量较大，信息丰富。为抵消卫星在轨运动时产生的多普勒效应，PS-InSAR运行过程均在斜距-零多普勒坐标系中完成，最后再通过地理编码进行可视化处理，将坐标系统转换为地理坐标系，生成可比数据。对1 007 958个累计形变点的平均形变速度进行IDW（Inverse Distance Weighting）插值，计算本区综合形变速度场，如图2所示。由于累计形变点在空间上属于非均衡分布，导致相干条件差的区域提取到的PS点数量较少，

插值效果不佳，主要为和林县北西侧的集约化农垦区域。总体上看，县域大部分区域PS点密度较高，反演的年均形变速度场能如实反映区域地表形变的一般性规律。研究区综合形变速度场分布特征为东南侧隆升、西北侧沉降，构造隆升与侵蚀夷平并存。构造上主要受NE方向的块断构造控制（《1∶25万呼和浩特市幅（K49C004003）区域地质调查报告》），属于阴山-燕山板内造山带的西段，鄂尔多斯盆地东缘。本区自燕山运动以来处于阶段性隆升态势[24]，全新世至今山体不断隆升，形成了山前广为发育的河流阶地和冲洪积平原。其中，新近纪古夷平面的发现，证实了本区地表形变的主因是构造隆升与侵蚀夷平的交互驱动。由图2可知，沉降区主要分布在黄土丘陵及陆前平原，位于NE向高角度正断层的北西侧（该断层分割了山区与平原），也是该县的人口聚集区，沉降成因可能包含多种因素，诸如地质凹陷、人类活动、滑坡、泥石流、土壤侵蚀等，需进一步细分形变机制，定向提取土壤侵蚀形变。

2.2 累计形变点筛选

雷达时序差分干涉在县域内共计提取累计形变点1 007 958个，包含不同区位条件、不同驱动因素的地表形变。从1 007 958个累计形变点中筛选出土壤侵蚀形变点，需要对全部累计形变点进行再次筛选。本文设计了多尺度滤波流程对累计形变点进行筛选，包括空间域滤波和时间域滤波。空间域滤波用于消除人类扰动影响。结合2017土地利用数据掩膜人类活动形成的不透水面及工作面；再使用振幅标准差滤波排除茂密森林植被及活跃耕地，使监测范围缩小至天然裸露（半裸露）地表；然后设置分割阈值，排除无形变或形变极小的稳定地表，依据前人研究设置稳定地表阈值为±1 cm/a[25]。空间域滤波结果如图3所示。宏观上，累计形变点集中分布于县域东南部的山区丘陵地区，山区丘陵是地表形变的多发易发区域，同时相干条件较好，永久散射体密度较高。滤波区主要为县镇村的建成区、集约化农业区、林区、公路铁路及设施用地等。

时间域滤波用于分离地表自然形变。不同类型自然形变的物候特征不同，在时间序列上的差异性显著，具体表现在形变方向、形变区间、形变幅度、形变趋势差异等。形变方向分为抬升正向形变和沉降负向形变；形变区间分为全年形变、季节形变和突发形变；形变幅度分为平缓形变和剧烈形变；形变趋势分为单向渐进形变和双向波动形变。根据土壤侵蚀的物候特征，其形变特征表现为全年均有发生、在雨季骤增的渐近式负向形变。根据时间序列的形变差异，进行有监督的阈值分割，可以初步圈定侵蚀发生地。再与同步观测的Sentinel-2A高分光学遥感影像叠加，根据贝叶斯理论增加已知信息量[26]。检验形变点所处区域的光谱特征、区位特征、形态学特征和动态变化轨迹是否与典型侵蚀景观一致。若否，则返回上一步，迭代更新初始阈值，直至确定最优分离阈值。

经由时间域滤波及多源数据信息符合，确定土壤侵蚀形变的最优分离区间及约束条件，见公式（2）～（6）

式中Sc为年累计形变量（mm），Si为单月形变量（mm），i为月份，n为年份，R为区域降雨系数。时间域滤波的基本逻辑是，基于土壤侵蚀的物候特征及先验知识，对累计形变点的形变方向、形变区间、形变幅度及形变趋势加以限定，从而分离出符合侵蚀特征的形变点。经由空间域和时间域滤波，从1 007 958个累计形变点中，提取出土壤侵蚀形变点 10 596个，占累计形变点总数的1.05%。

2.3 区域土壤侵蚀特征

2.3.1 空间特征

前人研究表明土壤侵蚀与植被[27]、地形[28]、降雨[29]、土壤[30]等因素具有显著的相关性。将10 596个土壤侵蚀点与同步观测计算的动态数据（NDVI分布图、降雨量插值图）及静态数据（坡度图、土壤可蚀性）进行信息复合，评估本区土壤侵蚀的宏观分布特征（图4）。其中NDVI由2017年7月的Sentinel-2A数据计算归一化植被指数得出；降雨量由和林格尔县气象局提供的19个气象观测站的 2017—2018降雨数据，进行克里金插值得到；坡度分布由地理空间数据云下载的DEM数据计算得出；土壤可蚀性由中科院南京土壤所的土壤数据计算K因子。

图4a可见侵蚀点云多分布于低植被或无植被覆盖区域，植被群落以草本植物为主，且较为稀疏，光学遥感特征近似裸地。在NDVI介于0至0.5区间内（NDVI大于0.5区段被滤除），土壤侵蚀点主要分布于NDVI零值附近，土壤侵蚀多发于无植被保护区域。图4b可见县域降雨量分布极不平衡，整体呈东高西低分布，东部山区的年降雨量是西北部平原的 3倍。本区土壤侵蚀主要分布于300 mm等降雨量线以下的干旱半干旱区域。图4c可知平原与山区丘陵的地形分界线约为北东向 40°，线性构造带发育，从地质地貌角度属于土壤侵蚀多发、易发地带。土壤侵蚀点主要分布于坡度 0～35°范围内，平原区的侵蚀类型主要为风蚀，山区的侵蚀类型主要为水蚀。图4d可见土壤可蚀性与侵蚀点云一致性较高。侵蚀点云聚集区多位于土壤可蚀性较高区域，在土壤可蚀性较低区域则分布较少且稀疏。总体而言，本区土壤侵蚀主要分布于山区平原过渡带，侵蚀点云约呈北东向40°分布，贯穿县域南北，在县域西北、县域东南也呈局部聚集态势。

2.3.2 时间特征

对PS-InSAR观测结果进行月度回溯，对全部土壤侵蚀点的累计侵蚀量按月取平均值，同时计算月度累计降雨量平均值。如图5所示。

研究区土壤侵蚀总量逐月连续增加，即便在没有降雨的12月、1月、2月，土壤侵蚀仍在发生，说明风力是本区土壤侵蚀的重要营力。此外，在降雨集中的 6—9月，土壤侵蚀量显著增加，变化梯度增高，表明降雨对研究区的土壤侵蚀具有重要驱动作用。研究区土壤侵蚀的主要营力为风力水力混合，侵蚀类型为风蚀水蚀混合。平均侵蚀速率为15 mm/a，风蚀强度为中度。《土壤侵蚀分类分级标准SL190—2007》将本区划分为内蒙古高原草原中度风蚀水蚀区，本文监测结果与之一致。

2.4 监测精度评估

PS-InSAR方法是一种直接观测方法，比之光学遥感的米级分辨率，多时相SAR监测地表形变的精度理论上可达厘米级甚至毫米级[8,12]，监测量纲提升2个数量级。然而，受限于初始数据精度和研究对象的复杂程度，实践中可能难以达到毫米级的测量精度。

InSAR方法的精度验证多为与其他高精度GPS数据进行比对[31]、及实地调查形变径迹[32]等方式，本文同时采用上述两种方式验证精度。依据《水土保持工程调查与勘测标准》（GB/T 51297—2018）规定，课题组赴研究区开展了野外调查验证，共调查 103处土壤侵蚀点，采用侵蚀径迹调查、高精度GPS测量、实地拍照的方式，验证所提土壤侵蚀点与侵蚀径迹的对应准确度（如图6）。

调查中发现，绝大多数调查点位于侵蚀沟（图7a）、风蚀柱（图7b）、新鲜红黏土（图7c）等典型侵蚀径迹区域。经统计，土壤侵蚀点的提取准确度为86.4%。此外，收集了和林格尔县境内 142个四等水准点数据，依据最临近原则对四等水准点及其最近的PS点的绝对高程进行比对，同时与SRTM DEM（30 m）在水准点位置的绝对高程进行比对（图7d）。由于PS点与四等水准点在空间位置上不重复，故设置最大搜索距离为 10 m，相距超过10 m的点对不纳入统计。统计后发现，四等水准点与其最近的PS点的绝对高程值非常接近，且接近程度与距离反相关。PS点和水准点相距1 m以内的绝对高程相对偏差为7～8 cm；相距10 m以内的相对偏差介于几米至十几米之间。尤其在水准点1 m以内的监测精度优于SRTM DEM两个数量级（后者为米级监测精度）。总体来看，本文方法监测土壤侵蚀形变具有可靠性。

3 讨 论

国内外使用InSAR方法监测土壤侵蚀的研究较少，该领域整体研究程度不高，缺乏成熟可供借鉴的研究经验。本文旨在探索并发展PS-InSAR监测土壤侵蚀的理论基础和技术框架，并论证技术框架的可行性和有效性，为后续更深入的研究提供参考。方法运行结果表明，监测到的土壤侵蚀形变在空间上与植被、土壤、降雨、地形因子相关，在时间上与降雨季节相关，在实地验证中也发现了侵蚀径迹，且与四等水准点的一致性较高，说明 PS-InSAR方法监测土壤侵蚀具有可行性和可靠性。从精度验证结果来看，所述方法具备厘米级监测精度，可以识别中度以上风蚀和极强烈水蚀，适用于中国北方干旱半干旱地区的土壤侵蚀动态监测及分析评估。后续研究包括：1）采用分辨率更高的InSAR和DEM数据源，改进优化技术框架，提高监测精度至亚厘米级或毫米级，并开展“三同（同时、同地、同精度）”地面精密水准测量以验证监测精度；2）复杂地形地貌区域的土壤侵蚀监测及侵蚀模数反演；3）基于InSAR监测的土壤侵蚀模型构建或现行水土流失模型修正等。

地表形变场是多重驱动力下的混合形变效应，提取侵蚀形变的基本逻辑是在空间上圈定侵蚀易发地，在时间上符合侵蚀发生发展的物候规律。本文设计了多尺度滤波流程，对 1 007 958个累计形变点进行筛选，提取出侵蚀形变点10 596个，占累计形变点总数的1.05%。需要说明的是，本文的研究假设是各类形变在水准面上具有二维可分性（比照光学遥感分类），然而部分区域可能构造运动与土壤侵蚀并存，土壤侵蚀量的一部分被构造隆升量所抵消，从而低估实际土壤侵蚀量。本文的处理办法是，在发现构造隆升趋势后，返回到提取永久散射体步骤，重新选择了参考点，对每个参考点进行多次回归迭代，更新高程误差改正值，修正了最终形变结果。因此，在现今构造运动活跃区域开展土壤侵蚀监测，还需深入解析形变机制机理，叠加或分离现今构造运动修正。

研究区位于黄土高原北坡，马兰期风成黄土属于自重湿陷性黄土，易发生侵蚀形变[33]。本区土壤侵蚀的一大特征是与人类活动相关，侵蚀严重区域的人类活动较为密集，约10%的侵蚀点位于耕地范围。耕地范围的侵蚀点主要位于撂荒地、坡耕地及大量井灌引发的地下水超采区域。土壤侵蚀对耕地耕层质量退化作用表现为土壤性质恶化、土壤质量劣化、土地生产力衰退3个方面。应注意避免陡坡开垦，加大退耕还林还草力度。在县域西北极干旱地区，年均降雨量只有150 mm左右，却是该县的粮食主产区，且未见显著的河流、湖泊等地表水补给源，集约化农业大量依赖地下水，势必引发地下水超采问题。地下水超采会降低土壤含水率，导致对土壤墒情敏感的植被群落衰亡，从而使局地生态环境逆向演替。建议加强对该区域地下水位的动态监测，疏导产业结构向低耗水量行业转型。

4 结 论

本文利用土壤侵蚀主要特征是土层厚度剥蚀、且合成孔径雷达对高程信息敏感的特点，研究PS-InSAR技术监测土壤侵蚀的可行性。基于土壤侵蚀的物候规律对技术框架进行了针对性的改进。通过多源数据融合增加已知信息量，构建时空域多尺度滤波方法分离混合形变效应，并给出了时间域滤波模型及参考值，使监测对象锚定于土壤侵蚀形变。经验证，本文方法获得了与实际侵蚀发生地一致的监测结果，土壤侵蚀点提取准确度为86.4%，高程向总体监测精度为厘米级。侵蚀点在空间上表现出与植被、降雨、坡度、土壤的强相关性，在时间上表现为在雨季骤增的非线性负向形变。侵蚀外营力为风力水力混合，平均侵蚀速率为15 mm/a。本文方法或可为土壤侵蚀的宏观面上研究提供参考。