基于流域形态特征分析的北京北山泥石流易发性评价因子研究

黄来源

（北京市地质灾害防治研究所，北京 100120）

泥石流是北京山区最严重、最具破坏性的地质灾害，曾造成大量人员伤亡和重大财产损失，具有暴发突然、危害大等特点（赵忠海，2009；刘连刚，2015）。泥石流易发性研究是泥石流灾害防治领域重要组成部分，是泥石流发育机理研究的关键，是开展泥石流隐患早期识别和泥石流监测预警模型建立的前提条件（南赟等，2020；黄来源，2020）。泥石流系统的物源系统和沟道系统在空间上具有一定依附性和相关性，统一称为影响泥石流过程的本底因素。本文泥石流隐患的易发性分析评价主要是基于自然沟道的本底因素开展的。

谭炳炎教授1986年负责的泥石流沟严重程度的数量化评判研究主要是在西南地区铁路建设期间识别泥石流隐患工作过程中开展的，经过了大量实践检验，属于专家现场判别经验的归纳总结，最终确定影响一次泥石流过程的三大因素15个指标，是对泥石流过程中涉及非水源条件的各类因素的全面归纳总结，为后来的泥石流易发性规律研究奠定了基础（罗守敬等，2021；王海芝，2019），促进了泥石流易发性评价从定性评价到定量评价的发展。现行泥石流过程易发程度评价指标中存在不可量化指标、难以量化指标、可量化指标三类指标混合评价的问题。1）现状泥石流易发评价体系中存在不可量化定性指标，如沟口地貌指标；沟道横剖面特征指标、沟道堵塞程度指标；地质因素中岩石类型指标、不良地质现象指标等，这些指标无法进行单位性数量化，只是根据经验人为规定赋分分级边界，非量化指标赋分评价客观性不足，对泥石流整体易发性评价结论造成消极影响。2）现状泥石流易发评价体系中存在量化困难或难以量化指标，如河沟因素中补给长度指标、充淤变形指标；地质因素中松散物平均厚度指标、松散物储量等指标；这些指标难以通过现有的技术手段精确获取，实际操作中，多以现场调查目视估测为主得到指标数据，指标数据置信度低。

目前开展基于沟道地形因素影响泥石流易发性研究主要有：Yu Bin等（2011）通过统计方法来探讨地形因素和泥石流发生的关系；朱渊（2014）、禹磊等（2016）和朱云波等（2015）通过现场调查讨论了汇水区面积、主沟长度、形状系数和沟床坡度对泥石流易发性的影响；马煜等（2018）通过讨论地形因子对特定流域泥石流发生的影响，建立了地形因子。可见沟道地形因素与泥石流易发程度存在一定相关性规律。本文采用统计和调查法，在前人研究成果的基础上，提出基于可量化指标的易发性因子，通过选取位于曾发生过泥石流的典型流域的若干泥石流单沟隐患样本，分析沟道形态特征中可量化指标与泥石流易发性之间的规律，建立沟道形态特征参数指标与泥石流易发性指标的数学关系模型，在提出沟道特征综合因子基础上，提出基于沟道特征等可量化指标的泥石流易发性因子。

1 怀北镇大沙河流域地质环境条件

怀北镇主要泥石流隐患分布于境内大沙河流域上游，大沙河流域属海河流域上游四级支流，为境内最大的河流，从椴树岭至神山流经境内大水峪，流出山进入平原区。

大水峪水库以上的大沙河流域上游流域面积达到55.6 km2，大水峪水库面积0.7 km2，常年水面面积为52.2万m2，总库容1 460万m3。

大沙河流域属暖温带半湿润大陆性季风气候，四季分明，雨热同期，夏季暖热湿润，冬季寒冷少雪。年平均降水量635.3 mm，年平均降雨日数为93.4 d，最多达147 d（1990年），最少为63 d（1997年）。极端年最大雨量1120.8 mm（1969年），极端年最少雨量354.4 mm（1962年）。降雨集中在每年6月至8月，7月、8月最多。

大沙河流域中上游处于燕山山脉南麓，地貌上属于黑坨山、云蒙山、猴顶山主峰线南侧的中低山地貌单元，属于北京北部传统的降雨中心区域，区域性降水特征基本一致。区域内岩性主要是由中生代火山活动形成的花岗岩、闪长岩、安山岩组成；地质构造主要以北南方向的断裂构造为主，区域地质特征基本一致保证了区域沟道物源条件基本相同（冒建，2015）。

2 泥石流隐患沟道形态特征参数

2.1 沟道形态特征参数分析选择

根据前人研究成果，用于描述泥石流沟的沟道形态特征的因素主要包括反映沟道平面形态特征的流域面积、流域长度等直接参数和流域形态系数、流域不对称系数等间接参数；反映沟道垂直形态特征的主沟长度和相对高差等直接参数和沟道纵比降等间接参数。

描述沟道形态特征的主要指标见表1。

表1 泥石流沟道流域形态特征指标Tab.1 Characteristic index of shape characteristic factors of debris flow gully

2.2 研究区域沟道形态特征参数特征值

1）沟道形态特征参数特征值

根据《北京市怀柔区突发地质灾害详细调查报告（1/50 000）》数据（北京中色资源环境工程有限公司，2014），位于怀柔区怀北镇大沙河流域的29条典型泥石流隐患的沟道形态特征因素特征值见表2。

表2 大沙河流域泥石流隐患沟道形态特征因素特征值Tab.2 Characteristic value of shape characteristic factors of debris flow gully

2）沟道形态特征可量化指标参数特征值帕累托图分布特征

通过沟道形态特征可量化指标参数特征值帕累托图分析（图1），可以得出研究样本的沟道形态指标的数据分布特征，研究样本沟道流域面积在［0.06，1.66］区间的沟道数量占比80%；主沟长度在［0.55，2.95］区间的沟道数量占比90%；主沟平均宽度在［0.12，0.76］区间的沟道数量占比83%；流域相对高差在［181，901］区间的沟道数量占比99%；主沟纵坡在［8，233］区间的沟道数量占比83%；主沟纵比降在［0.14，0.34］区间的沟道数量占比83%；形态因子在［0.13，0.85］区间的沟道数量占比90%；主沟伸长比在［0.41，0.81］区间的沟道数量占比90%；而且通过形态因子和伸长比可以人工剔除研究样本中的3个离散度较大的样本（Kf=1.26、Kf=2.59、Kf=2.66）。

图1 反映流域形态特征指标数据分布特征Fig.1 Distribution characteristics of dimensionless watershed area index

3 泥石流隐患沟道形态特征参数之间关系

3.1 标准流域面积

沟道流域面积决定了集水区域的大小和表层风化物源的规模，因此流域面积与泥石流发生的规模是正向相关的。为了消除支沟差异的影响，对直接参数的流域面积进行无量纲化，提出标准流域面积Adl。

式中，Adl为标准流域面积，Ad为流域面积，Am为标准流域面积特征值。

通过统计怀北镇大沙河上游流域26个泥石流隐患沟道流域面积指标数据，研究样本沟道流域面积最大值6.3，最小值0.18，算术平均值1.39，平均偏差0.91，中位数0.96。研究样本流域面积Ad主要集中于［0.18，1.78］区间，占比80%（图2）。设定本次研究样本标准流域面积特征值Am=0.96 km2作为样本流域面积无量纲化的单位值。研究样本的标准流域面积见表3。

表3 大沙河流域泥石流隐患沟道标准流域面积统计表Tab.3 Statistics of standard drainage area of debris flow hidden trough in Dasha River Basin

图2 流域面积指标分布特征Fig.2 Distribution characteristics of dimensionless watershed area index Adl

3.2 形态因子

沟道形态因子与沟道流域内汇水快慢程度相关性高，该参数越接近于1，说明流域形状接近于正方形，利于流域内地表径流快速集中。该参数为无量纲参数。

式中，Kf为无量纲形状因子，Ad为流域面积，L为流域长度。

根据图3，形态因子和流域面积指标分布分散，经数据趋势拟合优化，形态因子Kf与标准流域面积Adl之间经幂函数拟合，决定系数R2=0.036 8，决定系数R2取值范围为0～1，R2等于1或接近1时，说明趋势线的估计值与对应的实际数据之间的拟合程度高，趋势线的可靠性就越高，反之则可靠性较低。根据拟合结果，形态因子与流域面积之间关系不明确。沟道形态因子不能完全独自代表沟道形态特征。

图3 形态因子Kf与标准流域面积Adl的关系Fig.3 Relationship between form factor Kf and watershed area index Adl

3.3 沟道伸长比

流域形状伸长比Er，面积等于流域面积的圆的直径与流域长度的比值称伸长比，伸长比越小，流域越趋于狭长形。

沟道形状伸长比是面积等于流域面积的圆的直径与流域长度的比值，为无量纲参数，根据水文学原理，狭长形沟道，主沟形成汇流时间相对扁形沟道要长，不利于地表径流的快速汇集。

式中，Er为无量纲沟道伸长比，A为流域面积，L为流域长度。

根据优化后的趋势线拟合结果，决定系数R2=0.036 8，相对较低，伸长比Er与标准流域面积指标Adl之间的相互关系不明确（图4）。

图4 伸长比Er与标准流域面积Adl的关系Fig.4 Relationship between form factor Er and watershed area index Adl

3.4 沟道纵比降

沟道纵比降是沟道相对高差与沟道长度的比值，为无量纲参数，根据能量守恒定律，该参数与泥石流过程中泥石流的流动速度和能量正向相关。该参数越大，泥石流造成的冲击力和破坏程度越大。

式中，Gv为无量纲沟道纵比降，H为相对高差，L为流域长度。

根据图5，沟道纵比降和标准流域面积指标分布基本呈现幂函数关系，出现了随着流域面积指标的增大，沟道纵比降指标出现降低的规律，纵比降Gv与流域面积Adl按照幂函数拟合趋势线呈现一定相关性，样本分布于曲线两侧。

图5 沟道纵比降Gv与标准流域面积Adl的关系Fig.5 Relationship between vertical gradient Gv and watershed area index Adl

4 沟道形态特征综合因子和泥石流易发性因子的建立

根据水文学原理，流域形态因子指标对于沟道汇水速度及物源启动产生主要影响，沟道纵比降与一次山洪泥石流过程的能量等级相关。沟道形态因子与沟道纵比降控制着泥石流过程启动和过程特征。根据前述分析成果，尝试探索通过引入沟道形态特征综合因子Fg，使其即包含沟道水平方向形态特征的形态因子Kf，又包括反映沟道垂直方向形态特征的沟道纵比降Gv，设定Fg=Kf*Gv，综合特征因子Fg与标准流域面积Adl的关系如图6所示。根据优化拟合结果，Fg=0.063 3Adl-0.202，形态综合因子Fg与标准流域面积指标Adl之间具有明显的相关性分布规律，样本分布于曲线两侧。标准流域面积越大，形态特征综合因子也相应变小。

图6 沟道形态特征综合因子Fg与标准流域面积Adl的关系Fig.6 Relationship between factors of morphological characteristics Fg and watershed area index Adl

根据泥石流定义和发生原理，流域面积决定着研究范围内汇水量的多少和物源的多少，流域形态因子指标对于沟道影响汇水速度及物源启动特征，沟道纵比降与山洪泥石流过程的能量等级相关。综上所述，沟道形态因子与沟道纵比降控制着泥石流过程的运动特征，沟道流域面积影响泥石流物质中的水量和岩土体物源量。

一次完整的泥石流过程原则上需要在时间轴上具备2个前置空间系统条件，也可以称作本底因素，即沟道系统和物源系统，物源系统又存在于沟道系统中，尝试探索通过引入泥石流易发性因子Is，使其既包含形态因子和沟道纵比降，又能反映流域面积影响物源的重要作用，初步根据因子权重相当原则，代表物源系统因子取Adl，设定泥石流易发性因子公式为Is=Fg*Adl和Is=Kf*Gv*Adl。

5 泥石流易发性因子与大沙河流域泥石流隐患特征的关系

研究区大沙河流域上游流域面积达到55.6 km2，大沙河流域中上游处于燕山山脉南麓，地貌上属于黑坨山、云蒙山南侧的中低山地貌单元，属于传统的北京北部降雨中心区域，区域性降水特征基本一致。区域内岩性主要是由中生代火山活动形成的花岗岩、闪长岩、安山岩组成；地质构造主要以北南方向的断裂构造为主，区域地质特征基本一致保证了区域沟道物源条件基本相同。

以2018年8月精细调查项目怀北镇该流域的29条泥石流沟为基础数据，根据泥石流易发性评价N值均高于60，说明这些隐患沟道均是泥石流隐患，处于轻度易发和易发等级。

根据同一区域的泥石流发生的年代数据，14条已明确发生过泥石流灾害的沟道作为校核研究样本，爆发泥石流的年代集中于1940年、1969年、2002年。这些具有明确发生历史的隐患沟道作为真实泥石流样本数据，其余15条泥石流隐患为定性评价确认的隐患。剔除形状因子异常数据3个，本流域26个研究样本的泥石流易发性因子与组成指标的关系如下：

泥石流易发性因子Is与标准流域面积Adl的关系如图7所示。根据优化拟合结果，Is=0.0633Adl0.7981，决定系数R2=0.800 2，泥石流易发性因子Is与标准流域面积指标Adl之间具有明显的相关性分布规律，标准流域面积越大，泥石流易发性因子指标也越大。

图7 易发性因子Is与标准流域面积Adl的关系Fig.7 Relationship between susceptibility factor Is and watershed area index Adl

泥石流易发性因子Is与形态因子Kf的关系如图8所示。根据优化拟合结果，Is=0.261 5Kf1.0119，决定系数R2=0.331 6，泥石流易发性因子Is与形态因子指标Kf之间相关性分布规律弱于其与流域面积之间的规律，形态因子指标越大，泥石流易发性因子指标也越大。

图8 易发性因子Is与形态因子Kf的关系Fig.8 Relationship between susceptibility factor Is and morphological factor Kf

泥石流易发性因子Is与沟道纵比降Gv的关系如图9所示。根据优化拟合结果，Is=0.019 9Gv-0.851，决定系数R2=0.215 7，石流易发性因子Is与沟道纵比降指标Gv之间分布规律相关性一般，与形态因子指标相当，沟道纵比降指标越大，石流易发性因子指标越小。这明显与前面提出的沟道纵比降指标与泥石流发生因素正相关的常识相悖（魏永明等，1998）。所以继续分析沟道纵比降Gv与标准流域面积Adl的关系。

图9 易发性因子Is与沟道纵比降Gv的关系Fig.9 Relationship between susceptibility factor Is and vertical gradient Gv

沟道纵比降Gv与标准流域面积的关系如图10所示。根据优化拟合结果，本流域研究样本沟道纵比降指标Gv与标准流域面积Adl之间呈现减函数的分布规律，沟道纵比降指标越大，标准流域面积指标越小。结论与调查实际情况吻合。

图10 沟道纵比降Gv与标准流域面积Adl的关系Fig.10 Relationship between vertical gradient Gv and watershed area index Adl

以上分析说明大沙河流域的泥石流隐患（包括实际发生过的和人为定性确认的）影响泥石流发生的本底因素中，相关性最强的因素是沟道流域面积，其次是沟道形态因子，最后才是沟道纵比降。沟道纵比降与流域面积负相关程度高于其与泥石流发生的正相关程度，造成沟道纵比降作为流域面积的反向因子与泥石流易发性负相关。通俗地讲，本地区泥石流隐患基本特征是流域面积较大、平面形态非狭长、沟道总坡度较小。或者说，沟道流域平面形态特征指标更大程度上影响着本地区泥石流的易发性程度。

6 结论和建议

6.1 结论

1）传统泥石流易发性评价指标中存在不可量化参数、难以量化参数和可量化参数。易发性评价标准分N值客观性不足，不同人员对不同沟道的评价N值主观性高。泥石流沟道特征因素具有基础数据客观、可量化、获取可操作性强的特点。在相同气象系统条件下，泥石流沟道特征因素与泥石流的发生过程存在特定的联系。

2）通过对大沙河流域泥石流隐患的沟道形态特征统计分析，总结了研究样本的总体形态特征；挑选出流域面积、形态因子、沟道伸长比、沟道纵比降等客观数据指标参与沟道形态特征综合因子的建立；提出标准流域面积指标，建立基于沟道形态因子、沟道纵比降的流域形态特征综合因子Fg，Fg*Adl0.202=0.063 3。

3）影响泥石流过程本底因素包括沟道系统和物源系统，在建立沟道系统特征参数Fg的基础上，尝试探索通过引入泥石流易发性因子Is，使其即包含形态因子和沟道纵比降，又能反映流域面积影响物源的重要作用，泥石流易发性因子公式为Is=Fg*Adl和Is=Kf*Gv*Adl。

4）通过对大沙河流域泥石流隐患样本的分析发现影响泥石流易发性的本底因素中，最主要因素是沟道流域面积、其次是沟道形态因子、最后才是沟道纵比降。沟道纵比降与流域面积负相关程度高于其与泥石流发生的正相关程度，造成沟道纵比降作为流域面积的反向因子与泥石流易发性负相关。

5）通过基于沟道流域形态特征的泥石流沟易发性评价研究，发现以沟道流域形态特征为代表的环境背景因素与泥石流过程的相关联系，建立沟道各异的沟道流域形态特征因子，为开展涵盖水源系统因子、物源系统因子和沟道系统因子三元结构的泥石流预警模型研究奠定基础。

6.2 建议

本次研究样本限定于特定的小流域，并且样本数量有限，结论存在空间局限性，建议开展更大范围研究。

由于物源系统条件数据的置信水平低，造成建立泥石流易发性因子中物源因素只是简单地用标准流域面积指标代表物源系统条件，从而造成泥石流易发性因子对流域面积的依赖度较高，建议开展物源系统条件指标的建立和获取研究。