0.   引　言

岩溶地貌一般指可溶岩在水流的溶蚀作用下形成的地貌^[1]。岩溶地貌在全球广泛分布，约占据陆地总面积的12%^[2]。中国是世界上岩溶地貌分布最广泛、发育最充分、类型最齐全的国家，其独特的地形条件不仅塑造出奇异多样的岩溶地貌景观，同时也深刻影响岩溶地区资源分布、生境系统和人文格局。岩溶地貌是岩溶形态的综合反映，它具有数字表达其形成过程和空间特征的特性，并能用计算机来建立岩溶地貌的空间分布模式。

近年来，随着航空航天测高技术的发展，数字高程模型（Digital Elevation Model，DEM）数据的精度、分辨率、覆盖性得到不断提高^[3]，为地貌研究融入了一系列的新思路。由于DEM中高程这一直接描述地表起伏形状的几何量，在地貌研究中非常广泛，特别是地貌要素的提取、地貌类型的分类、区域地貌的分析、水力作用造成的地貌演化模拟以及三维可视化表达等方面已有诸多应用^[4-6]。而基于各种空间分辨率的DEM数据进行相应研究尺度的数字地形建模，提取一系列地形参数，即数字地形分析（Digital Terrain Analysis，DTA）的过程，进而在岩溶地貌的地形分析、形态建模、分布格局及其演变规律等方面都取得了一系列进展。

本文针对前人基于DEM在岩溶地貌应用研究中所涉及的关键技术方法，系统总结出近些年来DEM空间分析技术在中国岩溶地貌研究中的新进展，从岩溶地貌识别的尺度效应、岩溶地貌的识别与分类、岩溶地貌的形态及格局分析、岩溶区生态环境变化等方面进行了总结和评述，并对未来DEM在岩溶地貌研究中的方向进行了展望，以期能够扩展DEM在岩溶地貌研究中更多的应用，为岩溶地貌从地形研究拓展到地形演变的过程与机理研究提供技术思路。

1.   基于DEM的岩溶地貌识别尺度效应

DEM数据的获取方式众多，针对不同的地貌研究尺度及研究细节常使用三类数据源：即高精度的无人机（UAV）摄影测量和激光雷达（LiDAR）、国家数字高程模型数据库、全球数字高程产品。无人机摄影测量和激光雷达数据主要用于精细化三维地形建模，获取成本较高；中国数字高程模型数据库有1∶100万、1∶25万、1∶5万及1∶1万4种不同比例尺DEM数据，其栅格格网分辨率分别为1 km、100 m、25 m和5 m。而全球数字高程产品由于其覆盖范围广、零成本、易获取，得到学界更广泛应用。目前，已公开发布了多种全球数字高程产品，如SRTM、ASTER-GDEM、AW3D30等^[7-8]，这为岩溶地貌的数字地形分析提供了数据支持。

应用DEM进行岩溶地貌识别前不得不面对数据模型的适用性问题，不同的精度和分辨率都会导致岩溶地貌形态度量与可视化的不确定性。精度代表数据记录地形信息的准确能力，其误差会引入到模型和地形分析结果中，直接影响实验结果的可靠性。由于DEM数据具有空间性，其误差分布亦与地貌地形条件有较强相关性，特别是在峰丛洼地和岩溶峡谷地区，DEM的精度相对在平原丘陵区有所下降。对西南岩溶区进行对比试验发现，SRTM、ASTER-GDEM等开放DEM数据在大、中起伏的岩溶山区的误差绝对值尽管变大，但变化幅度较小^[9]，进行地形因子提取的结果质量均可满足研究需要^[10]。岩溶地形对DEM精度的响应并不如其他地区（如黄土丘陵沟壑区）明显，岩溶地区采用1∶5万DEM信息源进行量化因子的提取即有较高的精度^[11-12]。DEM分辨率表征数据的识别能力，影响DEM尺度效应和分析模型的匹配问题。岩溶地貌的坡面尺度、流域尺度、区域尺度特征对DEM的分辨率要求不同是毋庸置疑的，精细尺度需要更高的分辨率来展现更多的细节特征。但对岩溶地貌识别而言，必须考虑区域地形的主体特征，尤其在特定区域（例如洼地或漏斗底部）内坡度变化不大的情况下，空间分辨率过高，会导致像素的过度异质性，混淆地形的主体特征，产生较多的识别噪声。研究者通过使用目前流行的全球数字高程数据集进行岩溶地貌目标的识别精度验证，发现90 m分辨率的SRTM相较30 m分辨率的ASTER-GDEM和AW3D30等数据能够更准确识别出目标^[13]，这也说明了更高的空间分辨率并不一定能保证更好的岩溶自动探测效果。

由于局部地形之间的相关性会受到边缘效应的影响，DEM及地形分析还具有明显的尺度依赖性。在计算或提取局域地形参数时，为能恰好反映地貌的完整性，使计算结果在一定区域内具有较强的代表性，需确定适宜分析窗口（图1）。地貌类型、研究范围、数据类型和数据分辨率都会直接影响适宜分析窗口的大小。目前定量选择适宜分析窗口的基本方式是通过多窗口局部地表形态刻画法获得尺度效应曲线^[14]，在均值变点处确定适宜分析窗口大小。研究者在应用DEM进行岩溶地貌分析时，尝试不同条件下确定适宜分析窗口大小，结果表明西南岩溶地区使用ASTER GDEM提取地表起伏度的最佳统计单元为17$\times$17像元，即面积为2.34 km²的格网窗口^[16]。赵宇鸾等^[17]考虑到西南岩溶山区涉及到乌蒙山、武陵山脉等大的山地，较小的分析窗口只能反映局部地形的起伏，不能涵盖大型山脉的山顶、山脚，其提取二次均值变点，选取6.50 km²为黔桂岩溶区地形起伏度最佳统计单元。孟欣^[18]统计出不同岩溶地貌样区的形态特征和地形因子，结果表明最佳分析半径（最佳分析窗口）与山顶点密度之间存在明显的负相关关系，对样点做线性拟合，运用最小二乘法求解，得到两者间的相关系数。

图  1  尺度效应的定量刻画曲线示例 ^[15]

Figure  1.  Examples of curves for quantifying the scale effects in DTA

(a. scale effect with resolution change; b. scale effect with neighborhood size change)

下载: 全尺寸图片 幻灯片

2.   基于DEM的岩溶地貌识别与分类

地貌分类体系的研究目前已取得较大进展，已完成全国数字地貌五级的划分，建立全国地貌区划数据库^[19]。由于统一在全国基本地貌类型的框架下，岩溶地貌仅作为一种成因类型进行约束，虽从区域地形起伏的角度，进行了岩溶分布区基本地貌的划分^[15-20]，但难以体现岩溶发育对区域地貌的控制程度。因此，还需从岩溶地貌单元尺度出发，辨别岩溶个体和组合形态，目前利用DEM数据自动提取岩溶地貌单元的方法主要有四种：基于多因子复合的方法、基于地形特征点提取的方法、基于等高线树的提取方法和基于正负地形划分的方法。

（1）基于多因子复合的方法是将研究区域进行一系列的网格划分，计算每个单元网格的地形属性（表1），如高程、坡度、平均曲率、高程半变异系数、粗糙度等，构建包括峰林、峰丛和其他岩溶地貌单元的识别模型。实验证明，该方法能够对峰林、峰丛和常态山进行有效分类，总体精度在80%以上^[21]。然而，网格式的分类忽略了地貌单元边界，更重要的是该方法在两种不同地貌类型之间的边界上分类结果非常不稳定。相较于以网格为划分单元，流域单元在地表形态和地貌演变方面具有明确的地理意义，且能契合地貌和水文之间紧密的关系，因而通过构建面向小流域单元的地貌识别模型也有很好的效果^[22-23]，该方法的分类结果能进行岩溶地区总体地貌格局的统计与分析，但由于边界的模糊性，难以进一步进行形态计量学分析等精细化研究。

表  1  峰林、峰丛岩溶地貌的形态指数特征

Table  1.  Morphological characteristics of peak forest and peak cluster

测度	形态指数	峰林、峰丛地貌与一般非岩溶地貌对比
垂向测度	高程	峰林高度较低，一般几十到两百米；峰丛高度通常在两三百米以上，最大可达六百米以上
	起伏度	峰丛地表粗糙，高程起伏相对变化大
水平测度	邻近性	峰丛具均匀的峰洼分布，峰林零星分布
	隔离性	峰丛具孤立的洼地和互相连通的峰体
形状测度	地形表面积	峰丛具较大的地表面积，峰林具较小的地表面积
	形状指数	峰丛岩溶地区由高频率的简单形态山体构成，形状指数小；非岩溶地貌一般由低频率的复杂形态山体构成，形状指数大
垂向和水平测度	坡度	峰林、峰丛地貌较一般丘陵坡度更大，峰丛坡度一般大于30°，峰林坡度一般大于45°
特征要素	山顶点/洼地点	峰丛山顶点和洼地点密集且相间均匀分布
	鞍部点	峰林无明显鞍部点，峰丛有鞍部点且在一定范围内鞍部点围绕洼地点
	山脊线和山谷线	非岩溶地貌通常由山脊和山谷组成，山脊的海拔会逐渐下降到山谷，如果断面恰好沿着山脊线或谷底，则减少得更慢

下载: 导出CSV 

| 显示表格

（2）基于地形特征点提取的方法是利用DEM数据提取地形特征点，根据典型地貌单元的形态特征，构建不同地貌单元地形特征点的空间关系模型，进而自动提取不同的地貌单元。Yang 等^[24]认为洼地是峰丛间鞍部点围成曲面在水平面上的投影，可根据山脊线和山谷线的交叉点提取鞍部点，对鞍部点高程插值生成空间趋势面，原始地形与趋势面叠加计算即得到洼地的分布。基于这种思想，也可根据山顶点、鞍部和凹陷点之间的空间关系区分峰林和峰丛单元^[25]。首先通过山顶点和凹陷点构建泰森多边形来控制正负地形区域，如若控制区域内有且仅有一个山顶点没有鞍部点，该区域是一个峰林单元；如若峰顶之间存在明显的鞍部，且鞍部点围绕一个凹陷点形成封闭区域，那么这些峰顶控制的区域组成一个峰丛单元；如若鞍部存在于峰顶之间，无法形成封闭区域且不存在凹陷点，那么这些峰顶控制的区域构成链状峰林。该方法从岩溶地貌发育的规律出发，峰林地貌是通过外源水溶蚀形成地表小型冲积平原，峰丛地貌是通过内源水逐步向下侵蚀形成封闭洼地，因而，两者的主要区别是有无圈闭洼地。该方法的准确性依赖于区域特征点的提取，山顶点、凹陷点在复杂地形条件下自动提取还不够精确，然而特征点的提取精度会直接决定后续的识别效果。

（3）基于等高线树的提取方法通过建立等高线识别模型，提取地貌轮廓线，从而达到分类目的^[26-27]。其思路是首先提取最外层的闭合等高线作为地形实体的边界，通过实体轮廓内闭合等高线的空间关系及其高程的变化建立等高线树。峰林、峰丛、洼地地貌单元及其组合的等高线树结构特征模式不同，进而确定每种岩溶地貌实体的类型。如图2所示，简单地形（如洼地、孤峰）等高线树具有单支结构，从根节点到叶节点，高程逐渐升高，则判断为孤峰，若高程逐渐降低，则判断为岩溶洼地。若等高线树有分支结构，将其分解为多个单支树进行求算，分支树均为孤峰类型，则该多支树为链状峰林；分支树均为洼地类型，则该多支树为复合洼地；分支树由洼地和峰体组成，则该多支树为峰丛洼地。以往测绘特定地形实体主要通过人工检查地形图上的等高线来实现，而等高线树的方法符合人的认知，实际上是将传统人工目视的方法总结为逻辑模型，实现了计算机的自动化提取。更重要的是，该方法能够提取出岩溶地貌单元清晰的边界，将连续的岩溶地貌划分为离散的地貌实体，这对于后期定量描述岩溶地貌形态特征，进行岩溶地貌形态结构的空间分异性评价非常关键。

图  2  典型岩溶地貌的等高线树模型示意图

Figure  2.  Contour tree model of typical karst landscape

下载: 全尺寸图片 幻灯片

（4）基于正负地形划分的方法，将复杂多样的地形降维简化成正负地形，有效地突出形态差异，再根据峰林、峰丛等地貌单元的特征逐步划分^[28]。传统填洼处理的水文分析法在选择阈值的过程中，采用的是个体经验和试错法，一定程度上不是定量的。Yokoyama等^[29]提出了地形开度指数用以表达一个位置在不规则表面上的开放或封闭程度，它的特点是摒弃了传统无差异性的规则矩形的计算单元，而是计算八个方向的平均最小天顶角和平均最小天底角，这能在一定程度上自动调整和适应地形起伏变化情况，更好地顾及地貌空间分布的多尺度特性^[30-31]。因而，地形开度是划分岩溶地区正负地形的理想指标。

3.   基于DEM的岩溶地貌形态及格局分析

地形因子的提取是数字地形分析的核心，DEM可派生出多个地形因子，目前已有的地形因子高达上百种，基本形成系统性和层次结构的地形因子分类框架，从不同维度表征地貌单元形态及空间格局。根据描述尺度的范围，可将地形指标体系分为微观地形因子和宏观地形因子，宏观地形因子又可进一步分为坡面尺度、流域尺度、区域尺度；根据各因子的功能又可将地形因子划分为形态因子和过程与机理因子；此外，还可依据所描述的不同的点、线、面对象将地形因子进行分类。

在岩溶流域中，影响地貌形态发育的因素错综复杂，岩溶地貌形态特征并非单一指标所能描述，需构建有效反映岩溶地貌形态特征的组合指标^[32]。目前使用较多的方法流程是将统一量纲后的地形因子进行多元回归分析，采用雪式熵值法剔除相关系数过大的地形因子，进而选取最佳地形因子组合^[33]。国内外学者通过计算洼地的多边形边数、半径、投影面积、周长、峰洼高差、洼地标高，峰林、峰丛的坡度、体积、峰顶高程等指标，挖掘指标之间的相关关系，进行岩溶地貌个体形态特征统计与对比分析^[34]。通过聚类分析、邻近分析法统计岩溶地貌的空间随机分布特性，如空间密度、峰洼个数等，研究不同空间尺度下峰林、峰丛地貌形态的群体组合特征，进一步探寻岩溶地貌的空间格局。杨先武^[35]提供了一套具有实际应用意义的岩溶地貌空间分异特征指标，分别从坡度、剖面凹凸度、圆度、紧度、偏心率、延长率等描述个体形态，从地形粗糙比、点蚀指数、莫兰指数等描述群体形态，全面构建峰林、峰丛形态定量评价指标体系（表2）。

表  2  可用于表征岩溶地貌空间格局的典型指标参数^{[35, 39-40]}

Table  2.  Typical index parameters used to characterize the spatial pattern of karst landscape

类别	量化因子	计算公式	备注	地学意义
形态统计特征	坡度	$\beta = \arctan \sqrt {f_x^2 + f_y^2}$	式中：fx是X方向高程变化率；fy是Y方向高程变化率	反映地表面在该点的倾斜程度
	粗糙度	${\text{R}} = {S_s}/{S_p} = 1/\cos (\tan \beta )$	计算地表的曲面面积Ss与其在水平面上的投影面积Sp之比；tanβ为DEM栅格单元坡度	反映地表的起伏变化和侵蚀程度的指标
	复合地形指数	${{CTI} } = In(\alpha /\tan \beta )$	式中：α 表示单位等高线长度的汇水面积；tanβ为该处的坡度；CTI又称地形湿度指数TWI	对径流路径长度、产流面积等的定量描述，也可反映地形的复杂性
	分形维数	${{\text{F}}_{\text{d}}} = - \log N(\varepsilon )/\log \varepsilon$	式中：ε为栅格格网边长大小；N（ε）为栅格总数	表征不同地貌类型下峰体形态的自组织程度
空间展布特征	形状指数	$S = \displaystyle\sum\limits_{i = 1}^N {{W_i}} \dfrac{{{P_i}}}{{2\sqrt {\pi {A_i}} }}$	式中：N为该地区总斑块个数；Wi为第i个斑块的面积权重；Pi为第i个斑块周长；Ai为第i个斑块面积	反映地貌单元景观斑块在空间结构上的不规则程度
	邻近指数	$PI = \displaystyle\sum {({{{a_j}} / {h_{ij}^2}})}$	式中：aj表示斑块面积；hij表示斑块ij到同类型斑块的最近距离	可表征峰洼之间的隔离趋势
	莫兰指数	$I = \dfrac{{\displaystyle\sum\limits_{i = 1}^n {\displaystyle\sum\limits_{j = 1}^m {[ {( {{x_i} - {x_m}} )( {{x_j} - {x_m}} )} ]} } }}{{\displaystyle\sum\limits_{i = 1}^n {{{( {{x_i} - {x_m}} )}^2}} }}$	式中：xi和xj分别为在位置i、j的测量值；xm是所在所有i、j位置点测量值的均值；n为所有测量点的数目	反映地貌景观斑块在空间上的集聚程度
地貌发育 演化特征	峰洼密度	$D = \dfrac{{{N_p} + {N_s}}}{A}$	式中：Np为样区内峰顶个数；Ns为样区内洼地个数；A为样区的面积	表示峰体洼地的聚集程度，侧面反映岩溶发育程度
	面积－积分值	${E_i} = \dfrac{{\displaystyle\int_0^H {adh} }}{{HA}} = \displaystyle\int_0^1 {xdy}$	式中：a表示水平断面面积；h是等高线的相对高程值；H是样区的高差；A是样区的面积	通过构建不同等高线上的面积和相对高差之间的函数关系来评价地貌演化阶段和侵蚀动力差异

下载: 导出CSV 

| 显示表格

此外，研究者发展了地面纹理、坡谱信息熵、地形信息图谱、分形指数等用于特定领域的研究，进一步认识地貌的空间差异性规律。例如，利用岩溶区地形具有自相似分形的特征^[36]，基于高分辨率数字高程模型对岩溶地貌系统进行分形分析，刻画地表形态的复杂纹理特征^[37]。在地貌学中，空间换时间来获得地形演变的信息是常用的策略。通过构建Strahler高程—面积积分曲线可量化岩溶侵蚀地貌演化阶段，结合区域地质地貌发育证据，能大致划分出戴维斯地理循环论的青年期、壮年期和老年期地貌发育阶段^[38]。

4.   基于DEM的岩溶区生态环境变化分析

岩溶地貌具有脆弱的生态环境，不同地形条件使地表的光、热、水、气、土壤、植被等发生再分配作用，导致自然生态环境的差异^[41]，因而基于DEM可用于岩溶区气象、水文、土壤、生境等对象的地形统计、空间格局、空间估算和三维可视化，进而分析岩溶区生态环境变化（图3）。近年来专家学者也逐步关注岩溶地形要素与土壤侵蚀、土壤养分、植物物种多样性、甚至区域贫困格局之间的关联关系，利用地形因子对已建立的各种评价预测模型进行修正，以适用于岩溶特殊环境，这对于更加精准评价、规划、治理、修复岩溶地区生态环境具有重要意义。

图  3  DEM在岩溶区生态环境中的主要应用模式

Figure  3.  Main application modes of DEM in the ecological environment of karst area

下载: 全尺寸图片 幻灯片

土壤侵蚀是岩溶地区石漠化的重要原因。研究表明，土壤侵蚀随坡度和坡长的增大而增加^[42]。通用土壤流失方程（USLE）或其修正形式（RUSLE）是一种广泛使用的经验统计模型，该模型常用于流域或更大尺度的土壤侵蚀预测和评价，在国内外的岩溶地貌中得到应用^[43-45]。USLE和RUSLE的应用表明, 计算流域内各点的LS因子是区域土壤侵蚀评价的重要基础，为此学者们针对区域尺度和流域尺度提出了LS的计算方法。为弥补LS因子在区域土壤侵蚀状况判定方面的不足，也有学者从坡谱信息熵等其他角度来揭示土壤侵蚀程度^[46]。然而岩溶地区严重的石漠化和浅层土壤，RUSLE模型往往高估了土壤侵蚀速率，需要改进以提高模拟精度。确定岩溶土壤侵蚀的空间梯度分布及其主要影响因素，修正RUSLE模型，以适应岩溶地区的评价^[47]。土壤深度也是评价岩溶石漠化的重要指标之一。DEM的空间分辨率不同，土壤性质与地形属性之间的相关性也不同，特别是在复杂的岩溶山区如峰丛洼地更为明显^[45]。在确定实现精确土壤—景观模块所需的DEM最佳空间分辨率下，使用多元线性逐步回归分析结果表明，坡度、湿度指数和高程是预测土壤深度的最重要变量^[48]。

地形也是影响区域植物群落组成、分布及物种多样性的主要环境因子之一，利用相关分析法研究岩溶山地不同地形下乔木、灌木和草本植物物种多样性差异及地形因子对其分布格局的影响，结果表明岩溶山区植物物种多样性及其分布格局受以海拔、坡度为主，坡位和坡向为辅的多个因子共同影响^[49]。区域贫困格局与地形要素也有一定空间吻合性，主要表现为地形条件较差的地区农民收入水平较低^[50]，利用数字地形分析辅助扶贫开发模式策略，提高精准扶贫、乡村振兴服务工作的成效。

5.   总结与展望

随着数字高程模型（DEM）数据的不断丰富，以及数字地形分析（DTA）方法和技术的不断发展，岩溶地貌的数字地形分析研究已取得一些重要成果，且DEM在岩溶地貌研究中得到充分应用，深化了对岩溶地貌形态、岩溶地貌单元空间分布格局、峰林峰丛单元自动提取、岩溶生态环境空间格局与地形统计等方面的认识。然而，受限于现今的岩溶地貌分类体系以抽象的语义表达、定性描述为主，难以具体量化岩溶地貌的空间结构信息、多尺度特征，其具有语义模糊性、范畴化，部分地貌类型还有交叠现象，岩溶地貌分类和边界不清晰，缺乏分类依据和指标体系，更难以适应数字化、智能化的需求。此外，岩溶地貌具有明显的区域性特征，热带岩溶发育以峰林平原、峰丛洼地为主，亚热带以溶丘谷地（洼地）、台原峡谷为代表，温带岩溶以干谷、浅平漏斗为代表，目前岩溶地貌的形态和定量表达方面，对于发育强烈的典型岩溶地貌如峰林、峰丛、洼地应用广泛，模型识别效果也较好，但对非典型性发育的岩溶地貌，数字地形分析鲜有尝试。

岩溶地貌及其生态系统的数字地形分析研究包含岩溶学、形态计量学、生态学、地理信息科学等多学科，具有较强的现实意义和应用价值。针对现有研究存在的问题，围绕岩溶地貌环境中数字分类体系构建、智能识别模型开发、多尺度形态表达、空间规律发现、地形综合属性挖掘等方面，提出未来的发展方向主要有：（1）对岩溶地貌和地貌要素构建更加系统科学高效的量化分类体系，以适应数字化、智能化的需求。一方面，基于地貌本体论研究，构建岩溶地貌知识组织模型，突破以往的通用地形元素分类模式，探索地形特征关联和融合的方式，借助机器学习技术，开发基于语义的岩溶地貌特征识别算法。另一方面，建立岩溶地貌数字指标体系，利用地貌形态参数化、空间纹理信息、景观格局指数分析方法，表征岩溶地貌单元及组合的类型特征、空间特征、地学属性特征的空间分布、关系、变化、规律；（2）研究不同岩溶地貌类型与地形因子的内在联系，分析不同岩溶区地貌参数的相似性与差异性，讨论数字高程模型与数字地形分析在不同岩溶区的适用性和应用示范，为岩溶区地形地貌特征的定量表达提供研究基础；（3）岩溶地貌形态类型具有典型的尺度和层级特点，需从区域、流域、坡面多级尺度下提取更加全面的地貌信息。在岩溶地貌群组空间分布格局上，从地形现状研究拓展到地形演变的过程与机理研究，以揭示各指标信息与岩溶地貌发育的内在联系，并进行发育机理与演化规律等研究。在岩溶地貌单元形态与表达上，随着Lidar、DOM等高精度数据的获取，能更加精细地刻画细节特征，有望深入分析岩溶地貌更加微观的形态，如石牙发育、落水洞识别等；（4）以DEM数据为基础，加强与遥感等多源数据的融合，丰富数字高程模型的地学属性内涵，挖掘出土壤、植被、地质构造等领域更深层次的信息，提高数字地貌服务于岩溶区生态环境治理和国土空间的适宜性评价。