基于SPEI的广西干旱综合特征及农业旱灾风险研究

杨星星, 杨云川,2,3, 邓思敏, 廖丽萍,2,3, 谢鑫昌, 田 忆

(1.广西大学 土木建筑工程学院, 南宁 530004; 2.广西大学 工程防灾与结构安全教育部重点实验室,南宁 530004; 3.广西防灾减灾与工程安全重点实验室, 南宁 530004)

干旱是最常见的自然灾害之一，其历时长，影响范围广，对自然生态系统破坏大，严重影响着人类社会的发展[1]。干旱不仅发生在降水匮乏、气候干燥的地区，在降水丰沛、气候湿润的地区也常有发生，且由于物候影响，湿润地区植物生长过程中耗水量较大，干旱爆发对其生态系统的破坏也更加严重[2]。2010年，中国西南地区大旱便是典型的湿润地区干旱事件，该次干旱事件造成上千万人饮水困难，30多万hm2农作物受灾，直接经济损失近200亿元[3]。而在全球气候变暖的大背景下，干旱的强度及频次均显著增多[4]，使得抗旱减灾工作日益严峻。

干旱的巨大危害性也吸引了众多学者的关注[5-7]。张强等[2]分析了中国南方地区干旱受灾面积的时间特征，并以干旱强度为致灾因子评估了区域旱灾风险。刘小刚等[8]基于SPI和SPEI指数研究了云南省过去及未来的干旱频率特征。何娇楠等[9]则以SPI值作为干旱致灾因子采用灾害风险理论评估了云南省旱灾风险。贾艳青等[10]基于日SPEI指数分析了西南地区极端干旱事件的趋势变化及空间分布等特征。柳媛普等[6]以CI指数作为干旱强度评估了气候变暖背景下西南地区干旱风险。康永辉等[11]研究了广西大石山区的干旱特性并分析了其成因，结果表明岩溶地质特性以及水土流失严重等因素是当地旱灾加剧的主要原因。综上可知，目前对于典型湿润地区——西南地区的干旱特征以及旱灾风险等方面的研究已有很多，但相关研究也存在一些不足：(1) 大多数研究仅采用干旱指数值进行特征分析及风险评估，忽略了干旱的多维度特征，而干旱的致灾危险性是由干旱强度以及历时等多种特征共同决定；(2) 大多数研究仅对干旱的趋势或频率等一两种特性进行分析，缺乏对干旱特征的全面综合分析研究。

广西是西南地区的重要组成部分，作为典型的亚热带季风湿润区[12]，其年降水量丰沛但时空分布不均，且岩溶发育的石山丘陵、洼地广泛分布，土层浅薄且保水能力差，产汇流过程复杂[13]，致使多年来区域洪涝、干旱灾害频繁，损失巨大[14]，区域干旱特性及旱灾风险研究迫在眉睫。籍此，本文基于1961—2016年广西降水及气温数据计算出SPEI值，采用游程理论提取出干旱的强度及历时两个维度信息，在此基础上采用多种方法分析对比这两个维度的趋势、突变、周期、频率等多种特征，最后又以干旱强度及历时共同作为致灾因子并结合社会经济发展数据评估广西地区农业旱灾风险，以期为广西乃至西南地区的干旱监测、预警等工作提供科学参考。

1 研究数据与方法

1.1 数 据

文中主要用到广西地区1961—2016年降雨、气温数据，以及2015—2016年的广西地区社会发展统计数据。其中降雨数据和气温数据分别来源于中国气象数据网(http:∥data.cma.cn/)的《中国地面降水月值0.5°×0.5°格点数据集(V2.0)》和《中国地面气温月值0.5°×0.5°格点数据集(V2.0)》。两数据集均基于2 472个国家气象站地面观测数据通过薄盘样条法插值而来，数据集经严格交叉验证、误差分析，质量状况良好[30]。研究中通过Python从数据集中批量提取广西地区0.5°×0.5°格点(87个)的1961—2016年逐月降水和气温数据用于特征及风险分析(图1)。而计算各旱灾风险评估指标所需的广西社会发展统计数据，包括耕地面积、行政区面积、播种面积、农业生产总值、耕地面积等(表1)则主要来源于2015—2016年广西统计年鉴(http:∥tjj.gxzf.gov.cn/)。

图1 广西0.5度格点分布

1.2 干旱监测指数

标准化蒸散发指数(standardized precipitation evapotranspiration index,SPEI)是由Vicenteserrano等[15]提出的基于水量平衡原理的干旱监测指数，其主要通过计算净降水量偏离平均状态的程度来描述干旱强度[8]。其主要计算步骤包括[1]：(1) 计算区域潜在蒸散发。(使用Thornthwaite 公式或Penman-Monteith 公式);(2) 计算逐月净降水量。(即降雨与潜在蒸散发的差值);(3) 计算不同时间尺度的水分盈亏累积序列，当选取时间尺度为1月时，计算结果即为SPEI1，同理可得SPEI3，SPEI6，SPEI12;(4) 对水分盈亏累积序列进行正态化(采用log-logistic函数);(5) 进行标准正态分布转换，计算得SPEI值。SPEI的干旱强度分级表与SPI相同，可参考相关文献。

表1 旱灾风险评估指标

1.3 干旱特征分析

1.3.1 干旱多维度信息提取方法 研究中采用游程理论提取干旱强度及干旱历时两个维度信息，其基本操作步骤如下[16]：

(1) 确定R0，R1，R2。本研究中根据侯陈瑶等[17]以及芦佳玉等[18]的研究经验，选取R0=0，R1=-0.3，R2=-0.5。

(2) 确定潜在干旱事件。若某月SPEI1≤R1,则将该月视作潜在干旱月，对于连续潜在干旱月视作一次潜在干旱事件。如图2中有a，b，c1，c2，d这4次潜在干旱事件。

(3) 从潜在事件中筛选出干旱事件。对于历时仅有1月的潜在干旱事件，若其SPEI≥R2则将其剔除，否则视为一次干旱事件；若两次潜在干旱事件间隔一个月，且该月SPEI≤R0，两者合并为一次干旱事件；其余历时大于一个月的潜在事件均视为一次干旱事件。如图2中a被剔除，b为干旱事件，c1和c2合并为一次干旱，d为一次干旱事件。

图2 基于游程理论的干旱事件识别及干旱特征变量的确定

1.3.2 干旱时间序列特征分析方法 为较为清晰地分析干旱强度以及历时序列的各种特征，研究中采用线性倾向率[19]描述干旱强度及历时的年变化量，采用Mann-kendall检验干旱强度及历时变化趋势的显著程度，采用Hurst方法[20]评估两者变化趋势的未来延续性，采用Pettitt算法[21]检验两者的变异性，采用EEMD方法[22]分析两者的周期性。以上方法均是常用的水文气象时间序列分析方法，其详细计算过程可参考相关文献。

1.3.3 干旱频率分析方法 Copula函数是常用的多维干旱特征联合分布函数，其定义域为[0,1]，其不要求变量同分布，可以将符合不同形式边缘分布的干旱特征联合起来[17,23]。由Copula建立的联合分布函数形式如下[23]：

F(d,s)=P(D≤d,S≤s)=C[FD(d),FS(s)]

(1)

式中：D为干旱历时；S为强度，FD(d)，FS(s)分别为历时、强度的边缘分布函数；C为联合分布。

其具体构建步骤如下：

(1) 求出干旱历时及强度的边缘分布函数。采用Weibull，Gamma，Normal，Birnbaum-Saunders，Logistic，Log-Normal，Poisson，Exponential这8种分布函数[23]分别拟合基于游程理论提取的干旱历时及干旱强度，并采用Kolmogorov-Smirnov(KS)检验[24]拟合效果的显著性，最后根据均方差(OLS)最小原则[25]从通过显著性检验的备选函数中找到干旱历时及强度的最优边缘分布。

(2) 选取Normal，T，Frank，Clayton，Gumbel共5种常用Copula函数[26]作为干旱强度及历时的备选联合分布函数，采用极大似然法计算出每种Copula函数的参数，并采用Bootstrap Cramér-von Mises(CM)方法[24]进行假设检验，最后从通过检验的Copula函数中选择均方差(OLS)最小的作为最优Copula。

1.4 农业旱灾风险评估

1.4.1 旱灾风险评估指标选取 研究中采用基于旱灾损失风险构成要素的风险评估方法，即从灾害的危险性(H)、暴露性(E)、敏感性(V)以及抗旱能力(RE)等方面选取一定的指标进行旱灾风险评估[27]。并根据计算公式：

(2)

计算出区域旱灾风险值。对于危险性等各个维度指标的选取则从现有研究经验[27-28]以及资料的易获取性两个方面考虑，最终选取的指标见表1。每个功能层指标与权重乘积的加和即为该功能层的值。

1.4.2 旱灾风险评估指标权重确定方法 对于旱灾风险评估中各指标权重的确定则采用投影寻踪模型。投影寻踪法主要原理是将高维数据投影到低维空间，然后建立投影指标函数，当函数达到最优时，其对应的即为最佳投影向量[29]，即可用于指标权重的确定。对于投影函数的最优解求解则常采用遗传算法[30]。

2 结果与分析

2.1 干旱时间变化特征

2.1.1 干旱序列趋势变化特征 采用游程理论从广西SPEI1序列中提取出干旱强度以及干旱历时，并统计出整个区域每年干旱事件的平均强度以及平均历时(图3)，并计算出两序列的Cv值，倾向率、Mann-kendall趋势检验Z值，以及Hurst值(表2)。由图3可知，1961—2016年，广西地区干旱强度与历时的波动均较为平缓，两者Cv值均在0.5以下。在变化趋势上，干旱强度与干旱历时则呈现出较大的差异性，1961—2016年，干旱强度呈现出微弱的下降趋势，其下降率小于-0.001/a，而干旱历时则呈现出微弱的上升趋势，其上升率为0.007月/a，不过两者变化的趋势均不显著，其Mann-kendall检验值均小于1.28，未通过90%显著性检验。从Hurst值来看，干旱历时的Hurst值为0.729，有较大的可能会继续保持微弱上升趋势，而干旱强度的Hurst仅有0.591，其变化趋势的未来延续性不显著。

图3 广西1961-2016年干旱强度、历时序列变化曲线

表2 广西1961-2016年干旱强度及历时序列特征值

2.1.2 干旱序列突变及周期变化特征 为分析1961—2016年广西干旱的突变特征，采用Pettitt算法，分别检测干旱强度以及干旱历时序列的突变点。由图4可知，广西地区1961—2016年的干旱强度序列以及干旱历时序列均存在突变点，前者突变点在2007年，后者突变点在1985年，两突变点均通过了95%显著性检验。突变后的平均干旱强度相比与突变前增加了12.26%，而平均干旱历时则相对于突变前则增加了14.56%。

图4 广西1961-2016年干旱强度及历时序列突变检验

为分析1961—2016年广西干旱的周期性特征，采用EEMD算法分别对干旱强度以及干旱历时进行分解。如图5可知，干旱强度和历时序列均可被分解为4个IMF分量，以及一个Res余量。各分量间差异显著，但总体呈现出由IMF1到IMF4，振动频率以及振动幅度均逐渐减小的变化特性。而Res余量则体现出干旱强度以及干旱历时在不同时期的变化趋势，干旱强度在1961—1991年呈现出逐渐下降的趋势，随后在1991—2016年又呈现出逐渐上升的变化趋势；而干旱历时序列则是在1961—1981年逐渐下降，1981—2001年逐渐上升，2001—2016年又逐渐下降。为确定两序列主要周期，对各序列IMF分量进行显著性检验[31]。由图6可知，干旱强度及干旱历时序列各分量中仅有IMF1分量通过显著性检验，其中前者比后者更加显著，干旱强度IMF1落在90%置信区以上，而干旱历时IMF1落在80%置信区以上，其余分量中，干旱强度的IMF3分量，以及干旱历时的IMF3，IMF4分量也均接近80%显著性检验线，但均未通过检验。对强度和历时的IMF1序列进行Hilbert-Huang变换即可得到两者对应的中心频率分别为0.30，0.32，频率的倒数即为周期，则可得两者对应的主周期为3.38，3.12 a。

图5 广西1961-2016年干旱强度、历时序列EEMD分解

图6 干旱强度、历时各IMF分量显著性检验

2.1.3 干旱频率特征 为分析广西地区1961—2016年干旱的单变量频率特征，采用1.3.3中方法拟合出干旱历时及强度的分布函数。最终选取Log-Normal(μ=0.510，σ=0.564)作为干旱历时的最优分布，选取Birnbaum-Saunders(β=1.598，γ=0.656)作为干旱强度的最优分布，两者经验概率分布与理论概率分布曲线如图7所示，Log-Normal和Birnbaum-Saunders分布对干旱历时及干旱强度的拟合较好，两者经验分布曲线与理论概率分布曲线基本重合。广西地区大多数干旱事件的历时都在7个月以内，历时在1～4月的干旱事件最多(超过80%)；强度大多小于8，强度值在0～4的干旱事件超过80%。

图7 干旱历时及强度的经验分布和理论分布

为分析基于干旱强度及历时的双变量干旱频率特征，采用1.3.3中方法最终选取Gumbel Copula(θ=3.185)函数建立干旱强度及历时的二维联合分布。图8为基于干旱历时及强度的联合概率分布图，由图8可查得不同干旱强度及历时下的干旱频率，且由图可知广西地区短历时低强度干旱居多，极少有干旱事件强度大于12且历时大于12个月，但干旱历时不大于4个月且干旱强度不大于4的干旱事件约有60%；随着干旱历时和干旱强度值的不断增大，二者的联合累计概率值也不断增大，但联合累计概率值越大，其增长速度也越慢，这说明大多数干旱事件以短历时低强度为主，而强度越高且历时越长的干旱事件越少。

2.2 干旱空间分布特征

使用游程理论识别出每个格点1961—2016年所有干旱事件的强度及历时，并计算出其平均强度及平均历时空间插值。由图9可知，广西干旱平均强度为1.20～2.05，其中桂东南的贵港、玉林、北海等地较高，在1.65以上，其次在桂北的桂林、柳州以及桂西河池西部、崇左等地为1.50～1.60，最低的则是桂北的河池东部以及桂南的钦州、防城港等地在1.50以下。广西干旱历时则集中在4.4～4.9个月，且呈现出，桂北的桂林北部—河池，以及桂东南的梧州—防城港两个带状区域较低(4.4～4.7个月)，贺州—崇左的带状区域较高(4.7～4.9个月)。由干旱强度以及干旱历时的空间分布对比则可以看出，干旱强度及干旱历时的空间分布在整体上没有表现出显著的相关性，但在部分地区，例如桂北的河池东部以及桂南的钦州西部两者均偏小，而在桂中的来宾、贵港等地两者均偏大。

计算出每个格点年干旱强度及干旱历时序列的年变化率，并计算出每个序列的Mann-kendall趋势检验值和Hurst值(图10)。由图10可知，干旱强度及历时呈现出显著的空间分布特征。从干旱强度来看，桂北呈现出降低趋势，而桂南则呈现出上升趋势，但仅有桂北及桂东南极小部分区域的变化趋势较为显著，其余大部分区域的Mann-kendall检验值均为-1.28～1.28(显著性检验小于90%)，而从干旱历时来看，全区大部均呈现出上升趋势，仅有桂北小部分区域呈现出下降趋势，但除桂北小部分区域的增长趋势较为显著，其余地区均未通过90%显著性检验。在趋势未来延续性方面，干旱历时序列在全区域大部分地区的Hurst值均在0.50左右，未来延续性并不明确，而干旱强度序列的Husrt值在桂西以及桂南等小部分区域在0.50附近，其余区域Hurst值均在0.60以上，呈现出一定的趋势延续性，尤其是在桂中等区域，其Hurst值在0.70以上，其趋势延续性较强。

2.3 农业旱灾风险

计算出风险评估各指标值，进行归一化处理以后建立投影寻踪模型，其中各市SPEI强度及历时由各市区域内格点的强度及历时的多年算术平均求得。根据1.3.2内容建立投影寻踪模型并得到最佳投影方向A={0.152,0.292,0.401,0.189,0.470,0.145,0.341,0.098,0.282,0.404}，A中每一个元素顺序对应表1中a1—a12。通过投影向量可以计算出各市的危险性、暴露性、敏感性、抗旱能力值，并根据旱灾风险综合计算公式可得各市干旱风险(图11)。

图8 干旱强度及历时的联合概率分布

图9 广西干旱强度及历时空间分布

图10 广西干旱强度及历时的变化趋势空间分布

图11 广西干旱危险性、暴露性、脆弱性、抗旱能力及农业旱灾风险空间分布

由图11可知,广西地区的干旱危险性在贺州、贵港以及百色等地较高，在桂南的钦州、北海、防城港等地较低，其余地区则相差不大；而暴露性则在桂中地区较高，尤其是崇左、来宾等地，在桂西的百色以及桂东的贺州、梧州等地偏低，这与屈艳萍等[32]对全国旱灾暴露性研究结果较为相似(崇左等地较高，桂东地区偏低)；脆弱性则呈现出由桂东向桂西逐渐增加的变化趋势，这与康永辉等[11]对桂西北地区的脆弱性研究结果较为相似(桂西脆弱性较高)；抗旱能力则是在桂北的桂林以及桂东等区域较低，桂中的次之，桂西等地则较高，但最高的则是桂中的来宾。综合各项指标计算得广西地区各市干旱风险值，其中干旱风险最高的是桂西南的崇左，以及桂中的来宾，其余各区域中，桂西北的旱灾风险要明显高于桂东南，这与张强等[2]在研究华南地区旱灾风险的结论有一定的相同之处(桂西较高，桂东北及桂东较低)。

3 结 论

(1) 1961—2016年，广西干旱强度呈现出微弱的下降趋势，而历时则微弱上升，但两者的变化趋势均不显著，不过后者表现出更强的未来延续性；干旱强度序列在2007年存在突变点，而历时则在1985年；干旱频次和历时的周期均在3.20 a左右。

(2) Log-Normal和Birnbaum-Saunders函数可以分别较好地拟合干旱历时和干旱强度的概率分布；而Gumbel Copula函数则为两者最优的联合概率分布。就单变量来说，超过80%的干旱历时在1～4个月，超过80%的干旱强度在0～4；就双变量来说，短历时低强度干旱居多，约有60%的干旱事件强度及历时均小于4。

(3) 干旱强度在桂东南较大，桂西次之，桂南最低；干旱历时则在贺州—来宾—崇左的带状区域内较高，在桂东南较低。就变化趋势来说，强度在桂北以减小为主，在桂南以增长为主；历时在全区大部分区域均以增长为主，且桂中较高；但全区大部分区域强度及历时的变化趋势均不显著。

(4) 广西地区干旱风险最高的是桂东南的崇左，以及桂中的来宾，其余各区域中，桂西北区域的干旱风险要明显高于桂东南。