赵 瑛 王显彪 陈菊林 刘阳阳 杜文学 彭 亮
(1.青海省水文地质工程地质环境地质调查院,西宁 810008;2.青海省水文地质及地热地质重点实验室,西宁 810008;3.青海中煤地质工程有限责任公司,西宁 810008)
在地球地表形态演化过程中,形成了一定的古滑坡,受后期不利扰动因素影响,易在局部诱发次级滑坡活动,潜在威害较大,对其研究具有重要意义[1-3];滑坡稳定性及变形潜势分析可为灾害防治提供较强的理论指导,此两方面已成为滑坡灾害的研究热点.在滑坡稳定性研究方面,多利用数值模拟进行滑坡稳定性分析[4-5],忽略了滑坡土体的不均匀性,普适性存在不足,且考虑到滑坡稳定性影响因素相对较多,利用层次分析法构建其稳定性评价体系具有较强的可操作性[6],值得深入研究.在滑坡变形潜势分析方面,考虑到滑坡变形是其稳定性的直观体现,进而利用变形预测来实现滑坡变形潜势分析是可行的;受监测条件限制,滑坡变形数据含有一定的误差信息,且Kalman滤波能有效分解滑坡变形数据[7],适用性已被广泛验证,加之相关向量机[8]及马尔科夫链[9]也被广泛应用于岩土领域,因此,基于前述三类方法构建滑坡变形预测模型也是可行的.综上,以王家坡滑坡为工程背景,首先,利用层次分析法和德尔菲法构建滑坡稳定性评价体系,开展其稳定性评价;其次,基于滑坡变形监测成果,先利用Kalman滤波实现滑坡变形数据的去噪分解处理,将变形数据分解为趋势项和误差项,再利用优化相关向量机进行趋势项预测,并利用马尔科夫链实现误差修正预测,以达到高精度预测的目的.
由于滑坡稳定性影响因素相对较多,故而其稳定性评价模型具多样性和全面性;考虑到层次分析法能合理、系统地分析滑坡稳定性影响指标,因此,利用其构建滑坡稳定性评价模型是可行的.同时,1-9标度法已被广泛应用于评价指标的权值求解,确定其为评价指标的权值求解方法;德尔菲法能归纳、总结、反馈专家经验,并能实现定量评价,进而利用其确定评价指标的隶属度.故将稳定性评价模型的构建流程总结为:滑坡稳定性评价体系的构建;评价指标的权值求解;评价指标的隶属度求解;稳定性定量评价.
1)滑坡稳定性评价体系的构建
结合相关文献研究成果及实例,利用层次分析法构建滑坡稳定性评价体系,将滑坡稳定性评价体系确定为3层结构,其中,目标层为“滑坡稳定性评价”,下设两层评价指标,一级指标有4个,分别为滑坡几何因素、外部因素、内部因素和人类活动因素,二级指标共计13个,具体如图1所示.
图1 滑坡稳定性评价体系
2)评价指标的权值求解
1-9标度法的权值求解过程又可细分为3 个步骤,即判断矩阵的构建、一致性检验和权值计算.
①判断矩阵的构建.通过两两对比评价指标间的相对重要性(评价标准见表1),构建出判断矩阵,且判断矩阵属互反矩阵,矩阵元素间相对独立.
表1 评价指标间的重要性判别标准
②一致性检验.由于判断矩阵是由评价指标两两对比得出,难以保证所有指标间的相对重要性具一致性,需对判断矩阵进行一致性检验.先计算判断矩阵的最大特征值λmax,再计算一致性比率值CR:
式中:RI为误差控制值;n为评价指标数.
利用CR值即可评价判断矩阵的一致性,判据为:若CR值<0.1,则判断矩阵满足一致性检验;反之,不满足一致性检验,需对判断矩阵进行重新构建,直至满足一致性检验.
③权值计算.当判断矩阵通过一致性检验后,对最大特征值λmax对应的特征向量进行归一化处理,所得归一化值即为对应的权值.
3)评价指标的隶属度求解
德尔菲法是一种反馈匿名函询法,已被广泛应用于评价指标的隶属度求解过程中,其大致流程为:组建专家组,并向专家组提交问题背景材料及待回答问题;专家结合问题及背景材料,对问题进行书面答复,答复过程具相对独立性;对专家答复进行汇总统计,并将汇总结果返回专家,专家结合汇总结果进行第二次修改;最后,结合二次修改回复进行统计.
值得指出的是,由于不同专家的工作经验存在一定差异,则其回复结果的可信度也应区别对待,进而在专家回复意见统计过程中,对相应专家的评价结果应进行折减处理,具体标准见表2.
表2 专家回复结果的可信度折减标准
4)稳定性定量评价
以前述评价指标的权值和隶属度为基础,利用P×C分级法实现滑坡稳定性的定量评价;同时,在稳定性评价过程中,结合滑坡稳定性状态划分标准,将滑坡稳定性等级划分为4级,具体划分标准见表3.
表3 滑坡稳定性等级划分标准
滑坡变形潜势指的是滑坡体内部潜在应变能的强弱,因此,提出利用变形预测来实现滑坡变形潜势分析.预测过程可细分为3个阶段:其一,是利用Kalman滤波实现变形数据的去噪分解处理;其二,利用CSO-RVM 模型实现趋势项预测;其三,利用MC 模型实现误差修正预测.
1.2.1 变形去噪模型的构建
据文献[7],在滑坡位移监测过程中,受不确定因素影响,滑坡变形监测值含有一定的误差信息,即滑坡变形数据可分解为趋势项和误差项,即:
式中:yt为变形监测值;γt为趋势项;εt为误差项.
Kalman滤波是一种动态数据处理方法,可在最小误差准则基础上实现最优估计,适用于滑坡数据的去噪分解处理;由于文献[7]已详述Kalman滤波的基本原理,限于篇幅,不再赘述.同时,伴随理论与计算机技术发展,Kalman滤波也衍生出了多种类型,如标准型、自适应型、抗差自适应型和半参数型,各类型Kalman滤波的适用性存在一定差异,为保证滑坡变形数据的去噪分解效果,提出对4类Kalman滤波均进行去噪效果分解处理,效果最佳者即为滑坡变形数据的去噪分解方法.
同时,在变形数据的去噪分解过程中,提出以信噪比、均方根误差为基础,构建综合评价指标g:
式中:g1为信噪比的归一化值;g2为均方根误差的归一化值.
1.2.2 趋势项预测模型的构建
相关向量机(Relevance Vector Machine,RVM)在保留核心特征向量的同时,还不受Mercer定理约束,使之具有高效的预测能力,因此,提出利用其构建趋势项预测模型.
由于RVM 模型的基本原理已在文献[10-11]中详述,限于篇幅不再赘述其基本原理;但是,在RVM模型的应用过程中,核函数对预测效果具有直接影响,如核函数类型较多,不同核函数的适用性存在一定差异;同时,核函数宽度也直接影响预测效果.为保证预测精度,对上述两问题均应进行优化处理,即:
①核函数类型的优化处理
RVM 模型的常用核函数有RBF核、Sigmoid核、高斯核和ERBF核,4类核函数的基本原理是有差异的,为实现其优化处理,提出对4类核函数均进行试算处理,预测效果最优者即为RVM 模型的核函数.
②核函数宽度的优化处理
鸡群优化算法(Chicken Swarm Optimization,CSO)是于2014年提出的一种新型智能算法,其主要模拟鸡群等级制度和觅食行为,具体优化流程为:先设置一定规模的鸡群,包括公鸡、母鸡和雏鸡,并以单只公鸡为中心,将鸡群划分为若干子群;由公鸡引领,带领子群进行觅食行为,以实现寻优处理.
1.2.3 误差修正模型的构建
前述优化处理虽有效保证了趋势项的预测精度,但也难以完全刻画趋势项的变化特征,加之前述变形数据去噪分解产生的误差项,有必要进行误差修正预测.具体过程为:将趋势项的预测误差与前述数据去噪分解产生的误差项叠加,组成新的误差序列,并利用马尔科夫链(Markov Chain,MC)实现误差修正预测.根据MC模型的基本原理,可将误差修正流程详述如下:
①状态划分.对误差序列的分布范围进行统计,并进行状态划分,结合该文实例,划分5个分析状态.
②状态转移概率计算.以其中一个状态为基础,计算该状态与另一状态间的转移次数,并进一步计算该状态在总样本数中的出现频次,前者与后者相除,即可得到该状态的转移概率.
③构建构造矩阵.以转移概率为基础,构建出概率矩阵,并设置转移步长,以求得相应状态条件下的预测结果,达到误差修正预测的目的.
最后,将趋势项预测结果与误差修正预测结果叠加,所得值即为滑坡变形的最终预测值,并以其结果判断滑坡变形潜势,即:若滑坡变形呈持续增加趋势,则其变形潜势趋于不利方向发展;反之,若滑坡变形趋于稳定,则其变形潜势趋于有利方向发展.
王家坡滑坡位于西安市灞桥区,具白鹿塬塬边斜坡地貌,后缘高程为740.6 m,前缘高程为453.8 m,斜坡坡度介于15°~35°,高差相对较大;据钻孔资料,滑坡区上覆土层以崩坡积碎石土为主,下覆基岩以新近纪砂、泥岩为主;滑坡区未见明显地表水体,但受地形条件限制,在强降雨条件下,易于地表形成径流,对坡体水土流失具有一定影响,地下水则主要划分为孔隙水和裂隙水两类,主要接受大气降雨补给.
王家坡滑坡亦属古滑坡活动诱发产生的次级滑坡,其成因主要为:古滑坡形成后,重力促使其向下滑动,并形成反倾洼地,利于水流汇集,对坡脚浸泡,加之地表形态改变,促使地下水位发生变化,由此长期累积,诱发了王家坡滑坡次级滑坡的形成.
据现场调查,王家坡滑坡纵向长约300 m,宽约150 m,面积约3.7万m2,平均厚度约18 m,体积为66.6万m3,属中型土质滑坡.王家坡滑坡威胁17户64人,人口密度相对一般,且区内人类活动强度主要是以道路建设、房屋建设、耕种为主,受地形条件限制,活动强度一般.
近年,受不利因素诱发影响,滑坡出现了显著变形(如图2所示),如:①滑坡中部道路产生裂缝和不均匀沉降,且具持续变形特征;同时,道路外侧土体出现局部剥落和滑动迹象.据现场调查,变形裂缝的延伸长度约7 m,宽度间于4~25 cm,下错高度3~5 cm,可见深度约10 cm.②居民房屋出现明显裂缝,且裂缝具贯穿迹象,墙面变形明显,裂缝延伸长度约5 m,裂缝宽度约10 cm,可见深度3~7 cm.
图2 滑坡变形迹象照片
为切实掌握滑坡变形特征,在滑坡体上布设了变形监测点,其中,2号监测点位于滑坡中部居民区;3号监测点位于居民区前缘.监测时间由2018年10月17日至2019年8月15日,共计监测301 d,由于监测数据采集过程具非等距特征,不利于后续预测分析,为解决该问题,利用Matlab拟合工具箱对变形数据进行3 次样条拟合,并按7d一个周期,进行等距处理,共计得到43个周期的变形数据,具体变形曲线如图3所示.
图3 滑坡位移变形曲线
以前述图1滑坡稳定性评价体系为基础,先利用1-9标度法实现滑坡稳定性评价指标的权值求解,限于篇幅,以一级指标为例,详述权值求解过程.通过两两对比A1~A4指标的相对重要性,得一级指标的判断矩阵见表4.
表4 一级指标的判断矩阵
通过计算,得其λmax1=4.015,对应特征向量b1=[0.151,0.269,0.829,0.467],再计算得到CR=CI/RI=0.015/0.89=0.017<0.1,满足一致性检验.
对b1进行归一化处理,得A1~A4指标的权值为:[0.088,0.157,0.483,0.272].
再对其余二级指标进行权值求解,得到所有评价指标的权值见表5.
表5 滑坡稳定性评价指标的权值结果
其次,利用德尔菲法实现滑坡稳定性评价指标的隶属度求解,共计统计27位专家的反馈结果,所得隶属度结果见表6.
表6 滑坡稳定性评价指标的隶属度结果
最后,结合评价指标的权值和隶属度求解结果,进行滑坡各级评价指标的稳定程度评价.
①二级指标的稳定性评价
结合表3中的稳定性权值和表6中的隶属度,求得各二级指标的稳定性得分见表7.由表7知,各二级指标的稳定性得分存在一定差异,说明其稳定性程度是不同的,其中,70分以上的评价指标有9个,所占比例为69.23%;60~70分的评价指标有1个,所占比例为7.69%;低于60分的评价指标有3个,所占比例为23.08%;通过对比,得二级指标以70分以上得分为主.同时,按稳定性等级评价,Ⅲ级指标有10个,所占比例为76.92%,Ⅱ级指标有3个,所占比例为23.08%,得二级指标的稳定性等级以Ⅲ级为主.
表7 二级指标的稳定性评价结果
②一级指标的稳定性评价
对一级指标的稳定性得分进行计算,结果见表8.由表8知,4个一级指标的稳定性得分也存在一定差异,其中,以外部因素A2的稳定性得分相对最高,其次是滑坡几何因素A1、内部因素A3和人口活动因素A4,4者均为Ⅲ级.
表8 一级指标的稳定性评价结果
③滑坡稳定性的整体评价
同理,进一步开展王家坡滑坡的整体稳定性评价,通过计算得其稳定性得分为69.51分,属Ⅲ级,基本稳定.
最后,结合文献[12]的研究成果,通过强度折减法计算得到的滑坡稳定性系数约为1.24~1.25,与本文稳定性判定结果较为一致,且本文评价相对保守,安全储备相对更高,侧面验证了本文评价模型的合理性.
在变形预测过程中,以1~38周期为训练样本,39~43周期为验证样本,外推预测4个周期;同时,为验证去噪过程及优化过程的有效性,以2号监测点为例,详述去噪过程和优化过程.
首先,开展滑坡变形数据的去噪分解处理,得到不同类型Kalman滤波的去噪分解效果见表9.由表9知,不同类型Kalman滤波的分解效果存在一定差异,其中,半参数型的综合指标g相对最大,说明其分解效果相对最优,其次是自适应型、抗差自适应型和标准型,因此,确定滑坡变形数据的去噪分解模型为半参数型Kalman滤波.
表9 不同类型Kalman滤波的去噪分解结果
为进一步验证Kalman滤波在滑坡变形数据去噪分解处理过程中的适用性,先统计上述Kalman滤波综合指标g的均值,再对若干cof小波和sym 小波的分解效果进行类似统计,得到三者的分解效果如图4所示.由图4可知,3类去噪模型的分解效果存在一定差异,Kalman滤波的分解效果明显优于其余两类小波的去噪分解效果,因此,验证了Kalman滤波适用于滑坡变形数据的去噪分解处理.
图4 不同模型的去噪分解效果对比
其次,进行滑坡变形预测,以分析其变形潜势,其预测过程如下:
①核函数优化处理.对四类核函数的预测结果进行统计,结果见表10.由表10可知,四类核函数的预测效果存在明显差异,其中,ERBF 核的相对误差均值为2.58%,相对最小,其次是Sigmoid核、RBF 核和高斯核,说明ERBF 核的预测效果相对最优,确定其为RVM 模型的核函数.
表10 核函数的优化筛选结果
②核函数宽度优化处理.在核函数类型优化处理的基础上,再利用CSO 算法优化核函数宽度,经优化后的预测结果见表11.由表11可知,经CSO算法的优化处理,2号监测点趋势项的预测误差介于2.17%~2.53%,相对误差均值为2.40%,与优化前相比,预测精度有一定的提高,验证了CSO 算法的优化处理效果.
表11 核函数宽度的优化预测结果
③误差修正预测.前述已实现2号监测点的趋势项优化预测,再利用MC 模型实现其误差修正预测,经修正预测后的预测结果见表12.由表12知,在2号监测点的预测结果中,最大、最小相对误差分别为2.11%和1.86%,相对误差均值为1.97%,相较误差修正前,预测精度明显提高,充分说明了利用MC 模型进行误差修正预测的必要性和有效性.
表12 2号监测点的最终预测结果
④外推预测及变形潜势分析.前述已利用2号监测点详述了数据去噪分解效果和各阶段的优化处理效果,再对其进行外推预测,并对3号监测点进行类似预测分析,得两者预测结果见表13.由表13可知,3号监测点预测结果的最大、最小相对误差分别为2.05%和1.84%,相对误差均值为1.93%,预测效果与2号监测点的1.97%相当,初步验证了本文预测思路具有较高的预测精度.同时,由两监测点的外推预测结果可知,两者变形仍会进一步增加,且无收敛趋势,说明滑坡变形潜势趋于不利方向发展,需尽快进行灾害防治处理,避免成灾损失.
表13 最终变形预测结果
通过前述分析,王家坡滑坡的稳定影响因素主要属Ⅲ级稳定性,具基本稳定状态;少数影响因素属Ⅱ级稳定性,具欠稳定状态;综合评价滑坡整体稳定性为Ⅲ级,暂属基本稳定状态.但是,经变形潜势分析,滑坡变形仍会进一步增加,趋于不利方向发展,建议加强综合防治.
通过对王家坡滑坡的稳定性评价和变形潜势分析,主要得如下结论:
1)通过整体稳定性评价,得滑坡稳定性得分为69.51分,属Ⅲ级,具基本稳定状态.值得指出的是,不同评价指标的稳定性得分及等级存在一定差异,以Ⅲ级稳定性为主,说明不同评价指标对滑坡稳定性的贡献程度是有差异的,侧面验证了通过层次分析法构建滑坡稳定性评价体系是可行的.
2)通过变形预测,得CSO-RVM-MC 模型具有较高预测精度,预测结果的平均相对误差小于2%,验证了其在滑坡变形预测中的适用性,且外推预测显示其变形仍会进一步增加,无收敛趋势,得其变形潜势趋于不利方向发展.
3)限于篇幅,本文仅从变形预测角度开展滑坡变形潜势分析,建议在条件允许前提下,可进一步引入变形趋势判断方法或增加分析监测点数,以进一步佐证分析结果的准确性.