曹 青,刘 康 龙,郝 振 纯

(1.水利部水文气象灾害机理与预警重点实验室/协同创新中心,江苏 南京 210044; 2.南京信息工程大学 水文与水资源工程学院,江苏 南京 210044; 3.河海大学 水文水资源与水利工程科学国家重点实验室,江苏 南京 210098)

0 引 言

干旱是灾害性最大的常发自然灾害之一,对农业、水资源供给等人类活动构成了巨大威胁[1]。近年来,全球气象学者提出一系列气象干旱评价指标来研究不同地区的干旱演化特征,如标准化降水蒸散指数(SPEI)[2-3]、标准化降水指数(SPI)[4-5]、帕尔默干旱强度指数[6]等,并与农业、环境、水利水电工程等学科交叉结合以研究干旱对不同学科现象的影响特征[7-8]。李华伟[9]通过Copula函数开发出的气象-农业模型模拟得出,新疆谷物产量相应较敏感的干旱时间尺度较长,而棉花所响应尺度较短;Cuartas等[10]基于多种干旱指数评估了巴西的水资源储量情况,发现SPEI可以用于评估测站数量较少的流域干旱事件,且对巴西干旱形势下水力水资源管理提出警醒与理论依据。与此同时,研究学者尝试寻求地区干旱的年或年代际等多种时间尺度演变趋势的遥相关因素和动力机制[11-12],刘仲藜等[13]、张剑明等[14]研究发现,洞庭湖流域以及湖南大部分区域干旱事件皆与西太平洋副热带高压的形势有着显著的相关影响。更进一步地,部分学者通过构建数学模型并对其不断修正以期对未来干旱现象做出准确预警和处理[15-17]。Dikshit等[18]采用标准降水指数(SPI)和SHAP预测算法阐述了气候变量在月尺度和年尺度不同范围内的重要性,以此证明了在预测模型中添加气候变量指数作为预测因子的必要性。Malik等[19]通过有效干旱指数(EDI)研究得出GWO与SHO等元启发式算法对气象干旱预测具有更好的适用性和优越性,且混合SVR-GWO模型更优于SVR-SHO模型。

松花江流域是中国十大流域之一,作为中国东北地区重要的农业产区和重工业基地,其耕地面积大、土壤肥沃等优势使得东北地区农业发展十分发达,因此研究松花江流域的干旱特征对农业灾害防治预警有着重大意义[20-21]。赵兰兰等[22]通过VIC模型模拟土壤含水率距平百分率指数,证明了基于该指数的干旱等级划分标准具有较好的适用性。Faiz等[23]研究发现基于参考蒸散量的干旱指数能够捕捉到更多的干旱事件,且多个干旱指数在RCP8.5下皆显示出最高的干旱发生率。通过采用WEAP模型,预测气候变化条件下未来河流径流可能会减少24%,且轻、中度干旱发生概率增加[24]。Song 等[25]探究了松花江干旱事件与极端降水的演变特征,得出松花江流域春季与秋季干旱发生率较高。

综上所述,松花江流域已有的气象干旱研究主要集中在干旱指数的时空演变和多模态特征研究,以及与农业、植被、环境等其他领域的相互影响和干旱频率强度风险分析等方面,而在干旱现象年代际突变点识别以及突变前后变化差异等方面仍存在不足。在20世纪80～90年代之后,由于降水减少与全球变暖的因素,中国北方呈现“暖干化”趋势[26],松花江流域作为中国重要的商品粮基地,厘清干旱指数是否存在突变及其突变前后与遥相关因子的响应关系,有利于加深对松花江流域旱涝灾害风险趋势的理解。因而本研究结合多种非参数检验方法识别标准化降水蒸散指数时间序列突变点,并运用交叉小波变换与小波相干分析,探究突变前后标准化降水蒸散指数与遥相关因子相互作用关系的差异,具象化流域干旱突变特性的空间分布模式。研究成果可为松花江流域抗旱减灾、生态水资源管理、保障国家粮食安全等工作提供技术支撑。

1 研究区概况与数据来源

松花江是黑龙江在中国境内最大支流,为中国七大河流之一,总长度约达1 927 km,流域总面积达到55.68万km2,横跨内蒙古、吉林、黑龙江3个省份;流域多年平均气温为3～5 ℃,年内日平均温度最高月份为7月,可达25 ℃,1月最低,约为-20 ℃。多年平均降水量约为500 mm。本研究所用气温与降雨数据为中国气象数据网中的中国地面气候资料日值数据集V3.0;北极涛动Arctic Oscillation(AO)指数与太平洋年代际振荡Pacific Decadal Oscillation(PDO)数据来自NOAA网站(http:∥www.cpc.noaa.gov);El Nio-Southern Oscillation(ENSO)指数的计算采用全球海表温度再分析数据,来自英国气象局哈德利中心观测数据集(HadiSST1)。松花江流域地理信息和气象站点如图1所示,选取松花江流域52个气象站点1960～2020年共61 a逐日降水与气温数据,气象站点数量较多且在流域范围内分散均匀,且本研究选取的时间序列跨度较长,确保了研究分析与结论的可靠性。

图1 松花江流域气象站分布Fig.1 Meteorological station distribution in Songhua River basin

2 研究方法

2.1 遥相关因子

遥相关研究对于流域气象旱涝要素变化与扰动机制具有重大意义,研究显示中国东北地区干旱变化的第一二主周期与厄尔尼诺-南方涛动事件具有一致的振荡准周期[27-28]。北极涛动呈负位相分布、太平洋年代际振荡指数较高时,松花江流域西部易雨涝[29-30]。本文选取北极涛动、太平洋年代际涛动、厄尔尼诺-南方涛动作为遥相关研究因子。

2.2 标准化降水蒸散指数

SPEI是在标准化降水指数(SPI)的基础上提出的标准化蒸散指数,具备多时间尺度特征,计算简便的同时还引入了影响干旱的潜在蒸发项,能够更好地刻画中国区域的干旱事件[31]。沈国强等[32]、梁丰等[33]、徐一丹等[34]还通过干旱灾情数据、土壤水分检测值与标准化降水蒸散指数的显著关联性等多个角度论证了标准化降水蒸散指数在东北地区干旱定量化研究的适用性。

标准化蒸散指数首先采用Thornthwaite方法计算逐月潜在蒸散发量,将逐月降水量与潜在蒸散量的差值序列进行正态化处理,构建成不同时间尺度Log-logistic概率分布的累计函数,对时间序列进行拟合并标准化,得到不同时间尺度的SPEI指数序列。本文取12个月尺度的SPEI指数(SPEI-12),SPEI干旱等级分类如表1所列。

表1 SPEI-12干旱等级分类Tab.1 SPEI-12 drought classification

2.3 突变检验方法

(1) Yamamoto法。该方法是从信噪比的角度来分析突变问题,检验两子序列均值的差异是否显著,被广泛应用于识别气象要素年代际突变点[35]。

(2) Lepage法。该方法为无分布双样本的非参数检验方法,将序列中的两个子序列视为两个独立主体,若统计检验显示两个子序列有显著差异,则认定在两个子序列的划分基准线点时刻发生了突变[36]。

2.4 交叉小波变换和小波相干变换

交叉小波变换(XWT)结合了小波变换与交叉谱分析,可以识别出两个时间序列在不同时频域相互作用的显著性,并且能识别出序列间的相关一致性及在时频空间中的相位关系。设WX(s),WY(s)为两个时间序列X与Y的交叉小波变换,则定义其交叉小波功率谱为[37]

(1)

式中:左项绝对值即对应交叉功率谱密度,绝对值越大,表明两者高能量区相关越显著。

(2)

式中:σX,σY分别为时间序列X,Y的标准差;v为Morlet小波变换自由度,取2,当式中左项超过红色噪音功率谱的95%的置信限上界,则认为通过了显著性水平α=0.05的红色噪音标准谱的检验[37];Zv(P)是与概率P有关的置信度,在显著性水平α=0.05下,Z2(95%)=3.999。

小波相干变换(WTC)弥补了交叉小波在低能量区识别序列间相关关系的不足。定义序列X与Y的小波相干谱为

(3)

本文使用MATLAB中的小波工具箱(http:∥grinsted.github.io/wavelet-coherence/)来进行小波变换相关计算,将处理后的数据导入工具箱,从而分析得出交叉小波与小波相干谱图。

3 结果分析

3.1 干旱指数突变的时空分布特征

结合松花江流域52个气象站点1960～2020年共61 a逐日降水数据,计算逐年SPEI-12指数,使用Yamamoto法、Lepage法识别并检验指数年序列突变点,线性拟合突变前后变化趋势直线,并对该年序列取10 a滑动平均,结果如表2、图2所示。

表2 SPEI-12指数突变识别与检验Tab.2 Identification and test of the mutation of the SPEI-12

图2 SPEI-12与遥相关因子指数变化Fig.2 Variation of SPEI-12 and teleconnection factor index

表2表明,两种检验方法都显示SPEI-12指数序列最为显著突变点为1984年,本研究取1984年为松花江流域SPEI-12指数的突变节点。

由图2可知,SPEI-12指数突变前后趋势变化明显,突变前拟合斜率约为-0.017,突变后拟合斜率约为-0.046,突变后下降趋势增强。

利用空间插值方法处理所取突变点前后的SPEI-12指数,绘制出地理空间分布图,如图3所示。由图3(a)、(b)可知松花江流域SPEI-12突变前空间上北部与南部少部分区域数值在正常范围之内,其他大部分区域数值处在轻涝级别;突变后转变为除黑河市南部与齐齐哈尔市西北部小部分区域有轻旱现象外,其他绝大部分区域数值处在正常范围内,无明显旱涝现象。图3(c)则显示出松花江流域整体SPEI-12指数趋于减小,齐齐哈尔市以西北呼伦贝尔以东南地区、黑河市周边地区和吉林北部地区变化最为显著。

图3 SPEI-12突变前后空间分布Fig.3 Spatial distribution of SPEI-12 before and after mutation

3.2 突变前后干旱指数遥相关因子分析

3.2.1交叉小波与小波相干谱分析

本文对松花江流域SPEI-12指数突变前后序列与3个遥相关因子指数对应年际变化序列进行交叉小波分析及小波相干分析,探究松花江流域降水集中指数突变前后与遥相关因子的共振周期、相位关系与显著时段等遥相关特征的区别与变化,绘制突变前后交叉小波能量谱与小波相干谱图,如图4、图5所示。

注:(a)(c)(e)为突变前,(b)(d)(f)为突变后。图4 PDO、AO、ENSO与SPEI-12交叉小波功率谱Fig.4 Cross wavelet power spectrum between SPEI-12 and PDO,AO,ENSO

注:(a)(c)(e)为突变前,(b)(d)(f)为突变后。图5 PDO、AO、ENSO与SPEI-12小波相干谱Fig.5 Wavelet coherent spectrum between SPEI-12 and PDO,AO,ENSO

由图4～5可知:松花江流域SPEI-12指数在高能量区突变前与PDO无明显共振周期(见图4(a)),突变后1993～1997年存在位相角为90°、尺度为3～4 a的共振周期,表示在此期间高能量区SPEI-12变化滞后于PDO变化1/4个共振周期(见图4(b));在低能量区,突变前与PDO无明显共振周期(见图5(a)),突变后存在位相角为90°、尺度为2～3 a的共振周期(见图5(b)),即低能量区SPEI-12变化滞后于PDO变化1/4个共振周期。

SPEI-12与AO在高能量区突变前1966～1971年存在位相角0°、尺度2～4 a的共振周期(见图4(c)),即SPEI-12与AO负相关变化,突变后无明显共振周期(见图4(d));在低能量区突变前存在位相角约270°、尺度2～5a的共振周期(见图5(c)),即SPEI-12变化提前于AO变化1/4个周期,突变后无显著共振周期(见图5(d))。

SPEI-12与ENSO高能量区突变前1966～1975年显示有位相角180°、尺度为1.5～4.0 a的共振周期(见图4(e)),即两者同步变化。突变后在1994～2000年存在位相角0°、尺度为2.5～4.0 a的共振周期(见图4(f)),即两者负相关变化;低能量区1962～1972年存在尺度为2.5～5.0 a的同步变化共振周期(见图5(e)),突变后无显著共振周期(见图5(f))。

3.2.2时滞相关分析

交叉小波与小波相干分析显示,突变后PDO和ENSO与SPEI-12共振现象更加显著,因此对PDO、ENSO与SPEI-12指数序列分别进行突变前后的同期及滞后1～12个月的Pearson相关系数计算,结果如表3、表4所列。

表3 PDO与SPEI-12指数突变前后时滞相关分析Tab.3 Simultaneous and lag correlation analysis between PDO and SPEI-12 index before and after mutation

表4 ENSO与SPEI-12指数突变前后时滞相关分析Tab.4 Simultaneous and lag correlation analysis between ENSO and SPEI-12 index before and after mutation

由表3～4可知:SPEI-12突变前与PDO同期与滞后1～3个月具有逐渐减弱的显著正相关性,滞后9～12个月转变为显著的负相关性;突变后同期与滞后1～12个月具有极显著的正相关性,相关系数随滞后尺度增大逐渐减小。SPEI-12突变前与ENSO同期与滞后1～10个月具有显著的正相关性,其中滞后4～5个月相关性更为显著;突变后仅滞后6～8个月具有较显著的正相关性。

3.3 讨 论

综合上述,松花江流域SPEI-12在1984年左右突变后,在时间、空间上皆整体呈现减小趋势,齐齐哈尔以西北呼伦贝尔以东南地区、黑河市周边地区和吉林北部地区变化最为显著。这与吴燕锋等[38]在分析松花江水文干旱与气象干旱关系中气象干旱在1982年后趋于减轻的演变结论相符。曾鹏[39]也在其对中国未来干旱趋势模拟研究中以RCP4.5与RCP8.5两种模拟情景模式模拟出松花江流域的干燥趋势显著,SPEI指数整体呈现负值,与本文的研究结果一致。在PDO、AO、ENSO与SPEI-12的交叉小波与小波相干分析中,松花江流域突变前后SPEI-12指数与遥相关因子存在共振周期,其中PDO和ENSO与突变后SPEI-12共振相互作用更为显著,表明遥相关因子可能是导致松花江流域SPEI-12指数降低的重要因素之一。

而在进一步的相关分析中,PDO和ENSO与SPEI-12突变前后相关性质的变化与小波分析部分结论相符,即PDO和ENSO是导致松花江流域SPEI-12指数降低的重要的遥相关影响因素。ENSO表现为赤道中东部太平洋海面温度异常增高或降低,暖相与冷相皆会引起沃克环流的变化,Ren等[40]也在研究中证明了松花江流域SPEI指数与不同ENSO指数皆存在显著的相关关系;PDO体现为北太平洋西部和中部海温的异常。两者皆起源于太平洋与大气的相互作用,影响海洋大气水平-垂直相互作用,给太平洋沿岸大陆乃至全球气候带来热量及水汽输送的异常变化。但除了PDO与ENSO的影响外,松花江流域干旱的成因机制还与季风、人类活动和地形地质等有关。并且由于数据的限制,本文选择使用了Thornthwaite方法计算潜在蒸散发,该方法只考虑了温度的影响,忽略了其他如风速、辐射等气候因子的影响,因此松花江旱涝状况的评估与遥相关因子之间的耦合机理仍需进一步研究。

4 结 论

本文结合松花江流域52个气象站点1960～2020年共61a逐日降水与气温数据,计算12个月尺度的标准化降水蒸散指数,识别其时间序列突变点,并运用交叉小波变换与小波相干分析,得出以下结论:

(1) 松花江流域SPEI-12指数突变节点约为1984年,空间上呈现中部及中部偏东地区突变前轻涝、突变后轻旱的变化特征,流域整体SPEI-12指数趋于减小,中部及中部偏东地区变化最为明显。

(2) 交叉小波功率谱与小波相干谱显示,SPEI-12与PDO、AO、ENSO存在共振周期,其中PDO和ENSO突变前后差异更明显,对突变后SPEI-12共振相互作用更加显著。

(3) 同期与滞后相关性分析显示,突变前后SPEI-12与ENSO相关时滞尺度有较大变化。SPEI-12与PDO随时滞尺度变化突变前相关性质具有较大差异,突变后转变为极显著的正相关性,表明PDO与ENSO是导致松花江流域SPEI-12指数突变后数值减小的重要的遥相关影响因素。