基于maxEnt模型的哈尼梯田核心区滑坡易发性评价
2020-07-07赵冬梅角媛梅邱应美刘澄静徐秋娥
赵冬梅, 角媛梅, 邱应美, 刘澄静, 徐秋娥, 张 娟
(云南师范大学 旅游与地理科学学院, 昆明 650500)
滑坡是我国乃至世界范围内发生最频繁、最具危险的地质灾害之一[1]。云南省地处板块交界处,构造运动强烈,地形陡峻,切割破碎,加之局地强降水频繁,滑坡、泥石流等灾害易发高发,是造成生命财产损失、制约区域经济发展的主要灾种。被列为联合国教科文组织世界遗产名录的哈尼梯田作为云南省旅游胜地,其地质环境脆弱,人类活动影响强烈,易引发滑坡等地质灾害。该区因受连续强降雨的影响,遗产地内老虎嘴梯田片区在2018年6月26日发生了一起重大的山体滑坡,造成约11.5 hm2梯田受损,严重威胁哈尼梯田世界遗产地景观的可持续发展。因此,开展世界遗产地滑坡易发性对区域的防灾减灾及可持续发展具有重要意义。
滑坡易发性评估始于20世纪70年代中期[2],以滑坡灾害理论、形成机制等定性描述为主。近年来,涌现出大量的定量评估方法,可归纳为两类:一是以滑坡灾害机制为基础“白箱”型方法,二是以灾害监测历史数据为基础“黑箱”型方法。“白箱”型方法主要包括斜坡力学[3]和专家经验模型[4]。斜坡力学模型侧重单个滑坡体,其可靠性较高,但所需水文、岩土力学等物理参数难以获得,进而难以大范围推广[5]。然而专家经验模型主观性较强,使得模型精度具有很大不确定性。“黑箱”型方法主要包括概率统计、机器学习模型等。概率统计模型评价方法在中国、印度、意大利和土耳其等地运用广泛,如逻辑回归[6]、频率比、判别分析[7]、数据叠加、证据权重[8]、信息量模型[9]等,虽然一定程度上减少模型主观性且表现性能相对较好,但仍然存在:高精度数据难以获得;滑坡和地质环境资料需由专业学者初步处理;权重受研究区间划分的影响[10]等缺点。相比之下,最近学术界更偏向于机器学习模型,它是一组强大的数据驱动工具,使用算法上考虑了滑坡与滑坡影响因子之间的非线性关系[11],有效克服了概率统计模型的局限性。如交叉决策树、随机森林[12]、贝叶斯[13]、支持向量模型[14]与人工神经网络[15]等应用广泛,但截至目前,依然缺乏最佳模型评价滑坡易发性。
为探究最佳模型以准确评价滑坡易发性,本文将最大熵模型(maxent entropy Model,简称maxEnt)引入到滑坡易发性评价中。maxEnt模型是由Phillips等[16]开发的生态位模型,作者将其应用物种生境适宜性评价。主要是根据物种的不完整分布数据和环境特征计算其分布概率,从而达到预测该物种潜在分布的目的。与其他模型相比,模型具有操作简便、计算效率高、预测结果精确性和可信度高,可避免模型过拟合等优势,因而被生态学者广泛应用,现已成为生态学研究的热点工具。物种分布模式与滑坡易发性制图非常相似,都是利用已知事件对目标空间分布进行建模,建模过程涉及多个环境变量。国外已有学者尝试将maxEnt模型应用于滑坡易发性的研究[17-19],但国内引用较为鲜见。
鉴于此,研究以哈尼梯田遗产地为例,综合研究区滑坡灾害发育的地质环境、滑坡空间分布和发育特征等数据,选取气象、地形、地质、植被、人类活动等15个指标因子,应用maxEnt模型进行滑坡易发性评估,通过计算历史滑坡点落在较高危险性区域内的比例和接收器工作曲线下的面积AUC值定量检验。探究maxEnt模型在哈尼梯田核心区滑坡易发性评价中的适用性,并对其预测性能进行评估。根据统计结果确定预测危险性等级评价标准,由此对预测结果进行危险性等级划分,这对哈尼梯田心区的可持续发展、土地利用规划及灾害防治有一定的参考价值。
1 研究区概况
研究区位于我国西南部云南省元阳县,属红河哈尼梯田世界文化景观遗产核心区,面积为166.12 km2。其地理坐标22°49′—23°19′N,120°27′—103°13′E(图1)。该区地处哀牢山南段,为高山峡谷切割中山地貌,沟壑密布,切割强烈,地形陡峻,向源侵蚀较为严重。区内平均海拔2 000 m以上,梯田多分布在海拔400~1 900 m和>25°的陡坡上,且山脊走向为北东—南西向,整个地势南高北低,自南向北倾斜。气候属亚热带山地季风气候,素有“一山有四季,隔里不同天”之称。全年日照时数1 770.2 h,年均温14.2 ℃,雨量充沛,年降水量1 353.8 mm,年均湿度90.3%。土壤主要以黄棕壤、黄壤为主,土壤剖面发育完整。境内地层岩性主要为变质岩、砂岩、板岩等,在水作用下易于软化、崩解,且受长期复活的红河大断裂和哀牢山深大断裂的影响,新构造运动强烈,滑坡、崩塌等地质灾害频繁发生,滑坡是该区最为突出的地质灾害。此外,伴随着近年来哈尼梯田申遗成功,旅游业得以巨大发展,人类活动愈加频繁,加剧了滑坡生成,严重影响了哈尼梯田景观的稳定性及区域旅游的发展。
2 数据来源及研究方法
2.1 数据来源
2.1.1 易发性评价流程 滑坡易发性评价研究主要由3步骤组成。首先构建滑坡易发性指标体系,研究共选取了地质、地形、植被、气象和人类活动5类共15种环境变量;其次对所有环境变量进行预处理,统一转换格式、坐标系和像元值等,并将研究区内111个滑坡历史数据随机分成两部分,75%为训练数据,25%为验证数据;最后,将maxEnt模型引入哈尼梯田滑坡易发性评价,利用受试者工作特征曲线下面积ROC-AUC和滑坡密度指标分别对模型精度及分区结果进行检验。
图1 云南省红河哈尼梯田核心区位置及滑坡点分布
2.1.2 滑坡数据 基于哈尼梯田遗产核心区地质灾害实地调查资料并结合Google地球上获取的遥感影像数据,研究区共确定滑坡111处。其中,70%以上滑坡厚度小于5 m,80%以上属浅层滑坡。区内滑坡主要分布于道路两旁,且多发生于雨季和人为工程活动强烈的区域,最后将滑坡数据随机分为两组数据,75%的数据即84个滑坡点用做训练样本进行建模,剩余的25%即27个滑坡点作为测试样本对模型进行验证。设模型运行参数中的迭代次数为500次,取500次模拟结果的平均值作为最终模拟结果。
2.1.3 指标因子构建及处理 选取了气象、地形、植被、地质、人类活动等共15个影响因子作为易发性评价指标体系。各因子来源如下:(1) 地形因子使用地理空间数据云上的数字高程模型DEM栅格数据,精度为30 m,并从中提取了7个地形因子:海拔、坡度、坡向、地表粗糙度、平面和剖面曲率和地形湿度指数(topographic wetness index,TWI)(附图6-Ⅰ A—G)。(2) 地质数据经云南省地质环境监测院提供的1∶5万地质图数字化获得,包括岩性、断层(附图6-Ⅰ L,K),岩性主要包括六类:1为泥质灰岩,2为泥质粉砂岩、灰岩、白云质灰岩、白云岩,3为片麻岩及角闪斜长片麻岩互组,4为粉砂岩、石英砂岩、夹板岩,5为灰岩,6为长石石英砂岩。(3) 距道路、水系、断层距离3个指标(附图6-Ⅰ I—K)在ArcGIS 10.4中通过欧式距离计算得到这些因子的距离图层。(4) 多年平均降水数据(AMR)由黑河中心提供的2008—2014年的降雨量数据经空间分析克里金插值得到(附图6-Ⅱ O)。(5) 植被归一化指数(NDVI)由分辨率为30 m的地理空间数据云Landsat8遥感影像数据经过大气矫正、波段运算获得(附图6-Ⅱ M)。土地利用数据是云南省地理信息中心提供的0.6 m的高分辨率遥感数据采用监督分类的方法提取的,根据国家最新发布的《土地利用现状分类》(GB/T21010—2017),将土地利用划分为8类,分别是旱地、林地、水田、草地、园地、水域、建筑用地和其他用地(附图6-Ⅱ N)。(6) 居民点密度由云南省第一次地理国情普查数据库矢量图经ArcGIS软件的核密度分析获得(附图6-Ⅰ H)。本研究所有环境变量均在ESRI公司的ArcGIS 10.4软件中完成,统一到WGS_1984_UTM_Zone_48N投影坐标系,并以空间分辨率10 m×10 m的ASCⅡ grid格式的形式输入软件。将所有的环境变量数据通过ArcGIS 10.4将格式转换为ASC格式,并结合滑坡数据导入maxEnt Version 3.3.3 K软件中,采用jackknife刀切法计算各气候因子贡献率,其余选项采用模型默认设置。
2.2 研究方法
2.2.1 maxEnt模型 maxEnt模型是基于已知的滑坡灾害分布数据以及地质环境条件,通过特定算法求出区域内滑坡灾害发生在特定位置上的可能分布π,进行未知区域的滑坡分布概率预测[16,20]。具体是将研究区划分为有限个像元集X,设x表示研究区域上的随机点x∈X,π(x)是每个像元π发生滑坡的一个非负概率分布值P(x),且所有像元的概率值之和为。将任意像元的滑坡观测结果视为响应变量y。若像元内“发生”滑坡则定义为y=1,“不发生”滑坡时则定义为y=0。由此通过使用贝叶斯规则得到条件分布概率P(y|x)[16,21]。因此,滑坡发生的概率可表示为:
(1)
式中:P(y=1|x)为在特定x点发生山体滑坡的概率;P(y|x=1)为滑坡分布条件下特定点x发生滑坡的可能性,亦是π(x);P(y=1)为滑坡总体发生率;P(x)为任意点发生x滑坡的概率。由于P(x)在研究区内所有像元X的任意点x中等于1/|X|,因此上述等式可以改写为:
P(y=1|x)=π(x)P(y=1)|X|
(2)
可通过边缘化联合概率分布来计算P(x):
(3)
考虑滑坡发生和不发生的概率相等[P(y=0)=P(y=1) =0.5],可简化方程为:
(4)
maxEnt模型的应用直接取决于条件概率P(y=1|x),条件概率值越大,滑坡发生的可能性越大。事件发生的数据π(x)可用于建模,代替直接估算P(y=1|x)。最大熵原理估计的π(x)值等于由指数表示的吉布斯概率分布。如果考虑n个特征(f1,i=1,2,…,n),则吉布斯的概率分布可定义为:
(5)
式中:Zλ是一个归一化常数可确保qλ(x)和1,而是分配给特征的权重向量。在qλ(x)的估计中,模型试图利用正则化I2找到最接近约束条件下的分布,以避免过度拟合。因此,maxEnt模型的目的是最大限度地处理对数似然。如果研究区内事件出现m次,则对数似然与正则化之间的差异应为最大化,表示为:
(6)
式中:βj是jth特性fj的正则化参数。最大熵模型发现吉布斯分布不仅符合赋存数据,而且具有较好的推广意义。
2.2.2 模型精确性验证 模型模拟的精度验证是保证预测结果准确度的必要步骤。一般预测模型通常会产生两类误差:一是过低估算,将实际有滑坡发生的区域预测为低敏感区,为假阴性;二是过高估算,将实际不发生滑坡的区域预测为高敏感区,为假阳性,这两类误差都与阈值有关。受试者工作特征曲线下的面积(ROC-AUC)分析法是目前学者们认可度较高的模型诊断试验评价指标。该分析方法以灵敏度(sensitivity,也称作真阳性,1-遗漏率)为纵坐标,以特异度(specificity,也称作假阳性)为横坐标,绘制而成的曲线称之为ROC曲线。通常AUC-ROC值越大,表示指标因子与预测模型结果之间相关性越大,预测效果也就越好,其能很好地说明模型模拟的准确性。其评价标准为:当AUC-ROC的值为0.5~0.6时预测效果失败(fail),值在0.6~0.7预测效果较差(poor),0.7~0.8预测效果一般(moderate),0.8~0.9预测效果较好(good),0.9~1.0预测效果非常好(Excellent)[17]。
2.2.3 环境变量重要性评价 利用刀切法(Jacknife)对模型的因子贡献进行了估计。在这种方法中,每个环境变量都被有意排除在外,并且使用剩余的因素构建了一个模型。然后将使用所有环境变量创建的模型与使用依次排除一个变量构建的模型的预测性能进行比较。因此,可以审查排除的环境变量贡献值。还利用响应曲线推导了各因果因子与预测模型之间的关系。
3 结果与分析
3.1 预测结果检验
研究所用的训练样本和测试样本生成的ROC曲线如图2所示,模型运算结果表明训练样本下的ROC-AUC值达0.895,验证样本下的ROC-AUC值达0.739,显然两者值均大于0.7,远大于随机测试的ROC曲线下面积值(ROC-AUC=0.5),且ROC-AUC标准差为0.047,非正规训练增益为1.201 5,未经规范的测试增益为0.319。根据分类标准,此模型在滑坡易发性中的预测精度达到优秀水平,表明模型模拟效果具有很高的可信度和准确度,可用于滑坡易发性预测。
图2 模型精度检验ROC曲线
3.2 环境变量因子对模型的贡献度
本研究采用Jackknife检验模型中15个环境变量对模型预测的贡献程度表明(表1),所有环境变量的贡献率都大于0,贡献率大于1%的因子9个,排名前7的环境因子依次是距断层距离(26.80%)、多年平均降雨量(17.33%)、距道路距离(16.42%)、土地利用(14.91%)、居民点密度(4.80%)、坡向(4.67%)和岩性(4.5%),累计贡献率达到了89%,这7个因子对滑坡分布具有重要影响,其次是距水系距离、地表粗糙度、海拔、NDVI、坡度和TWI,平面和剖面曲率的影响最小。
3.3 环境变量对maxEnt模型预测结果的影响
图3为15个环境变量响应曲线,各环境因子的易发性分析表明,滑坡易发性在1 500~2 500 m出现峰值,在2 200 m时达到最大,之后随海拔的增加而减小,这与当地居民居住的海拔范围相同,另外居民点密度显示人类活动对滑坡的影响巨大。对坡度、地表粗糙度和坡向进行易发性分析发现随着坡度的增加,灵敏度值越大,表明滑坡发生的可能性和坡度呈正比,在坡度较大的地段,土壤与岩石之间十分不牢固,梯田灌水后将增大下滑动力,促使滑坡发生。坡向统计发现坡向为东北、西北、东坡时对滑坡灾害的影响最大,而北坡影响最小。从曲率值看,滑坡的发生可能性随着平面、剖面曲率值的增加而减少,呈负相关关系。地形粗糙度来看值为1.0时对滑坡的影响最大。
表1 环境变量的贡献率及对应环境变量的ROC-AUC值
从TWI来看,滑坡在值0~5呈增长趋势,值为5时达到峰值,表明此时对滑坡影响最大。居民点密度越大,表明越容易发生滑坡。从距离因子上看,距离道路越近时滑坡易发性越高,特别是在0~500 m影响最大。距水系0~300 m也表现出一定的距离效应,但当距离大于300 m后滑坡影响逐渐减少。当断层距离在4 000~6 000 m时对滑坡的影响最大。从岩性的不同类别上看,片麻岩与角闪斜长片麻岩对滑坡的影响最大,这类岩石在水作用下极易风化、崩解,形成粘性的红壤、黄壤等,这类土壤土质紧密,遇水泥泞,粘性很强,但透水性差,雨水不容易下渗,常导致滑坡产生。NDVI在-0.05~0.22时滑坡随植被指数的增加呈现缓慢的增加趋势,而大于0.22时又随植被指数的增加而减少,直至为0表现为几乎无影响的状态。从土地利用类型上看,林地对于滑坡易发性影响较大,水田影响相对较小,一定程度上梯田确实减缓了地表径流速度,减小了水土流失量,但其分布于山腰居多,其山体的坡度普遍较陡,会造成田埂垮塌,从而引起滑坡的发生。从AMR上看,滑坡随降雨量的先增加影响也在增加,在920 mm时出现峰值,后较为平稳,小于820 mm时影响较小,长时间降雨会导致水压力过大,将有增大滑坡灾害发生的可能性。据调查,研究区内近80%的滑坡发生在雨季,降雨成为诱发滑坡的主要因子。
3.4 滑坡易发性分区
模型预测结果值越大,滑坡发生的可能性越大。根据正态分布理论与专家经验法,将研究区模拟的滑坡易发性结果划分为5个等级:极低易发区、低易发区、中易发区、高易发区和极高易发区。maxEnt模型模拟的滑坡易发性指数值为0~1,其分级标准为:极低(0~0.13)、低(0.14~0.31)、中(0.32~0.51)、高(0.52~0.73)和极高(0.73~1)。
统计显示(附图7,图4)极低易发区占比最高,达到37.24%,低易发区占23.21%,其他中、高和极高易发区分别占17.78%,12.79%,8.96%。易发性区划图空间差异明显,极高和高易发区的占比不大。实际存在滑坡点在易发区类别上看,高和极高易发区的实际占比分别为27.92%,54.05%,极高和高敏感区主要分布核心区中部道路两旁,即人类居住集中区域,以中部偏东地区最为严重。据调查,中东部的道路是在哈尼梯田申遗后期修筑的,修建道路,致使边坡失稳,加大了边坡发生滑坡的可能性。研究区西南部也有零散分布,集中于距断裂较近处,该区由于受哀牢山的大断裂强烈影响,断裂带沿线岩体破碎,岩层节理裂隙发育,进而造成滑坡呈条带状的格局。
此外,研究还利用滑坡密度指标对滑坡易发性分区结果进行检验,它是通过滑坡点占比与各易发性等级占比相除计算得出,前人研究显示滑坡易发性等级越高滑坡密度越大,说明模型越可靠[22]。因此,根据统计得到滑坡易发性分区下对应的5个滑坡密度值,可以发现滑坡密度随着滑坡易发性的等级的升高而增加,且对应易发性区域值分别为0.1,0.27,0.46,2.18,6.03,符合以上事实,这再次证明在滑坡数据量较少的情况下,maxEnt模型模拟结果仍可信。
图3 环境变量的响应曲线
图4 基于自然断点法模型的滑坡易发性评价区划
4 讨 论
目前,滑坡易发性的研究已涌进了大量的评估方法。受调查程度、地区条件、使用者水平等的影响,不同地区方法不同时,预测效果不同;同一地区方法不同,预测性能也有所差别,学术界针对不同的模型对其方法进行了大量的比较分析并确定最佳模型,但在模型精度上并没有达成一致意见。这些模型不仅仅运用在地质灾害的评估中,大多也被运用到生物物种的生境适宜性分析、潜在矿物分布、地面沉降测绘及医学传染病等研究,注重不同学科模型及方法的交叉及多重使用,特别是物种分布模型与滑坡易发性研究的原理研究类似,都是利用已知的事物来模拟未来事物的分布。因此,研究引用了生态位模型最为常用的maxEnt模型,它也称为人工挖掘模型来评价哈尼梯田区滑坡易发性。
该模型与滑坡易发性相关的研究国内几乎没有,国际上较为少见。最早将maxEnt模型引入此研究是2013年,Felicísimo等[18]选择了20个环境变量和1 102个滑坡点将多变量回归模型、多元自适应回归、分类回归树和maxEnt模型运用在西班牙吉普祖科亚省德巴河谷,结果显示maxEnt模型模拟精度高于其他模型。Kornejady等[19]运用maxEnt模型模拟不同样点情景下的滑坡易发性,结果表明马氏距离样点方案下模拟结果精度达0.906,高于随机情景下模拟精度。以及Jiao等[23]将其与常用的信息量模型进行对比研究发现maxEnt模型的精度达到了优秀水平。相较于前人研究而言,本研究运用maxEnt模型的模拟精度相差不大,ROC-AUC值达到0.895,模拟效果较好,研究结果可信。研究再次证明了maxEnt模型具有较高的预测精度这一说法。
同时,在影响因子的选择上看,本文不仅考虑了气候条件对滑坡的影响,该研究还将人类活动因子作为当前和未来发生滑坡的重要驱动因素,从而增加了关于气候与人为干扰的讨论。前人滑坡易发性研究中常见的变量有地形、距离、岩性、土地利用等因子,对居民点密度这个变量考虑较少,但此类变量在本研究的贡献程度为4.8%,影响程度较大(表1),是哈尼梯田遗产核心区滑坡灾害产生的主要控制因子,缺少此因子时模拟结果精度ROC-AUC为0.69,较坡向、地表粗糙度、平面和剖面曲率、距水系距离、坡度、NDVI和TWI高。因此,在人类活动强度较高的区域应当充分考虑居民点对滑坡灾害的影响。此外,研究区内AMR单个因子的影响程度最大(表1),降雨成为该区滑坡发生的一个主要因子。有研究表明受区域持续强降雨的影响,滑坡等地质灾害产生的概率将会提高,因此,在未来的研究过程中,结合当地实际情况,充分利用遥感高分辨数据及ArcGIS等技术手段,多时段、多时间序列的观测区域的环境,将自然和人类活动等因子考虑到研究领域中,以便更加精确地预测滑坡等地质灾害。本研究存在一些不足之处,如使用的滑坡数据点通过高分辨率的Google遥感影像目视解译的,可能存在一定的遗漏,难免会对模拟结果产生一定的误差,但整体上模拟的结果与实际存在的数据大致相符,因此模型预测结果本身具有一定的参考价值,得到的结果更能体现由人为活动贡献较强的情况下滑坡发生的概率,所以本研究既具有学科意义,也对哈尼梯田遗产核心区的土地利用规划及灾害防治及减灾有一定的指导作用。
5 结 论
(1) 以哈尼梯田世界遗产地为研究对象,选取海拔等15个因子建立滑坡易发性评价体系。基于GIS和maxEnt模型进行区域滑坡易发性评价,所得结果与实际滑坡分布情况较为一致,ROC-AUC值达到了0.895,表明模型模拟效果具有很高的可信度,在滑坡易发性评价中适用。且通过滑坡密度可以看出模型具有较强的准确度。
(2) 根据Jackknife方法对模型预测的贡献程度表明距断层距离、多年平均降雨量、距道路距离和居民点密度是影响遗产地滑坡发育的主要影响因子。环境响应曲线也揭示了各变量的等级状态下滑坡分布的统计规律,对分类及连续型变量具有很好的适用性。
(3) 使用正态分布理论与专家经验法将滑坡易发性分成5类,滑坡极高和高易发区主要分布核心区中部道路两旁,即人类居住较为集中区域,中部偏东最为严重。而研究区内西南部集中分布于距断裂较近处,呈条带状分布格局的特点。