考虑融雪量的多变量综合干旱指数的构建与评估

摘 要：

【目的】为解决大多数干旱评估中雪因素常被忽略的问题，需要开发适用于高寒融雪流域的干旱指数。【方法】根据潜在蒸散量、融雪量、降雨量计算水分增益/损失，然后利用Copula函数计算水分增加/损失与土壤水分的联合分布概率，从而构建了一个基于降雨、融雪量、潜在蒸散量和土壤水分新的标准化融雪综合干旱指数(SSCDI),并比较和评估了3月和6月时间尺度的SSCDI在黄河源区的干旱监测性能。【结果】结果表明，Fischer-Hinzmann Copula是构建SSCDI的最佳连接函数；新构建的SSCDI与SPEI和MSDI相比，同SSI的拟合程度更好(R2高于0.7);SSCDI比SPEI和MSDI更符合观测到的干旱变化。【结论】特别是在高寒流域，对雪因素的考虑提高了SSCDI的干旱监测能力，并且SSCDI能够同时监测出SPEI和MSDI所监测到的干旱情况，总体而言，SSCDI可用于识别和监测高寒区不同时间尺度的气象和农业综合干旱。

关键词：

干旱指数;融雪量;Copula函数;黄河源区;

作者简介：

胡莹莹(1990—),女，博士研究生，主要从事生态水文研究。

*周毓彦(1991—),男，高级工程师，博士，主要从事水文水资源研究。

基金：

2022年度科技智库青年人才计划(20220615ZZ07110156);

青海省中央引导地方科技发展资金项目(2022ZY020);

第二次青藏高原综合科学考察研究(2019QZKK0207-02);

国家自然青年科学基金项目(51909275);

中国水利水电科学研究院基本科研业务费专项项目(WR110145B0052021);

流域水循环模拟与调控国家重点实验室开放研究基金项目(IWHR-SKL-KF202204);

引用：

胡莹莹, 周毓彦, 鲁帆, 等. 考虑融雪量的多变量综合干旱指数的构建与评估[J]. 水利水电技术(中英文), 2023, 54(8): 79- 90.

HU Yingying, ZHOU Yuyan, LU Fan , et al. Development and evaluation of a multivariate composite drought index accounting for snow [J]. Water Resources and Hydropower Engineering, 2023, 54(8): 79- 90.

0 引 言

干旱是常见的自然灾害之一，由于其持续时间长、影响范围广，往往造成巨大的经济损失。随着全球平均气温升高，极端天气事件发生频率增加，导致气候系统的稳定性下降，全球大部分地区干旱化趋势日益显著，随之而来的干旱灾害增加，严重威胁到农业和生态安全。2021年政府间气候变化专门委员会 (IPCC)第六次评估报告(AR6)指出，气候变化预测表明许多地区将经历更频繁和更严重的干旱。因此，有效的干旱监测、预测和早期预警对于决策和减轻负面影响的缓解措施至关重要。温度和大气蒸发需求是干旱变化的重要驱动因素，它们的增加预计会减少未来的降雪量,降水将更有可能以雨而不是雪的形式落下，这将对以融雪为主的流域水文产生重大影响。由雪变化引起的长期地表水失衡和水分亏缺可能导致更频繁的水文和农业干旱，进一步威胁水安全。特别是对于中高海拔地区，未来降雪量的减少可能会导致更严重的干旱。考虑到雪在干旱演变中的重要性和复杂性，有必要将雪因素纳入干旱监测。

因为干旱的定义复杂，许多干旱指数被逐渐发展出来表征干旱。干旱指数是研究干旱的基础和关键，只有合理的干旱指数才能准确监测和量化干旱，反映出干旱规律和趋势的变化。然而，仍然很难建立一个能够监测和识别各种干旱的通用干旱指数。AR6 将干旱分为气象干旱、农业干旱、水文干旱和生态干旱。干旱类型的划分是相对的，不同类型之间的干旱在一定条件下可以相互转化。为了量化不同类型的干旱，过去已做出一些努力来开发各种干旱指标，依靠这些干旱指标，干旱评估已经有了长足的进步：标准化径流指数(SRI)、归一化植被指数(NDVI)、植被状况指数(VCI)、标准化土壤湿度指数(SSMI)、标准化降水指数(SPI)、帕尔默干旱指数(PDSI)、标准降水蒸散指数(SPEI)已经在世界范围内广泛使用。SPI和SPEI通常用于表征气象干旱，SRI用于表征水文干旱，NDVI和VCI用于表征生态干旱，SSMI用于量化农业干旱，PDSI用于评估基于降水、蒸散量和径流的综合干旱。然而，他们通常忽略了雪对干旱变异性的影响。这些方法忽略了雪因素可能会降低它们监测和识别干旱的能力。因此，干旱指数需要考虑积雪和/或融雪。虽然STAUDINGER等人提出了标准化融雪和雨水指数，该指数考虑了影响干旱的雨水和融雪赤字，但它只评估了水文干旱，而且是利用温度和降水来模拟雪，并不是直接采用雪数据；ZHANG最近建立了结合雪动力学的标准化水分异常指数(SZIsnow),但它具有计算复杂和困难等局限性。

HAO和AGHAKOUCHAK等提出了一种多变量多指标干旱模拟方法—多变量标准化干旱指数(MSDI),该方法利用降水和土壤水分的联合分布函数来综合这两个变量的干旱信息。Javed等比较了SPI、标准土壤水分指数(SSI)和 MSDI,以确定最适合中国的农业干旱指数，结果表明MSDI在监测中国大陆农业干旱方面表现更好。然而，由于MSDI没有考虑蒸散量和雪因素，因此其在干旱监测中存在一定的缺陷。VICENTE SERRANO等人在SPI的基础上开发了SPEI,由于它物理意义和简单的计算方法、多标量和多类别特征，SPEI已经广泛应用于表征气候变化下干旱监测。SPEI考虑了供水和需求，并利用降水和潜在蒸散量之间的差异来评估水分亏缺或盈余。然而，潜在蒸散量并不是影响需水量和干旱条件的唯一因素，因为土壤水分在估算需水量方面也起着重要作用。

黄河源区位于青藏高原和黄土高原的过渡地带，是典型的旱地高寒草地生态系统，对气候变化响应最为敏感，也对生态环境变化最为脆弱。近年来全球变暖的趋势对黄河源区生态环境产生了重要影响。因此，本研究以黄河源区为研究对象，综合考虑了MSDI和SPEI的局限性，主要目标是： (1) 将融雪量纳入干旱评估，利用Copula函数建立一种新的干旱指标SSCDI,它基于MSDI 的联合概率方法，同时具有MSDI和SPEI的优点，并且可用于多时间尺度干旱评估；(2)评价新构建的干旱指数SSCDI在黄河源区干旱监测中的性能。

1 研究区概况

黄河源区是指为以唐乃亥水文站为流域出口的汇水区域。如图1所示，黄河源区位于东经95°50′—103°30′,北纬32°10′ —36°15′,面积为1.22×105 km2,海拔2 569～6 072 m。黄河源区流域年均径流量为 2.05×1011 m3,占全流域的35%,是黄河流域重要的水源产流区。草地是黄河源区最重要的土地覆被，占整个源区的80%,主要类型为高寒草甸、高寒沼泽草甸和高寒草原。黄河源区属典型的高原大陆性气候，兼具高寒气候和干旱气候的特点。黄河源区长期年平均降水量为 502.72 mm, 时空变化极不均匀,其中西北部较为干燥，年降水量不足350 mm, 而东南部气候较为湿润，年降水量约750 mm。气候变暖对黄河源区产生了重大的生态水文影响。特别是黄河源区的平均气温自1980年代以来以0.48 ℃/10 a 的速度升高。近40年来，黄河源区降雪量以2.43 mm/10 a的速率增加，降雪率(降雪量占降水量的比例)多年平均约为26%,黄河源区积雪期集中在11月到次年4月,融雪期主要为每年的4—6月和9—10月。

图1 研究区位置

2 研究方法和数据

2.1 SPEI、SSI 和 MSDI 的计算

SPEI是一种广泛使用的气象干旱指数，用于说明大气湿度与平均状态的偏差。SPEI计算的基本步骤为：(1)计算降水量与潜在蒸散量之差；(2)通过拟合一定时间尺度下差异的累积值来求出概率分布；(3)归一化累积概率密度函数以获得 SPEI 值，本研究采用Python计算gamma分布的SPEI。SSI通过充分考虑土壤水分作为输入来描述干旱状况,SSI 是根据正常降水百分比和 Palmer Z指数的概念以及SPI的统计结构开发的。MSDI是HAO等开发的基于SPI的指数的扩展，进一步将SPI扩展为基于降水和土壤水分的双变量模型。MSDI可以根据降水和土壤水分的联合概率计算得出。

2.2 SSCDI的构建

2.2.1 SSCDI的构建框架

开发和评估SSCDI的框架如图2所示。首先对NASA Global Land Data Assimilation System (GLDAS)数据进行预处理，以提取降水、降雨、潜在蒸散量、融雪和土壤水分。其次，根据潜在蒸散量、融雪量和降雨计算水分增减，然后用Copula函数计算水分增减与土壤水分之间的联合分布概率。在此基础上，本文选择了5个Copula函数来构建联合分布函数，通过赤池信息准则(AIC)、贝叶斯信息准则(BIC)、均方根误差 (RMSE) 和纳什-萨特克利夫效率 (NSE)对每个候选Copula函数的性能进行评估，以选择最佳拟合的分布函数。最后，通过对联合分布概率的归一化得到SSCDI。同时计算了SPEI、MSDI、SSI干旱指数，并将其用于SSCDI的可靠性评估。此外，根据从黄河源区实际干旱统计资料中选取的典型干旱年份，对SPEI、MSDI和SSCDI监测的干旱结果进行比较，验证SSCDI在干旱监测中的适用性。

图2 SSCDI的构建流程

2.2.2 SSCDI的构建方法

本研究首先根据潜在蒸散量、融雪量、降雨量、土壤水分计算水分增益/损失，然后采用Copula函数计算水分增加/损失与土壤水分的联合分布概率，Copula函数可以导出两个或多个变量的联合分布函数而不管它们的原始边际分布如何。计算公式为

式中，C 是 Copula函数；F(Di)和G(Si)分别是随机变量Di和Si的边际累积分布函数；p表示两个变量的联合分布概率；d和s表示需要满足的特定值；Di为水分增益/损失；Si是土壤水分；Mi是融雪量；Ri是降雨量；PETi是潜在蒸散量；i代表月。

式中，φ为标准正态分布函数。通过将联合概率标准化，可以将其与不同的干旱指数进行交叉对比，与MSDI相似，它可以用于监测不同的时间尺度(1个月、2个月、3个月、6个月、9个月等),不同的时间尺度可以反映不同的干旱类型，3月的尺度更能准确监测到黄河源区的季节性干旱事件，6月的时间尺度能较好地反映出黄河源地区干季和雨季的变化。因此，在本研究中使用3月和6月尺度的SSCDI表征黄河源区的干旱状况。如表1所列，SSCDI负值代表干旱，而正值表示潮湿。

2.2.3 Copula模型

Copula近年来大量应用于水文研究中。利用边缘分布和相关框架，Copula函数模型可以构建多维联合分布。Copula函数能够将多个随机变量与任何给定的依赖结构耦合，并构造一个包含独立指数的所有关键特征的因变量，Copula 定义为

式中，F1(x1),…,Fn(xn)是指F(x1,…,xn)的单变量边际分布；C(u1,…,un)是Copula的密度函数。

Copula函数种类较多，一般来说可分为四种类型，包括阿基米德型(Archimedean)、椭圆型(Elliptical)、极值型(Extreme value)和其他类型。Copula用于描述随机变量之间的依赖结构。所有的Copula都可以捕捉到几种依赖结构，如正、负的对称性、不对称性和尾部依赖性。然而单参数Copula只能捕捉一种类型的依赖关系, 而双变量，如BB1(Clayton-Gumbel),BB8(Joe-Frank)可以捕获随机变量之间的一个以上的依赖类型。本文选取5种Copula函数来构造联合分布函数，并采用最佳拟合表示的Copula 进行最终评估。5个不同 Copula 函数和参数范围如表2所列。

2.2.4 拟合度优选方法

我们使用了由SADEGH等人开发的MATLAB工具箱—多元Copula分析工具箱(MvCAT)来建模二元概率分布。MvCAT包括各种不同复杂度的Copula族，并使用基于残差的高斯似然函数的贝叶斯框架来推断Copula参数。MvCAT的主要优点是引入了一种混合进化马尔可夫链蒙特卡罗方法，用于耦合参数的后验分布的数值估计，并相对不确定性进行评估。在本研究中，使用AIC、BIC、RMSE)和NSE来评估不同Copula函数的性能。

2.3 相关分析

在构造多变量的联合概率分布之前，首先需要对变量间的相关性结构进行探究，本文采用Pearson、Kendall、Spearman三种相关系数来度量变量间的相关性。Pearson相关系数描述变量间的线性关系，Kendall秩相关系数和Spearman秩相关系数描述变量间的非线性关系。设(xi, yi)为取自随机变量(X,Y)的样本，三种相关系数的计算公式为

(1)Pearson相关系数 ρp

式中，x¯,y¯分别为样本xi和yi的均值。

(2)Kendall 秩相关系数 τ

式中，n代表数据系列的总长度。

(3)Spearman′s 相关系数 ρs

式中，Ri和Qi分别是样本值； xi, yi在序列中取值等级的顺序；R,Q分别是是样本Ri和Qi的平均等级。

2.4 GLDAS陆地表面数据集

GLDAS-2是通过使用全球气象强迫数据集作为输入，驱动Noah陆面过程模型(LSM)来模拟全球陆面过程而获得各地表数据。GLDAS-2能够较好地捕捉到土壤水分的时间变化，与观测土壤水分的相关系数在0.5以上,同时，GLDAS-2中的融雪量数据也是根据Noah陆面模式的算法计算出来的，该数据为推进与雪相关的水文效应研究提供了独特的机会，因为该数据集提供了高质量、长期、全球覆盖和多变量的地表和地表条件估计,它准确再现雪过程的能力是我们干旱指数发展的重要基石。本研究中采用1948年1月至2014年12月 的GLDAS-2.0 Noah LSM 月度数据和2000年1月至2021年12月的GLDAS-2.1 Noah LSM月度数据，它们的分辨率均为0.25度。GLDAS-2.0和GLDAS-2.1数据集都可以通过戈达德地球科学数据和信息服务中心获得 (http: //disc.sci.gsfc.nasa.gov/hydrology/data-holdings)。分析中使用的各种参数包括降水量、降雨量、潜在蒸散量、融雪量和土壤水分， 对其进行距离加权平均重映射进行插值，插值精度为0.05°,该插值方法是对所有输入字段的四个最近的邻居值进行距离加权平均重映射。GLDAS提供四个土层的水分数据，分别为0～10 cm、10～40 cm、40～100 cm和100～200 cm, 在本研究中我们采用顶层(即0～10 cm)的土壤水分数据。

2.5 气象灾害统计数据

青海省历年干旱灾害情况数据来源于《1945—1995中国灾情报告》《中国气象灾害统计年鉴》《气象干旱年鉴》和《中国气象灾害大典·青海卷》。这些数据用于反映现实世界干旱的发生情况，因此用于验证本研究中提出的SSCDI 在黄河源区干旱监测中的性能。《中国气象灾害大典-青海卷》数据时间范围为1950—2000年，下载自国家青藏高原科学数据中心(https: //data.tpdc.ac.cn/zh-hans/data/197ce7c7-31c4-445a-b3b1-953d91b7c54a)。其余三种资料均下载自中国知网，其中《中国气象灾害年鉴》时间范围为2004—2020年，《气象干旱年鉴》时间范围为2011—2017年。

3 研究结果与分析

3.1 最佳拟合Copula优选结果

验证变量之间的相依性结构是正确进行Copula分析的关键，如表3所列，3月和6月尺度的Di 和 Si的Pearson, Kendall和Spearman的相关系数值高达0.79,0.63和0.83,这表明Di 和 Si具有较强的相关关系，同时可以看出，3月尺度的三个相关系数均大于6月尺度的相关系数。

5种Copula函数用于推导Di 和 Si的联合概率分布的拟合优度检验结果如表4所列。从表中可以看出，5种候选Copula函数的3月和6月尺度RMSE,NSE,AIC,BIC检验结果非常接近，各时间尺度的平均值分别为0.22,0.99,-6713,-6704,表明所选Copula函数对变量的拟合效果都较好，具体来看，Fischer-Hinzmann Copula 在各时间尺度上的RMSE值、BIC值和AIC值最小，且NSE值最大，因此Fischer-Hinzmann Copula 被选为最佳拟合分布，故而本研究构建SSCDI时使用Fischer-Hinzmann Copula作为连接函数。

3.2 SSCDI与其他干旱指数的对比验证

为了解SSCDI在量化干旱事件方面的表现，我们将其与文献中其他已经建立的干旱指数进行比较。采用了SSI来初步评估SSCDI,过去虽然已经开发了多种能够衡量农业干旱的指数，例如PALMER开发的作物湿度指数(CMI),CAMMALLERI等人开发的土壤水分干旱严重程度指数(DSI)等等，然而这些指数只能采用固定时间尺度进行分析，而SSI可以在多时间尺度上衡量农业干旱，因此可以在不同尺度上与SSCDI进行对比。基于Copula函数的优选结果，以3月和6月尺度为例，计算了黄河源区SPEI、SSI、MSDI和SSCDI的指标值，并记为SPEI-3、SSI-3、MSDI-3和SSCDI-3,以及 SPEI-6、SSI-6、MSDI-6和SSCDI-6。计算一定时间尺度下干旱指数之间的线性回归拟合程度的目的是评价在相同条件下，哪种干旱指数在干旱识别方面表现最好，不同时间尺度的SSI与SPEI、MSDI和SSCDI的散点关系图如图3所示，其中3-月SSI-SSCDI的R2值为0.711,高于SSI-SPEI(0.525)和SSI-MSDI(0.630),6月尺度SSI-SSCDI的R2值为0.732,同样高于SSI-SPEI(0.694)和SSI-MSDI(0.571),并且两种时间尺度上SSI-SSCDI的RMSE值都最小。因此，SSCDI在各时间尺度上和空间位置上的误差是最低的，表明SSCDI具有更好的干旱监测性能。

图3 SPEI、MSDI、SSCDI与SSI在3月和6月时间尺度上的散点图

贵南县和玛多县SPEI、MSDI和SSCDI在3月尺度的时间序列如图4和图5所示。其中，玛多县平均海拔4 300 m, 多年平均融雪量15.83 mm, 而贵南县平均海拔3 100 m, 多年平均融雪量3.16 mm, 远远低于玛多县，因此本文称玛多县为融雪区域，贵南县为非融雪区域，将两者进行对比分析。如图4所示，非融雪区域的各干旱指数走势和变化规律具有一致性，SSCDI-3折线与大多数MSDI-3折线重叠，SPEI-3、MSDI-3和SSCDI-3的干旱开始和结束时间相似，这意味着SPEI-3、MSDI-3和SSCDI-3的干旱监测性能相似，原因在于非融雪区域雪对干旱表征的影响较小。融雪区域中三种干旱指数显示的干湿状况并不完全一致(见图5),根据《气象干旱年鉴》中2016年中国各省干旱事件统计表可知，2016年8月下旬至9月上旬青海省发生了特级干旱，如图5(b)所示，SSCDI指数值小于-2,监测出了该特干旱事件，而SPEI和MSDI仅监测出了中等干旱。说明融雪量的加入对干旱的评估具有一定的改善作用。

图4 贵南县1948—2021年 SPEI-3、 MSDI-3和 SSCDI-3时间序列干旱指数比较

图5 玛多县1948—2021年 SPEI-3、 MSDI-3和 SSCDI-3时间序列干旱指数的比较

3.3 典型干旱年份的比较分析

为了进一步测试 SSCDI 的干旱监测性能，选择了典型的干旱年份作为案例研究。据《中国气象灾害大典-青海卷》和《中国气象灾害年鉴》记载，青海省1960—1961年全省发生持续干旱。2013年3—4月，青海省遭遇严重干旱，全省重旱、特旱面积居全国之首。因此，本研究中，以1960年和2013年为案例，评估SPEI、MSDI和SSCDI在3月和6月时间尺度上的干旱监测性能，因为研究区域在这些年份具有典型的干旱状态。

1960年7月和2013年4月SPEI、MSDI和SSCDI在3月和6月时间尺度上的干旱监测结果如图6和图7所示。总体而言，与气象干旱相比，农业干旱影响的区域更大。SSCDI成功监测出了气象干旱和农业干旱，此外，SSCDI的干旱空间格局与MSDI相似。图中蓝色的干旱等级为D1,表示正常或偏湿润，颜色越深代表干旱等级越高。从图6中紫色圆圈标记的区域可以看出，1960年对于MSDI-3未能监测出干旱而SPEI-3监测出D3干旱的区域，SSCDI-3同样监测出了D3干旱，2013年对于MSDI-3未能监测到干旱而SPEI-3监测出D2干旱的区域，SSCDI-3同样监测出了D2干旱。从图7中紫色圆圈标记的区域可以看出，1960年对于SPEI-6未能监测出干旱，而MSDI-6监测出干旱的区域，SSCDI-6同样监测出了干旱情况，2013年对于MSDI-6未能监测出干旱，而SPEI-6监测出D2干旱的区域，SSCDI-6监测出了D3干旱。

图6 1960年7月和2013年4月SPEI、MSDI 和SSCDI 在3个月时间尺度上的空间格局

图7 1960年7月和2013年4月SPEI、MSDI 和SSCDI 在6个月时间尺度上的空间格局

因此可以得出，SSCDI能够同时监测出SPEI或MSDI监测到的干旱，且SSCDI量化的干旱等级大于或等于MSDI,可能原因在于MSDI没有考虑蒸散发因素。由此我们可以得出结论，SSCDI在黄河源区干旱监测中表现良好。

4 结 论

基于Copula理论，本文将融雪量纳入干旱评估体系，提出了一种新的综合干旱指数SSCDI。虽然本研究是基于GLDAS数据分析的，但SSCDI同样适用于站点数据，它是多变量的，可以在多个时间尺度上进行计算，并且计算过程简单清晰。该指数综合考虑了降雨、融雪量、潜在蒸散量和土壤水分，为干旱量化(特别是高纬度地区的干旱量化)带来了更广泛的视角。

(1)SSCDI是在SPEI和MSDI的基础上制定的，并且对5种Copula函数进行拟合度优选，得出了构建SSCDI时使用Fischer-Hinzmann Copula作为连接函数。

(2)通过与SSI进行线性回归及误差分析，结果显示在3月和6月时间尺度上SSCDI与SSI拟合度更高(R2分别达到了0.711和0.732)。SSCDI在融雪较多的区域表现更佳，而在非融雪区域与SPEI和MSDI具有较高的一致性，因此，一个流域受雪的影响越多，就越值得用SSCDI来评价其干旱情况。

(3)在黄河源区将这三者的结果与实际统计资料的干旱情况进行了对比，验证了SSCDI在干旱监测中的适用性。以1960年7月和2013年4月为例，分析了SSCDI、SPEI、MSDI在3月和6月尺度上的干旱监测性能，结果表明，SSCDI能够同时监测出SPEI和MSDI监测到的干旱，纳入融雪量使SSCDI在高寒区能够产生更可靠的干旱监测结果。

综上所述，构建的干旱指数SSCDI在高寒地区的干旱量化方面具有良好的监测效果，适用于气象和农业干旱的综合评估。然而，尽管构建的干旱指数很新颖，但仍存在一些局限性，本文仅使用GLDAS数据对SSCDI进行了验证和分析，未来还需要利用多种遥感数据和实测数据进行研究，进一步验证SSCDI的可靠性。

水利水电技术（中英文）

水利部《水利水电技术（中英文）》杂志是中国水利水电行业的综合性技术期刊（月刊），为全国中文核心期刊，面向国内外公开发行。本刊以介绍大陆水资源的开发、利用、治理、配置、节约和保护，以及水利水电工程的勘测、设计、施工、运行管理和科学研究等方面的技术经验为主，同时也报道国外的先进技术。期刊主要栏目有:水文水资源、水工建筑、工程施工、工程基础、水力学、机电技术、泥沙研究、水环境与水生态、运行管理、试验研究、工程地质、金属结构、水利经济、水利规划、防汛抗旱、建设管理、新能源、城市水利、农村水利、水土保持、水库移民、水利现代化、国际水利等。