工业革命以来, 化石能源的广泛应用导致大量温室气体激增, 生态系统在这一背景下逐渐失衡。植被作为区域生态系统中的重要组成部分, 在保持水土流失、调节气候及温室气体等方面起到了至关重要的作用^{[ 1]}。植被能通过光合作用捕获并封存CO₂来减少其排放, 从而达到植被固碳(carbon sequestration, CS)的目的^{[ 2]}, 研究植被固碳对于可持续发展战略的实施具有重要意义。目前世界各国主要通过增加地表植被数量来提高生态系统的自然固碳能力, 从而来降低大气中的CO₂浓度。近年来国内外学者尝试用植被总初级生产力^{[ 3]}(Gross Primary Productivity, GPP)、净初级生产力^{[ 4]}(Net Primary Productivity, NPP)等植被群落生产能力参数反映其固碳能力, 而NPP作为判定生态系统碳源碳汇的主要指标, 不仅能表征陆地生态系统的质量状况, 还能反映植被在自然环境下的生产能力。利用遥感估算NPP并结合生态系统服务功能实物量评估模型使植被固碳量的定量计算成为可能^{[ 5- 6]}, 而且已经有不少学者在各区域尺度^{[ 7- 10]}证实了NPP时空分异与气候因子的相关性^{[ 11]}, 此外还可能受到地形、植被类型等其他自然因子^{[ 12]}和人为因子的影响, 不同区域NPP的影响因素也不相同。因此, 研究区域尺度NPP变化对丰富碳源/碳汇研究具有重要意义^{[ 13]}。

黄土高原近年来受频繁的人类活动影响, 是中国陆地生态系统最脆弱和敏感的区域之一^{[ 14]}。黄土高原的人类活动是影响植被NPP变化的重要因素^{[ 15]}。土地利用/土地覆被变化(LUCC)不仅给黄土高原带来明显的地表结构变化, 而且影响以土地为载体的生态系统服务功能, 乃至区域生态环境和生态过程^{[ 16]}。此外, 黄土高原沟壑纵横, 其地形分异特征带来的植被垂直分带给不同植被类型NPP也造成很大差异, 而海拔高度是地形中最基础、最直接的信息, 除人为因素外海拔高度能直接或间接影响植被空间生长和分布差异^{[ 12]}。因此地表覆被变化以及综合海拔高度响应的研究对于黄土高原生态系统碳循环及周边地区环境保护和生态建设具有重要的指示意义^{[ 11]}。目前在黄土高原针对NPP的研究已有不少, 但前人或局限于某一亚区^{[ 17- 18]}, 或研究单一植被类型NPP^{[ 14, 19]}, 或其气候影响^{[ 20- 21]}, 而对黄土高原全区域范围的NPP与人为因素造成的土地利用转变之间的响应机制探讨较少, 且针对黄土高原随长时间序列的固碳变化情况, 以及土地利用转变与地形因子协同作用下的NPP变化研究也少有涉及。因此, 合理选用黄土高原遥感数据及适宜的模型定量估算NPP, 并分析黄土高原植被碳源碳汇时空格局变化具有重要的理论和实践意义。

本文从自然(海拔高度)和人为(土地利用转变)因子出发, 利用CASA(Carnegie Ames Stanford Approach)模型估算黄土高原NPP, 研究黄土高原在不同土地利用方式及地形分级中多时空尺度植被固碳变化及其影响, 并进一步加入土地利用变化以及海拔高度的协同作用定量评估其对NPP的影响, 可以更好地了解区域尺度生态系统对自然和社会环境变化的反馈能力, 以期为黄土高原的生态环境建设提供科学依据。

1 研究区概况

黄土高原(100°54′E-114°33′E, 33°43′N-41°16′N)包括中国青海省东北部、甘肃省中东部、宁夏回族自治区大部、山西省大部、陕西省中北部、河南省西部以及内蒙古自治区省区( 图 1), 东西长约1300 km, 南北长约750 km, 总面积约62.4万km²。黄土高原地处平原向高原的过渡带, 具有典型的暖温带半干旱大陆性季风气候特征, 且地带性差异特征明显, 自东南到西北依次为半湿润、半干旱、干旱气候, 年平均气温由15℃减小到3℃, 降水量由800 mm向200 mm递减, 太阳辐射总量由5000 MJ/m²向6300 MJ/m²递增, 土壤类型变化顺序为褐土、黑垆土、黄绵土、灰钙土、风沙土^{[ 22]}。

2 研究方法 2.1 数据来源与预处理

本研究使用中科院土地利用数据( 表 1), 并重分类成六种类型( 表 2)。GIMMS(global inventory modeling and mapping studies)NDVI每半月数据采用最大值合成(maximum value composite, MVC)月NDVI数据；并将研究区及周边MODIS NDVI影像拼接, 本研究所获得的两种数据集的NDVI数据都经过几何纠正、辐射校正、大气校正等预处理。两种土地利用和NDVI数据的重叠时间和重叠区域的数据都通过一致性精度检验, 满足本研究使用^{[ 20, 23- 26]}。气象数据来自中国气象信息中心, 选取黄土高原及周边80个站点地面观测数据, 并用反距离加权法^{[ 27]}得到气温、降水、辐射栅格数据。将所有栅格数据集的空间参考系统设置为正轴等积割圆锥(Albers)投影, 并将栅格重采样为1 km, 按研究区边界裁切。

2.2 计算与分析方法 2.2.1 NPP估算

本研究采用朱文泉改进的CASA模型^{[ 28]}, NPP估算主要由植被吸收的光合有效辐射和植物光能利用率完成。

(1)

式中, SOL(x, t)即月太阳总辐射量数据集上像元x位置在t月的数值, 单位为MJ m^-2月^-2, FPAR(x, t)根据某一植被类型NDVI以及SR的最大值和最小值线性推导而出^{[ 29]}, 常数γ是由植被在0.4-0.7 μm波长下所能利用的太阳有效辐射占太阳总辐射的比例确定的, 本实验常数γ设定为0.5；温度和水分胁迫系数T_ε1(x, t)、T_ε2(x, t)、W_ε(x, t)即由月平均温度、月累计降水量数据集上像元数值所决定；植被类型的最大光能利用率ε_max参数根据文献资料收集得到^{[ 30]}。

2.2.2 植被固碳估算

根据中国林业科学研究院森林生态环境与保护研究所起草、国家林业局发布的《LY/T 1721-2008森林生态系统服务功能评估规范》^{[ 31]}中森林生态系统服务功能实物量评估公式及参数设置一表知：

(2)

式中, G_植被固碳为植被年固碳量, 单位t/a；1.63代表森林每积累1 t干物质, 平均吸收CO₂ 1.63 t；R_碳为CO₂中碳的含量, 为27.27%；A为林分面积, 单位km², B_年为林分净初级生产力, 单位t km^-2 a^-1。

2.2.3 相关性分析

(3)

式中, Cor为相关系数, X、Y为X、Y的均值, 相关系数R介于[-1, 1]之间, 其绝对值越大, 表明变量X与变量Y的相关性越高；反之越低。

3 结果与分析 3.1 土地利用时空变化特征

从区域尺度, 黄土高原地表覆盖以农田、草地和林地为主, 分别占总面积37.34%、36.54%和22.75%, 其余三类用地占总面积3.37%。1990-2015年黄土高原的农田面积在时间尺度上剧烈减少, 林地持续增长, 草地先剧烈减少后有少许增长, 未利用地呈下降趋势, 水体和居民地均有不同程度增长；空间尺度上呈现东北-西南条带状分布, 土地利用类型从东南向西北依次为林地、农田、草地、未利用地, 居民地和水体零星分布在耕地周围( 图 2)。近30年未利用地减少8890 km², 其中69.8%转化为林草地, 该转变主要分布在研究区北部。空间尺度上植被由东南逐渐向西北覆盖, 其中林地面积增加115191 km², 增幅173%, 其次, 农田面积减少68036 km², 减幅25%, 未利用地、草地面积减幅为51%和18%。说明黄土高原用地类型变化明显, 且变化主要发生在耕、林草地之间( 表 2)。

3.2 NPP时空分布特征 3.2.1 不同用地类型NPP时空变化

黄土高原不同用地类型的NPP差异较大( 表 3)。林地年均NPP(619.5 gC m^-2 a^-1), 远高于农田(466.8 gC m^-2 a^-1)和草地(293.9 gC m^-2 a^-1)。1990年, 农田NPP总量最多(113.4 TgC), 草地次之(76.1 TgC), 林地较小(47.7 TgC)；2000年, 林地年总NPP最大(98.9 TgC), 农田次之(96.0 TgC), 草地较小(53.3 TgC), 这与退耕还林还草政策的实施密切相关；2015年, 植被覆盖持续上升, 年总NPP前三位仍是林地(119.3 TgC)、农田(112.9 TgC)和草地(71.8 TgC)。整体上, 年均NPP高值区落在林地、草地、农田植被分布区, 且NPP纹理与土地利用格局总体相似( 图 2)。

3.2.2 不同海拔高度下NPP时空变化

黄土高原年均NPP随高程的增加先降低后升高, 年总NPP则随高程变化趋势相反( 表 4)。近30年研究区在1至6级高程分级下的年均NPP分别为：582.49 gC m^-2 a^-1、533.89 gC m^-2 a^-1、385.92 gC m^-2 a^-1、399.63 gC m^-2 a^-1、451.80 gC m^-2 a^-1和573.07 gC m^-2 a^-1, 年总NPP分别为16.50 TgC、41.80 TgC、124.21 TgC、47.83 TgC、25.10 TgC以及11.51 TgC。

3.3 土地利用变化对NPP及固碳的影响

1990-2000年黄土高原土地利用变化面积125940 km², 占研究区面积的20.19%, NPP总和6.02 TgC, 植被累计固碳1.79 TgC。1990-2000年NPP均值增长排列为( 表 5)：耕地转为林地>未利用地转为林地>耕地转为草地, NPP减少程度排序：林地转为耕地>林地转为草地>耕地转为未利用地。1990-2000年数据表明, 林地转化为居民建设用地、未利用地, NPP数值为负, 未利用地转化为林地、农田, NPP数值为正, 即一般来说, 高等级向低等级用地类型转化会使NPP减少, 反之则增加。1990-2000年转变类型面积由大到小为：草地变为农田>草地变为林地>农田变为林地>农田变为草地, 而植被固碳总和最大的是农田变为林地(2.42 TgC)以及草地变为林地(0.31 TgC)。从NPP均值以及固碳量的角度看, 退耕还林还草政策的固碳效果良好, 但鉴于初期植被覆盖(尤其草地的不稳定性)以及固化的经济思想等人为因素都会造成固碳成效不明显。2000-2015年, 土地利用转变程度明显减弱, 尤其林地得到很好保持, 用地类型转变面积30815 km², NPP总和2.68 TgC, 研究区总计固碳1.19 TgC, 用地类型主要变化在农田转变为林地(7447.8 km²)、草地(5058.5 km²)上, 固碳变化主要发生在农田与林地(0.51 TgC)、草地(0.21 TgC)上。

3.4 土地利用变化与海拔高度协同作用下NPP的变化特征

前文统计了NPP在不同土地利用类型、海拔高度及土地利用变化中的变化特征, 但在何种海拔高度下适宜进行哪种土地利用转变尚不明确, 为探究土地利用转变与海拔高度的协同作用对NPP的影响, 将两因子叠加后结果显示( 图 3, 图中坐标轴1-6高程分级同 表 4)：耕、林草地未发生转变部分的年均NPP随海拔的增加逐级降低, 固碳量随海拔的增加先升高后降低, 且在同一高程分级下, 耕地转为林地的年均NPP最大, 如海拔低于500 m, 耕地转为林地的NPP(197.71 gC m^-2 a^-1)>草地转为耕地(113.35 gC m^-2 a^-1)>草地转为林地(67.64 gC m^-2 a^-1)>耕地转为草地(54.69 gC m^-2 a^-1)>林地转为耕地(31.84 gC m^-2 a^-1)>林地转为草地(-59.07 gC m^-2 a^-1), 故相同高程分级中, 退耕还林植被固碳效果最优, 而林地变为耕地或草地的年均NPP值为负, 均不能达到固碳目的。此外, 耕地转为林地的年均NPP随海拔的增加先降低后升高, 高海拔(>3000 m) 年均NPP高于其前两级海拔(1500-3000 m)NPP, 而固碳量随海拔增加先升高后降低, 证明高海拔虽然较适宜林地种植, 但种植条件差、范围小等原因导致固碳效果不理想, 因此黄土高原应充分利用3000 m以上的区域退耕还林, 其效果优于在海拔1500-3000 m退耕还林。在所研究的六种典型转变方式中, 除了退耕还林、变林为耕以及变林为草的年均NPP呈现明显正负趋势外, 另外3种转变(退耕还草、变草为林和变草为耕)年均NPP不稳定, 海拔低于1500 m, 变草为耕的NPP高于退耕还草和变草为林；海拔高于1500 m, 变草为耕的NPP剧烈降低, 海拔高于2000 m, 变草为耕的NPP值为负, 远低于退耕还草和变草为林, 海拔高于3000 m, 草地转为林地的NPP(78.90 gC m^-2 a^-1)远高于退耕还草(26.62 gC m^-2 a^-1)和变草为耕(-35.11 gC m^-2 a^-1)；而退耕还草的NPP除了在研究区低于500 m和高于3000 m海拔处低于变草为林, 在500-3000 m海拔, 退耕还草效果优于变草为林。

4 讨论 4.1 数据适用性分析

本研究模拟结果与朱文泉^{[ 28]}、苑全志^{[ 32]}、韩艳飞^{[ 33]}、石志华等^{[ 34]}学者在全国区域, 尤其是黄土高原得出一致结果, 说明研究选用数据和模型在黄土高原NPP研究中具有可靠性。为探讨其在黄土高原进一步分析的适用性, 本研究基于7个省(自治州)级行政区将它们在1990、2000和2015年的植被NPP与各土地利用类型面积占比、土地利用变化、海拔高度以及两者的协同作用做相关性分析^{[ 9]}。结果表明, 3个时期NPP与耕地和林地的相关性最高( 表 6), 呈显著正相关, 表明耕地、林地对NPP响应显著；其次1990年草地和居民建筑用地呈显著的正相关。而其余时期的草地、建设用地、水域和未利用地与NPP的无明显相关性, 说明本研究得出结果在阐释黄土高原在不同土地利用类型角度的NPP时空变化上具有可信度。同时, 植被NPP与土地利用转变因子的相关系数最高达0.74, 与海拔高度的相关系数达0.87, 植被NPP与土地利用变化和海拔高度协同作用的相关性为0.94, 可以看出黄土高原植被NPP与双因子协同作用相关性要高于各单因子, 且在低海拔比高海拔有更高的相关性( 表 6), 可能是因为植被在低海拔的覆盖率更高。

此外, 本研究结果异于谢宝妮^{[ 11]}、高志强^{[ 35]}等得出的结果, 表现为整体上略高于GLO-PEM遥感模型和BIOME-BGC、CEVST过程模型, 可能主要与模型选取有关。CASA光能利用率模型考虑光、温、水、营养物质等的影响, 而CEVST模型仅考虑气候因素变化的影响, GLO-PEM包括各种非气候因素如病虫害和火灾影响, BIOME-BGC所需参数包括多种复杂气象资料, 以及土壤和植被中的碳、氮、水等状态参数及分配比率等且忽略了空间异质性因素的影响, 这类模型仅在小范围空间, 以及空间斑块均质的条件下才能近似适用^{[ 36]}。此外, 模型参数、时间序列、研究区范围及数据质量的差异都会影响研究结果。

4.2 土地利用变化与海拔高度对植被NPP的影响

黄土高原作为我国生态环境最为脆弱的地区之一, 除了历史及自然气候环境的变异导致植被遭大量破坏以外, 最主要的原因还是在经济快速发展中生产生活建设的人为因素^{[ 37]}。人为因素中以开荒种田、开采森林和过度放牧等原因造成的土地利用变化对植被NPP以及固碳的扰动很大, 其中退耕还林还草的固碳效果显著, 反之释碳较多。本研究结果显示, 近30年黄土高原农田面积占比一直较大, 年均固碳47.7 TgC, 森林年均固碳39.4 TgC, 草地年均固碳29.8 TgC, 且植被覆盖度逐年增大(1990年植被覆盖度26%, 2000年28%, 2015年35%), NPP呈逐年增加趋势, 植被固碳也增多。同时, 研究区退耕还林达到的固碳效果较好, 年均固碳79.55 TgC。土地利用转变的确影响生态系统的碳源碳汇, 更有学者指出完全可以通过土地利用变化来实现《京都议定书》中设定的每年减少200亿吨碳排放的目标^{[ 38]}。

海拔高度对NPP的响应表现为：随时间序列逐年上升, 年均NPP随高程的增加先降低后升高, 年总NPP随高程增加变化趋势相反, 这与前人研究^{[ 10, 39]}结果一致, 黄土高原NPP受高程因子变化较大, 其原因在于：其气候等环境条件受地形的分异特征变化显著。山地影响陆间水循环, 进而控制地表水源分布, 使得该区域地表环境变化对地形的依赖程度较高, 因此, 海拔高度与NPP分级特征明显的主要原因是高程通过控制土壤以及水热条件影响植被垂直分布与气候变化, 进而影响植被区域分布^{[ 4, 10, 39- 41]}。

土地利用转变与海拔高度的协同作用对黄土高原植被NPP分异特征明显^{[ 17]}。研究区两因子叠加后, 仍以退耕还林固碳最多( 图 3), 由于林地的冠层相比耕地和草地对降水有阻截作用, 枯枝落叶层具有节流作用, 能很好地阻抗降水对地表的冲刷, 对地表土壤有一定的稳固作用, 导致林地在研究区具有较高的土壤保持效应^{[ 42]}。但耕地和草地之间的转变在两因子的协同作用下变得不稳定：海拔低于1500 m, 宜变草地为耕地, 可能是因为低海拔处人为耕作周期更长, 大部分农作物长势均比草地好, 且农作物更受到人为培育, 而草地借势生长但远不如耕地的固碳效果, 因此低海拔农田固碳效果优于草地；而海拔大于1500 m, 中海拔处频繁的农业活动会干扰碳循环, 短周期的耕作加速了土壤中的生物活性, 导致土壤中CO₂含量升高, 在影响碳酸盐的溶解和沉淀的同时影响地上植被的光合作用, 进而造成固碳能力的降低；而在中海拔草地受人为活动影响不大, 草地群落长期保持自然生长状态, 从而导致耕地的固碳能力不如草地, 这也与Hui, A A.等^{[ 43]}研究结果一致；但随着高程逐级增大, 耕地和草地的固碳效果均远不如林地, 对气候以及热量等环境条件的需求程度逐渐提高, 在高海拔处不易积蓄水分, 造成土壤水分含量低, 可耕作类型数量下降, 不宜于农作物和草地生长, 而在高海拔林地生态系统其固碳能力的鲁棒性较高^{[ 44]}, 因此, 海拔高的地区耕地和草地的固碳效果不如林地。

综上, 植被NPP对自然与社会因子的响应关系错综复杂, 但显著相关性仍存在于其时空变化上, 同时因研究区、尺度和方法等各种因素的不同也会给NPP的反演估算带来差异性和不确定性。本研究也仅从自然和社会因子中选取海拔高度、土地利用变化以及两者协同作用的角度分析了黄土高原植被NPP的时空分异特征, 在NPP响应因子的研究中对定量方法与权衡指标的选取值得进一步探索。

5 结论与展望

本文研究了黄土高原1990-2015年在土地利用转变与海拔高度协同作用中植被固碳的变化情况, 利用CASA模型进行NPP遥感估算, 在此基础上统计并分析土地利用转变及海拔高度的协同作用对植被固碳的综合影响。总体来看：(1)黄土高原年均NPP逐年上升, 年均NPP增速2.74 gC m^-2 a^-1, 年均固碳增速1.13 TgC/a, 研究区林地固碳效果理想(年均NPP为619.5 gC m^-2 a^-1)。(2)土地的开发利用, 即自然生态用地转为商服设施用地, 释碳增多；恢复自然生态用地, 固碳增多。(3)黄土高原年均NPP随高程的增加先降低后升高, 年总NPP和固碳量随高程增加变化趋势相反；(3)研究区退耕还林的植被固碳效果最好；而林地变为耕地或草地均不能达到固碳目的, 此外, 研究区低海拔和高海拔相比中海拔地区更适宜退耕还林, 海拔低于1500 m宜变草为耕, 海拔高于1500 m宜退耕还草, 海拔高于3000 m宜变耕、草为林。

不同遥感传感器提供的数据源可能会影响NPP估算模型的结果, 而现行方法针对不同传感器获取的遥感影像NDVI时间序列插补并不完善, 采用逐像元回归插补或是区域综合插补还需有待研究。关于长时间序列的区域尺度固碳研究, 从NPP遥感估算的角度来说, 更高精度、更长时间序列、更快解算速度是未来需要进一步研究的方向。