基于S型增长曲线组合预测模型的滑坡变形预测研究

2022-05-06贺小黑贺鑫焱彭必建

江西科学 2022年2期

彭鑫，贺小黑，贺鑫焱，彭必建

(1. 东华理工大学水资源与环境工程学院，330013，南昌；2. 中国人民武装警察部队研究院工程设计研究所，100012，北京；3.北京国信华源科技有限公司，100055, 北京；4.云南地质工程勘察设计研究院，650041，昆明)

0 引言

斜坡岩土体在构造、地应力等内在因素控制作用及降雨、人类工程活动等外在诱发因素共同影响下，沿着某一软弱结构面产生成块的滑移现象，称为滑坡。一旦下滑的滑体达到一定规模时，就可能造成其周边影响范围内的重大经济损失以及危害人员生命安全[1]。为了尽可能避免滑坡的危害，根据监测信息进行滑坡的预测预报成为滑坡防治工程的重要一环。

滑坡的预测预报研究包括时间尺度上的失稳时间预报和空间尺度上的变形趋势预测等内容[2]。学者们基于不同的理论对滑坡变形监测数据分析采用了多种方法，其中，在滑坡变形预测方面，学者们在基于灰理论、生长曲线理论的统计型模型以及依据突变理论、协同理论等非线性理论的非线性模型投入了非常大的精力，取得了丰硕的研究成果，比较典型的变形预测模型有基于生物生长曲线的S型生长曲线拟合模型(Verhulst、Pearl)[3-4]；基于灰理论的GM(1,1)模型[5-6]；基于机器学理论的神经网络模型[7-8]、支持向量机等模型[9-10]等。

然而，有些模型方法过程复杂，可能需要较深的数学理论知识储备(如协同理论[2,11]、混沌理论[12])，当模型的掌握不够深时容易造成模型预测精度不可控。Verhulst和Pearl模型在滑坡的变形预测和时间预报方面都有应用，其数学理论相对简单，有高等数学基础的人就能熟练运用。因此，针对Verhulst和Pearl单一预测模型预测精度的不足，本文基于加权组合预测的思想[13]，采用最优加权组合Verhulst和Pearl子模型对其进行优势互补，提高滑坡预测模型的变形预测精度，并通过已有滑坡实例对其进行验证。

1 滑坡的变形演化特征

通过总结大量的监测位移-时间曲线的形态特征，许强[14]等人将滑坡分为：稳定型滑坡、渐进型滑坡、突发型滑坡(图1)。对于突变型滑坡，有学者研究发现其一般具有较深的滑动面，且一般会发生液化现象，如修德皓[15]、许强[14]等人通过对甘肃黑方台滑坡研究发现，由于黄土固有的湿陷性，在降雨、地下水、灌水等因素下发生液化流动现象，是典型的突变型滑坡。

大量的研究表明，滑坡大多都有一定的变形

图1 不同滑坡演进类型的位移随时间的变化

演化过程[16-17]，现有的预测预报模型也大多针对渐进型这一类滑坡。渐进型滑坡的变形演化特征可分为缓慢变形、等速变形、加速变形及失稳破坏4个阶段[4]。晏同珍[18]、孙景恒[19]等学者认为滑坡的孕育、生长、成熟及消亡过程具有生物生长曲线类似的机制：滑坡的缓慢变形及等速变形阶段，滑坡主要以蠕滑为主，相当于生物生长的孕育阶段，加速变形阶段，滑带土的内摩擦角、抗剪强度不断弱化，曲线斜率逐渐增大，相当于生长模型曲线的生长阶段，滑坡进一步发展，变形急剧加速，曲线斜率呈陡崖式发展，然后趋于稳定，相当于生长模型曲线的成熟阶段[4,19]。用生物增长曲线取拟合滑坡历史监测数据对滑坡进行预测预报是可行的。

2 基于滑坡演化过程的预测模型

2.1 灰色Verhulst模型基本原理

Verhulst模型是1987年德国物理学家发现的一种生物生长模型，由于滑坡的位移特征与生物的生长规律类似，很多学者用Verhulst 模型对滑坡进行预报研究[3]。Verhuslt模型的白化微分方程形式为：

(1)

式中：a1，b1为系数。

对于一组等时距非负增量位移监测序列x(1)=[x(0)(1),x(0)(2),…,x(0)(n)]，经一次累加后得到原始累计位移监测序列x(1)=[x(1)(1),x(1)(2),…,x(1)(n)]。

对累加序列x(1)作紧邻均值生成序列z(1)，

式中

(2)

构造数据矩阵，用最小二乘法求解系数a1，b1，计算表达式为

(3)

(4)

2.2 以累积位移为参量的Pearl预测预报模型

Pearl预测预报模型最早由孙景恒提出并应用在新滩滑坡和意大利Vaiont滑坡监测数据上，取得了好的成果，少量学者利用Pearl模型开展了滑坡变形预测预报的研究，均取得不错的进展。Pearl生长曲线的数学模型表达式为：

(5)

式中：k为常数；f(t)为自变量t的多项式。

f(t)=a0+a1t+…

(6)

一般多项式阶数取1，则其数学表达式转换为一般logistic函数的表达式：

(7)

式中：a、b、c为待拟合系数；y为t时刻的位移拟合值。

因此，只需确定好a、b、c3个参数即可对滑坡进行拟合回归预测。采取非线性拟合求参的方法[4]，通过matlab中的Levenberg-Marquardt算法进行非线性拟合来确定a、b、c3个参数。

2.3 基于Verhulst增长曲线和Pearl增长曲线的组合预测模型

由于这些曲线各有不同，对某一滑坡的预测结果可能偏高或者偏低。因此，本文引入组合预测的思想，基于最小二乘原理将最优权重与前述2种单一模型组合在一起，通过对单一模型的取长补短以提高预测的准确性。

基于最优加权算法的组合模型计算流程如下。

令一组原始监测数据表示为xt,可用m个模型来拟合，拟合结果表示为xit，其中i=1,2,…,m；t=1,2,…,N。单个模型的权重可表示为wi，其满足下述要求：

(8)

组合预测模型的拟合值可表示为：

(9)

令eit为预测模型i在时刻t的拟合残差，则组合预测模型的拟合残差可表示为：

(10)

基于残差平方和最小的原则，通过最小二乘方法，最优权重wi可在下列约束条件下获取。

(11)

3 实例分析

3.1 模型输入数据的选取

依据文献[21]提供的卧龙寺新滑坡位移监测资料进行预测预报，监测数据见表1。新滑坡发生时间为1971年5月5日，因此，取最后2 d监测数据为验证数据用以预测，选取3月15日至5月3日经过平滑处理提取趋势项位移，其具有类似斋藤室内实验模拟滑坡累计位移的“三段式”特征。通过正态检验及相关性分析方法[5]综合判断得出1971年4月22号为滑坡体从等速变形阶段进入加速变形阶段的临界点，故选取1971年4月22日至1971年5月3日的监测数据进行预测预报。

3.2 单模型的预测分析

利用Matlab将前述Verhulst模型拟合选取的监测数据，通过式(2)～(3)得出模型参数a1=-0.204 6、b1=-0.009 5；将其代入式(4)，发现

表1 卧龙寺新滑坡监测数据

当i=2时，其相对误差最小。利用Matlab将前述Pearl模型拟合选取的监测数据，依据拟合优度最优原则得出模型参数a=0.004 974、b=0.398 9、c=-0.349 9，计算结果见表2，拟合曲线与监测数据曲线见图2。

表2 单预测模型拟合结果

图2 单模型拟合曲线与原始监测曲线的关系

图2虚线左侧是根据选取的建模数据得到的拟合值，虚线右侧是根据建模数据计算的预测值。从图2可以看出，根据选取的监测数据(39—50 d)进行拟合得到曲线与实测曲线相差不大，Pearl模型拟合的曲线较为精确，但两个模型在外推最后2 d得到的预测结果相差较大(一个较实测值偏大而另一个则偏小)，可以看出，往外预测期数越多，误差越大。由第52天预测结果可知，Verhulst模型的预测结果可能会起到一个提前触发临滑预警预报的效果，而Pearl模型则可能不会触发临界失稳预警预报。

3.3 组合预测模型预测分析

在得到前述单预测模型的拟合数据的基础上，根据式(8)～(12)构造基于Verhulst子模型和Pearl子模型的组合预测模型，由于本文只组合2个子模型，可令Verhust子模型的权重为w，则Pearl子模型的权重为1-w，计算得到的单预测模型权重见表3，单模型、组合模型的预测值及原始累计位移值之间的关系见图3。

表3 S型增长曲线子模型的权重

图3 组合模型与单模型预测结果对比

从图3可以看出，组合模型通过最优权重组合各子模型的优势，使得拟合值最大限度“逼近”实测值，尤其从第50期数据之后(模型外推预测值)更为明显，其预测值相较单模型预测值更加符合实际。

为了检验单模型及组合模型的预测精度，采取拟合优度R2及均方根误差RMSE 2个指标进行分析。分析结果见表4。

表4中的数据从拟合区精度、预测区精度以及综合精度3个方面对比分析，发现在选取建模数据基础上Pearl模型的拟合效果最好(拟合区)，但是在预测区域根据建模数据建立的拟合回归方程计算的预测值效果最差，而Verhulst模型在拟合区的拟合效果比Pearl模型差，但在预测区域却比Pearl模型的要好，虽然它的均方根误差也好大(表4，预测区90.456 1)，但是从图3中可以看出Verhulst模型能起到提前预警预报的作用(假设第51天或第52天的数据值为滑坡失稳时设定的预警阈值)。这表明单一预测模型根据滑坡监测历史数据拟合得到的回归方程可能不能很好地反映未来的滑坡变形趋势，有些模型可能会造成预警设备漏报或不报而滑坡已经发生的现象(如本文中得到的Pearl预测模型)。而通过最优加权组合单一模型的方法可以避免这一不足，拟合优度R2从单一模型的0.934、0.923提高到0.985，均方根误差RMSE从单一模型的4.69、4.72降低到1.29(表4中综合精度数据)。表明该组合预测模型优于上述的单一预测模型，基于最优加权的组合预测模型改进效果非常明显。