张 曦,陈丽霞,徐 勇,连志鹏
(1.中国地质大学(武汉)地球物理与空间信息学院,湖北 武汉 430074; 2.中国地质调查局武汉地质调查中心,湖北 武汉 430205)
斜坡单元划分是开展大比例尺区域滑坡灾害易发性、危险性、风险评价的首要工作和重要基础,其划分过程繁琐、耗时且复杂,难有统一、客观的流程与方法。基于斜坡单元的灾害评价工作已广泛开展,但少有研究讨论不同划分方式对斜坡单元形态、面积大小、数量以及滑坡灾害易发性评价结果的影响。
区域滑坡灾害评价单元一般分为5种类型[1]:栅格单元、地貌单元、特定条件单元、斜坡单元和地形学单元。其中,斜坡单元能较好地体现1∶10 000尺度区域斜坡的整体性,相较于栅格单元能更有效地应用于评价区域防灾减灾工程与风险管理实践,发挥灾害区划理论研究成果的科学指导作用。国内外学者一直较为重视基于斜坡单元的滑坡灾害评价研究,如武雪玲等[2]基于斜坡单元对三峡库区滑坡易发性的评价研究,指出采用斜坡单元作为滑坡灾害易发性评价的基本单元,能够综合体现滑坡孕灾环境和诱发因子的作用,提高滑坡空间预测的精度和效率;邱丹丹等[3]基于斜坡单元对芦山地区滑坡危险性进行了评价,斜坡单元相比栅格单元能更好地保留斜坡的整体性,且评价结果更为精确;薛强等[4]基于斜坡单元对延安宝塔区滑坡易发性进行了评价,认为以斜坡作为单元进行滑坡易发性评价可以提高与实际地形地貌的吻合度,能够更好地体现区域中滑坡的实际发育状况;Erener等[5]指出斜坡单元能有效地消除栅格单元中无法考虑单元间联系和地面形态元素的缺点。
目前,关于斜坡单元的划分方法研究已有一些研究成果,如Philip等[6]提出了遥感影像光谱与地形剖面曲率相结合的划分方法,但受限于遥感影像和数字高程模型(DEM)的精度,该方法对1∶10 000图幅的划分不易发挥较好的效果;Xie等[7]利用水文分析法提取分水岭与汇水网划分斜坡单元的方法已被广泛接受并应用于实践中,但划分过程中会出现较多的细小破碎面和不合理长条状面,后期需要进行大量手动修改工作;颜阁[8]利用分水岭分析地表曲率划分斜坡单元,其步骤简洁,单元面积与形状均匀,但缺乏对地形、地质界线和人为活动等要素的考虑。
在滑坡灾害易发性评价中,斜坡单元的划分也需要考虑研究区历史灾害的发育规模,以确定控制单元大小的阈值,使评价结果更符合实际。但是,由于地形、地质、水文环境等条件的差异性,各研究区内斜坡的分布也呈现不同的空间分布特征,如何在这些差异中寻找共性,并结合滑坡灾害实际发育规模和评价尺度,通过滑坡易发性评价效果选用合适的划分方式,是本文研究的重点。为此,本文将以我国鄂西山区咸丰县高乐山镇为研究区,对比分析两种不同斜坡单元划分方法(水文分析法和地表曲率分水岭法)的优劣和滑坡灾害易发性评价的效果。
斜坡单元的划分主要采用ARCGIS平台,使用数字高程模型(DEM)、地质图和遥感影像来实现,其划分过程依次为:去除平地单元、提取山谷山脊线、添加相关界线。其中,山谷山脊线的提取是斜坡单元划分中最复杂、最重要的步骤,本文对比分析了水文分析法和地表曲率分水岭法两种方法提取山谷山脊线的具体流程。
基于水文分析法提取山谷山脊线的流程(见图1)如下:
图1 水文分析法提取山谷山脊线的流程图Fig.1 Procedure of HM for extraction of valley and ridge line
(1) 获取正DEM和反DEM。反DEM即用研究区的最高高程减去原DEM得到。
(2) 分别对正、反两个DEM进行填挖,防止后期水文分析时出现方向混乱与积水区域。
(3) 对填挖后的DEM提取流向。通过提取研究区内的地表水流动方向能够快速地区分出地表形态极速变化的区域,也是后续水文分析的必备步骤。
(4) 建立河网的连接。在流向的基础上提取流量、汇水点,连接各级水系得到河网。其中,通过正DEM得到的河网大致相当于地形中的山谷;通过反DEM得到的河网大致相当于地形中的山脊。
(5) 将正、反两个DEM得到的河网汇总,即可得到山谷线和山脊线。
基于地表曲率分水岭法(以下简称曲率分水岭法)提取山谷山脊线的流程(见图2)如下:
图2 曲率分水岭法提取山谷山脊线的流程图Fig.2 Procedure of CWM for extraction of valley and ridge line
(1) 获取整个研究区的地表曲率。在得到的地表曲率图层中,数值为正表示此处的地形向上凸出,数值越大表示凸出程度越明显;数值为负则表示此处的地形向下凹陷,数值越小表示凹陷程度越明显。
(2) 以地表曲率为基础提取整个研究区的流向。虚拟水流的方向即为高值至低值,山脊至山谷。
(3) 计算洼地。由于水流方向基本确定,积水处(洼地)会落在山谷(地表曲率图层中的低值)附近。
(4) 基于流向与洼地计算分水岭。得到的分水岭回溯到DEM中即为山脊线。
(5) 反转地表曲率图层每个栅格的正负值,重复上述步骤可以得到山谷线,最终汇总得到山谷线和山脊线。
斜坡单元面积的大小和形状是影响斜坡单元划分结果的两个主要要素,本文引入形状指数F[9]来进行衡量,即
F=L2/(4πS)
(1)
式中:L为斜坡单元周长(m);S为斜坡单元面积(m2)。
当F值趋近于1时,表示斜坡单元形状趋近于圆形;F值越大,则表示斜坡单元的形状越显狭长。
本文选取湖北省鄂西南咸丰县高乐山镇为研究区,分别采用水文分析法和曲率分水岭法进行了斜坡单元提取和滑坡灾害易发性评价,并对比分析了两种斜坡单元划分方法对滑坡灾害易发性评价的效果。
研究区面积约23 km2,总体地势为南北高、中部较低,高程分布于700~1 100 m之间,坡度分布主要集中于0°~40°之间,高陡坡较为少见,坡向分布均匀,坡面形态多为凸型坡。
在地质条件方面,咸丰县高乐山镇地层岩性中除缺少下泥盆统、上下石炭统以及上三叠统地层外,从上古生界寒武系到新生代第四系地层均有出露。该地区大地构造位置处于扬子准地台上扬子台坪八面山台褶带恩施台褶束,构造体系归属新华夏系一级构造鄂西隆起带西部,组成褶束的褶皱以北东和北北东向为主,断裂构造则是沿背斜陡翼形成的北东向压性或压扭性断裂。
滑坡是咸丰县高乐山镇最主要的地质灾害类型,研究区内共有10处已发生的滑坡灾害,主要分布在研究区的东北至西南的中部条带区域内,形态多为舌形和半圆形,研究区滑坡分布及其形状指数和面积统计,详见图3和表1。本次研究涉及到的资料数据有咸丰县高乐山镇遥感影像、数字高程模型(3 m×3 m DEM)、断层、河流、地层岩性和历史滑坡,其中断层、河流、地层岩性、历史滑坡为线状矢量文件。数据来源于历史资料收集、无人机航测、遥感影像解译、区域地质环境详查、工程地质测绘等。
图3 咸丰县高乐山镇滑坡分布图Fig.3 Landslide inventory map of Gaolle Town,Xianfeng County,Hubei Province,China
滑坡编号形状指数F面积/(m2)滑坡编号形状指数F面积/(m2)XF11.0997249XF61.1118254XF21.17615411XF71.17110100XF31.07913156XF81.11010923XF41.07822490XF91.1383791XF51.3087754XF101.1857764
2.2.1 平地单元去除
依据DEM提取坡度、坡向、山体阴影(见图4),叠加分析坡度图层中0°~5°区域、坡向图层中平地(深绿色部分)、阴影图层中光影分割线,并综合遥感影像划分得到平地单元。若DEM精度较高导致坡向图层中无深绿色部分,则可降低DEM精度再次获取。
图4 平地单元提取辅助数据Fig.4 Auxiliary maps for removing flat areas
本次划分得到的平地面积共3.9 km2,占据研究区总面积约17.5%,主要分布于河流两岸,局部分布于山间、山顶平台。
2.2.2 水文分析法和曲率分水岭法提取山谷山脊线
基于水文分析法提取的斜坡单元数量为1 462个,斜坡单元平均面积约15 245 m2,当地势变化明显时,提取效果好,但当地面较平坦时,出现细小且狭长的破碎面和冗余的分界[见图5(a)],这是因为在水文分析法分析中地势平坦的区域水流无明确的方向。
基于曲率分水岭法提取出的地形单元分布均匀、大小一致,但提取的斜坡单元数量较多(8 500个),斜坡单元平均面积较小(约3 000 m2),这是因为在曲率分水岭法分析中,由于使用过高精度的DEM(3 m×3 m),使得地面细小的起伏均被识别为洼地,过多的洼地导致了斜坡单元数量的偏多和斜坡单元面积的细小[见图5(b)]。
图5 基于水文分析法和曲率分水岭法所得初步斜坡单元分布图(局部)Fig.5 Preliminary results of slope units by HM and CWM(part)
为了解决曲率分水岭法提取斜坡单元数量过多和斜坡单元面积小的问题,本文降低DEM和曲率图层精度,平滑时均采用了圆形均值滤波器,斜坡单位数量随平滑半径的变化见图6。
图6 斜坡单元数量随平滑半径的变化Fig.6 Variation of slope units number with radius of DEM and curvature resolution
由图6可见,当DEM平滑半径到10个栅格(30 m,0.9×103m2)后,斜坡单元数量快速上升;当曲率平滑半径到15个栅格(45 m,2.025×103m2)后,斜坡单元数量下降缓慢。当进一步增加滤波器半径时,不仅斜坡单元数量上升,且山脊线和山谷线出现较大程度偏移。综合以上两个特征,在对DEM平滑10个栅格提取曲率图的基础上,采用15个栅格平滑该曲率图并用于提取分水岭,斜坡单元数量降低为4 320个,斜坡单元面积平均增大1倍。经多次平滑降低精度,划分得到的斜坡单元分布及其形状指数和面积统计,详见图7和表2。
图7 多次平滑后斜坡单元分布图(局部)Fig.7 Distribution map of slope units by CWM after multiple smoothing (part)
形状指数F水文分析法占比曲率分水岭法占比1~20.45760.85812~30.24970.08913~40.12310.02714~50.05880.01255~60.03210.00356~70.02940.00197~80.01500.0016>80.03420.0062面积/m2水文分析法占比曲率分水岭法占比0~1030.43640.3269103~1040.18600.5042104~1050.36590.1690>1050.01160.0000
由表1可见,咸丰县高乐山镇发生的历史滑坡灾害面积主要集中在104m2范围内,形状指数F主要集中在1~1.5之间。综合表2可见,曲率分水岭法斜坡单元的初步划分结果在相应面积和形状指数区域分别达到85.81%和83.11%,均占比更大,划分效果更加优异。
最终,通过添加断层、水系、交通线路等明显分界线,即可完成一次斜坡单元的划分。
目前滑坡易发性评价体系已十分成熟,众多学者已经将信息量[10-11]、证据权[12-13]、支持向量机(SVM)[14-15]、神经网络[16]等模型或方法应用到实际问题中。本文选用证据权模型,结合咸丰县高乐山镇历史滑坡数据,将上述两种斜坡单元划分方法提取的斜坡单元用于滑坡灾害易发性评价。根据咸丰县高乐山镇实际情况与实测数据,最终选取坡度、坡向、高程、坡面形态、地层岩性、坡体结构、土地利用类型、道路、断层共9个因子进入评价体系,各评价因子证据权值统计结果见表3。
表3 研究区滑坡灾害易发性评价因子及其证据权值
两种斜坡单元划分方法得到的研究区滑坡灾害易发性分区结果见图8,两种斜坡单元划分方法的成功率(AUC值,即曲线线下面积)分别为71.69%和78.43%(见图9),其对应的不同等级滑坡灾害易发性分区面积占比统计见表4。
图8 两种斜坡单元划分方法滑坡灾害易发性分区图Fig.8 Landslide susceptibility map by HM(A) and CWM(B)
图9 两种斜坡单元划分方法滑坡易发性ROC曲线Fig.9 Receiver Operating Characteristic for landslide susceptibility by the two methods
由图9和表4可见,曲率分水岭法比水文分析法所得研究区滑坡灾害易发性分区精度略高6.7%,理论结果更可信;但水文分析法所得的滑坡灾害高、中易发区的面积达到44.9%,比曲率分水岭法略高7%,其评价结果更为保守。
表4 两种斜坡单元划分方法对应的不同等级滑坡灾害 易发性分区面积占比统计
本文针对大比例尺滑坡灾害区划工作中,斜坡单元划分方法不统一、使用效果不明确的问题,以我国鄂西地区咸丰县高乐山镇为研究区,对比分析了水文分析法和曲率分水岭法用于斜坡单元划分和滑坡灾害易发性分区的效果,得出以下结论:
(1) 基于曲率分水岭法划分得到的斜坡单元形状更接近历史滑坡形状(85%)和面积大小(83%),划分效果更优;但当DEM精度一定时,水文分析法在地势变化明显的区域斜坡单元提取效果较好。
(2) 曲率分水岭法划分斜坡单元数量分别随DEM和曲率精度下降而呈现先递减后增(DEM)或稳定(曲率)的规律,基于该规律,可确定咸丰县高乐山镇(1∶10 000尺度)斜坡单元划分的DEM平滑大小(30 m×30 m)和曲率平滑大小(45 m×45 m)。
(3) 曲率分水岭法比水文分析法所得滑坡灾害易发性分区精度略高6.7%;但水文分析法分区结果更为保守,滑坡灾害高、中易发区面积比曲率分水岭法略高7%。
[1] Guzzetti F,Carrara A,Cardinali M,et al.Landslide hazard evaluation:A review of current techniques and their application in a multi-scale study,Central Italy[J].Geomorphology,1999,31(1/2/3/4):181-216.
[2] 武雪玲,任福,牛瑞卿,等.斜坡单元支持下的滑坡易发性评价支持向量机模型[J].武汉大学学报(信息科学版),2013(12):1499-1503.
[3] 邱丹丹,牛瑞卿,赵艳南,等.斜坡单元支持下地震滑坡危险性区划——以芦山地震为例[J].吉林大学学报(地球科学版),2015(5):1470-1478.
[4] 薛强,张茂省,李林.基于斜坡单元与信息量法结合的宝塔区黄土滑坡易发性评价[J].地质通报,2015,34(11):2108-2115.
[5] Erener A,Düzgün H S B.Landslide susceptibility assessment:What are the effects of mapping unit and mapping method?[J].EnvironmentalEarthSciences,2012,66(3):859-877.
[6] Philip T G,Franklin S E.An automated approach to the classification of the slope units using digital data[J].Geomorphology,1998,21(3):251-264.
[7] Xie M,Esaki T,Zhou G.GIS-based probabilistic mapping of landslide hazard using a three-dimensional deterministic model[J].NaturalHazards,2004,33(2):265-282.
[8] 颜阁.华池县滑坡易发性制图[D].兰州:兰州大学,2016.
[9] 胡小芳,曾文雄,吴成宝,等.颗粒表面分维与其形状指数关系研究[J].中国粉体技术,2007(2):14-17.
[10]周超,殷坤龙,向章波,等.基于GIS的淳安县滑坡易发性定量评价[J].安全与环境工程,2015,22(1):45-50,55.
[11]朱吉祥,张礼中,梁国玲,等.基于GIS的信息熵理论在滑坡危险性评价中的应用——以四川省青川县滑坡为例[J].安全与环境工程,2012,19(3):15-19.
[12]许冲,戴福初,徐锡伟,等.基于GIS平台与证据权的地震滑坡易发性评价[J].地球科学(中国地质大学学报),2011,36(6):1155-1164.
[13]Lee S,Chol J.Landslide susceptibility mapping using GIS and the weight-of-evidence model[J].InternationalJournalofGeographicalInformation,2004,18(8):148-789.
[14]Marjanovic'm,Kovacevic' M,Bajat B,et al.Landslide susceptibility assessment using SVM machine learning algorithm[J].EngineeringGeology,2011,123(3):225-234.
[15]姜琪文,许强,何政伟.基于SVM多类分类的滑坡区域危险性评价方法研究[J].地质灾害与环境保护,2005,16(3):328-330.
[16]许冲,徐锡伟.基于GIS与ANN模型的地震滑坡易发性区划[J].地质科技情报,2012,31(3):116-121.