基于多元非线性回归的极震区泥石流物源量估算模型

任志刚,裴向军，顾文韬

(1.中国电建集团成都勘测设计研究院有限公司，四川 成都　610072;2.成都理工大学 地质灾害防治与地质环境保护国家重点实验室,四川 成都　610059)



基于多元非线性回归的极震区泥石流物源量估算模型

任志刚1,裴向军2，顾文韬2

(1.中国电建集团成都勘测设计研究院有限公司，四川 成都610072;2.成都理工大学 地质灾害防治与地质环境保护国家重点实验室,四川 成都610059)

汶川地震极震区高川乡、清平乡、映秀镇的28组泥石流样本，考虑沟域面积、相对高差、发震断裂距离、沟域岩性及岩体结构4个物源量影响因子，基于控制变量法找出单个因子与泥石流物源储量的相关关系，通过MATLAB多元非线性回归方法建立了综合因素影响下的极震区泥石流物源量估算模型。最后，利用非拟合数据以外的极震区6条沟对模型进行了验证，结果表明模型对于汶川地震极震区泥石流物源储量的估算具有较强的准确度及适宜性。

极震区；泥石流；物源；储量

0　引　　言

“5.12”汶川地震具有震级大、震源浅、持时长、烈度高等特点。震后，在强大的地震力作用下，直接触发了约15000处地质灾害[1]，沿龙门山主断裂带展布的狭长极震区(地震烈度大于Ⅸ地区)山体产生了大量崩滑灾害。大量崩滑堆积体及潜在震裂山体在泥石流沟域的剧增，使得降雨条件下极震区泥石流表现出明显的“灾后效应”，给灾区人民生命财产造成了严重的二次损害。

目前，国内外关于泥石流物源特征的研究方法也不拘一格，主要是在定性描述物源分布规律的基础上通过一定的方法对泥石流物源量进行估算。对于泥石流物源量估算的主要方法有形态评价法[3]、比例统计法[4]、灰色理论法[5]、遥感解译法[6]等，但这些方法极少与地震效应进行结合。因此，研究震后泥石流的物源量估算方法对于了解该类泥石流发展趋势、评价危害程度及科学治理具有重要意义。

1　研究区泥石流概况

由于汶川地震极震区震后泥石流表现出规模大、危害强、反复发作等特征，本文主要选取汶川地震极震区震后泥石流物源量较多、启动条件较低、危害性较大的地区作为研究区域进行特大地震后泥石流物源特征研究。主要选取的研究区为地震烈度为Ⅺ度的汶川县映秀镇岷江流域、地震烈度为Ⅹ度的安县高川乡高川河流域、地震烈度为Ⅹ度的绵竹市清平乡绵远河流域。其中映秀镇岷江流域研究区选取11条泥石流沟作为研究样本，高川乡高川河流域研究区选取9条泥石流沟作为研究样本，清平乡绵远河流域研究区选取8条泥石流沟作为研究样本。

2　极震区泥石流物源类型及发育特征

2.1物源类型

(1)地震崩滑堆积型。强震时地震波在坡体内产生的反复拉-压-剪切作用，导致斜坡以不同的动力破坏形式失稳，崩滑体进入沟道后以直接参与或间接补给的方式参与泥石流运动，是极震区震后泥石流物源量暴增的主要途径。

(2)震裂山体失稳型。地震作用下山体震裂损伤主要表现在两个方面：一是在地震波反复作用下岩体劈裂、结构面扩张而导致岩体结构劣化[7]，震裂山体在后期降雨作用及时效变形的作用下发生崩滑灾害；二是地震作用下坡表覆盖层产生震动松弛现象[8]，使得密实度降低，后期暴雨条件下产生了大量土体滑坡。

(3)震裂坡表侵蚀型。震后，坡面的侵蚀现象主要表现在两个方面：一是地震作用下坡表植被覆盖率降低，土体密实度减小，暴雨条件下覆盖层转换成坡面侵蚀物源；二是震后松散的崩滑堆积体在强降雨作用下产生坡面侵蚀现象，成为坡面侵蚀物源。该类物源的产生与带动多与泥石流水源的汇集同步发生。

2.2发育特征

极震区泥石流物源的发育特征是泥石流沟域斜坡自身发育特征以及地震作用、降雨作用等诱发因素的综合叠加作用结果。

(1)坡体内因。极震区泥石流沟道两侧的地形、侵蚀坡度及临空条件是斜坡失稳的重要因素，坡度的陡缓以及临空程度的高低控制坡面的张力范围以及坡脚的应力集中带[9]；斜坡岩土体的性质及组合方式，是斜坡应力集中、边界发展贯通以及坡体失稳的重要因素；泥石流沟道两侧坡体相对于河谷的高差大小将直接决定坡体地震加速度与振幅放大效应的程度。

(2)地震作用。地震发生后，地震波在介质中传播，能量呈对数衰减趋势衰减，根据震后地质灾害调查，极震区地质灾害沿主断裂带表现出明显的“距离效应”[10]；地震时逆冲发震断裂的地震加速度峰值上盘衰减较慢，而下盘衰减较快，研究表明汶川地震所诱发的上盘的地质灾害不仅在密度分布上，而且在规模上都大于下盘；发震断裂的局部“锁固段”在地震过程中被剪断、破裂，将会释放出较大的能量，产生强烈震动，形成“次级震源”，从而加强地表的地质灾害的发生。

3　极震区泥石流物源量影响因子

在考虑地震效应、岩土体结构特征、沟域规模、地质环境条件等情况下，将极震区泥石流物源量主控影响因子概括为：沟域面积、相对高差、离发震断裂的距离、岩性及结构特征。通过控制变量法找出各因子与物源储量的相关关系后，进行因子的多元非线性拟合。

3.1沟域面积

泥石流沟域面积对物源的影响主要表现在两个方面：首先，沟域面积概化了沟的规模，同等条件下沟域面积与物源储量呈线性增长(见图1)；其次，沟域面积还影响静储量—动储量的转化，极震区泥石流沟多年储量变化统计表明，沟域越大其汇水面积就越大，暴雨时坡表侵蚀量及沟道带动量就越多，静储量—动储量的转化量就越高。

图1　控制变量条件下沟域面积与总储量关系

3.2相对高差

研究表明，地震过程中斜坡端部地震峰值加速度(PGA)总体随高程增大[11]。为摒除泥石流沟底原始高差的影响，本文选择沟顶与沟底的相对高差作为研究依据。高川乡与映秀镇研究区分段高差与震后崩滑堆积体方量的统计数据表明：映秀镇和高川乡相对高差分别在1 200 m和1 000 m后高程端部放大效应开始突显，崩滑体平均体积与总体积呈现陡增的趋势，分段高差与崩滑体总体积总体呈指数关系增长(见图2)。清平乡研究区各崩滑灾害点的统计数据表明：崩滑体方量在相对高差800～900 m处出现陡增现象，崩滑体方量与相对高差总体呈指数关系增长(见图3)。

图2映秀镇研究区崩滑体方量与高差关系图3高川乡研究区崩滑体方量与高差关系

3.3发震断裂距离

对高川乡、清平乡各沟单位面积物源储量与断裂距离的统计结果表明：控制变量条件下，随断裂距离的增加，单位面积物源储量总体呈对数关系衰减(见图4)。

3.4岩性及岩体结构

根据沟域岩体结构特征采用各沟控制变量条件下的岩石坚硬程度、结构面张开度、结构面连通率、岩体裂隙度四个指标评价沟域岩体质量等级。

图4　高川乡、清平乡研究区距离效应拟合关系

(1)沟域岩体质量等级评级系统。极震区研究区各沟域岩体质量评价模型的建立主要运用层次分析法(AHP)[14]。

本文岩体质量评价模型所选择的因子为：沟域岩体坚硬程度、沟域岩体结构面张开度(e)、沟域岩体结构面连通率(k)、沟域岩体节理平均间距(d)四个因子。

将归一化后的矩阵进行求行和后，再进行列归一化处理，得到Aij的特征向量w：

在评价因子和因子权重值的基础上，结合地质力学系统的理论方法(地质力学系统，也称岩体权重系统，即“岩体评分”)，建立岩体质量等级评价系统(见表1)。

根据表1，沟域岩体的质量等级指数N为：

aibi

式中a——为评价因子影响强度评分值；

b——为评价因子的权重值。

将各沟域岩体质量等级评价因子的参数值代入岩体质量等级评价系统，计算得出研究区各沟域岩体的质量等级指数N值见表2。

表2　各沟域岩体质量等级指数

注：N值越大则岩体质量越差，越容易发生崩滑灾害。

(2)沟域岩体质量与物源储量的相关关系。本文对上述10条沟单位面积物源储量进行统计，结合各沟域岩体质量等级评价指数，得出沟域岩体质量与单位面积物源储量整体呈线性关系(见图5)。

图5　沟域岩体质量与单位面积储量关系

4　极震区泥石流物源量估算模型

泥石流物源储量决定了泥石流的规模与危险等级，是研究震后泥石流的基础。本文利用极震区三大研究区28个沟道随机样本(见表3)，在已分别建立的沟域面积、相对高差、断裂距离、沟域岩体质量与总储量相关关系(见第3章)的前提下，利用多元非线性回归方法建立四个因子与总储量的关系，再建立总储量与动储量关系，得到极震区泥石流动储量估算模型。

4.1基于MATLAB的多元非线性回归模型

4.1.1回归模型的建立

根据以上对控制变量条件下极震区物源储量与沟域面积、相对高差、断裂距离和沟域岩体质量的关系分析可知：

表3　汶川地震极震区三大片区泥石流沟样本物源参数总汇

(1)在控制其余变量条件下，沟域面积与总储量呈线性正相关(见图1)；(2)控制变量条件下，相对高差与总储量呈指数正相关(见图2、3)；(3)控制变量条件下，断裂距离与总储量呈对数负相关(见图4)；(4)控制变量条件下，沟域岩体质量与总储量呈线性正相关(见图5)。

根据以上分析及控制变量条件下的各因子拟合结果，极震区泥石流物源总储量(V总)回归模型可概化为：

V总=a1X1+a2ea3X2+a4ln(X3)+a5X4+a6

式中X1——为沟域面积,km2;

X2——为相对高差,m;

X3——为断裂距离,km;

X4——为沟域岩体质量指数(无量纲);

a1、a2、a3、a4、a5、a6——为特征系数。

4.1.2特征系数的确定

通过MATLAB多元非线性拟合运算程序，将表5中Xi与Yi的28组数据输入，将ai赋入初值，并进行非线性回归计算后，得到系数矩阵ai为：

ai={42.300 60.039 20.002 4

-0.026 90.458-1.986 2}T

将系数代入回归模型，可得极震区泥石流总储量计算模型：

V总=42.300 6X1+0.039 2e0.002 4X2-

0.026 9ln(X3)+0.458X4-1.986 2

上式计算模型回归检验结果如表4，r>r0.05(28)，F>F0.05(4，28)，检验结果表明回归模型相关性较好，拟合效果较好[16]。

4.1.3误差分析

将表3中各沟域参数代入上式，可得回归模型

表4　非线性回归结果及检验

注：r表示拟合相关性，F表示回归效果显著性检验的统计量值，r0.05(28)、F0.05(4，28)检验临界值。

计算值，与实测值相比较，得出各点计算值与实测值的误差关系如图6。

图6　实测值与回归值关系

绝对误差模型绝对误差与相对误差见图7、8。

从以上误差结果统计中可知，除样本6、8、12外，其余样本误差均较小，符合相关性标准。其中，样本6为高川乡三叉沟，样本8为高川乡黄洞子沟，样本12为清平乡文家沟。根据许强等对汶川地震大型滑坡发育规律的分析[17]，清平—高川段为地震主断裂的转折和错裂末端部位，具有“锁固段效应”，因此该区域为断裂压应力高度集中区，地震时会释放出更多能量。而研究区三叉沟、黄洞子沟和文家沟正好位于映秀—北川断裂高川—清平段锁固段内，震时三沟内分别诱发了火石沟滑坡、大光包滑坡和文家沟滑坡三个巨型滑坡，从而大大增加了此3个样本的物源储量。

图7样本绝对误差关系图8样本相对误差关系

4.2极震区泥石流物源动储量估算模型

调查及实际研究表明，极震区泥石流物源总储量与动储量存在一定的相关性。将三大研究区28个泥石流样本中实测物源总储量与动储量进行相关性拟合[18]，得到相关关系如图9。

图9　极震区泥石流物源总储量与动储量拟合关系

得到极震区三大片区泥石流物源总储量V总与动储量V0拟合关系如下：

V0=0.295V总-5.192 (R2=0.855)

将上式代入极震区泥石流物源总储量计算模型，可得极震区泥石流物源动储量估算模型：

V0=12.479A+0.011 56e0.002 4H-

0.007 936lnL+0.135 11N-5.777 93

式中V0——为极震区泥石流沟域物源动储量;

A——为泥石流沟域面积;

H——为泥石流沟域相对河流高差;

L——为沟域离发震断裂的垂直距离;

N——为沟域岩体质量等级指数(可根据沟域岩体质量等级评级系统求得)。

4.3模型准确性及适宜性检验

为检验模型的适宜性，随机选择三大研究区以外的极震区6条泥石流沟作为适宜性检验样本，根据所得的极震区泥石流动储量预测模型，代入相应参数(见表5)进行计算及误差分析。

表5　极震区泥石流物源预测模型适宜性检验样本参数汇总

将适宜性检验样本参数代入评价模型，得出泥石流物源动储量计算值，与实测值比较，其误差关系见图10。

图10　极震区泥石流物源预测模型适宜性检验误差关系

由图10可知，除北川县城以外，其余样本相对误差及绝对误差均较小。北川县西山坡沟沟内发育有一个110万m3的大型滑坡[19]，具有难以预测的袭扰因素。摒除袭扰因素，本模型适宜于评价、估算极震区未发育有巨型滑坡泥石流沟域的物源动储量。

5　结论与建议

本文主要大量现场调查资料及统计数据对汶川地震后极震区泥石流物源的类型、分布及发育特征、影响因素及估算方法等方面展开论述研究。结合前人研究，初步概括了汶川地震极震区泥石流物源的分布和发育特征，初步研究了泥石流物源量的影响因子，并建立了极震区物源估算模型对极震区泥石流物源量进行预测研究。主要概括为以下几点：

(1)根据前人研究成果分析了极震区物源的形成机制，通过现场调查总结了极震区震后泥石流物源类型，并从坡体内因以及地震断裂效应两方面分析了震后泥石流物源的形成机制及发育特征；

(2)在物源发育特征分析的基础上总结了影响物源发育的主控因子，其中包括沟域面积、相对高差、发震断裂距离、沟域岩性及岩体结构等。通过大量统计数据分析了在控制变量作用下各因子与物源量的函数关系，在此过程中建立了沟域岩体质量等级评价模型；

(3)在各因子显著性分析的基础上，将以上建立的各因子的函数关系概化为多元非线性方程，通过MATLAB程序语言将极震区三大片区泥石流物源因子及物源量的28组数据代入方程进行回归分析，从而确定回归方程的系数，最终通过相关性验证及误差分析；

(4)通过极震区拟合数据以外的另外6条沟对模型适宜性进行检验，检验结果显示模型适宜性较强。

[1]殷跃平.汶川八级地震地质灾害研究[J].工程地质学报,2008,16(4): 433-444.

[2]冯文凯,黄润秋,许强.地震波效应与山体斜坡震裂机理深入分析[J]. 西北地震学报.2011,33(1):20-25.

[3]周必凡.泥石流防治指南[M].北京:科学出版社,1991.

[4]蒋忠信.泥石流固体物质储量变化的定量预测[J].山地研究,1994,12(3):155-162.

[5]乔建平，黄栋，杨宗佶，孟华君. 汶川地震极震区泥石流物源动储量统计方法讨论[J]. 中国地质灾害与防治学报，2012，6(2):1-6.

[6]Chuan Tang, Zhu Jing. Approach for urban debris flow risk assessment [J]. Advances in Water Scient ,2006,17(3):67-71.

[7]唐川，朱静.城市泥石流灾害风险评价探讨[J].水科学进展，2006,17(3):67-71.

[8]Newmark N M. Effects of earthquakes on dams and embankments [J].Geotechnique 1965，15(2)：139-160.

[9]张倬元,王士天,王兰生,黄润秋,许强,陶连金.工程地质分析原理(第三版)[M].北京:地质出版社,2008.

[10]黄润秋,李为乐.汶川大地震触发地质灾害的断层效应分析[J].工程地质学报,2009,17(1):19-27.

[11]邹威.强震作用下均质斜坡动力响应的大型振动台试验研究[D].成都理工大学,2011.

[12]于海英,王栋,杨永强,卢大伟,解全才,张明宇,周宝峰,江纹乡,程翔,杨剑.汶川8.0级地震强震动特征初步分析[J].震灾防御技术.2008,3(4):321-336.

[13]Huang Runqiu ,Li Weile.Fault effect analysis of geohazards triggered by Wenchuan Earthquake ,China[J]. Bull. Of Engineering Geology,2009, 17(1):19-88.

[14]许树柏.实用决策方法-层次分析法原理[M].天津:天津大学出版社.1988.

[15]袁进科.斜坡震裂岩体结构特征与震后崩塌识别体系研究[D].成都理工大学,2012.

[16]崔圣华,裴向军,黄润秋.直线型斜坡滚石运动速度特征研究[J].工程地质学报.2013,21(6):912-919.

[17]许强,李为乐.汶川地震诱发大型滑坡分布规律研究[J].工程地质学报.2008,3(4):818-826.

[18]顾文韬,裴向军,裴钻,张雄,李天涛.极震区泥石流物源特征研究[J].自然灾害学报.2015,24(2):107-114.

[19]唐川,梁京涛.汶川震区北川9.24暴雨泥石流特征研究[J].工程地质学报,2008,12(6):751-758.

2016-03-07

任志刚(1972-)，男，四川仁寿人，高级工程师，从事水力水电工程勘测设计工作。

P642.23

B

1003-9805(2016)03-0056-06