1. 引言
近年来,对泥石流进行危险性评价已成为泥石流灾害的预防和应对研究的重点[1]。泥石流易发性评价旨在识别评价区域内泥石流可能发生的地点,但不涉及泥石流发生的时间或规模[2]。这些方法通常依赖于静态环境因子,如坡度、坡向、植被覆盖度以及与河流的距离等,而未能充分考虑时间因素对泥石流发生条件的影响。在现实中,地理环境并不是静态的,泥石流的动态特征往往随地理环境的时序变化而演变。但这一关键变量往往被现有泥石流易发性评估模型所忽略,导致对活动泥石流的识别不够精确,进而影响了评估结果的可靠性。
随着遥感技术进步,干涉合成孔径雷达(InSAR)技术已成为监测区域地表形变、早期识别泥石流现象和可靠的泥石流预测的有效工具[3] [4]。InSAR技术能够捕捉到地表微小的形变信息,小基线集雷达干涉技术(SBAS-InSAR)能克服传统InSAR方法在时间和空间分辨率上的局限性[5]。SBAS-InSAR技术通过减少大气和轨道误差的影响,能够提供更加连续和精确的地表形变数据。这些数据不仅能够揭示泥石流的动态形变特征,而且对于提高泥石流易发性评估的精度具有重要意义。
因此,本文引入InSAR形变因子作为动态变化特征,进行泥石流易发性评估。目前的泥石流易发性评价,除了有传统的多因子综合评价法、信息量法等方法,还大量运用了机器学习算法模型,如随机森林、SVM和BP神经网络等模型。由于不同模型运用在不同目标区域和采取不同因子会得到不同的精确度,因此本文构建信息量模型、随机森林、SVM和BP神经网络进行泥石流易发性评估,采用不同精度评价指标对其结果进行对比分析,以获得研究区最佳评价算法及结果。
2. 研究区概况和数据来源
2.1. 研究区概况
本文的研究区域为云南省昆明市东川区的小江流域,见图1,处于103˚0'~103˚30'E,26˚0'~26˚30'N,是一个典型的暴雨型滑坡泥石流区域[6]。小江是东川的自然分界线,其东侧是乌蒙山系,最高峰牯牛寨海拔高4017.3 m;小江西侧为拱王山系,最高峰雪岭是昆明市海拔最高处,海拔4344.1 m;该区域的最低点海拔为695米,整体高差达3649.1米。该区域地形起伏显著,地处小江断裂带,地质条件复杂且地震多发,亚热带季风气候雨季暴雨频繁,土壤质地松软,水土流失问题尤为严重。这些因素共同作用,使得该区域成为地质灾害的高发地带。
Figure 1. Location of study area
图1. 研究区位置
2.2. 数据来源
泥石流是一种由多种因素共同作用而引发的地质灾害,其形成机制复杂且多变。黄发明等认为影响地质灾害的主要指标可分为四大类:地形、地质、水文环境以及人类活动[7]。本研究进一步细化并深化了这些基本指标,选择了地形、地质、水文、人类活动和InSAR形变五个影响因子大类,其中地形、地质、水文和人类活动为静态因子,InSAR形变为动态因子,具体为10个影响因子。这些因子共同构成了对泥石流易发性动态评估的数据基础,涵盖了从自然条件到人为因素的全方位影响因素。重采样为30 m,处理后得到泥石流影响因子信息,见表1。
Table 1. Mudslide conditioning factors
表1. 泥石流影响因子
  
Figure 2. Spatial distribution of static impact factors
图2. 静态影响因子空间分布
在地形因子方面,本文选取数字高程模型(DEM)、坡度、坡向、地形起伏度和地形粗糙度5个指标,分别见图2(a)、图2(b)、图2(c)、图2(d)、图2(e)。地质因子着重考虑了距断层的距离。由于多条区域性断层穿过本区,断层活动易造成岩体破碎,可为泥石流提供充足物源,且与断层距离越近,情况越严重[8],因此与地质断层距度是泥石流发育的重要影响因素,见图2(f)。水文环境因素是激发泥石流发生的决定性因素[9],包括距河流的距离和降雨量,分别见图2(g)、图2(h)。为了定量评估人类活动对地表影响的程度,引入归一化植被指数(NDVI),见图2(i)。
(5) InSAR形变因子
InSAR技术自问世以来,其在地表形变监测领域的应用不断深化,技术日益成熟[10]。作为一种先进的遥感技术,InSAR能够捕捉到地表微小的形变信息,对于地质灾害预警和环境变化研究具有重要意义。本文中,我们采用了SBAS-InSAR技术,这是一种特别适合于区域性形变监测的方法,因其能够有效降低大气延迟的影响,提高测量精度[11]。我们获取小江流域在2023年至2024年4月期间的形变数据,数据处理流程包括数据预处理、基线估计、相位解缠、形变时间序列分析等关键步骤,见图3。
Figure 3. Flow chart of data processing by SBAS-InSAR technology
图3. SBAS-InSAR技术处理数据流程图
(a)                                                       (b)
Figure 4. Average deformation rate graphs
图4. 平均形变速率图
为了深入分析小江流域的形变特征,根据图4(a),我们利用离散形变速率数据,运用克里金插值法对这些数据进行连续化和平滑处理。克里金插值[12]能够根据已知样本点的空间分布和变异特性对未知区域进行估计,从而生成连续的形变速率场。通过这种方法生成小江流域的平均形变速率等值线,以完善研究区域内的形变速率数据,见图4(b)。
3. 易发性评价模型
3.1. 信息量模型
信息量模型(IM)是一种基于历史泥石流灾害数据统计分析的综合性评估手段。该方法通过将一系列影响泥石流影响因子的实际观测值转化为量化的信息价值量,以精确度量这些因子对泥石流孕育与发生过程的贡献程度。该模型旨在揭示不同泥石流影响因子对泥石流灾害潜在易发性的影响力大小,即所谓的“信息量”。以下是该评估方法的数学表达式[13]。
                               (1)
式中:I为评价栅格单元中的总信息量,n为影响因子的数量,
为评价因子xi对泥石流灾害发生(H)所贡献的信息量,Ni为分布在评价因子xi内的泥石流灾害面积,N为整个研究区内泥石流灾害的总面积,Si为研究区内含有评价因子xi的面积,S为整个研究区的总面积。信息价值评估法的核心在于,通过计算得到的总信息量I,可以直观地评估泥石流灾害的潜在易发性。具体而言,当总信息量I较高时,表明该评价单元内受多种有利因子共同影响,泥石流灾害的发生概率显著增加,即灾害更易发;反之,若总信息量I较低,则意味着该区域泥石流灾害的发生条件相对不利,灾害相对不易发生[14]。
3.2. 随机森林
随机森林(RF)作为一种高度灵活且效能卓越的多分类器集成策略,其核心在于构建并整合多个决策树模型[15]。这些决策树并非孤立存在,而是作为一个集体共同工作,以实现更为全面和稳健的分类性能。每一棵决策树在随机森林的构建过程中均独立地接受训练[16],这一过程确保了模型的多样性和鲁棒性,因为不同的树可能会基于数据集的不同子集或特征组合进行学习,从而捕捉到数据中更为丰富的模式和结构。在构建每棵决策树时,若原始数据集D中共有n个样本,随机森林算法会对训练集进行Bootstrap采样,从总训练集中抽取x (x < n)个样本,即有放回地随机抽取样本来生成不同的训练数据集。在单个决策树的训练过程中,随机森林通常采用Gini系数进行节点的不纯度评价,公式为:
                                    (2)
其中,代表第N颗决策树中的第j个决策节点,K代表数据的总类别数量,代表第j个决策节点的分类结果中,类别i的数量占总类别的比例,系数越小代表该决策节点分出的样本纯度越高,分类效果越好。对N颗决策树重复上述步骤,将训练好的决策树结合起来即生成随机森林,其输出结果为:
                                  (3)
其中x为输入的测试数据,且代表第n颗决策树的分类结果,代表所有决策树的平均投票结果。这种随机性能够增加不同决策树彼此之间的差异性,避免了RF训练过程中的模型过拟合问题。
3.3. 支持向量机
支持向量机(SVM)是一种强大的机器学习算法,最初由Cortes等学者于1995年引入并正式提出[17]。该算法广泛应用于各类分类任务中,其核心思想在于通过一种映射机制,将原始的低维输入空间数据转换至一个高维的潜在特征空间中,使得原本在低维空间中难以线性分割的数据集能够在这个新构建的高维空间内实现有效的线性分离。其带入拉格朗日函数后的线性拟合函数为:
                                (4)
                                   (5)
其中,
为确定超平面方向的权重向量,
为偏置,
为拉格朗日乘子。本文采用的核函数为径向基核函数(RBF),其表达式为:
                                 (6)
其中,
为可调参数,控制RBF的宽度,
为数据点之间的欧式距离。
3.4. BP神经网络
BP神经网络(BPNN)由保罗·保姆哈特和大卫·鲍曼于1986年提出,广泛应用于诸多领域是一种按误差反向传播算法训练的多层前馈网络,分别由信息的正向传播和误差的反向传播两个过程组成[18]。该算法的学习过程包括多层前馈网络和反向误差修正两个阶段。多层前馈数学模型为:
                             (7)
其中,
为第
层第
个节点的输出值,
为第
层第
个节点的激活值,
为第
层中第
个节点到第
个节点的权重值,
为第
层第
个节点的阈值,
为总层数。
在优化过程中,反向误差调整阶段采用了梯度下降策略,该策略精细地调整神经网络各层级间连接的权重参数,旨在引导系统总误差沿其最小化方向逐步收敛,从而实现模型性能的显著提升。其表达式为:
                           (8)
4. 小江流域泥石流易发性评价
4.1. 影响因子共线性及重要性分析
利用SPSS26软件中的方差膨胀因子(VIF)计算了上述10个泥石流易发性影响因子的多重共线性值。当VIF ≥ 10或容差小于0.1时因子存在共线性问题[19]。10个选定的滑坡影响因素VIF值均大于1且小于5,容差均大于0.4小于1。VIF值及容差表明,10个泥石流影响因素相互独立,不存在共线性问题。这表明本文研究的10个泥石流影响因素可用于模型的学习和评估,以确保数据的有效性和准确性。
为明晰每一个因子的影响程度,进一步研究了每一个滑坡因子。随机森林模型可以在分析数据时提供可变重要性分数[20],利用随机森林模型的基尼重要性可以计算不同影响因子对泥石流发生的重要性得分,这些结果从高到低为DEM > 降雨量 > InSAR形变 > 地形起伏度 > 距河流距度 > NDVI > 地质断层距度 > 坡向 > 坡度 > 地形粗糙度。结果表明,InSAR形变因子在模型评估中的作用是十分重要的,其重要性超过了所选取的7个静态影响因子。因此,本文共使用4种模型进行泥石流易发性评价,以对比哪个模型评价效果更佳;每种模型分为两组进行对照比较,一组加入InSAR形变因子,另一组去除InSAR形变因子,以对比InSAR形变因子对模型的影响程度。
4.2. 信息量模型影响因子分析
Table 2. Impact factor status informativeness
表2. 影响因子状态信息量
 
  
    | 影响因子 | 等级划分 | 信息量 | 排序 | 影响因子 | 等级划分 | 信息量 | 排序 | 
  
    | DEM/m | 25~43 | 1.1505 | 5 | 地质断层距度/m | 0~2602 | 0.2616 | 14 | 
  
    | 43~58 | 1.1005 | 6 | 2602~5811 | −0.275 | 44 | 
  
    | 58~72 | 0.0253 | 21 | 5811~9714 | −0.454 | 46 | 
  
    | 72~85 | −1.88 | 58 | 9714~14485 | 0.0786 | 19 | 
  
    | 85~97 | −0.633 | 49 | 14,485~22,204 | −0.002 | 32 | 
  
    | 97~108 | −0.63 | 48 | 距河流距度/m | 0~579 | 0.7982 | 8 | 
  
    | 108~119 | −0.022 | 34 | 579~1190 | −0.265 | 41 | 
  
    | 119~134 | −0.273 | 43 | 1190~1844 | −0.915 | 55 | 
  
    | 134~167 | 0 | 23 | 1844~2605 | −0.722 | 50 | 
  
    | 坡度/(˚) | 0~0.6 | −0.169 | 40 | 2605~5123 | 0 | 23 | 
  
    | 0.6~1.4 | −0.081 | 37 | 降雨量/mm | 163.4~174.4 | 1.2021 | 3 | 
  
    | 1.4~2.1 | 0.0678 | 20 | 174.4~182.1 | 1.1598 | 4 | 
  
    | 2.1~3.2 | 0.4423 | 10 | 182.1~188.3 | −0.051 | 35 | 
  
    | 3.2~12.3 | 0 | 23 | 188.3~193.3 | −0.897 | 54 | 
  
    | 坡向(˚) | −1~58 | −0.347 | 45 | 193.3~201.9 | −0.879 | 53 | 
  
    | 58~140 | −0.832 | 51 | DNVI/% | 0~15 | 0.2008 | 16 | 
  
    | 140~213 | 0.3745 | 11 | 15~34 | 0.2618 | 13 | 
  
    | 213~283 | 0.4901 | 9 | 34~55 | −0.267 | 42 | 
  
    | 283~358 | 0.148 | 18 | 55~80 | −0.538 | 47 | 
  
    | 地形起伏度/m | 0~4 | −0.119 | 38 | 80~100 | −0.162 | 39 | 
  
    | 4~8 | 0.264 | 12 | InSAR形变/m | −118~−32 | 0 | 23 | 
  
    | 8~12 | 0.0151 | 22 | −32~−19 | −0.005 | 33 | 
  
    | 12~16 | −0.838 | 52 | −19~−12 | −1.628 | 56 | 
  
    | 16~36 | 0 | 23 | −12~−7 | −1.809 | 57 | 
  
    | 地表粗糙度 | 1~1.000235 | −0.076 | 36 | −7~−2 | 0.1522 | 17 | 
  
    | 1.000345~1.000885 | 0.2229 | 15 | −2~3 | 0.8512 | 7 | 
  
    | 1.000885~1.001825 | 0 | 23 | 3~9 | 1.5392 | 2 | 
  
    | 1.001825~1.004212 | 0 | 23 | 9~21 | 1.5914 | 1 | 
  
    | 1.004212~1.023367 | 0 | 23 | 21~96 | 0 | 23 | 
 结合影响因子共线性及重要性分析结果选取DEM、坡度、坡向、地形起伏度、地表粗糙度、地质断层距度、距河流距度、降雨量、NDVI、InSAR形变10个指标作为信息量模型的评价指标体系,将各指标因子的信息量值进行统计并进行排序,见表2,排序前5位的分别为:InSAR形变[9 m, 21 m) (信息量值1.5914)、InSAR形变[3 m, 9 m) (信息1.5392)、降雨量[163.4 mm, 174.4 mm) (信息量值1.2021)、降雨量[174.4 mm, 182.1 mm) (信息量值1.1598)、系距DEM [25 m, 43 m) (信息量值1.1505)。这表明小江InSAR形变、降雨量、DEM是影响泥石流发育的主要影响因子[21]。
4.3. 易发性评价结果对比分析
采用自然断点法将信息量值以及易发性概率值划分为5个等级:极低易发区、低易发区、中易发区、高易发区和极高易发区[22],以此得到了4种模型的泥石流易发性评价结果,见图5(a)、图5(c)、图5(e)、图5(g),并与原始静态评价因子模型进行对比,见图5(b)、图5(d)、图5(f)、图5(h)。
 (a) InSAR_IM                         (b) IM
 (c) InSAR_RF                             (d) RF
 (e) InSAR_SVM                          (f) SVM
 (g) InSAR_BP                             (h) BP
Figure 5. Results of mudslide susceptibility in the study area
图5. 研究区泥石流易发性结果
结果表明,4种算法评价得出的研究区泥石流易发性结果在空间位置分布上存在一定的相似性和差异性,见图5。具体而言,在4种算法得到的易发性结果:极高易发区集中分布在小江沿岸的大沟、老凹沟和蒋家沟等地区,说明这一带相对其他区域发生泥石流的可能性较大,和历史泥石流灾害点一致;极低易发区集中分布在离小江水系较远的区域,如拱王山和牯牛寨山,这些区域离小江主干远且植被覆盖度相对较高,不利于泥石流的发生,因此被赋予较低的易发性。四种算法都与实际情况相吻合,但四者的差异,主要是在小江东岸和西岸中低易发区的划分,其中信息量模型与其他三种机器学习算法差异最大。从是否联合InSAR形变因子来看,两者划分大体相似,略有不同。
通过定性的比较,难以比较模型的优劣。因此,我们对模型评价结果准确性进行量化分析,利用Python分别计算模型的准确率(ACC),精确率(PRE),AUC值,见表3,并画出其ROC曲线,见图6。
Table 3. Evaluation model accuracy comparison
表3. 评价模型精度对比
 
  
    | 评价模型 | ACC | PRE | AUC | 评价模型 | ACC | PRE | AUC | 
  
    | InSAR-IM | 0.67 | 0.67 | 0.777 | IM | 0.61 | 0.61 | 0.761 | 
  
    | InSAR-RF | 0.75 | 0.75 | 0.846 | RF | 0.73 | 0.73 | 0.837 | 
  
    | InSAR-SVM | 0.77 | 0.77 | 0.849 | SVM | 0.72 | 0.72 | 0.840 | 
  
    | InSAR-BP | 0.79 | 0.79 | 0.855 | BP | 0.76 | 0.76 | 0.842 | 
 
Figure 6. Comparison of model ROC curves
图6. 模型ROC曲线对比图
从总体看,InSAR-BP神经网络的ACC值、PRE值和AUC值最高分别为0.79、0.79、0.855,InSAR-SVM模型次之为0.77、0.77和0.849,IM最低为0.61、0.61和0.761,见表3。从横向对比来看,联合InSAR形变因子的模型指标均优于原始模型。其中从AUC值来看,四种模型分别提高了1.6%,0.9%,0.9%和1.3%。从纵向对比来看,BP神经网络的所有指标均优于其他三种算法。其中,InSAR-BP神经网络的AUC值较InSAR-IM、InSAR-RF和InSAR-SVM分别提升了7.8%、0.9%和0.6%,BP神经网络的AUC值较IM、RF和SVM分别提升了8.1%、0.5%和0.2% (图6),表明BP神经网络在小江流域具有更高的泥石流预测能力。
此外,为了进一步验证模型的精度和检验模型易发性分区的科学合理性,实验统计了各模型易发性分区的检验结果,见表4。结果表明,所有模型随着易发性等级的增加,泥石流的数量和密度均逐渐递增,表明本次研究结果有效且合理[23]。采用历史泥石流灾害点落在极高易发区、高易发区的比例表征模型预测的精度[24]。结果表明,从总体看,InSAR-BP神经网络的精度最高为89.5%,IM最低为63%。从横向对比来看,联合InSAR形变因子的模型精度均优于原始模型。四种模型分别提高了5.5%、1%、1%和1.5%。从纵向对比来看,BP神经网络的所有指标均优于其他三种算法,其中,InSAR-BP神经网络的精度较InSAR-IM、InSAR-RF和InSAR-SVM分别提升了21%、8%和2%,BP神经网络的精度较IM、RF和SVM分别提升了25%、7.5%和1.5%,检验结果与上述模型评价指标吻合,验证了模型的精度,表明了BP神经网络结果的精度和合理性皆优于其他模型,能够为该区域的防灾减灾提供可靠参考。
Table 4. Model susceptibility partitioning test results
表4. 模型易发性分区检验结果
 
  
    | 模型 | 易发性分区 | 泥石流数量 | 分级面积 | 面积比 | 密度 | 模型 | 易发性分区 | 泥石流数量 | 分级面积 | 面积比 | 密度 | 
  
    | InSAR-IM | 低 | 3 | 584.2 | 14.7% | 0.5% | IM | 低 | 9 | 685.2 | 17.2% | 1.3% | 
  
    | 较低 | 21 | 1150.4 | 28.9% | 1.8% | 较低 | 28 | 1313.4 | 33.0% | 2.1% | 
  
    | 中 | 39 | 1164.4 | 29.2% | 3.3% | 中 | 41 | 1007.3 | 25.3% | 4.1% | 
  
    | 高 | 42 | 575.3 | 14.4% | 7.3% | 高 | 42 | 509.0 | 12.8% | 8.3% | 
  
    | 极高 | 95 | 509.2 | 12.8% | 18.7% | 极高 | 84 | 468.7 | 11.8% | 17.9% | 
  
    | InSAR-RF | 低 | 2 | 969.2 | 24.3% | 0.2% | RF | 低 | 1 | 862.9 | 21.7% | 0.1% | 
  
    | 较低 | 9 | 998.1 | 25.1% | 0.9% | 较低 | 11 | 1054.8 | 26.5% | 1.0% | 
  
    | 中 | 26 | 786.2 | 19.7% | 3.3% | 中 | 27 | 901.1 | 22.6% | 3.0% | 
  
    | 高 | 63 | 797.6 | 20.0% | 7.9% | 高 | 60 | 702.4 | 17.6% | 8.5% | 
  
    | 极高 | 100 | 431.9 | 10.8% | 23.2% | 极高 | 101 | 461.9 | 11.6% | 21.9% | 
  
    | InSAR-SVM | 低 | 1 | 1093.5 | 27.5% | 0.1% | SVM | 低 | 2 | 1100.4 | 27.6% | 0.2% | 
  
    | 较低 | 8 | 780.0 | 19.6% | 1.0% | 较低 | 8 | 778.8 | 19.6% | 1.0% | 
  
    | 中 | 16 | 669.7 | 16.8% | 2.4% | 中 | 17 | 662.3 | 16.6% | 2.6% | 
  
    | 高 | 35 | 650.1 | 16.3% | 5.4% | 高 | 35 | 650.6 | 16.3% | 5.4% | 
  
    | 极高 | 140 | 789.7 | 19.8% | 17.7% | 极高 | 138 | 791.0 | 19.9% | 17.4% | 
  
    | InSAR-BP | 低 | 1 | 1241.1 | 31.2% | 0.1% | BP | 低 | 1 | 1249.7 | 31.4% | 0.1% | 
  
    | 较低 | 8 | 675.7 | 17.0% | 1.2% | 较低 | 7 | 652.1 | 16.4% | 1.1% | 
  
    | 中 | 12 | 585.8 | 14.7% | 2.0% | 中 | 16 | 549.6 | 13.8% | 2.9% | 
  
    | 高 | 34 | 599.0 | 15.0% | 5.7% | 高 | 30 | 599.2 | 15.0% | 5.0% | 
  
    | 极高 | 145 | 881.4 | 22.1% | 16.5% | 极高 | 146 | 932.4 | 23.4% | 15.7% | 
 5. 结论
(1) 4种算法的易发性结果在空间分布上具有较高的一致性,极高易发区分布在小江沿岸的大沟、老凹沟和蒋家沟等地区,极低易发区集中分布在离小江水系较远的区域,如拱王山和牯牛寨山。
(2) 联合InSAR形变因子的模型AUC较原始模型分别提高了1.6%、0.9%、0.9%和1.3%,表明联合InSAR形变因子可以提升泥石流易发性评估的准确度。InSAR-BP神经网络较InSAR-IM、InSAR-RF和InSAR-SVM分别提升7.8%、0.9%和0.6%,BP神经网络的AUC值较IM、RF和SVM分别提升了8.1%、0.5%和0.2%,表明BP神经网络在小江流域具有更高的泥石流预测能力。
(3) 联合InSAR形变因子的模型精度均优于原始模型。四种模型分别提高了5.5%、1%、1%和1.5%。BP神经网络的所有指标均优于其他三种算法,其中,InSAR-BP神经网络的精度较InSAR-IM、InSAR-RF和InSAR-SVM分别提升了21%、8%和2%,BP神经网络的精度较IM、RF和SVM分别提升了25%、7.5%和1.5%,表明了BP神经网络的精度和合理性皆优于其他模型。
NOTES
*第一作者。
#通讯作者。