沙漠作为地球上生态环境十分脆弱的地区之一，其植被变化对于沙漠发展至关重要。特别是在沙漠边缘地区，植被对抑制风沙活动和沙漠化具有举足轻重的作用，植被生长状况直接影响着环境的优劣^[1~3]。同时，植被作为连接土壤、大气和水分的自然“纽带”，在全球变化研究中起到“指示器”的作用^[4]。归一化植被指数(Normalized Difference Vegetation Index，简称NDVI)作为评价植被覆盖度和植被状态的有效参数，是监测全球植被和生态环境变化的有效指标^[5]，在植被-气候响应^[6]、地表植被动态变化^[7]和生物量估算^[8]等方面已经被广泛应用。利用遥感影像可以定期、及时获取研究区的NDVI信息，但由于受到遥感影像资料年限的制约，20世纪80年代之前的NDVI数据无法通过这种手段直接获取，使得长时间序列的植被变化研究也受到一定程度的限制。

树木年轮作为一种分布广泛、分辨率高、时间序列长的代用指标，已被应用于重建过去气候^[9~13]、水文^[14~16]、积雪^[17]和冰川^[18]变化过程。随着树轮年代学的发展，有不少学者通过分析树轮宽度与NDVI的关系，将树轮资料成功应用于长时间序列NDVI的重建^[19~21]。腾格里沙漠南缘山地虽临近沙漠，但分布着罕见的天然针叶林。近10年来，在本区已经陆续开展了树轮气候学研究，建立了多条降水^[22~23]、温度^[24]和Palmer干旱指数(Palmer Drought Severity Index，简称PDSI)^[25~26]序列。这些研究丰富了对本区历史时期气候变化过程的认识。腾格里沙漠南缘作为我国土地沙漠化较为严重的地区^[27]，其植被变化对于本区沙尘天气的发生和发展具有重要影响，但之前对本区历史时期植被生长状况研究还相对较少。本文利用采自腾格里沙漠南缘哈思山的油松树轮资料，建立树轮宽度年表，通过分析树轮宽度年表与NDVI的关系，重建研究区NDVI变化序列，并深入分析本区植被长期变化特征，揭示其与沙尘天气的关系，为本区植被恢复和生态环境改善提供重要的科学依据。

哈思山(图 1，36°58.33′~37°02.67′N，104°18.67′~104°35.00′E)位于甘肃省靖远县境内，属于祁连山系。山体位于贺兰山以西的腾格里沙漠南缘地区，主要由变质岩、花岗岩、闪长岩组成，地表多为第四纪风成黄土覆盖，部分地区岩石裸露，平均海拔2400~2700 m，最高峰3017.8 m^[28~29]。本区位于中国北方干旱半干旱区，受地形的影响，年降水量高于其他地区，为350~410 mm，降水多集中在7~9月，占全年降水量的75%以上，年平均气温3~10 ℃，平均最低气温-12.3 ℃，出现在1月，平均最高气温28.7 ℃，出现在7月^[29~30]，盛行西北、西西北和北西北风^[31]。哈思山林区植被覆盖较好，但森林植被类型较少，且层次结构简单、树种组成单纯^[28]，生长着甘肃中东部干旱区罕有的天然油松纯林，被誉为“沙漠里的绿洲”^[29]。

图 1(Fig. 1)

图 1 树轮采样点位置图 Fig. 1 Location of tree-ring sampling site

哈思山地区油松(Pinus tabulaeformis)天然林垂直分布在海拔2200~2700 m的陡峭阴坡及半阴坡面，在半阳坡多以小片状与山地草原呈复合状分布，其分布处于我国温带针叶林西北部及其邻近的蒙新高原半荒漠干旱草原区^[29]。2010年8月对本区油松进行采样，每棵树取两个钻芯作为复本，共采集55棵树，114芯。

采集的样品经过实验室晾干、固定、打磨预处理后，在北京师范大学环境演变与自然灾害教育部重点实验室利用LINTAB树轮宽度仪对处理好的样本进行宽度测量，并利用骨架示意图交叉定年，最后用COFECHA程序^[32]对交叉定年结果进行检查。在建立树轮宽度指数序列前，剔除那些年代序列短、与主序列相关较差的样芯。利用ARSTAN软件^[33]以负指数或线性函数拟合树木生长趋势，然后将实测树轮宽度值除以每年生长的期望值进行年轮宽度指数化，再以双权重平均法合成平均树轮宽度年表。

本文所采用的气象数据来自距离哈思山最近的景泰气象站(37°10.8′N，104°03′E)，时间区间为1957~2007年。所采用的NDVI数据来自美国国家航天航空局的GIMMS(global inventory modeling and mapping studies)NDVI数据集。此数据集包括从1981年7月到2006年12月的全球植被指数变化，时间分辨率是15天，空间分辨率0.22°×0.22°。从该数据集中提取覆盖研究区(35.99°~37.08°N，104.08°~104.95°E)的逐月平均NDVI值来表示该区的植被覆盖情况，建立本区1981年7月~ 2006年12月的月NDVI数据序列。

哈思山地区多年月平均NDVI变化见图 2，从中可以看出，其多年月平均NDVI变化表现为明显的单峰型，变化范围为0.10~0.22，低于祁连山东部山区平均水平^[34]，说明本区植被相对较为稀疏，植被覆盖度总体较低。这与哈思山脱离祁连山山脉主体，相对孤立，同时又靠近沙漠有关。从全年来看，本区NDVI在12~3月处于低值期，从4月开始增加，8月达到最大值，9月开始降低，到11月NDVI值仅为0.12。春季是植物生长季开始的季节，能直接表现区域植被覆盖的分布变化，而夏季的NDVI，则反应了植被生长状态的优差^[35]。从NDVI的季节变化情况可以看出，本区植被从4月份开始生长，随着植物的生长，地表植被覆盖度逐渐增加，在7~8月植物生长最为繁茂，而从9月开始植物开始枯萎，10月以后基本停止生长。将NDVI与景泰气象站记录的多年月平均气温(图 2a)和降水(图 2b)数据对比可以发现，NDVI月变化滞后于温度变化，而与降水的变化更为一致。研究区4~5月气温迅速升高，但由于此时植物处于生长期初期，植被生长缓慢，NDVI值仍处于较低水平。7~8月水热条件最为适宜，NDVI也达到一年中的最大值。

图 2(Fig. 2)

图 2 NDVI与温度(a)和降水(b)的月平均变化 Fig. 2 Monthly variations of mean temperature (a) and precipitation (b) with NDVI

树轮宽度变化不仅与当年气候条件有关，而且还受到生长前期气候的影响^[36]。本文分析了树轮宽度标准年表(STD)与前一年9月至当年12月降水和温度的关系，相关分析结果显示(图 3)，哈思山油松树轮宽度变化仅与前一年9月降水呈显著正相关，当年任何月份降水对油松径向生长的限制作用均不明显；与前一年9月、前一年11月和当年1~3月、当年11月温度呈显著负相关，说明秋末和冬季温度对油松径向生长有较强的限制作用。前一年9月丰富的降水使树木体内可以储存足够的水分供来年生长所需^[37]。秋末和冬季过高的温度，会使树木呼吸作用和代谢作用增强，加速体内储存物质的消耗，对植物的次年生长不利^[38]。

图 3(Fig. 3)

图 3 树轮宽度与温度和降水的相关系数 -9~-12指前一年9~12月，虚线表示通过0.05的显著性检验 Fig. 3 Correlation coefficients between tree ring chronology and monthly mean temperature and monthly total precipitation. -9~-12 indicates the September to December of previous year, the dotted lines represent statistically significant level(P=0.05)

利用相关分析法分析了STD和差值年表(RES)与当年各月份NDVI的关系，发现树轮差值年表与NDVI的相关性更好(表 1)，这与王文志等^[21]在整个祁连山山区的研究结果一致，因此，在本文的后续研究中选择差值年表进行分析。相关分析结果显示，树轮宽度年表与9月NDVI相关性最强，相关系数为0.58，与5~11月NDVI呈正相关关系，且与6~9月NDVI相关性在0.01置信区间显著相关，而与冬季(12~2月)的NDVI呈负相关关系。

表 1(Table 1)

表 1 树轮宽度年表与各月NDVI相关系数 Table 1 Correlation coefficients between tree ring chronology and NDVI 月份 1 2 3 4 5 6 7 8 9 10 11 12 6~9 STD 0.208 0.047 0.078 -0.108 0.405* 0.544** 0.487* 0.430* 0.444* 0.287 0.199 0.260 0.562** RES -0.099 -0.085 -0.015 -0.227 0.380 0.562** 0.502** 0.573** 0.583** 0.326 0.067 -0.033 0.670**   *通过0.05的置信水平   **通过0.01的置信水平

表 1 树轮宽度年表与各月NDVI相关系数 Table 1 Correlation coefficients between tree ring chronology and NDVI

基于树轮差值年表与NDVI的相关分析结果，建立了两者的定量关系，拟合线性回归方程如下：

其中，NDVI_6~9是6~9月平均NDVI值，RES为差值年表，方程方差解释量为44.9%，调整后方差解释量为42.7%。图 4为1981~2006年NDVI实际值与重建值的对比，可以发现虽然重建序列离散程度要小于实际值，但多数年份重建值与实际值比较接近，两者变化趋势一致。为保证重建序列的可靠性，对重建序列采用多种方法进行检验，包括符号检验法、乘积平均值检验和误差缩减值检验^[37]。低频符号检验(5-/21+)、高频符号检验(7-/19+)和乘积平均值(t=2.18)检验都达到了95%的置信水平，误差缩减值(RE=0.45)大于0，这些检验结果都验证了重建序列的可靠性。基于该线性回归方程，我们重建了1804年以来腾格里沙漠南缘6~9月的NDVI序列(图 5)。

图 4(Fig. 4)

图 4 NDVI实际值与重建值比较 Fig. 4 Comparison between the actual and reconstructed NDVI

图 5(Fig. 5)

图 5 1804年以来6~9月NDVI重建值和5年滑动平均序列 Fig. 5 Reconstructed June-September mean NDVI since 1804 and a 5-year moving average of NDVI reconstruction

重建NDVI序列多年平均值为0.21，标准差为0.02。NDVI最低值出现在1966年，仅为0.16，最高值出现在1804年，为0.27。将NDVI值介于Mean(平均值)±σ(标准差)之间的年份作为植被生长正常的年份，而低于Mean-σ的年份作为植被生长较差的年份，高于Mean+σ的年份作为植被生长较好的年份^[39]。依据这个标准可以发现，从1804年到2009年植被生长较好的年份共24年，占总年份的11.65%，植被生长较差的年份共27年，占总年份的13.11%。为了更好地分析重建NDVI序列的低频变化信息，对NDVI序列进行Mann-Kendall趋势分析和5年滑动平均分析，结果显示，1804~2006年腾格里沙漠南缘NDVI虽然经历了多次波动，但植被没有显著的增加趋势(Z=0.91，介于± 1.96之间)，其中植被生长较好的时期为1826~1829年、1849~1852年、1885~1895年、1934~1963年和1969~1978年，植被较差的时期为1815~1823年、1831~1847年、1897~1909年、1924~1932年、1988~1993年和2001~2006年。

在沙漠及其邻近的干旱与半干旱区，由于其特定的地理环境和下垫面条件，往往会发生沙尘暴^[3]。以往的研究表明，西北地区植被覆盖率和沙尘暴发生频率之间存在负相关关系^[40]，且这种相关关系在不同地貌类型区和不同季节有所差异^[41]。利用中国气象数据网提供的中国强沙尘暴序列，我们还分析了1954~2007年腾格里沙漠南缘景泰和靖远两个地区的沙尘暴发生日数与本区NDVI序列的关系，可以发现NDVI极低值年份与沙尘暴多发期有良好的对应关系(表 2)，特别是20世纪70年代之前，在NDVI值极低的年份，景泰和靖远两个邻近地区的沙尘暴发生日数都高于多年平均值；70年代之后，即使在NDVI极低值年份，沙尘暴发生频率也明显减少，这主要与区域气候背景有关。全林生等^[42]通过分析中国北方地区83个站平均的沙尘天气(沙尘暴)频次与年降水量、冬季月平均气温及春季气旋频次的逐年变化发现，受到春季气旋活动减少和冬季气温变暖的影响，20世纪70年代以后中国沙尘暴发生日数呈现明显减少的趋势。景泰和靖远两个站点记录的1985~2007年沙尘暴平均发生日数减少也为3天和1天。1990年、1991年和2000年这3年沙尘暴发生日数虽然较前期显著降低，但与后期其他年份相比仍为沙尘暴多发期。在干旱半干旱区，地表植被变化不仅能反映区域温度和降水的变化，也会引起区域土壤条件和地表粗糙度的改变，因而成为影响沙尘天气的一个重要因子。已有研究发现，植被覆盖状况是影响中国北方沙尘暴频次变化的重要因素^[43]。但由于沙尘天气的发生还受到其他因素的影响，如天气系统^[42]、下垫面性质^[44]等，而地表植被作为影响因子之一，其对沙尘暴的发生所产生的影响是有限的。魏桂英等^[45]通过分析腾格里沙漠南缘沙尘暴气候变化特征，发现年大风日数对沙尘暴的影响程度最强，而植被条件次之。矫梅燕等^[46]也认为大气的动力条件决定着沙尘天气的强度和发生范围，而沙源区地表特征的变化对沙尘天气的发生起着重要作用。

表 2(Table 2)

表 2 NDVI低值年份与景泰、靖远沙尘暴日数对比 Table 2 Periods with low NDVI values and sandstorm days in Jingtai and Jingyuan 年份 NDVI 景泰沙尘暴日数 靖远沙尘暴日数 1957 0.19 33 10 1966 0.16 25 11 1973 0.17 52 18 1981 0.17 16 6 1982 0.19 19 6 1990 0.17 1 3 1991 0.18 4 1 2000 0.17 4 2 多年平均值 0.21 13 6

表 2 NDVI低值年份与景泰、靖远沙尘暴日数对比 Table 2 Periods with low NDVI values and sandstorm days in Jingtai and Jingyuan

植被变化受多种因素控制，而气候作为主要的影响因子，NDVI对其变化的响应十分敏感^[47]。在干旱半干旱区由于降水稀少，水分条件往往成为地表植被生长的制约因素^[48]。本文将重建的NDVI序列与兰州旱涝等级^[49](图 6c)和本区重建的PDSI^[50](图 6d)分别进行相关分析。图 6c中纵坐标的旱涝等级序列采用1、2、3、4、5分别代表涝、偏涝、正常、偏旱和旱5个等级；而PDSI作为表征干旱强度的一个重要指标^[51]，正值表示湿润的气候条件，负值指示干旱的气候条件。结果显示(图 6c)，NDVI与兰州旱涝等级在0.01置信水平上有显著的负相关关系，相关系数为-0.29，这与王亚军等^[52]在邻近的祁连山东部余脉兴隆山地区研究结果一致，与PDSI有显著的正相关关系，相关系数为0.38。这些都说明水分条件对本区植被变化有重要影响。据《西北灾荒史》^[53]记载，1821~1847年西北地区断续发生大旱共13次，灾区主要集中在甘宁青地区，而重建的NDVI序列显示，NDVI值虽然在1826~1829年偏高，但之后的一段时间内NDVI值几乎都低于多年平均水平。1877~1878年在北方多省发生持续的干旱，史称“丁戊奇荒”^[53]，这次罕见的旱情对本区植被也产生了重要影响，NDVI值由1876年的0.24急剧下降到1877年的0.17(图 6a)。20世纪20年代中国北方普遍发生干旱事件，其中1928年、1929年和1930年这3年是毁灭性大旱，为西北历史上成灾最重的大旱灾^[53]，重建NDVI序列显示，1924~1932年也是腾格里沙漠南缘植被生长较差的一个时期。王文志等^[21]在整个祁连山地区的研究也揭示这段时期因为降水减少和持续干旱，植被生长最差(图 6b)，这也从一个方面说明哈思山虽与祁连山主体分离，山体较小，但其植被变化对区域性气候变化的响应也是敏感的。

图 6(Fig. 6)

图 6 本文重建NDVI序列与其他重建序列对比 (a)腾格里沙漠南缘1804~2009年6~9月NDVI变化曲线；(b)祁连山地区1843~2003年6~8月NDVI变化曲线[21]；(c)兰州旱涝等级[49]；(d)腾格里沙漠南缘地区重建的6~8月PDSI格点数据[33] Fig. 6 Comparison of reconstruction series in different areas. (a)The reconstructed June-September NDVI from 1804 to 2009 in the southern margin of the Tengger Desert, (b)the reconstructed June-August NDVI from 1843 to 2003 in the Qilian Mountains[21], (c)the dryness/wetness index in Lanzhou[49] and (d) the reconstructed June-August PDSI from same grid derived from the Monsoon Asia Drought Atlas[33]

(1) 气候对哈思山油松径向生长表现出明显的滞后性，树轮宽度变化主要受到前一年9月降水和前一年秋末和冬季温度的影响。

(2) 从月际变化来看，腾格里沙漠南缘NDVI变化呈现明显的单峰特征，NDVI最高值对应水热条件最好的7~8月份。从年际变化来看，Mann-Kendall趋势分析和5年滑动平均结果显示，1804~2009年腾格里沙漠南缘NDVI没有显著的增加趋势，但存在多次波动，在1826~1829年、1849~1852年、1885~1895年、1934~1963年和1969~1978年5个时期植被生长较好，在1815~1823年、1831~1847年、1897~1909年、1924~1932年、1988~1993年和2001~2006年6个阶段植被较差。

(3) 腾格里沙漠南缘NDVI变化与当地沙尘天气有一定联系，NDVI极低值年份与景泰和靖远两地的沙尘暴多发期有很好地一致性。在20世纪70年代之后，当地沙尘暴发生日数显著减少，而NDVI却没有明显的增加趋势，这说明本区植被变化对沙尘暴发生所产生的影响是有限的。

(4) 腾格里沙漠南缘NDVI序列与兰州旱涝等级有显著的负相关关系，与PDSI有显著正相关关系，1821~1847年、1877~1878年和1924~1932年发生的极端旱灾对本区植被也产生很大影响，植被生长较差，而且同一时期在整个祁连山地区也存在类似的情况，这说明本区地表植被变化受到水分条件的制约，对区域气候变化响应敏感。

致谢: 感谢唐清亮、尚媛在样品采集和处理时的帮助；真诚地感谢审稿专家和编辑部老师建设性的修改意见。