高级检索

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

基于混合效应的长白落叶松一级枝条基径预估模型

靳晓娟 孙玉军 潘磊

靳晓娟, 孙玉军, 潘磊. 基于混合效应的长白落叶松一级枝条基径预估模型[J]. 北京林业大学学报. doi: 10.12171/j.1000-1522.20200133
引用本文: 靳晓娟, 孙玉军, 潘磊. 基于混合效应的长白落叶松一级枝条基径预估模型[J]. 北京林业大学学报. doi: 10.12171/j.1000-1522.20200133
Jin Xiaojuan, Sun Yujun, Pan Lei. Prediction model of base diameter of primary branch for Larix olgensis based on mixed effects[J]. Journal of Beijing Forestry University. doi: 10.12171/j.1000-1522.20200133
Citation: Jin Xiaojuan, Sun Yujun, Pan Lei. Prediction model of base diameter of primary branch for Larix olgensis based on mixed effects[J]. Journal of Beijing Forestry University. doi: 10.12171/j.1000-1522.20200133

基于混合效应的长白落叶松一级枝条基径预估模型

doi: 10.12171/j.1000-1522.20200133
基金项目: 国家自然科学基金项目“基于树木生长过程的长白落叶松树冠模型”(31870620)
详细信息
    作者简介:

    靳晓娟。主要研究方向:森林资源监测与模型模拟。Email:1072733074@qq.com 地址:100083 北京市海淀区清华东路 35 号北京林业大学林学院

    通讯作者:

    孙玉军,教授,博士生导师。主要研究方向:森林资源监测与模型。Email:sunyj@bjfu.edu.cn 地址:同上

  • 中图分类号: S757

Prediction model of base diameter of primary branch for Larix olgensis based on mixed effects

  • 摘要:   目的  利用非线性混合效应建模方法构建龄组−单木两水平长白落叶松一级枝条基径模型,为探索不同龄组下枝条基径的生长特点及差异提供理论依据。  方法  对4个基础模型进行改进,通过调整系数($R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$),均方根误差(RMSE)选择备选模型,在此基础上构建长白落叶松枝条基径非线性混合效应模型。利用独立数据验证模型拟合结果,用平均绝对误差(MAE)、平均相对误差绝对值(MRAE)评价模型预测能力。并对基础模型与混合模型的预测值进行比较,利用龄组水平的随机参数模拟各龄组枝条基径的分布。  结果  以改进后的Gompertz方程为基础模型,当龄组随机效应作用于参数b、单木随机效应同时作用于参数bcd上,随机效应的方差协方差结构为广义正定矩阵,异方差结构为幂函数时,模型的拟合效果最优。混合模型的$ R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $有所提升,RMSE、MAE和MRAE都明显降低。最终模型的$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$、RMSE、MAE和MRAE分别为0.699 8、4.768 4 mm、3.705 8 mm和0.391 6 mm。混合模型的预测值可体现单木间的差异,枝条基径的分布范围随着龄组的增大逐渐增大,各龄组的枝径生长均有差异。  结论  考虑龄组与单木水平所构建的枝条基径混合效应模型能提高模型的预测精度。利用龄组水平的随机效应参数模拟枝条基径的生长可以反映其规律和差异,也符合树木生理学意义。因此基于龄组与单木两水平所构建的混合效应模型可对不同年龄长白落叶松一级枝枝条基径的生长进行合理的预测。
  • 图  1  枝条基径随着枝深度、胸径的分布

    Figure  1.  Distribution of branch base diameter with depth into crown and DBH

    图  2  基础模型(a)与混合模型包括广义正定矩阵(b)、幂函数(c)、指数函数(d)残差分布

    Figure  2.  Residual distribution of basic model (a) and mixed model with generalized positive definite matrix (b), mixed model with power function (c), and mixed model with exponential function

    图  3  枝条基径基础与混合模型预测值分布

    1 ~ 3为幼龄,4 ~ 12为中龄,13 ~ 16为近熟,17 ~ 18为成熟。1−3 are young age, 4−12 are middle age, 13−16 are near mature, 17−18 are mature.

    Figure  3.  Distribution of predicted values of branch base diameter by basic and mixed models

    图  4  不同龄组枝条基径分布

    Figure  4.  Distribution of branch base diameter of different age groups

    表  1  长白落叶松样地、解析木及枝条因子统计

    Table  1.   Statistics of sample plots, analytic trees and branch attributes for Larix olgensis

    因子 Factor建模数据 Fitting data检验数据 Validating data
    最大值 Max.最小值 Min.均值 Mean最大值 Max.最小值 Min.均值 Mean
    样地、样木数 Number of sample plot and sample tree18(幼龄3、中龄9、近熟4、成熟2)
    18 (young age 3, middle age 9,
    near mature 4, mature 2)
    6(幼龄2、中龄2、近熟1、成熟1)
    6 (young age 2, middle age 2,
    near mature 1, mature 1)
    坡度 Slope/(°)15081004
    海拔 Altitude/m3760304456266380
    立地指数 Site index231720221720
    林分平均胸径 Average DBH of stand/cm29.04.016.521.28.115.8
    林分密度/(株·hm−2) Stand density/(plant·ha−1)2 8004338582 6175661 105
    枝条数 Branch number2 089700
    胸径 DBH/cm28.24.516.820.98.215.1
    年龄/a Age/year55827441124
    冠长 Crown length/m13.54.08.410.35.36.6
    着枝深度 Depth into crown/m11.500.023.6010.300.103.40
    枝条基径 Branch base diameter/mm52.01.414.445.02.012.4
    枝条长度 Branch length/m4.00.11.414.20.11.2
    下载: 导出CSV

    表  2  长白落叶松枝条基径基础及改进模型

    Table  2.   Basic and improved models of branch base diameter for Larix olgensis

    模型名称
    Model name
    模型编号
    Model No.
    基础模型形式
    Basic model form
    模型编号
    Model No.
    改进模型形式
    Improved model form
    Gompertz1${\mathrm{BD}}=a{\rm{e}}^{(-b{\rm{e}}^{-c{\mathrm{DINC}}})} $5${\mathrm{BD}}=a{\rm{e}}^{(-b{\rm{e}}^{-c{\mathrm{DINC}}})}{\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)} $
    Mitscherlich2${\mathrm{BD}}=a(1-b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $6${\mathrm{BD}}=a(1-b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}){\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)} $
    Logistic3${\mathrm{BD}}=a/(1+b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $7${\mathrm{BD}}=a{\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)}/(1+b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $
    单分子式 Monomolecular formula4${\mathrm{BD}}=a(1-{\rm{e}}^{\left(-b{\mathrm{DINC}}\right)}) $8$ {\mathrm{BD}}=a(1-{\rm{e}}^{\left(-b{\mathrm{DINC}}\right)}){\rm{e}}^{\left({{\mathrm{DBH}}}^{c}\right)} $
    注:BD和DINC分别是一级枝条基径和着枝深度;DBH是单木胸径;abcd是待拟合的参数。Notes: BD and DINC are base diameter and depth into crown of primary branch, respectively; DBH is individual tree DBH; a, b, c and d are parameters to be fitted.
    下载: 导出CSV

    表  3  基础及改进模型拟合与检验结果

    Table  3.   Fitting and validating results of basic and improved models

    项目
    Item
    模型编号
    Model No.
    拟合优度 Goodness of fitting检验指标 Validating index
    均方根误差 RMSE$ R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $平均相对误差 MAE平均相对误差绝对值 MRAE
    基础
    Basic model
    16.080 70.512 04.660 60.533 8
    26.095 50.509 74.675 30.536 1
    36.076 10.512 94.663 50.536 0
    46.155 40.500 24.687 30.521 2
    改进
    Improved model
    55.518 80.597 93.727 40.402 4
    65.534 90.595 53.746 60.405 9
    75.520 80.597 63.732 20.405 4
    85.552 70.593 13.745 10.397 7
    下载: 导出CSV

    表  4  不同随机参数最优混合模型拟合结果比较

    Table  4.   Comparison in fitting results of optimal mixed models with different random parameters

    随机参数个数
    Number of random
    parameter
    模型编号
    Model No.
    随机参数 Random parameter统计指标 Statistical index似然比检验 Likelihood ratio (LRT) test
    龄组水平
    Age group level
    单木水平
    Individual tree level
    AICBICLoglik似然比 LRTP
    25-1cb12 728.3512 767.86−6 357.175
    35-2bcd12 645.0512 695.85−6 313.52387.303 6 < 0.000 1
    45-3bbcd12 621.6712 689.40−6 298.83329.380 4 < 0.000 1
    55-4b、db、c、d12 625.6712 704.69−6 298.8340.001 30.999 4
    65-5a、c、db、c、d12 633.6412 729.59−6 299.8191.970 90.578 5
    75-6a、b、c、da、b、d12 665.6312 784.16−6 311.81323.988 70.000 1
    85-7a、b、c、da、b、c、d12 701.6612 842.77−6 325.831212.709 5 < 0.000 1
    下载: 导出CSV

    表  5  不同方差协方差结构拟合结果比较

    Table  5.   Comparison in fitting results of different variance-covariance structures

    方差−协方差结构
    Variance-covariance structure
    AICBICLoglik似然比
    LRT
    P
    广义正定矩阵 Generalized positive definite matrix12 621.6712 689.4−6 298.83
    对角矩阵 Diagonal matrix12 625.0112 687.1−6 301.505.341 00.020 8
    下载: 导出CSV

    表  6  不同异方差结构混合模型拟合结果比较

    Table  6.   Comparison in fitting results of mixed model with different heteroscedasticity structures

    异方差结构
    Heteroscedasticity structure
    AICBICLoglik似然比
    LRT
    P
    无 None12 621.6712 689.40 −6 298.830
    幂函数 Power function12 151.5312 230.55−6 061.764474.14 < 0.000 1
    指数函数 Exponential function12 220.7212 294.10−6 097.35971.19 < 0.000 1
    下载: 导出CSV

    表  7  基础与混合模型拟合结果及检验

    Table  7.   Fitting and validating results of basic and mixed models

    模型形式 Model form参数 Parameter (标准差 SD)拟合优度 Goodness of fitting评价指标 Evaluating index
    abcdRMSE$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $MAEMRAE
    基础模型(5)
    Basic model (5)
    3.172 82.017 40.384 90.259 35.518 80.597 93.727 40.402 4
    (0.190 4)(0.069 4)(0.007 6)(0.023 4)
    混合模型(5-3)
    Mixed model (5-3)
    2.088 11.882 30.379 60.328 74.768 40.699 83.705 80.391 6
    (0.194 6)(0.164 5)(0.040 4)(0.014 3)
    异方差(幂函数)
    Heteroscedasticity (power function)
    $ \alpha $ = 0.63
    模型方差
    Model variance
    $ {\sigma }^{2} $ = 6.597 8
    方差协方差结构
    Variance-covariance structure
    ${{D} }_1=\left|3.82 \times 10^{-34}\right|$$ {{D}}_2=\left|\begin{array}{ccc}0.450 \; 9& 0.257 \; 0& -0.143 \; 0\\ 0.257 \; 0& 0.021 \; 4& -0.752 \; 0\\ -0.143 \; 0& -0.752 \; 0& 0.000 \; 4\end{array}\right| $
    注:D1D2分别是龄组、单木水平的方差协方差结构。Notes: D1 and D2 are variance-covariance structures of age group and individual tree level, respectively.
    下载: 导出CSV

    表  8  固定效应与随机效应检验结果的比较

    Table  8.   Comparison of validating results of fixed effect and random effect

    单木编号
    Individual tree No.
    固定效应 Fixed effect随机效应 Random effect
    平均绝对误差 MAE平均相对误差绝对值 MRAE平均绝对误差 MAE平均相对误差绝对值 MRAE
    12.691 10.457 02.309 90.364 1
    22.772 30.296 42.477 10.229 6
    35.490 90.264 64.368 70.227 2
    43.579 90.274 13.420 80.267 7
    53.913 00.469 63.753 10.340 8
    64.746 10.687 52.902 90.359 0
    下载: 导出CSV
  • [1] Brown P L, Doley D, Keenan R J. Stem and crown dimensions as predictors of thinning responses in a crowded tropical rainforest plantation of Flindersia brayleyana F. Muell[J]. Forest Ecology and Management, 2004, 196(2): 379−392.
    [2] Maguire D A, Kershaw J A, Hann D W. Predicting the effects of silvicultural regime on branch size and crown wood core in Douglas-fir[J]. Forest Science, 1991, 37(5): 1409−1428.
    [3] 王小青, 刘杏娥, 任海青. 树冠特征对小黑杨木材性质和生长量的影响研究[J]. 林业科学研究, 2007, 20(6):801−806. doi:  10.3321/j.issn:1001-1498.2007.06.011.

    Wang X Q, Liu X E, Ren H Q. Effects of crown attributes on wood characteristics and increments of Populus × xiaohei[J]. Forest Research, 2007, 20(6): 801−806. doi:  10.3321/j.issn:1001-1498.2007.06.011.
    [4] Groot A, Schneider R. Predicting maximum branch diameter from crown dimensions, stand characteristics and tree species[J]. Forestry Chronicle, 2011, 87(4): 542−551.
    [5] 孙晓, 李春友, 贺红月, 等. ‘107杨’一级枝条基径和长度的动态模拟研究[J]. 林业与生态科学, 2018, 33(4):364−372.

    Sun X, Li C Y, He H Y, et al. Dynamic simulation of length and diameter of first branches of ‘107 poplar’[J]. Forestry and Ecological Sciences, 2018, 33(4): 364−372.
    [6] 李春明. 基于两层次线性混合效应模型的杉木林单木胸径生长量模型[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.
    [7] Yang Y, Huang S. Allometric modelling of crown width for white spruce by fixed- and mixed-effects models[J]. The Forestry Chronicle, 2017, 93(2): 138−147.
    [8] 朱万才, 贾炜玮. 基于随机效应的红松人工林一级枝条动态生长模型[J]. 森林工程, 2015, 31(4):26−32. doi:  10.3969/j.issn.1001-005X.2015.04.007.

    Zhu W C, Jia W W. Dynamic growth model for First-hierarchy branch of Korean pine plantation based on mixed effect model[J]. Forest Engineering, 2015, 31(4): 26−32. doi:  10.3969/j.issn.1001-005X.2015.04.007.
    [9] 姜立春, 李凤日, 张锐. 基于线性混合模型的落叶松枝条基径模型[J]. 林业科学研究, 2012, 25(4):464−469. doi:  10.3969/j.issn.1001-1498.2012.04.009.

    Jiang L C, Li F R, Zhang R. Modeling branch diameter with linear mixed effects for Dahurian larch[J]. Forest Research, 2012, 25(4): 464−469. doi:  10.3969/j.issn.1001-1498.2012.04.009.
    [10] Chun-Sheng W, Cheng T, Sebastian H, et al. Branch development of five-year-old Betula alnoides plantations in response to planting density[J/OL]. Forests, 2018, 9(1): 42 [2019−11−12]. https://doi.org/10.3390/f9010042.
    [11] Weiskittel A R, Maguire D A, Monserud R A. Response of branch growth and mortality to silvicultural treatments in coastal Douglas-fir plantations: implications for predicting tree growth[J]. Forest Ecology and Management, 2007, 251(3): 182−194.
    [12] 姜立春, 潘莹, 李耀翔. 兴安落叶松枝条特征联立方程组模型及树冠形状模拟[J]. 北京林业大学学报, 2016, 38(6):1−7.

    Jiang L C, Pan Y, Li Y X. Model systems of branch characteristics and crown profile simulation for Larix gmelinii[J]. Journal of Beijing Forestry University, 2016, 38(6): 1−7.
    [13] 陈东升, 孙晓梅, 李凤日. 落叶松人工林枝条直径和长度的非线性混合模型[J]. 南京林业大学学报(自然科学版), 2015, 39(6):74−80.

    Chen D S, Sun X M, Li F R. Nonlinear mixed models of branch diameter and length in larch plantation[J]. Journal of Nanjing Forestry University (Natural Sciences Edition), 2015, 39(6): 74−80.
    [14] Lemay A, Pamerleau-Couture É, Krause C. Maximum branch diameter in black spruce following partial cutting and clearcutting[J/OL]. Forests 2019, 10(10): 913 [2020−01−03]. https://www.mdpi.com/1999-4907/10/10/913.
    [15] Jia W, Chen D. Nonlinear mixed-effects height to crown base and crown length dynamic models using the branch mortality technique for a Korean larch (Larix olgensis) plantations in northeast China[J]. Journal of Forestry Research, 2019, 30(6): 2095−2109.
    [16] 贾炜玮, 孙守强, 李凤日, 等. 基于混合模型的红松人工林枝条动态研究[J]. 植物研究, 2015, 35(3):425−430. doi:  10.7525/j.issn.1673-5102.2015.03.016.

    Jia W W, Sun S Q, Li F R, et al. Branch dynamics for the Korean pine plantation based on linear mixed model[J]. Bulletin of Botanical Research, 2015, 35(3): 425−430. doi:  10.7525/j.issn.1673-5102.2015.03.016.
    [17] 王春红, 李凤日, 贾炜玮, 等. 基于非线性混合模型的红松人工林枝条生长[J]. 应用生态学报, 2013, 24(7):1945−1952.

    Wang C H, Li F R, Jia W W, et al. Branch growth of Korean pine plantation based on nonlinear mixed model[J]. Chinese Journal of Applied Ecology, 2013, 24(7): 1945−1952.
    [18] Sharma R P, Vacek Z, Vacek S. Individual tree crown width models for Norway spruce and European beech in Czech Republic[J]. Forest Ecology and Management, 2016, 366: 208−220.
    [19] 高慧淋, 董利虎, 李凤日. 基于混合效应的人工落叶松树冠轮廓模型[J]. 林业科学, 2017, 53(3):84−93. doi:  10.11707/j.1001-7488.20170310.

    Gao H L, Dong L H, Li F R. Crown shape model for Larix olgensis plantation based on mixed effect[J]. Forest Research, 2017, 53(3): 84−93. doi:  10.11707/j.1001-7488.20170310.
    [20] Fu L, Sun H, Sharma R P, et al. Nonlinear mixed-effects crown width models for individual trees of Chinese fir (Cunninghamia lanceolata) in south-central China[J]. Forest Ecology and Management, 2013, 302: 210−220.
    [21] 罗恒春, 张超, 魏安超, 等. 云南松林分平均高生长模型及模型参数环境解释[J]. 北京林业大学学报, 2018, 40(4):67−75.

    Luo H C, Zhang C, Wei A C, et al. Stand average height growth model and environmental interpretation in model parameter of Pinus yunnanensis[J]. Journal of Beijing Forestry University, 2018, 40(4): 67−75.
    [22] 洪奕丰, 陈东升, 申佳朋, 等. 长白落叶松人工林单木和林分水平的相容性生物量模型研究[J]. 林业科学研究, 2019, 32(4):33−40.

    Hong Y F, Chen D S, Shen J P, et al. Compatible biomass models forLarix olgensis plantation based on tree-level and stand-level[J]. Forest Research, 2019, 32(4): 33−40.
    [23] 罗恒春, 张超, 魏安超, 等. 云南松林分平均胸径生长模型及模型参数环境解释[J]. 浙江农林大学学报, 2018, 35(6):1079−1087. doi:  10.11833/j.issn.2095-0756.2018.06.011.

    Luo H C, Zhang C, Wei A C, et al. Average DBH growth model of a stand with environmental parameters for Pinus yunnanensis in central Yunnan, China[J]. Journal of Zhejiang A&F University, 2018, 35(6): 1079−1087. doi:  10.11833/j.issn.2095-0756.2018.06.011.
    [24] Huff S, Poudel K P, Ritchie M, et al. Quantifying aboveground biomass for common shrubs in northeastern California using nonlinear mixed effect models[J]. Forest Ecology and Management, 2018, 424: 154−163.
    [25] 姜立春, 张锐, 李凤日. 基于线性混合模型的落叶松枝条长度和角度模型[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.
  • [1] 牛亦龙, 董利虎, 李凤日.  基于广义代数差分法的长白落叶松人工林地位指数模型 . 北京林业大学学报, doi: 10.12171/j.1000-1522.20190036
    [2] 贾炜玮, 冯万举, 李凤日.  基于节子剖析数据的长白落叶松人工林枝条丢失年轮数研究 . 北京林业大学学报, doi: 10.12171/j.1000-1522.20190038
    [3] 徐奇刚, 雷相东, 国红, 李海奎, 李玉堂.  基于多层感知机的长白落叶松人工林林分生物量模型 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20190035
    [4] 王涛, 董利虎, 李凤日.  基于混合效应的杂种落叶松人工幼龄林单木枯损模型 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20170437
    [5] 王烁, 董利虎, 李凤日.  人工长白落叶松枝条存活模型 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20170203
    [6] 王冬至, 张冬燕, 张志东, 黄选瑞.  塞罕坝华北落叶松人工林断面积预测模型 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20170072
    [7] 邵英男, 田松岩, 刘延坤, 陈瑶, 孙志虎.  密度调控对长白落叶松人工林土壤呼吸的影响 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20170029
    [8] 王曼霖, 董利虎, 李凤日.  基于Possion回归混合效应模型的长白落叶松一级枝数量模拟 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20170204
    [9] 姜立春, 潘莹, 李耀翔.  兴安落叶松枝条特征联立方程组模型及树冠形状模拟 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20150339
    [10] 刘宇, 郭建斌, 王彦辉, 刘泽彬, 邓秀秀, 张桐, 熊伟, 左海军.  宁夏六盘山不同密度华北落叶松人工林枯落物水文效应 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20160007
    [11] 臧颢, 雷相东, 张会儒, 李春明, 卢军.  红松树高-胸径的非线性混合效应模型研究 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20160008
    [12] 欧光龙, 胥辉, 王俊峰, 肖义发, 陈科屹, 郑海妹.  思茅松天然林林分生物量混合效应模型构建 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20140316
    [13] 毛斌, 彭立群, 李乐, 徐程扬.  侧柏风景林美景度的林内色彩斑块非线性模型研究 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20140481
    [14] 王健铭, 张天汉, 于琳倩, 董芳宇, 李景文, 李俊清, 赵秀海.  基于生态保护格局的长白山北部森林景观特征分析 . 北京林业大学学报, doi: 10.13332/j.1000-1522.20140413
    [15] 张建亮, 崔国发, 黄祥童, 郭子良, 周海城.  国家一级保护植物长白松种群结构与动态预测 . 北京林业大学学报, doi: 10.13332/j.cnki.jbfu.2014.03.004
    [16] 李耀翔, 姜立春.  基于非线性混合模型的落叶松木材管胞长度模拟 . 北京林业大学学报,
    [17]
    孙志虎, 毕永娟, 牟长城, 蔡体久
    基于FORECAST模型的长白落叶松人工林经营措施对长期生产力的影响 . 北京林业大学学报,
    [18] 李春明, 李利学.  基于非线性混合模型的栓皮栎树高与胸径关系研究 . 北京林业大学学报,
    [19] 王雄宾, 余新晓, 徐成立, 谷建才, 周彬, 范敏锐, 贾国栋, 吕锡芝, .  间伐对华北落叶松人工林边缘效应的影响 . 北京林业大学学报,
    [20] 李春明.  利用非线性混合模型进行杉木林分断面积生长模拟研究 . 北京林业大学学报,
  • 加载中
图(4) / 表 (8)
计量
  • 文章访问数:  99
  • HTML全文浏览量:  39
  • PDF下载量:  21
  • 被引次数: 0
出版历程
  • 收稿日期:  2020-05-05
  • 修回日期:  2020-07-07
  • 网络出版日期:  2020-10-09

基于混合效应的长白落叶松一级枝条基径预估模型

doi: 10.12171/j.1000-1522.20200133
    基金项目:  国家自然科学基金项目“基于树木生长过程的长白落叶松树冠模型”(31870620)
    作者简介:

    靳晓娟。主要研究方向:森林资源监测与模型模拟。Email:1072733074@qq.com 地址:100083 北京市海淀区清华东路 35 号北京林业大学林学院

    通讯作者: 孙玉军,教授,博士生导师。主要研究方向:森林资源监测与模型。Email:sunyj@bjfu.edu.cn 地址:同上
  • 中图分类号: S757

摘要:   目的  利用非线性混合效应建模方法构建龄组−单木两水平长白落叶松一级枝条基径模型,为探索不同龄组下枝条基径的生长特点及差异提供理论依据。  方法  对4个基础模型进行改进,通过调整系数($R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$),均方根误差(RMSE)选择备选模型,在此基础上构建长白落叶松枝条基径非线性混合效应模型。利用独立数据验证模型拟合结果,用平均绝对误差(MAE)、平均相对误差绝对值(MRAE)评价模型预测能力。并对基础模型与混合模型的预测值进行比较,利用龄组水平的随机参数模拟各龄组枝条基径的分布。  结果  以改进后的Gompertz方程为基础模型,当龄组随机效应作用于参数b、单木随机效应同时作用于参数bcd上,随机效应的方差协方差结构为广义正定矩阵,异方差结构为幂函数时,模型的拟合效果最优。混合模型的$ R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $有所提升,RMSE、MAE和MRAE都明显降低。最终模型的$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$、RMSE、MAE和MRAE分别为0.699 8、4.768 4 mm、3.705 8 mm和0.391 6 mm。混合模型的预测值可体现单木间的差异,枝条基径的分布范围随着龄组的增大逐渐增大,各龄组的枝径生长均有差异。  结论  考虑龄组与单木水平所构建的枝条基径混合效应模型能提高模型的预测精度。利用龄组水平的随机效应参数模拟枝条基径的生长可以反映其规律和差异,也符合树木生理学意义。因此基于龄组与单木两水平所构建的混合效应模型可对不同年龄长白落叶松一级枝枝条基径的生长进行合理的预测。

English Abstract

靳晓娟, 孙玉军, 潘磊. 基于混合效应的长白落叶松一级枝条基径预估模型[J]. 北京林业大学学报. doi: 10.12171/j.1000-1522.20200133
引用本文: 靳晓娟, 孙玉军, 潘磊. 基于混合效应的长白落叶松一级枝条基径预估模型[J]. 北京林业大学学报. doi: 10.12171/j.1000-1522.20200133
Jin Xiaojuan, Sun Yujun, Pan Lei. Prediction model of base diameter of primary branch for Larix olgensis based on mixed effects[J]. Journal of Beijing Forestry University. doi: 10.12171/j.1000-1522.20200133
Citation: Jin Xiaojuan, Sun Yujun, Pan Lei. Prediction model of base diameter of primary branch for Larix olgensis based on mixed effects[J]. Journal of Beijing Forestry University. doi: 10.12171/j.1000-1522.20200133
  • 长白落叶松(Larix olgensis)是我国东北地区栽培面积较大的人工林种,生长迅速,轮伐期相对较短,是重要的用材树种之一。树木的枝条分布构成树冠的基本形态,一级枝作为树木骨架,是叶片的重要支撑结构,它的长短、大小及其空间排列直接影响了树木的光合、蒸腾等生理过程,不仅影响树木生长及生产力分配,还决定着树木的木材质量[1-3]。因此研究枝条的生长规律是研究树冠结构和功能以及树木有机物产量的基础,对林分生长动态具有重要的作用。

    多数学者以模型为基础,加入相关的林分变量来研究树木一级枝基径生长的规律,如Groot等利用木材易测因子(树高、胸径、林分断面积、冠幅半径)建立了加拿大主要用材树种的最大枝条直径线性模型,得出最大枝直径与冠幅半径、树高显著正相关,与林分断面积呈负相关[4]。孙晓等以对数线性模型为基础,利用胸径树高子模型嵌套,构建出自变量为树龄的‘107杨’(Populus × euramericanna ‘74/76’)一级枝枝径模型[5]。但是在实际调查中,枝条数据取自不同的林木和样地,形成了样地−单木−枝条的数据层次结构,如果统一拟合固定模型,模型会出现较大误差,混合效应模型可以解释这种情况,其在树木生长模型的研究中已广泛使用。如李春明考虑样地、区域两水平效应,建立了江西省杉木( Cunninghamia lanceolata)人工林单木胸径生长模型[6],Yang等使用固定效应模型和混合效应模型建立了白云杉(Picea glauca)冠幅异速生长模型,经过比较发现混合效应模型可以显著提高模型预测精度[7]。朱万才等通过非线性混合效应模型,分别考虑样地效应、单木效应和枝条效应,建立了红松(Pinus koraiensis)人工林一级枝条动态生长模型,它比基础模型更能描述单木差异[8]。姜立春等采用线性混合效应拟合长白落叶松枝条基径模型,结果表明林木胸径和着枝深度能很好地反映枝条基径的变化[9]

    目前已有很多基于样地、单木水平或者两水平嵌套的混合效应模型来描述长白落叶松枝条生长规律,但是不同年龄下树木的生长特点有较大差异,而结合龄组水平对枝条基径生长的研究较少,本文以黑龙江省东折棱河林场的长白落叶松人工林为研究对象,结合基础理论方程,加入变量,构建龄组和单木水平嵌套的非线性混合效应模型,并对各参数意义作出解释,模拟并分析不同龄组枝条基径的差异,为林场的经营管理决策提供依据。

    • 研究区位于黑龙江省伊春市朗乡林业局东折棱河林场,该林场处于小兴安岭南部低山丘陵区域,地理坐标为128°55′30″ ~ 129°15′21″E、46°31′58″ ~ 46°49′38″N,海拔200 ~ 970 m之间,林场气候特征属北温带大陆湿润性气候,年平均相对湿度65% ~ 75%,年平均气温1.9 ℃,降雨量620 mm。森林土壤类型为暗棕壤,土层厚度大约30 ~ 60 cm。林场总面积7 039 hm2,森林植被面积可达到95%,长白落叶松为该地区优势树种,其人工林占比37.8%。主要分布在海拔300 ~ 500 m之间,属于低山区,坡度较小,平均坡度7°。

    • 以东折棱河林场长白落叶松为研究对象,于2018年7月在该林场设置7块面积为20 m × 30 m的临时样地,测定样地坡向、坡度、海拔,实测胸径大于5 cm的全部乔木,记录树木胸径、树高、冠幅(东西、南北)、枝下高。依据$ {D}_{\mathrm{g}}=\sqrt{\dfrac{1}{N}\displaystyle\sum _{i=1}^{n}{{d}_{i}}^{2}} $计算林分平均胸径,选取3株胸径与$ {D}_{\mathrm{g}} $相近的“平均木”,求其树高均值作为林分平均高,据此选择一株平均木进行树干解析,测量每一段着生一级枝的枝长、枝径、枝角,着枝深度(一级枝着生部位到树顶的距离)。结合历史数据,共有来自24块样地的24株解析木(2 789组枝条),将其按照3∶1的比例分类,18株(2 089组枝条)作为建模样本,6株(700组枝条)用于检验(表1)。

      表 1  长白落叶松样地、解析木及枝条因子统计

      Table 1.  Statistics of sample plots, analytic trees and branch attributes for Larix olgensis

      因子 Factor建模数据 Fitting data检验数据 Validating data
      最大值 Max.最小值 Min.均值 Mean最大值 Max.最小值 Min.均值 Mean
      样地、样木数 Number of sample plot and sample tree18(幼龄3、中龄9、近熟4、成熟2)
      18 (young age 3, middle age 9,
      near mature 4, mature 2)
      6(幼龄2、中龄2、近熟1、成熟1)
      6 (young age 2, middle age 2,
      near mature 1, mature 1)
      坡度 Slope/(°)15081004
      海拔 Altitude/m3760304456266380
      立地指数 Site index231720221720
      林分平均胸径 Average DBH of stand/cm29.04.016.521.28.115.8
      林分密度/(株·hm−2) Stand density/(plant·ha−1)2 8004338582 6175661 105
      枝条数 Branch number2 089700
      胸径 DBH/cm28.24.516.820.98.215.1
      年龄/a Age/year55827441124
      冠长 Crown length/m13.54.08.410.35.36.6
      着枝深度 Depth into crown/m11.500.023.6010.300.103.40
      枝条基径 Branch base diameter/mm52.01.414.445.02.012.4
      枝条长度 Branch length/m4.00.11.414.20.11.2
    • 在构建枝条直径生长模型时,部分学者通过逐步回归法,选择相关性较高的因子,并且利用显著性指标、合理的生态学解释和方差膨胀因子(VIF)不大于3来选择模型的预测变量[10-12],或者基于林业上常用的生长方程,向模型中增加相关变量,如林分株数密度、胸径、树高、树冠宽度等[13-15],然后对比相应的预测指标选择合适的模型[16-17]。本研究基于各调查变量与枝条基径的相关性,最终选择4个基础模型:Gompertz、Mitscherlich、Logistic、单分子式模型进行拟合。并添加相关的林分变量,将基础模型再参数化,基础模型及改进式见表2

      表 2  长白落叶松枝条基径基础及改进模型

      Table 2.  Basic and improved models of branch base diameter for Larix olgensis

      模型名称
      Model name
      模型编号
      Model No.
      基础模型形式
      Basic model form
      模型编号
      Model No.
      改进模型形式
      Improved model form
      Gompertz1${\mathrm{BD}}=a{\rm{e}}^{(-b{\rm{e}}^{-c{\mathrm{DINC}}})} $5${\mathrm{BD}}=a{\rm{e}}^{(-b{\rm{e}}^{-c{\mathrm{DINC}}})}{\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)} $
      Mitscherlich2${\mathrm{BD}}=a(1-b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $6${\mathrm{BD}}=a(1-b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}){\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)} $
      Logistic3${\mathrm{BD}}=a/(1+b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $7${\mathrm{BD}}=a{\rm{e}}^{\left({{\mathrm{DBH}}}^{d}\right)}/(1+b{\rm{e}}^{\left(-c{\mathrm{DINC}}\right)}) $
      单分子式 Monomolecular formula4${\mathrm{BD}}=a(1-{\rm{e}}^{\left(-b{\mathrm{DINC}}\right)}) $8$ {\mathrm{BD}}=a(1-{\rm{e}}^{\left(-b{\mathrm{DINC}}\right)}){\rm{e}}^{\left({{\mathrm{DBH}}}^{c}\right)} $
      注:BD和DINC分别是一级枝条基径和着枝深度;DBH是单木胸径;abcd是待拟合的参数。Notes: BD and DINC are base diameter and depth into crown of primary branch, respectively; DBH is individual tree DBH; a, b, c and d are parameters to be fitted.
    • 不同年龄、不同树木的枝条基径在生长过程上可能存在差异性。为了减少此问题对拟合精度的影响,本研究在最优模型的基础上,建立龄组与单木两水平非线性混合效应模型,模型表达式如下:

      $$ \left\{\begin{array}{l}{y}_{ijk}=f\left({\varnothing }_{ijk,}{\vartheta }_{ijk}\right)+{\varepsilon }_{ijk},i=1,\cdots,M,\\ j=1,\cdots,{M}_{i},k=1,\cdots,{n}_{ij},{\varepsilon}_{ijk}\sim N(0,R)\\ {\varnothing }_{ijk}={{{A}}}_{ijk}\beta +{{{B}}}_{i,jk}{{{u}}}_{i}+{{{B}}}_{ij,k}{{{u}}}_{ij}, \\ {{{u}}}_{i}\sim N\left(0,{{D}}_{1}\right),{{{u}}}_{ij}\sim N\left(0,{{D}}_{2}\right)\end{array}\right. $$ (1)

      式中:M为第1水平龄组数量;$ j $为第$ i $龄组水平内的树木数量;$ k $为第$ i $个龄组中第$ j $棵树包含的枝条个数;$ {y}_{ijk} $$ {\vartheta }_{ijk} $分别为第$ i $个龄组中第$ j $棵树对应的第$ k $个一级枝枝条基径和着枝深度的观测值;$ f $是关于参数向量$ {\varnothing }_{ijk},\;{\vartheta }_{ijk} $的非线性函数;$ \;\beta $$ p\times 1 $维固定效应向量;$ {u}_{i} $$ q\times 1 $维的第1水平随机效应向量,假定服从期望为0、方差为${{D}}_{1} $的正态分布;$ {{{u}}}_{ij} $$ {q}_{2}\times 1 $维的第2水平随机效应向量,假定服从期望为0、方差为${{D}}_{2} $的正态分布;$ {{{A}}}_{ijk} $$\; \beta $的固定效应矩阵,$ {{{B}}}_{i,jk} $$ {{{B}}}_{ij,k} $分别为$ {{{u}}}_{i} $$ {{{u}}}_{ij} $的随机效应设计矩阵;$ {\varepsilon }_{ijk} $为随机误差项,假定服从期望为0、方差为R的正态分布,并假定$ {{{u}}}_{i} $$ {{{u}}}_{ij} $和误差项$ {\varepsilon }_{ijk} $之间相互独立。

    • 以改进后的最优模型为基础,首先将所有参数都作为随机参数进行拟合,通过比较不同随机参数组合的赤池信息准则(AIC)和对数似然值Loglik确定最优的混合效应模型[18]。AIC、BIC越小,Loglik值越大,模型拟合效果越好。对不同参数个数的模型进行似然比检验,以确定最佳参数个数[19]。在AIC、BIC越小,Loglik值越大的前提下,选择参数个数较少的模型。

    • 随机效应的方差协方差矩阵反映了各水平间的差异,选择林业上常用的广义正定矩阵、对角矩阵2种结构进行拟合[9]。以包含3个随机参数的方差协方差结构为例:

      $$ {\text{广义正定矩阵:}}{{D}}=\left[\begin{array}{ccc}{\sigma }_{{b}_{1}}^{2}& {\sigma }_{{b}_{1}{b}_{2}}& {\sigma }_{{b}_{1}{b}_{3}}\\ {\sigma }_{{b}_{1}{b}_{2}}& {\sigma }_{{b}_{2}}^{2}& {\sigma }_{{b}_{2}{b}_{3}}\\ {\sigma }_{{b}_{1}{b}_{3}}& {\sigma }_{{b}_{2}{b}_{3}}& {\sigma }_{{b}_{3}}^{2}\end{array}\right] $$ (2)
      $$ {\text{对角矩阵:}}{{D}}=\left[\begin{array}{ccc}{\sigma }_{{b}_{1}}^{2}& 0& 0\\ 0& {\sigma }_{{b}_{2}}^{2}& 0\\ 0& 0& {\sigma }_{{b}_{3}}^{2}\end{array}\right] $$ (3)

      式中:$ {\sigma }_{{b}_{1}}^{2} $$ {\sigma }_{{b}_{2}}^{2} $$ {\sigma }_{{b}_{3}}^{2} $分别是随机参数$ b_{1} $$ b_{2} $$ b_{3} $的方差,$ {\sigma }_{{b}_{1}{b}_{2}} $为随机参数$ b_{1} $$ b_{2} $的协方差,$ {\sigma }_{{b}_{2}{b}_{3}} $为随机参数$ b_{2} $$ b_{3} $的协方差,$ {\sigma }_{b_{1}{b}_{3}} $为随机参数$ b_{1} $$ b_{3} $的协方差。

    • 在林业上常采用下式描述组内方差协方差结构[9]

      $$ {{{R}}}_{i}={\sigma }^{2}{{{G}}}_{i}^{0.5}{{\varGamma}}_{i}{{{G}}}_{i}^{0.5} $$ (1)

      式中:$ {{{R}}}_{i} $为树木内误差方差协方差矩阵,$ {\sigma }^{2} $为模型的误差方差,$ {{\varGamma }}_{i} $为树木内误差相关性结构,$ {{{G}}}_{i}^{0.5} $为描述方差异质性的对角矩阵。

      为了确定树木内的方差协方差结构,必须解决自相关性和方差异质性的问题,采用指数函数和幂函数矫正异方差现象[20]。比较他们的拟合结果是否差异显著,并通过残差图分析模型改进后的优劣。

      $$ g\left({\mu}_{ij},\alpha \right)={\left|{\mu }_{ij}\right|}^{\alpha } $$ (2)
      $$ g\left({\mu}_{ij},\beta \right)=\mathrm{exp}\left(\beta {\mu}_{ij}\right) $$ (3)

      式中:$ {\mu }_{ij} $是利用固定参数拟合的预测值,$ \alpha $$ \;\beta $是待拟合的参数,g是混合模型的方差。

    • 对混合模型进行检验时,固定效应部分与传统模型的检验方法相同,而随机效应检验需要二次抽样计算随机参数值,本研究采用下式计算随机参数值[6]

      $$ \hat{b}_{i}\approx{\hat{{{D}}}{\hat{{{Z}}}}_{i}^{{\rm{T}}}\left({\hat{{{Z}}}}_{i}\hat{{{D}}}{\hat{{{Z}}}}_{i}^{{\rm{T}}}+{\hat{{{R}}}}_{i}\right)}^{-1}{\hat{e}}_{i} $$ (4)

      式中:$ \hat{{{D}}} $为随机效应参数的方差协方差矩阵,$ {\hat{{{R}}}}_{i} $为树木内方差协方差结构,$ {\hat{{{Z}}}}_{i} $为关于随机参数的偏导数矩阵,$ {\hat{{{Z}}}}_{i}^{{\rm{T}}} $$ {\hat{{{Z}}}}_{i} $的转置矩阵,$ {\hat{e}}_{i} $为实际值减去用固定效应参数计算的预测值的差。

      通过均方根误差(RMSE)、和调整系数($R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$)对基础模型和改进模型精度进行比较,RMSE越小,$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$越接近1,则表明所构建模型的精度越高。通过平均绝对误差MAE、平均相对误差绝对值MRAE检验模型效果[21-22]

      $${\mathrm{RMSE}}=\sqrt{\dfrac{\displaystyle\sum _{i=1}^{n}{\left({y}_{i}-{\hat{y}}_{i}\right)}^{2}}{n}} $$ (5)
      $$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}=1-\dfrac{\displaystyle\sum _{i=1}^{n}{({y}_{i}-{\hat{y}}_{i})}^{2}}{\displaystyle\sum _{i=1}^{n}{({y}_{i}-{\bar{y}}_{i})}^{2}}\left(\frac{n-1}{n-p-1}\right) $$ (6)
      $${\mathrm{AIC}}=-2\mathrm{l}\mathrm{n}\left(L\right)+2p $$ (7)
      $${\mathrm{BIC}}=-2\mathrm{ln}\left(L\right)+p \ln\left(n\right) $$ (8)
      $${\mathrm{MAE}}=\sum _{i=1}^{n}\left|\frac{{y}_{i}-{\hat{y}}_{i}}{n}\right| $$ (9)
      $${\mathrm{MRAE}}=\dfrac{1}{n}\sum _{i=1}^{n}\dfrac{\left|{y}_{i}-{\hat{y}}_{i}\right|}{{y}_{i}} $$ (10)

      式中:$ {y}_{i} $是枝条基径的观测值,$ {\hat{y}}_{i} $是模型的预测值,n为样本个数,L是模型的对数似然值,p为模型参数个数,${\bar{y}}_i $枝条基径观测值的均值。

    • 为了确定枝条基径与林分、单木变量的关系,首先对林分及单木变量进行相关性检验,结果显示枝条基径与胸径的相关系数为0.35,而与树高、冠长的相关系数为0.31、0.32。然后通过图形分析其关系,分别绘制枝条基径与着枝深度、胸径的散点图(图1),可以看出枝条基径与枝条的着枝深度高度相关,随着着枝深度的增加,枝条基径逐渐增大,达到峰值之后,略有减小。当树木胸径越大时,对应单木的枝条基径分布范围越大,对枝条基径的最大值影响较高。由于林木树高、胸径、年龄之间存在高度相关性,着枝深度可通过树高与枝下高进行计算,为了避免同时加入会出现共线性,因此只考虑加入林木胸径进行改进,将林木胸径以指数函数形式加入基础方程。

      图  1  枝条基径随着枝深度、胸径的分布

      Figure 1.  Distribution of branch base diameter with depth into crown and DBH

      对基础及改进模型式进行拟合检验,在基础方程中,除了单分子式(模型4)的调整系数较低之外,其余模型调整系数值均在0.51左右,其中Logistic方程(模型3)的调整系数($R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2}$)最高,均方根误差(RMSE)最小,Gompertz方程次之,Mitscherlich方程最小。说明这3个基础模型都能较好的解释着枝深度与枝条基径的关系。但检验结果却显示在这3个模型中,Gompertz方程(模型1)的平均绝对误差MAE以及平均相对误差绝对值MRAE最小,说明模型(1)相比于模型(2)、(3)能较好的预测枝条基径生长。因此向3个模型中添加变量,拟合结果显示加入胸径之后各模型的调整系数均有提高,其中模型(5)的调整系数最大,RMSE最小,分别为0.597 9和5.518 8 mm,说明以Gompertz方程作为基础式的改进模型能较好地描述胸径与着枝深度对枝条基径生长的影响。检验结果显示模型(5)的MAE为3.727 4 mm,MRAE为0.402 4 mm,在3个改进模型中也达到最优(表3)。并且所有模型的参数均显著(P < 0.05)。且模型(1)与模型(5)方差分析的结果表明两者差异显著(F = 446.26,P < 0.05)。因此将改进后的的Gompertz方程(模型5)作为最优模型,用于构建龄组与样地两水平的非线性混合模型,并进行进一步的比较。

      表 3  基础及改进模型拟合与检验结果

      Table 3.  Fitting and validating results of basic and improved models

      项目
      Item
      模型编号
      Model No.
      拟合优度 Goodness of fitting检验指标 Validating index
      均方根误差 RMSE$ R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $平均相对误差 MAE平均相对误差绝对值 MRAE
      基础
      Basic model
      16.080 70.512 04.660 60.533 8
      26.095 50.509 74.675 30.536 1
      36.076 10.512 94.663 50.536 0
      46.155 40.500 24.687 30.521 2
      改进
      Improved model
      55.518 80.597 93.727 40.402 4
      65.534 90.595 53.746 60.405 9
      75.520 80.597 63.732 20.405 4
      85.552 70.593 13.745 10.397 7
    • 确定最优基础模型后,利用Rstudio软件以及其中的nlme包进行混合模型的计算,考虑随机效应参数所在位置,模型(5)共有210种随机效应参数的组合,其中有146种形式可以收敛。表4列出了不同随机参数个数中最优模型的拟合结果,当考虑2个随机参数时,龄组效应作用于参数c,单木效应作用于参数b,AIC值最小,Loglik值最大(模型5-1)。当考虑3个随机参数,龄组效应作用于参数b,单木效应同时作用于参数c、d时,拟合效果最好(模型5-2)。对于4个随机效应参数,选择龄组效应作用于参数b,单木效应同时作用于参数b、c、d为组合的模型效果最好(模型5-3),此时AIC值为12 621.67,Loglik值最大,为−6 298.833。当模型的随机参数个数分别为5、6、7、8时,模型的AIC、BIC值均高于模型(5-3)。对参数个数不同的模型进行似然比检验,检验结果表明随机参数个数为5、6的混合模型差异不显著(P > 0.05),其余模型之间差异均显著。综合来看,模型(5-3)的AIC、BIC值最小,且模型的随机参数个数为4,因此将模型(5-3)作为最优枝条基径模型。模型形式为

      表 4  不同随机参数最优混合模型拟合结果比较

      Table 4.  Comparison in fitting results of optimal mixed models with different random parameters

      随机参数个数
      Number of random
      parameter
      模型编号
      Model No.
      随机参数 Random parameter统计指标 Statistical index似然比检验 Likelihood ratio (LRT) test
      龄组水平
      Age group level
      单木水平
      Individual tree level
      AICBICLoglik似然比 LRTP
      25-1cb12 728.3512 767.86−6 357.175
      35-2bcd12 645.0512 695.85−6 313.52387.303 6 < 0.000 1
      45-3bbcd12 621.6712 689.40−6 298.83329.380 4 < 0.000 1
      55-4b、db、c、d12 625.6712 704.69−6 298.8340.001 30.999 4
      65-5a、c、db、c、d12 633.6412 729.59−6 299.8191.970 90.578 5
      75-6a、b、c、da、b、d12 665.6312 784.16−6 311.81323.988 70.000 1
      85-7a、b、c、da、b、c、d12 701.6612 842.77−6 325.831212.709 5 < 0.000 1
      $${\mathrm{BD}}=a{{\rm{e}}}^{\left(-\left(b+{u}_{i1}+{u}_{ij1}\right){{\rm{e}}}^{-(c+{u}_{ij2}){\mathrm{DINC}}}\right)}{{\rm{e}}}^{\left({{\mathrm{DBH}}}^{(d+{u}_{ij3})}\right)} $$ (11)

      式中:a、b、c、d为固定参数,$ {u}_{i1} $为龄组水平的随机参数。$ {u}_{ij1}$${u}_{ij2} $${u}_{ij3} $为单木水平的随机参数。

    • 在构建混合模型时,最关键的就是解决模型的随机效应方差−协方差结构以及异方差问题,对广义正定矩阵、对角矩阵2种结构进行拟合,并进行似然比检验(表5)。拟合结果表明广义正定矩阵与对角矩阵有明显差异(P < 0.05),但广义正定矩阵的AIC值更小。因此选择广义正定矩阵作为随机效应的方差-协方差结构。

      表 5  不同方差协方差结构拟合结果比较

      Table 5.  Comparison in fitting results of different variance-covariance structures

      方差−协方差结构
      Variance-covariance structure
      AICBICLoglik似然比
      LRT
      P
      广义正定矩阵 Generalized positive definite matrix12 621.6712 689.4−6 298.83
      对角矩阵 Diagonal matrix12 625.0112 687.1−6 301.505.341 00.020 8
    • 由于本研究中的枝条数据不是连续观测数据,因此不考虑误差的自相关性,在确定随机参数的方差协方差结构之后,通过模型的残差图检验是否存在异方差性。绘制基础与改进模型的残差图(图2),图2a是基础模型的残差分布,图2b是选择广义正定矩阵进行改进后的混合模型的残差分布,可以看出确定随机效应的方差协方差结构之后,混合模型的残差在横轴方向分布更加离散,说明应用混合效应模型使模型预测精度得到改善,但是仍呈现喇叭状分布,这表示模型仍存在异方差问题,利用常用的幂函数、指数函数方差结构进行矫正。

      图  2  基础模型(a)与混合模型包括广义正定矩阵(b)、幂函数(c)、指数函数(d)残差分布

      Figure 2.  Residual distribution of basic model (a) and mixed model with generalized positive definite matrix (b), mixed model with power function (c), and mixed model with exponential function

      图2c和图2d分别是幂函数和指数函数作为异方差结构时模型的残差分布,可以看出,以幂函数作为模型的异方差结构时,残差均匀分布于0值附近,说明模型的异方差性得到明显改善。而指数函数为异方差结构时模型的拟合效果略差。对这两种结构的模型进行似然比检验,结果表明,加入幂函数与指数函数矫正后的模型与未修正的模型有显著差异,都能够减小模型的AIC值。而幂函数作为异方差结构时AIC比使用指数函数修正异方差时小,Loglik更大(表6)。因此,选择幂函数作为混合模型的异方差结构。

      表 6  不同异方差结构混合模型拟合结果比较

      Table 6.  Comparison in fitting results of mixed model with different heteroscedasticity structures

      异方差结构
      Heteroscedasticity structure
      AICBICLoglik似然比
      LRT
      P
      无 None12 621.6712 689.40 −6 298.830
      幂函数 Power function12 151.5312 230.55−6 061.764474.14 < 0.000 1
      指数函数 Exponential function12 220.7212 294.10−6 097.35971.19 < 0.000 1
    • 表7列出了基础模型与混合模型的固定参数值,随机参数的方差协方差结构,幂函数估计值,模型方差值,结果表明,相比于基础模型,混合模型的调整系数值由0.597 9增加到0.699 8,RMSE值由5.518 8 mm减小为4.768 4 mm。利用式(7)计算检验数据中各单木的随机参数值,分别对混合模型的固定效应与随机效应进行检验,表明考虑随机效应之后,无论是MAE还是MRAE,都明显低于只考虑固定效应参数的模型,表明加入随机效应之后模型的预测精度得到了提升(表8)。

      表 7  基础与混合模型拟合结果及检验

      Table 7.  Fitting and validating results of basic and mixed models

      模型形式 Model form参数 Parameter (标准差 SD)拟合优度 Goodness of fitting评价指标 Evaluating index
      abcdRMSE$R_{\mathrm{a}\mathrm{d}\mathrm{j} }^{2} $MAEMRAE
      基础模型(5)
      Basic model (5)
      3.172 82.017 40.384 90.259 35.518 80.597 93.727 40.402 4
      (0.190 4)(0.069 4)(0.007 6)(0.023 4)
      混合模型(5-3)
      Mixed model (5-3)
      2.088 11.882 30.379 60.328 74.768 40.699 83.705 80.391 6
      (0.194 6)(0.164 5)(0.040 4)(0.014 3)
      异方差(幂函数)
      Heteroscedasticity (power function)
      $ \alpha $ = 0.63
      模型方差
      Model variance
      $ {\sigma }^{2} $ = 6.597 8
      方差协方差结构
      Variance-covariance structure
      ${{D} }_1=\left|3.82 \times 10^{-34}\right|$$ {{D}}_2=\left|\begin{array}{ccc}0.450 \; 9& 0.257 \; 0& -0.143 \; 0\\ 0.257 \; 0& 0.021 \; 4& -0.752 \; 0\\ -0.143 \; 0& -0.752 \; 0& 0.000 \; 4\end{array}\right| $
      注:D1D2分别是龄组、单木水平的方差协方差结构。Notes: D1 and D2 are variance-covariance structures of age group and individual tree level, respectively.

      表 8  固定效应与随机效应检验结果的比较

      Table 8.  Comparison of validating results of fixed effect and random effect

      单木编号
      Individual tree No.
      固定效应 Fixed effect随机效应 Random effect
      平均绝对误差 MAE平均相对误差绝对值 MRAE平均绝对误差 MAE平均相对误差绝对值 MRAE
      12.691 10.457 02.309 90.364 1
      22.772 30.296 42.477 10.229 6
      35.490 90.264 64.368 70.227 2
      43.579 90.274 13.420 80.267 7
      53.913 00.469 63.753 10.340 8
      64.746 10.687 52.902 90.359 0
    • 选择改进后的Gompertz方程作为基础模型,这与罗恒春等的研究结果一致,他在构建云南松(Pinus yunnanensis)林分平均胸径生长模型时选择Gompertz方程,表明Gompertz方程比起Logistic方程更适用于描述树木生长[23]。混合模型对分组数据的适用性使它在森林生长建模中更具说服力和实用性[24],其随机参数值反映了对应的变量对枝条基径的影响程度。参数d作用于单木胸径,它与固定参数a结合,决定了该单木枝条基径的最大值,反映了胸径对枝条基径生长的影响。参数b作为生长参数[23],同时含有龄组与单木随机效应,说明不同龄组间、不同单木间枝条基径的生长存在差异。参数c则表明单木不同,着枝深度对枝条基径的影响也不同。基于基础和混合模型结果,加入各单木的随机参数对单木进行预测(图3),可以看出基础模型(图3a)展现了枝条基径的总体变化趋势,但一定程度上将实际值平均化。而应用混合模型(图3b)不会过度增大或降低枝条基径的分布范围,能够反映不同单木间枝条基径的差异。考虑龄组水平,模拟各龄组下枝条基径在不同着枝深度时的分布(图4),模型中包含胸径和着枝深度两个变量,将胸径值固定为对应龄组胸径的平均值,参数同时考虑固定效应和各龄组水平对应的随机参数,可以看出,(1)枝条基径与着枝深度呈正相关。(2)随着年龄的增大,同一着枝深度对应的基径逐渐增大且最大枝条基径出现的位置越靠近树干下部。表明应用混合模型模拟枝条基径能够体现不同层次的差异,得到更精确的拟合结果。这与多位学者应用混合模型探讨树木建模的研究结果相似[6, 8, 17, 19]。姜立春等采用线性混合效应所构建的落叶松枝条基径模型R2值为0.663 2,与本文构建的模型结果进行比较,说明以非线性为基础描述枝条基径的生长更为合理[9, 25]。陈东升等所构建的枝条基径模型以单分子式为基础,R2为0.783 7[13],高于本文所构建的模型,产生这种问题的原因可能是样本数据的范围不同,选取枝条样本时他采取的方式是在每一段选取长势良好的2 ~ 3个标准枝,这种取样方式在拟合时更容易获得较好的拟合优度,我认为应当考虑所有一级枝进行建模更能阐述枝条基径的生长情况。且其研究中落叶松的胸径值分布范围为15.4 ~ 35.7 cm,本文中落叶松的胸径最小值为4.5 cm,涵盖了4个龄组,所以本文所构建的模型更适用于对不同年龄长白落叶松进行预测。

      图  3  枝条基径基础与混合模型预测值分布

      Figure 3.  Distribution of predicted values of branch base diameter by basic and mixed models

      图  4  不同龄组枝条基径分布

      Figure 4.  Distribution of branch base diameter of different age groups

      对不同龄组枝条基径的生长进行分析,由分布图可以看出(图4),随着龄组的增大,枝条基径最大值出现的位置越靠近树干下部,分龄组进行讨论,①幼龄:枝条基径随着枝深度的变化程度较小,在7 ~ 8 m时枝条基径达到峰值,约为14 mm,由于这一时期该地区长白落叶松树高多分布于5 ~ 12 m之间,而且幼龄林的郁闭度相对较小,所以枝条的垂直分布范围广。②中龄:枝条基径多分布于6 ~ 23 mm之间,在枝径值较小时分布较多,在着枝深度为7 ~ 8 m时趋于稳定,平均可达到为22 mm。结合树木生长过程来看,此时树木进入快速生长期,树干中上部新生枝条增多。枝条的生长速度大致保持一致。③近熟,当着枝深度为8 ~ 9 m时枝条基径达到峰值,为32 mm。相比幼龄到中龄阶段,枝条基径生长幅度明显增大,在垂直分布上较为集中,分布于树干下部的枝条开始减少。④成熟:这一时期枝条基径在着枝深度为9 ~ 10 m时达到最大,约为38 mm,枝径范围在10 ~ 38 mm之间,树高在这一时期达到20 m,而由近熟到成熟相较于中龄到近熟枝径增长的幅度较小,此时树木的活枝下高更高,表明这一时期枝条生长速度减缓,下层枝条出现衰亡,这种规律也符合树木生长的特点。这表明通过龄组水平的随机参数可阐述不同年龄枝条基径的分布规律。不需要分别建立各龄组的模型。

      由于数据的有限性,本文主要基于龄组与单木来构建枝条基径混合效应模型,以求研究不同龄组枝条生长的不同之处。除此之外,枝条生长还可能受到林木种植密度与立地质量的影响[10],因此还应该进一步考虑不同林分密度、立地条件等因素对枝条基径生长的影响,构建更适用于该地区的最优枝条基径模型。

    • 研究基于龄组与单木构建两水平长白落叶松一级枝枝条基径非线性混合效应模型,并对不同龄组枝条分布进行模拟。得出以下结论:(1)以改进后的Gompertz方程为基础模型,广义正定矩阵作为随机效应的方差协方差结构,幂函数作为模型的异方差结构。所构建的两水平混合效应模型能很好地解释组间与组内的差异,有效提升预测精度。(2)结合龄组水平的随机参数能合理的表述各龄组枝条基径的生长规律及差异。无论树木处于什么年龄,枝条基径与林木胸径和枝条着枝深度都呈正相关。但枝条基径大小的分布范围为幼龄 < 中龄 < 近熟 < 成熟。幼龄枝条基径垂直分布范围最广,中龄与近熟林枝条基径生长速率较高。到成熟时期生长减缓。总之,此枝条基径模型能够为东折棱河林场不同年龄长白落叶松人工林的经营管理提供理论支持。

参考文献 (25)

目录

    /

    返回文章
    返回