气候系统的综合观测和多项指标表明，全球变暖趋势仍在进一步持续。2016年，全球表面平均温度再创新高，成为有气象观测记录以来的最暖年份，比工业化前水平高出约1.1 ℃；亚洲陆地表面年平均气温是1901年以来的第2高值；中国年平均气温亦属明显的偏暖年份^[1]。这种全球增暖带来的时空上非均匀的降水强度和分布的变化，必定对气候的干湿状况产生重要影响^[2-3]。当前，气候干湿状况及演变规律一直是科学家研究的热点问题之一。

湿润指数(humidity index，H_i)作为热量与水分的综合因子，控制陆地生态系统与大气之间能量和物质交换，可以指示某一地区能量和湿度从地面到大气的转换情况，并直接影响到植被生产及其对水分需求^[4]。近年来，学者们基于H_i对中国^{[3, 5-6]}及其各区域^[7-11]干湿变化的研究取得了大量成果，并就其主导影响因子进行了广泛的讨论^{[3, 7-8, 10-11]}。青藏高原作为“第三极”和“亚洲水塔”，众多学者从干湿状况的时空分布特征、变化趋势及其影响因子等方面做了大量的研究^{[4, 12-15]}。自1990年以来青藏高原正处于加速升温阶段^[16-17]，干湿状况发生了明显变化，这种态势在进入21世纪后表现的尤为突出。羌塘高原约占青藏高原总面积1/5，是世界上海拔最高、气候条件最恶劣的高原。羌塘国家级自然保护区(以下简称自然保护区)就位于羌塘高原，这里气候异常寒冷，属高原寒带季风干旱气候，具有独特的高寒荒漠生态系统，生态环境十分脆弱。作为气候变化的指示因子和植被生态系统的影响因子，研究干湿指数具有重要的科学意义。目前有关青藏高原地表干湿状况的研究，时段大多截至21世纪初，而近10年最新变化的研究尚未见报道，尤其是针对羌塘自然保护区湿润指数的变化更是鲜有报道。鉴于此，笔者利用自然保护区边缘5个气象站点1971—2018年逐日气象资料，基于Penman-Monteith模型计算潜在蒸散和湿润指数H_i，系统地分析近48年自然保护区干湿状况的时空变化规律，并对干湿变化的成因进行讨论，明确全球气候变暖下羌塘高原的干湿分布和变化趋势，旨为认识生态环境现状与过去的差异，为羌塘国家级自然保护区生态环境保护和恢复提供依据。

1 研究区概况

自然保护区位于西藏自治区西北部，昆仑山、可可西里山以南，冈底斯山和念青唐古拉山以北，地处E 79°59′~90°26′、N 32°12′~36°29′，平均海拔5 000 m以上，被称为“世界屋脊的屋脊”。自然保护区于1993年经西藏自治区人民政府批准成立，2000年4月4日经国务院批准晋升为国家级自然保护区。保护区面积约29万8 000 km²，是我国面积最大的自然保护区，仅次于格陵兰国家公园的世界第2大陆地自然保护区，也是平均海拔最高的自然保护区。行政隶属西藏那曲市西部3县(安多、尼玛、双湖)和阿里地区北部3县(改则、日土、革吉)^[18-19]。

自然保护区属于羌塘高原湖盆高寒草原区^[20]，羌塘高原是青藏高原的主体，处于昆仑山和冈底斯—念青唐古拉山之间，东自内外流水系的分水岭，西以公珠错—革吉—多玛一线与阿里西部山地为界。整个地势南北高、中间低，高原面较完好，北部海拔4 900 m左右，南部约4 500 m。羌塘高原是世界上海拔最高的内流区，流域集水面积小，多季节性河流。该区是高原湖泊集中分布区域，是著名的高海拔湖泊群区，湖泊大多为咸水湖和盐湖，淡水湖极少。这里气候寒冷，气温的年日变化大。最暖月平均气温6~10 ℃。局地可达12 ℃；最冷月平均气温在-10 ℃以下。年降水量约100~300 mm，自东南向西北递减。冬春多大风。

2 资料与方法

自然保护区内无长时间序列观测资料的气象站，当前只能在其周边选取同属羌塘高原湖盆高寒草原区的狮泉河、改则、申扎、班戈和安多5个气象站(图 1)1971—2018年逐日地面气象资料(包括平均气压、平均气温、平均最高气温、平均最低气温、平均相对湿度、平均水汽压、平均风速、日照时间、日照比例和降水量等)。为了保证数据的可靠性，对所有站点的数据均进行了质量控制。

2.1 地表湿润指数计算

地表湿润指数关系表达式为

$ {H_i} = P/{E_0} 。$

(1)

式中：H_i为湿润指数，量纲为1；P为降水量，mm/d；E₀为潜在蒸散量，mm/d，应用1998年FAO推荐并修订的Penman-Monteith模型计算E₀^[21]。

2.2 分析方法

利用5个站气象资料，建立自然保护区4季(3—5月为春季，6—8月为夏季，9—11月为秋季，12月—翌年2月为冬季)以及年潜在蒸散、H_i和气象因子等序列，以分析自然保护区地表H_i的气候变化特征。气候基准期为1981—2010年，多年平均值采用此时段30年的平均值。年代际变化距平为各年代平均值与1981—2010年多年平均值的差，主要用于气候要素年代之间的比较。

要素的变化趋势可用下式估计：

$ Y=a_{0}+a_{1} t_{\circ} $

(2)

式中：Y为要素(如：H_i；潜在蒸散，mm；平均气温，℃；气温日较差，℃；日照时间，h；平均风速，m/s；平均相对湿度，%；平均水汽压，hPa；降水量，mm)；t为时间，a；a₀为常数项；a₁为线性趋势项，表示要素每年的变化趋势。采用F检验法对拟合的线性回归方程进行显著性检验(P<0.10、P<0.05、P<0.01和P<0.001)。

另外，采用Mann-Kendall(后称M-K)法^[22-24]对年、季H_i进行突变检验。利用主成分回归方法^[25]讨论影响H_i的主导因子。

3 结果与分析 3.1 湿润指数的时空分布特征

在多年平均状态下，采用年H_i作为干湿区域划分标准^[2]，分析结果表明：自然保护区大部分地方H_i＜0.5，属于半干旱和干旱区，其中西部狮泉河和改则H_i≤0.2，为干旱区；中部申扎与班戈0.2＜H_i≤0.5，为半干旱区；东部安多，H_i＞0.5，为半湿润区。从自然保护区H_i月际变化(图 2)来看，呈单峰型变化，其中3月份最低，为0.03；8月份最高，为0.80；1—5月降水稀少，H_i较低，气候极为干燥，属于干旱期；进入6月后，降水迅速增加，但仍小于潜在蒸散，H_i为0.45，气候处在半干旱期；7—9月降水量猛增，为全年最为集中期，占年降水量的69%，大于此同期的潜在蒸散，H_i＞0.5，气候处在半湿润期。10—12月降水明显减少，H_i＜0.5，气候处在干旱期。在季尺度上，春、冬2季均为干旱期，夏季为半湿润期，秋季为半干旱期。

3.2 湿润指数变化 3.2.1 年际变化

从保护区年、季H_i的变化来看，近48年(1971—2018年)年H_i呈显著的增加趋势(图 3a)，但增幅不大，仅为0.001/a(P<0.05)，主要表现在春、夏2季，增幅分别为0.001/a(P<0.05)和0.003/a(P<0.05)；尤其是近28年(1991—2018年)夏季H_i增幅更明显(图 3b)，增幅达0.007/a(P<0.05)。秋季和冬季H_i在过去48年里的变化趋势不明显，但近28年2季H_i表现为较明显的减小趋势(图 3c和d)，幅度分别为-0.004/a(未通过显著性检验)和-0.001/a(P<0.01)。

表 1为自然保护区各站H_i的变化趋势，在年尺度上，近48年除狮泉河H_i无变化外，其他各站呈增加趋势，增幅为0.002~0.003/a，以申扎增幅最大(P<0.01)；进入20世纪90年代后，除安多年H_i趋于减小(-0.002/a)外，其他各地H_i均呈增加趋势，尤其是改则增幅最为明显，增幅达0.004/a(P<0.01)。在季尺度上，近48年春季H_i除狮泉河无变化外，其他各地为增加趋势，增幅为0.001~0.003/a，以安多增幅最大，这种趋势在近28年表现的更明显。夏季H_i变化趋势的分布与春季相同，不过增幅略大一些，为0.002~0.005/a(申扎最大)；其中，近28年改则夏季H_i增幅突出，达0.010/a(P<0.01)。秋季，近48年大部分站点H_i呈现为弱的减少趋势，但在近28年里H_i减小趋势明显，以申扎最突出，减幅达-0.008/a(P<0.01)。冬季，近48年大部分站点H_i无变化，不过安多站在近28年里H_i变小趋势极为显著，为-0.007/a(P<0.001)。总体来看，自然保护区春、夏2季H_i趋于增大，秋季H_i为减小趋势，而冬季H_i无变化。

3.2.2 年代际变化

在年代际变化度上(图 4)，自然保护区年H_i在1970s—1990s为负距平，其中80年代H_i最低；进入21世纪后转为正距平，以2000s年代最大。H_i的季节分布来看，1970s除夏季外H_i均为负距平；1980s 4季H_i明显偏低，以夏季最为明显，对区域1980s的H_i减小起到了最大的贡献；1990s春、夏2季为负距平，而秋、冬2季为正距平；到2000s—2010s，春、夏2季H_i回升为正距平，尤其是夏季(2000s偏高0.09)，而秋季和冬季H_i为负距平，以秋季偏低最为明显，其中2010s偏低0.06。

3.2.3 突变特征

通过对近48年自然保护区年和4季H_i进行M-K突变检验发现(图 5)，年H_i的UF曲线在1971—1994年期间呈振动变化，从1995年开始快速上升，与UB曲线在1996年相交，且在临界线±1.96之间，确定在1996年发生了突变，这与汪步惟等^[15]的研究基本相同。其中2002—2018年H_i增加趋势超过了95%置信度的临界线，表明年H_i增加趋势是十分显著的；突变前后H_i平均值分别为0.282和0.336，突变后较突变前增加了19.1%。同样，春季H_i也在1995年发生转折，由偏低跃变为偏高态势，突变点前后H_i分别为0.073和0.118 mm，后者较前期偏高61.6%；夏季H_i突变较晚，出现在2000年，突变后较突变前H_i增加21.6%；秋、冬2季H_i未发生突变。

3.3 潜在蒸散及其影响因子的变化

由于影响潜在蒸散变化的因素众多，不同因素之间也相互影响，所以潜在蒸散的变化成因十分复杂。黄会平等^[25]认为影响潜在蒸散量的因子中，第1主成分为热力学因子(气温、水汽压和纬度)，第2主成分为水分因子(相对湿度、降水量)和辐射因子(日照时间)，第3主成分为地理因子(经度)和空气动力学因子(平均风速)，第4主成分为高程因子。以此为据，因该区域气象站点稀少，本文不考虑地理因子，选取7个主要气候因子(平均气温、气温日较差、日照时间、平均风速、平均相对湿度、平均水汽压和降水量)，来分析潜在蒸散与各气候因子的变化。

根据1971—2018年自然保护区年潜在蒸散的变化趋势(表 2)分析，近48年潜在蒸散呈减少趋势，平均每年减少0.27 mm(未通过显著性检验水平)，不过近28年(1991—2018年)潜在蒸散却呈现出显著的增加趋势，增幅为2.26 mm/a(P<0.05)。就季节而言，近48年4季潜在蒸散的线性变化趋势都不是很明显，秋季略为增加，其他3季趋于减少；而近28年秋、冬2季潜在蒸散趋于增加，增幅分别为0.85 mm/a(P<0.05)和1.36 mm/a(P<0.001)，春、夏2季潜在蒸散基本无变化。

近48年，自然保护区年平均气温表现为极显著的升高趋势，升幅为0.046 ℃/a(P<0.001)，4季平均气温也以0.038~0.060 ℃/a的速度显著升高，以冬季升温率最大，尤其是近28年升温率达0.11 ℃/a。年气温日较差趋于变小，为-0.023 ℃/a(P<0.001)，4季气温日较差亦呈减小趋势，其中夏季减幅最大。年日照时间倾向于减少态势，平均每年减少1.92 h(P<0.05)，以夏季减幅最为明显，为-1.71 h/a(P<0.01)。一年4季平均风速均表现为显著减小趋势，为-0.022~-0.039 (m/s)/a(P<0.001)，以春、冬2季尤为突出。年、季平均相对湿度波动大，线性变化趋势不明显，均在±0.10%/a之内。年、季平均水汽压均呈增加趋势，其中夏季增幅最大，为0.012 hPa/a(P<0.05)。年降水量以1.18 mm/a(P<0.05)的速度显著增加，主要表现在春季和夏季。

综上所述，平均气温、降水量和平均水汽压趋于上升(增加)，气温日较差、日照时间、平均风速呈下降趋势，平均相对湿度变化不大，这些气候因素变化的综合影响最终造成了近48年自然保护区潜在蒸散的减少。

3.4 湿润指数与其影响因子的主成分分析 3.4.1 影响湿润指数因子的主成分分析

为了进一步分析各气候因子对年H_i的影响，利用数据处理系统(data processing system, DPS)^[26]的主成分分析方法对7个气候因子进行分析(表 3和表 4)。表 3是各影响因子的总方差分解，可以看出前4个主成分的特征值占总方差的比例达到91.54%，即前4个主成分已经对7个因子所涵盖的信息进行了概括，选用前4个主成分代表潜在蒸散影响因子，其中第1主成分的方差贡献率最大，达到48.06%。表 4是旋转后的主成分载荷矩阵，从表 4中可知第1主成分在平均水汽压和气温日较差2个因子上有较大的载荷量，平均水汽压与第1成分呈正相关，而气温日较差与第1成分呈负相关。第2主成分方差贡献率为21.25%，是次重要的影响因子，该主成分在平均相对湿度、平均气温上的载荷量较大，前者与第2成分呈正相关，后者与第2成分呈负相关。第3主成分方差贡献率为13.6%，该主成分在日照时间、平均风速上的载荷量较大，与日照时间呈正相关，与平均风速呈负相关。第4成分中气温日较差的载荷量最大，与其呈正相关。

3.4.2 主成分回归分析

利用DPS^[26]中的主成分回归方法，为了尽可能少的损失信息，提取前4个主成分(p₁、p₂、p₃、p₄)，累积贡献率达到91.54%，根据主成分矩阵(表 4)，得到主成分方程：

$ \begin{gathered} p_{1}=0.285 X_{1}-0.421 X_{2}-0.383 X_{3}-0.418 X_{4}+ \\ 0.249 X_{5}+0.438 X_{6}+0.408 X_{7}；\end{gathered} $

(3)

$ \begin{gathered} {{p_2} = - 0.571{X_1} + 0.128{X_2} + 0.061{X_3} + 0.127{X_4} + }\\ {0.715{X_5} + 0.344{X_6} - 0.089{X_7}}；\end{gathered} $

(4)

$ \begin{gathered} p_{3}=0.329 X_{1}+0.008 X_{2}+0.592 X_{3}-0.504 X_{4}+ \\ 0.146 X_{5}+0.193 X_{6}-0.479 X_{7} ；\end{gathered} $

(5)

$\begin{gathered} p_{4}=0.428 X_{1}+0.714 X_{2}+0.019 X_{3}+0.243 X_{4}- \\ 0.008 X_{5}+0.433 X_{6}+0.246 X_{7} 。\end{gathered} $

(6)

式中：p_i为主成分；X₁~X₇表述见表 4，将以上4个主成分作为自变量，年H_i作为因变量建立回归模型，即：

$ {H_i} = 0.026\;7{p_1} + 0.003 \; 2{p_2} - 0.024\;4{p_3} + 0.009\;8{p_4} $

(7)

其中修正R²=0.814(P<0.001)，将式(3)~(6)代入式(7)，最终得到主成分标准化回归系数方程：

$ \begin{gathered} H_{i}=0.001 \; 9 X_{1}-0.004 \; 0 X_{2}-0.024 \; 3 X_{3}+0.004 \; 0 X_{4}+ \\ 0.005 \; 3 X_{5}+0.012 \; 3 X_{6}+0.024 \; 7 X_{7} 。\end{gathered} $

(8)

从式(8)标准化回归系数来看，降水量X₇的影响程度最大，表现为正相关；日照时间X₃的影响次之，呈负相关；随后是水汽压X₆(正相关)，其他因子影响相对较小。同样，也分析了4季H_i的主成分回归方程，得到影响程度最大的气候因子都是降水量，次要因子有所不同：春、秋2季是日照时间(负贡献)、夏季是气温日较差(负贡献)、冬季是平均气温(负贡献)。

国内学者采用相关系数分析^{[6-8, 13-14]}、敏感性系数^[5]、主成分分析等^{[25, 27]}方法，选取整个研究时段探讨了不同气象要素对潜在蒸散变化的影响，而汪步惟等^[15]和曹雯等^[28]则分析了潜在蒸散突变前后不同阶段的变化趋势及其影响因素。鉴于此，笔者利用主成分回归分析方法分析了自然保护区年H_i突变前后不同的主导影响因子，结果发现：1996年之前降水量是主导因子(正贡献)，气温日较差是次要因子(负贡献)；1996年之后主导因子也是降水量(正贡献)、次要因子是日照时间(负贡献)。也就是说降水量少、气温日较差大是1971—1995年保护区年H_i偏低的主要因素，而降水量增多、日照时间减少则是1996—2018年保护区年H_i明显偏大的重要因素。从季节看，春季H_i在1995年突变前后，主导因子没变，均是降水量(正贡献)，不同是次要因子，突变前后分别是气温日较差(负贡献)和日照时间(负贡献)；而夏季H_i在2000年突变前，主导因子和次要因子都没有变化，降水量是主导因子(正贡献)，气温日较差(负贡献)是次要因子。

以上分析表明，降水量的多少是主导自然保护区干湿状况的决定因子，随着降水增多H_i增大，在气温升高的背景下，近48年该区域呈明显暖湿化趋势，十分有利于生态环境的改善。

4 讨论

1) 尽管自然保护区地表干湿状况的变化直接影响到自身及其周边区域冰冻圈和生态系统的变化。但由于保护区无气象站，观测资料严重匮乏，当前只能利用周边与其同属羌塘高原湖盆高寒草原区的站点资料来反映自然保护区的气候变化特征。对于自然保护区核心区而言，国内学者目前多以短期观测、遥感反演和融合产品为主获取资料^[19]，来分析羌塘高原环境变化。若需准确认识自然保护区地表干湿状况的变化，从长期来看，在自然保护区合理布设综合气象观测站，实时、有效地进行观测；短期来看，应加大卫星遥感反演数据、再分析资料在羌塘高原的适用性研究，以弥补数据不足的缺憾。

2) 1961—2014年华东、华南以及华中中东部、新疆西北部、东北中部、青海部分地区H_i呈增加趋势，H_i趋于减小的地区主要发生在半干旱半湿润气候区，其中山东东部、河北西部、山西、甘肃和四川盆地、云贵大部分地区呈显著变干的趋势^[3]。2001—2010年青藏高原有25%的区域在逐渐变干，主要集中在高原南部，特别是高原腹地、西藏东部及甘孜州南部等部分区域，柴达木盆地及青海湖区域也有变干倾向^[29]。高原西部、高原北部及高原东北部^[29]、羊卓雍湖流域^[14]、藏北地区^[30]、大通河流域^[31]地表H_i亦呈增大趋势，与自然保护区趋同存异。但与其邻近的黄河源区^[32]、若尔盖湿地^[33]年H_i却呈减小趋势。以上也说明不同地区干湿状况变化差异很大。

3) 从H_i的定义来看，它的变化取决于降水和潜在蒸散2个分量，气候变干或变湿则取决于这2个分量的变化速率，但更与潜在蒸散变化有关，因为潜在蒸散受辐射、温度、相对湿度和风速等气象要素之间相互作用的综合影响。目前有关青藏高原潜在蒸散变化归因研究的文献较多，如Chen等^[4]认为，1961—2000年潜在蒸散减少的主要因子是风速和相对湿度；Zhang等^[34]指出，1971—2004年风速、日照时间减小分别是青藏高原北部、东南部参考蒸散显著下降的主导因子；尹云鹤等^[35]认为，风速减小是影响青藏高原年潜在蒸散降低的主要原因；汪步惟等^[15]分析表明，风速变小是1971—1996年高原年参考蒸散减少的主要因素，相对湿度降低则极大地促进了1997—2014年高原主体(除高原北缘外)参考蒸散的显著增加；周秉荣等^[36]研究表明，最高气温上升、总辐射增加和相对湿度降低是三江源地区年潜在蒸散呈增加趋势的主要原因，风速的影响较小。以上研究也说明，对影响潜在蒸散变化的主导因子尚存在分歧。而针对青藏高原H_i归因研究也有报道，杜军等^[30]和Du等^[37]利用相关分析得到，西藏H_i对降水量、相对湿度和气温日较差的响应最为敏感；王敏等^[29]指出降水是导致高原区域干湿气候空间格局差异的主要因素；笔者也认为降水量是自然保护区1年4季H_i变化的主导因子，不过降水量趋于明显增加，而H_i并明显增大。因此，影响H_i变化的原因仍是个尚未解决的重大科学问题，它的影响机理需要进一步研究。

4) 近48年自然保护区气候暖湿化愈发明显，植被NDVI值趋于增加^[38-39]，这也佐证了地表H_i有增大的趋势。而未来40年西藏H_i趋于增加，环境水热要素相对提高，干旱化程度在逐步减小，较利于生态环境的改善，尤其对干旱的羌塘自然保护区草地荒漠化过程有抑制作用^[40]。因此，分析研究羌塘自然保护区H_i的变化，可为当地生态环境监测提供基础性决策依据。

5 结论

1) 在多年平均状态下，自然保护区大部分地方H_i＜0.5，属于半干旱和干旱区。H_i月际变化呈单峰型，其中3月份最低，为0.03；8月份最高，为0.80。在季尺度上，春、冬2季均为干旱期，夏季为半湿润期，秋季为半干旱期。

2) 近48年(1971—2018)自然保护区年平均气温以0.046 ℃/a的速度显著升高，年降水量以1.18 mm/a的变化趋势明显增加，暖湿化气候特征明显，但自然保护区年H_i却没有表现出明显的增大趋势，只是略有增大，仅为0.001/a，不过夏季H_i增大较明显(0.003/a，P<0.05)，特别是近28年(1991—2018年)增幅尤为突出，达0.007/a(P<0.05)。

3) 年H_i在20世纪70年代到90年代为负距平，气候偏干，其中80年代H_i最低；进入21世纪后H_i转为正距平，以2000s年代最大，气候相对较为湿润。在时间突变点上，年H_i在1996年是一个由相对偏干期转为相对偏湿期的突变点，同样春季和夏季H_i分别在1995年、2000年发生了由偏干转偏湿的突变，而秋、冬2季H_i未出现突变。

4) 主成分回归分析显示，降水量是年H_i变化的主导因子(正贡献)，日照时间起着次要作用(负贡献)。就季节而言，影响H_i变化的主导因子都是降水量。从转折前后来看，1996年之前降水量是主导因子(正贡献)、次要因子(负贡献)是气温日较差；1996年之后主导因子仍是降水量(正贡献)、次要因子却是日照时间(负贡献)。

感谢中国气象科学研究院郭建平研究员对本文英文摘要的修改和润色。