张金萍,张 航
(1.郑州大学水利科学与工程学院,郑州 450001;2.郑州大学黄河生态保护与区域协调发展研究院,郑州 450001;3.郑州市水资源与水环境重点实验室,郑州 450001)
由于中上游水土流失、人类活动影响及气候变化等因素导致黄河水生态愈发脆弱[1,2]。降水作为黄河流域水资源的主要制约因素,受地形地貌及季风气候的影响,在黄河流域内部呈现明显的空间异质性。如何从已发生变异的降水变量中分析出水文特征规律[3-5],对于准确把握黄河流域水资源特征,促进流域内工农业生产,制定生态保护措施具有重要意义。
在黄河流域已有的降水研究中,杨特群[6]根据黄河流域1951-2007年降水资料认为其汛期降水量和年降水量变化规律大致相同。王雁[7]则利用黄河流域1951-2008年月降水资料对流域上、中、下游的降水趋势进行了分析。王国庆[8]、何金梅[9]分别从不同角度对兰州上游区间的降水进行了研究分析。周帅[10]基于标准化降水指数(SPI)分析了黄河流域干旱时空演变规律。陈磊[11]基于黄河流域106 个代表气象站1960-2010年逐日降水资料,对季节降水的变化趋势和变异情况进行了相关研究。黄建平[12]则利用欧洲中期天气预报中心ERA5再分析资料对黄河流域过去40年的降水变化情况进行了研究。绍晓梅[13]、李占杰[14]、鱼京善[15]等人采用小波分析或最大熵谱分析的方法对黄河流域2003年之前的降水周期性规律进行了研究。
以上研究主要对黄河流域的降水趋势进行了分析,并采用小波分析和最大熵谱分析对其年降水周期变化规律进行了揭示,但受限于小波窗口和最优截止矩的确定,具有一定的局限性,而对其年降水和汛期降水相关关系的研究则更加少见。为了对黄河流域降水在不同时间尺度的演化规律进行更深层次挖掘,本研究在黄河流域年降水序列趋势检验的基础上,运用信息熵理论计算黄河流域年、汛期降水在不同阶段的时空信息熵,从而对不同时间尺度上的降水复杂性进行识别;使用CEEMDAN 方法对流域年、汛期降水序列进行多时间尺度分解,进而提取原始序列蕴含的周期性规律;最后应用Copula 理论建立年降水和汛期降水的联合分布,从而更精确地把握两者的相依关系。本研究不仅揭示了黄河流域降水的变化趋势,而且从时间和空间维度对其演变规律进行了细致剖析,相关研究结果可为黄河流域生态保护和高质量发展提供支撑作用。
黄河流域贯穿中国半干旱半湿润地区,大部分区域属温带大陆性季风气候。本研究将黄河流域分为8 个二级水资源区,分别为龙羊峡以上(简称龙上区间)、龙羊峡至兰州(简称龙-兰区间)、兰州至河口镇(简称兰-河区间)、内流区、河口至龙门(简称河-龙区间)、龙门至三门峡(简称龙-三区间)、三门峡至花园口(简称三-花区间)、花园口以下(简称花下区间)。
流域内部分气象站点由于迁移、新建原因导致数据资料序列较短且存在残缺,若使用插补处理会对精度产生较大影响。因此以中国国家青藏高原数据中心(http://data.tpdc.ac.cn/)空间分辨率为0.008 333 3°(约1 km)的中国逐月降水数据为基础进行研究分析。该数据集是通过Delta 空间降尺度方案在中国地区降尺度生成,其结果已使用496 个独立气象观测站点数据进行验证,验证结果可信。为了进行相关指标的计算,在研究中以黄河流域内321 个代表气象站点基础(如图1),提取各个气象站点所在栅格1960-2017年的月降水数据作为站点资料进行研究。
1.2.1 信息熵
熵最初来源于热力学,1997年,Singh[16]较为详细地介绍了信息熵在水文水资源等方面的应用,国内也有学者[17,18]使用信息熵对部分区域的降雨规律进行了研究。熵的计算如公式(1)所示:
式中:pi为变量xi所对应的概率。
Ebrahimi 等人[19]对信息熵的计算公式进行了改进,提出了在总体分布未知情况下的熵值计算公式,如公式(2)所示:
其中:
式中:yi是原始数据升序后的序列。当i≤m时,yi-m=y1;i≥nm时,yi+m=yn。
1.2.2 CEEMDAN分解
CEEMDAN 算法通过在原始信号上加入有限方差约束的多组独立同分布自适应白噪声实现对非平稳原始信号不同频率成分的提取[20]。设x(t)为待分解序列,w(t)表示均值为0,方差为1 的高斯白噪声,r(t)为残余差值,IMFk为第k阶模态分量,CEEMDAN分解过程如下:
(1)向x(t)中加入高斯白噪声,并进行I次实验,第i次实验信号可表示为:
(2)对第i次加入白噪声后的信号xi(t)进行EMD 分解并进行I次实验平均,得到第一阶模态分量IMF1,然后计算x(t)与IMFi分量的残余差值r1(t)。
(3)再次向残余差值r1(t)中加入白噪声,进行I次实验平均,第i个残余差值可表示为:
(4)对第i次加入白噪声的信号ri1(t)进行EMD 分解,第i次为IMFi2分量,通过I次实验平均,得到第二阶IMF分量为:
(5)此时,分解得到的残余差值为r2(t) =r1(t) -IMF2,返回步骤(4)、(5)计算下一阶模态,则第k阶模态为:
重复第3 步至第5 步直到差值rk(t)满足以下条件之一时终止分解:①不能被EMD 进一步分解;②IMF分量满足相应条件;③局部极值点的个数小于3 个。最终,原始信号x(t)可分解为若干个IMF分量和一个趋势项rK(t):
1.2.3 Copula函数
Copula 函数具有将联合分布和各自的边缘分布联结起来的作用,从而实现变量间随机性和相依性的分离,广泛应用于降雨、洪水、干旱和枯水等多特征属性频率分析研究中[21-23]。应用最广泛的为Archimedean Copula,其类型主要包括Gumbel Copula、Clayton Copula 和Frank Copula 三种(见表1)。表中θ为参数,τ为Kendall秩相关系数。
表1 常用Archimedean Copula函数Tab.1 Three commonly used Archimedean Copula functions
1960-2017年间,黄河流域年均降水量500.6 mm,龙上、龙-兰、兰-河、内流区、河-龙、龙-三、三-花、花下区间各水资源区的年均降水量分别为543.4、458.5、250.3、349.1、443.3、563.6、619.6、654.5 mm。整体来看,黄河流域河源区及下游降水较为丰富,高于流域平均水平。
MK 检验结果(图2)表明,黄河流域年降水序列自20 世纪70年代以来开始呈现一种明显的下降趋势,在20 世纪初,这种下降趋势超过了0.05 的显著性水平(u0.05= 1.96)。这和王飞等人[24]采用SPEI指数得到黄河流域干旱呈显著增加趋势的结论是一致的。河-龙、龙-三、三-花和花下4 个水资源区变化趋势与流域整体比较接近;流域上游的龙上和龙-兰水资源区年降水序列具有缓慢增加的趋势,王远见等人[25]通过研究同样发现黄河源区年降雨量呈增加趋势,只不过在其结论中认为这种上升趋势是显著的;兰-河、内流区水资源区年降水序列主要呈下降趋势,但相对于流域整体而言,其下降趋势并不明显。
此外,从图2中UF统计量可以初步判断黄河流域各水资源区年降水序列在2000年左右变化最为显著,为进一步确定突变时间,结合滑动T检验综合确定2002年为显著突变时间点,并以此划分突变前(1960-2001年)和突变后(2002-2017年)时段。
图2 年降水序列MK检验结果Fig.2 MK test results of annual precipitation series
2.2.1 降水时间熵
表2、3 分别为各水资源区年、汛期降水量突变前后两阶段的时间信息熵。由表2 可知,黄河流域年降水序列时间熵在突变后相比突变前降低5.8%,即2002年以后,黄河流域降水量在年际间的波动减弱。各水资源区中,龙-兰水资源区时间熵减小幅度最大,降水量在年际间的波动相比突变前出现了显著降低;三-花、花下水资源区年降水时间熵变化幅度较小,降水在年际间波动较为稳定,其余水资源区降水在年际间的波动位于两者之间。
表2 年降水序列突变前后时间信息熵Tab.2 Temporal information entropy value of annual rainfall in different stages
由表3 可知,兰-河、黄河内流区水资源区所对应的汛期降水时间熵突变后相比突变前降幅均超过10%,说明2002年以后汛期降水量在年际间的波动显著减弱;三-花、花下水资源区汛期时间熵变化较小,说明其汛期降水量在整个研究时段内的波动比较平稳。
表3 汛期降水突变前后时间信息熵Tab.3 Temporal information entropy value of flood season rainfall in different stages
此外,对比年和汛期时间信息熵可得到如下结论:
(1)无论年尺度还是汛期尺度,突变后的时间信息熵相比突变前均出现了降低,说明近些年的降水量虽然呈下降趋势,但降水在黄河流域和各水资源区内部沿时程的复杂程度也出现了下降,不确定性降低,有利于对降水规律进行研究。
(2)降水时间信息熵在不同水资源分区内具有空间差异性。首先表现为突变前后各水资源区在汛期和年尺度上的时间信息熵值大小不同,且发生变异后下降幅度不一。其原因主要是不同水资源分区在地理位置、局部气候、地势地貌等方面不同所导致。其次,位于龙-三水资源分区内的部分气象站点汛期降水时间熵在突变后相比突变前出现了增加,但对龙-三水资源分区而言,突变后的信息熵相较于突变前是下降的,因此存在局部降水复杂性上升、不确定性增加的现象。
2.2.2 降水空间熵
由黄河流域年降水空间熵结果(见图3)可以发现:
图3 年降水空间信息熵Fig.3 Spatial information entropy value of annual rainfall
(1)2001年以后,黄河流域年降水空间熵出现小幅度降低,说明从流域整体角度而言,突变后年降水量在空间上的不确定性出现了下降,但这种趋势并不明显;
(2)河-龙及龙-三水资源分区空间熵出现了升高,其余水资源分区的空间熵下降,说明突变以后各水资源区降水复杂性出现了不同方向的改变,但这种改变以复杂性降低为主。
通过对空间熵进行MK 检验发现,在研究时段内,黄河流域和龙上、龙-兰、兰-河、内流区、河-龙、花下水资源区在1990年后空间熵呈现下降趋势,龙-三、三-花水资源区的空间熵则呈现一种微弱的上升趋势。说明黄河流域内部降水在空间上的复杂性主要是降低的,即降水的空间分布更加均匀。
进一步对黄河流域年降水空间熵(Hs)和年降水变异系数(Cv)进行相关分析(图4)可以看出,降水空间熵和变异系数之间存在明显的正相关。变异系数代表着流域内部各站点降雨的不均匀程度,变异系数越大,说明各流域内部降水越不均匀。由此可见,空间差异性蕴藏在空间信息熵所蕴含的复杂性之中,即差异性是复杂性的一种体现。此外,黄河流域年降水空间熵>汛期降水空间熵>9月降水空间熵,说明时间跨度越长,所包含的不确定性越大,这与信息熵所代表的物理含义是一致的。就黄河流域而言,其汛期空间熵突变后是降低的,但9月份降水空间熵在突变后是增大的。由此可见,黄河流域9月份降水所包含的不确定性是增大的,汛期和年的变化相一致。因此在今后的降雨观测过程中,应密切注意9月份降水的变化情况。
图4 Hs与Cv相关关系Fig.4 Correlation between Hs and Cv
图5(a)为黄河流域年降水序列,具有较强的非线性和非平稳特性。图5(b)、(c)为多时间尺度分解结果,突变前的原始序列最终被分解为5 层,包括4 个IMF分量和1 个趋势项(RES余量);突变后的原始序列被分解为4层,包括3个IMF分量和1个趋势项。分解出的IMF分量较好地体现了原始信号序列在不同周期上的波动规律和变化特征。
图5 黄河流域年降雨序列多时间尺度分解Fig.5 Multi-time scale decomposition results of annual rainfall series
突变前后各IMF分量周期波动和振幅变化特征如表4所示。
由表4可知,突变前:从IMF1分量到IMF4分量,周期增大,振幅减小;IMF1 分量的准周期在2~3 a 范围内,说明其变化频率较高,且其最大振幅和平均振幅都处于各分量的首位;IMF1~IMF4 和RES分量的方差贡献率分别为74.41%、11.9%、2.64%、1.23%和9.81%,由此可以看出,具有短周期、高振幅的高频分量包含了原始序列较多的信息量,在原始序列中起着主要作用;IMF1分量的平均振幅较大,IMF2分量的平均振幅位居其次,IMF3、IMF4 分量的平均振幅相差不大,这与最大振幅的变化规律是相一致的。从突变前的RES余量可以看出,黄河流域年降雨量在突变前呈现一种下降趋势。
表4 年降水序列多时间尺度分解结果Tab.4 Multi-time scale decomposition results of annual rainfall series
突变后:原始序列分解出3 个IMF分量和一个趋势项。IMF1~IMF3 分量周期逐渐增大,IMF1 分量包含2 a 和3 a 两个准周期,IMF2 分量为4 a 的准周期,IMF3 分量为9 a 的准周期。各分量的贡献率分别为84.23%、7.08%、6.10%和2.59%;RES分量呈现先下降后上升的变化趋势,前半程的下降趋势对应着突变前RES分量的下降趋势,是相一致的,后半程的上升趋势表明从全局而言,黄河流域降水在最近几年逐渐开始恢复,这和原始序列的变化趋势也是相一致的。
对各二级水资源区进行CEEMDAN分解,鉴于篇幅原因,图6 仅列出龙上和花下区间年降水序列多周期分解结果,并得到以下结论:
图6 各水资源区年降水多周期分解结果(以龙上和花下区间为例)Fig.6 Multi-time scale decomposition results of annual rainfall in secondary water resources area
(1)各水资源区IMF1 分量准周期均为2~4 a,即占主要贡献的高频分量周期基本一致,表明黄河流域水资源区包含着相似的短周期变化规律。因此,无论是对黄河流域还是各水资源分区进行研究时都应关注短周期尺度上的变化规律,即对2~4 a的周期进行重点研究。
(2)突变后序列分解出的IMF分量个数均减少,造成该结果的原因可能有两点:①突变后黄河流域内部降水的复杂程度和不确定性是降低的,包含的信息量减少,导致分解出的IMF分量减少;②突变后原始数据的时间序列较短。
(3)由RES分量可知,突变前各子水资源区降水量变化趋势和黄河流域整体基本一致,呈下降趋势;突变后,龙-兰、兰-河、内流区、河-龙、龙-三5 个水资源区和流域整体在研究时段末呈现一种上升趋势,龙上、三-花和花下3 个水资源区在研究时段后期则呈现一种下降趋势。表明突变前各水资源分区和流域整体的变化趋势是一致的,都呈现一种下降趋势,但突变造成各水资源分区存在一定差异性。
通过对各IMF分量进行归一化处理,计算得到多时间尺度熵(表5),其代表了原始序列在不同时间尺度上的变化特征:
表5 黄河流域年降水序列突变前后多时间尺度熵Tab.5 Multi-time scale entropy of rainfall in different stages
(1)突变后各IMF分量熵值相比突变前均出现不同程度的下降,说明变异导致降水要素在不同尺度上的复杂程度下降,随机性降低,这与原始序列时空信息熵的变化规律是一致的。
(2)无论突变前还是突变后,IMF1分量的熵值最大,表明高频短周期分量所包含的信息量比较丰富,最能体现原始序列的变化规律。因此对于黄河流域而言,对降水的研究应集中在2~4年的短周期上。
黄河流域年降水和汛期降水关系密切,具有较强的相关关系。本研究使用参数法确定各变量边缘分布,并利用Kolmogorov-Smirnov(KS)方法进行检验,检验结果显示,突变前年降水、汛期降水分别以伽马分布和韦伯分布进行拟合时效果最优,见图7。
图7 突变前年降水序列经验频率和理论频率拟合图Fig.7 Fitting results of empirical frequency and theoretical frequency of annual rainfall before mutation
采用相关性指标法对Copula 函数进行参数估计,并以AIC信息准则、OLS规划原则作为拟合优度评价指标,结果(表6)表明Gumbel-Copula 函数在构建突变前年降水量和汛期降水量之间的联合分布时,拟合效果最优。
表6 突变前拟合优度检验结果Tab.6 Goodness of fit test results before mutation
突变前最终的Copula 函数如式(10)所示,其分布函数与等高线图如图8所示:
图8 年雨量-汛期雨量联合分布Fig.8 Joint distribution of annual and flood season precipitation
同理,求出突变后最优Copula函数如式(11)所示:
以37.5%、62.5%分别作为丰枯分界,计算突变前后年、汛期降水丰枯遭遇概率,计算结果如表7、8所示,可以发现:
表7 突变前年-汛期降水联合概率Tab.7 Joint probability of annual and flood season precipitation before mutation
(1)不同阶段(突变前后),年和汛期降水联合分布丰枯同步频率大于丰枯异步概率,且同步频率主要集中在同丰或同枯状态;
(2)降水要素发生突变后,同步频率增大,从突变前的0.655 增加到突变后的0.722。主要体现在同丰和同平状态的增加。
上述两个结论表明2002年以后黄河流域降水在年和汛期尺度上更加一致,即年内降水更集中在汛期,这对于防洪是不利的,需要有关部门提前做好防汛措施。
以黄河流域1960-2017年月降水数据为基础,基于信息熵、多时间尺度和Copula 函数分析了黄河流域降水在时间和空间上的变化规律,主要得出以下结论。
(1)过去50年间,黄河流域降水整体表现出一种下降趋势,但源区降水在缓慢增加,且降水在黄河流域内部呈现明显的空间异质性;
(2)突变后(2002年后),流域降水时空信息熵均出现了降低,其中龙羊峡-兰州水资源区的时间信息熵减小幅度最大,超过10%,空间熵则在黄河流域内部呈现不同方向的变化;
(3)通过CEEMDAN 方法将突变前后黄河流域年降水序列分别分解为4层(3个IMF分量和1个RES余量)和3层(2个IMF分量和1 个RES余量)。结果表明对黄河流域降水的观测和研究应主要集中在2~4 a 的短周期上,且RES余量反映了原始序列整体下降但在近几年逐渐恢复的趋势;
表8 突变后年-汛期联合降水概率Tab.8 Joint probability of annual and flood season precipitation after mutation
(4)通过Gumbel-Copula 函数发现在不同阶段黄河流域年、汛期降水联合分布以丰枯同频为主,同丰和同枯的重现期分别在3~4 a、4~5 a之间。此外,2002年以后,同步频率增加,年内降水更集中在汛期,对于防洪提出了更高要求。
本研究主要基于信息熵、CEEMDAN 和Copula 理论对黄河流域年、汛期降水的时空变化规律进行了研究,而黄河流域极值暴雨同样是影响流域内社会经济活动和生态稳定的重要因素。受限于资料原因,在此并未对其进行研究分析,未来在获取更短时间步长降水数据的基础上,对极值暴雨事件开展研究,对于指导流域内生产生活具有更加实际的意义。 □