基于随机森林模型和地理探测器的赣南兴国县崩岗易发性评价

郭 飞 吴 迪 王秀娟 代巧怡 赖 鹏 陈 洋 夏 栋

(1.湖北长江三峡滑坡国家野外科学观测研究站, 湖北 宜昌 443002;2.三峡大学 土木与建筑学院, 湖北宜昌 443002;3.中国科学院 南京土壤研究所, 南京 210000;4.三峡大学 水利与环境学院, 湖北 宜昌443002)

崩岗侵蚀是发育于中国南方山地丘陵区一种受水力与重力复合作用而产生的山坡岩土分解、崩塌的地带性侵蚀现象[1-2],它是中国南方一种特殊的水土流失类型,国际上同类地貌分布少,仅马达加斯加的lavaka侵蚀沟谷、巴西的vocorocas地貌以及日本的“崩坏”地貌类似[3].尽管崩岗在地质灾害分布和侵蚀面积中所占比例不大,但其具有侵蚀量大、爆发力强、发展速度快等特点,易诱发泥石流和山体滑坡等地质灾害,还严重制约生态恢复、农业生产和区域社会经济可持续发展[4-5].因此,开展崩岗空间预测研究对崩岗预警防控、国土空间规划具有重要意义.

崩岗易发性,即基础环境背景下发生崩岗侵蚀的概率[6],是土壤侵蚀研究领域所关注的重点.随着人工智能技术的进步,机器学习方法迅速发展成为崩岗易发性评价模型的核心.

崩岗的发育是一个复杂的动态过程,受到岩性[7]、土壤类型与性质[8]、植被[9]、土地利用[10]、地形[11]等因素的显著影响,这使其发生具有强烈的复杂性和不确定性,故明确主导因素、筛选代表性的影响因素对崩岗易发性评价具有重要意义.国内外研究中常见的指标筛选方法可概括为线性选择方法和非线性选择方法.最普遍的线性选择方法包括主成分分析[12]、共线性分析[13]、逐步回归[14],方差显著性检验[15]等.Yu等[12]构建共线性分析排除高度相关的影响因素,对福建省崩岗进行易发性建模,取得了良好的精度.此外,灰色关联度、递归特征消除、随机森林等非线性选择方法也被应用[16-17].廖凯涛等[15]构建灰色关联度分析赣州市崩岗侵蚀关键影响因子,认为其侵蚀主要受到地形和水力因素的影响.

尽管上述方法均取得不错的预测精度,适用于崩岗易发性评价,但均没有考虑到地理空间内不同要素对崩岗发育过程的影响.地理探测器(Geodetector,GD)[18]是一种空间统计工具,该工具凭借合理分析各种现象的驱动力和影响因子的能力,已经被广泛应用到自然等领域,并取得符合现实且预测精度高的结果.如张若婧等[19]运用参数最优地理探测器分析江西省山洪灾害的主要驱动因素及作用机制,探究其空间分异特征,为山洪灾害防治提供科学依据.Yang等[20]将GD 与逻辑回归相结合,与传统的逻辑回归模型相比,精度提高了11.9%.利用GD 对崩岗进行空间分异性分析,有助于崩岗研究.

近些年来,崩岗易发性评价多选用数据驱动模型,最常用的有信息量模型[21]、逻辑回归模型[22]、多项式逻辑回归模型[23]、随机森林模型(Randomforest,RF)[24]、多层感知机模型[22]、人工神经网络模型[22]、支持向量机模型[25]等.亦有学者将多个模型进行对比分析,较一致地认为,随机森林模型可处理高维度、大数据量的数据集,对评价单元的限制小,无需考虑数据尺度和数据分布,超参数调试代价相比神经网络、深度学习等明显较低,被广泛应用于易发性评价模型中[26].

综上,本文以赣南兴国县花岗岩区为典型案例,利用地理探测器筛选指标体系,构建随机森林模型对该地区崩岗易发性展开评价,以期为崩岗预警和防控工作提供参考.地理探测器-随机森林(GeoDetector-RandomForest,GD-RF)模型主要包含4 步(如图1所示):①基于频率比分析,选取崩岗相关的环境因子;②创建非崩岗点,构建样本集;③使用地理探测器筛选主导因子,构建崩岗易发性评价指标体系;④将样本与评价指标因子输入RF模型中,进行崩岗易发性指数计算.

图1 GD-RF模型流程

1 研究区概况与环境因子数据

1.1 研究区概况

研究区位于江西省中南部兴国县境内(115°01'～115°51'E,26°03'～26°41'N),总面积3 215 km2.地处南岭东西向复杂构造带东段北侧、赣江支流的平江流域,形成低山、丘陵、平原等破碎折叠地貌,总地势自东北西边缘逐渐向中南部倾向,海拔114～1 186 m(如图2所示).属亚热带季风气候,气候温暖湿润、雨热同季,多年平均气温为18.8℃,降雨量达1 515.6 mm,且主要集中于5～9月.地带性土壤为第四纪红黏土,母岩以花岗岩、红砂岩、页岩为主,土壤含沙量大、持水性差、易侵蚀.由于近现代以来强烈人为活动(砍伐、开垦、工程开发等)造成严重植被破坏,并诱发了持续性土壤侵蚀,因此被为“江南红色沙漠”.

图2 研究区位置图及崩岗分布图

据2015年江西省水土保持规划崩岗调查数据显示,兴国县现存崩岗点2 933个,分布密度为0.91个/km2,其中约80%(2 460)分布于花岗岩区(图2).

1.2 崩岗环境因子数据

崩岗容易受到地形、气象、土壤和植被等方面的影响.因此,本文依据前人研究中提出的各类崩岗诱发因素,并结合现有可利用的数据源,选取10个环境因子(如图3所示)作为崩岗潜在影响因素.数据的来源及采集时间见表1,各因子是基于SAGAGIS软件提取的[26],将其重采样为30 m×30 m 作为栅格单元大小.并用频率比法[27]分析崩岗与环境因子间的非线性关系,其统计结果见表2.

表1 实验数据来源

表2 各影响因子的频率比值

图3 研究区各环境因子

1)地形因素

据观察(表2),崩岗分布集中在朝南方向和坡度小于15°地区,即主要发生在中坡度与缓坡上的阳坡和半阳坡.此外,崩岗发生率随地形湿度指数值增大而增大,表明含水率越高,更利于崩岗的发生.就坡面形态而言,崩岗发生率随剖面曲率、平面曲率值呈现先上升后下降的趋势,在数值接近零是达到峰值,79.70%和78.17%的崩岗集中分布在平面曲率和剖面曲率在-0.5～0.5之间的地区,进一步说明在平缓区更易孕育崩岗.

2)土壤因素

崩岗在黏土含量5%～30%之间,频率比值均大于1,说明在这个范围内崩岗较易发生侵蚀;砂含量在大于50%时,其频率比值最大,达到2.01,对崩岗发生起到正向作用.

3)气象因素

降雨量低于1 150 mm 的条件下,崩岗的发生占总量的78.65%;降雨在1 130 mm 以下时,频率比值2.77为最大值.年均降雨侵蚀力在不断增大的条件下,崩岗密度及发生率都在不停的波动,在年均降雨侵蚀力小于9 700 MJ·mm/(hm2·h·a)时,崩岗的发生率为47.04%,在9 700～10 100 MJ·mm/(hm2·h·a)条件下,频率比值最大.可见崩岗的出现主要是受低降雨量影响.

4)植被因素

频率比值随着植被覆盖度值的增加呈现先增加后减小的趋势,可见,并不是植被越茂盛崩岗就一定不会发生.

结合以上环境制约因素,赣南兴国县崩岗地貌主要集中在坡度平缓至水平,半阳坡至阳坡,植被覆盖稀疏至中等的花岗岩区.总体而言,所选取的环境因素与崩岗存在非线性关系,可作为崩岗的影响因素.

2 构建指标体系

2.1 地理探测器

地理探测器[18]无线性假设,能够探测崩岗的空间分异性,分析不同分层内影响因子对崩岗发生的解释力度.作为空间数据探索和分析十分可靠的工具之一,主要被用来分析各种现象的驱动力和影响因素,定量化筛选出贡献较大的因素.地理探测器由因子探测、交互探测、风险探测、生态探测构成.在定量化探测某因子对属性的空间分异解释程度时,属因子探测,用q值[28]来度量,其表达式为:

式中:q为崩岗空间分异影响力探测力指标;Nh和N分别为分层h和全区的样本单元数;L为变量Y(崩岗)或因子X(影响因子)的分层,即分类或分区;σ2h和σ2分别是分层h和全区的Y值的方差.q的值域为[0,1],值越大说明崩岗的空间分异性越明显;如果分层是由自变量X生成的,则q值越大表示自变量X对属性Y的解释力越强,反之则越弱.q值表示X解释了100×q%的Y.

2.2 评价指标选取

地理探测器无法处理连续数据,故将连续数据重分类.在花岗岩区随机选取与已发生崩岗数目相同的非崩岗作为负样本,构建30 m×30 m 的渔网,空间连接各因子属性值及样本集.样本集为Y(是否为已发生崩岗),各环境因子为X,代入地理探测器,探究10个环境因子对崩岗发生的贡献(如图4所示).

图4 因子探测结果

结果显示,地理分区解释力最强的是年均降雨侵蚀力,能够对崩岗解释29%,年均降雨量这类水文气象指标对崩岗解释力次之,能够解释25%的崩岗.其次是植被因素和土壤指标,说明这两类指标对于崩岗的发生具有一定的影响意义.此外,除坡度之外的其他地形因子对崩岗的解释力与其他因子相比相差少了1个数量级,对崩岗发生的解释力弱.

最后,选取年均降雨侵蚀力、年均降雨量、植被覆盖度、坡度、砂含量和黏土含量作为崩岗的主导因素,构建GD-RF模型进行易发性评价.

3 评价与检验模型

3.1 随机森林算法

随机森林是一种监督学习算法,所构建的“森林”是众多决策树的集成,以“bagging”方法来对其进行训练,用交叉验证方法选取模型样本以获得更准确和稳定的预测[29],被认为是一种集成学习方法,根据自举技术产生的一些弱分类树的聚合结果对未知样本进行分类.RF主要思想是有放回地抽取样本以及随机选取不同的输入特征数量来构建不同的训练集,使产生的决策树多样化[30].通过计算评价因子k在节点分割时基尼指数的减少值DGk;将森林中所有节点的DGk求和后对所有树取平均,即为评价因子k的重要性.以评价因子平均基尼减小值占所有因子平均基尼减少值总和的百分比度量评价因子的重要程度.按式(2)计算:

式中:m、n、t分别是评价因子总数、分类树棵数和单棵树节点数;DGkhj为第k个评价因子在第h棵树的第j个节点的基尼指数减少值;Pk为第k个评价因子在所有评价因子中的重要程度.

3.2 GD-RF模型精度检验

受试者工作特征曲线(receiver operating characteristic curve,ROC)可用于评价地质灾害预测模型的准确性[31-32].本研究将容易发生崩岗的预测数据与已有的崩岗灾害点进行比较.ROC 曲线显示了真阳性率(RTP)与假阳性率(RFP)的函数,突出观测数据和预测数据之间的一致性.具体来说,RTP描述了一组崩岗像素中被正确分类为崩岗的比例,而RFP是被错误分类为崩岗的非崩岗像素的比例,即

式中:真正(mTP)和真负(mTN)为正确分类的像素个数,假正(mFP)和假负(mFN)为错误分类的像素个数.

4 结果与分析

4.1 崩岗易发性模型性能评价

ROC与预测准确度之间的定量-定性关系范围为[0,1][33].基于Rstuido软件,得到兴国县易发性预测评价结果的ROC 曲线(如图5 所示),其对应的AUC值为0.864 4,表示评估精度良好.这说明基于地理探测器的随机森林模型所构建的崩岗易发性评价方法符合崩岗发生规律,具有客观性、可行性和现实可操作性,为兴国县花岗岩区的地质灾害防治起到一定的指导作用.

图5 崩岗易发性结果检验曲线

4.2 崩岗易发性空间分布及统计特征

将崩岗发生点位与易发性分区图进行叠加分析,得到崩岗易发性等级与崩岗分布预测结果,见表3.

表3 崩岗易发性等级与崩岗分布预测结果图

结果表明,极高-高易发区面积占兴国县花岗岩区域的41.49%,其中崩岗已发生点落入极高-高易发区共有1 999 处,占花岗岩区总崩岗发生数的81.26%,由此可知,极高-高易发区崩岗相对点密度为1.96.同时,极低易发区、低易发区分别占兴国县花岗岩区域的32.13%、15.91%,相应的崩岗数量密度分别为3.90%、7.28%.此外,每个易发分区类的面积百分比与每个类别中发生的现象百分比的比率被称为SCAI(Seed Cell Area Index)[34],它背后的逻辑是,在非常保守的区域范围内对是否容易产生崩岗进行正确的分类.Süzen等认为[35]在极高、高易发区应有较小的SCAI值,在极低、低易发区应具有较高的SCAI值,这是较为理想的SCAI值分布.本研究中,当易发性类别由低到高时,SCAI值也由较高的8.23递减到较低的0.45,可见该研究在随机森林模型下得到的结果(如图6所示)是合理的.

图6 随机森林下的崩岗易发性图

5 结 论

1)兴国县崩岗易发性评价指标贡献排序为年均降雨侵蚀力＞年均降雨量＞植被覆盖度＞坡度＞砂含量＞黏土含量＞平面曲率＞坡向＞地形湿度指数＞剖面曲率.

2)采用随机森林模型AUC 值为0.864 4,相应SCAI值由极高易发区到极低易发区不停增大分别为0.45、0.83、1.39、2.19和8.23,与实际情况较吻合,表明GD-RF模型适用于该区崩岗易发性评价.极高-高易发区主要集中在永丰镇、古龙岗镇、江背镇和鼎龙乡等地区.