Development of Stand Density Management Diagram for Carbon Storage Management for Masson Pine Plantations

湖南省位于中国中部、长江中游，省境绝大部分在洞庭湖以南，故称湖南；湘江贯穿省境南北，故简称湘。湖南省位于东经108°47′～114°15′，北纬24°38′～30°08′之间，为大陆性亚热带季风湿润气候，冬寒冷而夏酷热，春温多变，秋温陡降，春夏多雨，秋冬干旱^[21]。气候的年际变化也较大。年均气温一般为16-19摄氏度，热量充足，雨水集中，日平均气温稳定，无霜期253-311天，年均降水量在1200-1700mm之间。湖南省是我国林业大省，属于亚热带常绿阔叶林带，森林资源丰富，2020年底，湖南省林地面积为1.298 × 10⁷ hm²，森林面积为1.002 × 10⁷ hm²，全省森林覆盖率大于59%^[22]，森林蓄积量达5.88 × 10⁸ hm³，湖南省天然林的主要组成树种为马尾松、杉木、栎类、樟树（Cinnamomum camphora Nees ex Wall）、木荷（Schima superba Gardner & Champ.）、枫香（Liquidamba formosana Hance）等；人工林的主要组成树种为马尾松、杉木、柏木（Cupressus funebris Endl.）、木荷、桉树（Eucalyptus robusta Smith）等^[23]。

本研究数据来源于湖南省第9期森林资源连续清查数据，包括样地信息和样地每木检尺数据。每木检尺数据针对样地内胸径≥ 5 厘米（DBH ≥ 5 cm）的乔木，测定并记录了样地号、样木号、树种、胸径等信息。样地信息记录样地号、样地类别海拔、坡向、坡位、坡度等信息。调查样地是面积为0.06公顷的方形样地，调查间隔期为5年。本研究对数据进行筛选，共计835块马尾松人工林纯林样地用于密度控制图研建。835块样地统计信息如表1所示。从表1中可知，林分密度N的平均值为358株·hm⁻²（15～2850株·hm⁻²），公顷断面积BA为5.74 m²·hm⁻²（0.03～31.9 m²·hm⁻²），林分平均胸径QMD为16.6cm（5～59.8 cm），公顷蓄积量VOL为35.474 m³·hm⁻²（0.094～229.304 m³·hm⁻²），碳储量CS为15.37 Mg·hm⁻²（0.0458～100.93 Mg·hm⁻²）。将生物量乘以生物量-碳的转换系数（也就是碳含量比例）从而得到碳储量，其中生物量-碳转化系数（即含碳率）取0.5 ^[24]，马尾松生物量方程如表2。

完满立木度林分样地的筛选采用Solomon^[26]提出的确定方法，在双对数坐标的条件下ln(N)-ln(QMD)的关系中，斜率为-1.605恒定值，之后利用最大ln(N)和ln(QMD)组合的样地数据，计算出截距a；然后用得出的自然稀疏模型计算每块纯林样地的最大林分密度（Nmax）；根据N/Nmax求出每块样地的相对密度指数（RD），当RD>0.7就可以视为完满立木度样地，用于构建自然稀疏模型。

本研究采用Reineke自然稀疏模型^[27]，该模型描述了林分平方平均胸径和每公顷株数的关系，模型的方程形式如下：

式中：N表示每公顷林木株数，QMD表示林分平方平均直径（cm），$ {\alpha }_{1} $为模型的截距系数。

普通最小二乘法（OLS）是线性回归模型最基本的参数估计方法，也常用于树木生长模型的估计^{[28, 29]}。然而在生态数据中的大多数变量存在多重共线性、离异值等缺点，很难满足这些要求^[30]；且普通最小二乘法（OLS）拟合出的通常是一条体现数据点分布的中心趋势线而非最大密度线，因此普通最小二乘法（OLS）拟合出来的模型参数有其局限性。分位数回归法（QR）对于异常点较为敏感，尤其是林分自疏线上限附近的异常点，且进行统计推断比较困难^[31]；部分学者认为RMA 法比OLS法更客观科学，更适合被用在拟合最大密度线上，其能够阐明植物个体大小与密度之间的真实关系^[32-34]。因此，本研究采用RMA对Reineke自然稀疏模型进行拟合。

为了能更准确反映马尾松人工林的生长状况，碳储量密度控制图除自然稀疏模型外，还有优势高模型、蓄积模型、碳储量模型及立地指数等4个辅助模型构成。本研究采用非线性最小二乘法（nonlinear least squares，NLS）拟合优势高模型、蓄积模型、碳储量模型。相关辅助模型的方程式如下所示：

式中：$ QMD $为林分平方平均直径（cm），$ N $为每公顷株数（株·hm⁻²），$ DH $为林分优势高（m），$ V $指材积（m⁻³）,$ CS $为碳储量（Mg·hm⁻²）,$ t $为林龄，$ a\mathrm{、}b\mathrm{、}c\mathrm{、}d $均为模型参数。

优势木树高根据Tang^[14]和 Castedo^[35] 以每公顷100株最大的优势木平均高来确定，本研究的样地为一亩地，所以样地内取7株胸径最大的优势木的平均高来计算优势木树高^[36]。由于一类数据没有观测林木树高，本研究采用杜志^[37]建立的湖南省马尾松树高方程计算林木树高，具体方程式如下：

式中：H为林木树高（m），DBH为林木胸径（cm）。

本研究利用刘燕英等^[38]建立的马尾松人工林立地指数模型，模型具体形式如下：

式中：Hd为优势木高（m），t为林龄。

本研究采用均方根误差（$ RMSE $）、绝对偏差（Absolute bias）、决定系数（$ {R}_{adj}^{2} $）3个统计量对模型的拟合优度进行评价，相关统计量公式如下所示：

式中，$ L $是对数似然函数的最大值；$ p $是模型参数个数；$ n $是样本单元的数量；$ {y}_{i} $和$ {\widehat{y}}_{i} $分别是因变量的观察值和预测值，$ \stackrel{-}{y} $是观察值的均值。

密度控制图在坐标轴上进行绘制，其中以林分平均胸径QMD为横坐标，以每公顷株数密度N为纵坐标。Vacchiano^[8]建议将自然稀疏模型的截距向上平移到其95%置信区间的上限，作为最大密度线。相关学者^{[13, 39, 40]}指出，一般将最大密度线的60%作为最优密度的上限，保证了林木不会因为密度过高而出现激烈竞争，导致林木个体死亡；最大密度线的30%作为最优密度的下限，保证林木能够具有充足的生长空间，林地空间得到有效利用；最大密度线的20%作为林分郁闭值，保证林分的郁闭状态。将上述3条密度线、等优势高线、等蓄积量线、等碳储量线叠加到坐标系上，生成最终的碳储量密度控制图。本研究采用R语言的ggplot2包进行密度控制图绘制。

当林分密度位于合理密度区间范围内（最大密度线的30%-60%），此时的林分中的林木个体处于互利竞争状态，即不会发生因为对自然资源的恶性竞争，所导致的自然稀疏现象。此时的林分处于一种合理、稳定、健康的状态；当林分密度高于最优密度上限时（最大密度线的60%），此时的林木个体之间会发生对温度、光照、养分等自然资源的互斥性竞争，导致林木死亡现象的发生，造成了碳释放；当林分密度低于最优密度下限时（最大密度线的30%），会造成林地空间资源的浪费。本研究以湖南省慈利县马尾松林为例，基于二类调查数据马尾松小班林分密度，将位于最大密度线30%-60%的林分，定义为合理密度林分，将位于最大密度线的30%以下的林分定义为低密度林分，将位于最大密度线的60%以上的林分定义为高密度林分，最终构建慈利县马尾松人工林密度状态分布诊断图。

基于马尾松人工纯林835个样地数据，采用Solomon提出的完满立木度样地确定的方法，计算每一块样地的相对密度指数（RD），当相对密度指数（RD）大于等于0.70，就属于完满立木度样地。本研究最终从835块马尾松相对纯林中筛选出54块完满立木度样地（如图1）。采用RMA回归，拟合得到湖南马尾松纯林Reineke自然稀疏模型：

Figure 1.  Natural sparse line and perfect standing degree plot.

式中R²为0.926 ，RMSE为0.127。

基于马尾松人工纯林835个样地数据，对优势高模型（式2）、蓄积量模型（式3）、碳储量模型（式4）进行了拟合，模型参数估计结果与拟合优度统计量如表3、表4所示。辅助模型的参数估计值均表现出显著性（p < 0.05），决定系数R²_adj均在0.90以上，均方根误差RMSE和绝对偏差Bias较小（表4），说明马尾松人工纯林数据对三个辅助模型方程的拟合效果较好。

马尾松蓄积模型：

马尾松优势高模型：

马尾松碳储量模型：

将自然稀疏线方程（式10）水平上移到截距的95%置信区间的上限，得到最大密度线。将最大密度线，分别乘以60%，30%，20%，得到最优密度上限、最优密度下限及林分郁闭度线，其方程式如下所示：

基于上述密度线模型，将最大密度线、最优密度线上限、最优密度线下限、林分郁闭线添加到对数轴上，密度平行线都采用红色线条。此外，基于辅助模型，添加等优势高线、等蓄积量线、等碳储量线。蓝色线条表示12、15、20、25 m的林分等树高线，绿色线条表示50、100、200、300 m³·hm⁻²的等蓄积量线，紫色线条表示10、25、50、100 Mg·hm⁻²的等碳储量线。马尾松人工林碳储量密度控制图如图2所示。

Figure 2.  Density management diagram for carbon sequestration of Pinus massoniana plantation.

将林分平均胸径值带入最大密度线、最优密度线上限、最优密度线下限、林分郁闭线方程式，得到马尾松人工纯林不同林分平均胸径下的株数密度阈值（表5）。例如，当林分平均胸径为10 cm时，林分最大密度为3231 株·hm⁻²，最优密度区间（最大密度线的35%～60%）为1130-1938 株·hm⁻²，林分郁闭时（最大密度线的20%）的密度值为645 株·hm⁻²。

基于慈利县二类调查数据，采用所构建的马尾松碳储量密度控制图，对慈利县马尾松人工林密度合理性进行诊断，生成湖南慈利县马尾松人工林密度诊断图（图3）。淡黄色区域代表非马尾松林地，其余 3 个颜色是马尾松林地。其中，绿色表示最优密度区间的马尾松人工林，即林分的密度在最优密度上、下限之间；蓝色表示马尾松人工林的密度较大，高于最优密度上限；红色表示马尾松人工林的密度较小，低于最优密度下限。从马尾松密度诊断图可以看出，高密度的马尾松人工林所占比例占压倒性优势，急需开展密度调控。

Figure 3.  Density warning diagram of Pinus massoniana

本研究基于湖南省第9期森林资源连续清查数据，构建了湖南马尾松人工林碳储量密度控制图，根据构建出的密度控制图，利用慈利县二类调查样地数据，研建了慈利县马尾松人工林的密度诊断图，通过密度诊断图可以对慈利县马尾松人工林密度情况进行诊断。结果显示：（1）绝大部分马尾松人工林的密度高于合理密度区间，林木就会对有限的资源产生竞争，甚至有些林分会因为竞争产生严重的自然枯死现象；（2）少部分林分位于合理密度区间范围内，此时的林分能够充分利用光照、温度、水分等自然资源条件，同时林木个体之间处于互相促进的状态；（3）少部分林分密度低于最优密度下限，此时的林分不能充分利用林地空间，造成林地资源的浪费。通过抚育间伐、林下补植等措施，将林分密度调整到合理的密度是促进人工林正向发展的关键措施^[41]。本研究针对慈利县马尾松人工林密度问题，提出如下建议：

（1）高密度林分（高于最优密度的上限）。高密度林分蓄积量集中在100 m³·hm⁻²-300 m³·hm⁻²之间，其碳储量集中在50 Mg·hm⁻²-100 Mg·hm⁻²之间。对于高密度的中幼龄林，采用目标树作业体系，选择目标树，伐除目标树周围的干扰树，为目标树做自由冠，同时伐除林内IV、V级立木，降低林分密度。此外，在产生的林窗下补植红锥、香樟、桢楠等乡土阔叶树种，营造近自然异龄复层混交林。Andrew^[42]研究发现楠木、青冈、黄檀、香樟等慢生树种的固碳能力远高于速生树种；对于高密度的近熟林，如维持同龄林作业体系，则可以开展渐伐作业，在一个龄级期内分批次收获主林层林木，并同时完成林下造林。此外，也可以采用目标树作业体系，将当前林分改造为异龄复层混交林；对于高密度过成熟林，则急需开展皆伐作业，避免出现因森林自然成熟而导致的碳排放。森林自然成熟是当林分或树木生长到开始枯萎阶段时的状态，出现负生长，是森林收获的底线，自然成熟的不会再储存额外的碳^{[43, 44]}。此外，值得说明的是抚育间伐使活立木的碳转移到木材产品中去，收获的林木虽然减少了森林生态系统碳储量，但是只要木材产品不被燃烧，则不被视为森林碳排放^[45]。

（2）位于合理密度区间范围内的林分（高于最优密度的下限且低于最优密度上限）。合理密度林分蓄积量多集中在50 m³· hm⁻²-200 m³· hm⁻²之间，其碳储量集中在25 Mg·hm⁻²-100 Mg·hm⁻²之间。如果继续维持当前的皆伐作业体系，此时的林分不需要采取任何林分密度调控措施，但是需要根据林分密度控制图，及时监测未来森林密度情况，当超出最优密度范围，则需要及时开展森林抚育，将林分密度维持在合理区间范围内；如果计划将当前单层马尾松纯林导向异龄复层混交林，则需要选取目标树，低密度采伐影响目标树生长的干扰树，在产生的林窗下补植红锥、香樟、桢楠等乡土树种，将当前林分逐步导向以乡土树种为主体的异龄复层混交林。

（3）低密度林分（低于最优密度的下限）。对于低密度林分，首先，建议在林下及时补植红锥、香樟、桢楠等乡土树种，必要时辅助于扩堰、松土等措施，促进其尽快进入主林层。其次，在主林层选择干形通直、生活力旺盛的林木作为目标树，对目标树进行适当修枝等作业措施，促进目标树的快速生长。

本研究所构建的马尾松碳储量密度控制图是典型的全林分模型，仅从林分密度方面给出森林经营决策信息。因此，适应于马尾松人工纯林，而对于以马尾松为优势树种的多树种混交林并不适应。国内外众多研究表明，生物多样性与森林固碳能力呈现正相关^[46-51]例如，Mensah^[49]通过研究南非云雾带森林生物多样性和生态系统功能之间的关系发现，树种多样性（丰富度）与地上碳储量呈正相关。Huang^[51]发现，物种多样性对亚热带森林的碳汇具有积极的影响，并强调了混交林在缓解气候变化和维护生物保护方面的重要性。此外，异龄复层混交林在生物多样性保护^[52]、水源涵养^{[53, 54]}、森林游憩^[55]等生态服务功能上，对于人工纯林具有不可比拟的优势。因此，构建多树种混交林密度控制图显得尤为必要。国外不乏多树种混交林密度控制图的案例。例如，Edwin^[56]构建了云杉和冷杉混交林的密度控制图，主要目的是在混交林中如何实施有效的经营管理工作，为我国类似地区的森林经营管理提供参考依据；Mehari^[57]研制了非洲山地冷杉和罗汉松混交天然同龄林的密度控制图，在实施疏伐管理时，既要提高森林的生产力，也要维护生态平衡和生物多样性。

此外，本研究密度控制图的功能仅聚焦于碳储量、蓄积等目标，而密度控制图可以解决林冠火规避^[58]、动物栖息生境构建^[59]、水源涵养^[15]以及防风灾能力提升^[60]等问题，例如Barrio^[58]通过使用密度控制图，可以有效地监测和控制森林密度，降低林冠火的发生概率；John^[59]发现啄木鸟种群的恢复和整个共同分布区的长叶松树林的恢复之间有着重要的联系，构建了为维持濒危啄木鸟的栖息地的密度控制图。因此，建议将森林的其他功能的等密度线加入到密度控制图，构建多目标密度控制图，有助于更全面地评估森林的可持续发展状况，为我国森林资源的经营管理提供科学依据。

林分密度控制图是森林可持续经营管理的一个有效工具，本研究通过林分密度状况与固碳能力之间的量化关系，合理调控湖南省马尾松人工林的经营密度，使其向最优密度范围内发展，有利于提升森林的固碳能力，发挥森林多重功能与效益。对于森林经营管理者来说，利用林分密度控制图确定合理密度区间，对于低密度林分，及时补植固碳能力强的树种，对于高密度林分，确定合理的采伐强度促进碳储量的增长，精准确定采伐强度，对于提高不同林分的碳汇功能具有重要意义，为其经营管理森林提供决策支持依据。除此之外，随着人工林的多目标经营，单一目标的密度控制图向多功能密度控制图转变，根据不同地区不同树种的具体情况来构建才能满足当地经营目标的需求，能够更精准的预测林分生长状态，为森林经营者提供更优的抚育方案，对提高森林结构多样性和功能具有重要的意义。