植被作为陆地生态系统的主体, 也是气候系统中重要的、动态的组成部分^[1-2]。陆地植被可以通过改变地表反照率及陆地表面净辐射在感热、潜热与地面热通量之间的分配^[3]；或通过改变表面粗糙度^[4]来调节大气动力特征, 进而调节地-气间的水分和能量循环, 影响地表面能量平衡, 从而影响气候特征。目前在局地和区域的尺度探寻植被对气候的影响过程与机制方面已经开展了众多研究^[5-6]。但是, 由于地球表面物理、化学和生物特性存在明显区域差异, 植被对气候要素影响的机制也可能会因区域不同而存在差异。例如, 砍伐森林可能导致热带地区的气温升高, 但在寒冷地区气温反而下降^[7]。因此, 为了进一步认识植被与气候因子之间的关系, 有必要再在不同地区和不同时空尺度下分析植被分布和变化过程、地表参数及气候对其的响应。

黄土高原是中国气候变化的敏感区域之一, 其独特的自然气候状况以及生态环境变化受到中外学者广泛关注^[8-9]。为缓解中国北方地区的水土流失和土地退化等问题, 自1999年开始在黄土高原开展了大规模的退耕还林(草)等生态工程, 植被增长效果显著^[10]。目前对黄土高原退耕还林(草)政策下的农户特征变化^[11]、经济效益^[12]研究较多, 而对于该政策区域环境效应, 尤其是植被恢复后其与大气之间的相互作用的探讨较少。郭建侠等使用中尺度气象模式耦合陆面植被模式模拟退耕还林对局地环境要素的影响^[13]；王平等基于CLM模型对黄土高原植被覆盖变化对气温和降水的影响作了模拟研究^[14]；梁玲等使用中尺度气候模式模拟了黄土高原植被改变对夏季气候影响^[15]；张静辉等使用区域气候模式研究黄河源区植被覆盖度对区域气候影响^[16]。但以上研究大多基于模型模拟的方法, 较少使用观测数据证明局地或区域尺度上植被恢复对气候的反馈作用。对于黄土高原来说, 基于观测数据的退耕还林(草)工程是否可能导致区域气温变化尚未清楚。本文利用遥感卫星反演的植被覆盖数据, 结合54个位于黄土高原还林(草)工程区的气象站资料及EIN-Interim热通量数据, 采用空间分析交叉验证及地表热量平衡分析的方法, 从站点尺度探寻三者间的可能关系, 分析植被变化对气温影响的基本过程与其物理机制。

1 研究区与数据

黄土高原地区位于黄河中上游地区, 占地面积大约62.4×10⁴ km², 介于33°43′—41°16′N, 100°54′—114°33′E之间, 包括太行山以西, 青海日月山以东, 秦岭以北阴山以南广大地区^[17](图 1)。该区土壤深厚, 主要土壤类型为褐土、黑垆土和黄绵土, 这些土壤毛管孔过大, 易于蒸发, 保水性较差^[18]；而且该区为典型大陆性季风气候, 大部分地区平均年降水量400 mm左右, 空间分异较大, 而且季节性分布明显, 约有69.33%降雨量以暴雨形式集中在6—9月。黄土高原的独特地形地貌结构及其土质疏松和暴雨频发等特征, 导致该地水土流失十分严重, 是我国生态脆弱区和水土保持重点区。

本文主要采用基于地面气象站的气温观测数据、基于遥感卫星反演的植被NDVI数据与EIN-Interim再分析资料3类研究数据。所用地面实际观测资料来自中国气象数据共享网1998—2000年、2008—2010年黄土高原站点的气候资料月值数据集。

在黄土高原及周边地区选取了54个气象站点, 选取原则为：(1)站点气温数据连续, 在1998—2000年、2008—2010年5月至9月无缺失情况；(2)站点位于还林还草工程的覆盖区, 保证站点周围在退耕还林工程后期主要土地类型为草地或林地(绿地占比大于33.33%)。选取的站点分布如图 1所示, 气象站10 km绿地比重是根据黄土高原2005年1:25万土地利用数据^[19]在ArcGIS空间分析平台中计算得到。

NDVI数据集由SPOT卫星搭载的植被传感器SPOT VGT反演得到, 空间分辨率为1 km, 使用每月10 d最大值合成的NDVI产品(VGT-S10)。1998—2008年VGT-S10数据来源于“黑河计划数据管理中心”(http://westdc.westgis.ac.cn)；2008—2010年使用比利时VITO研究所图像处理与存档中心提供的VGT-S10数据产品(http://www.vito-eodata.be)。本文选取1998—2000年、2008—2010年ERA-Interim月平均地表热通量资料。ERA-Interim数据集的水平分辨率为0.125°×0.125°, 垂直分辨率60层, ERA-Interim再分析资料同化方案、数据源等参照相关文献^[20]。

2 方法 2.1 数据处理

气温数据为植被生长季(5—9月)54个气象站的月平均气温(T_mean, Mean temperature)、月最高气温(T_max, Maximum temperature)与月最低气温(T_min, Minimum temperature)的均值。NDVI月值为基于VGT-S10产品的DN值转换的NDVI数据的月最大值合成值。选取ERA-Interim数据集的植被生长季感热通量(SH, Sensible heat flux)、潜热通量(LH, Latent heat flux)与大气下行长波辐射(LD, Atmospheric downward longwave radiation)3个热通量值分析地表热量平衡状况。计算生长季热通量均值并使用克里金插值的方法估计54个气象站的辐射量值。所有的站点气温值、NDVI值及热通量采取站点周围9个像元平均值。

2.2 气温、NDVI的时空变化特征

对1998—2000年和2008—2010年两个时段的植被生长季NDVI(NDVI in growing season, 以下记为NDVIg)和气温变化的空间格局进行了分析。同时利用ArcGIS 10.0的空间分析工具计算NDVI和气温变量的最大值、最小值、平均值和标准差, 分析其时间变化特征, 并分析了NDVIg和气温的变异系数(CV, Coefficient of variation)。

2.3 植被变化与气温、热通量变化空间交叉验证

为了量化生长季(5—9月)黄土高原植被生长变化对气温的影响, 本文选取两个典型时段进行差值分析：1998—2000年(黄土高原退耕还林(草)工程初期, 以下简称S1)与2008—2010年(工程开展10年后, 以下简称S2)。首先计算S1、S2两时段间平均年NDVIg、气温变量及热通量的变化量, 然后通过ArcGIS 10.0的空间分析工具对变化量进行叠加分析, 以检测植被覆盖度与气温和地表热通量之间的关系, 并探寻其影响机制。

3 结果与讨论 3.1 植被时空变化

在退耕还林(草)开始的10年间, 黄土高原NDVIg整体增加(表 1), 其中S1和S2期时段NDVIg均值分别为0.38和0.45, 增加幅度高达18.42%；同期, NDVIg最大值也增加了13.16%。而且, S1时段NDVIg的CV值比S2高, 表明S1时段区域内单元的NDVIg差异大, 分布不均, 而S2时段的区域内差异较小。

两个时段NDVIg及两时段变化值的空间分布差异明显(图 2), 植被恢复区域明显增加。无植被覆盖或植被较少覆盖区域(0＜NDVIg≤0.3)从S1时段的总面积占比39%下降至S2时段的25.10%；植被高覆盖区域(NDVIg≥0.6)的面积从5%增加至10.32%(图 2a、图 2b)。黄土高原中部、东部、东南部的绿地覆盖面积增加较快, NDVIg增加的区域占到全区的92.37%；强烈的恢复过程伴随着有小部分区域的植被退化过程(图 2c), 在兰州以北较为明显, NDVIg下降最低值达到-0.28。退化区域可能受到自然状况与人类活动的共同作用, 该区域地处半干旱偏旱区, 土壤瘠薄, 水资源短缺, 生态治理难度较大^[21]；可能存在的过度放牧也会对干旱区风沙草原植被退化有着重要的驱动作用^[22]。

3.2 站点尺度植被覆盖与气温变化

统计得到54个站点S1、S2时段NDVIg与气温变量统计值如表 2所示。站点尺度的植被覆盖变化与像元尺度的结果一致(表 2、图 2), NDVIg各统计值在S2时段均出现上升；同时各气温变量最大值与均值下降。其中, 各站点的NDVIg最大值增加幅度最大, 达0.07(表 2)。研究时段内, 生长季各气温变量T_mean、T_max、T_min的均值与最大值都下降, 表明降温是该区的一致趋势, 尤以T_max最为明显。与之相反的是, S1至S2时段, 统计的气温变量的最小值都升高。这说明黄土高原气温状况仍受到大尺度全球变暖背景的影响, 使得气温变量最小值出现上升；但其作用不如局地强迫明显, 黄土高原气温总体上呈下降趋势。而且, S1时段NDVIg与各气温变量的CV值比S2高, 表明S1时段区域内单元的NDVI差异大, 分布不均, 而S2时段的区域内差异较小。进一步对两个时段NDVI以及气温指标的均值统计量进行99%置信区间的配对样本T检验, NDVI、T_mean、T_max的均值的在P=0.01水平上存在显著差异, T_min值在两时段差异不显著。

研究时段内NDVIg变化(ΔNDVIg)、T_mean变化(ΔT_mean)、T_max变化(ΔT_max)与T_min变化(ΔT_min)的空间分布如图 3所示。图 3表明气温变化的空间格局与NDVIg变化大致呈相反模式, 植被覆盖增加(减少)的站点气温下降(上升)。在退耕还林(草)开始的10年间, 黄土高原中部偏南靠近黄河主干道的站点NDVIg增长幅度远大于其他地区站点, 最大的在离石站, 增幅达到0.19。研究时段内黄土高原植被覆盖整体呈增加趋势, 少数站点的NDVIg出现下降, 最大值出现在黄土高原东北部的集宁站, 平均NDVIg下降0.04。黄土高原西北部站点的NDVIg值变化较小。

植被覆盖增加地区的气温下降现象比较普遍(图 3), 各站点T_mean值在中部偏南靠近黄河主干道下降情况较明显, 在阳泉站达到最低-1.1 ℃；西南部有少数变暖现象, 变暖幅度为0—1 ℃, 最高的增温出现在佛坪站0.79 ℃(图 3b)。ΔT_max值分布与ΔT_mean值大致相同, 西南部升温中部降温, 分别在玛曲站(0.82 ℃)与榆林站(-1.95 ℃)达到最高(图 3c)。ΔT_min空间分布模式与ΔT_mean与ΔT_max分布模式存在差异, 其在中部地区有少数变暖站点, 幅度在0—0.5℃, 变暖站点仍集中分布于西南地区, 最高的增温出现在佛坪站达1.25℃；降温最低的站点西移至固原北部的盐池站, 达到-1.54 ℃(图 3d)。黄土高原及周边植被覆盖站点的日最低气温空间分布变异性高于日平均气温与日最高气温。

对研究时段内平均NDVIg、气温变化的交叉验证分析表明, 植被覆盖增加与气温变量降低在空间上呈正相关, 且植被对气温的影响在日间较夜间大(图 3、图 4)。在研究时段内T_mean、T_max的第Ⅳ象限站点数最多, 占黄土高原的大部分区域, 分别约占64.81%、72.22%；T_min第Ⅰ、Ⅳ象限站点数大致相同, 分别约占44.44%、46.30%。这表明植被对气温的影响在日间较夜间大, 可能是因为植被与大气间的大部分能量交换活动都在有光照的日间。站点数在第Ⅲ象限空间上呈极少数分布, T_mean、T_max与T_min站点数分别占5.56%、3.70%和5.56%。

3.3 地表热量平衡分析

上述结果表明, 植被覆盖增加可能会对黄土高原的局地热量状况造成影响, 带来一定的降温效应。植被可以通过改变下垫面的反照率、粗糙度及蒸散等方式影响大气, 从而影响地面环流和能量平衡, 引起局地气温的变化^[2-3]。为进一步验证研究结果, 使用热量平衡分析的方法对植被覆盖变化对气温影响的基本过程及其可能的物理机制作一些探讨, 选择影响地表能量平衡关系的物理量^[23]：地-气之间的感热(ΔSH)、潜热通量(ΔLH)与大气下行长波辐射(ΔLD)。另外, 还对研究时段内生长季植被覆盖度变化分级为3类：高绿化站点(ΔNDVIg＞0.1)、中绿化站点(0.05＜ΔNDVIg≤0.1)与低绿化站点(0＜ΔNDVIg≤0.03), 定量统计其对应气温与热通量变化值如表 3所示。

各气温变化与ΔNDVIg变化高度一致(表 3), 气温变化对应NDVIg绿化区间呈梯度分布, 即随着绿化程度的增加, 各气温变量的降温程度也相应增加。T_max变量在中绿化与低绿化站点出现少量增温。除ΔLD在低绿化站点下降幅度较中绿化站点高外, 热通量变化与ΔNDVIg变化相一致, 随着绿化程度增加, ΔLH增加的幅度增加, ΔSH下降的幅度增加, ΔLD下降的幅度增加。

为探讨ΔNDVIg与气温变化之间的物理关联, 交叉验证了热通量变化与NDVIg变化之间的关系。由图 5可得, 在退耕还林(草)开始的10年间, 黄土高原植被覆盖增加与潜热通量增加、感热通量与下行长波辐射下降在空间上呈正相关关系(表 3, 图 5)。植被覆盖增加后, 地表潜热通量由于区域土壤含水量的增加、地表粗糙率及叶面积指数的增大及蒸散发的增加而出现了显著增加, 感热通量减少。如图 5所示, 约有64.81%的站点植被覆盖与潜热通量同时增加, 增加幅度为0.46—44.77 W/m²；约为74.07%的站点植被覆盖增加但感热通量减少, 减少幅度为-0.22—30.53 W/m²。

长波辐射变化可以从植被地面反照率方面进行解释, 一般而言, 植被的反照率比荒漠地要低得多^[24]。黄土高原地区植被覆盖增加后引起地表反照率的改变, 会对地表吸收短波太阳辐射和发射长波辐射造成影响。许多研究证明, 植被覆盖的增加, 可以通过减少地表反照率, 增加地表面吸收的太阳辐射量, 从而产生对气温的正反馈^[25]。然而, 这种正反照率反馈可能部分被更高的蒸散发的冷却作用所抵消^[26]。这种蒸散发带来大气中水汽增加, 使得云量增加, 会造成下行长波辐射的明显减少, 如图 5所示, 减少的站点分别占黄土高原的79.63%, 幅度在-0.02—13.91 W/m²。

从热量平衡的角度分析, 植被增加使得土壤含水量增加, 表面阻抗增加, 植被的蒸腾增加, 以及粗糙度增加而造成的湍流垂直扩散增强, 会造成潜热通量增加, 感热通量的减少。感热通量减小幅度小于潜热的增加幅度, 故该区的有效热通量(感热+感热)是增加的, 数值在0.24—14.24 W/m²。再加上反照率减少, 吸收的净短波辐射增加。通过分析热量平衡分量变化, 由于地表净吸收的辐射增加, 且有效热通量的值增加, 要在更低的气温下才能使地表吸收的短波、发射的长波、输送的有效热通量达到新的平衡。通过植被蒸散作用加强并可能导致降温的结果与梁玲等^[15]、苟娇娇等^[27]在黄土高原进行的植被变化对气候影响试验数值模拟研究一致。这与北极地区植被增长反馈作用分析^[28]及中国的土地覆盖变化敏感性试验研究结果^[29]不同, 他们认为植被增长使得地表反照率减少, 增加了表面吸收的太阳辐射量, 进而导致植被对气温的正反馈。然而, 从局地尺度上来看, 黄土高原气候冬季严寒、夏季暖热, 且大量覆盖草地覆被。植被的蒸散发作用在该地区有着重要的地位, 其可以在植被增加的情况下削弱由于反照率减少导致的变暖, 植被需要蒸腾更多的水分, 从而在生长季持续降温反馈。植被恢复可能会对黄土高原气候变暖现象有减缓作用。

4 结论

本文基于气象观测数据与卫星遥感NDVI数据, 对黄土高原生态工程实施前后10年间(S1：1998—2000年, 工程初期；S2：2008—2010年, 工程10年后)的植被增长状况与气温变量时空变化状况进行交叉验证, 探讨黄土高原植被恢复对气温的反馈作用；同时基于EIN-Interim再分析资料提取的地面观测热通量值变化, 分析植被变化对气温影响的基本过程及其可能的物理机制, 探讨植被与大气间相互作用。通过研究得出以下结论：

(1) 黄土高原植被恢复比较明显, 生长季植被增加显著；与之对应的是生长季3个气温变量(T_mean、T_max、T_min)的均值与最大值在S2时段都下降, 尤以T_max最为明显。

(2) 植被覆盖增加与气温变量降低在空间上呈正相关。同时, 植被覆盖增加与潜热通量增加、感热通量与下行长波辐射下降在空间上也呈正相关关系, 植被恢复通过加剧植被蒸散作用导致近地面降温。

研究还存在一定不足, 数据的不确定性和模型缺陷会阻碍对植被活动增加引起的气温反馈大小的定量结论。未来应考虑加大观测数据的测量时空精度、在模型模拟中引入其他参数(如叶面积指数、表面反照率)使陆面过程参数化更符合实际状况。这些改进预期能够更好地了解植被在气候系统中的作用。