-
为保持天然林的多元结构和功能，天然林的收获通常采用择伐的方式[1]，但是由于天然林异龄混交的复杂性，择伐木的选择并无统一定论[2-4]。近年来一些学者以调整林分空间结构为采伐目标对此问题进行了研究[5-7]。陈亚南[5]通过构建林木的空间结构多样性综合指数评定采伐木。惠刚盈等[6]也提出依据结构化经营理论进行林木格局调整确定采伐木。Alexander通过量化林木的生命活力和成长空间选取了采伐木[7]。这些研究主要是从发挥林分的生态效益和林分健康的角度出发制定采伐方案，而没有将林分的生产力指标纳入到经营方案的设计中。
林地作为一种有限的资源，使其发挥较高的生产力并保持长期稳定应该是其经营目标之一[1]，因此采伐木的选取也应以此作为筛选标准之一。人工林常以林木的数量成熟，即连年生长曲线与平均生长曲线的交点作为采伐依据，因为在这个时间点采伐、更新，即可获得最大木材收获量，使每株林木的位置都得到最大化的利用。但是在天然林中，不同树种占据的营养空间不同，不再适宜用单位时间的生长量代表其对林分空间的利用效率。
林木的生长主要来源于光合作用，而树冠是光合作用的主要场所，所以林木在地上部分对营养的竞争主要是树冠占据的林分空间[8-9]，因此可以以树冠的大小来衡量林木所占据的营养空间。Arney[10]认为林木的最大生长空间面积等于具有相同胸径的自由树的树冠面积，并据此提出了竞争压力指数。其他学者也纷纷将树冠因子加入到竞争指标中预测林木生长情况，拟合结果均得到显著提高[11-13]。
本文用林木生长量和树冠构建林木空间利用率指标，以长白山地区天然云冷杉针阔混交林为研究对象，建立林木空间利用率线性模型。同时，为消除林业数据中广泛存在的自相关异方差等问题[14]，提高模型对混交林的适应性，本文在线性模型的基础上进一步建立按树种分组的混合效应模型，由此拟合各树种的空间利用率，从而为天然林单木成熟的判定提供依据。
-
研究地区位于吉林省汪清林业局东北部的金沟岭林场（130°10′E、43°22′N），海拔300 ~ 1 200 m，地貌属低山丘陵，坡度多在5° ~ 25°。全年平均气温约3.9 ℃，积温2 144 ℃；年降水量600 ~ 700 mm，生长期为120 d。土壤多为针叶林灰棕壤，粒状结构，湿松，根系多，平均厚度在40 cm左右。植被丰富，立地较好。
研究区主要为天然针阔叶混交过伐林。主要树种有：云杉（Picea koraiensis）、冷杉（Abies holophylla）、红松（Pinus koraiensis）、枫桦（Betula costata）、椴树（Tilia amurensis）等，其他树种还有：色木槭（Acer mono）、水曲柳（Fraxinus mandschurica）、胡桃楸（Juglans mandshurica）、黄菠萝（Phellodendron amurense）、白桦（Betula platyphylla）、青楷槭（Acer tegmentosum）、花楷槭（Acer ukurunduense）等。
本研究采用了20块云冷杉针阔混交林永久样地2006—2010年间的调查数据。该区域主要是云冷杉林群落，是原始林型在不同强度的择伐后[15]，经过20多年的恢复形成的林分，林分蓄积量在150 ~ 400 m3/hm2（原始林为350 ~ 400 m3/hm2），林分株数分布以小径木为多。样地设置于1986/1987年，面积0.12 ~ 0.2 hm2不等，记录了坡度坡向等立地因子，并进行每木检尺及定位，此后每隔2到3年复测1次。研究选用了一个调查间隔期期初和期末的2 268株样木数据。随机抽取2/3用于建模，其余用于模型检验。由于天然混交林中阔叶树树种多而株数少，本研究将阔叶树分为中阔（白桦、榆树和杂木）和慢阔（色木槭、水曲柳、椴树和枫桦）两个树种组[16]。样木统计信息如表1所示。
表 1 建模数据和检验数据主要指标统计信息
Table 1. Summary statistics for modeling and validation data sets
数据
Data树种组
Tree species group株树
Plant number统计指标
Statistic index胸径
DBH/cm树高
Tree height/m平均冠幅
Average crown width/m样地蓄积/（m3·hm− 2）
Sample plot volume/ (m3·ha− 1)针叶树比例
Conifer percentage建模数据 Modeling data 红松
Pinus koraiensis219 平均值 Mean 24.4 18.4 2.2 296.5 0.91 标准差 SD 8.4 3.9 0.7 73.7 0.06 云杉
Picea koraiensis271 平均值 Mean 22.1 17.9 2.1 288.3 0.88 标准差 SD 9.7 5.2 0.7 74.1 0.07 冷杉
Abies holophylla717 平均值 Mean 21.5 18.7 1.8 297.6 0.90 标准差 SD 7.9 4.4 0.6 76.9 0.06 落叶松
Larix olgensis80 平均值 Mean 22.0 20.7 2.2 263.6 0.90 标准差 SD 7.4 4.5 0.7 62.3 0.07 中阔
Medium-speed broadleaf trees57 平均值 Mean 19.5 19.7 2.0 333.0 0.87 标准差 SD 6.6 4.5 0.9 53.2 0.05 慢阔
Slow broadleaf trees168 平均值 Mean 17.2 16.8 2.2 298.3 0.87 标准差 SD 6.6 4.7 0.8 71.4 0.06 检验数据
Validation data红松
Pinus koraiensis122 平均值 Mean 24.2 18.9 2.1 294.5 0.91 标准差 SD 7.8 4.0 0.7 71.7 0.06 云杉
Picea koraiensis141 平均值 Mean 21.6 17.9 2.1 290.8 0.88 标准差 SD 8.6 5.0 0.7 76.0 0.07 冷杉
Abies holophylla357 平均值 Mean 22.8 19.3 1.9 305.4 0.90 标准差 SD 7.8 4.4 0.6 72.4 0.06 落叶松
Larix olgensis26 平均值 Mean 23.1 21.4 2.5 246.2 0.88 标准差 SD 6.7 3.9 0.8 51.4 0.08 中阔
Medium-speed broadleaf trees21 平均值 Mean 21.9 21.0 2.3 329.7 0.88 标准差 SD 6.0 2.9 0.9 58.4 0.06 慢阔
Slow broadleaf trees89 平均值 Mean 17.0 17.0 2.2 297.7 0.86 标准差 SD 6.5 4.6 0.8 78.7 0.06 注：SD为标准差。Note: SD means standard deviation. -
为评价单株林木对单位营养空间的利用效率，本研究用定期生长量与营养空间之比来定义林木的生产力指标——林木空间利用率。
胸高断面积增长量和材积增长量是衡量林木生长量的常用指标。树冠占据的林分空间主要是二维的水平空间和三维的立体空间，也就是树冠的投影面积和树冠体积。因此分别考察胸高断面积增长量和材积生长量与树冠投影面积和树冠体积的比值。具体计算方式如下：
$$g = {\text{π}} \cdot {\rm{C}}{{\rm{W}}^2}$$ (1) $$ {V}_{\mathrm{c}}={f}_{\mathrm{c}}g\mathrm{C}\mathrm{L} $$ (2) $$ {B}_{\mathrm{g}}={B}_{\mathrm{z}}/g $$ (3) $$ {B}_{\mathrm{v}}={B}_{\mathrm{z}}/{V}_{\mathrm{c}} $$ (4) $$ {V}_{\mathrm{g}}={V}_{\mathrm{z}}/g $$ (5) $$ \begin{aligned} &\; \\ &{V}_{\mathrm{v}}={V}_{\mathrm{z}}/{V}_{\mathrm{c}} \end{aligned} $$ (6) 式中：CW为东西南北4个方向的平均冠幅；
$ g $ 为树冠投影面积；$ {f}_{\mathrm{c}} $ 为冠型指数，针叶树近似为圆锥体，冠型指数为1/3，阔叶树近似为抛物线体，其冠型指数为1/2[17]；CL为冠长；$ {V}_{\mathrm{c}}$ 为树冠体积；Bz、Vz分别为林木胸高断面积、蓄积年生长量；Bg、Bv、Vg、Vv是4个待定的林木空间利用率指标。胸径是林木最常用也最易测得的指标。为了空间利用率指标的拟合和应用，分别考察4个待定指标与胸径的关系以确定最佳指标。
-
林木生长主要受4个方面影响，即林木大小、林分特征、立地条件和竞争状态。林木空间利用率是由林木材积生长量和冠幅构成的指标，二者的生长皆受这4个方面的影响，故林木空间利用率指标也以这4组因子作为候选自变量拟合。
林木大小主要包括期初胸径DBH、胸高断面积（直接表示为胸径的平方DBH2）、树高H、平均冠幅CW及冠长CL。林分特征包括样地蓄积S、样地胸高断面积之和BA、林分密度指数SDI[18]、针阔比PF（本文以针叶树所占比例表示）。立地条件包括样地平均高Hm、优势木高Ha、坡度和坡向。坡度坡向利用Stage提出的转化公式[19]为：
$$ {\rm{SLC}} = {\rm{SL}} \cdot {\rm{cosASP}} $$ (7) 式中：SL为坡度；ASP为坡向。
竞争指标则选用了CIS[20]，选取高于对象木且距离对象木小于6 m的林木为竞争木，计算公式为：
$${{\rm{CIS}}_{i}} = \mathop \sum \limits_{j = 1}^{{r}} \left. {\left[ {\frac{{{D_j}}}{{{D_i}}} \frac{{{D_j}/{d_{ij}}}}{{ \mathop \sum \limits_{j = 1}^{{r}} {D_j}/{d_{ij}}}}}\; \right]}\right/\frac{{\left( {{Z_i} + 1} \right)}}{{{\rm{CW}}}}$$ (8) 式中：Di为对象木i期初的胸径， Dj 为竞争木 j 期初的胸径，r为竞争木个数，dij为对象木i与竞争木j之间的距离，CW为对象木i期初的冠幅，Zi为对象木i的胸径生长量。
为使因变量林木空间利用率指标更符合正态分布，模型更加精确，建模前先用Box-Cox变换对其进行转化[21]：
$$f\left( {x,{\rm{\lambda }}} \right) = \left\{ {\begin{aligned}& {({x^\lambda } - 1)/\lambda,\lambda \ne 0}\\& {\log x,\lambda = 0} \end{aligned}} \right.$$ (9) $ \lambda $ 的值通过R软件forecast包的BoxCox.lambda函数计算得到。 -
为了增强模型的适应性，本文以树种分组，加入随机效应及误差项方差协方差结构，建立林木空间利用率混合效应模型。
线性混合效应模型的形式[22]为：
$$ \left\{{\begin{aligned}&{{{y}}_{{i}}={{X}}_{i}{{\beta}} +{{Z}}_{i}{{b}}_{i}+{{\varepsilon }}_{i}, i={1,2},\cdots,m}\\& {{{b}}_{i} \sim N\left(0,{D}\right),{{\varepsilon }}_{i} \sim N\left(0,{{R}}_{i}\right)}\end{aligned}}\right. $$ (10) 式中：
$ {y}_i $ 是第i个树种组中$ \left({n}_{i}\times 1\right)$ 维的样木观测值，$ {n}_{i} $ 是第 i 个树种组的观测株树；m是树种组个数；Xi是$ \left({n}_{i}\times p\right)$ 维已知设计矩阵，p为固定参数个数；β是$ \left(p\times 1\right)$ 维固定参数向量；Zi是$ \left({n}_{i}\times q\right)$ 维随机效应的设计矩阵，q为随机参数个数；bi是$ \left(q\times 1\right)$ 维随机参数向量；D是随机参数的方差协方差矩阵；$ {\varepsilon }_i $ 是$ \left({n}_{i}\times 1\right)$ 维模型的误差项；Ri是树种组内方差协方差结构。目前常用的误差方差协方差结构为：
$$ {{{R}}}_{i}{={\sigma }_{i}^{2}{{\psi}} }_{i}^{0.5} {{{\varGamma }}}_{i} {{{\psi}} }_{i}^{0.5} $$ (11) 式中：
$ {\sigma }_{i}^{2} $ 为误差扩散的比例因子，$ {{{\varGamma}} }_{i} $ 为树种组内误差相关性结构；$ {{{\psi}}}_{i} $ 为描述方差异质性的对角矩阵。由于本次研究数据均为一个调查期内，不存在典型的自相关问题，故相关性结构$ {{{\varGamma}}}_{i} $ 采用单位矩阵[8]。 -
基础模型的拟合和检验结果通过均方根误差（RMSE）和调整决定系数（
$ {R}_{{\rm{a}}}^{2} $ ）评价[14]。混合模型的拟合和检验结果通过赤池信息量准则（AIC）、贝叶斯信息准则（BIC）两个统计量评价[23]。$$ \mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}=\sqrt{\frac{\displaystyle\sum \limits_{i=1}^m{\displaystyle\sum \limits _{j=1}^{n_i}}{\left({y}_{ij}-{\hat{y}}_{ij}\right)}^{2}}{n-p-1}} $$ (12) $$ {R}_{{\rm{a}}}^{2}=1-\frac{\displaystyle \sum \limits_{i=1}^m{\displaystyle \sum \limits _{j=1}^{n_i} }{\left({y}_{ij}-{\hat{y}}_{ij}\right)}^{2}}{\displaystyle \sum \limits_{i=1}^m{\displaystyle \sum \limits _{j=1}^{n_i}}{\left({y}_{ij}-\bar y\right)}^{2}} \frac{n-1}{n-p-1} $$ (13) $$ \mathrm{A}\mathrm{I}\mathrm{C}=-2\mathrm{l}\mathrm{o}\mathrm{g}\mathrm{l}\mathrm{i}\mathrm{k}\mathrm{e}\mathrm{l}\mathrm{i}\mathrm{h}\mathrm{o}\mathrm{o}\mathrm{d}+2p $$ (14) $$ \mathrm{B}\mathrm{I}\mathrm{C}=-2\mathrm{l}\mathrm{o}\mathrm{g}\mathrm{l}\mathrm{i}\mathrm{k}\mathrm{e}\mathrm{l}\mathrm{i}\mathrm{h}\mathrm{o}\mathrm{o}\mathrm{d}+p\mathrm{l}\mathrm{o}\mathrm{g}n $$ (15) 式中：
$ {y}_{ij} $ 为观测值；$ {\hat{y}}_{ij} $ 为预测值；$ \bar y $ 为观测值的平均值；$ {n}_{i} $ 为树种组内的观察株数；m为树种组个数；n为样本数；p为模型中参数的个数；loglikelihood为模型的极大似然函数。模型的拟合和检验通过R软件的nlme包完成。
-
公式（3） ~ （6）的4个指标与林木胸径的相关系数分别为0.351、0.039、0.467、0.109，即相关性最高的为材积增长量与树冠投影面积之比Vg（见图1）。另外，考虑到材积增长也正是森林经营的最终目标，所以选定Vg为研究对象，即下文所指的林木空间利用率。
图 1 林木空间利用率与胸径的关系
Figure 1. Relationship between forest space utilization and DBH
-
计算得到式（9）中λ = 0.35。因此，因变量为：
$$ {V_{\rm{t}}} = \left( {{V_{\rm{g}}}^{0.35} - 1} \right)/0.35 $$ (16) 对Vt和建模的4组因子进行逐步回归，最终得到模型的基础形式如下：
$$\begin{split} {V_{\rm{t}}} =& {b_0} + {b_1}{\rm{DBH}} + {b_2}{\rm{DB}}{{\rm{H}}^2} + {b_3}H + {b_4}{\rm{CW}} +\\& {b_5}S + {b_6}{\rm{PF}} + {b_7}{\rm{SLC}} + {b_8}{\rm{CIS}} + {\rm{\varepsilon }}\end{split} $$ (17) 式中：DBH为胸径，H为树高，CW为冠幅，S为样地每公顷蓄积，PF为针叶树比例，SLC为坡向坡度转化式，CIS为竞争指标，ε为估计误差。
经检验，所有保留变量的方差膨胀因子均小于10，且回归系数显著（P < 0.01）。拟合结果如表2所示。
表 2 林木空间利用率基础模型模拟结果
Table 2. Simulating results of basic model for forest space utilization
变量
Variable估计值
Estimating value方差
Variancet P (Pr > |t|) 截距 Intercept (INT) − 2.767 04 0.019 55 − 141.542 < 0.000 1*** 胸径 DBH 0.013 76 0.000 76 18.184 < 0.000 1*** 胸径平方 Square of the DBH (DBH2) − 0.000 18 0.000 01 − 13.335 < 0.000 1*** 树高 Tree height (H) 0.001 26 0.000 35 3.553 0.000 39*** 冠幅 Crown width (CW) − 0.051 17 0.001 71 − 29.975 < 0.000 1*** 蓄积 Stock volume (S) − 0.000 14 0.000 01 − 9.249 < 0.000 1*** 针阔比 Proportion of conifer to broadleaved trees (PF) 0.065 63 0.016 59 3.955 < 0.000 1*** 坡度坡向 Slope degree and slope aspect (SLC) 0.001 36 0.000 42 3.245 0.001 19 ** 竞争指数 Competition index (CIS) − 0.005 58 0.000 94 − 5.928 < 0.000 1*** 注：***表示在P < 0.001水平下差异显著，**表示在P < 0.01水平下差异显著。Notes: *** represents significant difference at P < 0.001 level, ** represents significant difference at P < 0.01 level. -
用前述基本模型（公式（17））分别拟合各树种组数据，得到各树种组固定效应系数的置信区间，如图2所示。从图2结果可以看出，树种组间DBH、DBH2、H、SLC 4个变量参数的置信区间没有完全重叠的部分，所以为这4个变量增加随机效应以扩展模型的适用性[24]。加入4个变量及截距不同组合形式的随机效应进行拟合，收敛的结果如表3所示。比较评价指标可以看出，在截距（INT）和DBH2两项上加入随机效应的模型VTLME5拟合结果最优。
图 2 各树种组线性模型拟合参数的95%置信区间
Figure 2. 95% confidence intervals on the linear model parameters for each tree species group
表 3 不同随机效应拟合结果统计检验
Table 3. Goodness-of-fit statistics for mixed model of different random effects
模型编号
Model No.随机参数
Random parameterAIC BIC loglikelihood VTLME1 INT − 5 142.105 − 5 083.638 2 582.053 VTLME2 INT DBH − 5 161.782 − 5 092.684 2 593.891 VTLME3 INT H − 5 145.099 − 5 076.001 2 585.550 VTLME4 INT SLC − 5 150.207 − 5 081.109 2 588.104 VTLME5 INT DBH2 − 5 165.259 − 5 096.161 2 595.630 VTLME6 INT DBH2 SLC − 5 172.721 − 5 087.678 2 602.360 注：VTLME1、VTLME2、VTLME3、VTLME4、VTLME5、VTLME6分别为加入相应随机参数的混合模型编号。Notes: VTLME1，VTLME2，VTLME3，VTLME4，VTLME5，VTLME6 are the mixed model codes with corresponding random parameters. 在模型VTLME5中加入不同形式的残差方差结构，结果如表4所示。模型评价指标显示幂函数形式的残差方差结构拟合结果最优。
表 4 不同残差方差结构拟合结果统计检验
Table 4. Goodness-of-fit statistics for mixed model of different residual variance structure
模型编号
Model No.残差方差模型
Residual variance modelAIC BIC loglikelihood VTLME5 无 None − 5 165.259 − 5 096.161 2 595.630 VTLME51 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma {\mathrm{D}\mathrm{B}\mathrm{H}}^{2}\right)$ − 5 167.432 − 5 093.019 2 597.716 VTLME52 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma \mathrm{C}\mathrm{W}\right)$ − 5 175.094 − 5 100.681 2 601.547 VTLME53 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma \mathrm{C}\mathrm{I}\mathrm{S}\right)$ − 5 174.558 − 5 100.145 2 601.279 VTLME54 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma f\right({x}_{i}\left)\right)$ − 5 185.268 − 5 110.855 2 606.634 VTLME55 $\scriptstyle {\sigma }^{2}{\mathrm{D}\mathrm{B}\mathrm{H} }^{2 \times2\gamma } $ − 5 170.525 − 5 096.112 2 599.263 VTLME56 $\scriptstyle {\sigma }^{2}{\mathrm{C}\mathrm{W}}^{2\gamma } $ − 5 175.849 − 5 101.436 2 601.925 VTLME57 $\scriptstyle {\sigma }^{2}{\mathrm{C}\mathrm{I}\mathrm{S}}^{2\gamma } $ 不收敛 Not converge VTLME58 $\scriptstyle {\sigma }^{2}{f\left({x}_{i}\right)}^{2\gamma } $ − 5 185.182 − 5 110.769 2 606.591 注：VTLME5、VTLME51、VTLME52、VTLME53、VTLME54、VTLME55、VTLME56、VTLME57、VTLME58分别为加入相应残差方差模型的混合模型编号；$\scriptstyle f\left({x}_{i}\right)$为固定效应及随机效应的预测值。Notes: VTLME5,VTLME51, VTLME52, VTLME53, VTLME54, VTLME55, VTLME56, VTLME57, VTLME58 are the mixed model codes with corresponding residual variance models; $\scriptstyle f\left({x}_{i}\right)$ represents the predicted value of the fixed effect and the random effect. 最终得到的混合效应模型如下：
$$ \left\{ \begin{array}{l} V_{{\rm{t}}_{ij}} = - 2.778 \;4{\rm{}} + {b_{0i}} + 0.015\;3{\rm{DB}}{{\rm{H}}_{ij}} + \\ \quad\quad\;\; \left( { - 0.000\;2 + {b_{2i}}} \right){\rm{DBH}}_{ij}^2 + 0.000\;9{H_{ij}} - \\ \quad\quad\;\; 0.047\;0{\rm{C}}{{\rm{W}}_{ij}} - 0.000\;1{S_{ij}} + 0.044\;5{\rm{P}}{{\rm{F}}_{ij}} + \\ \quad\quad\;\; 0.001\;6{\rm{SL}}{{\rm{C}}_{ij}} - 0.004\;7{\rm{CI}}{{\rm{S}}_{ij}} + {\varepsilon _{ij}} \quad\quad\\ {\left( {{b_{0i}} {b_{2i}}} \right)^{\rm{T}}} \sim N\left( {0,{{D}}} \right),{{D}} = \left({\begin{array}{*{20}{c}} {0.012\;4 \;\; \; - 0.005}\\ {\; - 0.005 \;\;\;\; 0.000\;02} \end{array}} \right) \\ {\varepsilon _{ij}} \sim N\left( {0,{{{R}}_i}} \right),{{{R}}_i} = 9.22_{}^2 {{\psi}} _i^{0.5} \times {\bf{I}} \times {{\psi}} _i^{0.5}, \\ {{{\psi}} _i} = {\rm{exp}}\left( {2 \times 2.045 \times f\left( {{{{x}}_i}} \right)} \right) \end{array} \right.$$ 式中：i为树种组编号；j为样木编号；
$ {b}_{0i} $ 、$ {b}_{2i} $ 为样地层次随机效应参数；I为单位矩阵；$ {{{x}}}_{i} $ 为自变量，$ f\left({{{x}}}_{i}\right)$ 为固定效应及随机效应的预测值。分别拟合建模数据和检验数据，结果均显示混合效应模型略优于基础模型（见表5）。
表 5 混合效应模型和基本线性模型拟合结果统计检验
Table 5. Goodness-of-fit statistics for mixed model and general linear model
项目 Item 建模数据 Modeling data 检验数据 Validation data Ra2 RMSE Ra2 RMSE 基础线性模型 General linear model 0.580 0.043 0.583 0.043 混合效应模型 Mixed effect model 0.599 0.042 0.598 0.042 -
用前述得到的混合效应模型预测各树种林木空间利用率，结果如图3所示。
图 3 整体预测（固定参数）、分组预测（树种组）与森林空间利用率观察值对比图
Figure 3. Population predictions (fixed), within-group predictions (species group), and the observed forest space utilization versus DBH
由图3可知，除慢阔组外，各树种空间利用率相近，其中：云杉的空间利用率明显高于树种平均水平，并且维持时间长，在胸径40 cm左右时达到峰值；冷杉、落叶松、红松、中阔组与平均水平接近，约在37 cm左右时达到峰值；慢阔组生长显著低于平均水平，约在32 cm时空间利用率最高，之后各自开始下降。
-
在天然混交林中，不同树种占据的营养空间不同，故直接比较不同树种的生长速率是不合理的。林木空间利用率以单位营养空间的生长量作为指标，使得不同树种间的比较更为合理。
经过比较，材积生长量与树冠投影面积之比与胸径相关系数最高，且材积生长是林业经营的重要目标之一，所以以此定义为本文的林木空间利用率。通过逐步回归，最终选定林木胸径，胸径平方（代表胸高断面积），树高，冠幅，样地蓄积，针阔比，坡向坡度，竞争指标作为自变量建立林木空间利用率的线性模型，并以此为基础模型，加入按树种分组的随机变量及误差方差结构，提高模型在不同树种间的适应性。
利用本文建立的林木空间利用率混合模型拟合研究数据，可以得到不同树种空间利用率的生长过程。除慢阔组外，各树种空间利用率相近，其中，云杉的空间利用率明显高于树种平均水平，并且维持时间长，在胸径40 cm左右时达到峰值；冷杉、落叶松、红松、中阔组与平均水平接近，约在37 cm左右时达到峰值；慢阔组生长显著低于平均水平，约在32 cm时空间利用率最高。在此之后，林木空间利用率下降，到达成熟收获点。
与空间结构指标相比，空间利用率指标更注重林地的生产功能，前者通过调整林分空间结构改善林内竞争条件，促进林木生长，后者是根据林木生长特性为不同林木配置适宜的生长空间，实现林地的最大化利用，如果将二者相结合则可以更好的对林分进行优化。
本文所定义的林木成熟仅是根据林木自身生长规律及立地条件所作的判断，在实际生产中，采伐计划的制订还要综合考虑采伐后林下幼苗幼树的生长及保持林分的生态稳定性。因此，在下一步研究中会考虑加入林分空间结构指数[5]等其他林分生长的相关因素，以长周期全林分林木空间利用率之和最大为经营目标制定采伐计划，在保证林分生态稳定性的同时提高林地的生产力。
本文所用数据集中在胸径为10 ~ 40 cm的林木，有一定局限性，计算所得数量成熟林木均未超过40 cm可能也是受此影响，在今后的研究中会尝试将该指标应用于更多大径级木，以期得到更合理的结果。另外，因为本文所采用的数据中，优势树种云杉冷杉占绝对优势，其他树种比例很小，所以以树种组分组的混合模型与一般线性模型的残差区别很小，混合效应优势不明显。但是计算结果显示各树种组空间利用率的生长过程接近，也从侧面验证了林木的生长量与冠幅的关系紧密。
Mixed model of forest space utilization in spruce-fir coniferous and broadleaved mixed forest of Changbai Mountains, northeastern China
-
摘要:
目的建立林木空间利用率模型，为天然混交林中不同树种间生产力的比较提供依据，为单木成熟的判断提供参考。 方法利用林木生长量与树冠大小比值定义林木空间利用率，以长白山地区云冷杉针阔混交林为研究对象，基于20块标准地的2 268株单木数据，建立林木空间利用率混合模型，拟合各树种的空间利用率。 结果（1）备选指标中蓄积生长量和树冠投影面积之比与胸径相关系数最高，适宜作为计算指标。（2）通过逐步回归，最终选定林木胸径、胸径平方（代表胸高断面积）、树高、冠幅、样地蓄积、针阔比、坡向坡度、竞争指标作为林木空间利用率基础模型的自变量。（3）确定按树种分组，包含胸径平方及截距随机效应参数、指数函数异方差结构的混合模型，经检验，混合模型在建模数据及检验数据中的表现均略优于一般线性模型。（4）利用所构建的混合效应模型，对研究数据进行拟合预测，各树种空间利用率最高时期的胸径分别为云杉约40 cm，冷杉、落叶松、红松、中阔组约37 cm，慢阔组约32 cm。 结论林木利用率模型得到的数量成熟是以单位营养空间的生产力为基础的，使不同树种间的比较更为合理，且计算结果符合一般林学规律，可以作为该地区云冷杉针阔混交林判定单木成熟及优化林分结构的参考依据。 Abstract:ObjectiveIn this paper, the forest space utilization model was established to provide the basis for the comparison of productivity among different tree species in the natural mixed forest, and provide reference for the judgement of tree maturity. MethodThe combination of tree growth and canopy size was used to define the forest space utilization. Based on the 2 268 individuals of 20 standard sites, which located in spruce-fir and broadleaved mixed forests in Changbai Mountain area of northeastern China, a mixed-effects model of forest space utilization was established to fit the spatial utilization of different species. Result(1) The correlation coefficient between DBH and the ratio of volume growth to the crown projection area was the highest, so it is suitable to be the calculation index. (2) Through the stepwise regression, the tree DBH, DBH square (represent area), tree height, crown width, plot accumulation, proportion of coniferous trees, slope gradient and competition index were selected as the independent variables of forest space utilization. (3) Mixed effects model was developed, including the DBH square and intercept random effect parameters by species groups and the exponential function heteroscedasticity structure. The results showed that the mixed effects model was better than the general linear model both in the modeling data and the test data. (4) Using the mixed effect model to predict the research area data, the DBH of the highest spatial utilization of each tree species was obtained: spruce (Picea koraiensis) was about 40 cm; fir (Abies holophylla), Korean pine (Pinus koraiensis), larch (Larix olgensis), and medium-speed broadleaf trees were about 37 cm; slow broadleaf trees were about 32 cm. ConclusionThe mature DBH obtained from the forest utilization model is based on the productivity of the unit nutrient space, so that the comparison between different tree species is more reasonable. The calculation results are in accordance with the general forest rule, so it can be used as a reference to judge the tree maturity, optimize the stand structure in the coniferous and broadleaved mixed forest in that area. -
图 1 林木空间利用率与胸径的关系
Bg、Bv、Vg、Vv是待定的林木空间利用率指标。Bg, Bv, Vg and Vv are the undetermined indexes of forest spatial utilization.
Figure 1. Relationship between forest space utilization and DBH
图 2 各树种组线性模型拟合参数的95%置信区间
Figure 2. 95% confidence intervals on the linear model parameters for each tree species group
图 3 整体预测（固定参数）、分组预测（树种组）与森林空间利用率观察值对比图
Vt是正态化后的林木空间利用率。Vt is normalization of forest space utilization.
Figure 3. Population predictions (fixed), within-group predictions (species group), and the observed forest space utilization versus DBH
表 1 建模数据和检验数据主要指标统计信息
Table 1. Summary statistics for modeling and validation data sets
数据
Data树种组
Tree species group株树
Plant number统计指标
Statistic index胸径
DBH/cm树高
Tree height/m平均冠幅
Average crown width/m样地蓄积/（m3·hm− 2）
Sample plot volume/ (m3·ha− 1)针叶树比例
Conifer percentage建模数据 Modeling data 红松
Pinus koraiensis219 平均值 Mean 24.4 18.4 2.2 296.5 0.91 标准差 SD 8.4 3.9 0.7 73.7 0.06 云杉
Picea koraiensis271 平均值 Mean 22.1 17.9 2.1 288.3 0.88 标准差 SD 9.7 5.2 0.7 74.1 0.07 冷杉
Abies holophylla717 平均值 Mean 21.5 18.7 1.8 297.6 0.90 标准差 SD 7.9 4.4 0.6 76.9 0.06 落叶松
Larix olgensis80 平均值 Mean 22.0 20.7 2.2 263.6 0.90 标准差 SD 7.4 4.5 0.7 62.3 0.07 中阔
Medium-speed broadleaf trees57 平均值 Mean 19.5 19.7 2.0 333.0 0.87 标准差 SD 6.6 4.5 0.9 53.2 0.05 慢阔
Slow broadleaf trees168 平均值 Mean 17.2 16.8 2.2 298.3 0.87 标准差 SD 6.6 4.7 0.8 71.4 0.06 检验数据
Validation data红松
Pinus koraiensis122 平均值 Mean 24.2 18.9 2.1 294.5 0.91 标准差 SD 7.8 4.0 0.7 71.7 0.06 云杉
Picea koraiensis141 平均值 Mean 21.6 17.9 2.1 290.8 0.88 标准差 SD 8.6 5.0 0.7 76.0 0.07 冷杉
Abies holophylla357 平均值 Mean 22.8 19.3 1.9 305.4 0.90 标准差 SD 7.8 4.4 0.6 72.4 0.06 落叶松
Larix olgensis26 平均值 Mean 23.1 21.4 2.5 246.2 0.88 标准差 SD 6.7 3.9 0.8 51.4 0.08 中阔
Medium-speed broadleaf trees21 平均值 Mean 21.9 21.0 2.3 329.7 0.88 标准差 SD 6.0 2.9 0.9 58.4 0.06 慢阔
Slow broadleaf trees89 平均值 Mean 17.0 17.0 2.2 297.7 0.86 标准差 SD 6.5 4.6 0.8 78.7 0.06 注：SD为标准差。Note: SD means standard deviation. 
下载: 导出CSV
表 2 林木空间利用率基础模型模拟结果
Table 2. Simulating results of basic model for forest space utilization
变量
Variable估计值
Estimating value方差
Variancet P (Pr > |t|) 截距 Intercept (INT) − 2.767 04 0.019 55 − 141.542 < 0.000 1*** 胸径 DBH 0.013 76 0.000 76 18.184 < 0.000 1*** 胸径平方 Square of the DBH (DBH2) − 0.000 18 0.000 01 − 13.335 < 0.000 1*** 树高 Tree height (H) 0.001 26 0.000 35 3.553 0.000 39*** 冠幅 Crown width (CW) − 0.051 17 0.001 71 − 29.975 < 0.000 1*** 蓄积 Stock volume (S) − 0.000 14 0.000 01 − 9.249 < 0.000 1*** 针阔比 Proportion of conifer to broadleaved trees (PF) 0.065 63 0.016 59 3.955 < 0.000 1*** 坡度坡向 Slope degree and slope aspect (SLC) 0.001 36 0.000 42 3.245 0.001 19 ** 竞争指数 Competition index (CIS) − 0.005 58 0.000 94 − 5.928 < 0.000 1*** 注：***表示在P < 0.001水平下差异显著，**表示在P < 0.01水平下差异显著。Notes: *** represents significant difference at P < 0.001 level, ** represents significant difference at P < 0.01 level.
下载: 导出CSV
表 3 不同随机效应拟合结果统计检验
Table 3. Goodness-of-fit statistics for mixed model of different random effects
模型编号
Model No.随机参数
Random parameterAIC BIC loglikelihood VTLME1 INT − 5 142.105 − 5 083.638 2 582.053 VTLME2 INT DBH − 5 161.782 − 5 092.684 2 593.891 VTLME3 INT H − 5 145.099 − 5 076.001 2 585.550 VTLME4 INT SLC − 5 150.207 − 5 081.109 2 588.104 VTLME5 INT DBH2 − 5 165.259 − 5 096.161 2 595.630 VTLME6 INT DBH2 SLC − 5 172.721 − 5 087.678 2 602.360 注：VTLME1、VTLME2、VTLME3、VTLME4、VTLME5、VTLME6分别为加入相应随机参数的混合模型编号。Notes: VTLME1，VTLME2，VTLME3，VTLME4，VTLME5，VTLME6 are the mixed model codes with corresponding random parameters.
下载: 导出CSV
表 4 不同残差方差结构拟合结果统计检验
Table 4. Goodness-of-fit statistics for mixed model of different residual variance structure
模型编号
Model No.残差方差模型
Residual variance modelAIC BIC loglikelihood VTLME5 无 None − 5 165.259 − 5 096.161 2 595.630 VTLME51 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma {\mathrm{D}\mathrm{B}\mathrm{H}}^{2}\right)$ − 5 167.432 − 5 093.019 2 597.716 VTLME52 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma \mathrm{C}\mathrm{W}\right)$ − 5 175.094 − 5 100.681 2 601.547 VTLME53 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma \mathrm{C}\mathrm{I}\mathrm{S}\right)$ − 5 174.558 − 5 100.145 2 601.279 VTLME54 $\scriptstyle {\sigma }^{2}\mathrm{e}\mathrm{x}\mathrm{p}\left(2\gamma f\right({x}_{i}\left)\right)$ − 5 185.268 − 5 110.855 2 606.634 VTLME55 $\scriptstyle {\sigma }^{2}{\mathrm{D}\mathrm{B}\mathrm{H} }^{2 \times2\gamma } $ − 5 170.525 − 5 096.112 2 599.263 VTLME56 $\scriptstyle {\sigma }^{2}{\mathrm{C}\mathrm{W}}^{2\gamma } $ − 5 175.849 − 5 101.436 2 601.925 VTLME57 $\scriptstyle {\sigma }^{2}{\mathrm{C}\mathrm{I}\mathrm{S}}^{2\gamma } $ 不收敛 Not converge VTLME58 $\scriptstyle {\sigma }^{2}{f\left({x}_{i}\right)}^{2\gamma } $ − 5 185.182 − 5 110.769 2 606.591 注：VTLME5、VTLME51、VTLME52、VTLME53、VTLME54、VTLME55、VTLME56、VTLME57、VTLME58分别为加入相应残差方差模型的混合模型编号；$\scriptstyle f\left({x}_{i}\right)$为固定效应及随机效应的预测值。Notes: VTLME5,VTLME51, VTLME52, VTLME53, VTLME54, VTLME55, VTLME56, VTLME57, VTLME58 are the mixed model codes with corresponding residual variance models; $\scriptstyle f\left({x}_{i}\right)$ represents the predicted value of the fixed effect and the random effect.
下载: 导出CSV
表 5 混合效应模型和基本线性模型拟合结果统计检验
Table 5. Goodness-of-fit statistics for mixed model and general linear model
项目 Item 建模数据 Modeling data 检验数据 Validation data Ra2 RMSE Ra2 RMSE 基础线性模型 General linear model 0.580 0.043 0.583 0.043 混合效应模型 Mixed effect model 0.599 0.042 0.598 0.042
下载: 导出CSV
-
[1] 亢新刚. 森林经理学[M]. 4版. 北京: 中国林业出版社, 2011: 152. Kang X G. Forest management[M]. 4th ed. Beijing: China Forestry Publishing House, 2011: 152. [2] 方旭东. 异龄林林木成熟确定的方法[J]. 林业科学, 2007, 43(6):83−87. doi: 10.3321/j.issn:1001-7488.2007.06.015 Fang X D. Research of the confirmation method of mature level in uneven-aged forest[J]. Scientia Silvae Sinicae, 2007, 43(6): 83−87. doi: 10.3321/j.issn:1001-7488.2007.06.015 [3] 田猛. 天然异龄林主伐关键技术研究[D]. 北京: 北京林业大学, 2015. Tian M. The key technologies of the natural uneven-aged forest cutting[D]. Beijing: Beijing Forestry University, 2015. [4] Van Galen L G, Jordan G J, Musk R A, et al. Quantifying floristic and structural forest maturity: an attribute-based method for wet eucalypt forests[J]. Journal of Applied Ecology, 2018, 55(4): 1668−1681. doi: 10.1111/1365-2664.13133 [5] 陈亚南. 落叶松云冷杉林空间结构多样性及择伐模型的研究[D]. 北京: 北京林业大学, 2016. Chen Y N. Study on spatial structure diversity and selective cutting decision model of larch-sprue-fir stands[D]. Beijing: Beijing Forestry University, 2016. [6] 惠刚盈, 赵中华, 胡艳波, 等. 我国西北主要天然林经营模式设计[J]. 林业科学研究, 2016, 29(2):155−161. doi: 10.3969/j.issn.1001-1498.2016.02.001 Hui G Y, Zhao Z H, Hu Y B, et al. Management models for natural forests in northwestern China[J]. Forest Research, 2016, 29(2): 155−161. doi: 10.3969/j.issn.1001-1498.2016.02.001 [7] Gradel A, Ammer C, Ganbaatar B, et al. On the effect of thinning on tree growth and stand structure of white birch (Betula platyphylla Sukaczev) and Siberian larch (Larix sibirica Ledeb.) in Mongolia[J]. Forests, 2017, 8(4): 105. doi: 10.3390/f8040105 [8] 符利勇, 孙华. 基于混合效应模型的杉木单木冠幅预测模型[J]. 林业科学, 2013, 49(8):65−74. doi: 10.11707/j.1001-7488.20130810 Fu L Y, Sun H. Individual crown diameter prediction for Cunninghamia lanceolata forests based on mixed effects models[J]. Scientia Silvae Sinicae, 2013, 49(8): 65−74. doi: 10.11707/j.1001-7488.20130810 [9] Lee S, Lee J, Kim S, et al. Forest structure and carbon dynamics of an intact lowland mixed dipterocarp forest in Brunei Darussalam[J]. Journal of Forestry Research, 2018, 29(1): 199−203. doi: 10.1007/s11676-017-0410-y [10] Arney J D. Tables for quantifying competitive stress on individual trees[R]. Victoria: Pacific Forestry Centre, 1973. [11] 吕飞舟. 基于CSI的蒙古栎林木竞争与分级研究[D]. 长沙: 中南林业科技大学, 2016. Lü F Z. Study on the forest competition and classification of Quercus mongolica based on CSI[D]. Changsha: Central South University of Forestry & Technology, 2016. [12] 张晔珵, 张怀清, 陈永富, 等. 基于树冠因子的林木竞争指数研究[J]. 林业科学研究, 2016, 29(1):80−84. doi: 10.3969/j.issn.1001-1498.2016.01.012 Zhang Y C, Zhang H Q, Chen Y F, et al. Study of tree competition index based on crown feature[J]. Forest Research, 2016, 29(1): 80−84. doi: 10.3969/j.issn.1001-1498.2016.01.012 [13] 刘强, 李凤日, 董利虎. 基于树冠竞争因子的落叶松人工林单木生长模型[J]. 植物研究, 2014, 34(4):547−553. doi: 10.7525/j.issn.1673-5102.2014.04.020 Liu Q, Li F R, Dong L H. Individual tree growth model of larch plantation based on crown competition factors[J]. Bulletin of Botanical Research, 2014, 34(4): 547−553. doi: 10.7525/j.issn.1673-5102.2014.04.020 [14] 臧颢, 雷相东, 张会儒, 等. 红松树高−胸径的非线性混合效应模型研究[J]. 北京林业大学学报, 2016, 38(6):8−16. Zang H, Lei X D, Zhang H R, et al. Nonlinear mixed-effects height-diameter model of Pinus koraiensis[J]. Journal of Beijing Forestry University, 2016, 38(6): 8−16. [15] 孔雷. 金沟岭林场三种林型最优林分结构的研究[D]. 北京: 北京林业大学, 2013. Kong L. Optimal structures for three ecological forests in Jingouling Forest Farm[D]. Beijing: Beijing Forestry University, 2013. [16] 雷相东, 李永慈, 向玮. 基于混合模型的单木断面积生长模型[J]. 林业科学, 2009, 45(1):74−80. doi: 10.3321/j.issn:1001-7488.2009.01.014 Lei X D, Li Y C, Xiang W. Individual basal area growth model using multi-level linear mixed model with repeated measures[J]. Scientia Silvae Sinicae, 2009, 45(1): 74−80. doi: 10.3321/j.issn:1001-7488.2009.01.014 [17] 辛营营, 韦新良, 王鑫. 利用均值类非完备型竞争指标分析林木竞争[J]. 西南林业大学学报, 2011, 31(4):29−34. doi: 10.3969/j.issn.2095-1914.2011.04.006 Xin Y Y, Wei X L, Wang X. Application of incomplete competition mean index to analysis of forest plant competition[J]. Journal of Southwest Forestry University, 2011, 31(4): 29−34. doi: 10.3969/j.issn.2095-1914.2011.04.006 [18] 赵俊卉. 长白山云冷杉混交林生长模型的研究[D]. 北京: 北京林业大学, 2010. Zhao J H. Growth modeling for spruce-fir forest in Changbai Mountains[D]. Beijing: Beijing Forestry University, 2010. [19] 李春明. 基于两层次线性混合效应模型的杉木林单木胸径生长量模型[J]. 林业科学, 2012, 48(3):66−73. doi: 10.11707/j.1001-7488.20120311 Li C M. Individual tree diameter increment model for Chinese fir plantation based on two-level linear mixed effects models[J]. Scientia Silvae Sinicae, 2012, 48(3): 66−73. doi: 10.11707/j.1001-7488.20120311 [20] 孔雷, 亢新刚, 赵浩彦, 等. 长白山云冷杉针阔混交林最优竞争指标的研究[J]. 西北农林科技大学学报(自然科学版), 2014, 42(4):88−94, 101. Kong L, Kang X G, Zhao H Y, et al. Optimal competition indexes of spruce-fir conifer and broadleaf mixed stands in the Changbai Mountains[J]. Journal of Northwest A & F University (Natural Science Edition), 2014, 42(4): 88−94, 101. [21] Osborne J W. Improving your data transformations: applying the Box-Cox transformation[J]. Practical Assessment Research & Evaluation, 2010, 15(12): 1−9. [22] 姜立春, 张锐, 李凤日. 基于线性混合模型的落叶松枝条长度和角度模型[J]. 林业科学, 2012, 48(5):53−60. doi: 10.11707/j.1001-7488.20120508 Jiang L C, Zhang R, Li F R. Modeling branch length and branch angle with linear mixed effects for Dahurian larch[J]. Scientia Silvae Sinicae, 2012, 48(5): 53−60. doi: 10.11707/j.1001-7488.20120508 [23] Burnham K P, Anderson D R. Multimodel inference: understanding AIC and BIC in model selection[J]. Sociological Methods & Research, 2004, 33: 261−304. [24] Pinheiro J C, Bates D M. Mixed-effects models in S and S-plus[M]. New York: Springer-Verlag, 2000: 143. -
点击查看大图
图(3) / 表 (5)
计量
- 文章访问数: 1144
- HTML全文浏览量: 168
- PDF下载量: 36
- 被引次数: 0
