长江上游流域生长季气象干旱分异特征

2022-03-14杨少康董晓华刘艳丽丛方杰

水土保持研究 2022年2期

杨少康, 刘冀, 魏榕, 董晓华, 刘艳丽, 丛方杰

(1.三峡大学水利与环境学院，湖北宜昌 443002; 2.三峡库区生态环境教育部工程研究中心，湖北宜昌 443002;3.水资源安全保障湖北省协同创新中心, 武汉 430072; 4.南京水利科学研究院水文水资源与水利工程科学国家重点实验室,南京 210098; 5.水利部应对气候变化研究中心, 南京 210029; 6.烟台大学网络与教育技术中心, 山东烟台 264005)

干旱是一种缓慢的自然灾害，作为一种经常出现的极端气候事件，其影响范围可以跨域多个国家和地区[1]。及时开展区域干旱监测及特征研究有利于清晰认识干旱演变规律并科学有效开展抗旱减灾措施。因此，随着对干旱深入研究，国内外学者陆续提出不同的干旱指数来对干旱进行量化分析，以便于能够更好认识干旱变化规律，其中包括气象干旱指数、农业干旱指数、水文干旱指数及社会经济干旱指数[4]。

IPCC气候评估报告指出全球气温将不断持续升高，且这种变化在北半球更为明显[5]。在全球变暖的背景下，长江上游流域正表现出暖干化趋势，尤其是2005年之后干旱频次、程度、历时增加显著[6]。因此，许多学者对长江上游流域干旱变化的时空特征进行了研究，沈浒英等[7]基于PAI指出，长江上游地区发生大范围的严重干旱都是在特定的气候背景及大气环流异常情况下发生；洪兴骏等[8]认为对长江上游干旱指标评估干旱等级时，应该考虑不确定性的影响，以此更为准确有效地为抗旱提供决策；刘志雄等[9]通过对长江上游流域进行旱涝分级，指出对于范围比较大的区域可以通过划分为更小的区域开展研究。不过，多数研究是以长江上游流域干旱的整体变化趋势作为研究对象，且多是基于月、季、年尺度上，不能很好地反映植物生长季所需的气象要素的特征，且长江上游流域地理跨度大，地形复杂，气候条件和敏感性在不同的区域差异很大[10]，将长江上游流域作为单一的研究区域，无法揭示其流域内部的干旱变化。因此为更好了解长江上游流域干旱演变特征，有必要对长江上游流域进行区域空间细分。

基于上述原因，本文选取Vicente-Serrano等[11]构建的气象干旱指数SPEI，SPEI不仅结合气温和降水，更综合比较潜在蒸散量对干旱的影响，是目前国内研究干旱应用最为广泛的指数之一[12-14]。对长江上游流域植物生长季SPEI利用REOF进行计算，完成长江上游流域空间上的细分，进而选取Yevjevich于1967年提出的游程理论对长江上游流域各分区干旱特征进行分析[15]，以期对长江上游流域农业生产、生态调度以及水资源安全提供决策支持。

1 数据及研究方法

1.1 研究区域和数据来源

长江发源于青藏高原各拉丹东，至湖北宜昌以上为上游流域，总体跨越6个省市区，干流全长4 511 km，流域面积100万km2，介于24°46′—35°91′N，90°47′—113°31′E。长江上游流域自然地理条件分异明显，地跨我国第一、二级阶梯，流域地形复杂，不仅受西南季风和东南季风的影响，又受“中华水塔”青藏高原的影响，流域地理位置及气象站分布如图1所示。本文采用中国气象科学数据共享服务网(https:∥data.cma.cn/)提供的长江上游流域气象站1961—2019年气象数据。综合考虑站点缺测、数据起止时间及空间代表性等情况，最终选取流域内及周边共67个气象站点，对于少数站点缺测数据采用临近站进行插补。

图1 长江上游流域气象站点分布

1.2 研究方法

1.2.1 标准化降水蒸散指数标准化降水蒸散指数(SPEI)是在标准化降水指数(SPI)的基础上发展而来，综合考虑了温度、降水及蒸散作用的影响，是较为理想的气象干旱指标。同SPI一样，SPEI可以计算多尺度，可分为月尺度(SPEI-1)、季节尺度(SPEI-3)、半年尺度(SPEI-6)及年尺度(SPEI-12)[16]。由于SPEI-3更适用于流域的连续性干旱事件识别，因此本文主要根据长江上游流域植物生长季(4—10月)各站点数据计算得到站点SPEI-3，运用反距离平方法进行空间插值，从而得到长江上游流域各像元SPEI-3值，并对各像元进行干旱识别。SPEI的具体计算方法如下：首先计算逐月降水与潜在蒸散量的水汽平衡Di。

Di=Pi-PETi

(1)

式中：i为月份；P为月降水量(mm)；PET为月潜在蒸散量(mm)。其中PET采用联合国粮食及农业组织修订的Penman-Monteith模型。

(2)

式中：Rn为净辐射[MJ/(m2·d)]；G为土壤热通量[MJ/(m2·d)]；T为日均气温(℃)；u2为2 m高风速(m/s)；es为饱和水汽压(kPa)；ea为实际水汽压(kPa)；Δ为饱和水汽压—温度曲线斜率(kPa/℃)；γ为干湿表系数(kPa/℃)。

然后采用三参数log-logistic概率分布函数计算概率密度函数f(x)。

(3)

式中：α，β，γ分别为尺度参数、形状参数和位置参数，这些参数主要通过线性矩法估计得到。因此，可以进一步得到三参数log-logistic概率分布函数F(x)。

(4)

参数α，β，γ分别采用线性矩法获得。

(5)

式中：Γ为阶乘函数；ω0，ω1，ω2为原始数据Di的概率加权矩。

对概率分布函数F(x)进行标准化处理，即可得到SPEI。

(6)

(7)

式中：c0=2.515517，c1=0.802853，c2=0.010328，d1=1.432788，d2=0.189269，d3=0.001308。在本文中，选用SPEI作为定义干旱的指标。根据国家气象等级标准[17]，将干旱分为5个等级，具体分类情况见表1。

表1 SPEI值干旱等级分类标准

1.2.2 旋转正交分解经验正交分解EOF的功能是把随时间变化的向量场分解为不随时间变化的空间函数部分及依赖时间变化的时间函数部分，对包含空间点的场随时间变化进行分解，将变量场时间和空间变化特征进行分离，展开得到前几个特征向量，用尽可能少的模态表达主要的时空变化特征。针对EOF存在过于强调整体结构而忽视了重要的局部特性的问题，可使用旋转正交分解REOF，即对EOF载荷矩阵进行极大方差旋转，将特征向量高载荷集中于某一较小区域上，突出气候要素在空间分布的局域特征。REOF具体计算步骤可参考文献[18]。

1.2.3 集合经验模态分解经验模态分解EMD是一种自适应的时频局部化分析方法，该法可从离散无序的原信号序列中分解出有限本征函数IMFs分量，进而体现出原信号序列不同时间尺度的局部特征，常适用于非线性、非平稳的时间序列。集合经验模态分解EEMD在EMD的基础上通过添加人为白噪声和求取集合平均方法，有效解决了采用EMD分解时间序列时可能出现的模态混叠现象。具体计算步骤详见文献[19]。

1.2.4 Mann-Kendall检验及Sen趋势度 Mann-Kendall趋势检验常用于气象领域的时间序列趋势分析，主要是根据检验统计量Z判断序列的趋势显著性[20]，通常结合Sen趋势度使用，Sen趋势度检验是由Sen[21]提出的研究长时序变量的方法，趋势度β计算公式如下：

(8)

式中：xj，xi分别为不同时刻下的变量值；median为中位数函数。当β>0时表示变量的时间序列呈上升趋势，反之当β<0表示变量的时间序列呈下降趋势。

1.2.5 干旱特征指标目前一般基于游程理论来进行干旱特征的识别，即先选定干旱指标阈值，若某个像元的干旱指标超过阈值就被认定属于干旱[22]，本文将干旱阈值指标选为-1，即判定SPEI小于-1的像元点为干旱像元。干旱指标识别过程如图2所示，当某个像元关于SPEI的时间序列在一个时间或多个时间内连续大于-1时，出现正游程，反之则出现负游程。在干旱特征指标的研究中，负游程的长度被称为干旱历时d，干旱强度s为干旱历时与干旱阈值水平包围的面积，干旱面积则是像元SPEI值低于阈值的干旱像元面积总数。

图2 干旱特征识别过程

1.2.6 干旱离散度干旱离散度(AD)是一种反映干旱区域干旱强度聚集度的指标，其值越小，说明区域干旱发生越集中，反之则越发散。计算AD，首先需要计算各区域的干旱集中点。本文根据不同分区每年4—10月份的像元数据，得到7个逐月干旱集中点，再利用这7个干旱集中点计算出当年干旱集中点，依次方法可以得出不同分区1961—2019年59 a的逐年干旱集中点，利用不同分区59个干旱集中点，算出各分区59 a平均干旱集中点，最后算出各分区59个逐年干旱集中点与该分区59 a平均干旱集中点的相对比值，得到该分区AD，具体计算方法见公式(9)—(12)。

首先计算不同分区每年4—10月份逐月干旱集中点：

(9)

式中：si为第i个干旱像元的强度值；xi，yi为第i个干旱像元的经度、纬度坐标；xnm，ynm为第n年中第m个月干旱逐月集中点的经度、纬度坐标。

根据每年4—10月算出来的7个逐月干旱集中点，得到当年干旱集中点：

(10)

式中：sm为第m个月逐月干旱集中点强度值；xm，ym为第m个月干旱集中点的经度、纬度坐标；xn，yn为第n年逐年干旱集中点的经度、纬度坐标。

再根据逐年干旱集中点，算出平均干旱集中点：

(11)

式中：sn为第n年逐年干旱集中点强度值；xn，yn为第n年干旱集中点的经度、纬度坐标；xa，ya为第a分区平均干旱集中点的经度、纬度坐标。

最后计算各分区干旱离散度：

(12)

式中：x，y为该分区经度、纬度的长度差；AD为各分区离散度。

2 结果与分析

2.1 干旱指数空间分布特征

为探究长江上游流域不同分区的干旱演变特征，本文开展了基于EOF、REOF的长江上游流域生长季SPEI-3特征分区的对比研究。由表2可知，由于长江上游流域较大，地形复杂，导致EOF分解的各个模态的方差贡献率不高，且呈现逐渐减小趋势，前8个模态的累积方差贡献率达到67.71%，表明分解出来的模态函数可解释长江上游流域大部分地区干旱变化。将前8个模态进行极大方差旋转，旋转后得到的各个特征向量方差贡献率均匀分布。模态中的载荷值反映了空间相关程度，是气候分区的主要依据，可将载荷值高的地区视为同一区域。图3给出了长江上游流域各空间模态分布图。

表2 长江上游流域生长季SPEI前8个模态EOF/REOF方差贡献率 %

将前8个模态以空间载荷值绝对值大于0.15进行划分，可将长江上游流域划分为8个子区域，高载荷区域将长江上游流域绝大多数地区涵盖，存在个别相邻分区在低值区域少量重叠现象，但整体未影响分区结果。具体分区结果如图4所示。

图3 长江上游流域基于REOF的8个空间模态分布

图4 基于REOF的长江上游流域SPEI-3指数分区

2.2 各分区干旱变化的趋势特征

为分析长江上游流域各分区SPEI-3的周期特征，利用EEMD方法对8个分区的SPEI-3时间序列进行分解。由图5可知，4个分量IMFs及1个趋势项RES，IMFs分量依次代表SPEI-3序列分解之后从高频到低频的变化，纵轴没有实际意义，趋势项RES则用来解释SPEI-3的总体变化趋势。从RES变化可以看出，SPEI-3呈明显下降的分区为Ⅴ区，表示该区逐渐趋于干旱。Ⅱ区、Ⅳ区及Ⅷ区的趋势呈现U型，其拐点分别在1992年、1997年及1982年，这三区干旱呈现“干旱—湿润”态势，区别在于Ⅱ区及Ⅷ区变化较为剧烈，Ⅳ区变化则较为平缓。Ⅰ区、Ⅲ区及Ⅶ区呈现上升趋势，表明这3区趋于湿润，且Ⅶ区的湿润趋势比其他两区较为明显。Ⅵ区趋势呈现波浪线型，“湿润—干旱”态势交替出现，整体趋势变化不显著。

各分区IMFs分量的周期及贡献率见表3。各区IMF1分量的周期普遍为3 a左右，且各分区1～10 a周期占主导地位。比较长江上游流域不同分区周期贡献率，年代际及趋势项的贡献率远远低于年际贡献率，但各分区之间又互相有所差异。长江上游整体流域的干旱周期以年际变化为主，Ⅰ和Ⅲ区的年代际贡献率相比其他区域较高，趋势变化贡献较高的区域则主要是Ⅴ区及Ⅷ区，其他分区的年代际贡献率高于趋势变化。

为确定各分区SPEI-3的主要周期，对各分区的IMFs分量进行显著性检验，具体结果如图6所示。长江上游流域8个分区的IMF1分量除Ⅲ区及Ⅳ区落在80%的置信区间以上，其余各分区IMF1分量均接近80%的显著性检验线，但未通过检验。在关于IMF2分量的显著性检验中，仅有Ⅱ区、Ⅴ区及Ⅶ区通过显著性检验，且Ⅶ区IMF2分量落在95%置信区间以上。IMF3分量的显著性检验中，除Ⅶ区外，其余各区均未通过显著性检验。Ⅰ区、Ⅱ区、Ⅴ区及Ⅶ区的IMF4分量均落在80%的置信区间以上，均通过显著性检验，其中Ⅰ区及Ⅴ区通过了95%的显著性水平检验。由此可知，以IMF1分量周期为主要周期的为Ⅲ区、Ⅳ区、Ⅵ区及Ⅷ区，Ⅱ区及Ⅶ区主要以IMF2分量周期为主要周期，Ⅰ区及Ⅴ区则以IMF4分量周期为主要周期。整体而言，长江上游流域以IMF1及IMF2分量的年际周期为主要周期。

2.3 各分区干旱特征

由表4可知，长江上游流域各区生长季的平均干旱历时均在1.4月以上，其中Ⅷ区干旱历时最长为1.59月，最短的为Ⅶ区的1.41月；Cv值及σ值用来描述系列的离散程度，各区干旱历时的Cv值均大于0.3，σ值均大于0.5，其中Ⅶ区干旱历时的离散程度最大，Cv值达到0.71，σ值达到1，Ⅲ区的Cv值和σ值最小，说明该区干旱历时离散程度较低；Cs值主要用来分析系列不对称的程度，除Ⅰ区、Ⅳ区及Ⅵ区外，其余5区均表现为正偏；在变化趋势分析中，仅有Ⅶ区干旱历时的Z值为-2.05，β值为-0.01，说明该区通过95%置信水平检验，干旱历时呈显著下降趋势，其余各区均未通过检验。通过分析长江上游流域各区植物生长季的干旱强度可知，各区的平均干旱强度差别不大，且各区β值均为0，说明长江上游流域整体干旱强度变化不大。对干旱面积进行分析，Ⅶ区的干旱面积最大，达到30.36%，干旱面积最小的为Ⅷ区，仅占19.18%，各区的Cs值均表现为正偏，Ⅰ区和Ⅶ区的干旱面积均通过95%的置信水平检验，呈显著下降趋势。整体而言，Ⅶ区的干旱历时和干旱强度最小，但其干旱面积最大。Ⅷ区的干旱面积最小，干旱历时和干旱强度最大，且历时和强度的离散程度较低，说明该区干旱发生范围较小，但其强度大。