基于加权信息量法的黄土滑坡易发性评价
——以1∶5万天水市麦积幅为例
2022-06-24孟晓捷张新社曾庆铭王冬
孟晓捷,张新社,曾庆铭,王冬
(自然资源部黄土地质灾害重点实验室,中国地质调查局西安地质调查中心,陕西 西安 710054)
据统计,中国的滑坡约1/3发生在黄土高原地区,甘肃省即存在各类规模的滑坡约4万余处(吴玮江等,2006)。黄土滑坡作为山区常见的突发性地质现象,常常埋没村庄、破坏农田、毁坏道路、摧毁工厂和堵塞江河等,还会进一步诱发泥石流等次生危害;长期以来威胁着黄土高原广大人民群众的生命财产安全,制约着国民经济发展,也为近年来的黄河流域生态保护战略计划和国土空间规划管控工作带来了不小的挑战。因此,针对黄土滑坡地质灾害进行深入调查和易发性评价,可为区域地质灾害防治及政府规划管理等提供科学依据,具有重要的现实意义。
滑坡易发区指具备滑坡发生的地质构造、地形地貌和气候条件,容易或者可能发生地质灾害的区域,主要依据地质环境条件,参考地质灾害现状,考虑潜在的滑坡而划定(唐亚明,2012)。近二十年来,随着3S技术(GIS、RS、GPS)的迅速发展,国内外学者将众多数学模型(证据权重模型、确定系数法模型、层次分析法模型、信息量模型、模糊评判模型、BP神经网络模型、逻辑回归模型等)引入滑坡易发性评价。基于GIS技术广泛应用于中国各个地质灾害高发区域,极大地提高了评价效率和精度,取得了丰硕的成果(沈芳等,1999;向喜琼等,2000;Ohlmacher,2003;赵海卿等,2004;王哲,2009;Guo C B, et al.,2015;Zhang,2017;张晓东等,2018;吴常润等,2020)。近年来,由于各类评价方法自身局限性及对滑坡灾害调查数据的精确性要求不同,研究者逐渐开始重视将不同的模型进行组合使用,选取符合研究区的最优评价模型(牛全福等,2017;张晓东等,2018;樊芷吟等,2018;邱维蓉等,2020;王高峰等,2021)。因此,笔者在天水市麦积区1∶5万图幅黄土滑坡地质灾害调查的基础上,基于ArcGIS平台,针对麦积区黄土滑坡地质灾害,采用层次分析法求得各评价因子的权重值,随后利用加权信息量模型,按照网格划分进行易发性评价,将评价结果进行验证,并划分天水市麦积区地质灾害易发性分区。旨在为该地区进行地质灾害防治规划、重大工程规划建设提供科学依据,同时也为实现区域地质灾害风险评价奠定良好基础。
1 研究区与研究方法
1.1 研究区地质概况
天水市麦积区位于青藏高原东北缘,大地构造上位于南北活动构造带与西秦岭构造带的交汇部位,周边发育多条强活动断裂(图1),历史上地震相当频繁。另外,在地貌上位于秦岭山区、陇山山区、陕北-陇东黄土高原及青藏高原的过渡地区。
图1 研究区位置图
麦积区幅位于天水市东侧,地理坐标为东经105°45′—106°00′,北纬34°30′—34°40′,总面积约为438.73 km2;属温带大陆性季风气候,多年平均气温为9.0 ℃;多年平均降水量为516.9 mm,降水量主要集中在5月—9月,占全年降水量的80%;区内河流属于黄河流域水系,发育有渭河、籍河、牛头河、东柯河和颍川河等多条河流(图2)。区内地势最高点位于西南基岩山区,海拔1 923 m,最低点位于渭河东缘最下游,海拔1 061 m,相对高差862 m。区内地貌单元可分为河谷地貌、低山地貌及黄土丘陵地貌(图2),其中纵横起伏、千沟万壑的黄土丘陵地貌总面积为363.84 km2,占全区面积的82.93%,加上区内多条活动断裂发育(图2),地质构造复杂,新构造运动强烈,为地质灾害的高发易发创造了客观地质条件。因此,天水市自古以来饱受滑坡之害,是中国滑坡灾害多发的城市之一。
F1-1.陈家湾-王家崖断裂;F1-2.李家湾-杨家湾断裂;F2-1.三阳川断裂;F2-2.凤凰山-麦积断裂;F2-3.水眼寨断裂;F3-1.孙家山断裂;F4.曹家沟断裂;F5.孙家沟断裂
项目组在前人调查的基础上(杨为民等,2016),进一步对麦积区幅进行地质灾害详细调查、核实,区内滑坡共计460处,其中有458处发生于黄土丘陵地区,占99.56%,绝大多数为老滑坡和新滑坡,所占比例分别为61.5%、37.4%。滑坡多为中-深层牵引式滑坡,体积多为中-大型,滑坡面多位于上覆黄土与下伏泥岩的接触带。大多数滑坡目前基本稳定,活动迹象不明显,但在将来的地震或强降雨等因子的诱发下,可能会再次复活。
1.2 研究方法
“信息量法”最早由美国数学家Shannon提出,并在概率论和逻辑推理基础上,推导出信息量法计算公式(Shannon et al.,1948)。20世纪80年代以来,信息量法被引入至滑坡灾害预测研究中,结合GIS技术,可快速动态地对研究区做出分析评价。其核心思想为:通过研究已发生或已变形的地质灾害地区实际情况,分析相关影响因子,并研究各影响因子的信息数量和质量,选取对地质灾害易发性影响最大的相关组合因子,反映不同因子对地质灾害的形成发生的贡献大小。其模型计算公式(这里省略推导过程)为:
(1)
式中Ii指各影响因子对滑坡提供的信息量;N指研究区内已发生滑坡的总个数;S指研究区总面积;Ni指分布在各因子类别中发生滑坡的个数;Si指研究区内某评价因子的面积。最终计算出的Ii为该影响因子的总信息量值,其值越大表示越容易引起滑坡的发生,即滑坡的易发性越高(孟庆华,2011)。
各学者将信息量模型广泛运用于中国诸多地质灾害易发高发区域,取得了显著的成果,其优势在于能够结合GIS技术快速计算出各影响因子的总信息量值。但该模型存在一定的局限性,地质灾害发生原因复杂,各影响因子在不同地质环境中的灾害发生过程里起到的作用大小不一,该模型只反映了因子特定类别在组合情况下对灾害发生的影响,不能充分反映各因子影响程度的差异(杜国梁等,2016)。
因此,笔者拟采用层次分析法(这里不再赘述该原理)求得各影响因子的不同权重,体现各因子的重要程度,结合加权信息量模型来进行易发性评价,其计算公式为:
(2)
其中Wi为层次分析法计算出的某影响因子的权重。
1.3 数据来源
本研究所用的地质灾害基础调查数据源于中国地质科学院地质力学研究所于2016年所完成的《西秦岭北缘断裂带天水段地震地质灾害调查评价成果报告》(杨为民等,2016)。随后由中国地质调查局下达的“关中平原城市群综合地质调查”二级项目于天水市麦积区进行核查,其统计结果包括滑坡灾害点460处。此外,笔者进行滑坡调查及易发性评价过程中所使用的基础数据包括DEM数据(25 m分辨率)和1∶5万地形图基础数据(等高线、高程点、水系、道路等)源自自然资源部信息中心;植被指数通过中国科学院资源环境科学数据中心数据注册与出版系统下载获取。
2 评价因子的选择与分级
根据分析各个因子对研究区地质灾害发生的控制结果,共确定7个指标因子(高程、坡度、地层岩性、植被归一化指数、距断裂距离、距水系距离和距道路距离)组成本次滑坡灾害易发生评价的指标体系。另外,笔者需要补充的是,经现场调查发现,几乎所有滑坡发生前后其坡形有着较大的变化,现今多为凹型坡或复合型坡,难以推测判断原始坡形。经对滑坡分布进行分析,该区460个滑坡在各方向分布上无明显数量差异,阳坡和阴坡上分布较为平均。另外,由于地形起伏度难以获取,在ARCGIS软件中随着计算邻域半径的变化,求出的地形起伏度有着较大的变化,不利于准确进行易发性评价。本研究区地震动峰值加速度均为0.3 g,无法体现不同等级地震动峰值加速度的影响。研究区内及周边有天水市、麦积区、清水县和张家川县等4个气象站点,降雨量数据难以完整获取。因此,前人常用的“坡形或曲率”、“坡向”、“地形起伏度”、“地震动峰值加速度”和“降雨量”等因子在这里不做考虑。
2.1 高程
麦积区黄土丘陵及基岩山区沟壑纵横,切割强烈,地形相对高差较大,对黄土滑坡的发生有着较为重要的影响;研究区内地形高程范围为1 060.78~1 923.21 m,将其分成5个区间,分别为≤1 200 m,1 200~1 350 m,1 350~1 500 m,1 500~1 650 m,>1 650 m。据统计,其中435处滑坡,即94.57%的滑坡发生于高程1 200~1 650 m区间(图3a)。
2.2 坡度
斜坡坡度控制着斜坡的稳定性,对滑坡的破坏方式和运移有着显著影响。在其他地质背景相似的条件下,通常斜坡坡度越陡,坡体越倾向于失稳,在一定角度范围区间,破坏以滑坡为主,超过一定角度,破坏以崩塌为主。研究区内坡度范围为0°~61°,将其分为5个区间,分别为0°~10°、10°~20°、20°~30°、30°~40°、>40°。据统计,其中有370个滑坡,即80.43%的滑坡发生于10°~30°区间(图3b)。
2.3 工程地质岩组
研究区内的工程地质岩组按照区域分布面积由大到小依次为:松散岩组>软弱岩组>较硬岩组>坚硬岩组>较软岩组(图3c)。松散岩组广泛分布于全区,岩性主要以第四系黄土和冲洪积物为主;软弱岩组常下伏于黄土沉积底部作为基底,岩性以新近系杂色泥岩为主,在较高斜坡处常有出露,整体产状近水平;较软岩组主要分布于牛头河北段西岸,岩性主要为古近系和侏罗系砂砾岩,整体产状近水平;较硬岩组主要分布于牛头河和渭河两岸岸坡,岩性主要由下古生界牛头河群的大理岩、片麻岩组成;坚硬岩组少量分布于渭河两岸局部地区,大部分分布于牛头河北段两岸,岩性主要由侵入岩、下古生界牛头河群石英片岩、角闪片岩、花岗岩和闪长岩组成。据统计,其中有450处滑坡,即97.83%的滑坡发生于松散岩组和软弱岩组地区,其中松散岩组地区有410处滑坡,软弱岩组地区有40处。
a.滑坡分布与高程;b.滑坡分布与坡度;c.滑坡分布与工程地质岩组;d.滑坡分布与植被归一化指数;e.滑坡分布与距断裂距离;f.滑坡分布与距水系距离;g.滑坡分布与距公路距离
2.4 植被归一化指数
植被能起到护坡和防止水土流失的作用,有利于斜坡的稳定。植被覆盖程度可用植被指数来表述(薛强等,2015),本文选取标准化差值植被指数(NDVI)作为分析指标,使用Spot5遥感数据进行标准化差值计算,求取植被指数NDVI,随后按自然断点法划分为0.08~0.22、0.22~0.31、0.31~0.37、0.37~0.42、0.42~0.49、0.49~0.66等6个区间(图3d)。前2个区间NDVI值较低,主要分布于盆地内城市区域及河谷两岸主要交通干道。据统计,其中有424处滑坡,即约92.17%的滑坡发生于NDVI值为0.31~0.66的区间。
2.5 距断裂距离
断裂破碎带岩体破碎,裂隙发育,易风化剥蚀,沿断裂带常会出现负地形,处于斜坡地带的断裂带会形成大量的松散层堆积,为滑坡形成创造了有利条件。其次,断裂带附近构造应力集中,在黄土丘陵区的黄土中易产生构造裂隙,在降雨等条件影响下,发育成较大的软弱面,容易发生大规模的破坏作用。再次,构造活动区地壳沿断裂带发生区域性抬升和下降,当地壳抬升时,河流下切作用增强,导致斜坡前缘临空面发育,加大了崩塌、滑坡的易发性(杨为民等,2016)。
研究区内断裂较为发育,按其走向可分为北西西向、北西向、北东向3组。其中,北西西向为西秦岭北缘断裂带的分支断裂,断裂特征明显,活动性较强,控制了天水渭河盆地的发育,主要包括2条断裂:F1为天水-宝鸡断裂,F2为凤凰山-天水断裂(图2)。其中,F1断裂由2条次级断裂组成,F2断裂由3条次级断裂组成,均为左旋走滑断裂。北西向为研究区西南角的F3(孙家山断裂),在遥感影像上有清楚显示,区域上属于左行走滑断裂。北东向为F4(曹家沟断裂)和F5(孙家沟断裂),沿线多为村落和庄稼地,村内常见墙体、地面开裂。曹家沟、孙家沟两侧斜坡及盘山公路上大量发育滑坡和地裂缝,可见活动性较为显著。将研究区滑坡距断裂的距离分为<1 500 m、1 500~3 000 m、3 000~4 500 m、4 500~6 000 m、>6 000 m等5个区间(图3e)。据统计,发育在距断裂1 500 m以内的滑坡数量为277处,占全滑坡的60.22%,随着距断裂距离的增加,滑坡数量显著减少,可见断裂对滑坡发育的控制作用。
2.6 距水系距离
河流侵蚀与水对坡脚岩土体的浸泡软化及河流下切形成的斜坡临空面也会对斜坡岩土体稳定性造成影响。将研究区滑坡距水系的距离划分为<250 m、250~500 m、500~750 m、750~1 000 m、1 000~1 250 m、1 250~1 500 m等6个区间(图3f)。据统计,距水系1 000 m以内的滑坡数量为434处,占全滑坡的94.35%。
2.7 距公路距离
人类工程活动中的修建公路对斜坡稳定性影响较大,修路对边坡的开挖破坏了原始斜坡的自然平衡状态,是地质灾害重要的诱发因子之一。对研究区内主干公路两侧500 m范围内进行距离分析,其中盆地内河谷区无地质灾害分布,不纳入统计范围。另外,黄土丘陵区内部没有位于斜坡地带、未改造边坡形态的乡间土路也不纳入统计。将区内滑坡距道路的距离分为<100 m、100~200 m、200~300 m、300~400 m、400~500 m、>500 m等6个区间(图3g)。据统计,距道路500 m以内,滑坡数量为389处,占全部滑坡的84.57%。
3 基于GIS的滑坡易发性评价
3.1 单因子信息量计算
基于上一章节选择的各个评价因子及分区进行信息量计算,得到各评价因子信息量值(表1)。
3.2 权重值计算
本研究拟采用层次分析法(AHP)求得上述7个评价因子指标的权重值;该方法可以解决多指标定性、定量组合决策的复杂问题。其具体步骤通常包括建立目标层-指标层-方案层的结构模型、构建模型成对比较矩阵、权向量计算和一致性检验。参考前人研究(石菊松等,2007;张春山等,2008;孙炜锋等,2008;唐亚明等,2011),对指标因子之间的相互重要性进行打分,结合研究区滑坡易发性因子级别分级结果,构建判断矩阵(表2);针对高程、坡度、岩组、NDVI、断裂、水系、公路分析得到特征向量为0.794,0.986,2.083,0.515,1.593,0.515,0.515;最大特征根为7.178,CI值为0.03,查明一致性指标RI为1.36。随后根据检验性指标CR=CI/RI=0.03/1.36=0.022<1,判断该矩阵具有较好的一致性,最终归一化处理,得到各影响因子的权重Wi(表2)。权重计算结果表明,影响研究区内滑坡稳定性的因子其重要性依次为:岩组、断裂、坡度、高程、NDVI、水系和公路。在信息量表(表1)中,高程和距公路距离排序靠前,就有很大的滑坡灾害发生的可能性。如果仅使用信息量模型,而无各因子的权重对此进行约束,则对将来潜在地质灾害发展趋势的判定就局限于过往的数据统计中,做出的易发性分区图极有可能过度强调高程和公路的作用。
表1 各评价因子状态信息量表
表2 因子权重配对比较矩阵表
3.3 滑坡易发性评价
运用ArcGIS重分类功能对各因子已分级的栅格图层赋予各自信息量值,再运用栅格计算器按照公式(2)叠加各单因子信息量与权重之积,得到研究区各评价单元的加权总信息量值,范围为-1.49~0.56(四舍五入保留小数点后两位,下同),最后对叠加后的栅格图层进行重分类,采用自然间断点分级法将其分为4个区间(图4),从大到小依次分为极高易发区(0.09~0.56)、高易发区(-0.14~0.09)、中易发区(-0.37~-0.14)和低易发区(-1.49~-0.37)。上述分区的面积依次为97.29 km2、124.71 km2、94.99 km2和121.73 km2,对应发生滑坡点数为254、150、37、19个。
图4 麦积区滑坡灾害易发性分区图
4 评价结果分析
在地质灾害易发性评价过程中,常使用ROC(Receiver Operating Characteristic Curve)曲线对分区的结果进行检验。通过计算AUC(Area Under Curve)值来评价区划结果的准确性(Chung C J F et al.,2003;Jadda M et al.,2011)。
ROC曲线指受试者工作特征曲线/接收器操作特性曲线(receiver operating characteristic curve),是反映敏感性和特异性连续变量的综合指标。在地质灾害易发性或风险性评价中,ROC曲线的横轴一般为信息量值由高到低的累计面积百分比,纵轴一般为对应信息量叠加值范围内地质灾害点个数累计百分比。一般认为,AUC为0.5~0.7时,准确性较差;AUC为0.7~0.9时,具有较好的准确性;AUC为0.9以上时,准确性极好;AUC越接近1,模型精度越好。笔者根据加权信息量法所评价的结果,绘制ROC曲线(图5)。计算加权信息量法曲线AUC值为75.04%,介于0.7~0.9,其评价精度具有较好的准确性。
图5 麦积区滑坡灾害易发性分区评价结果:ROC曲线图
5 讨论
(1)以天水市麦积区滑坡地质灾害为研究对象,在野外调查基础上,经权衡现实条件及资料获取程度,选择高程、坡度、工程地质岩组类型、植被归一化指数、距断裂距离、距水系距离及距公路距离等因子作为易发性评价的考量因子。基于ArcGIS平台,采用加权信息量法进行了易发性评价,按自然资源部中国地质调查局地质调查技术标准(DD 2019-08)—地质灾害调查技术要求(1∶50 000),将研究区分为极高易发区、高易发区、中易发区及低易发区。经ROC曲线验证AUC值为75.04%,评价结果具有一定的准确性,可为天水地区滑坡灾害防治规划提供科学依据。
(2)采用加权信息量法对研究区易发性评价的分类显示,极高易发区面积为97.29 km2,区内对应发生滑坡数量为254个;高易发区面积为124.71 km2,区内发生滑坡数量为150个;中易发区面积为94.99 km2,区内发生滑坡数量为37个;低易发区面积为121.73 km2,区内发生滑坡数量为19个。
极高和高易发区主要分布于渭河两岸,距离断裂带1.5 km范围内,地貌类型主要为黄土丘陵,岩性主要为黄土、冲积沙砾石;黄土一般覆盖厚度较大,且坡度为10°~30°,冲沟切割较深,灾害多分布在深切冲沟两岸。中易发区和低易发区主要分布于渭河、牛头河、葫芦河河谷地区及两岸支流的冲积阶地及西南基岩山区。
(3)对评价因子的信息量及权重分析表明,信息量值排名靠前的评价因子不一定在地灾发生过程里起到最重要的主导作用,也不一定象征着今后发展趋势是更容易发生在这些评价因子排名靠前的分类区间当中。如果单纯依据调查结果或资料搜集中所获得的编录数据来进行易发性评价,则对将来潜在的地质灾害发展趋势的判定就局限于过往的数据统计中,无法体现专家在宏观和微观上对地灾易发性的专业考虑。因此,在该类评价中要注重结合各评价因子间的相互重要性来进行综合评价。
(4)本文由于受限于研究水平及数据资料的限制,其深度、广度及精度仍存在不足之处。例如,降雨数据无法完整获取,地形起伏度数据计算量庞大。另外,由于地形复杂,为了计算方便,将评价单元按照网格划分,未按斜坡单元划分,没有很好的体现原始斜坡特征和地貌的吻合度。这些都需要再进一步的细致研究。