当前碳循环已经成为评价全球气候变化影响的重要指标^[1], 而植被净初级生产力(NPP)作为陆地碳循环的重要组成部分, 可有效的量化陆地生态系统的生产能力, 是全球碳循环研究的热点问题^[2]。随着气候变化加剧, 区域NPP的定量化研究受到越来越多研究者的关注, 已成为众多国际间合作研究计划(如国际地圈-生物圈计划(IGBP)、京都协定(KP))的核心内容之一^[3]。而生态脆弱区碳循环对全球气候变化响应十分敏感, 定量评估NPP可为生态脆弱区实现“碳达峰”、“碳中和”提供理论支撑, 为应对气候变化、区域生态环境保护作出重要指导。

传统的NPP估算方法包含直接观测法和间接估算法^[4]两种, 其中直接观测法精度较高, 但因其适用区间较小, 难以大范围推广^[5]。间接估算法包括模型法(统计模型、参数模型)和遥感反演法, 模型法由于其精度不高、涉及参数多、数据获取困难等缺点, 误差往往难以控制。遥感反演法集遥感数据、地理信息系统(GIS)、模型模拟三者优势于一体, 在全球及区域NPP的估算中得到了广泛应用。其中之一的CASA(Carnegie-Ames-Stanford Approach)模型具有模型简单、参数较少、多时相数据获取便捷等优势^[6], 是目前区域NPP估算中应用较广的方法。利用CASA模型, 部分研究者开展了诸多典型区域的NPP定量化研究, 如欧洲地中海流域^[7]、非洲萨赫勒草原^[8]、北美平原区^[9]、澳大利亚^[10]、亚洲的伊朗半干旱牧区^[11]以及国内的东北平原^[12]、西北干旱区^[13]、青藏高寒区^[14]、陕西甘肃地区^[15]等区域, 获得了许多区域NPP的重要信息, 为区域碳循环及气候变化研究提供了重要支撑。此外, 亦有研究者对传统的CASA模型进行改进, 通过修正原模型的植被光能利用率^[16]、水分胁迫系数^[17]、下垫面参数^[18]等参数以提高模型在各区域NPP估算过程中的适用性, 取得了众多研究成果。

吕梁连片贫困区地处黄河中游, 是我国半干旱和半湿润区域的过渡地带, 生态环境极为脆弱^[19]。自改革开放之后, 该区域对煤、铁、铝等资源的大规模无序开发导致了植被破坏, 水土流失日趋严重, 生态环境日益恶化。2000年以来, 此区域开展了一系列生态恢复工程(国家实施的三北防护林工程(4、5期)、退耕还林还草工程以及山西省实施的吕梁山生态脆弱区林业生态建设工程等), 在此影响下区域植被覆盖度得到了一定提高^[20], 但上述措施是否带来了区域NPP的变化, 目前尚不清晰, 因此亟需开展吕梁连片贫困区在大规模生态恢复措施背景下的NPP定量评估。基于此, 本文选取吕梁连片贫困区生态恢复措施实施以来(2000年)的多期MODIS遥感数据、气象数据和土地利用类型数据, 利用CASA模型定量评估了区域NPP时空分布变化特征, 系统分析了不同土地利用类型下NPP的变化状况。研究结果将有利于增加气候变化下对区域碳循环过程的认识, 同时也为区域生态环境恢复、吕梁山区精准扶贫、黄河流域高效开发利用等提供重要的理论依据。

1 研究区概况

如图 1所示, 本文研究区地处黄土高原东部, 西隔黄河与陕西相望, 位于东经110°22′-112°37′, 北纬35°53′-39°24′之间, 总面积约为3.3万km², 海拔在400-2800 m之间, 包括神池县、五寨县、岢岚县、静乐县、兴县、临县、石楼县、岚县、吉县、大宁县、隰县、永和县、汾西县等13个国家级贫困县^[21]。研究区气候属于温带大陆性季风气候, 雨热同期, 降水较少(300-700 mm)且集中于7-8月份, 年内气温为-6-24℃, 潜在蒸散发为1000-1400 mm。地形以山地、丘陵为主, 地表沟壑纵横分布, 崎岖不平。区域内部水土流失严重, 植被稀少, 天然植被以侧柏、刺槐及小灌木为主, 是我国半湿润半干旱地区过渡地带以及农牧交错带^[22]; 境内水系主要为黄河水系, 重要的河流包括汾河、昕水河、三川河、岚漪河、蔚汾河、湫水河等, 区域水资源相对匮乏, 生态系统脆弱, 环境承载力低下^[23]。作为我国革命老区, 该区域贫困程度较深、脱贫攻坚难度较大。区域交通不便, 城镇化水平处于起步发展阶段, 经济社会发展较为落后。区域日益恶化的生态环境引起了政府部门的高度重视。

自2000年以来, 为解决区域水土流失及生态环境退化等问题, 此地先后开展了一系列生态恢复工程。国家级工程有三北防护林工程(4、5期)、天然林保护、退耕还林还草等林业重点工程; 省级工程有2013年开始实施的吕梁山生态脆弱区林业生态建设工程、2018年开始实施的太行山吕梁山生态系统保护和修复重大工程; 在地市级层面, 临汾市采取多种措施创建生态经济型三北防护林工程, 具体有吉县留村800 hm²的黄土高原综合治理林业示范建设项目工程、乡宁县堡子村933.33 hm²的三北造林示范工程等, 吕梁市也实施了离石区可视范围山体绿化工程建设项目等生态恢复措施。

2 数据来源与分析方法 2.1 数据来源

本文主要应用气象数据(下载自中国气象局国家气象信息中心(http://data.cma.cn/site/index.html))、土地利用类型数据(下载自欧洲航天局(http://maps.elie.ucl.ac.be/CCI/viewer/))及遥感影像数据(来源于美国LAADS网站(https://ladsweb.modaps.eosdis.nasa.gov/)发布的空间分辨率为1 km的MODIS数据MOD13A3数据集)三种数据(表 1)。气象数据包含该区11个站点的月均温度(℃)、月均最低最高温度(℃)、月降水量(mm)、月均日照时长(h)、月均风速(m/s)、月均水汽压(kPa)和相对湿度(%)等。土地利用类型数据为研究时段内空间分辨率为1 km的土地利用类型数据, 本文将其分为耕地、林地、草地、居民地、水域和未利用土地6个一级类型。遥感影像数据主要为5个时段内60景的月合成MODIS数据, 并使用MRT软件将HDF格式文件转换为TIFF文件, 投影转换为WGS84/Albers Equal Area Conic投影, 最终经掩膜处理获得研究时段内逐月标准化NDVI数据。

2.2 研究方法

本研究选用基于光能利用率的CASA模型估算吕梁连片贫困区2000-2018年NPP, 其结构如下所示:

(1)

式中, NPP(x, t)表示t月份单位像元x内的植被净初级生产力(gC m^-2), APAR表示植物所吸收的光合有效辐射(MJ/m²), ε表示光能利用效率, t表示时间, x表示空间位置^[24]。

其中APAR参数由下式求得:

(2)

式中, SOL(x, t)表示t月份像元x处的太阳总辐射量; FPAR(x, t)表示植被层对入射光合有效辐射(PAR)的吸收比例; 常数0.5表示植被所能利用的太阳有效辐射(波长为0.38-0.71 μm)占太阳总辐射的比例。

SOL采用经验公式进行计算, 由大气上界太阳辐射量和日照百分率计算得出^[25]。

由前人研究可知, FPAR与归一化植被指数(NDVI)和比值植被指数(SR)都存在一定的线性关系, 结合两者进行计算可提高计算精度, 故参考陈福军等人^[26-27]的研究计算像元FPAR值。

模型中光能利用效率ε是指植物把吸收的光合有效辐射(PAR)转化为有机碳的效率。由下式求得:

(3)

式中, T_ε1(x, t)表示在低温和高温环境下植物自身生理作用对光合作用产生限制而降低其净初级生产力的程度; T_ε2(x, t)表示植物在生长最适温度向高温和低温转变时光能利用率逐渐变小的趋势; W_ε(x, t)表示水分状况对光能利用效率的限制程度; ε^*表示理想条件下的最大光能利用率, 本研究采用前人研究成果, 取0.389 gC/MJ。

T_ε1(x, t)和T_ε2(x, t)的计算过程见参考文献^[28]。水分胁迫系数W_ε(x, t)在环境影响下取值在0.5(极干旱)到1(极湿润)之间波动, 计算公式如下:

(4)

式中, E(x, t)为区域实际蒸散发, E_p(x, t)为区域潜在蒸散发。潜在蒸散发采用Penman-Monteith模型^[25]进行计算; 实际蒸散发由于数据难以获取, 故采用周广胜张新时模型^[29]进行模拟。

3 结果与分析 3.1 研究区NPP变化特征分析 3.1.1 区域NPP的年际变化

图 2为2000-2018年间研究区NPP年均值空间分布特征图, 该区年均NPP值介于241.24-331.70 gC m^-2 a^-1之间, 与黄土高原其他区域^[30-31]相比处于中等水平, 表明该区域植被整体生产能力一般, 生态环境条件不是十分优越。研究区NPP空间分布整体呈“中部高、四周低”的变化趋势, 中部高值区为吕梁山山地, 西部为黄河河谷, 东部为汾河谷地, 中部NPP年均值较高与该区域植被覆盖较好有一定关联性。

2000-2018年间区域NPP年均值整体呈波动上升趋势, 其变化大致可分为“上升-下降-上升”三个阶段。其中2000-2010年研究区NPP年均值呈稳定上升趋势, 区域年均NPP值由2000年的241.24 gC m^-2 a^-1上升至2010年的303.81 gC m^-2 a^-1, 年增速为6.26 gC m^-2 a^-1。2010-2015年间, 研究区NPP年均值呈整体下降趋势, 由2010年的303.81 gC m^-2 a^-1降低至2015年的258.49 gC m^-2 a^-1, 年降速为9.1 gC m^-2 a^-1。2015-2018年间, 区域NPP年均值再次快速上升, 由2015年的258.49 gC m^-2 a^-1升至2018年的331.7 gC m^-2 a^-1, 增长速率为24.4 gC m^-2 a^-1。

3.1.2 NPP的月际变化

NPP的月际变化反映了植物在年内净初级生产力的变化状况, 吕梁连片贫困区NPP逐月(月均)NPP分布特征如图 3所示:

由图 3可知, 各年内月均NPP值均呈现先上升后下降的整体趋势, 春季NPP多年平均值为48.38 gC m^-2季^-1, 夏季NPP多年平均值为158.1 gC m^-2季^-1, 秋季NPP多年平均值为62.85 gC m^-2季^-1, 冬季NPP多年平均值为10.94 gC m^-2季^-1。各月NPP均值以7、8两月为最高, 多年平均值分别为53.26 gC m^-2月^-1和65.41 gC m^-2月^-1; 冬季的1、12两月的NPP多年均值最低, 分别为3.65 gC m^-2月^-1和3.59 gC m^-2月^-1。各月间NPP增长速率最高时段为6-7月, 增长速率达到28.86 gC m^-2月^-1; 下降速率最高时段为8-9月, 下降速率达到52.39 gC m^-2月^-1。这些变化与研究区绿色植被的生长规律一致, 冬季地表植被光合作用微弱, 此时区域NPP处于全年的低值时段, 自春季开始, 作物萌芽生长, 5月开始植物生理过程变得活跃, 进入快速生长季, NPP开始显著增加, 区域NPP的变化也是该区域物候特征直接反映。研究区地处温带大陆性季风气候区, 雨热同期且集中于夏季, 水热组合条件为一年中最适宜植物生长阶段, 故NPP值在8月达到峰值; 秋季开始, 植物生长速度减缓, 进入凋落季, 农作物亦大量收割, NPP下降迅速; 至冬季植物枯萎殆尽, 生命活动基本停止, 故冬季为一年中NPP最低时期。

3.2 研究区NPP变化率空间分布特征

为了清晰地揭示研究区NPP值变化的空间分布特征, 根据以往研究^[32], 本研究将5个时段的NPP变化状况以显著减少(＜-20%)、微弱减少(-20%-0)、稳定增长(0-20%)和显著增长(＞20%)进行了划分。结果如图 4所示:

结果表明: 2000-2018年区域NPP整体呈上升趋势, 只有吕梁山东南部呈轻微下降趋势, NPP呈增长状态的区域占研究区的93.46%, 呈减少状态的仅为6.54%。2000-2005年区域北部NPP整体呈微弱下降趋势, 兴县、临县交界处以及岢岚县和五寨县的局部区域是发生下降的主要地区, 中南部呈上升趋势, NPP呈增长状态的区域占总面积的66.06%。2005-2010年间, 区域NPP值在研究区中北部发生显著变化, 由轻微下降变为显著上升, 其中有86.5%的区域NPP呈增长状态, 呈减少状态的仅为13.5%。在2010-2015年间, 研究区NPP整体下降, 这与其它几个研究时段截然不同, 其中以中南部下降最为显著; 这一时期, NPP呈减少状态的区域达到了88.95%, 生态环境在这一时期出现退化。近几年来(2015-2018年), 区域NPP整体呈显著上升趋势, 说明区域生态环境整体改善, 植被覆被度以及NPP稳定提高, 这一时期NPP呈增长状态的区域占比为99.65%。

3.3 不同土地利用下的区域NPP时空分布特征

根据以往研究, LUCC变化对区域NPP的空间分布规律具有显著影响^[33], 2000年以来在人类活动影响下研究区土地利用和覆被类型发生显著变化。为揭示不同土地利用方式下NPP的变化特征, 本文对比分析了不同时间段四种土地利用类型(耕地、林地、草地、城镇)下的NPP时空变化特征, 如图 5所示:

研究区耕地区域NPP值时空分布如图 5所示, 耕地NPP均值整体呈波动上升趋势, 多年平均值为243.31 gC m^-2 a^-1, 最大值出现在2018年, 为305.89 gC m^-2 a^-1, 最小值出现在2000年, 为199.74 gC m^-2 a^-1。耕地土地利用类型变动广泛存在于各县市之中, 在研究时段内其面积缩减了2.42%, 主要转化为林地和草地, 总体保持平稳。耕地转出为林地的区域主要集中在吕梁山中北部(方山县、临县、兴县交界处, 中阳县)以及大宁县、吉县交界处, 耕地转化为草地的区域主要集中在研究区的东北角、西南角以及中部区域。

过去18年来, 研究区林地区域NPP变化整体呈波动上升趋势(图 5), 多年平均值为425.09 gC m^-2 a^-1, 最大值出现在2005年, 为451.74 gC m^-2 a^-1, 最小值出现在2015年, 为379.51 gC m^-2 a^-1。林地土地利用类型变化以转入为主, 面积在研究时段内增长了6.84%, 整体变化呈现稳定上升趋势。变化以中北部最为明显, 增长源主要来自于耕地以及草地。草地转化为林地的主要集中于吕梁山中部的离石区、中阳县以及方山县、临县、兴县交界处。

研究区草地NPP变化整体呈波动上升趋势(图 5), 其中研究区中部变化最为显著, 其多年平均值为272.55 gC m^-2 a^-1, 最大值出现在2018年, 为326.1 gC m^-2 a^-1, 最小值出现在2000年, 为226.71 gC m^-2 a^-1。草地利用类型在研究时段内整体变动较大, 转入和转出都较为明显, 草地转出为耕地的区域主要集中在岢岚县、五寨县交界处、石楼县、交口县交界处以及蒲县, 以蒲县最为明显; 草地转出为城镇用地的区域主要集中在各县市城区, 其中以离石区、柳林县、中阳县最为明显。草地面积整体较为平稳, 在研究时段内下降了0.7%。

区域城镇化一直是影响地表NPP变化的重要因素, 随着近年来吕梁连片贫困区城镇化迅速推进, 区域城镇建筑面积增加显著(图 5), 城镇土地利用类型变化明显的区域主要位于各县市城区地带, 面积在研究时段内增加了609.14%, 其中47.42%来源于耕地, 52.2%来源于草地。而城镇区域的NPP值呈波动上升趋势, 多年平均值为179.78 gC m^-2 a^-1, 最大值出现在2018年, 为229.26 gC m^-2 a^-1, 最小值出现在2005年, 为149.43 gC m^-2 a^-1。

3.4 区域NPP变化与主要气象因子的相关性分析

先前研究表明气象因子是影响区域NPP值的重要因素^[34-37], 为进一步分析研究区NPP与主要气象因子(气温、降水、潜在蒸散发)的关系, 本研究将月均NPP值与气象因子月均值进行Pearson双侧相关检验以及偏相关检验, 结果如表 2所示。由表可知, 吕梁连片贫困区月均NPP值与月均气温、月均降水量、月均潜在蒸散发的平均相关系数分别为0.897, 0.962, 0.897, 均在0.01的置信水平上达到显著水平; 平均偏相关系数为-0.387, 0.887, 0.523。这说明吕梁连片贫困区NPP变化受到气温、降水以及潜在蒸散发共同控制且影响程度均达到极显著水平, 其中NPP与降水量相关程度最高, 表明研究区内降水量的变化是影响区域NPP变化的最主要因素; 在进行偏相关分析以消除另外两个气象因子的影响时, 温度与NPP呈负相关关系但相关程度不高, 潜在蒸散发与NPP呈正相关关系但相关程度不高, 只有降水量依然与NPP显著相关, 表明在吕梁连片贫困区降水量为区域NPP在月尺度变化的主要影响因子。

4 讨论

2000-2010年是“三北防护林”第二阶段第四期工程实施阶段, 吕梁连片贫困区与吕梁山生态脆弱区高度重合, 是该工程的重点实施区域, 期间该区域NPP呈现显著的上升趋势, 这表明生态恢复工程对植被净初级生产力具有明显的促进作用。2010-2015年间NPP出现下降趋势, 可能与区域城镇化快速推进带来的负面影响有关, 根据山西省统计年鉴(http://www.shanxi.gov.cn/sj/tjnj/)可知, 研究区城镇人口占总人口比重由2000年的14.31%上升至2018年的47.79%, 增长了234.05%。建筑的增加及路面硬化致使绿地(其中以草地与耕地为主)面积减少, 间接导致了区域水土流失及生态退化, 环境条件变得脆弱, 是2015年NPP出现低值的原因之一; 加之2015年发生了历史上第二强的厄尔尼诺气候事件, 导致了明显的气候异常。据研究表明, 华北地区雨季开始晚、结束早, 降水量为近13年来次少, NPP累积最重要的季节--夏季降水异常偏少。根据表 2可知, 研究区NPP变化与降水量变化呈正相关关系, 降水的减少直接导致了NPP的降低, 在图 3上有较为直观的反映。人类活动和气候变化共同导致了2015年的NPP低值, 致使这一时段NPP年均值呈现下降趋势。2015-2018年间, 区域NPP年均值再次快速上升, 由2015年的258.49 gC m^-2 a^-1升至2018年的331.7 gC m^-2 a^-1, 增长速率为24.4 gC m^-2 a^-1。2015年起吕梁山生态脆弱区林业生态建设工程被列入国家级工程, 区域生态环境恢复得到国家的高度重视以及政策的倾斜支持, 加之2018年气象条件较2015年好, 促使这一时期区域NPP呈上升趋势。

近年来, 部分学者研究了许多典型区域NPP的动态变化特征、探讨了其与各影响因子(降水量、温度、干旱指数、潜在蒸散发等)间的相关性等^[34-37], 如张镱锂^[34]等发现青藏高原高寒草地NPP在空间上呈现由东南向西北逐渐递减趋势, 且在研究时段内呈波动上升趋势。秦景秀^[35]等发现新疆地区的NPP也呈现明显的上升趋势, 人类活动是影响新疆植被恢复和退化的主要原因。陈杰等^[36]发现川西地区大面积区域NPP呈显著上升, 且不同土地利用及植被类型的NPP差异较大, 海拔与研究区NPP相关性非常显著, 人类活动对地区NPP变化的负干扰明显。这些研究结果对于区域碳循环研究、环境保护、经济可持续发展等具有重要科学价值。在黄土高原地区, 谢宝妮^[37]等证实黄土高原中东部NPP在过去20年亦呈现增加趋势, 与本研究的结论较为相似。人类活动对于吕梁连片贫困区NPP的干扰也十分显著, 2000年以来该区域开展生态治理恢复工程, 三北防护林布设, 天然林保护、退耕还林还草等措施的实施取得了显著成果, 区域绿地面积不断增加, 生态环境治理效果显著。但同时人类活动的负干扰在该地区NPP的变化中也得到了明显的体现, 如2010-2015年间区域NPP整体呈下降趋势, 这与区域该时期人口的增长、城镇化水平的快速上升有直接关系, 这种负干扰作用与前人在川西生态脆弱区的研究结果较为相似^[36], 警示了人们在进行生产生活的同时必须注重对生态环境进行保护。

5 结论

本文以2000-2018年为研究时段, 基于MODIS遥感数据、土地利用类型数据以及气象数据等, 利用CASA模型对吕梁连片贫困区NPP值进行估算并分析得出以下结论：

(1) 自2000年以来, 吕梁连片贫困区NPP年均值呈波动上升趋势, 93.46%的区域NPP呈增长状态, 仅有6.54%的区域呈减少状态。

(2) 研究时段内耕地面积缩减, 草地面积基本保持稳定, 林地、城镇面积增加且城镇面积扩张迅速; 各土地利用类型NPP年均值增长以耕地最为迅速, 达5.9 gC m^-2 a^-1, 林地最为平稳, 为1.32 gC m^-2 a^-1。

(3) 影响区域NPP变化的气象因子中, 与NPP的相关关系从大到小依次为降水量、潜在蒸散发、气温, 降水量为影响NPP变化的最主要气象因子。