气候变化背景下中巴经济走廊崩滑灾害危险性动态评价及预测

余文秀， 李秀珍， 姚杰， 郑玲静

(1.中国科学院、水利部成都山地灾害与环境研究所， 成都 610041； 2.中国科学院大学工程科学学院， 北京 100049)

崩塌滑坡是山区最为常见的地质灾害，世界银行数据显示，全球约370万km2的陆地处于崩塌滑坡的高易发区内，崩塌滑坡灾害威胁着3亿多人的生命财产安全[1]。近年来，全球气候变暖直接引发了冰川冻土退化、地表蒸发量增加、极端降雨次数增多，也间接导致世界各地发生崩滑灾害的数量及破坏程度显著上升。这不仅严重威胁着人类的生命和财产安全，也造成了严重的环境破坏[2-3]。中巴经济走廊以中巴公路为导向，纵贯帕米尔高原，北接“丝绸之路经济带”，南连“21 世纪海上丝绸之路”，是中国经巴基斯坦通往南亚和中东地区的重要通道[4-5]。走廊沿线地质条件复杂、生态环境脆弱、区域稳定性极差，是全球崩滑灾害高发区之一，也是气候变化最为敏感和复杂的地区。随着全球气温的升高，致使走廊内冰雪融水剧增，一方面导致河流径流量增多，加剧了河流对于岸坡的侵蚀；另一方面，融雪水的入渗会降低岩土体的抗剪强度，促进崩滑灾害的发育，使中巴经济走廊的建设和民生安全受到严重威胁。为了保证走廊的顺利建设，降低因发生崩滑灾害而造成的人员伤亡和经济损失，急需对该区域崩滑灾害进行系统科学的危险性评价及预测。

常用的危险性评价方法分为定性评价及定量评价两类[6]。定性评价主观性较强，是通过专家打分对主要因素进行评判和确定，常见的有层次分析法、专家评分法等[7]；相比而言，定量评价更为客观和准确，常见的有信息量法、逻辑回归法、支持向量机法、频率比法等。其中信息量及逻辑回归两种模型由于操作简单、准确率高，在区域崩滑灾害危险性评价中应用尤为广泛，如文献[8-10]利用信息量模型及ArcGIS软件完成了区域崩滑灾害的风险评价；文献[11-13]采用逻辑回归分析分别对汶川地震滑坡、巫山县滑坡及思南县滑坡进行了危险性评估，证实了模型的可靠性。上述信息量及逻辑回归模型虽在崩滑灾害危险性评价中应用广泛，却各有优缺点，如信息量模型仅考虑了同一因子不同等级对灾害发生的影响，忽略了因子间的相关性。而逻辑回归模型虽无法考虑同一因子不同等级对灾害的影响，却可以得到不同影响因子的权重系数。二者联合形成的信息量-逻辑回归耦合模型优势互补，能明显提高预测精度。

在因子选取方面，已有研究大多只考虑了相对稳定的地质环境因素，得到的也是静态不变的评价结果[14]。然而，随着全球气候变暖，降雨、气温等诱发因素的动态变化必然会导致崩滑灾害的危险性随之变化。科学评估和预测崩滑灾害的动态危险性，对精准防范山地灾害具有重要意义[15]。

现充分考虑气候变化及冰雪融水对崩滑灾害危险性的动态影响，选取与崩滑灾害发生密切相关的评价因子，基于信息量-逻辑回归耦合模型对中巴经济走廊北部山区崩滑灾害的危险性进行静态和动态评价及预测，以期为中巴经济走廊崩滑灾害的防灾减灾提供科学决策和指导。

1 研究区自然地质环境及灾害概况

中巴经济走廊北部山区位于33°N～40°N、71°E～78°E，起于中国喀什，终于巴基斯坦伊斯兰堡北部(图1)。北接“丝绸之路经济带”、南连“二十一世纪海上丝绸之路”，是贯通南北丝路的关键枢纽，是一条包括铁路、公路、光缆和油气通道的贸易走廊，也是“一带一路”的重要组成部分[16]。

图1 研究区地理位置图Fig.1 Geographical location map of the study area

1.1 自然地质环境概况

地形地貌是形成地质灾害的内在因素，也是影响崩滑灾害发育的重要因素，对地质灾害的发生起着控制性作用[17]。中巴经济走廊地形高差悬殊、地势险峻、群山耸立、侵蚀切割强烈、冲沟发育密集、切割深度大、峡谷地形众多，其特殊的地形条件为崩滑灾害的发生提供了良好的孕灾环境。

研究区地质条件复杂，寒武纪到第四纪各时代的地层均有出露(图2)。研究区东北部以花岗岩、闪长岩及镁铁矿石为主，南部以灰岩及石英砂岩为主，整体稳定，崩滑灾害发育较少；中西部以泥岩、砂岩、粉砂岩、膏泥岩、石英片岩等易泥化或软化的岩石以及第四系松散堆积体为主，稳定性极差，崩滑灾害最为密集。本研究将区内岩体依据硬度及完整程度划分为5个岩组：Ⅰ第四系松散堆积物；Ⅱ泥片岩；Ⅲ砂岩、泥质砂岩、粉砂岩、泥岩、石英片岩；Ⅳ灰岩、白云岩、石英砂岩、镁铁质矿；Ⅴ花岗岩、闪长岩等。

中巴经济走廊地处喜马拉雅山、喀喇昆仑山和兴都库什山三大山脉交汇处，自新生代以来，欧亚板块与印度板块相互碰撞导致该区成为大陆内部新构造运动最为活跃的地区之一，活动断裂发育，地震活动频繁。据地震局1902年至今存档数据，区内共发生地震17 996次，主要为3～5级的小震，占总数的89.26%(图3)。频繁的地震及密集发育的活动断裂使得研究区内岩体破碎，区域稳定性显著降低，为崩滑灾害的发育提供了有利条件[18-19]。

研究区属暖温带大陆性干旱气候类型，垂直分带明显，高山冰川与干热河谷气候截然不同。中、低海拔为典型的山地气候，夏季炎热湿润，冬季寒冷干燥。海拔5 000 m以上的群峰被永久积雪和冰川覆盖，常年积雪。区内水系发达，支流众多，崩滑灾害大多沿河流两岸发育。主要河流如印度河、洪扎河、吉尔吉特河等均以高山冰雪融水补给为主，以地下水和降雨补给为辅。通过分析近十几年来的气候数据，发现中巴经济走廊内年均气温及降雨量呈现明显的上升趋势。融雪水及降雨量的逐年增加，加剧了河流对岸坡的侵蚀，也降低了岩土体内部的稳定性，致使区内崩滑灾害的发育愈加剧烈。

图2 地层与断裂带分布图Fig.2 Distribution map of strata and fault zones

1.2 崩滑灾害发育特征

中巴经济走廊内典型的高山陡坡地貌、强烈的构造运动、脆弱的生态环境及敏感的气候条件致使崩滑灾害十分发育。现通过文献资料收集、遥感解疑和现场调查查明走廊北部山区共发育崩滑灾害1 062处，其中滑坡401处，崩塌661处。中巴公路(Karakoram Highway，KKH)沿线及重要路段或工程影响区以野外调查及资料收集为主、遥感解疑为辅，其他区域以遥感解译及资料收集为主、野外调查为辅。统计分析灾害规模及数量，区内共发育大型、特大型崩塌137处，大型、特大型滑坡111处，中小型崩塌524处，中小型滑坡290处(图4)。

图4 研究区灾害点空间分布图Fig.4 Distribution map of geological disasters

2 评价模型和评价因子

2.1 评价模型

灾害的发生受各种地形及地质环境因子的影响，不同的因子对于灾害发生所起的作用大小不同[20]，对于同一影响因子，其值的变化也会直接影响灾害发生概率的大小。信息量法认为事件的发生受到多种因素共同作用，且同一因素的不同等级对于事件发生与否的影响存在差异，其理论公式为

(1)

式(1)中：I=(Y,x1,x2,…,xn)为影响因素x1,x2,…,xn对灾害提供的信息量；P(Y,x1,x2,…,xn)为影响因素x1,x2,…,xn共同作用下灾害发生的概率；P(Y)为灾害发生的先验概率。

当涉及的因子较多，且同一因子不同等级对灾害的影响也不同时，采用原式计算较为复杂，故先计算单一因素的信息量，再将各因素的信息量值进行求和，得到综合信息量[21]，原理论模型简化为

(2)

式(2)中：I为评价单元内的综合信息量；Ni为在影响因素Xi内的灾害面积；N为研究区内灾害分布总面积；Si表示影响因素Xi内的研究区面积；S为研究区总面积；H表示滑坡事件。

计算所得信息量值介于-1～1，值越大表示发生灾害的可能性越大；反之，值越小，越不易触发灾害；值为0时表示无法判断灾害是否发生。

逻辑回归模型侧重于分析不同影响因子与灾害发生的相关性大小[22]，属于多变量模型，其最明显的优势在于自变量无需满足正态分布，且自变量可以是多个，但因变量只能有两种情况。理论公式为

(3)

式(3)中：z=w1x1+w2x2+…+wmxm+b是一个加权线性组合模型，其中b为回归函数截距，wm为回归系数，xm为自变量；F(z)为事件发生的概率，其值介于0～1，值越大发生事件的概率越大。

信息量法由于意义明确、操作简单、实用性强，在地质灾害评价领域得到了广泛应用，但信息量模型只能反应同一因子不同等级对灾害发生的影响程度，忽略了因子间的相关性。逻辑回归模型以灾害是否发生作为因变量，影响因子值作为自变量，可以得到影响因子的权重大小。为了充分利用两种模型的优势，建立信息量-逻辑回归耦合模型，以影响因子归一化后的信息量值为自变量，以灾害发生或者不发生作为因变量，既能反映不同因子间的相对重要性，又能反映同一因子不同等级对于灾害发生的影响[23]。其理论公式为

(4)

式(4)中：Z为危险性值；Ii为第i个因子的分级指标提供的信息量；Ii，max为第i个因子的信息量最大值；Ii，min为第i个因子的信息量最小值；βi为第i个因子的逻辑回归系数；β0为逻辑回归所得常数项。

2.2 评价因子选取

科学、合理地选择评价因子是保证评价结果可靠性的前提。中巴经济走廊内发育的崩滑灾害与其特殊的地貌、地质、水文、生态、气候环境密不可分。本研究选取相对高差、坡度、坡向、曲率、岩组、断裂带密度、河网密度、距冰川距离、地震峰值加速度(peak ground acceleration,PGA)、土地利用10个因子(图5)作为静态危险性评价因子体系，在验证静态危险性评价结果可靠性的基础上通过引入动态降雨及气温与静态评价因子构成动态危险性评价因子体系。

3 崩滑灾害静态危险性评价

3.1 信息量计算

根据式(2)计算出各评价因子的信息量值I，如表1所示。

表1 信息量计算表Table 1 Calculation results of information content

信息量值有正有负，正值越大发生灾害的概率越大；负值越小发生灾害的概率越小，由此得到中巴经济走廊北部山区对崩滑灾害发生影响较大的因子类别：岩组为Ⅱ，坡度为50°～60°，PGA为0.46～0.5，断裂带密度为0.065～0.098 km/km2，曲率为-7～0.5，距冰川6～12 km，坡向为南东及南西，相对高差为1 130～1 410 m，河网密度为0.077～0.1 km/km2，土地利用类型为草地和灌木地。以上几种因子类别中发生崩滑灾害的概率较大。

3.2 基于耦合模型的静态危险性评价

中巴经济走廊北部山区共发育崩滑灾害1 062处，通过Arcgis随机生成等量的非灾害点并结合遥感影像进行调整。以各点信息量为自变量，以灾害发生与否为因变量，导入SPSS软件进行逻辑回归分析。将表2中回归系数代入式(4)中得

Z=0.195+0.531ITD+1.123IPGA+0.834IBC+0.563IPD+0.572IHW+0.475IPX+0.74IQL+0.493IGC-0.087IDLD+0.825IYZ

(5)

式(5)中:Z表示危险性系数；ITD、IPGA、IBC、IPD、IHW、IPX、IQL、IGC、IDLD、IYZ分别表示归一化后的土地利用、PGA、距冰川距离、坡度、河网密度、坡向、曲率、相对高差、断裂带密度、岩组的信息量。

表2 逻辑回归系数计算结果Table 2 Calculation results of logistic regression coefficients

结合前人危险性分级的经验，利用Arcgis自然断点法将危险性系数分为1～5级，分别对应极低、低、中等、高、极高危险性。统计分析各级危险区面积及灾害个数，结果显示危险等级越高，崩滑灾害数量及灾害点密度越大，约有89.92%的灾害点分布于高、极高危险区，在低、极低危险区内的灾害点仅占2.32%(表3)，评价结果与实际吻合，表明耦合模型用于该区崩滑灾害危险性评价具有较高的可信度。危险性区划图(图6)显示，中巴经济走廊北部山区崩滑灾害的高危险区主要集中在中西部地区，特别是红其拉甫至哈维连一带的坡度陡峭山区和地形起伏大的河谷两岸。

4 崩滑灾害动态危险性评价及预测

崩滑灾害的发生与气候变化密切相关，据历史记载，崩滑灾害的发生往往紧随极端降雨，降雨是诱发灾害最重要的因素之一。研究区内气温升高引起区域范围内融雪水剧增，成为了诱发崩滑灾害的另一重要因素，故本研究同时引入动态降雨和气温两个因子进行动态危险性的计算。

表3 危险性等级与崩滑灾害分布对比Table 3 Comparison of risk grade and disaster distribution

图6 危险性区划图 Fig.6 Hazard zoning map

4.1 气候数据来源及因子选取

气候数据下载自中国科学数据库，包括日均降雨、日最高温及日最低温三种类型，时间跨度为1961—2015年，覆盖了整个中巴经济走廊，是基于中巴经济走廊及其周边地区的65 个气象站点数据，以数字高程模型(digital elevation model, DEM)为协变量，采用 ANUSPLIN 软件进行空间插值得到的结果。 已有研究大多选取平均气温或温差来探索地质灾害与气温变化的关系，本文需选取能直接反应冰川融化程度的因子，因此引入了新的概念：一年内最低温高于0°的天数。当温度高于0°时，冰川消融；温度低于0°时，冰川积累[24]。本次通过ARCGIS软件处理得到各年最低温高于零(Tmin>0)的天数分布图，如图7所示。

图7 部分年份降雨及Tmin>0分布图Fig.7 Distribution map of rainfall and Tmin>0 in some years

4.2 历史危险性动态评价

以10年为间隔，分别选取1985年、1995年、2005年及2015年进行历史危险性动态分析。首先计算气候因子的信息量(表4)。

表4 气候因子信息量Table 4 Information quantity of climate factors

再将所有动态评价因子进行逻辑回归分析，得到新的回归系数与动态危险性计算公式为

H=0.189+0.533ITD+1.148IPGA+0.842IBC+0.555IPD+0.583IHW+0.471IPX+0.74IQL+0.5IGC-0.116IDLD+0.836IYZ+0.812IJY,i+0.791IQW,i

(6)

式(6)中：H为动态危险性值；ITD、IPGA、IBC、IPD、IHW、IPX、IQL、IGC、IDLD、IYZ分别表示归一化后的土地利用、PGA、距冰川距离、坡度、河网密度、坡向、曲率、相对高差、断裂带密度、岩组各个因子的信息量值；IJY,i为第i年归一化后的降雨信息量值；IQW,i为第i年归一化后的气温信息量值。

通过式(6)计算出1985年、1995年、2005年及2015年的危险性系数，采用自然断点法对危险性系数分级得到危险性区划图(图8)。动态危险性区划图显示，研究区不同年份的崩滑灾害危险区分布范围随降雨及气温变化而出现明显差异。如 2005年吉尔吉特以西出现极端强降雨，导致西部高危险区面积显著增加。动态危险性结果表明，崩滑灾害的危险性随降雨量增加、气温升高而增加。

图8 1985—2015年动态危险性分布图Fig.8 Dynamic risk distribution map from 1985 to 2015

4.3 危险性预测

基于雨量站点及随机点构成的提取点提取出1961—2015年多年的气候数据，通过SPSS软件进行时间序列预测分析，得到未来不同年份的降雨及气温值，再利用ArcGIS中的IDW插值处理得到预测降雨和气温图层(图9)。再将预测降雨及气温代入式(6)计算得到预测危险性区划图(图10)。

图9 部分年份预测降雨、气温及提取点分布图Fig.9 Distribution map of predicted rainfall, temperature and extraction points in some years

统计分析1985—2035年各年的较高危险区(极高、高危险区)面积占比并绘制变化趋势图(图11)。

图10 2035年危险性区划图(预测)Fig.10 Hazard zoning map for 2035 (forecast)

图11 1985—2035年较高危险区占比分布图Fig.11 Distribution of proportion of high-risk areas from 1985—2035

结果表明：研究区崩滑灾害较高危险区占比随气候变化呈现波动式上升趋势，预示未来中巴经济走廊内崩滑灾害的发生将会更加频繁，崩滑灾害发生的风险也将会随着气候变暖趋势的增加而逐步增大。

5 结论

以中巴经济走廊北部山区为研究区，采用信息量-逻辑回归耦合模型进行区域崩滑灾害危险性评价，得到以下结论。

(1)选取相对高差、坡度、坡向、曲率、岩组、断裂带密度、河网密度、PGA、土地利用、距冰川距离10个因子构成静态危险性评价因子体系。通过耦合模型计算得到研究区静态危险性分布图，分布图显示研究区较高危险区集中在中、西部地区，特别是红其拉甫至哈维连段的坡度陡峭的山区和地形起伏大的河谷两岸。统计分析灾害点在各危险等级内的占比，结果显示有89.92%的灾害点分布于较高危险区，且灾害点密度随危险等级升高而增加，表明模型具有较高的可信度。

(2)引入1985—2015年的年均降雨及一年内Tmin>0天数两个动态因子进行动态危险性评价。结果表明危险区分布范围随降雨及气温变化而出现明显差异，崩滑灾害的危险性随降雨量增加、气温升高而增加。

(3)将历史气候数据进行时间序列预测分析即得到预测气候因子图层，将其代入耦合模型公式计算出预测危险性分布图。预测危险性分布图显示研究区崩滑灾害较高危险区面积占比随气候变化呈波动上升趋势，预示未来中巴经济走廊内崩滑灾害的发生将会更加频繁。