广东省针叶树种蓄积量和生物量生长模型研究

黄金金，刘晓彤，张逸如，李海奎＊

(1. 中国林业科学研究院资源信息研究所，北京 100091；2. 国家林业和草原局森林经营与生长模拟重点实验室，北京 100091)

森林蓄积量和生物量是反映森林资源水平的基本指标，精准估测森林蓄积量和生物量，可以为林业资源管理决策提供科学依据[1]。国内外关于蓄积量和生物量模型的研建主要集中于单木水平，而林分水平的较少，尚未形成我国林业行业标准[2]。林分水平的生长模型能直接预测林分蓄积量和生物量，减少单木的误差积累[3-4]。目前估计区域尺度蓄积量和生物量的建模方法有多种，应用较多的一是以异速生长方程为基础，根据林分总胸高断面积、林分优势高、林分平均高等估计林分以及区域的地上生物量[5-6]；二是生物量转换因子法，张会儒等[7]将全国森林类型按照林分优势树种归并为9 类，构建了全国主要树种组林木与材积相容的生物量回归模型；三是通过单木模型推算，傅煜等[8]通过Monte Carlo 法反复模拟由单木生物量模型推算区域尺度地上生物量，估计了江西省杉木地上总生物量，秦立厚等[9]分析量化单木生物量模型中各种不确定性来源，估计区域地上生物量。

大量研究表明，林木胸径与林木年龄之间存在紧密的生长关系，为建立有关的胸径生长模型奠定了基础[10-14]。曹磊[15-16]和龙时胜等[17]采用Richards生长方程拟合树木直径与年龄的生长关系，模型的拟合优度较高，其中龙时胜等[17]利用模型估计林分平均年龄时，拟合年龄与真实年龄的绝对误差最大值为2 年，相对误差最大值为8.2%，两者间的差异性较小，说明利用林木多期直径估计异龄林林分平均年龄具有较高的可行性。大量研究[18-20]认为，林分蓄积量、生物量生长主要受树木自身大小、竞争因子和立地条件的影响。

本研究以广东省3 个针叶树种155 块样地为对象，利用5 期连清数据，基于林分平均胸径与年龄的生长关系估计样地林龄，进而建立林分水平的蓄积量和生物量生长模型，以体现生长潜力的参数a分级和与生长速度有关的参数b分级、考虑林分竞争与否、分步建模（一元非线性回归）和联合建模（非线性误差变量联立方程组）共8 种组合，构建模型系，评价模型的拟合优度和预测效果，筛选最优模型，以期实现对广东省针叶林蓄积量和生物量的准确估计。

1 研究区概况与数据来源

1.1 研究区概况

广东省位于20°13′～25°31′ N 和109°39′～117°19′ E，处于东亚季风区，包括了热带、中亚和南亚热带3 种气候，年均气温22.3℃，月降水量约313 mm，年日照时数约1 746 h，太阳辐射总量约4 200～5 400 MJ·m−2，生长期水热充足，利于植被生长。地貌类型复杂多样，以低山、丘陵为主，地势北高南低。

1.2 数据来源

建模数据：从广东省第5—9 次连清数据中，选择马尾松（Pinus massonianaLamb.）、湿地松（Pinus elliottiiEngelm.）、杉木（Cunninghamia lanceolata(Lamb.) Hook.）3 个针叶树种，筛选5 次清查中皆存在且保留木大于5 株的样地，共155 块（马尾松林91 块、湿地松林16 块、杉木林48 块），样地分布见图1。固定样地为0.066 7 hm2，样木调查因子包括样木号、胸径，单株材积按广东省地方一元材积表[21]计算获得，单株生物量按行业标准生物量模型[22-24]获得。为避免株数差异造成样地蓄积量和生物量的差异，本研究中样地蓄积量和生物量为样地内5 期保留木的蓄积量和生物量算数平均值，统计量如表1 所示。

表1 建模数据统计量Table 1 Modeling data statistics

图1 建模数据固定样地分布示意图Fig. 1 Distribution map of fixed sample plots of modeling data

样地胸径采用林分平均胸径：

式中，Dg为林分平均胸径，n为样木株数，di为样木胸径。

验证数据：选取2002 年达到起测胸径、存在于4 次连清中，且保留木大于5 株的样地，共126块，其中马尾松59 块、湿地松22 块、杉木45 块。各统计量如表2 所示。

表2 验证数据统计量Table 2 Verification data statistics

2 方法

2.1 模型构建

采用Mitscherlich 理论生长方程，拟合5 期林分平均胸径和相应时间间隔之间的关系；以Richards理论生长方程，拟合样地平均蓄积量和生物量与相应时间间隔之间的关系。

式中：Dgij、Vij、BIOij分别为第i个样地第j期的样地平均胸径、蓄积量和生物量，agei是第i个样地的初期林龄，adj为第j期与初期（1997 年）的时间间隔，a1、a2、a3、b1、b2、b3、c2、c3为参数。

基于上述胸径、蓄积量和生物量3 个生长方程，分别以体现生长潜力的参数a分级和与生长速度有关的参数b分级、考虑竞争与否、分步或联合建模等8 种组合，构建模型系（表3）。

表3 8 种模型系Table 3 Eight model systems

2.1.1 含林木竞争的胸径生长方程 对于林分水平模型来说, 与树冠有关的竞争指数能较好地反映树木间的竞争关系，因此选用基于交角且对象木与竞争木距离小于10 m 的竞争指数ju10，详见惠刚盈等[25]。由于林木竞争对单木胸径的生长具有较大的影响，但对林分平均胸径的影响并不算太大，为解释林木竞争对林分生长的影响，故将竞争指数引入了胸径生长方程：

式中：b0为竞争指数参数，ju10ij为第i个样地第j期的竞争指数（样地竞争指数为样地内所有保留木竞争指数的算数平均值），其余同上。

2.1.2 参数分级 为了兼顾区域尺度上立地质量的多样性以提高模型在大尺度范围内的适用性，对参数a或b分级代替立地质量的差异，以参数a分级为例：

将式（2）作为立地质量划分等级的分级方程，已知林分平均胸径Dgij和时间间隔adj，以样地号和树种号构成的唯一编号作为参数a和age的哑变量，树种作为参数b的哑变量，拟合方程（6）。将参数a1i从低到高排序，并按照有序样本聚类分为6 类。

式中：a1i(i=1,2,...,155)是第i个样地的参数，b1k(k=1,2,3时分别代表马尾松、湿地松和杉木)是第k个树种的参数，其余同上。

2.1.3 分步建模 以Model 5 为例，基于式（2），以参数a区别树种，以参数b区别立地质量差异，以唯一编号作为初期林龄age（未知）的哑变量，根据5 期林分平均胸径Dgij和时间间隔adj之间的关系来估计样地初期林龄agei，见式（7）。将agei作为式（3）和（4）的自变量，以参数c区别树种，其余参数处理同上，分别拟合样地平均蓄积量Vij和生物量BIOij与林龄agei的关系，即式（8）和（9）。

式（7）～（9）中：a1k、a2k、a3k(k=1,2,3)分别为胸径、蓄积量和生物量生长方程中的第k个树种的参数，b1m、b2m、b3m(m=1,2,…,6)分别为胸径、蓄积量和生物量生长方程中的第m个立地质量等级的参数，c2k、c3k(k=1,2,3)分别为蓄积量和生物量生长方程中的第k个树种的参数，其余同上。

2.1.4 联合建模 为了比较建模方法对生长模型拟合优度和估测精度的差异，除分步建模法外，本研究还采用了联合建模法，以Model 7 为例。由于林龄未知，且3 个方程中年龄参数agei是一致的，基于式（2）～（4），参数处理同分步建模法，根据林分平均胸径、蓄积量和生物量与时间间隔之间的关系，采用非线性误差变量联立方程组法同时拟合3 个方程，即式（10）。

2.2 模型评价

考虑各模型参数不同，为了检验8 种模型系的拟合效果，采用以下4 项模型评价指标：确定系数R2、估计值标准误差SEE、平均预估误差MPE、总体相对误差TRE，公式详见文献[26]。

2.3 模型验证

（1）确定立地质量等级和林龄（以参数b 分级为例）

根据胸径生长方程反推林龄计算公式，如下：

式中，pageij为第i个样地第j期的样地林龄，当不引入竞争指数时b0=0，其余同上。

由于验证样本所在样地的立地质量等级m和初期林龄皆未知，需先确定立地质量等级，再由式（11）求出林龄，才可估计蓄积量和生物量。当立地质量分级时，至少需要两期胸径才能确定样地的立地质量等级，薛春泉等[27]发现长周期确定的等级比短周期稳定，因此本研究采用期初（2002年）和期末（2017 年）2 期胸径来确定样地的立地质量等级，确定方法如下：

式中，pagei1和pagei4为第i个样地期初和期末的样地林龄，该方法的直观解释就是使推算的2002 年和2017 年林龄间隔与实际林龄间隔（15 年）相差最小的立地质量等级，即为该样地的立地质量等级。

根据式（12）确定样地的立地质量等级后可知参数bm，再由式（11）计算初期（2002 年）林龄age，则2007 年、2012 年、2017 年时样地的林龄分别由初期林龄相应加上5、10、15 年。

（2）模型估测效果评价指标

由蓄积量和生物量生长方程估计样地平均蓄积量和生物量，样地总蓄积量和总生物量可由平均蓄积量和生物量乘以株数获得，按式（13）～（14）计算各树种4 个时期下区域尺度蓄积量和生物量估计的总相对误差EV和EBIO。

式中：n为样地数，Ni为第i个样地的保留木株数，PVi、Vi分别为第i个样地平均蓄积量的估计值和真实值，PBIOi、BIOi分别为第i个样地平均生物量的估计值和真实值。

2.4 生长量函数

连年生长量函数Z(t)为总生长过程曲线y(t)的一阶导数，为树木在某年的实际生长速度，即连年生长量随年龄t变化的函数。Richards 生长方程的连年生长量方程为：

平均生长量函数θ(t)用来表示树木在某一时刻的平均生长速度，总生长过程曲线y(t)除以年龄t则可得平均生长量依年龄变化的函数。Richards生长方程的平均生长量方程为：

3 结果

3.1 模型拟合

表4 列出了8 种模型系的拟合评价指标。整体上看，8 种模型系的3 个生长方程拟合确定系数R2皆达0.980 以上，总相对误差TRE在±0.5%以内，平均预估误差MPE在±1%以内，拟合效果良好，但各模型间的4 项指标差异不大。当控制建模方法和立地质量分级的参数相同、仅比较竞争影响时，引入竞争指数对模型拟合优度的提升没有明显的作用，如Model 5 和Model 6 拟合指标间差异甚微、Model 2 较Model 1 的蓄积量和生物量生长方程的拟合效果甚至出现了一定程度的下降；类似地，可以发现参数b分级模型的拟合优度较参数a分级模型稍好；当控制竞争影响和分级参数相同时，相对于分步建模，采用含度量误差的联立方程组时，蓄积量和生物量生长方程的拟合效果有所提高，但是胸径生长方程的拟合效果稍差。

表4 模型拟合评价指标Table 4 Model fitting evaluation index

3.2 模型验证

8 种蓄积量和生物量生长方程的模型验证指标分别如图2 和图3 所示。总体而言，各模型的预测效果都较好，从图2 可知，区域尺度蓄积量估计误差最大的仅有17.62%，说明蓄积量估计精度可达82.38%以上；从图3 可知，区域尺度生物量估计误差最大的仅有12.51%，则生物量估计精度可达87.49%以上。8 个模型中以基于参数b分级的分步模型（即Model 5 和Model 6）预测效果最优，其对3 个树种4 期蓄积量和生物量的估计误差较其余模型稍小，且在不同树种和不同时期的表现都较为稳定，说明模型适用性较好。以Model 5 为例，3 个树种中对杉木的估计效果稍好，估计误差最大为−5.25%，而估计马尾松和湿地松时误差最大为10.00%和10.36%。8 个模型对3 个针叶树种的4 期区域尺度蓄积量和生物量估计的误差，大都表现出中期误差高于首末2 期，初期估计精度较高，但明显出现低估的现象，伴随着样地林龄的增大，中期出现高估且精度稍低，但到末期（2017 年）时估计误差降低，精度有所提高。

图2 8 种蓄积量生长模型的区域尺度估计误差Fig. 2 Regional-scale estimation errors of eight volume growth models

图3 8 种生物量生长模型的区域尺度估计误差Fig. 3 Regional-scale estimation errors of eight biomass growth models

3.3 最优模型

综合8 个模型的拟合效果（表4）和验证效果（图2、图3），Model 5 最优，Model 6 次之。表5 和表6 列出了基于参数b分级的分步模型的模型参数估计值。由表5 可知，马尾松胸径、蓄积量和生物量的生长潜力均为最大，杉木最小，其中马尾松蓄积量和生物量的生长极限值可达湿地松和杉木的2 倍以上。由表6 可知，立地质量等级越高，树木的生长速度越快，但是各等级下的林分蓄积量生长速度均为最快，其次为生物量，胸径生长速度最慢。

表5 树种哑变量参数和竞争指数参数Table 5 Dummy variable parameters of tree species and competition index parameters

表6 基于参数b分级的立地质量等级哑变量参数Table 6 Site quality grade dummy variable parameters based on parameterbclassification

3.4 蓄积量和生物量生长量函数

基于参数b分级不考虑竞争的分步建模（即Model 5）最优，为节约篇幅，选择样地立地质量等级占比较多的等级2、3 和4，绘制Model 5 下蓄积量、生物量的连年生长量和平均生长量的关系曲线图（图4）。蓄积量和生物量的连年生长量、平均生长量均显示出立地质量越好，生长量能达到的极限值越大，且所需年限更短，各相邻等级间的差异随着等级的增加而逐渐加大。对于不同的树种来说，同一立地质量等级下的生长量均为马尾松>湿地松>杉木，生长量达到最大时的林龄为杉木>湿地松>马尾松，各树种间生长量差异随着等级的增加而逐渐加大。立地质量等级越高，平均生长量和连年生长量曲线相交时的样地年龄越小。

4 讨论

处理多树种问题是建立林分模型的一个难点，通常是将组成树种进行分组或合并[28]。曾伟生等[2]建立了东北林区林分水平的三储量模型，筛选出10 种主要森林类型并引入哑变量体现差异后，蓄积量模型、生物量模型和碳储量模型的确定系数分别从0.945、0.805、0.839 提高到0.959、0.949、0.951。本研究对存在多个组成树种的样地按优势树种分属不同的样地，并构建哑变量区别广东省的3 个主要针叶树种，从表5 和图4 可知，各树种间的生长潜力和生长量有较大的差异，表现为马尾松>湿地松>杉木，其中马尾松远超2 倍左右，这是因为建模数据（表1）中马尾松保留木大多为天然林，而杉木的起源大多为人工林，Zeng 等[29]研究结果表明同一直径的单木地上生物量估计值在不同起源上存在一定差异，天然林的生物量大于人工林。此外，在同一立地质量下，不同树种间生长量的差异主要是由于树种特性的差异引起的，不同树种对同一区域的生长适应程度并不一致，马尾松可能更适宜广东的立地条件。各树种间区域尺度蓄积量和生物量的预测精度也有较大的差异，对杉木的估计效果较马尾松和湿地松稍好，可能原因是不同树种的建模样本量和检验样本量不一致，另一个可能的原因是不同的树种有不同最优的基础生长模型，如吴宏炜等[30]研究发现，Schumacher 模型是湿地松林分蓄积量模型的最优基础模型。对于估计误差，大都表现出中期高于首末两期的现象，一方面可能是由于采用首末两期胸径数据来确定验证样本所在样地的立地质量等级，相较于中期，推算的各样地立地质量等级更贴合期初和期末时的实际情况，因此用于林龄确定、估计蓄积量和生物量时，3 个生长方程中反映立地质量差异的参数更为稳定，模型的预测效果会更好。另一方面随着林龄增大，蓄积量和生物量逐渐积累，公式（13）和（14）分母变大，这可能使后期的估计误差降低，估计效果更好。

非线性独立回归估计法的估计效果优于含度量误差变量联立方程组法，因为分步建模时，3 个方程都是独立拟合给出参数的最优估计，而采用联立方程组法联立求解时，由于方程中林龄一致的约束，为保证方程之间的联系，得到更优化的参数估计，各方程较独立拟合时都有了不同的变化。从拟合评价指标（表4）可以看出，对比分步建模法，采用联立方程组法可以提高蓄积量和生物量生长方程的拟合效果，但是却降低了胸径生长方程的拟合优度，导致林分年龄的估计精度有所下降，使用验证数据集检验模型的适用性时，由于林分年龄未知，确定首期林龄成为关键的一步，所以胸径生长方程的拟合效果直接影响了验证时的估计精度。由于林分大部分为混交林，具有年龄不一、树种结构关系复杂等特点，因此探究林木间竞争作用对林分结构和生长变化的影响至关重要。一般而言，模型中的参数个数越多，拟合效果会更优，但是引入竞争指数后对于模型的拟合效果（表4）和估计精度（图2、3）没有明显的提高。虽然从理论上讲包含树木空间位置的竞争指数对预测林木生长更有效, 但Biging 等[31]研究发现与距离无关的竞争指数比与距离有关的竞争指数预测单木生长更好，臧灏等[32]探究人工杉木林胸径生长对竞争的响应时也发现对比于与距离有关的单木竞争指标，采用与距离无关的单木竞争指标的生长模型估计精度更高。竞争指数表现的预测能力依赖于树种、立地条件等因素，这些因素的变化会导致不同类型的竞争指数表现不一致，且不具有决定性[33-34]。本研究引入竞争指数对林分模型的适用性没有明显影响，可能与选取的竞争指数有关，也有可能和立地的限制有关，因为在较大范围内密度对林分蓄积量、生物量的影响不大[35]。

模型的应用范围和模型的预测精度相关联，当模型应用于大尺度的森林资源测算时，若忽略立地质量的差异必然导致模型估计结果的不准确。立地质量是指在某一立地上既定森林或其他植被类型的生产潜力[35]，既定立地上的不同树种和相同的树种在不同的立地上的立地质量都是有差异的。本研究为体现广东省样地间的立地质量多样性，分别对与生长潜力有关的参数a和与生长速度有关的参数b进行分级，结果表明以参数b分级体现立地质量的差异优于参数a分级。薛春泉等[26]研究发现生长速度分级的胸径生长方程与地上生物量的异速生长方程构成的模型系，不仅拟合效果达到了0.95以上，区域尺度生物量的估计误差也下降到了8%以内，这与本研究参数b分级模型预测效果较好具有一致性。另外，林分蓄积量和生物量的生长量随着立地质量的提高而增加，说明立地质量等级越高，森林蓄积量和生物量密度越大。立地质量与森林生产力相关，不同的立地质量代表着不同的生产力，立地质量越高，林分的生长速率越快，达到最大生长量所需的年限也更短。赵菡等[36]用优势木树高等级代替立地等级估计江西省不同立地质量等级的马尾松林生物量，研究发现不同立地等级的地上生物量均值估计结果随着立地质量等级的升高而增大，这与本研究结果一致。

5 结论

本研究以广东省3 种针叶优势树种为研究对象，利用1997—2017 年的5 期森林资源连续清查数据，基于林分平均胸径、蓄积量和生物量与相应调查时间间隔的关系建立林分生长模型，实现了对广东省蓄积量和生物量的准确估计。本研究结论如下：（1）基于与生长速度有关的参数b分级不考虑竞争的分步模型最优，引入竞争指数与否差异不大。（2）最优模型的区域尺度蓄积量和生物量估计误差最大分别为10.36%和10.22%，模型的适用性较好，4 期的估计误差表现为中期高于首末两期。（3）不同树种中以马尾松的生长潜力最大，杉木最小，但是模型对杉木的估计效果较马尾松和湿地松稍好。（4）立地质量等级越高，生长量极值越大，达到最高峰所需年限也更短；同一立地质量下马尾松的生长量最大，其次为湿地松，杉木最小。