1 引言

多年冻土与冰川、 积雪是青藏高原冰冻圈的重要组成部分，在全球气候变化的背景下，冰冻圈的响应对区域环境和生态系统必将产生重要影响^{[1, 2]}。青藏高原多年冻土面积约1.05×10⁶-1.5×10⁶km^2[3]，分布广、 厚度薄、 稳定性差^[4]。在过去几十年的气候变暖背景下，青藏高原多年冻土的分布格局正在变化，气候变暖显著和人类活动日益增加，加剧了冻土退化，已经影响到了高原多年冻土环境，引发了一系列生态和环境问题，亟待系统的长期观测、 试验和模型研究^{[5, 6]}。青藏高原地区多年冻土的特征和演化趋势与全球其他地区的多年冻土有一定相似性。使用现代气候数据模拟得出的模型对多年冻土发展趋势进行预测，结果表明在西伯利亚^{[7, 8, 9, 10]}、 北极圈内^[11]、 阿拉斯加^[12]、 加拿大北部^[10]、 北欧甚至整个北半球^[13]的多年冻土随着温度上升都有面积逐渐减小、 活动层厚度增大的趋势，尤其是北极圈内多年冻土对温度变化最为敏感^[14]。针对北半球^[15]、 北美^[16]、 欧亚大陆^[17]和西北欧^[18]分别模拟重建末次冰盛期多年冻土分布，显示这一时期各地的多年冻土面积、 厚度较现代增大，并且在降温条件下有扩张趋势。Renssen和Isarin^[19]重建了欧洲新仙女木时期的气候并发现在降温条件下多年冻土面积明显增加并有扩张趋势，Wetterich等^[20]针对俄罗斯晚第四纪气候变化的研究则发现在升温时期如全新世大暖期多年冻土面积大幅减小。

青藏高原的多年冻土是第四纪地质时期高原隆升过程与气候变化叠加的产物。随着第四纪冰期-间冰期循环，多年冻土周期性扩张和收缩^[21]。冰期时多年冻土扩张的范围界线、 间冰期时多年冻土收缩的范围界线以及扩张-收缩过程中的地理表现是青藏高原环境变化研究的重要内容^{[4, 22, 23, 24]}。本文根据气象和冻土观测等资料绘制了青藏高原地区现代气温等温线和多年冻土分类边界，在此基础上对各类多年冻土边界分布与等温线的对应关系进行了分析。应用现代气温等温线和古气候研究资料对末次盛冰期(Last Glacial Maximum，简称LGM)和全新世大暖期(Holocene Megathermal)的气温等温线进行了重建，对多年冻土区内地表过程响应机制展开讨论，为探讨末次冰期以来青藏高原环境变化提供依据。

本文使用数据主要来自:《中国冰川冻土沙漠图》^[25]； 1981-2010年新疆、 西藏、 甘肃、 青海、 云南及四川等省气象站点月平均气温资料(气象科学数据共享网http://cdc.cma.gov.cn)； ASTER GDEM 90m分辨率高程数据(NASA)，采用WGS1984坐标系。

由于研究区的气象站点东部多，西部少( 图1)，应用气象观测数据资料和插值方法绘制研究区等温线受到气象站点数量有限、 分布不均匀的限制。为了弥补这一缺陷，需要在研究区进行插值补点，实现更高的科学性，使等温线更加客观。采用SPSS和EXCEL工具对数据进行回归分析，从而尝试给出已有站点温度和海拔、 纬度间的关系，通过温度随纬度、 海拔的直减率，计算出补点的温度。利用ArcGis软件进行插值和图像处理。采用回归分析法绘制现代平均气温等温线，在此基础上重建末次盛冰期和全新世大暖期的温度分布和多年冻土边界。

图1(Fig.1)

图1 自然邻域法绘制的多年平均气温等温线图 Fig.1 Isotherm of Qinghai-Tibet Plateau by natural neighbor interpolation

本文使用的基础数据包含青藏高原及邻近地区147个气象站点的月均温数据，取1980-2010年共30年的资料进行分析。气象站点在高原地区分布不均匀，中西部地区站点非常稀疏； 研究区高差大、 面积大，跨越的经度和纬度范围较大。在绘制研究区年均温度分布图时，必须考虑到上述因素的影响，进行必要的插值补点。对上述资料进行了质量控制检查，并对所选台站个别缺测资料用与其相关最好的几个站的资料通过多元回归方法进行插补。

本文中补点方法分两种情况: 1)传统方法。此方法主要用于青藏高原地区东部和南部。选取相邻两个站点并在两点间连线，沿线寻找海拔变化 ±1000m 的点并沿其地貌走向补3-5个点，此方法能使控制点覆盖站点之间无法观测的高海拔和低谷地区，并使用回归方法补全温度变化较大的点从而得到更详细的等温线。2)网格方法。主要用于西北无站点地区，该区域平均海拔高站点分布极少，需要间隔经纬度绘制20×20网格进行补点。在缺少观测数据的海拔变化幅度小的地区使用此方法可有效预测补充温度数据空白。

已有研究结果表明^{[26, 27]}，样本本身的空间分布是影响插值精度的重要因素，合理的采样设计是必要的前提； 对于台站稀少的地区，必须把随机插值方法和确定性方法结合起来估计气候变量的空间分布。

普通Kriging插值能够反映出青藏高原气温分布的一定的空间结构，但结果并不好； 由于考虑了高度变量，Co-kriging插值表现出一定程度的性能改进，但因为台站海拔高度偏低，插值结果依然不理想^[28]。由于青藏高原上尤其是高原西北地区气象台站稀少( 图1)，有必要根据已有站点数据，采用综合方法分析气温随海拔和纬度的变化，将其订正、 丰富后，再对它们进行插值分析，其结果较为正确地反映了青藏高原气温空间分布的特征，误差远小于其他方法。

利用Excel对研究区内各站点的30年内平均气温与海拔、 纬度之间的关系进行了线性多元回归分析。可以得到:

其中L为某地北纬纬度(°N)，H为某地海拔(m a.s.l.)，y为某地的多年平均气温(℃)。

该回归结果的R²=0.883269，且通过α为0.05的F值检验，所以结果是比较可信的。这样只要知道海拔和纬度就可经计算出某地的平均温度。海拔和纬度可以从DEM图读出。

线性回归分析后即可得到补点处的推算温度结果，以补点和站点的温度数据为基础，利用空间插值可提供每个计算栅格的温度值。在ArcGIS中关联到dbs中后生成shp文件，使用Kriging方法和自然邻域法分别进行了空间插值，然后分别对两个栅格文件进行了表面分析中的等值线分析，得到青藏高原现代年平均气温分布图( 图1)。通过与实际气象数据相比对，发现通过自然领域法插值得到的等温线图与实际情况更为符合，因此采用该方法进行进一步分析。

青藏高原现存的多年冻土层主体是末次冰期所形成的^[21]，后经全新世以来的多次气候波动，演变成目前的状况。多年冻土在冰期扩张以及间冰期退缩留下的证据广泛分布在青藏高原地区，成为古气候重建的重要依据^{[4, 5, 29, 30, 31, 32]}。

基于前人相关研究总结，末次盛冰期的最低温时间始自22ka左右^[33]，根据前人^{[34, 35, 36, 37, 38]}的相关研究数据，本文选取22-15ka作为末次盛冰期时间。关于青藏高原末次盛冰期的降温值，本文采用Thompson等^[39]、 姚檀栋等^[40]和施雅风^[41]的相关数据，选择降温值在6-11℃，平均7℃。末次盛冰期的降温幅度研究在国际上也是备受关注，应用材料不同、 模型不同，得出的结果有差异。Deimling等^[42]应用模型推测末次盛冰期全球降温幅度在4.3-9.8℃，两极地区降温幅度大，热带小； 全球平均降温5.8℃。Shakun和Carlson^[43]的研究报告指出末次盛冰期全球降温至少是4.9℃； 末次盛冰期时，由两极到赤道，降温幅度有很大差异，北半球格陵兰降温幅度最大，达到17℃，南半球南极大陆降温幅度在8-9℃，中纬度地区变化较大，在5-7℃，低纬到赤道，仅有1-4℃降温。以上所有资料，缺东亚和青藏高原，考虑到高山高原的地形效应，青藏高原的变温幅度要大，因此，施雅风先生^[33]指出高原在末次盛冰期降温幅度在7℃是可信的。

全新世大暖期的时代，国内学者多采用9-3ka，也有取较窄的8-4ka^[44]； 国外多用9-5ka^{[45, 46]}。基于前人的研究^{[47, 48, 49]}本文取全新世大暖期的起始时间为9-3ka。综合植物孢粉、 湖泊、 冰缘和泥炭等证据^[50]，喜马拉雅山地区全新世大暖期孢粉化石所在地记录的林线较现代林线高出500-800m，按照每100m升温0.6℃计算，喜马拉雅山地区全新世大暖期平均气温较现代高3-5℃^[51]。考虑到高山地区对温度变化的敏感性，可推测青藏高原地区全新世大暖期气温较现代普遍高出3-5℃，结合施雅风^[33]的《中国第四纪冰川与环境变化》一书中所绘制的“全新世大暖期盛期年平均温度与现今差值的等值线图”，本文按照平均升温4℃重建高原多年冻土分布界线。

现代多年冻土的分布范围资料来源于《中国冰川冻土沙漠图》^[25]和Ran等^[3]，将其图中多年冻土进行数字化，建立坐标系统，进行地理配准，分大片多年冻土、 岛状多年冻土和高山多年冻土进行绘制，得到青藏高原现代多年冻土界线的DEM图( 图2)。

图2(Fig.2)

图2 青藏高原现代多年冻土分布图 Fig.2 Existing permafrost of Qinghai-Tibet Plateau

对得到的青藏高原现代多年冻土界线( 图2)和现代平均气温分布图( 图1)进行匹配，发现现代大片多年冻土界线处的年平均气温为-2-0℃； 岛状多年冻土界线处的年平均气温比较复杂，界线西南部的年平均气温约为3℃，界线东南部约为-1-3℃，界线东北部约为0-2℃； 高山多年冻土主要存在祁连山和喜马拉雅山地区，祁连山区高山多年冻土界线处的年平均气温约为1-2℃，喜马拉雅山地区界线约为3-4℃。通过对比，现代多年冻土界线与现代等温线间的一致性总体来说还是比较好的^{[33, 52, 53]}。

在青藏高原现代多年冻土与温度制图的基础上，参照米德生^[54]对各类型多年冻土分布的温度与海拔条件，以年均温下降7℃和升高4℃在青藏高原现代多年冻土界线的基础上，对末次盛冰期和全新世大暖期的三类多年冻土(大片多年冻土、 岛状多年冻土和高山多年冻土)界线进行重建，结果如 图3和 图4。本文中的末次盛冰期和全新世大暖期都是范围较长的时间段，期间都发生过程度不同的冷暖波动事件，导致两个时期内多年冻土的面积范围也是波动的。这里选取了两时期具有代表性的温度变化数值对多年冻土边界范围重建，因此代表了两时期内多年冻土在理想的平均温度条件下的稳定状态。

图3(Fig.3)

图3 青藏高原末次冰盛期(22-15ka)多年冻土分布图 Fig.3 Permafrost extent of the Qinghai-Tibet Plateau during the Last Glacial Maximum(22-15ka)

图4(Fig.4)

图4 青藏高原全新世大暖期(9-3ka)多年冻土分布图 Fig.4 Permafrost extent of the Qinghai-Tibet Plateau during the Holocene Megathermal(9-3ka)

在全球变暖的大背景下，青藏高原现代多年冻土整体呈退化趋势。 首先，多年冻土区地温升高，冻土厚度减小，活动层厚度加大。大片多年冻土区的昆仑山口的钻孔信息显示，在1985-2000年，冻土层6-15m深度的平均温度上升了0.2-0.4℃，15m以下部分地温变化微弱^[4]，整体厚度变化介于80-100m，比较稳定，羌塘盆地活动层厚度最小(0.87-1.78m)，以羌塘盆地为中心向四周活动层厚度不断增加^[55]； 岛状多年冻土区北界的西大滩(35.72°N，94.14°E) 自1983年以来多年冻土平均地温升高0.2-0.3℃，冻土厚度缩小至少4m，部分地区多年冻土已消失，活动层厚度增大明显，甚至出现双向退缩^[4]； 高山多年冻土主要受海拔控制，其地温略有增长但不显著且多年冻土和活动层厚度较为稳定。

其次，多年冻土区广泛发育热喀斯特地貌。截至2009年9月青藏高原地区热喀斯特湖总积水面积约3.5×10⁵m^2[56]。高含冰量多年冻土区发育的热喀斯特湖占总数量的80 ％ ，约占总面积的70 ％ ，低含冰量多年冻土区发育的热喀斯特湖占总数的20 ％ ^[57]。

多年冻土退化对生态系统也有很大影响。随着多年冻土的持续退化，土壤温度逐渐升高，土壤含水量逐渐下降，有机质含量降低。多年冻土的退化不仅不利于高原多年冻土生态系统功能的发挥，而且不利于草地畜牧业的发展，最严重的还会造成沙漠化^[58]。

多年冻土的消融引起冻胀和沉陷，并直接影响工程建筑的稳定性。据吴青柏和牛富俊^[59]报道，青藏公路的破坏路段中80 ％ 以上是由于这一原因所导致。

(1)末次盛冰期 青藏高原广泛分布着第四纪冰缘地貌和沉积证据，成为恢复古多年冻土范围和古环境重建的重要证据，周幼吾等^[21]对1998年以前的文献进行了全面梳理评述。青藏高原有可靠年代证据的古多年冻土出现在“倒二冰期”(MIS 6，即191-130ka)、 末次间冰期(MIS 5，即130-71ka)、 末次冰期(MIS 4-MIS 2，即71-14ka)，而全新世古多年冻土的证据就更多了^{[60, 61]}。

末次盛冰期时高原大片多年冻土界线在现代年平均气温5-7℃附近，范围向南部和东部大为扩展； 现代为季节冻土区的柴达木盆地和共和盆地等地区在末次盛冰期时年平均气温降到-2℃以下，基本被大片多年冻土所占据( 图3)。同时根据古冰缘遗迹的分布情况( 图3)，末次盛冰期时柴达木盆地冷湖、 青海湖周围及共和盆地以及若尔盖盆地北缘等地的多年冻土为连续分布，即大片多年冻土区。

本文的重建和分析表明，大片多年冻土北界在高原西北部的帕米尔高原和喀喇昆仑山地区(38.47°-38.79°N)，多年冻土下界约为2800m，比现代低1500m以上； 到高原北部扩展到塔里木盆地南缘(37°-39°N)，下界海拔2500-3000m； 到祁连山北麓(39°N)，多年冻土下界约为2000-2300m，比现代低1500m以上； 到高原的东北部，多年冻土北界(36°-39°N) 下界海拔在 36°-37°N约为2900-3100m、 37°-39°N约为3700-3900m。Yang等^[62]在昆仑山北麓-塔里木盆地南缘的河流阶地中发现大量末次盛冰期的冰楔和冰碛物，证明该地区在末次盛冰期温度较现代低10℃以上。

根据潘保田和陈发虎^[22]发现的砂楔遗迹，大片多年冻土东界(34°-36°N) 向东扩展至柴达木盆地东缘(98°E) 。我国境内的大片多年冻土西界到达帕米尔高原西缘(74°E) 。李树德和贺益贤^[63]、 张维信^[64]、 王绍令和李树德^[65]、 徐叔鹰等^[66]、 常晓丽等^[67]在末次冰盛期大片多年冻土区发现大量该时期冰楔、 土楔、 砂楔和冻融褶皱遗迹，证明末次盛冰期多年冻土在高原地区广泛发育。

末次盛冰期大片多年冻土的南界(34°-29°N) 向南扩展到 29°N附近，与岛状多年冻土北界重合； 在现代被岛状多年冻土区所包围的狭长的季节冻土区中，末次盛冰期大片多年冻土区达到 30°N左右，界线位于羊八井附近，下界海拔约为3900-4000m； 岛状多年冻土界线在南北两侧受到大片多年冻土和高山多年冻土的影响强烈收缩，喜马拉雅山北麓多年冻土退缩至 31°N。 但岛状多年冻土界线向东扩展至 100°E左右( 图3)。

喜马拉雅山地区高山多年冻土界线在平面上变化不大，喜马拉雅山北坡多年冻土下界在 29°N附近，海拔约为3900-4000m。祁连山北麓多年冻土向北扩张，下界在 40°N附近，海拔约为2000m； 向东扩张至 103°E，下界海拔约1800m； 向南扩张至 27°N的横断山脉附近，下界海拔约2000m( 图3)。笔者根据野外考察和年代数据，发现个别地区(如格尔木河)的末次盛冰期冰缘证据较现代多年冻土下界低1500m左右。潘保田和陈发虎^[22]在共和切吉发现的冰缘遗迹海拔较祁连山地区现代多年冻土下界^[4]低1000m 左右； 徐叔鹰等^[66]在祁连山地区发现了大量LGM时期的砂楔、 石海等遗迹。比较不同纬度处的多年冻土下界，符合现代多年冻土下界每降低1个纬度升高200-300m的变化规律^[66]。

根据重建的末次盛冰期多年冻土分布图，青藏高原多年冻土的面积扩大到2.71×10⁶km²，是现在的195 ％ 。

(2)全新世大暖期

全新世大暖期青藏高原气温高出现代3-4℃，部分地区达到5℃^{[49, 68]}，按照大暖期高原平均升温4℃的思路对这一时期青藏高原的多年冻土界线进行了恢复( 图4)。多年冻土大幅度退缩，多年冻土下界抬升。大片多年冻土区由于其海拔较高、 气温较低，多年冻土比较稳定； 岛状多年冻土区北部靠近高海拔地区比较稳定，南下界受温度影响较大，在增温趋势下向北剧烈收缩； 高山多年冻土的变化最为剧烈，喜马拉雅山地区多年冻土在升温4℃情况下面积骤减； 显示出高原多年冻土对气候变化的响应的敏感性在边缘山区比腹部地区更为敏感。

大暖期时，大片多年冻土界线在现代年平均气温-2--4℃附近，范围较盛冰期大为收缩。柴达木盆地和共和盆地等地区的气温达到0℃以上，多年冻土消退。大暖期时的大片多年冻土北界在昆仑山南麓和柴达木盆地南缘(36.01°-37.78°N)，南下界则在噶尔-改则-唐古拉山-昌都一线(31.67°-33.75°N)，东下界到达 34°N左右，海拔约为4400m，西下界与现代范围大致相当。王绍令和李树德^[65]在风火山南麓82道班房后发现了砂楔，证明多年冻土曾在该地区发育； 青藏公路沿线西大滩8号钻孔在4.4m深处腐殖质层年龄为 7530±300a B.P. ，纳赤台Ⅰ级阶地中的灰烬状泥炭质砂土年龄为 4910±100a B.P. ，昆仑河两岸纳赤台至西大滩多处发现人类用火的灰烬，说明这一带的气候和环境条件适宜人类居住^[5]。唐古拉山以南，109道班和120道班腐殖质-泥炭沉积结束的时间分别为 5058±443a B.P. 和 4363±178a B.P.^[5]，很好地指示了大片多年冻土在大暖期大面积退化。

岛状多年冻土区北部主要是松潘高原和横断山脉北部地区，大多数地区海拔都在4300m以上，在升温的条件下依然比较稳定。南部阿里高原、 念青唐古拉山山区受温度影响较大，在增温趋势下向北剧烈收缩，其南界平均在 30°N左右，多年冻土下界可达5000m，比现代升高1000m以上( 图4)。李世杰等^[69]根据若尔盖RH-IV孔冰缘证据和年代数据指出，若尔盖地区在 8000-8500a B.P. 和 2400a B.P. 之后出现冷期，之间为暖期。鄂陵湖北岸和黄河源区两岸的第二级阶地曾发现圆滑锅状的冰楔假型^[30]证明岛状多年冻土东下界达到 97°-99°E。 120道班、 当雄乌马曲^[5]、 羊八井七弄尕^[70]泥炭沉积结束的时间分别为4363±178a B.P.、 3575±80a B.P. 和 3050±120a B.P.，这些均为大暖期时段的产物，间接说明高原多年冻土已大面积退化。

本文的重建和分析表明，在升温4℃情况下，喜马拉雅山地区多年冻土面积骤减，只在 27.50°-29.30°N范围内零星分布，且平均海拔在5100m以上； 帕米尔高原的高山多年冻土退缩程度略小于喜马拉雅山地区，在 35.5°-38.2°N内少量分布； 祁连山地区多年冻土同样大幅减少，分布于 37.7°-38.2°N之间； 青海南山、 岷山地区(34.75°-37.80°N) 、 横断山脉在末次盛冰期的高山多年冻土退缩最严重，仅发育在少数高海拔地区。野马滩地区^[30]、 若尔盖盆地北缘^[22]发现的冰楔假型及共和盆地达连海砂楔^[66]的年龄分别为 5430±410a B.P. 、 8500±180a B.P. 和 7700±90a B.P. ，证明该时期后高山多年冻土剧烈退化。

根据重建的全新世大暖期多年冻土分布图，青藏高原多年冻土的面积扩大到1.02×10⁶km²，是现在的73 ％ 。

根据青藏高原气象站点的现代气温数据和前人研究资料绘制了青藏高原现代多年冻土分布界线，并重建了末次盛冰期和全新世大暖期两个时期多年冻土的分布，对比分析后得出以下结论:

末次盛冰期时期多年冻土较现代冻土面积扩大，达到现代多年冻土面积的195 ％ 左右。其中大片多年冻土南北界扩张明显，北界可达 37°-39°N，东西界分别达到 98°E和 74°E； 岛状多年冻土受到大片多年冻土和高山多年冻土扩张的影响导致发育空间被压缩，北界与大片多年冻土南界重合，东界扩张至 100°E； 喜马拉雅山地区多年冻土向北扩张至 29°N，祁连山地区多年冻土在北、 南、 东方向有大范围扩张趋势。

全新世大暖期多年冻土强烈收缩，仅是现代多年冻土面积的73 ％ 。大片多年冻土区由于其海拔较高气温较低，多年冻土比较稳定，其北界(36.01°-37.78°N) 向南移动但幅度较小，南界位于 31.67°-33.75°N； 岛状多年冻土北部靠近高海拔地区比较稳定，南界受温度影响较大北移至 30°N附近，多年冻土下界海拔5000m； 高山多年冻土面积大幅减小，仅在喜马拉雅山山脉、 祁连山山脉、 帕米尔高原和横断山脉等部分高海拔地区发育。

致谢 感谢北京大学城市与环境学院赵昕奕老师在气象资料及分析中提供的帮助和指导。特别感谢匿名评审人的修改意见，感谢杨美芳老师在文献和文字方面的帮助。