基于贝叶斯-粒子群算法的溜砂坡稳定性评价

娄超华，田荣燕，旺 久，孙威宇，罗 进

（西藏大学工学院，西藏 拉萨 850000）

0 引言

溜砂坡作为一类在我国西部地区常见的地质灾害，依据国际通用的DJ Varnes滑坡分类方法，属于广义滑坡中的一类特殊滑坡[1]。典型溜砂坡的形成一般都要经历砂砾石的产生，溜动及坡脚堆积三个过程[2]。作为一种在西藏地区常见的地质灾害，溜砂坡的发生常常会造成严重的后果[3]，因此对其稳定性评价分析是提前预防灾害发生的关键。目前在实际中对溜砂坡稳定性评判多通过力学计算，根据实地勘察所得主要影响因素及相关数据，采用静力平衡理论计算得出一个溜砂坡的稳定系数，根据计算的稳定系数来评价溜砂坡的稳定性[4]，本文尝试使用贝叶斯-粒子群算法分析溜砂坡的稳定性，并在此基础上将溜砂坡稳定等级划分为五个级别：稳定（Ⅰ）、较稳定（Ⅱ）、基本稳定（Ⅲ）、较不稳定（Ⅳ）和不稳定（Ⅴ）。

目前常用的溜砂坡稳定性评判方法有传统安全系数法、模糊综合评价法和层次分析法。传统安全系数法，即选取一评价因子作为溜砂坡的稳定状况，并将该参数与经验所得的安全警戒值相比较来确定其稳定状态[5]。但实践发现，该方法存在以下问题：（1）经验所得的安全警戒值具有较高的不准确性，这也使得到的最终结果不可靠；（2）该方法未考虑现实复杂环境下影响溜砂坡的多方面因素，因此所得出的结果仅仅代表某一特定情况下的稳定状况，不能反映溜砂坡整体是否安全[6]；（3）该方法并未考虑不同溜砂坡间存在的差异性，而是使用统一的安全系数来评价，结果自然不准确，将其运用到实际中必定存在很高的隐患。模糊综合评价法实际使用时存在计算繁琐、指标权重的确定有较强主观性的问题；更重要的是当其在权矢量和为1 的限制下当指标集个数U较大时，相对隶属度权系数往往偏小，导致权矢量与模糊矩阵R不匹配，最后结果出现超模糊现象，分辨效果很差，没有办法区分哪个的隶属度更高，甚至导致评判无效[7]。层次分析法存在当参数因子指标较多时数据统计量大，且权重很难确定的问题，由于层次分析法是模拟人们决策过程的思维方式的一种方法，因此必然带有较多的人为定性影响，导致客观定量数据少，定性成分多，结果令人信服程度低[8]。

因此选用一种更加科学准确的方法分析评判溜砂坡就显得很是重要。

1 贝叶斯-粒子群算法

1.1 贝叶斯网络

贝叶斯网络（Bayes Network/Bayesian Network，BN）又称为置信网络，是一种能够简单高效地处理不确定性问题的方法网络[9]。

BN 由三部分组成：条件分布概率组、结构和参数，其结构式是一个有向无环图，将变量作为其节点，结构就可以分成n个节点的集合m以及节点的边组成的集合。可得出其联合概率分布表达式：

式中：p(X1,···Xn)—— 表示一组随机变量 (X1,···Xn)的联合概率分布；

Xa—— 一个假设变量，且规定Xa是Xb的父节点，并用Pf(Xa)表示Xa的直接前驱节点联合。

BN 通过一个数据样本集进行参数和结构学习进而确定出网络的拓扑以及定义数值参数，从而形成一个能构建有效BN 的最优数据集，用贝叶斯评分 BIC函数（式2）来评价 B N的拓扑准确性，通过函数式（3）搜索全部可能的网络来优化 B N ，使用 A kaike信息判据确定惩罚项[10]，确保数据与结构的拟合准确度。

式中：J-表示训练数据样本集；

n-表示变量个数；

ma- 表示第a个 变量的取值数目；

ha- 表示第a个 变量的父节点的取值组合数目；

S abc- 表示第a个 变量取c时且其父节点取b时的样本数量；

S ab- 表示第a个 变量的父节点取b时的样本数量；

S-表示样本数量。

从以上公式可以看出，实际应用中n个 变量的网络连接数是十分庞大的，在如此多的数据变量中找到一个最优值是非常困难的，这也一直阻碍贝叶斯网络运用[11]的原因。

1.2 粒子群算法

粒子群算法（Particle Swarm Optimization, PSO）是一种进化计算技术，是基于对动物集群活动行为的观察基础上，根据动物群体的觅食行为利用群体中个体对信息的共享使整个群体运动在问题求解空间中产生从无序到有序的演化过程[12]，从而获得最优解，依据适应度判断解的质量。它没有遗传算法的“交叉”和“变异”操作，仅通过追寻当前所能搜索到的最优解来确定整体最优解。因其搜索速度快、效率高以及算法简单，被广泛运用在动态经济调度中。但其存在容易早熟收敛的缺点，并且容易陷入局部最优[13]。

粒子的位置和速度函数表达式：

式中：Xa、Ya、Zqa—— 分别为第a个 粒子在N维空间的位置、移动速度和最优位置；

Zq——全局最优位置；

w——惯性权重；

c1、c2——加速系数；

m1、m2——[0，1]区间内的随机值；

t——迭代时刻。

1.3 B N-PSO算法

2005年 DU 等[14]首次将 PSO应 用到 BN结构中，提出了BN-PSO算法[14]，采用 PSO来解决BN算法结构学习的问题，通过贝叶斯网络良好的推理能力来迅速构建起粒子节点和全局最优值之间的概率关系[15]。基于f1、f2、f3来更新粒子的离散速度和位置：

式中：m1、m2、m3——[0,1]区间内的随机值；

f1、f2、f3——分别对应式（5）的记忆、自身认知和群体认知三个部分；

f1out、f2out——分别对应式（7）和式（8）的输出结果；

M、C——变异和交叉运算，f2、f3中的运算逻辑与其相同（图1、图2）。

图1 M 运算逻辑图Fig.1 Logic diagram of M algorithm

粒子更新过程中，由于存在网络边缘粒子的插入和删除，可以将关联性差的粒子构成的回边去除，从而加速 B N 的学习[16]。惯性权重w可以确定留下来的粒子速度的数量，w越大，其搜索没有覆盖的区域能力相对越强，能够提高算法全局寻优能力和避免局部最优；w越小，其局部搜索能力相对越强，提高收敛速度。w的计算式为：

式中：w1、wn——分别表示初始权重和最终权重；

ev、evmax——分别表示当前迭代次数和最大迭代次数。

图2 C 运算逻辑图Fig.2 Logic diagram of C algorithm

同样，加速系数c1、c2如公式10、11 所示是线性变量，可以避免局部循环和早熟，且在后期能使寻找全局最优的进程加快。

式中：c11、c1n和c21、c2n—— 分别表示c1和c2的初始加速系数和最终加速系数。

根据以上算法绘出B N-PSO 的流程图见图3。

图3 BN-PSO 流程图Fig.3 BN-PSO flowchart

2 评价指标及评价权重的确定

2.1 评价指标的选择

影响溜砂坡稳定性的因素较多，综合起来可分为3 方面：即自然地质条件，外力作用和人为因素[17]。自然地质条件包括溜砂坡的坡度、坡高、植被覆盖率、地质活动和岩体风化等[18]。外力作用影响包括地震、暴雨冲刷、积雪融化、降雨和风作用等。人为因素包括不合理开挖坡脚、坡表植被破坏和动荷载震动等。但通过野外实地调研，同时结合拉萨本地的气候特点和地质环境综合分析，本文选取的拉萨本地溜砂坡灾害多为野外自然生成，人为因素影响较小，从地质条件和外力作用条件下选取主要评价指标，同时介于文章理论方法在溜砂坡稳定性评判方向上的研究尚处初级阶段，目前只选取了4 个主要因素即研究所取值：降雨量、坡高、坡度以及植被覆盖率。

经过三次实地调研测量和相关资料收集，得到所研究的各灾害点近期算例基础指标原始数据如表1所示。

2.2 评价权重的确定

为了定量计算主要影响因素大小顺序，利用信息熵[19-20]对其进行计算分析。将四个主要影响因素的评价指标分为5 级并赋值[21]，通过专家打分确定其权重，评价指标为[0, 1]区间内的值，0 代表影响作用最小，1 代表影响作用最大。溜砂坡评价指标等级划分及标准见表2。

评价指标通过打分已经全部变为无量纲值。再利用式（12）对随机抽取的12 个溜砂坡进行标准化处理。

表1 灾害点原始数据Table 1 Raw date of disaster points

表2 溜砂坡评价指标等级划分及标准Table 2 Classification and standardization of evaluation factors of debris slope stability

式中：Xa,b——表示第a（a=1, 2, 3, ···n）个样本的b（b=1, 2, 3, ···m）项评价指标的值；

Xmax(a)、Xmin(a)——分别为最大值和最小值。

然后根据下式求第b项评价指标的信息熵。

Pa,b——为评价指标b出现的频率。

再根据式（14）求出第b项评价指标的权重。

0≤wb≤1，wb越大，表示该指标在溜砂坡稳定性中的影响因素就越大。

各因素的信息熵和权重计算结果如表3所示。

由表3可知，各评价指标的权重从小到大为降雨量，坡度，坡高和植被覆盖率，同时能够看出：对溜砂坡稳定性影响的因素大小排序即为权重的排序。

表3 碎屑斜坡稳定性影响因素的信息熵及权重Table 3 Entropy and weight of evaluation factors of debris slope stability

3 实例分析

3.1 灾害点概况

文章使用了12 组主要位于途径曲水县、尼木县、仁布县、达孜县及墨竹工卡县等地的 G318国道沿线两侧以及拉萨至当雄县S109 省道沿线两侧的灾害点数据。分团队三次外出采集，调研地多属拉萨市辖区，各灾害点均具有典型溜砂坡灾害特征。G318 国道及S109省道拉萨市段主要为板岩、花岗岩、砂砾岩和玄武岩等，岩石在寒冻风化的内外共同强烈动力作用下破碎崩解，溜动堆积至坡脚，形成溜砂坡。拉萨市属于高原温带半干旱季风气候区，日照时间长，辐射较强，年变化气温较小，昼夜温差大，降雨量少且集中，干湿季节分明[22]。根据国家气象信息中心（http://cdc.cma.gov.cn）地面气候资料整编日值多年资料显示，拉萨市多年平均降雨量454.8 mm，多年平均最大降雨量607.2 mm，多年平均最小降雨量280.4 mm，6——9月的降雨量约占全年的90%。

3.2 评判值计算及稳定性等级划分

根据文献[23]中所使用的节点微网系统，文章采用节点微网建立适当的各微网数学模型列出相关数据见表1，将各溜砂坡灾害点相关数据代入算法程序，在BN-PSO 算法中，w1和wn分别为0.95 和0.4；c11和c1n分别为0.82 和0.5；c21和c2n分别为0.4 和0.83；evmax取10 000。

将上述参数代入MatlabR2016 中进行仿真分析，使用优化后的贝叶斯网络中的贝叶斯信息准则 BIC作为各稳定性系数指标评分函数，仿真分析得到各灾害点评分函数数据后，将各灾害点BIC 值按照从小到大的顺序排列，得到BIC 值与稳定等级关系成以顶点位置处BIC 值为0.400 0 的抛物线型分布，根据《滑坡防治工程勘察规范》（GBT 32864——2016）和《岩土工程勘察规范》（GB 50021——2018）将最大值与最小值间的间距等划分为五部分，即分为五类稳定性等级状态，其中当 BIC值远离（大于或小于）0.400 0 时，溜砂坡稳定性越差，当其BIC值越接近0.400 0 时，该灾害点的稳定性越高；在算法程序中，第一轮输出全局最优解 B IC值为灾害点2 所代表的0.408 9，即表明灾害点2 的稳定性最好，逐次循环进行排序，最后一轮循环所得 B IC值为灾害点6 所代表的0.682 3，表明灾害点6 的稳定性最差。求出各灾害点 BIC值进行划分等级，根据所得各灾害点BIC 值，可将灾害点分为五类稳定性状态，根据各灾害点 B IC值与0.400 0 的差值10%、11%、12%、13%以及超过14%为间距将边坡稳定程度分为五个级别：稳定（Ⅰ）、较稳定（Ⅱ）、基本稳定（Ⅲ）、较不稳定（Ⅳ）和不稳定（Ⅴ），即 BIC值在0.360 0～0.440 0 为稳定，0.350 0～0.360 0以及0.440 0～0.450 0 为较稳定，0.340 0～0.350 0 以及0.450 0～0.460 0 为 基 本 稳 定，0.330 0～0.340 0 以 及0.460 0～0.470 0 为较不稳定，小于0.330 0 以及大于0.470 0 为不稳定[24-27]。

3.3 稳定性评价结果

以上数据经计算和对比得到的处理结果如表4所示。

表4 灾害点数据处理结果Table 4 Date processing results of disaster points

由处理结果可知，12 个灾害点中处于基本稳定、较稳定和较不稳定状态的共有4 处，且有6 处灾害点处于不稳定，只有2 处是处于稳定状态，对比现场溜砂坡真实稳定状态，文章评价结果与实际稳定性状态相符合，这也证明了该算法运用在溜砂坡稳定性评价上的可行性[28-29]。拉萨地区特殊的地形地貌以及气候环境使得溜砂坡以较快的速度育和扩散，所以目前该地区溜砂坡灾害的发展趋势不容乐观，需引起重视。

4 结论

（1）文章首次将BN-PSO 算法运用在溜砂坡稳定性评判上，对比常见的模糊综合评价法、层次分析法，文章使用的BN-PSO 算法、信息熵赋权法在溜砂坡评价过程中无需确定隶属函数，评价因子的选择、评价因子的分层组合及判断矩阵的构建，因此结果评价计算、寻优都更加简单、高效，计算结果更为精确。

（2）证明了将粒子群算法引入贝叶斯网络进行溜砂坡稳定性评判问题的可行性，对溜砂坡的稳定性评判上有一定的参考价值。

（3）文章只针对一些典型的溜砂坡案例，阐述了引入粒子群算法在贝叶斯网络的寻优办法，将该方法运用在溜砂坡稳定性评价的研究尚还处在初级阶段，理论和方法还需不断完善，因此关于溜砂坡稳定性评价因素的选取和权重的确定都有待进一步探讨。