张玉红，闫 浩

（1.哈尔滨师范大学地理科学学院，黑龙江 哈尔滨 150025；2.自然资源部第二大地测量队，黑龙江 哈尔滨 150025）

引言

森林是地球生态系统的主体，在陆地生态系统中拥有最大的面积和最丰富的物质资源，约占陆地总面积的30%［1］。净初级生产力（Net primary productivity，简称NPP）是绿色植物在单位时间和面积内累积的有机物［2］，是表征地球表层植被覆盖生物量的重要指标，它反映了森林固定大气CO2的能力，对全球气候研究中起到重要的参考价值［3］。随着遥感技术（RS）与地理信息系统技术（GIS）的应用与发展，卫星遥感数据被大量应用到林业科学中来。通过遥感与地理信息系统技术，各类观测数据更容易获取，空间分析与模拟更为便捷准确，同时估测范围更加广泛。例如利用不同遥感平台的卫星影像进行森林火灾的监测与管理［4-6］，或者利用地理信息系统强大的空间分析功能进行各类森林火灾的模拟研究［7］。作为森林的监测与管理的重要指标，NPP的计算和模拟可以通过遥感技术解决传统定点估算NPP数值的局限性和长期性等问题，使NPP的估算获得更多的思路和方法［8］。同时通过遥感技术与不同的模型方法相结合，可以对森林植被的特征参数进行精确反演，获取更多的信息数据［9］。主要的NPP 估算模型分为统计模型、过程模型和参数模型三类［10］。统计模型主要考虑气候因子的影响，估算植被的NPP 主要考虑降水、气温、太阳辐射等因子与植被NPP的统计回归与相关性。例如Miami 模型［11］、Thornthwaite Memorial 模型［12］和Chikugo 模型［13］，它们输入参数简单，基础数据较容易取得，能够方便研究，可以对陆地生态系统净初级生产力进行估算和预测［14］。但相对的能够作为依据的生态理论缺少严密性，资料的准确性也有待提高。过程模型又被称机理模型，代表模型如TEM 模型［15］、CERTURY 模型［16］、BEPS 模型［17］。这些模型从一个整体系统的角度，对物质和能量的交换进行模拟并建立了相应的模型或模型库，从机理上模拟植被的光合作用和蒸腾蒸发等过程，从而推进了大空间尺度的植被NPP 研究。但是为了保证严谨的特点，过程模型的构造和机理相当繁琐，需要的参数如植物生态特征、土壤植被、气象变化数据等，数据需求大并且难以一一获取，在空间上较难定量化。应用较为广泛的则是参数模型，即光能利用率模型，其理论基础是植被光能利用率，将植被累积的生物量认定为是植被将太阳入射辐射吸收和转化等一系列光合作用过程得到的。这类模型的代表如CASA模型［18］、GLOPEM 模型［19］，C-FIX 模型［20］。与统计模型和参数模型相比，这类模型相对简单便利，它是通过使用遥感数据来计算光合有效辐射。因为容易获得周期性强和范围广泛的遥感数据，使得这类模型的时空分辨率较好，能够进行大至全球尺度范围的估测与观测。这些优势促进了这类模型的发展研究，成为现阶段NPP研究中的主要研究模型。

作为全球变化的敏感区域之一，大兴安岭地区位于北半球的中高纬度，受到了世界各国学者的关注。大兴安岭地区是重要的林业基地，也是我国唯一拥有大面积兴安落叶松（Larix gmelinii）的原始林区。该区高达70%以上的森林覆盖度和巨大的森林蓄积量，使大兴安岭森林地区有着重要的生态作用和社会效益［21-22］。1987 年5 月6 日，大兴安岭地区发生了严重的森林火灾。这场特大火灾使得森林资源遭受到了严重损失，植被覆盖以及森林功能都需要较长时间的恢复。文中基于遥感数据并结合森林调查数据对大兴安岭地区图强林业局的多年森林NPP的数量和空间分布进行估测，并分析影响森林NPP的主要自然因素以及对森林火灾之后的NPP 分布及恢复情况进行遥感监测，可以为林区管理提供基础信息，为林业规划和持续发展提供科学的依据。

1 研究区及数据源

1.1 研究区概况

黑龙江省大兴安岭地区图强林业局位于122°18'E-123°28'E，52°15'N-53°33'N 之间（图1），属于漠河县管辖范围，是我国最北的木材生产基地。整个林业局的总面积约为5 000 km2，其中林业用地约为4 970 km2，而森林覆盖度达到87%。该研究区属于寒温带大陆性气候，空气湿度小，气温年较差与月较差均较大。适应这种气候的植被种类较少，大约70%以上的树木都为兴安落叶松（Larix gmelinii），其余部分又以白桦（Betula platyphylla）和樟子松（Pinus sylvestris var.mongolica）为主，分别占整个林业局森林蓄积量的21%和7.5%，其他树种总共不足0.5%。该区地势较为平坦，但也有低矮丘陵，多处于东部，整体地势成南高北低的分布。1987 年5 月份大兴安岭地区发生了特大火灾，总过火面积超过10 000 km2，有林面积达到了70%左右。其中，图强林业局的森林资源受到了严重损失，过火区域超过一半以上。

图1 大兴安岭地区图强林业局的位置示意图及过火强度图Fig.1 Location and fire intensity map of Tuqiang Forestry Bureau，Daxing’anling Region

1.2 数据来源

1.2.1 遥感数据

文中所使用的遥感卫星数据均来自美国Landsat 系列卫星，轨道号为122/23。根据NDVI 估值与植被生长季的关系［23］，并考虑到7 月份植被的生长进入稳定阶段，因此选择时相为1987 年到2013 年中7 月份的遥感影像，云量皆为5%以下。经过影像数据的预处理之后，获取相关植被指数数据，并利用决策树进行植被信息提取。

1.2.2 气象数据

气象数据来源于英国东英吉利大学气候研究所（CRU）提供的1987年到2013年的逐月网格气温和降雨数据，以及附近气象站点记录的月平均温度、月总降水量和太阳辐射数据。

2 研究方法

2.1 CASE模型

1972 年，Monteith 在光能利用率的原理上提出了利用光合有效辐射的观测值和光能利用率来估算植被生产力的方法，即建立了Monteith 方程［24］。在此基础上，Potter 等学者相继提出了计算全球植被净初级生产力模型并建立了CASA 模型（Carnegie-Ames-Stanford Approach），完成了基于光能利用率方法的陆地净初级生产力全球估算，是目前对NPP估算研究中使用较多的模型［25］。该模型的主要优点在于：模型相对简单，需要的参数不多，因而便于计算和处理；由于部分参数来源于遥感数据，因而不受地面实测站点的条件限制，可以进行大区域上的NPP 估测；同时在遥感数据的基础上可以对植被进行分类来反映植被的变化状况，能够借用多时相遥感数据来对NPP进行动态监测等［26-27］。

CASA模型公式如下：

式中，APAR(x,t)代表的是t月中像元x吸收的光合有效辐射（gC·m-2），ε(x,t)代表的是t月中像元x的实际光能利用率（gC·M-1J）。

2.1.1 APAR的估算

通过遥感数据来估算植物叶片所吸收的光合有效辐射（APAR），其公式为：

式中，SOL(x,t)所代表的是像x处在t月的太阳总辐射量，FPAR(x,t)所代表的是植被层能够吸收的入射光合有效辐射所占比例，常数0.5所代表的是植被能够使用的太阳有效辐射所占的比例。在一定范围内，可以利用FPAR 与植被NDVI来确定FPAR 的最大值和最小值。同时也可以使用FPAR 与比值植被指数之间的线性关系来进行估算。前者估算的数值较高于实测值，后者则低于实测值但误差相对较小。因此将两个结果结合起来可以提高FPAR的估算精度。

2.1.2 光能利用率的估算

光能利用率一般是指植被在地表单位面积内，通过光合作用所产生的有机物中的化学潜能，与同一时间内该块面积上收到的有效太阳辐射之比。气温、土壤湿度以及大气湿度、气压等环境因子会通过影响植物的光合能力，进而调节植被的NPP数值。所以，在模型中光能利用率是估算NPP的重要因子。其中，最大光能利用率主要受温度和水分含量的影响，公式如下：

式中：Tε1(x,t)反映在低温和高温时植物内在的生化作用对光合作用的限制，从而降低净第一性生产力。Tε2(x,t)表示环境温度从最适温度Topt(x,t)向高温或低温变化时植物光能利用率逐渐变小的趋势。水分胁迫影响系数Wε(x,t)反映了植物所能利用的有效水分条件对光能利用率的影响。εmax是理想状态下的最大光能利用率（gC·MJ-1）［28-29］。因此，归纳其模型参数主要如表1所示。

表1 模型参数简表Table 1 Model parameter profile

2.2 模型参数获取

CASA 模型参数主要包括固定参数、植被类型、NDVI、月平均温度、月降水总量以及月太阳总辐射。NDVI通过预处理后的遥感数据在ENVI软件平台上计算得到；植被类型图是基于遥感影像并使用决策树方法获得；气温和降水数据均来自于CRU数据，通过空间插值处理对气象数据进行网格内插，得到研究区整体的气象栅格数据；太阳辐射数据则来自于气象站点和历史数据统计得到。通过以上数据的处理获得1987年至2013年共13期模型参数数据。

3 结果与分析

3.1 图强林业局森林NPP估算

3.1.1 图强林业局森林NPP的数量变化

将以上模型参数数据输入到模型中，获得13期7月份森林NPP估算结果。由于文中主要研究森林所带来的NPP产量，所以将非森林地区的数值刻意减少让它近于0，但这同时也降低了研究区全部整体的NPP平均值。因此在统计过程中，除了统计全部区域的NPP 平均值外，也统计了森林区域的NPP 平均值。图强林业局地区的1987-2013年7月的NPP产量统计值如图2所示。

图2 7月份森林NPP变化趋势Fig.2 Trend of forest NPP Change in July

从NPP 数值上看，1987 年7 月的总产值、平均值和最大值由于森林火灾都远远低于其他年份。各类统计值在2010 年达到最大，2005 年达到最小值。其他年份数据均在207 gC/m2～209 gC/m2范围之内，数值变动很小，平均值为208.03 gC/m2。从整体趋势上看，研究区7 月份森林NPP 的产值呈上升趋势。在火灾发生之后的第二年即1988 年7 月的数值就直接提升了近乎一倍。在之后的20 年里，数值增加相对较少。1993、1994、1995 三年7 月的NPP 产量平均值为4.817×1011gC，2000、2001、2002 三年7 月的NPP 产量平均值4.975×1011gC，NPP 产值增加了0.16×1011gC；2005年和2006年7月份的NPP 产值下降幅度很大，到了2008年之后开始回升；2009、2010、2011 三年7 月的NPP 产量平均值为5.711×1011gC，与之前相比NPP 产值增加了0.74×1011gC。2013年各统计值略有下降。

从以上结果可知，火灾后森林NPP 恢复很快，尤其是最大值。但是平均值的恢复就需要一定的时间。同时也可以看出，除了突发的火灾影响之外，森林NPP的变化还受到其它因素的影响。

3.1.2 图强林业局森林NPP的空间分布变化

文中选取图强林业局1987年-2013年之间几个典型的7月森林NPP空间分布图，包括1987年、1988年、1995年、2002年、2005年、2006年、2009年和2013年。

（1）火灾后恢复初期（1987-1988年）

从1987 年7 月的NPP 产量值分布图清晰的看出图强林业局受到森林大火的影响，在森林大火过后，整个研究区北部受到了严重的损失，很多地区的森林NPP产量值低于10 gC/m2，轻度过火区域的数值也在平均NPP产量值之下。未受火灾的地区仍有很多高NPP产值地区，但总值仍较往年减少很多（图3（a））。1988年研究区进入了初步恢复阶段，过火区域的NPP产量得到了显著提升，但与未过火区域相比，平均NPP产值仍相对较低。由于部分树木生长快速或恢复密集，在过火区域也出现了高NPP产值（图3（b））。

图3 图强林业局7月份NPP分布Fig.3 NPP distribution map of Tuqiang Forestry Bureau in July

（2）火灾后恢复中期（1988-2009年）

1995 年后数值呈缓慢上升趋势，森林的NPP 产值进入了平稳恢复阶段。这一期间，森林的恢复不像1988年一样迅速，仍可以从图中看出火烧痕迹，北部有较多地区数值偏低。但对比1988年的数据可以看出，与未过火区域相连的区域或轻度过火区域恢复进展良好，已经出现大面积的NPP高产值地区。以1988年最先恢复的地区为中心向外的区域的数值也达到了标准值以上（图3（c））。到2002年，7月份森林NPP 的产值相比7 年前有了整体的上升，数值较低的区域也减少很多，但仍能看出北部数值较南部数值偏低，重度过火区域仍然达不到森林植被NPP 产值的平均值。相比之下最北部的区域已经基本恢复，高数值区域也较1995 年有所增多（图3（d））。2005 年森林NPP 各统计值突然降低，无论南北均出现大量的低数值区域。2006年较2005年相比，虽然有部分区域正常，但也可以清晰地看出多处数值偏低区域。从历史资料来看，这两年并未发生重大事件。从分布图上看，区域整体数值变低，并非是局部区域减少。结合现有资料，初步判断是这两年7月的气候较往年不同而影响了森林植被的NPP产值。从图中还可以看出，虽然整体数值偏低，但是南北的差异已经并不明显（图3（e），图3（f））。

（3）火灾后恢复后期（2009-2013年）

到了2009 年，森林的NPP 产值进入到稳定阶段，之后的3 年时间里7 月的森林NPP 产值上下浮动不大，且整体数值都明显高于10 年前。从分布上很难看到曾经的过火区域，南北分布趋势也呈现一致，高产值区域也大面积出现（图3（g））。因此判断从2009 年开始，森林已经近乎恢复到原有的生态环境。2013 年的数据由于基础遥感影像有云覆盖的原因，导致南部底端部分森林没有数值，大大降低了NPP的数值，但当月的NPP最大值浮动并不大，总产值与前几年相比略有下降，但仍然保持在平均值以上（图3（h））。

从整个研究时段来看，研究区北部地区受到大火影响，多年总体数值一直较低，且变化明显呈增长趋势，到2009 年后则趋近平稳；南部未过火区域没有明显变化，多年来数值多由当月气候决定。整体上来看，非森林地区的周边数值变动较为明显。结合2011 年和2013 年的分布图进行对比，可以看出非森林周边多处的NPP 数值明显降低，距离远处数值则相对偏高。结合当地考察与专家咨询得知，图强林业局作为重要的木材储蓄地，要对森林木材进行采集，道路两边的树木方便采集，这会间接导致NPP 产值降低；采集过后会撒下树种，重新出现新生树木，而新生树木地区的森林覆盖度相对较低，使得该区域的NPP数值偏低。

3.2 影响因素分析

森林NPP 估算所利用的CASA 模型中考虑到的因子主要有气温、降水、太阳辐射、树种。文中选取的数据估算的均为7月份的森林NPP值，考虑到当地每年7月的太阳总辐射和温度变化都不大，因此文中重点分析降水和不同树种对森林NPP的影响。

3.2.1 降水对森林NPP的影响

图4为1987年至2013年期间研究区的温度和降水变化，结合前面的森林NPP估算结果可以发现几个特殊年份。如1993 年和2002 年的降水量远高于邻近年份，但NPP 值却十分相似。相反2005 年，2006 年，2008年降水相对较少，NPP 数值也较往年减少很多。为了更清楚的研究降水对NPP 生产量的影响，本文制作控制实验数据，即保留1987年、1993年、2002年、2005年、2006年和2008年的地面遥感数据、温度、太阳总辐射的数据，使用2009年的气象降雨数据，再根据模型估算出在此情况下的研究区NPP的生产量情况，与实际的估测值进行对比。

图4 1987-2013年7月气温降雨变化图Fig.4 The variation of temperature and total precipitation in July，1987-2013

结果表明，1993年和2002年结果与控制实验结果相比，虽然控制实验结果各项数值都有所减少，但变化不大，总产值、最大值与平均值均减少的数值都不到0.8 gC/m2。因而证明并不是降水的数值越多NPP 的数值就会越大。若水分过多超出植被和土壤能够吸收利用的数值，还有可能抑制植被NPP 的产值。1987 年、2005 年、2006 年、2008 年的控制实验结果数据则表明：随着月总降雨量的增加，NPP 的产值得到了大量的提升。

结合2005 年实际估测值与控制数据的区域分布图（图5）来看，其实在2005 年，就已经看不到曾经的过火区域，南北分布趋势也呈现一致，与2009年以后的分布图相似。判断影响2006年和2008年NPP产值的降低的主要原因为水分过低，没能达到较适宜的降水量。

图5 2005年实际模拟分布图与控制降水数据分布图对比Fig.5 Comparison of actual simulated distribution and control data in 2005

3.2.2 不同树种对NPP 产值的影响

植物的特性影响其光能利用率，如树叶的形状与大小、树干的高度、对阳光的喜好、水分的利用吸收等。一般情况下，阔叶林的光能利用率要高于针叶林。因此在模型中，不同的树种有着不同数值的最大光能利用率，光能利用率较高的树种的NPP产值会相对较高。文中研究区林种主要有白桦、落叶松和樟子松。表2为不同树种的NPP产值统计结果。

表2 研究区不同树种NPP产值数据Table 2 NPP value of different tree species in the study area

白桦、落叶松和樟子松的NPP数值增长趋势与森林整体的增长趋势保持一致。同样都是在1988年大幅度增长后逐渐缓慢，在2005年和2006年数值减少，到2008年数值增加并再度进入平稳。

白桦是喜光的树种，因为它生命力强、生长周期快等特征，使它在森林火烧后能够最先生长起来，形成大片的新生白桦林。白桦的产值为图强林业局总产值的一半，最多比例为51.41%，最少为49.17%。7 月产值最多的在2010 年为2.863×1011gC，平均值从1988 年的144.16 gC/m2提升到了2010 年的166.97 gC/m2，比落叶松的同月的NPP 平均值高出36.92 gC/m2，其主要原因是白桦树属于阔叶林，落叶松和樟子松属于针叶林，相比这2 种树种，白桦树本身的叶面积指数，植被冠层覆盖度都相对较高，白桦树叶内吸收的入射光合有效辐射以及光能利用率也要高于落叶松和樟子松，因此单位面积的NPP产量要高于落叶松和樟子松。落叶松也是喜光的树种，适应性强、耐寒，它的树干与树高一般高于白桦树，且全年生长。落叶松的NPP 多年均值大约占总产值的34.38%左右，落叶松最多比例为35.83%。樟子松也是适应力强的树种，但由于它树冠稀疏，针叶稀少，再加上樟子松在图强林业局的占地面积较少仅有7.5%，因此樟子松的总产值与其他两类树种相比数值相差很多，NPP 的平均值也比其它两类树种低很多，多年均值只有95.58 gC/m2，只在2009 年、2010年、2011 年3 年超过了100 gC/m2。从NPP 的总产值比例来看，相比于其它年份的樟子松比例大多数都低于10%，1987 年的15.33%就显得格外突出，这是由于1987 年落叶松大面积被烧毁，从而导致樟子松的NPP 产值比例相对上升。

可见，不同的树种对NPP 的数值有很大影响。在图强林业局，7 月的主要NPP 产量来源为白桦树，约占全部产值的50.13%，落叶松和樟子松约则分别约占34.38%和10.27%。

3.3 过火区域的森林恢复分析

研究区在1987 年的特大火灾中过火面积超过一半，其中绝大部分是重度过火区域，且均位于图强林业局总局较远的北部地区，能够投入的人力物力与南部地区相比较少，无法在种植树苗之后进行长期的人工抚育措施，栽种的新生树苗的生长环境得不到人为的改善，导致过火后北部的森林NPP 产值没能得到显著提升。

从森林NPP 总产值上来看，1987 年7 月的过火地区NPP 产值仅占当月总产值的3.97%，第二年7 月的产值大幅度提升，占当月总产值的21.92%，其主要原因是在火灾过后，林业局采取了人工造林的恢复措施，在过火地区大面积种植落叶松的幼苗，并采取人工更新等技术加快森林的恢复。1993 年到1995 年期间，过火区域的总产值并没有较大的提升，所占比例也略微有所下降。这是由于图强林业局只有6、7、8 三个月适宜植被的生长，其余月份温度都较低并且几乎无降雨，在这样的气候环境下，幼苗的成活率会大大降低，生长速度也会变得较为缓慢。2000年开始，过火地区的NPP产值上升趋势就变得较为明显。白桦和落叶松本来就是适应当地气候的树种，存活下来的幼树会逐渐进入高生长旺盛期，此时即使没有人工的抚育树木也能够正常的生长，因此过火地区的森林开始的自然恢复程度较好，产值也逐年提升。2005 年由于降雨量过低的原因，过火地区的总产值略微有所下降，但是占总产值的比例仍然提升到了25.02%，此时过火区域与非过火区域的差异已经逐渐减少。2009年以后占总产值的比例一直在26%以上。

结合研究区7 月份的过火森林区域NPP 平均值和全部森林区域NPP 平均值来看，可以明显看出其动态变化（图6）。1987 年二者相差较大，这一差距在1988 年大幅度缩减，到1995 年开始差距逐渐缩小，到2009年以后就趋近一致，说明过火区域的森林植被基本恢复。

图6 过火区和全部森林区域的NPP平均值对比Fig.6 Comparison of NPP mean value between the burned area and total forest area

4 结论与讨论

文中以黑龙江省大兴安岭地区图强林业局为研究区，充分利用GIS 和RS 技术，结合生态模型对火灾后的森林植被恢复情况进行模拟和监测。从森林植被净初级生产力NPP 入手，对研究区1988 年以来至2013年的森林恢复情况进行遥感估算，并通过控制实验对比分析影响森林NPP估算模型的主要因素和过火区域的植被恢复情况。研究表明：1987年森林大火造成的研究区森林NPP数值下降严重。但是在适当的气候条件下，森林NPP 逐年恢复。1988 年进行了森林人工种植，之后逐年恢复；1993 年之后过火区域的NPP 产值逐渐增加，到2009年与未过火区域的NPP产值接近，植被覆盖恢复到平均水平，恢复时间大约在20年左右。从空间分布来看，研究区北部地区受火灾影响较大，多年森林NPP总体数值变化明显呈增长趋势，到2009年后则趋近平稳；南部未过火区域没有明显变化。在降水量适中的情况下，森林NPP产值呈现增加的趋势，恢复较快；同时不同的树种的恢复程度有所差别。

可以看出，结合参数生态模型的遥感监测解决了传统定点估算NPP 数值的局限性与长期性等问题，同时能够更有效地计算光合有效辐射，使得大尺度的估测与监测得以实现。但是作为一种参数模型，CASA 生态模型存在一定的不足之处，例如缺乏生物学基础，对植物、大气与土壤之间的联系的解释尚不完整等。但是与传统的统计模型相比，利用容易获得的周期性较强且大范围的遥感数据来计算光合有效辐射，其时空分辨率更好，更容易进行大尺度的NPP 估测与观测，是目前NPP 研究的主要模型。文中由于遥感数据和其它资料的获取限制，缺乏研究时段内的逐年和逐月的森林NPP产值估算；同时，研究区在1987年特大火灾之后原有的火烧迹地又发生过多起火灾，对森林的恢复程度与恢复时间均产生了影响，因此不能全面分析而精准分析出森林的恢复情况；今后应尝试其他渠道获取研究区的遥感数据和环境数据，多方面多因素综合考虑以提高结果的准确性；同时要充分利用实测数据对模型估算的结果进行精度验证。