黄 森,崔素丽,辛 鹏,王 涛,梁昌玉
(1.大陆动力学国家重点实验室 西北大学地质学系,陕西 西安 710069; 2.中国地质科学院 地质力学研究所,北京 100081)
降雨滑坡是我国黄土高原地区发育最多的滑坡类型[1-4]。2013年7月,受极端降雨影响,延安及天水地区发生了大规模地质灾害,延安地区引发了8 000余处地质灾害,严重威胁当地群众的人身财产安全[5-6];天水地区的群发性灾害造成了24人死亡,1人失踪,9 052间房屋出现了不同程度的损坏,直接经济损失达到了82.75亿元[7-10]。因此,本次研究希望通过研究雨量变化过程与滑坡分布之间的关系,可以为减灾防灾提供良好的基础支撑作用。
诱发群发灾害的雨量变化过程可分为持续累积降雨、短期强降雨、小时降雨等。针对2013年7月天水极端降雨诱发大规模群发灾害事件,于国强等[8]分析其成因时,利用四次强降雨数据,分析研究得知此次灾害事件具有群发性、普遍性局地暴发性特征;王敏龙等[9]对四次强降雨与滑坡分布进行分析后认为,此次群发灾害的形成与降雨的周期性有一定关系,尤其是小型滑坡的发生与降水周期密切相关;郭富赟等[10]在累积降雨与四次强降雨的基础上,经研究发现极端降水对灾害形成的主要促进作用表现在前期降水丰富,持续时间长,当期降雨强度大;刘林通等[11]以秦州区教场坝沟为例,研究分析小时降雨与滑坡易发性关系,结果表明此次灾害事件的实际降雨强度与滑坡发生情况基本吻合。
对于天水地区2013年7月降雨型群发灾害事件,并没有学者开展有关极端降雨雨量变化过程与滑坡分布相关性研究,前人主要分析了此次极端降雨诱发群发性地质灾害的特征及成因,在成因分析的基础上利用强降雨及累积降雨数据分析了滑坡的分布特征。笔者在选取典型受灾区域的基础上,利用获取的小时精度降雨数据,对此次极端降雨事件降雨过程与群发灾害空间分布进行了相关性分析,其中小时降雨量最大为48.4mm,累积降雨最高为665.1mm。通过对不同降雨模式与滑坡分布进行相关性分析,选取相关性最优的时间段并统计分析,得出该地区降雨诱发群发性滑坡的临界降雨量,并对此次灾害的空间分布规律进行相应分析,为该地区突发强降雨或持续强降雨情况下的地质灾害监测预警工作提供相应的支持。
2013年7月25日,甘肃省天水市出现了大范围的强降雨天气。在极端降雨的影响下多地区出现大规模地质灾害,此次地质灾害分布广、数量多,且多以狭长、小型浅层黄土滑坡及泥流型黄土滑坡为主。灾害前后对比如图1所示,图示范围为研究区内柳林村,影像源自灾害前后高精度遥感影像图,滑坡泥流相互影响、
图1 灾害前后对比图Fig.1 The contrast before and after the disaster
相互促进,整体上表现出灾害的链式效应[12-13]。众多受灾区域以天水秦州区南部娘娘坝地区最为严重,因此选取该地区作为研究区域。
研究区地处西秦岭北缘,甘肃省天水市秦州区南部,地质环境复杂,植被茂密,属长江水系白家河流域,地理坐标为34.12°~34.41°N,105.62°~106°E,海拔为1366.1~2244.8 m,处于暖温带半湿润区,降水多,温差显著,年均降水量在600 mm以上,降雨多集中在6~9月份。区内地貌以侵蚀构造中低山为主,西北部红土-黄土丘陵地貌占比较少,中低山区山陡沟深,沟谷形态多呈深“V”型,出露地层岩性以泥盆系大草滩群紫红色砂岩及板岩互层夹砾岩、泥盆系舒家坝组变质砂岩和第四系风成黄土为主,上覆黄土层较薄,平均厚度约1 m,植被多以高大乔木为主;红土-黄土丘陵区地层以新近系红色泥岩及第四系黄土层为主,厚度普遍较大,植被覆盖较差。
在收集研究区受灾前后高精度遥感影像的基础上,通过前后对比进行滑坡的解译工作,以此确保解译工作的准确性。解译工作历时5个月,笔者完成的解译工作约占总数的90%,剩余10%工作在笔者野外实地调查期间由本项目其余小组成员进行解译。结果显示,在娘娘坝及周边3个乡镇约747.61 km2的区域内,共解译出滑坡45 446处,滑坡面积达到了16.89 km2,占研究区面积的2.26%。本研究选取的降雨数据源自娘娘坝地区雨量监测站,研究区内部及周边监测站点共22座,解译结果及站点位置如图2所示,研究区东南部为基岩山区,故监测站点主要集中在中西部及西北部,每个站点的雨量数据精确到了小时。
图2 滑坡解译及站点分布图Fig.2 The map of landslide interpretation and site distribution
对降雨数据进行整理可知,此次降雨过程从6月19日到7月26日共持续38天,累积降雨量整体表现为75~650 mm,钱家坝地区累积降雨量更是达到了665.1 mm。在此期间共经历了四次短期强降雨,四次强降雨降雨量分别为30~285 mm、0~150 mm、0~45 mm、0~180 mm,数据显示,降雨多集中在前两次及最后一次强降雨,第三次强降雨雨量相对较少,不同强降雨及累积降雨降雨量分布图如图3所示。
图3 降雨等值线分布图Fig.3 Contour map of rainfall
为了研究降雨区间与滑坡分布之间的对应关系,首先需要对降雨数据进行分区,即降雨量等值线图的制作,本文选取的制图方式是克里金插值法,插值结果如前文图3所示,方法如下:
克里金插值法(Kriging),又称空间局部估计或空间局部插值法[14],是依据协方差函数对随机过程/随机场进行空间建模和预测(插值)的回归算法,该方法是在变异函数理论和结构分析的基础上,在有限区域内对区域化变量进行无偏最优估计的一种办法[15-17]。
克里金法包括普通克里金法、泛克里金法和协同克里金法等方法,其中普通克里金插值法假设条件少、需求参数较为简单,是比较常用的一种统计方式,公式为:
(1)
式中:Z*(K0)为待计算位置的降水预测值;λi为插值过程中不同雨量站点对估算位置降水量的权重;Ki表示雨量站点的实际位置;Z(Ki)表示雨量站点的实测值。
考虑到山区降雨与地形海拔有一定的因素,因此文章采用协同克里金法(Co-Kriging)进行插值分析,该方法为克里金插值法的改进方法,在建模过程中不仅需要主变量,还需要引入协变量进行联合分析,协变量数量不限,但主变量与协变量必须具有一定的相关性。本文以雨量监测站点的实际降雨数据为主变量,以研究区高精度DEM高程数据为协变量进行降雨插值,两个变量的协同克里金插值公式如下[18-19]:
(2)
式中:Z*(K0)为待计算位置的降水预测值;λi为参与插值的站点对估算位置降水量的权重;Z(Ki)表示雨量站点的实测值;λj为参与插值的DEM高程数据对估算位置降水量的权重;Z(Kj)为DEM高程数据的实测值。
在降雨量等值线图的基础上,通过统计不同降雨区间的滑坡数可初步获取降雨数据与滑坡分布的对应关系,统计过程借助GIS中以表格显示分区统计法,统计结果可直接为后续的相关性分析计算提供数据支撑。
要实现两者之间的统计分析前提有两点:①需要对降雨量等值线图进行重分类;②按多边形提取解译滑坡的中心点,利用不同时段的降雨数据提取每一个滑坡点的降雨量。满足两个前提条件后对不同的降雨数据进行滑坡数量的分区统计,得到不同降雨区间的滑坡总数,选取的降雨数据包括四次强降雨、累积降雨以及强降雨期间具有代表性的小时降雨数据。
相关性分析是反映降雨数据与滑坡分布两者之间相关程度的高低,分为有影响和无影响两种情况。由于降雨数据分区后为有序变量,且两者关系不符合正态分布,因此选用斯皮尔曼(Spearman)相关系数进行计算分析,计算方法如下:
斯皮尔曼相关系数常用字母ρ表示,它是衡量两个变量的依赖性的非参数指标,该方法相比积差相关系数而言,不需要较为严格的数据条件,只需要两个成对的变量观测值即可[20]。具体公式如下:
(3)
在求得斯皮尔曼相关系数后,为了检测系数的可靠程度,需要对所得系数进行相应的检验,该检验方法称为t检验,公式如下:
(4)
式中:r代表相关系数,n代表所求样本的数据总量,则自由度为n-2。
若所求t值大于0.05,说明两个变量不存在相关性,若所求t值小于0.05,则表示两组变量存在显著相关性。相关性强弱取决于斯皮尔曼相关系数的大小,ρ的取值介于-1到1之间,系数为负表示负相关,系数为正表示正相关,
将降雨数据插值结果与滑坡分布图进行叠加,在相关性判别的基础上,剔除无相关性的降雨数据,在有显著相关性的降雨数据中选择相关程度最大且呈正相关的数据,利用该降雨数据与滑坡分布的叠加图进行降雨阈值的统计分析,最后利用累积降雨分析降雨与滑坡分布的空间分布规律。
不同降雨数据与滑坡分布的相关性及t值计算结果如表1、表2和表3所示。
由表1可知,10次有代表性质的小时降雨中有3次与滑坡分布具有显著相关性,其中两次计算结果小于0,表明两者之间呈负相关,仅6月20日2-3点相关性计算结果为0.636>0,故选取此次小时降雨过程分析诱发群发性浅层滑坡的小时降雨阈值。
表1 小时降雨与滑坡分布相关性分析结果Table 1 Results of correlation analysis between hourly rainfall and landslide distribution
由表2可知,四次强降雨过程中有两次与滑坡分布有显著相关性,两次计算结果分别为0.683>0和-0.643<0,故选取前者分析诱发群发性浅层滑坡的强降雨阈值,即第二次强降雨过程。
表2 强降雨与滑坡分布相关性分析结果Table 2 Results of correlation analysis between heavy rainfall and landslide distribution
由表3可知,总的累积降雨与滑坡分布有显著相关性,因此通过对累积降雨与滑坡分布进行统计分析可以得到该地区诱发大规模浅层滑坡的累积降雨阈值。
表3 累积降雨与滑坡分布相关性分析结果Table 3 Results of correlation analysis between cumulative rainfall and landslide distribution
通过相关系数的计算、筛选,不同降雨模式下最优降雨数据与滑坡分布叠加图及数据统计如图4、图5所示。
图4 降雨数据与滑坡分布叠加图Fig.4 Overlay of rainfall data and landslide distribution
图5 降雨数据与滑坡分布统计图Fig.5 Statistics of rainfall data and landslide distribution
图4(a)、图5(a)为最优小时降雨与滑坡分布对应情况,具体时间为2013年6月20日2~3点,从图中可以看出,群发滑坡多集中分布在降雨量15~20 mm之间,约占总数的48.86%,因此可以认为诱发群发性浅层滑坡的小时雨强最小为15 mm/h,最大为20 mm/h,小时降雨量小于15 mm或大于20 mm诱发群发滑坡的概率相对较小。
图4(b)、图5(b)为第二次强降雨阶段降雨量与滑坡分布对应情况,图中显示滑坡占比最多的区域集中在降雨量105 mm~120 mm之间,约为26.36%,在此之前,滑坡分布随降雨量增加整体呈上升趋势,高于该降雨区间后呈下降趋势,因此可以认为诱发群发性浅层滑坡的强降雨雨量最优值处于该降雨区间内。
从小时降雨及强降雨数据统计可以看出,在短时间内,随着降雨量的增加,滑坡数量呈现出先递增后减小的趋势,造成这种现象的原因在于,研究区内多为中低山地貌,山体坡度普遍较大,因此在雨水落至坡体后没有充足时间产生下渗作用就会在坡体上形成地表径流,特别是在暴雨期间,雨量大,来势猛,地表径流的形成时间会比雨量较小时提前,地表径流的提前形成降低了雨水在浅层土体中的入渗量,因此在降雨量超过一定界限时,滑坡的形成概率反而会减小。
此次持续降雨事件共持续了38天,在此期间部分地区累积降雨已超过了往年年平均累积降雨量,对累积降雨数据与滑坡分布进行统计分析得到结果如图6所示。
图6 累积降雨与滑坡分布统计图Fig.6 Statistics of cumulative rainfall and landslide distribution
从图中可以看到,最高点横坐标为275 mm~300 mm,此时滑坡数量为6179,约占总数的13.63%,除该点外,滑坡发育数量在累积降雨量250 mm~625 mm之间基本保持稳定,共占据总数的85.61%,经对比降雨等值线图可知图中局部突出是由于有较大的降雨量区域,综合考虑降雨时长及黄土层厚度可知,累积降雨量达到250 mm时研究区内上覆黄土层含水率已经饱和,此时黄土层与下伏基岩之间摩擦力降至最低,随着累积降雨的持续增加,降雨量已不是浅层黄土滑坡发生的主要影响因素,因此可认为,诱发群发性浅层滑坡的临界累积降雨量为250 mm。
相关性计算结果显示,滑坡分布与累积降雨有着显著相关性,两者叠加如图7所示。
图7 累积降雨与滑坡分布关系图Fig.7 Relationship between cumulative rainfall and landslide distribution
图中颜色从红色到蓝色表示累积降雨量的增长,通过与滑坡点进行叠加后发现,淡黄色到蓝色区域滑坡分布最多,以中部蓝色累积降雨量多的区域分布最为集中;红色、黄色区域累积降雨量少,对应的滑坡数量少,尤其是红色区域基本没有滑坡灾害的出现。
从图中可以反映出滑坡分布与降雨强度成正相关,降雨强度大的地区群发滑坡特征非常明显,研究区降雨强度从中部向两侧递减,最多降雨量集中在娘娘坝镇区域,方向在NNE向,从中部向东南向及西北向雨量逐渐减少,相对于东南向,西北向降雨递减速率大,从地形考虑主要是因为研究区东南部为基岩山区,西北部为黄土丘陵地带,群发滑坡的递减规律同降雨变化趋势保持一致,西北部强度衰减速率大于东南部。因此研究区滑坡空间分布整体表现为中部集中,东南部偏多,西北部较少,滑坡空间分布与降雨落区高度一致,降雨集中区域滑坡多,降雨少的区域滑坡少。
解译滑坡数据显示,此次群发灾害分布最密集区域为研究区中部位置,整体呈NNE向斜向分布,西北部分布少,东南部分布较多,联系滑坡解译图与四次强降雨降雨量等值线分布图,易知两者之间的对应关系,在此基础上,可对四次强降雨如何影响滑坡空间分布进行分析讨论,并明确每次强降雨在此次灾害事件中所承担的作用。
第一次强降雨的降雨多集中在研究区东部靠北区域,滑坡数据显示,滑坡多集中在研究区NNE向,首次强降雨雨量集中区域与滑坡集中分布区域不一致,从相关性分析结果也可确定,首次强降雨与滑坡分布无显著相关性。
第二次强降雨多集中在研究区中部娘娘坝镇地区,最高降雨量为150 mm,降雨量从大到小变化趋势整体上较为符合此次群发滑坡的分布情况,相关性分析结果显示,第二次强降雨与滑坡分布具有显著相关性,且相关性为正。
第三次强降雨降雨量多集中在研究区南部偏西娘娘坝镇与大门乡交界位置,面积相对较小,最高降雨量为45 mm,此次降雨在空间上表现出局部性,整体变化小,相关性计算结果显示,此次强降雨过程与滑坡分布之间无显著相关。
第四次强降雨降雨量最大为180 mm,位于研究区中部位置,该处为滑坡密集分布区域,但是此次降雨过程降雨量大小表现为中部>西南向>东北向,与滑坡点分布情况不符,且相关性计算结果显示为负,表示负相关。
通过对比四次强降雨与滑坡分布关系可知,第一次降雨与滑坡分布没有明显相关性,但首次强降雨降雨量是四次强降雨中最大的,因此可以认为此次降雨过程应为后期灾害形成做前期降雨的准备工作;第二次强降雨与滑坡分布相关性最好,由此可以推断,第二次强降雨应是此次群发灾害的转折点,在第一次强降雨的基础上,再一次的强降雨冲刷及雨水下渗作用,奠定了后期群发灾害的空间分布状态,但是由于降雨致灾具有一定的迟滞性,此次降雨并没有引发大规模灾害的发生;第三次强降雨降雨量是四次降雨过程最少的,对于后期灾害的形成作用应表现在对第二次强降雨结果的推动作用;最后一次强降雨过程降雨量大小基本与第二次持平,此次降雨过程应是整个群发灾害事件的最后助力,并最终形成了由第二次强降雨决定的群发滑坡空间形态。
本文以天水市2013年6、7月雨量监测站降雨数据及高清遥感影像为基础,利用协同克里金插值法与斯皮尔曼相关系数研究不同降雨模式与滑坡分布的相对关系,在此基础上通过统计分析得出不同模式下的降雨阈值及降雨与滑坡空间分布规律,结论如下:
(1)本次群发性浅层滑坡的形成主要是在前期四次强降雨的周期作用下形成的,做降雨阈值分析时选取的最优强降雨过程为第二次强降雨,四次强降雨过程对此次群发灾害事件空间分布的影响表现为首次降雨的准备工作,第二次降雨的决定作用,第三次降雨的推动作用及第四次降雨的成灾作用。
(2)通过研究小时降雨、强降雨及累积降雨与滑坡空间分布的相关性,最终认为研究区群发性浅层滑坡的小时降雨应以2013年6月20日2点—3点降雨模式为主,诱发群发性浅层滑坡的最小雨强为15 mm/h,小时降雨量在15 mm~20 mm为群发滑坡最易形成区域;强降雨过程应以2013年7月7日—7月8日降雨模式为主,诱发群发性浅层滑坡的最小降雨阈值为105 mm、降雨滑坡发育强度最大区域对应强降雨雨量为105 mm~120 mm;诱发群发性浅层滑坡的临界累积降雨量为250 mm。
(3)滑坡空间分布主要表现在与降雨强度的一致性,降雨集中区域滑坡多,降雨少的区域滑坡少,实际表征为中部多,逐级向两侧递减,东南部衰减速率小于西北部。