基于图像三维连通性识别方法的长江流域干旱事件特征

邓翠玲，佘敦先，张利平，张 琴，柳 鑫，王书霞

（1.武汉大学水资源与水电工程科学国家重点实验室，武汉 430072；2.海绵城市建设水系统科学湖北省重点实验室，武汉 430072）

0 引 言

干旱通常是指在一定时期内由降水不足等因素引发的区域水份供需不平衡的气候事件，具有持续时间长、范围波及广、发生频率高、产生危害大等特征[1]。全球约有45%的土地受到干旱威胁[2]，每年受旱灾影响的人口数超过其他任何一种气象灾害[3]。自20世纪70年代以来，干旱在中国造成的经济损失每年高达上千亿人民币[1]。干旱不仅会给农业生产带来巨大损失，影响区域经济可持续发展，还会造成区域水资源短缺、生态恶化、荒漠化加剧等诸多不利影响[4-5]。IPCC（2019）第六次评估报告《气候变化与土地》特别报告[6]指出1961—2013年间，受干旱影响的面积平均每年增加幅度高于1%。有研究表明1950年之后中国区域干旱发生频次和影响范围逐渐增加，造成的经济损失呈大幅上升趋势[7-8]。干旱的多发、频发已成为制约区域社会经济可持续发展不可忽视的因素。因此，明晰区域干旱发生发展演变规律具有重要的现实意义与价值。

科学辨识区域干旱事件是干旱演变规律分析的前提[9]。干旱事件的发生、发展和恢复不仅在时间上连续，在空间上也呈现动态扩展或收缩性[10-11]，具有明显的时空三维特征。不同干旱事件的发生、演变及时空结构差异显著[11]。早期关于区域干旱的辨识研究大多对干旱事件进行降维处理，侧重关注干旱在时间维度或空间维度的变化规律[5,7,12-13]，提出和发展了一系列干旱特征分析方法，如游程理论[14]、小波分析[11]等可以从时间维度识别干旱和提取干旱特征；主成分分析方法[15]、经验正交分解[16]、聚类分析[17-18]等方法被广泛应用于分析干旱在空间上的演变规律。这类分析方法不能有效、完整地呈现干旱事件的三维结构及其在时空发展上的关联性[9,19]。当前，综合考虑区域干旱事件在时空维度上的连续性以及对干旱动态规律分析是区域干旱研究的重要方向[9-10,20]。

近年来，越来越多的研究开始从时空三维视角对区域干旱进行提取和分析。其中，具有代表性的研究是将聚类算法引入干旱事件识别[21]，该方法基于干旱指标，首先在二维空间上识别干旱“斑块”，然后在时间维度上对空间重合面积超过一定阈值的干旱“斑块”进行合并，最终得到时空连续的区域干旱事件。在应用过程中学者对此方法进行了改进和发展。如基于不同干旱指数分析了干旱事件识别结果的差异[22]；为了使识别出的干旱事件更加符合区域实际情况，需要设置合理的面积阈值，学者探讨了不同区域最小干旱“斑块”面积、相邻时段“斑块”空间重合面积和场次干旱事件干旱“斑块”最大面积差等面积阈值对干旱事件识别的影响[9,11,19]；将此方法拓展到其他类型干旱事件识别过程，基于土壤水分数据，识别和提取了农业干旱事件[23]。基于聚类算法的干旱事件识别方法对每一个干旱“斑块”均需要与其相邻时段所有干旱“斑块”进行空间关联性判断，在大面积区域进行长时间尺度干旱事件识别时耗时较长，因此发展一种效率更高的干旱事件识别方法有助于区域干旱事件演变规律研究。

一场完整的干旱事件可以概化为在三维空间上相互连通的三维连通体，与三维连通的图像有相似之处。当前图像连通性识别技术发展比较成熟，形成了计算简单高效的系列算法[24-25]，因此本文将其引入干旱事件的识别过程，该方法基于干旱指标先识别出三维连通体，然后综合考虑区域干旱事件的发生和发展过程等因素，对连通体进行拆分，最后实现对区域干旱事件的识别和提取。与基于聚类算法的干旱事件识别方法相比，该方法只需判断干旱“斑块”与所属连通体内相邻时段所有干旱“斑块”是否属于同一场干旱事件，减少了干旱“斑块”空间关联判断次数。本文选择长江流域作为研究区域，应用提出的方法识别并提取了长江流域 1960－2015年的干旱事件，基于干旱强度、干旱面积、干旱历时，干旱逐月集中点等干旱事件特征变量描述和分析了长江流域干旱特征及其发展演变特征。

1 研究区和数据

1.1 研究区概况

长江发源于青藏高原的唐古拉山脉，向东流经中国11个省级行政区后注入东海，长江干流全长6 387 km，是世界上第三长的河流。长江流域位于 90°33'E～122°25'E，24°30'N～35°45'N（图1），总面积达180万km2，约占中国国土面积的18.8%，流域内人口占中国总人口的1/3[26]。其中，宜昌站以上、从宜昌到湖口、湖口至入海口分别为长江流域上游、中游和下游，流域控制面积分别约为100万、68万和12万km2。长江流域气候温暖，雨量较为丰沛，多年平均降水量为1 127 mm[27]，受水汽输送和地形因素的影响，降水量时空分配不均，干旱在各个季节均有发生，且影响面积较大[28]。有研究表明1961—2015年期间长江流域年均干旱发生次数为3.77次[29]，给区域内水资源管理、农业生产及社会经济可持续发展带来巨大压力[29-30]。因此，开展区域干旱事件的识别并研究其发展迁移规律对长江流域防旱抗旱和水资源优化配置等均具有重要现实意义。

1.2 数 据

本文选取长江流域及周边206个气象站点（图1）的逐日气象资料（降水、日最高气温、日最低气温、日平均气温、日照时数、风速、相对湿度、大气压）开展长江流域的气象干旱分析，数据来源于国家气象科学数据中心（http://data.cma.cn/）。潜在蒸散发由联合国粮农组织推荐的Penman-Monteith公式[31]计算得到。利用反距离加权（Inversed Distance Weighted, IDW）插值方法将计算得到的气象站点月降水和潜在蒸散发月数据插值成0.25°×0.25°的空间格点数据，最后根据SPEI指数计算方法[32]得到长江流域各格点1960—2015年SPEI3指数序列。

2 研究方法

2.1 干旱指数的选取

本文采用SPEI指数刻画气象干旱，该指标考虑了降水和潜在蒸散发对干旱的影响，具有计算简单、多时间尺度等优点[32]，在气象干旱评估中应用广泛。干旱具有较为明显的季节性特征，研究表明，3个月尺度的干旱指数的既不会因为时间太短造成干旱指数变化激烈，又不会因为时间太长而均化干旱指数序列变化过程[10]，被国内外学者广泛应用到不同地区的气象干旱事件识别[9-10,20,22]。因此，本文选取3个月尺度的SPEI（记为SPEI3）指数开展后续分析。参照SPEI干旱等级划分标准[32]和长江流域干旱有关研究[28,33]，当SPEI3<−1时，认为发生干旱；当 SPEI3≥−1时，认为未发生干旱。

2.2 基于图像三维连通性识别技术的干旱事件提取方法

干旱事件具有明显时空连续性和动态演变三维特征，发生干旱的格点与图像中的有效像素格点类似，从空间中识别三维图像连通体与从干旱指标中识别三维气象干旱事件具有相似性，因此可以将三维图像连通性识别技术引入到气象干旱事件的识别。图像三维连通性识别技术基于图像有效像素格点之间的距离判断图像之间是否相互连通，若距离为1则表示相互连通，大于1则不连通[34]。棋盘距离常被用来描述图像像素格点之间的距离，它是指两点间横纵（竖）坐标之差的最大值。二维D(a2,b2)和三维D(a3,b3)棋盘距离可表示为

式中a2(x1,y1) ,b2(x2,y2)和a3(x1,y1,z1) ,b3(x2,y2,z2)分别代表二维和三维图像上的两个像素格点。因此，在二维空间上，图像中心像素格点与其周围 8个像素格点相互连通；在三维空间上，图像中心像素格点与其周围 26和像素格点相互连通（图2）。

基于图像三维连通性识别技术的干旱事件识别方法首先在干旱指标空间体中识别出三维连通体，然后根据干旱发生和重合面积阈值将连通体逐层拆分，最终得到干旱事件。并从历史文献、相关记录和土壤水分数据等不同角度对识别出的干旱事件的合理性进行验证。研究表明干旱发生和重合面积阈值与研究区面积相关，为合理识别区域干旱事件，通常将其设为研究区面积的百分比。干旱事件识别的具体步骤如下：

1）干旱指标二值化。根据 SPEI指标的干旱等级划分标准，当 SPEI<−1时视为发生干旱，将该格点视为有效像素格点；否则该格点未发生干旱，视为无效格点。

2）空间三维连续体识别。基于二值化的干旱指标数据和三维图像连通性识别原理，识别并标记三维干旱连通体。

3）三维干旱连通体拆分。在实际研究中，小面积的干旱并不会引发严重的经济损失和不良后果，也不会被有关部门记载。在干旱事件提取过程中，为避免小面积干旱事件对分析区域干旱整体发生规律和演变特征的影响，当干旱“斑块”面积占比（干旱“斑块”面积占研究区面积比值）大于干旱发生面积阈值（A，%）时认为干旱事件发生。对第2）步中识别的持续时间超过一个月的三维干旱连通体，需要设置相邻时段空间重合面积阈值占比（B，%）判断相邻时段干旱“斑块”是否属于同一场干旱事件。三维干旱连通体拆分的具体步骤如下：

1）当三维干旱连通体持续时间为1个月时，干旱“斑块”面积占比大于A时认为此“斑块”有效并保留干旱标记，否则删除标记；

2）当三维干旱连通体持续时间大于1个月时，若相邻月份（如第i月和i+1月）干旱“斑块”空间重合面积占比大于B，则属于同一场干旱事件，保留原有干旱事件标记；否则，判断第i+1月干旱“斑块”面积占比是否大于A，若大于，则开始新的干旱事件标记，若小于等于，则删除该月干旱标记，直至拆分并标记完该三维干旱连通体，开始下一个三维干旱连通体的拆分和干旱事件提取。图3展示了三维干旱连通体拆分和干旱事件提取示例流程，若假设3、4、5月时，干旱“斑块”空间重合面积占比小于B，且3、4月时单个干旱“斑块”面积占比小于A，则该干旱连通体被拆分成两场干旱事件。

2.3 干旱事件特征变量

干旱事件特征变量可以用来分析干旱事件的特征及其发展规律，干旱强度、干旱历时和干旱面积常用于分析干旱事件的整体特征；干旱逐月集中点常用于分析干旱中心随时间的移动特征，可以表征干旱事件的中心的演变迁移路径[10]；干旱累积发生次数和年累积干旱发生频率能反应一定时间内干旱发生的频繁程度，其中，干旱事件累积发生次数为特定时间内干旱事件发生的次数之和，年累积干旱发生频率为每年发生干旱的月数除以12个月。因此本文选用上述6种干旱事件特征变量描述和分析干旱事件发展演变规律。上述有关特征变量的具体定义如下：

1）干旱历时（D，月）：干旱的持续时间，数值上等于干旱发生时间与干旱结束时间之间的时间间隔。

2）干旱面积（W，km2）：干旱所影响区域在平面上的平均投影面积，数值上等于干旱持续时间内逐月干旱面积之和除以干旱历时。

式中i表示一场干旱事件的第i个月，Wi表示一场干旱事件第i个月所影响的面积。

3）干旱强度（S）：干旱事件的干旱强度即事件的平均强度，等于逐月SPEI3绝对值在时空上的均值。

式中j为特定月份干旱事件的第j个格点，Sj为第j个格点的干旱强度，等于该点SPEI3指标的绝对值，J表示一场干旱事件在第i月发生干旱的格点数。

4）干旱逐月集中点[10]：一场干旱事件第i月的干旱强度加权质心（X,Y），假设该月干旱强度集中于该点[10]，干旱事件干旱逐月集中点的移动过程能反应干旱事件发展过程中干旱中心的移动方向。

式中Xj、Yj为某一月份干旱第j个格点的经度和纬度。

3 结果与分析

3.1 基于图像三维连通性识别方法的干旱事件识别

面积阈值的选取对干旱事件识别影响较大，不合理的面积阈值往往会提取出不符合实际情况的长历时干旱事件[35]。研究表明面积阈值设置为研究区面积的1.6%能有效识别研究区内的干旱事件[10,20,23,36]，考虑到不同研究区干旱发生演变特征不同，本文将面积阈值A、B设定为同一值，分析了不同面积阈值（从1.0%至2.0%）对长江流域干旱事件识别的影响（图4）。可以看出，当干旱面积阈值取值为研究区面积的1.5%时，在长江流域识别出的最长干旱事件历时和对应事件的干旱面积逐渐趋于平稳。考虑到长江流域属于亚热带季风气候，根据历史干旱事件发生记录，该区域很少发生20个月以上的长历时干旱事件[37]，因此，本文将面积阈值设置为研究区面积的1.5%（约为2.7万km2）。

基于1960－2015年长江流域SPEI3数据，采用基于图像三维连通性识别技术的干旱事件识别方法，共识别出281场干旱事件，表1列出了干旱面积排名前10的干旱事件的发生时间、结束时间和影响范围与相关文献和记录的对比情况。总体来说，这10场干旱事件与有关研究文献和官方报告记载相符合，每场干旱事件的开始时间、结束时间和影响范围与记载的一致性较高，并且累积干旱强度较大的区域与文献记载的干旱峰值出现区域重合度较高。个别干旱事件与相关文献记载存在一定的差异，如第4、5场干旱事件的影响范围大于相关文献记载，其原因可能是上述干旱事件在空间存在连续性，但由于某些区域干旱强度较小或因干旱对社会和自然环境的影响较小，没有被记录或者被记录成不同场次干旱事件。综合来看，该方法能有效地识别长江流域干旱事件。

表1 基于图像三维连通性识别技术识别出的干旱面积排名前10的干旱事件与历史干旱事件记录的对比Table 1 Comparison of the top ten largest drought events ranked by area between records and identified based on three-dimensional graph connectedness recognition method

为进一步验证本文提出的干旱事件识别方法的合理性，选取了长江流域2006年和2009年发生的两场典型干旱事件，分析了干旱发展演变过程与相关文献和记录的吻合情况，以及基于聚类算法识别出的典型干旱事件的差异（图5）。据2006年全国旱灾及抗旱减灾情况[38]，2006年干旱事件发生在长江中上游地区，为夏秋连旱，该场干旱始于2006年6月，汉江流域和长江流域西南地区等地最先出现干旱，6月下旬干旱有所缓解，7－8月下旬旱情发展到中国西南以及江南部分地区，9月旱情得到缓解，10月旱情再次加重，11月干旱逐渐消失。本文识别出的干旱事件表明，此次干旱最早于 6月份发生在长江流域西南东部地区，逐渐向西扩张，8－9月份干旱最为严重，云南、贵州和四川等地干旱强度较大，10月份以后，干旱情况得到缓和，11月份结束。基于聚类算法识别出的干旱事件发生时间相较于历史记录和本文识别的干旱事件开始时间略早，其影响范围小于本文识别的干旱事件（图5a）。2009－2010年的干旱事件是发生在中国西南地区的秋冬春连旱，据相关记载[39]，2009年入秋以后干旱开始发展，2010年3－4月初旱情发展到高峰，西南地区多数水利工程无水可用，5月份进入雨季后，干旱逐渐解除。识别出的干旱事件表明干旱始于2009年8月，持续时间长达10个月，影响范围从贵州等地逐步发展到西南五省部分区域，次年3月干旱达到高峰，云南地区遭受严重的干旱，5月以后西南地区干旱开始逐步解除。本文提出的干旱事件识别方法和基于聚类算法的干旱事件识别方法识别出的干旱事件的开始时间、结束时间和影响范围一致性较高，且干旱较严重的区域重合度较高。但本文提出的干旱事件识别方法识别出的干旱事件影响范围比基于聚类算法识别方法识别出的干旱事件稍大，可能是由于两种方法在计算相邻时段干旱“斑块”空间重合面积时所使用的方法不同造成的。根据与历史典型干旱事件相关记录和不同方法的识别结果可以看出，本文发展的干旱事件识别方法能较好的识别出干旱事件的发生、结束、重灾区以及整场干旱事件的发展演变过程。

气象干旱发生后在短时间内得不到缓解就影响土壤水分补给，甚至造成土壤水分亏缺，因此可以分析干旱事件影响范围内的土壤水分情况来验证识别出的干旱事件是否合理。考虑到土壤水分对气象干旱具有一定滞后性，本文对比分析了历时大于 1个月的干旱事件的干旱强度与对应事件影响范围内表层（0～10 cm）土壤水分亏缺程度之间的相关关系（图6）。图中土壤水分亏缺为标准化土壤水分的相反数，其值大于 0表示土壤水分不足，越大表明土壤水分亏缺越严重。识别出的干旱事件的干旱强度与其影响范围内的表层土壤水分亏缺具有显著正相关关系，干旱强度越大的干旱事件对应的土壤水亏缺越严重；长历时干旱事件的干旱强度与土壤水分亏缺程度的相关性更高。不同角度的验证结果表明，基于图像三维连通性识别技术的干旱事件识别方法能有效识别区域干旱事件。

3.2 长江流域干旱特征变化分析

图7展示了长江流域1960－2015年长历时干旱事件的时空发展过程。可以看出，时间上，强度高、影响面积广的干旱事件大多发生在1960－1962年、1966－1969年、1978－1979年、1992－1994年、1996－1999年、2002－2004年、2006－2007年、2009－2010年、2012－2014年等时间段；空间上，长江流域中部地区长历时干旱发生相对较为频繁。1980－1990年长历时干旱事件平均干旱面积相对较小，2000年之后长历时干旱事件发生次数明显增多，1960－1999年40 a间共发生43场长历时干旱事件，平均每年1.1次，平均干旱持续时间和面积分别为8.1个月和2.74×104km2；2000—2015年16 a间发生长历时干旱事件21场，平均每年1.5场，平均干旱持续时间和面积分别为7.6个月和3.50×104km2。

刘君龙等[28]的研究结果表明2000年之后长江流域干旱频率是1960—2018年的1.5倍，干旱影响面积逐渐增加；Zhang等[40]基于GRACE重力卫星的水储量数据的研究表明，长江流域2003年、2006－2008年、2008－2010年等年份发生了长历时干旱，与本文的研究结论一致。

长江流域干旱累积发生次数、累积干旱强度的空间分布和变化趋势空间分布如图8所示，长江流域西部的四川省，中部的陕西省、湖北省、湖南省和东部的江西省干旱累积发生次数较多，上述省份部分区域1960－2015年发生了45场以上干旱事件；流域东南部的云南省干旱事件发生次数相对较少。累积干旱强度也具有明显的空间异质性，四川省、重庆市、贵州省、湖南省和湖北省等地累积干旱强度较强，而长江流域源头区域、云南省和东南部地区干旱累积强度较小。对比图8a和8b，陕西省和江西省干旱事件发生次数较多，但累积干旱强度较小，表明这些区域轻旱发生次数较多，而在四川省和贵州省部分区域，干旱事件次数较少但累积干旱强度较大，表明这些区域重旱发生次数较多，这与历史统计资料分析结论上游发生重旱的频率明显高于中下游一致[37]。

中部区域和云南等地干旱发生频率呈现出增加趋势，其中四川盆地和云南省部分区域增加较大；源头区域、四川省西部和流域东部部分地区干旱事件发生频率呈下降趋势，其中四川省西部部分地区下降最明显（图8c）。与干旱事件发生频率趋势空间分布相似（图8d），年累积干旱强度变化在流域中部和云南省部分区域呈现出增加趋势，增加较大区域集中在四川盆地和云南省等地；流域西北部和东部大部分地区年累积干旱强度呈下降趋势，且下降趋势明显的区域多集中在流域西北部。

长江流域干旱逐月集中点迁移路径较为复杂，长历时干旱多数存在先向一个方向迁移后又折回的现象（图9）。不同年代干旱事件的主要迁移方向不同，1960－1969年干旱逐月集中点主要东南方向迁移，1970－1979年和20世纪初期主要向东南方向迁移，1970－1999年主要向西北方向迁移，2000年及之后主要向西北和东南方向迁移（表2）。整体上，长江流域长历时干旱事件逐月集中点主要向西北和东南方向迁移，约占长历时干旱事件总数的 83%，向东北和西南方向迁移的干旱事件较少，约占长历时干旱事件总数的 17%干旱事件迁移方向可能与区域降水和水汽输送等因素有关[41-42]。

表2 长江流域干旱逐月集中点向不同方向迁移的长历时干旱事件数量统计表Table 2 Statistical table of the number of long duration drought events which drought center migrate to different direction in the Yangtze River Basin

4 结 论

本文将图像三维连通性识别方法引入到气象干旱事件识别过程中，不同于基于聚类算法的干旱事件识别方法先识别干旱斑块再合并，该方法先识别三维连通体，再根据阈值条件对连通体进行拆分，从而得到干旱事件。应用此方法识别了 1960－2015年长江流域气象干旱事件，并分析了干旱事件发展演变特征。主要结论如下：

1）将本文提出的基于图像三维连通性识别方法的干旱事件识别方法所识别出的干旱事件与历史记录、干旱事件影响范围内的土壤水分亏缺情况和基于聚类算法识别出的干旱事件进行综合对比，结果表明本文提出的干旱事件识别方法能有效识别区域气象干旱事件。基于SPEI3数据，长江流域1960－2015年共识别出281场气象干旱事件，其中长历时干旱事件64场，；

2）长江流域中部区域长历时干旱事件发生次数较多，1980－1990年长历时干旱事件发生频率较低，2000年之后长历时干旱事件发生频率和影响面积较 1960－1990年增加，约83%的长历时干旱事件向西北和东南方向发展；

3）长江流域中部干旱事件发生频率和年累积干旱强度大于源头区域和流域东部，且存在明显增加趋势。

本文发展的三维干旱事件识别方法为基于事件角度分析干旱提供了新途径，对长江流域干旱事件演变规律分析有助于更好理解长江流域历史干旱事件发展变化特征。但本文只分析长江流域干旱事件整体的变化和迁移特征，下一步研究可以考虑分析影响干旱发生和干旱事件逐月集中点迁移的因素，以更好地从影响干旱发生机理角度分析长江流域干旱发生演变规律。