近10 a来青藏高原地表温度时空变化特征分析

朱 智，师春香，谷军霞，梁 晓

(国家气象信息中心，北京 100081)

青藏高原被称为地球的“第三极”，平均海拔高度在4 000 m以上，是世界上平均海拔最高的高原[1]，其强大的动力和热力作用，显著地影响着东亚气候格局、亚洲季风进程和北半球大气环流[2]。研究表明，青藏高原的气候近几十年来发生了明显变化，导致了该地区冰川、湖泊与冻土的改变[3-5]，因此青藏高原陆-气相互作用以及气候变化是全球气候变化研究的重点研究领域之一。地表是陆地和大气之间的交界，其温度能够影响陆-气间的能量和水分交换，通过研究分析青藏高原地表温度的时空变化特征对于进一步研究青藏高原地区的陆-气相互作用具有重要的意义。江灏等[6]基于ISCCP 地表温度反演产品，分析了青藏高原地表热状况的时空分布特征，认为由于海拔高度的影响，青藏高原地表温度与周边地区相比较低；李栋梁等[7]利用青藏高原地区地面观测站的地表温度观测资料，分析了高原异常变化的空间结构和时间演变趋势，结果表明高原地表温度在1970 年代前期呈现明显的下降趋势，而且1970 年代前期至 1990 年代在波动中呈现上升趋势；杨成松等[8]使用MODIS地表温度产品对青藏高原地表温度的空间分布和年际变化进行分析，认为整个青藏高原的地表温度在2003—2013年呈变暖趋势，年平均地表温度(LST)以每年0.015 K的速度升高，振幅以每年0.076 K的速度增长；黄芳芳等[9]利用青藏高原北部高寒气候环境观测研究站地面观测资料，对2001—2012年藏北高原地表温度的年际变化和年内变化规律及其对气候变化的响应特征进行了分析，认为藏北高原地表温度呈缓慢上升趋势，其中冬季上升最快且地表温度年内变化和日变化特征明显，且气温的变化相对于地表温度存在一定时间的滞后响应。

由于青藏高原地区的地表温度观测站点不仅数量较少，并且呈现东部密集西部稀疏、南部密集北部稀疏的空间分布特征，并且高原地区气候条件复杂，因此需要一套时空分布连续、适用性较好的地表温度模拟资料作为观测资料的补充，应用于气候变化研究。目前对于长时间序列地表温度的模拟主要依赖于再分析资料和陆面模式模拟，何冬燕等[10]利用青藏高原地区的地面气象台站地表温度观测数据，比较和分析了ERA-Interim、NCEP/NCAR、NCEP/DOE再分析地表温度资料在青藏高原的适用性，认为三种再分析地表温度地表均能够较好地描述了高原地表温度的季节变化和年际变化特征，但是对于地表温度的长期变化趋势估计不足；秦艳慧等[11]基于青藏高原高海拔地区地表温度观测数据，评估了ERA-Interim再分析地表温度的精度，结果表明ERA-Interim再分析地表温度与观测数资料的变化趋势一致，但是呈现显著低估现象。

中国的地表温度历史观测资料主要来源于中国气象局布设的地面气象观测站，历史上使用地面温度表(0 cm温度表)对地表温度进行定时人工观测(3～4次/日)，随着观测手段的不断改进，从2002年开始，地表温度观测逐步由人工定时观测变更为自动化逐小时观测；由于观测自动化是逐步推进的，因此在实际业务中，存在人工定时观测与自动化逐小时观测并存的现象。刘小宁等[12]对国家级基准站自动站观测与人工观测并存时期的地表温度观测资料进行了分析，认为观测方式的改变会对地表温度序列的均一性产生影响。因此在应用于气候变化研究时，需要对人工观测和自动观测的地表温度予以区分，而现有的研究均忽略了基于不同观测方式得到的地表温度观测资料存在的非均一性现象。同时随着第四代再分析资料的发布和陆面过程模拟技术的发展，可用于青藏高原地区陆-气相互作用研究的地表温度模拟资料不断丰富，有必要基于长时间序列观测资料对多种地表温度模拟资料在青藏高原地区的适用性进行评估分析。

所使用的观测数据来源于国家气象信息中心制作的地面基本气象要素观测数据集“中国国家级地面气象站基本气象要素定时值数据集(V3.0)”，并针对数据集中地表温度观测方式进行了分析，选取了2008—2017年青藏高原地区已转为自动化观测的133个地面气象观测站地表温度观测资料，对近10 a来青藏高原地区的地表温度时空变化特征进行了分析，并对五套地表温度模拟资料(ERA-5、ERA-Interim、GLDAS2.1-Noah、MERRA-2、JRA-55)在青藏高原地区的适用性进行评估，以期为青藏高原地区陆-气相互作用以及气候变化研究提供有益的参考。

1 数据与方法

1.1 观测数据

所使用的观测数据来源于国家气象信息中心制作的“中国国家级地面气象站基本气象要素定时值数据集(V3.0)”，在数据集的制作过程中，对可疑和错误的观测数据进行了人工核查与更正，并对数据缺测现象进行了补充，该数据集相对于以往发布的地面观测数据产品，在数据完整性和数据质量方面均有明显提高，也是目前内容较为完整、要素较为齐全、可用性较强、保留原始观测信息较多的地面观测数据集。随着城市的发展，地面气象站存在站点迁移的现象，有可能会对研究产生一定的影响。为了剔除站点迁移因素，选取了数据集中青藏高原地区2008—2017年连续观测且站点位置未发生变化的133个站的每日4次(02：00、08：00、14：00、20：00)的定时地表温度观测数据，站点空间分布如图1所示，其中海拔1 000～2 000 m的站点9个、海拔2 000～3 000 m的站点43个、海拔3 000～4 000 m的站点58个、海拔4 000～5 000 m的站点23个。

图1 青藏高原地区地表温度观测站点空间分布

1.2 地表温度模拟资料

1.2.1 陆面模式地表温度资料

所使用的陆面模式地表温度资料是GLDAS2.1-Noah地表温度资料。

全球陆面数据同化系统(GLDAS)产品是目前陆面水文过程研究常用的模拟产品[13-15]，目前发布的最新版本为GLDAS2.1。GLDAS2.1驱动数据来源于三个部分：①近地面气温、气压、湿度、风速驱动数据来源于美国国家海洋和大气管理局(National oceanic and atmospheric administration，NOAA)的全球数据同化系统(global data assimilation system，GDAS)；②辐射驱动数据来源于经过偏差订正后的美国空军气象局的AGRMET系统(the Air force weather agency’s AGRicultural METeorological modeling system)；③降水驱动数据来源于GPCP降水产品和全球数据同化系统的降水预报结果[16]。目前GLDAS2.1使用了Noah模式(Noah 3.3)对全球陆面状况进行模拟，时间分辨率是3 h，空间分辨率有两种(0.25°×0.25°和1°×1°)，使用的是空间分辨率0.25°×0.25°、时间分辨率3 h的GLDAS2.1-Noah地表温度资料。

1.2.2 再分析地表温度资料

使用的再分析地表温度资料包括：ERA-5再分析地表温度资料、ERA-Interim再分析地表温度资料、MERRA-2再分析地表温度资料、JRA-55再分析地表温度资料。

ERA-5再分析资料是欧洲中期天气预报中心(European Centre for Medium-range Weather Forecasts，ECMWF)研制的最新一代再分析资料[17]，使用Ensemble-4Dvar同化方法，也是首个能够提供不确定性信息的再分析资料，大气模式分析场的空间分辨率为0.25°×0.25°，时间分辨率是6 h，截至2018年12月发布了2005年至今的再分析产品；ERA-Interim再分析资料欧洲中期天气预报中心研制的第三代再分析资料[18-19]，空间分辨率是0.75°×0.75°，时间分辨率是6 h，使用了4Dvar同化方法，时间覆盖是1979年至今；MERRA-2是美国国家航空航天局(National Aeronautics and Space Administration，NASA)研制的第三代再分析资料[20]，空间分辨率是0.5°×0.625°，时间分辨率是1 h，在再分析资料制作过程中首次同化了地基和天基气溶胶光学厚度观测数据，使用了GSI-3Dvar同化方法，时间覆盖是1980年至今；JRA-55是日本气象厅研发的第三代再分析资料[21-23]，使用了4Dvar同化方法，空间分辨率T319L60，时间分辨率6 h，时间覆盖是1958年至今。

2 结果与分析

2.1 青藏高原地区地表温度气候态空间分布

图2给出了2008—2017年共10 a观测、ERA-5、ERA-Interim、GLDAS2.1-Noah、MERRA-2、JRA-55的年平均地表温度空间分布。从图2(a)可以看出在高原东南部的雅鲁藏布江流域和西北部的柴达木盆地存在两个高值中心，地表温度能够达到8 ℃以上，最高值出现在青藏高原最南端的云南省福贡地面气象观测站(98.8°E、26.9°N)，达到19 ℃以上；在高原的中部存在一个低值中心，年平均地表温度在0～4 ℃之间，最低值出现在海拔4 600 m以上的青海省五道梁地面气象观测站(93.1°E、35.2°N)，年平均地表温度接近0 ℃；同时在青藏高原地区，地表温度呈现边缘地区高、中心地区低，由边缘地区向中心地区递减的变化趋势，这与何冬燕等[10]得出的结论基本一致。

五套地表温度模拟资料(ERA-5、ERA-Interim、GLDAS2.1-Noah、MERRA-2、JRA-55)在青藏高原的空间分布均符合“两高一低”的空间分布特征[图2(b)～图2(f)]，能够描述出在雅鲁藏布江流域和柴达木盆地的两个地表温度高值中心，以及在青藏高原地区中部的低值中心。在青藏高原地区的柴达木盆地，ERA-5再分析地表温度资料和GLDAS2.1-Noah地表温度资料达到4～8 ℃，特别是在柴达木盆地的低海拔地区，达到8 ℃以上，能够与实际观测资料的空间分布状况较好地吻合，这是由于ERA-5再分析产品的研制过程中，使用了更新的数值预报模式以及更先进的同化分析方法，而在GLDAS2.1-Noah资料的研制过程中，对于驱动数据中的地面气温要素和辐射要素进行订正，因此两者能够更好地描述出柴达木盆地的地表温度空间分布特征；而其他三种地表温度模拟资料(ERA-Interim、MERRA-2、JRA-55)尽管能够描述出柴达木盆地的地表温度高值中心，但是相对于实际观测资料，出现了明显的低估，这表明三种产品所使用的数值模式和同化方法还存在一定的不足，需要进一步改进。在青藏高原地区的雅鲁藏布江流域和中西部地区，五套地表温度资料均出现了不同程度的低估现象，这可能是由于在地形复杂区，数值模式和陆面模式所使用的高程数据与实际观测站点的海拔高度存在一定的差异。

图2 2008—2017年年平均地表温度空间分布

图3 2008—2017年年平均地表温度观测资料与站点地理信息散点图

为了进一步研究高程差异对青藏高原地表温度模拟的影响，将五套地表温度资料的多年平均值和所使用的高程数据利用双线性插值方法插值到观测站点，如图3所示。根据图3(a)中给出的地表温度随高程的递减率(-0.34 ℃/100 m)，对地表温度插值结果进行了高程订正，并以地表温度观测资料为基准，对高程订正前后的五套地表温度模拟资料精度进行了定量的比较，如图4所示。

图4 高程订正前后的五套地表温度模拟资料与观测资料散点图

相对于地表温度观测资料，五套地表温度模拟资料均出现了低估现象[图4(a)～图4(d)]；从平均偏差来看，五套地表温度资料的平均偏差均超过-6 ℃，ERA-5再分析地表温度资料的平均偏差最大，达到-9.898 ℃，其次是ERA-Interim再分析地表温度资料，平均偏差达到-8.639 ℃；JRA-55和MERRA-2再分析地表温度资料的平均偏差相近，在-8～-7 ℃之间；平均偏差最小的是GLDAS2.1-Noah地表温度资料，但超过-6 ℃；从相关系数来看，五套地表温度资料与观测资料均在0.5以上，能够通过P<0.01的显著性检验，其中与观测资料相关性最好的是MERRA-2再分析地表温度资料，相关系数达到0.7以上，与观测资料相关性最差的是ERA-5再分析地表温度资料，相关系数为0.545；在经过高程订正后[图4(e)～图4(h)]，五套地表温度模拟资料的精度均有一定的改进，相关系数提升了0.2左右，平均偏差下降了2 ℃左右，高程订正后最接近观测资料的是GLDAS2.1-Noah地表温度资料和MERRA-2再分析地表温度资料，相关系数达到0.89以上，能够通过P<0.01的显著性检验，偏差在-5 ℃以内。从高程订正前后，五套地表温度资料与观测资料的对比来看，数值模式和陆面模式所使用的高程数据与实际情况之间存在的差异，是造成地表温度模拟存在负偏差的原因之一。

2.2 青藏高原地区地表温度时间变化特征

按照青藏高原地区的地表温度观测站点海拔高度，划分为4个层次(1 000～2 000、2 000～3 000、3 000～4 000、4 000～5 000 m)，各海拔层次的地表温度观测站空间分布见图1。由2008—2017年各海拔层次五套地表温度模拟资料与观测资料的变化距平(图5)可见，青藏高原地区的地表温度呈现增暖的趋势，且增暖幅度随着海拔高度层次的增加而更加平缓，海拔层次1 000～2 000 m的地表温度增暖幅度最为剧烈，达到0.0716 ℃/a，海拔层次4 000～5 000 m的地表温度增暖幅度最为平缓，达到0.0157 ℃/a；青藏高原地区的地表温度年际变化在±1.5 ℃以内，年际变化波动范围随海拔高度增加而减少，海拔层次1 000～2 000 m的地表温度年际变化波动范围最大，年平均地表温度最高值与最低值的差距超过1.3 ℃，海拔层次4 000～5 000 m的地表温度年际变化波动范围最小，年平均地表温度最高值与最低值的差距达到0.84 ℃，年平均地表温度最高值与最低值的差距超过1.3 ℃。在海拔高度1 000～2 000 m层次，2007—2012年的地表温度均低于平均值，而在2013年出现了明显的增长，增长幅度超过1.3 ℃，并在其后的2014—2017年均高于平均值；其他海拔高度层次的地表温度变化趋势与1 000～2 000 m海拔高度层次基本一致，当1 000～2 000 m海拔高度层次某一年地表温度较前一年升高(或降低)时，其他海拔高度层次地表温度表现出基本一样的变化特征，并且均在2007—2008年出现了变暖趋势，变化幅度随着海拔高度层次的增加而更加剧烈，特别是在4 000～5 000 m层次，变化幅度超过0.8 ℃，也是近10 a来该海拔高度层次地表温度变化最剧烈的一次。

图5 不同海拔高度层次五套地表温度模拟资料与观测资料距平的年际变化

在五套地表温度模拟资料中，ERA-5地表温度模拟资料年际变化幅度较为明显，而其他四套地表温度模拟资料年际变化较小。在海拔高度1 000～2 000 m层次，地表温度观测资料在2010—2011年出现了小幅度的增加(增加幅度为0.04 ℃)，五套地表温度模拟资料均出现了相反的变化特征(小幅度减小)；ERA-5、GLDAS2.1-Noah、MERRA-2、JRA-55均能够描述2013—2014年地表温度观测资料的小幅度降低现象，而ERA-Interim地表温度模拟资料出现了相反的变化趋势；GLDAS2.1-Noah地表温度模拟资料在2014—2017年出现了持续降低现象，与实际观测不符。在海拔高度2 000～3 000、3 000～4 000、4 000～5 000 m层次，ERA-5地表温度模拟的变化幅度较为剧烈，特别是在2008—2011年，出现了尽管变化趋势与观测资料一致，但变化幅度更为剧烈；其他四种地表温度模拟资料的变化特征和变化幅度与观测资料基本一致，当某一年地表温度观测资料较前一年升高(或降低)时，能够表现出相同的变化特征。

图6 不同海拔高度层次五套地表温度模拟资料与观测资料距平的月际变化

从青藏高原地区月平均地表温度距平的年内变化(图6)可以看出，各海拔高度层次地表温度的最高值均出现在7月，最低值均出现在1月，而且地表温度在1—7月单调递增，在从7—12月单调递减，呈现单峰形分布特征；各海拔层次月平均地表温度距平的变化范围基本一致，峰值与谷值的差距在25 ℃左右；夏季的地表温度距平月际变化较为平稳，这是由于随着夏季青藏高原地区气温的升高和降水的增加，植被覆盖度也开始增加[24-25]，而植被有调节地表温度的作用[26]；而在其他季节，月际变化更为明显，这样的变化特征与秦艳慧等[27]的研究结果较为一致。五套地表温度模拟资料(ERA-5、ERA-Interim、GLDAS2.1-Noah、MERRA-2、JRA-55)均能够表现出与观测资料相似的单峰形变化特征，但是相对于观测资料，ERA-5和ERA-Interim再分析地表温度资料的变化幅度更为剧烈，其他三种资料的变化幅度更加接近实际观测。

2.3 青藏高原多年冻土区地表温度变化特征分析

青藏高原是中纬度多年冻土分布面积最大、海拔最高的地区[28]，而地表温度作为陆-气能量交换中重要的物理量，是许多冻土分布模型和区域陆面模式重要的参数[27]，在多年冻土区的气候变化研究中起着重要的作用，因此对于青藏高原多年冻土区的地表温度研究是非常重要的工作。但是青藏高原多年冻土区的地表温度观测资料较为缺乏，仅有中国气象局布设的3个地面气象站(表1)能够提供长时间序列的多年冻土区地表温度观测资料，因此基于多年冻土区3个地面气象站2008—2017年的地表温度观测数据，对近10 a来青藏高原多年冻土区温度变化特征进行了分析，并检验五套地表温度模拟资料(ERA-5、ERA-Interim、GLDAS2.1-Noah、MERRA-2、JRA-55)在青藏高原多年冻土区的适用性，为地表温度模拟资料在观测资料较为缺乏的青藏高原多年冻土区的应用提供参考。

表1 青藏高原多年冻土区地面气象站信息

图7 青藏高原多年冻土区地面气象站月平均地表温度时间序列

从青藏高原多年冻土区3个地面气象站地表温度的月际变化来看(图7)，地表温度呈现明显的季节性变化特征，夏季达到一年中的最高值，冬季达到一年中的最低值；五套地表温度模拟资料均能够很好地模拟出地表温度的季节性变化，且能够与观测资料的变化趋势基本保持一致,相关系数均超过0.96，能够通过P<0.01的显著性检验。从偏差来看，相对于观测资料，五套地表温度模拟资料均呈现负偏差，出现了较为显著的低估现象；在五套地表温度模拟资料中，GLDAS2.1-Noah地表温度模拟资料的表现最好，偏差在-2.5～1.9 ℃之间，JRA-55地表温度模拟资料的表现略差于GLDAS2.1-Noah地表温度模拟资料，偏差分布在-4.5～3.3 ℃之间，其他三套地表温度模拟资料的偏差均大于-5 ℃，这样的现象与孙帅等[29]的研究结果较为一致。

在中国气象局布设的地面气象观测站中，主要使用钼电阻地温传感器进行自动化的地表温度观测；根据《地面气象观测规范》，地表温度观测仪器必须布设在疏松、平整、无草的裸地上[30]，这与地表温度模拟资料的所使用的地表覆盖产品存在一定的差异，而在模式模拟的过程中，植被覆盖能够影响大气与土壤之间的能量交换，进而降低地表温度，这可能是地表温度模拟资料均出现了低估的原因之一[31-32]。

从图7也可以看出，地表温度模拟资料的低估现象具有一定的季节性特征，在春季和秋季更为接近观测资料，而在夏季和冬季的低估更为明显；相对于其他四套地表温度模拟资料，ERA-5再分析地表温度资料在冬季出现了更为显著的低估现象，特别是在安多地面气象站[图7(b)]和沱沱河地面气象站[图7(c)]，出现了明显的负偏差，这可能与ERA-5再分析资料对于积雪观测资料的同化策略有关，在ERA-5再分析资料的生成过程中，仅同化了海拔1500 m以下的IMS(interactive multi-sensor snow and ice mapping system)积雪产品[33]，造成高海拔地区地表温度模拟缺少积雪观测信息的影响。

3 结论与讨论

基于2008—2017年青藏高原地区的地表温度观测数据和五套地表温度模拟资料，分析了近10 a来青藏高原地区地表温度变化趋势，并对地表温度模拟资料的适用性进行了评估，得到了以下主要结论：

(1)从空间分布看，青藏高原地区的地表温度在雅鲁藏布江流域和柴达木盆地存在两个地表温度高值中心，在高原中部地区存在一个地表温度低值中心；青藏高原地区的地表温度主要受海拔高度影响，地表温度随着海拔高度的增加而降低，同时服从纬度的地带性规律，随着纬度的增加，地表温度呈现下降趋势。

(2)五套地表温度模拟资料均能够描述出青藏高原地区地表温度的空间分布特征；在地表温度的低值区，相对于其他三种地表温度模拟资料，ERA-5再分析地表温度资料和GLDAS2.1-Noah地表温度资料能够与实际观测资料的空间分布更好地吻合；五套地表温度模拟资料均出现了不同程度的低估现象，主要是由于模拟资料所使用的高程数据与实际观测存在差异，在进行高程订正后，五套地表温度模拟资料均更加接近实际观测。

(3)在近10 a，青藏高原地区的地表温度年际变化在±1.5 ℃以内，呈现增暖的趋势，且增暖幅度随着海拔高度层次的增加而更加平缓；从年内变化看，青藏高原地区的地表温度呈现单峰形的分布特征，谷值出现在在冬季，峰值出现在夏季；五套地表温度模拟资料能够模拟出地表温度的年际变化和年内变化特征，但是相对于实际观测，ERA-5和ERA-Interim再分析地表温度资料的变化幅度更为剧烈。

(4)在青藏高原多年冻土区3个地面气象站，五套地表温度模拟资料的变化趋势能够与实际观测很好地吻合，相关系数均超过0.96，但是均出现了低估现象，且低估现象在夏季和冬季更为显著，具有明显的季节特征，这是因为根据观测业务要求，地表温度观测地段为无植被覆盖的裸土，与地表温度模拟资料所使用的地表覆盖静态数据存在一定的差异。