关中地区气象干旱时空变化特征分析

2023-10-31薛瑞华罗军刚李少轩

西安理工大学学报 2023年2期

薛瑞华, 罗军刚, 李少轩, 景鑫

(西安理工大学西北旱区生态水利国家重点实验室, 陕西西安 710048)

干旱是一种在世界各地频繁发生的慢性自然灾害,其发生与缓解都不如其他自然灾害那么迅速[1],且干旱的发生通常由地区长期异常降水造成[2],这导致与其他大多数自然灾害相比,干旱的影响更大、范围更广、持续时间更长。干旱通常分为四类:气象、农业、水文和社会经济干旱[3]。气象干旱主要表现为降雨量不足[4]。当一定时期内的降水量持续少于正常状态时,河流径流量会相应减少形成水文干旱,当干旱持续发展时,土壤含水量将继续减少直到不能满足生物生长的需水量而形成农业干旱,最终导致作物减产造成经济损失形成社会经济干旱[5]。因此,气象干旱是形成其他类型干旱的原因,监测气象干旱对减轻干旱危害至关重要[6]。

关中地区是陕西省最典型、最重要的农业区,受气候地形的影响,关中地区干旱频发。据统计,近百年来关中地区河流径流量严重减少,水资源开发难度大,但随着经济的发展以及人口的增长,对水的需求迅速增加,如今缺水率约为25%且仍在上升[7]。最近一次大型旱灾是2014年8月上旬,关中地区遭遇了57年以来最严重的大旱,渭河、洛河等河流的径流比往年同期偏枯一至四成,近20万人用水紧张,咸阳65万亩作物受旱,西安97万亩作物受旱、11万亩绝收[8]。因此,需密切关注干旱持续时间、严重程度等干旱特征,分析其变化趋势与重现期,并根据干旱特征采取合理的应对措施以减少旱灾造成的影响。

目前,关于关中地区的气象干旱研究已有许多。Zhang等[9]通过标准化降水蒸散指数分析了1955—2015年间关中盆地的干旱变化趋势;粟晓玲等[10]采用Gumbel Copula函数构建气象水文综合干旱指数,探讨了该指数的适用性,并使用该指数分析了关中地区的干旱周期与频率;Chen等[11]探讨了1928—1931年陕西关中地区极端干旱事件的社会影响。综上,已有的关中地区气象干旱研究,主要是基于干旱指数进行干旱频率、干旱强度等特征变量的单独分析,少有对多特征变量进行联合分析的。而干旱的发生与多个变量相关,其特征以多种方式体现,多特征变量联合分析可以更好地表征复杂的干旱条件及其影响[12],因此,分析干旱多特征变量的联合关系对研究干旱有着十分重要的意义。

Copula是最常用的多元分析方法之一,适用于多变量分析,已在极端降水、洪水频率分析和干旱特征分析等方面得到广泛应用[13]。刘瑞琳等[14]利用游程理论、单变量和多变量Copula方法,得到了横断山地区的干旱特征在不同重现期下的联合分布特征。李京芳等[15]运用Copula函数计算干旱特征变量之间的联合累积概率,分析了联合重现期和同现重现期以及汉江流域水文干旱风险及驱动因素。因此,基于Copula的干旱特征分析在单独考虑边际分布或多变量相关性建模方面具有优势,然而在使用Copula进行分析时,存在参数计算、边际函数选择和Copula函数选择等多个方面的不确定性[16]。

贝叶斯网络利用可信区间解决了参数的不确定性问题[17],因此本文采用包括基于贝叶斯框架的马尔可夫链蒙特卡罗方法(Markov chain Monte Carlo,MCMC)与局部优化算法在内的多种方法进行参数计算,优选最合适的参数,并通过比较广泛的边际函数与Copula函数来有效解决不确定性问题,为关中地区多特征变量联合分析提供有效方法。

本文基于SPI12时间序列,通过游程理论识别干旱特征,采用局部优化算法与基于贝叶斯框架下的马尔可夫链蒙特卡罗模拟两种算法进行Copula参数估计,得到最优的Copula函数及其参数;对干旱历时与干旱烈度进行联合分析得到关中地区的干旱分布特征,可为关中地区干旱风险评估与水资源管理规划提供依据。

1 研究区概况与数据

关中地区总面积5.55万km2,位于陕西省中部,属于渭河流域,由南向北从秦岭延伸至北山,由西向东从宝鸡延伸至黄河,处于温暖带半湿润与半干旱气候的过渡地带,属大陆性季风气候[18],西部、北部和南部海拔较高,东部和中部地区海拔较低。地貌类型主要由漫滩、阶地和黄土台地组成,行政区域包括杨凌示范区以及西安、宝鸡、咸阳、渭南、铜川5市,是陕西省经济最发达的地区。

本文收集到关中地区20个气象站1979—2018年的日降水资料(来源于国家青藏高原科学数据中心http://data.tpdc.ac.cn/zh-hans/),在对数据进行计算前,对缺测的数据采用链式方程多重插补法进行插补。研究区域与气象站的分布如图1所示。

图1 关中地区气象站点分布

2 研究方法

2.1 标准化降水指数

标准化降水指数(standardized precipitation index, SPI)是用来表示某时段降水量出现概率多少的指标,适用于不同地区的干旱检测与评估[19]。

(1)

(2)

(3)

其中,H(x)=q+(1-q)G(x),q是降雨量为0的概率,G(x)是服从伽马分布的累计概率;C0=2.515 517,C1=0.802 853,C2=0.010 328,d1=1.432 788,d2=0.189 29,d3=0.001 308。

SPI指数计算结果:SPI>-0.5表示无旱,-1.0

SPI指数用于度量多种时间尺度上的降水不足,关中地区干旱频发,几乎年年有旱,因此采用12个月时间尺度的SPI值(记为SPI12)反映长期降雨模式,分析关中地区气象干旱的年际变化。

2.2 游程理论

游程理论是由Yevjevich[20]提出的研究干旱特征最常用的方法之一[21]。考虑到实际中会突然出现较大降雨,导致一次严重干旱事件的过程分割,从而弱化干旱的实际严重程度,故本研究采用王晓峰等[22]提出的三截取游程理论,并在此基础上进行调整,得到适用于SPI指数的干旱事件识别游程理论。

1) 根据《气象干旱等级》(GB/T 20481—2017),当SPI值小于或等于-0.5时,有干旱发生。

2) 若只识别出单独一个月的事件,且这个月的SPI值大于-0.75,就认为这个月没有发生干旱。

3) 若两次干旱事件之间仅间隔一个月,且间隔月的SPI值小于0,则将两次干旱事件合并为一次干旱事件。

通过确定出的干旱事件,基于三截取游程理论可提取出持续时间和烈度等干旱事件的基本特征[23]。其中,干旱历时表示单次干旱过程的持续时间,干旱烈度表示单次干旱过程持续时间内SPI的累加值。

2.3 自适应Copula函数选择

自适应Copula函数选择是依据拟合优度的评价结果选取其中最优的方法作为Copula参数的估计方法。本文采用了两种算法进行Copula的参数估计,分别是局部优化算法与基于贝叶斯框架下的马尔可夫链蒙特卡罗模拟算法。

2.3.1边缘分布与相关性检验

在Matlab中采用MvCAT工具箱所包含的指数分布、伽马分布、威布尔分布等在内的16种函数进行单变量边缘分布函数拟合,并对拟合后的函数进行K-S检验,通过赤池信息准则(Akaike information criterion,AIC)进行拟合优度检验,根据AIC最小值原则选出合适的边缘分布函数。

采用Kendall秩相关系数法、Pearson相关系数法和Spearman秩相关系数法对干旱历时与干旱烈度进行关联性检验。

2.3.2Copula函数

Copula函数通常用于从单变量边际分布生成多变量分布,该技术的理论基础是由Sklar[24]定义的。根据Sklar定理,如果FXY是相依随机变量X和Y的二维联合分布函数,其边缘分布函数为FX和FY,则存在一个CopulaC,即:

P(X≤x,Y≤y)=FXY(x,y)=C(FX(x),FY(y))

(4)

上面定义的函数FXY是一个具有边缘分布函数FX和FY的二维函数。

本文采用了MvCAT工具箱中的Gaussian、t、Clayton、Frank、Gumbel等26种Copula函数进行计算[25]。

2.3.3局部优化算法

本文采用的局部优化算法是一种基于梯度的内点优化算法[26]。此算法通过一个稠密的拟牛顿近似方法(dense quasi-Newton approximation)估算Hessian矩阵,并使用内点法寻找可行解空间。

在求解参数的过程中,每次的参数估计将从随机的30个点出发寻找最优解,以此减少陷入局部最优的概率。

2.3.4贝叶斯框架

贝叶斯统计推断对于不确定性评估和参数估计是有效的,其中假定参数为具有概率分布的随机变量[16]。贝叶斯推理的基本定理是:

(5)

2.3.5马尔科夫链蒙特卡洛模拟

本研究使用一种新混合进化的MCMC方法,该方法采用自适应采样分布来描述贝叶斯环境中的后验参数区域。首先通过拉丁超立方抽样(Latin hypercube sampling, LHS)算法对整个先验空间进行随机搜索,通过样本的随机分配从中选择出可能性值最高的样本作为马尔科夫链的起点,这种智能起点的选择使不同吸引区域的点作为起始点成为可能,从而避免了退化。在进行后验参数空间搜索时,采用10%概率的斯诺克方向更新和90%概率的并行方向更新,在并行方向中有NAM条马尔科夫链采用自适应采样算法(adaptive metropolis, AM) 丰富跳跃方向,增强MCMC早期的搜索能力,其他链则采用差分进化算法(differential evolution, DE)加快收敛速度。

2.3.6拟合优度评价

在本研究中采用均方根误差(root mean square error,RMSE)、纳什效率系数(Nash-sutcliffe efficiency coefficient,NSE)以及AIC进行Copula函数拟合优度检验。

1) 均方根误差

(6)

(7)

(8)

2) 纳什效率系数

(9)

3) 信息准则

AIC=nln (MSE)+2

(10)

式中:Pei为联合分布经验频率;Pi为联合分布理论频率;n为样本容量;mi为干旱特征小于某一阈值的干旱事件发生的次数。

2.4 重现期计算

单变量重现期通常会高估或低估某一特定事件的风险率。因此,本研究计算了气象干旱特征的二变量重现期。该方法假设X1、X2代表与干旱事件相关的特征变量值,x1、x2代表与干旱事件相关的特征变量序列对应干旱特征的给定阈值。

二维联合回归周期可以表示为:

(11)

二维共现重现期可以表示为:

(12)

式中:N表示总观测次数;n表示样本容量。

如图2所示,本文主要包括干旱指数计算、干旱特征识别、边缘分布函数筛选、构建Copula联合分布以及关中地区干旱特征分析这五个部分。

图2 干旱特征分析流程图

3 结果分析

3.1 干旱特征

根据各测站的月累计降雨量数据,采用算数平均法计算出关中地区五个市的逐月累计降雨量,并根据逐月累计降雨量得到各市的SPI12时间序列。关中地区及五个市1980—2018年逐月累计降雨量与SPI12时间序列的变化如图3所示。从图中可以看出,关中地区整体的SPI12指数呈现下降趋势,宝鸡市与渭南市的SPI12指数有较明显的下降趋势,西安市和铜川市的SPI12指数呈现不明显的下降趋势,咸阳市的SPI12指数呈现不明显的上升趋势,表明在研究时段内关中地区整体的干旱愈发严重。

图3 关中地区月累计降水量与SPI12指数变化图

从图3可知,关中地区经常发生轻旱与中旱,但重旱与特旱的情况较少。关中地区在1997年间干旱情况十分严重,这与1994—1997年间降雨量偏少有着十分密切的关系,根据《中国气象灾害记录》记载,关中地区1995年全年降水严重偏少,其中上半年偏少达6～8成,9月初至12月底降水偏少达5～8成,1997年4～7月降水比常年偏少3～9成,8、9月降水偏少6～9成[27],与研究结果一致。从图3中还可知,近年来较为严重的干旱发生在2014年,这与报道的多个县区降水量偏少5成以上,出现大面积未灌或少灌区域,部分地区的猕猴桃出现枯死现象的结果一致。

通过游程理论,对各市SPI12时间序列进行识别,得到宝鸡市、咸阳市、西安市、铜川市以及渭南市的干旱次数分别为16次、14次、17次、21次以及18次。图4为不同地区干旱历时与干旱烈度的分布情况。关中地区整体的平均干旱历时为7.55月,其中宝鸡、咸阳、西安的平均干旱历时分别为8.94月、9.86月以及9.24月,高于关中地区的平均干旱历时;渭南的平均干旱历时为7.83月,与关中地区的平均干旱历时相接近;铜川的平均干旱历时为6.38月,低于关中地区的平均干旱历时。关中地区的平均干旱烈度为8.10,宝鸡、咸阳、西安、渭南的平均干旱烈度分别为9.90、11.06、9.70以及9.24,均高于关中地区的平均干旱烈度;铜川的平均干旱烈度为7.10,低于关中地区的平均值。宝鸡、咸阳以及西安的平均干旱历时与平均干旱烈度均较大,干旱情况较为严重;渭南市平均干旱历时较小,平均干旱烈度较大,且烈度的中位数远远小于平均值,表明渭南市多发小干旱事件,但也有极端干旱情况的出现;铜川市的平均干旱历时与平均干旱烈度均较小,表明铜川市的干旱情况较轻。

3.2 边缘分布

表1给出了20个测站识别出的干旱历时与干旱烈度的边缘分布函数与相关性检验。由表1可知,关中地区的干旱历时主要服从指数分布和广义帕累托分布,占所有测站的70%,其次是疲劳寿命分布以及逆高斯分布;干旱烈度主要服从指数分布和广义帕累托分布,占所有测站的80%,其次是疲劳寿命分布以及逆高斯分布。采用Kendall rank、Spearman’s rank-order以及Pearson product-moment三种相关系数进行检验,见表1的第4～6列。从表1中可以看出,除了扶风站,所有测站的Kendall值都在0.8以上,Pearson相关系数都在0.9以上,Spearman相关系数除长武站都在0.9以上,且所有站点的拟合结果均通过显著性为5%的K-S检验,表明所有测站的干旱历时与干旱烈度均有密切的相关关系。

表1 干旱历时与干旱烈度边缘分布情况与相关性检验

3.3 Copula联合分布

采用MvCAT工具箱内的局部优化算法与贝叶斯框架下的马尔科夫链蒙特卡洛模拟两种方法进行参数估计,并结合RMSE、NSE、AIC三种拟合优度评价指数选择出最优的Copula函数,结果见表2,其中θ1、θ2分别对应Copula函数中待定参数的最优值。从表2中可以看出,有30%的测站采用局部优化算法,有55%的测站采用马尔科夫链蒙特卡洛方法,有15%的测站两种方法的最优情况相同。适合关中地区干旱历时与干旱烈度联合分布的Copula函数主要为Gaussian、Clayton以及BB1这三种函数。其中对BB1函数的参数估计更适合使用局部优化算法,对Gaussian以及Clayton函数的参数估计更适合使用马尔科夫链蒙特卡洛模拟。

3.4 重现期分析

干旱重现期的计算结果结合了干旱发生的概率,是一种可用来评估干旱风险的简单实用的方法[28]。根据计算结果与克里金插值法绘制干旱历时与干旱烈度的联合重现期与同现重现期,如图5、图6所示。其中(a)、(b)、(c)、(d)4个图分别为根据关中地区总体干旱特征分布情况确定的小干旱特征事件(干旱历时为2个月,干旱烈度为1.5)、平均干旱特征事件(干旱历时为8个月,干旱烈度为8)、严重干旱特征事件(干旱历时为12个月,干旱烈度为12)以及极端干旱特征事件(干旱历时为26个月,干旱烈度为20)下的重现期。

图6 关中地区干旱历时与干旱烈度同现重现期分布

从图5中可以看出,对于小干旱事件,关中地区整体的联合重现期均在一年多;平均干旱事件的短历时联合重现期主要分布在宝鸡以及西安市东部地区,为1 a多,长历时联合重现期主要分布在铜川、咸阳东部、渭南西部以及西安北部地区,最长为4.7 a;严重干旱事件的联合重现期自西向东呈增加趋势,联合重现期最短出现在宝鸡的扶风,为2.19 a,联合重现期最长出现在铜川市的洛川测站,为5.3 a;极端干旱事件长历时联合重现期主要分布在铜川东部与渭南西北部,最长联合重现期出现在铜川市的洛川测站,为11 a,短历时联合重现期主要分布在宝鸡东南部与西安西部地区,最短联合重现期出现在西安市的武功测站,为3.7 a。

从图6中可以看出,长同现重现期主要分布在铜川市,短同现重现期主要分布在宝鸡市。小干旱事件的短同现重现期主要集中在铜川市以及咸阳市东部地区,最短同现重现期出现在宝鸡市扶风站,为1.1 a,长同现重现期主要集中在宝鸡市、西安市西部以及渭南市东部地区,最长同现重现期出现在铜川市的宜君站,为1.6 a;平均干旱事件的短同现重现期主要集中在宝鸡市以及渭南市东部,最短同现重现期出现在宝鸡市扶风站,为2 a,长同现重现期主要分布在铜川市,最长同现重现期出现在铜川的宜君站,为5.5 a;严重干旱事件的短同现重现期主要分布在宝鸡市、西安市西部、咸阳市西南部以及渭南市东部地区,最短同现重现期出现在宝鸡市的扶风站,为3.1 a,长同现重现期主要分布在铜川市以及渭南市西部,最长同现重现期出现在铜川市白水站,为10 a;极端干旱事件的同现重现期自西向东呈增加趋势,最短同现重现期出现在咸阳市的永寿站,为10.6 a,最长同现重现期出现在铜川市白水站,为138.9 a。

从图5、图6中可以看出,联合重现期与同现重现期有着相似的分布规律,由于同现重现期的出现有着更高的要求,因此同现重现期随着烈度历时的增大而快速增加。关中地区最频发的为小干旱特征事件,几乎每年都会出现,其他干旱事件自西向东重现期逐渐增加,长重现期主要分布在铜川市以及渭南市西北部地区,短重现期主要分布在宝鸡市、西安市西部以及渭南市东部少部分地区。

关中地区干旱特征呈现此种变化的原因与降水分布有关,通过冯晶等[29]对关中地区降雨特征的分析可知,除咸阳市外,其他城市的降雨量均呈下降趋势,这与本文的研究结果一致,同时冯晶还指出,降雨自西到东呈增加趋势,因此采用基于降水数据的SPI12指数计算出的重现期也呈现西边重现期短、东边重现期长的变化特征。Zhang等[9]通过计算干旱频率得到关中地区频发轻旱事件,这也与本文的研究结果一致。

本文通过重现期计算,得到关中西部、西安周边以及渭南市东部地区的重现期较短,干旱较为严重;乔亮[30]通过对关中地区旱灾累积程度的分析,认为宝鸡、西安、蓝田、韩城、潼关(宝鸡位于关中西部,西安和蓝田位于西安,韩城以及潼关位于渭南市东部地区)的旱灾累积程度最为严重,这也与本文的研究结果一致。