基于CF-LR组合模型的地质灾害易发性评价——以云南省宁洱哈尼族彝族自治县为例

秦 红 富，谈 树 成，2，施 旖 奇，李 红 梅，汪 旭

(1.云南大学 地球科学学院，云南 昆明 650500； 2.云南省地理研究所，云南 昆明 650223)

0 引 言

地质灾害作为一种破坏性的地质事件，制约着人类社会的可持续性发展，对人类的生命财产和生存环境构成严重威胁，引发了社会各界的关注，已成为灾害学中的研究热点[1]。地质灾害多发生于坡度较大、断裂带密集、河网密布、土质疏松、岩土体脆弱、植被稀疏的斜坡上，此外强降雨和过度的人类活动也是滑坡、泥石流等地质灾害形成的因素[2]。中国地质灾害主要分布在西南、西北地区，该地区区域地质地理环境复杂多变使得地质灾害具有类型多、数量多、频率高、分布广、危害大等特点，其中尤以四川省、云南省为代表。采用合适的评价模型或方法计算地质灾害发生的可能性，并对地质灾害易发性结果分区，对区域防灾减灾具有重要意义。

国内外有关区域地质灾害易发性评价的研究工作已取得大量的成果。常用的易发性分析方法有层次分析法[3]、随机森林法[4]、信息量法[5]、支持向量机[6]、逻辑回归[7]、确定性系数法[8]、最大熵模型[9]、人工神经网络和决策树[10]等。结合各方法的优缺点，进行模型组合，可以有效提高模型精度，采用不同方法结合进行地质灾害易发性研究已经成为研究的热点和趋势。与此同时，GIS技术的飞速发展为进行地质灾害易发性区划深入研究提供了一个卓有成效的技术平台与研究途径[11]。

因此，以云南省宁洱哈尼族彝族自治县为研究区(简称“宁洱县”)，GIS为技术平台，采用确定性系数和逻辑回归模型(CF&LR)相结合的方法开展以滑坡、泥石流、崩塌等为主体的地质灾害易发性评价，以期为宁洱县的防灾减灾、区域防治、地质灾害风险预测等工作提供参考。

1 研究区概况

宁洱县地处普洱市中部，东经100°42′～101°37′，北纬22°40′～23°36′。研究区东西最大横距91 km，南北最大纵距101 km，辖6镇3乡，国土总面积3 669.77 km2，地处无量山脉西南部，属红河水系支流李仙江与澜沧江的分水岭地区，雨量充沛，年平均气温18.2 ℃，年平均降雨量为1 414.9 mm，夏季降水集中而冬季降水较少。海拔高差最大达2 300 m，相对高差一般500～1 000 m，总体属中浅切割的中山-低中山地形。研究区山地面积占总面积的96.77%，公路建设和矿山开发等人类社会工程活动对生态环境有一定影响和破坏，诱发了一些地质灾害。宁洱县地质灾害类型有滑坡、不稳定斜坡、崩塌、泥石流、地面塌陷等，是云南省滑坡、泥石流等地质灾害较为严重的地区之一。

2 材料和方法

2.1 数据来源与预处理

灾害点数据来源于项目组成员根据2018年高分辨率遥感影像和Google影像联合解译得到，并于2019年3月前完成野外调查验证，其中滑坡36处、崩塌38处、泥石流5处(见图1)；DEM数据来源于地理空间数据云平台，用于提取高程、坡度、坡向；降雨量数据来源于中国科学院资源环境科学数据中心；断裂带和岩性数据来源于1∶50 000中国地质图。河流、道路及基础地理数据来自国家基础地理信息中心。利用ArcGIS工具将所有的栅格数据的坐标统一为WGS_1984_UTM_Zone_47N，并对各个评价指标进行重分类。

图1 宁洱县地理位置及主要灾害点分布Fig.1 The geographical location and distribution of major geological disaster points in Ninger County

2.2 评价方法与评价单元的确定

2.2.1评价方法

(1) 确定性系数模型。CF模型是由Shortliffe[12]提出的函数，用来表示概率。后来，Heckerman[13]完善了该模型。CF 模型计算过程如下：

(1)

式中：ppa是地质灾害事件在地质环境因子a中发生的概率[0，1]，即地质环境因子a在特定单元中地灾点个数与地质环境因子a总面积的比值。pps是地质灾害在整个研究区中发生的先验概率，为地灾点总数与研究区总面积之比[14]。CF取值区间为[-1，1]，当CF大于0时，代表发生地质灾害可能性高，值越接近于1，发生地质灾害得可能性越高；当CF小于0时，代表发生地质灾害的可能性较低，值越接近于-1，发生地质灾害得可能性越低。

(2) 逻辑回归模型。在地质灾害易发性评价中将评价指标看作独立变量，地质灾害是否发生看作二元因变量(0指代地质灾害发生，1指代地质灾害不发生)[15]。逻辑回归模型描述的是二元因变量和独立变量xi(i=1，2，3….n)的关系。其中，独立变量不必满足正态分布，LR模型对数据类型要求较低，可以减少难以量化的数据对研究结果的影响，同时可以将研究步骤化繁为简。函数如下式：

(2)

式中：P表示地质灾害发生的概率，范围为[0，1]，值越大发生地质灾害的概率越大；α为常数项；xi(i=1，2，3，…，n)为因子i中各分类级别的CF值；β为各xi对应的回归系数[16]。

2.2.2评价单元的确定

合理的评价单元，决定数据的预处理格式和评价结果精度。总的来说，地质灾害易发性评价中较常用的评价单元有行政单元、规则栅格单元、斜坡单元等[17]。一方面，研究区地质地貌类型复杂多样、河流众多，降雨量充沛；另一方面考虑到比例尺和数据精度，以及研究区地质灾害主要以滑坡、崩塌、泥石流为主，本文采用规则栅格单元作为地质灾害易发性评价单元，目前较为常用的基于专家经验[18]的计算公式为

Gt=7.49+0.0006t-2.0×10-9t2+

2.9×10-15t3

(3)

式中：Gt指合适的栅格单元大小；t指等高线精度的分母值。研究区所使用的DEM数据为1∶50 000，则由式(3)计算得到合适的栅格单元大致为30 m，利用ArcGIS工具将研究区划分为4 075 195个30 m×30 m的评价栅格单元。

2.2.3评价因子的选取和分级

(1) 降雨量。在诱发地质灾害的众多因素当中，降雨发挥着十分重要的作用，为地质灾害发育提供了外动力条件。本文共收集近10 a宁洱县降水量空间插值数据，利用ArcGIS中的栅格计算器工具得到宁洱县年平均降雨量数据，可以发现近10 a宁洱县平均降雨量最小值为1 232 mm，最大值为1 798 mm。然后利用自然间断法将计算结果重分类分为1 230～1 380，1 380～1 458，1 458～1 537，1 537～1 628，1 628～1 800 mm共5个级别(见图2(a))。将降雨量栅格数据转换为面文件与地灾点数据建立空间连接，进行不同分级范围内的地质灾害数量统计，绘制直方图(见图3(a))分析不同分级降雨量与灾害点相对密度的关系。结果显示，降雨量在1 458～1 537 mm范围内，单位面积灾害点比例最大，说明该区域发生地质灾害的可能性较高。

(2) 距断裂带距离。断裂构造一直被认为是地质灾害发生的重要影响因素。在断裂带附近，发育地质灾害的斜坡常较为集中的发生，一旦受到诱发因素(降雨、人类工程活动等)的作用就更容易发生地质灾害。本次研究在ArcGIS中对宁洱县断裂带数据建立缓冲区，分别构建0～500，500～1 000，1 000～1 500，1 500～2 000，2 000～2 500，>2 500 m共6个等级(见图2(b))，并对不同等级范围的地灾数量进行统计。然后绘制直方图(见图3(b))分析距断裂带不同距离与灾害点相对密度的关系。结果显示，在距断裂带0～500 m和500～1 000 m范围内灾害点相对密度较其它分级范围大得多，研究区地质灾害主要集中在距断裂带0～1 000 m范围内。

(3) 距道路距离。在诱发地质灾害发生的众多因素当中，人类活动也起着十分重要的作用，包括切坡建路、采石采矿、乡镇建设等。一般来说仅考虑单因素对地质灾害的影响时，距道路距离越近相对较容易发生地质灾害。基于宁洱县乡道、省道、国道、高速矢量数据，利用ArcGIS中多环缓冲区工具对道路数据建立缓冲区，得到0～500，500～1 000，1 000～1 500，1 500～2 000，2 000～2 500，>2 500 m共6个等级范围(见图2(c))，并对不同等级地灾数量进行统计，并绘制直方图(见图3(c))分析距河流不同距离与灾害点相对密度之间的关系。结果显示：在距道路距离0～500 m范围内灾害点相对密度为0.034，仅次于距道路距离2 000～2 500 m范围的0.04，进一步证实了地质灾害的发生不仅受到距道路距离的影响还受到其他因素的影响。

(4) 距离河流距离。河流对地质灾害的影响在于对河流两岸斜坡前缘的侵蚀作用，增加了斜坡的临空面，导致河流边缘日益陡峭，破坏其稳定性。利用ArcGIS中多环缓冲区工具对河流矢量数据建立缓冲区，得到0～200，200～400，400～600，600～800，800～1 000，>1 000 m共6个等级范围(见图2(d))。并对不同等级地灾数量进行统计，并绘制直方图(见图3(d))分析距道路不同距离与灾害点相对密度之间的关系。结果显示：研究区地质灾害主要发生在距离河流0～200 m范围内，随着距河流距离的增加，灾害点相对密度逐级递减，说明距离河流越近越容易发生地质灾害。

(5) 岩石硬度。岩土体是地质灾害的物质组成，也是地质灾害发生的重要因素之一。岩石的类型和软硬程度不同，其抗风化和抗腐蚀能力不同，导致其抗滑动能力也不同。基于宁洱县地层矢量和岩土体坚硬程度，按照一般岩石坚硬程度分类标准，将岩石划分为坚硬岩类、较坚硬岩类、较软岩类、软岩类共4类(见图2(e))，统计不同岩石硬度范围灾害点的数量并绘制直方图(见图3(e))，分析不同岩石硬度与灾害点相对密度之间的关系。结果显示，较坚硬岩类的灾害点相对密度最大，但是仅有包含9个灾害点，主要是因为已有地质灾害点多分布于坚硬若和较软岩类交界处导致的，较软岩类的相对密度次之，但包含了69个地灾点，说明该区域地灾多发育于相对密度第的较软岩类，而不是相对密度最大的较坚硬岩类。

图2 各评价因子分级Fig.2 Grading diagram of each evaluation factor

图3 各评价因子与灾害点相对密度关系Fig.3 The relative density diagram of each evaluation factor and geological disaster points

(6) 归一化植被覆盖指数(NDVI)。植被影响着土壤的抗侵蚀能力和地表径流，进而影响岩土体。因此，NDVI也是地质灾害发生的影响因子之一。基于2018年的Landsat8数据，利用ArcGIS中的栅格计算器构建NDVI计算模型计算得到宁洱县NDVI数据。利用ArcGIS中的重分类工具将计算结果分为-0.2～0，0～0.2，0.2～0.4，0.4～0.6共5类级别(见图2(f))，统计不同NDVI分级范围灾害点的数量并绘制直方图分析灾害点相对密度与不同NDVI分级范围之间的关系。结果显示：NDVI值处于0～0.2之间灾害点相对密度值最大，较容易发生地质灾害。

(7) 高程。高程对地质灾害的作用，一方面是控制岩土体应力值的大小，岩土体的应力值和势能会随着所处高程的增高而增大，易在诱发因素的影响下发生地质灾害；另一方面是不同海拔高度的地区、气候、植被类型、降雨量、人类活动均有所不同，可能提供不同的孕灾环境。宁洱县整体高程相差较大，最低高程为534 m，最高为2 830 m。基于ArcGIS对DEM数据进行处理后，将高程按照534～1 000，1 000～1 450，1 450～1 900，1 900～2 350，2 350～2 830 m进行分级(见图2(g))，统计不同高程分级范围内灾害点数量并绘制直方图分析灾害点相对密度与不同高程之间的关系。结果显示：高程范围在1 450～1 900 m之间灾害点相对密度值最大，高程范围在1 000～1 450 m之间灾害点相对密度值次之，地质灾害主要发生在高程范围为1 000～1 900 m范围内，该范围内灾害点个数占总灾害点个数的93.67%，主要是因为在该高程范围内交通线路密集、人类活动频繁。

(8) 坡度。坡度一般通过决定斜坡体应力大小和方向，影响地表径流和物质移动等进而控制斜坡岩土体的稳定性。利用ArcGIS中的坡度工具，基于宁洱县DEM数据计算研究区坡度，得到坡度数据。结果显示，研究区坡度最大值为69.2°。并利用重分类工具对坡度数据进行分级，划分为0°～15°，15°～30°，30°～45°，45°～60°，60°～75°共5个坡度等级范围(见图2(h))，统计不同坡度分级范围内灾害点数量并绘制直方图(见图3(h))分析灾害点相对密度与不同坡度范围之间的关系。结果显示：在坡度范围为45°～75°范围内没有地灾的发生，坡度范围在30°～45°灾害点相对密度值最大地灾灾害较易发生，在坡度范围15°～30°内灾害点数量最多。

(9) 坡向。坡向对地质灾害的影响主要体现在不同坡向的斜坡受到的太阳辐射强度和光照时长有所区别，因此，同一山体的局部气候会有所不同，也就是温度、降雨量有所不同，进而影响研究区的植被覆盖度，岩石风化速率、土壤湿度等。研究表明，在相同岩土体条件下，阳坡面发生地灾的概率大于阴坡面，因为阳坡面具有岩石风化速率快、降雨量多等特点。基于宁洱县DEM数据，利用ArcGIS计算坡向，得到坡向数据，并对其进行划分，平地(-1°)、东北向(22.5°～67.5°)、正东向(67.5°～112.5°)、东南向(112.5°～157.5°)、正南向(157.5°～202.5°)、西南向(202.5°～247.5°)、正西向(247.5°～292.5°)、西北向(292.5°～337.5°)、正北向(337.5°～22.5°)共9个坡向等级范围(见图2(i))。统计不同坡向灾害点数量并绘制直方图(见图3(i))分析灾害点相对密度与不同坡向的关系。结果显示：平地(-1)内没有地质灾害发生，灾害点相对密度值最大的坡向为西南向(202.5°～247.5°)，西南向较其它坡向更容易发生地质灾害。

3 易发性评价流程

基于各评价指标因子的栅格数据利用CF模型计算出坡度、坡向、降雨量、岩石硬度、距河流距离、距道路距离、距断裂带距离、NDVI、高程9个评价指标因子不同分类级别的确定性系数即CF值，然后将9个评价指标因子各分类级别的CF值看作独立变量，将地质灾害是否发生看作因变量(0表示不发生，1表示发生)。利用python中的逻辑回归分析方法进行逻辑回归，得到各评价指标因子的回归系数，看作各指标的权重，同时进行显著性检验，求得逻辑回归方程。最后利用多值提取至点工具对研究区4 075 195个格网坡度、坡向、降雨量、岩石硬度、距河流距离、距道路距离、距断裂带距离、归一化植被指数(NDVI)、高程数据进行赋值并代入逻辑回归方程，得到研究区地质灾害易发性P值，最后利用等间距法进行分区即可得到研究区地质灾害易发性分区结果。

3.1 计算各评价指标的CF值

在对各评价指标不同分级范围的地质灾害进行统计分析后，根据2.2节中确定性系数模型，计算得到各评价指标因子不同分级范围的CF值，CF值越大表明地质灾害发生的确定性越高，反之越低。各评价指标确定性系数计算结果见表1。

表1 各评价指标确定性系数

由表1可以看出，各评价指标的CF值除了与各分级范围内的地质灾害数量有关还与各分级范围的面积有关。

(1) 从降雨量5个分级范围的CF值来看，确定性系数最大的范围是1 458～1 537 mm，CF值为0.344，确定性系数最小的范围是1 628～1 800 mm，CF值为-0.586。表明降雨量为1 458～1 537 mm范围时发生地质灾害的确定性最高，为研究区地质灾害发生的最佳孕灾环境之一。

(2) 从距断裂带距离的6个分级范围的CF值来看，0～500 m范围CF值为0.331，500～1 000 m范围CF值次之，表明距离断层距离越近地质灾害更容易发生。

(3) 从距道路距离的6个分级范围的CF值来看，2 000～2 500 m范围内CF值最大为0.467，0～500 m范围内CF值次之为0.373，在0～500 m范围内发生的地质灾害的数量要多于2 000～2 500 m范围，但0～500 m分级范围的分级面积大于2 000～2 500 m分级范围的分级面积，导致0～500 m范围内CF值小于2 000～2 500 m范围内CF值。表明2 000～2 500 m范围地质灾害的发生不只受到距道路距离的影响，还受到其他因素的共同作用。

(4) 从距河流距离的6个分级范围的CF值来看，0～200 m范围内CF值最大为0.487，表明距离河流越近越容易发生地质灾害。

(5) 从岩石硬度的4个分级范围的CF值来看，较坚硬岩类的CF值最大为0.304，较软岩类CF值次之为-0.023，这是因为较软岩类的分级面积为3 278.03 km2远大于较坚硬岩类的293.49 km2，故即使较软岩类发生地质灾害的数量最多达到69处，但确定性系数却较低。

(6) 从NDVI的4个分级范围的CF值来看，0～0.2范围内CF值最大为0.645，其余分级范围CF值为负值，表明植被覆盖越少的地方越容易发生地质灾害，反之越不容易发生地质灾害。

(7) 从高程的5个分级范围的CF值来看，在1 450～1 900 m范围内CF值最大为0.098，1 000～1 450 m范围内CF值次之，1 450～1 900 m范围内较其它分级范围更容易发生地质灾害。

(8) 从坡度的5个分级范围的CF值来看，坡度范围在30°～45°的CF值最大为0.047，其它坡度范围CF值均为负值，30°～45°范围内较其他分级范围更容易发生地质灾害。

(9) 从坡向的分级来看，坡向西南向(202.5°～247.5°)的CF值最大为0.525，西南向较其他坡向更容易发生地质灾害。

3.2 评价指标独立性检验

采用相关性分析进行评价指标之间的独立性分析，剔除相关性较大的指标因子，使用python求得各因子之间相关系数矩阵见表2。

表2 评价指标间的相关系数矩阵

表2中X1，X2，X3，X4，X5，X6，X7，X8，X9分别代表岩石硬度、距断裂带距离、距河流距离、降雨量、距道路距离、NDVI、坡度、坡向、高程。结果显示：各因子之间的相关系数绝对值都小于0.3，表明所选因子之间的相关性较小，9个评价指标因子都可以代入评价模型。

3.3 计算各评价指标权重

研究区一共有79个灾害点，在计算得到各评价指标的CF值后，随机选取研究区内79个没有地质灾害发生的点，利用ArcGIS中的多值提取至点工具和python工具将各评价指标的CF值赋值给这158个样本点，然后对这158个样本点基于逻辑回归模型来进行逻辑回归分析，结果见表3。

表3 逻辑回归分析结果

结合式(2)和表3的回归系数可以得到以下逻辑回归方程组：

(4)

式中：P为发生地质灾害的可能性，值区间为[0，1]；x1为岩石硬度中各分级范围的CF值；x2为距断裂带距离中各分级范围的CF值；x3为距河流距离中各分级范围的CF值；x4为降雨量中各分级范围的CF值；x5为距道路距离中各分级范围的CF值；x6为NDVI中各分级范围的CF值；x7为坡度中各分级范围的CF值；x8为坡向硬度中各分级范围的CF值；x9为高程中各分级范围的CF值。

由各评价指标因子的回归系数可知，代入评价模型的9个评价指标因子对研究区地质灾害敏感性由高到低依次为距河流距离、距断裂带距离、坡向、高程、NDVI、降雨量、坡度、岩石硬度、距道路距离。

4 结果分析与验证

4.1 评价结果

宁洱县地质灾害易发性评价结果在0.045～0.861之间。基于ArcGIS使用等间距法将易发性分区评价结果分为：极高易发区(0.657～0.861)，高易发区(0.453～0.657)，中易发区(0.249～0.453)，低易发区(0.045～0.249)4个等级(见图4)。再结合图3可知极高易发区主要分布在断裂带密度较大、路网密集、人类活动强度较大、年平均降雨量在145～1 537 mm之间的西南和中部区域。高易发区主要分布于研究区西部区域和南部区域。中易发区主要分布于研究区的东部区域以及北部部分区域。低易发区主要分布于研究区的北部区域和东南区域。结合灾害点的空间分布来看，灾害点越密集的区域，发生地质灾害的风险就越大，整体上，研究区地质灾害呈现北少南多、东少西多的特征。

4.2 评价结果验证

得到宁洱县评价结果后，通过两个方面对评价分区结果进行验证，一是通过现状灾害点空间分布对分区结果进行合理性检验；二是利用ROC曲线对结果进行准确性检验。

基于ArcGIS工具统计现状灾害点的空间分布结果判断分区结果的合理性及可靠性，对采用CF&LR组合模型得到的研究区易发性分区结果进行统计可得表4。结果显示：极高易发区分布灾害点35个，占灾害点总数的44.30%，分区面积占研究区面积的8.23%，单位面积灾害点比例最大；高易发区和中易发区均分布灾害点22个，占灾害点总数的27.85%，但单位面积灾害点密度高易发区较中易发区大；低易发区内无现状灾害点分布。现状地质灾害在极高易发区强发育，高易发区中等发育，中等易发区弱发育，低易发区内较弱发育。从验证结果来看，本次评价分区合理。

ROC曲线是一种不受临界约束的结果评价方法，能有效的对评价结果的准确性进行检验。ROC曲线以假阳性率(1-specificity)为横坐标，真阳性率(Sensitivity)为纵坐标绘制的曲线。ROC曲线下的面积为AUC值，是衡量模型准确性的指标，AUC的取值区间为[0.5，1]，当AUC=0.5时，模型预测无效；当AUC范围在[0.5，0.7]时，模型准确性较低；当AUC范围在[0.7，0.9]时，模型准确性较高；当AUC>0.9时模型准确性特别高。利用python工具计算得到评价结果的ROC曲线 (见图5)。结果显示：AUC值为0.74，表明基于CF&LR组合模型的宁洱县易发性评价结果准确性较高，可较为客观地对宁洱县地质灾害易发性进行分区评价。

图5 评价结果ROC的曲线Fig.5 Curve of evaluation results for ROC

5 结 论

(1) 研究区9个评价指标对地质灾害发生的敏感性依次为距河流距离、距断裂带距离、坡向、高程、NDVI、降雨量、坡度、岩石硬度、距道路距离。

(2) 研究区最易发生地质灾害的条件为距河流距离0～200 m、距断裂带距离0～500 m、坡向西南向(202.5°～247.5°)、高程1 450～1 900 m、NDVI范围在0～0.2、降雨量范围1 458～1 537 mm、坡度30°～35°、岩石硬度为较软岩类、距道路距离0～500 m。

(3) 极高易发区、高易发区、中易发区、低易发区分区面积分别占研究区面积的8.23%，37.17%，43.39%，11.21%，研究区地质灾害呈现北少南多、东少西多的特征。

(4) 对易发性分区结果进行合理性和准确性检验，从合理性检验来看，基于CF&LR组合模型的分区结果基本合理。从准确性检验来看，准确性检验结果显示AUC值为0.74，表明基于CF&LR组合模型的宁洱县易发性评价结果准确性较高。