文章信息
- 陈东升, 李凤日, 孙晓梅, 贾炜玮
- Chen Dongsheng, Li Fengri, Sun Xiaomei, Jia Weiwei
- 基于线性混合模型的落叶松人工林节子大小预测模型
- Models to Predict Knot Size for Larch Plantation Using Linear Mixed Model
- 林业科学, 2011, 47(11): 121-128.
- Scientia Silvae Sinicae, 2011, 47(11): 121-128.
-
文章历史
- 收稿日期:2010-01-18
- 修回日期:2010-03-22
-
作者相关文章
2. 东北林业大学林学院 哈尔滨 150040
2. School of Forestry, Northeast Forestry University Harbin 150040
落叶松(Larix spp.)作为东北地区的主要造林树种,面积约400多万hm2,占东北人工林总面积的70%左右,如何提高其木材质量对于该地区人工林产业化发展具有重要意义。节子的大小和分布严重影响到木材质量,导致木材质量的降低(DIN-EN,1997; 中华人民共和国国家质量监督检验检疫总局,2004)。节子分为健全节和疏松节,健全节是枝条从形成到死亡时未与树干分离的部分; 疏松节是枝条死亡后被树干包藏起来的部分,也是容易脱落的部分,所以疏松节对板材的质量影响最大,会造成板材中出现空洞(Vestøl et al., 2000),同时节子还会增大木材受潮和真菌感染的概率(Shigo,1989)。因此,只有了解节子的各种属性(大小、角度、长度)及其形成机制才能找到合适的控制方法,从而提高木材质量。节子的形成过程取决于树木的生物学进程和外部的环境条件,如大雪或风力等对其有很大的影响,有研究表明,在好的立地条件、自然环境和林分条件下,树木疏松节长度要短一些(Lamsa et al., 1990;Hägg et al., 1995)。
目前一些学者已经对节子的形成机制及其大小特征进行了研究,但建立的节子预测模型都集中在采用传统线性和非线性模型(应用最小二乘法估计参数)上(Mäkinen et al., 2002; 贾炜玮,2006; Kantola et al., 2007; Hein et al., 2007a; 2007b; 2007c),这些模型都假定误差是服从独立同分布的,而现实中获取的节子观测数据很难满足这一条件,势必会对估计结果有一定的影响;另外,其主要说明研究对象的平均行为,忽略了个体或群体之间存在的序列相关性及差异性,例如以往节子大小模型很好地反映了节子大小平均变化规律,但却忽略了林木个体及样地之间的差异。而混合模型通过构造固定效应参数和随机效应参数(唐守正等,2002),消除了层次结构数据产生的误差,既能反映总体的平均变化趋势, 又可以提供数据方差、协方差等多种信息来反映个体之间的差异(李春明,2009a),在很大程度上弥补了以往模型估计方面的不足,提高了模型的预估精度。国外已经有学者将混合模型应用到树木节子研究中,如Moberg(2001)应用混合模型建立了欧洲云杉(Picea abies)节子内部特征模型,Hein(2008)采用41株山欧洲水青冈(Fagus sylvatica)解析木资料,建立了疏松节形成时间和所占节子比例模型,取得了很好的拟合效果。国内混合模型的应用主要集中在单木和林分生长模型方面(李永慈等,2004; 李春明,2009b; 雷向东等,2009),还没有对树木节子开展相关研究。因此,本文应用黑龙江省东北部地区孟家岗林场95株落叶松解析木、4 073个节子剖析数据,采用线性混合模型建立落叶松节子大小(直径、着生角度、健全节和疏松节的长度)预测模型,改善和提高模型的预估精度,从而为提高落叶松人工林集约化经营水平提供科学依据。
1 数据收集与整理研究区域位于黑龙江省佳木斯市孟家岗林场,其地理坐标为130°32′42″—130°52′36″E,46°20′—46°30′50″N。该地区地处完达山西麓余脉,以低山丘陵为主,坡度较为平缓,平均海拔250 m。属东亚大陆性季风气候,冬季漫长、寒冷且干燥; 夏季短促、温暖而湿润; 早春少雨、风大易干旱; 秋季降温迅速,常有冻害发生。年平均气温2.7 ℃,极端最高气温35.6 ℃,最低气温-34.7 ℃,年平均降水量为550 mm,无霜期120天左右。土壤以暗棕壤为主。
2007—2008年,在具有代表性的落叶松人工林中设置19块标准地(长25 m,宽40 m,面积0.1 hm2),在各标准地中测定每株树的胸径、树高、冠幅、枝下高和林木坐标(表 1)。每块标准地按照等断面积径级标准木法,按林木大小划分5个径级,将各径级林木的平均胸径、树高和冠长作为选取解析木的标准,并在每块标准地附近选伐样木。
本研究采用节子剖析技术获得各解析木节子生长数据,具体方法如下:在树冠基部以下、含有死枝或节子的轮枝附近,截取10~30 cm长的木段,并对每轮的死枝或节子进行编号,同时测量每个死枝或节子离地面的高度及方位角。在每个木段中选取2个最大的节子,沿树干方向用手持油锯通过树干髓心进行纵剖,得到节子的纵剖面(图 1)。
枝条生长发育过程中,可以分为枝条形成时的年龄、枝条停止形成年轮时的年龄、枝条死亡时的年龄和节子愈合(包藏)时的年龄。节子从枝条形成时到枝条死亡时的部分称为健全节(sound knot),而从节子死亡时到被完全包藏的部分称为疏松节(loose knot),疏松节的大小直接影响木材质量和等级。在节子纵剖面上测量节子的直径(DE)、角度(θK)、健全节长度(Lsound)、疏松节长度(Lloose)以及枝条各个时期的年龄等因子,各节子变量的统计值详见表 2。本研究获取95林解析木,其中13块标准地65株解析木的2 740个节子剖析数据作为模型拟合数据,剩余6块标准地30株解析木的1 333个节子剖析数据作为模型检验数据(表 2)。
随机效应的线性混合模型为:
(1) |
式中: Y为一个由来自不同样地、不同树木节子大小因子组成的n维矩阵,x为已知设计矩阵,β为模型中未知的固定参数向量; z为随机效应设计矩阵,u为模型中未知的随机参数向量,e是未知的随机误差向量。其中u~N(0,G),e~N(0,R),G与R二者之间不相关。对G和R必须选择其协方差结构,因此,它允许数据间存在相关性或异方差。因为样地和树木是随机选择的,所以随机效应包括样地间、树木间和树木个体内3个部分(唐守正等,2002)。
对于随机效应,包括2种情况: 1)随机截距模型(random intercept model),假定在不同时间的相关性相同,即只有一个随机效应,通过截距来体现; 2)随机截距和斜率模型(random coefficient model),即不同样地的截距和斜率都不相同。在本研究中,只考虑随机截距模型。
本研究中, 模型自变量选择过程用方差膨胀因子(VIF)来判断自变量间的多重共线性。一般认为,当VIF>10时,有严重的共线性。此时,剔除共线性较严重的自变量,保留共线性弱而对因变量贡献大的自变量。只有回归系数显著(P < 0.05)、且方差膨胀因子小于5的因子才能进入模型。
用混合模型对参数进行估计时,首先要确定合适的误差协方差结构。本研究通过比较不同协方差结构的AIC(akaike information criterion)指标来选择,AIC越小越好(Burnham et al., 2002)。采用最大似然法估计进行模型参数估计和方差求解(Verbeke et al., 2000)。
2.2 节子预测模型节子直径模型采用可转化为线性混合效应模型进行预测,所确定的最优模型为:
(2) |
式中: DKij为第i块样地第j株解析木节子直径(cm); HKij为第i块样地第j株解析木节子的高度(m); DBHij为第i块样地第j株解析木胸径(cm); ui,vij为样地(i=1~13)和样木的随机作用(i=1~13,j=1~65);eij为误差。
节子着生角度按照活枝着生角度的方法进行研究。枝条的分枝角度是衡量树木空间分布能力的重要指标,其向空间扩展的能力影响枝叶对光照、温度、CO2等的利用。本研究通过一个变化后的线性混合模型对节子的角度进行模拟,所确定的最佳模型为:
(3) |
式中: θKij为第i块样地第j株解析木节子与树干的倾斜角度。
预测节子长度的模型采用Mobery(2001)研究云杉所建立的节子长度预测模型,包括健全节长度和疏松节长度模型,健全节长度模型应用节子变量和树木变量建立,疏松节长度模型采用节子相对大小(DKrel)和胸径(DBH)为自变量,使用一个可转化为线性模型的自然对数模型进行预测。模型形式为:
(4) |
式中: Lsound ij为第i块样地第j株解析木节子健全节长度(cm); DKrel ij为第i块样地第j株解析木相对节子大小(DK/DKmean); DKmean为1 m区分段内节子的平均直径(cm)。
(5) |
式中: Lloose ij为第i块样地第j株解析木节子疏松节长度(cm)。
2.3 模型拟合优度指标和检验文中采用RSS,MSE,R2这3个拟合统计量作为评价模型拟合效果的指标。
剩余残差平方和(RSS):
(6) |
剩余均方差(MSE):
(7) |
相关指数(R2):
(8) |
应用表 2中模型的检验数据,采用以下5个统计量(Li et al., 2001)对所建模型进行独立性检验:
1) 平均偏差(mean error):
(9) |
2) 均绝对偏差(mean absolute error):
(10) |
3) 均相对偏差(mean percent error, M%E):
(11) |
4) 均相对偏差绝对值(mean absolute percent error, MA%E):
(12) |
5) 预估精度
(13) |
式中:
基于表 2中模型的拟合数据,首先进行模型自变量选择,采用方差膨胀因子(VIF)来判断自变量间是否具有多重共线性,判断后建立了上述节子大小预测模型。然后使用SAS 9.1统计软件对所建立的节子预估模型[方程(2),(3),(4),(5)]进行拟合(SAS 9.1,2004),参数估计值见表 3。从表 3中可以看出: 4个节子预测模型的参数估计值都通过了t检验,说明在建立节子预估模型中选取的自变量是合理的。节子直径预测模型的R2达到了0.53,MSE为0.012,模型的残差均匀分布在±0.6之间(图 2a); 节子角度模型的拟合精度最高,R2达到了0.84,MSE为0.018,残差均匀分布在±0.6之间(图 3a); 健全节长度模型的R2为0.57,MSE为1.25,残差均匀分布在±6之间(图 4a); 疏松节长度模型的R2为0.48,MSE为0.039,残差均匀分布在-0.8~+0.6之间(图 5a)。4个模型中拟合效果最差的为疏松节长度预估模型,因为疏松节的长度受到外界环境条件的影响很大,所以导致疏松节长度的预测比较困难,模型的预测精度较低,如果能在模型中加入外界的环境因子,可能会提高模型的精度。
通过对节子直径模型分析发现:节子直径随其着生高度的增大而增大,大约在树干的中部达到最大值,然后会随着生高度的增大而逐渐减小,在接近树冠的部位达到最小值(图 6)。这和活枝直径在树冠内的变化规律是相同的; 而且胸径越大的树木,节子直径越大,这是因为胸径大的树木光合能力强,对空间的占有能力也强,会产生大的枝条,枝条死亡后就会形成大的节子。节子角度变化规律和直径的变化规律相似,节子角度也随胸径的增大而增大,但是随着生高度的增大而逐渐减小(图 7),因为一般顶部的枝条着生角度都小于其下部枝条的着生角度,所以才会产生这样的结果。
通过分析健全节长度与预测因子之间的关系可知:健全节长度随节子着生高度增大而逐渐变大,也随节子相对直径的增大而增大。但是,随着节子相对直径增大健全节长度增大的程度有所不同,在节子相对直径小时,健全节长度的增加速度缓慢; 当节子相对直径大时,健全节长度增加速度明显加快(图 8)。这个变化规律可以表明:节子直径越大,其内部健全节的长度也越大,在同样高度的位置,直径小的枝将比直径大的枝死的早。在健全节拟合模型中包括了胸径,说明树木胸径的变化也会影响健全节的长度,树木越大,产生的健全节也越长。疏松节长度随着胸径和节子相对直径的增大而增大,在研究中发现,疏松节长度在节子小时变化趋势明显,但当节子变大后,疏松节长度的增大趋势逐渐减慢(图 9),这可能是由于节子越大会导致树干包藏这个节子就比较困难的原因。
通过方程的拟合优度指标对混合效应模型和固定效应模型的拟合效果进行比较(表 4),从表 4中可以看出: 4个模型中混合效应模型的R2都大于固定模型,而RSS和MSE值均小于固定模型,说明混合模型提高了拟合精度。模型的残差分布图中明显可以看出混合效应对预测结果有很大的影响,节子直径混合模型的残差分布比较均匀,集中在±0.6之间(图 2a),而不加混合作用的固定效应模型的残差集中在±1.0之间(图 2b); 节子角度混合模型的残差分布明显好于固定作用模型,不再出现残差聚集在一起的情况(图 3a,b); 健全节长度混合模型的残差分布和固定模型的残差分布相比,混合模型的残差沿着横坐标方向进行了拉伸,使得残差更加分散(图 4a,b); 疏松节长度混合模型的残差分布在-0.8~+0.6之间,而固定效应模型的残差分布在-2.0~+1.5之间(图 5a,b)。总的来说,混合模型的残差分布更趋于分散和均匀,说明混合模型可以对预测值进行调整,对模型的预估更加精确,充分表明样地和样木间的随机作用是明显的,因此随机效应的作用在节子大小的预测模型中是不可忽略的。
应用表 2中的模型独立样本检验数据,通过检验程序,对模型(2),(3),(4),(5)进行检验,检验结果见表 5。从表 5中可以看出:建立的节子大小预测模型的ME都很小,而预估精度都很高,达到了90%以上,说明这些模型可以用来描述节子的变化规律,对未来节子大小进行预测。
文中通过树木胸径和节子所在高度对节子直径进行预测,发现节子直径随胸径的增大而增大,并且会在树干中部达到最大值,然后又随着生高度的增大而逐渐减小; 节子的着生角度可以用节子直径、着生高度和树木胸径来描述,节子角度变化规律和直径的变化规律相似,节子角度也随胸径的增大而增大,但是随着生高度的增大而逐渐减小; 健全节长度受到很多因子的影响,包括胸径、节子直径和节子高度,胸径大的树健全节的长度在靠近地面部分随着生高度的增大有明显的增大趋势,但是在靠近活冠时开始缓慢地减小,胸径小的树健全节的长度在靠近地面时随着生高度的增大有缓慢的增大然后在靠近活冠时有显著的增大趋势,而且健全节长度和节子相对直径大小有关,在节子相对直径小时,健全节长度的增加速度缓慢,当节子相对直径大时,健全节长度增加速度明显加快。疏松节长度的模拟则采用可以转化为对数的线性混合模型进行拟合,疏松节的长度随胸径和节子相对直径的增大而增大。总的来说,这些节子的模拟研究可以作为分析木材质量的理论基础,为进一步提高落叶松木材质量的研究奠定基础。
本文采用混合模型对节子直径、着生角度以及长度进行了模拟。研究发现混合模型的预测精度高于固定效应模型(模型的R2都明显提高,剩余残差平方和、均方误差均显著减小),表明在节子预测模型中加入随机效用是可以提高模型精度的。因此说明广义的传统模型理论有所不足,仅能够有助于讨论固定参数的估计,无助于讨论随机参数的估计及推断。随机效应参数的引入使得混合模型对误差分布的描述更为精细,参数的估计结果更为精确。最后对建立的节子混合模型固定效应部分进行独立样本检验,结果表明模型的预测精度较高,可以很好地描述节子的变化规律。
贾炜玮. 2002. 落叶松人工林树冠构筑型及枝生长动态的研究. 东北林业大学硕士学位论文. http://cdmd.cnki.com.cn/Article/CDMD-10225-2002121972.htm
|
贾炜玮. 2006. 樟子松人工林枝条生长与节子大小预测模型的研究. 东北林业大学博士学位论文. http://cdmd.cnki.com.cn/article/cdmd-10225-2006111189.htm
|
雷相东, 李永慈, 向玮. 2009. 基于混合模型的单木断面积生长模型[J]. 林业科学, 45(1): 74-80. DOI:10.11707/j.1001-7488.20090113 |
李春明. 2009a. 利用非线性混合模型进行杉木林分断面积生长模拟研究[J]. 北京林业大学学报, 31(1): 44-49. |
李春明. 2009b. 混合效应模型在森林生长模型中的应用[J]. 林业科学, 45(4): 131-138. |
李凤日. 2004. 长白落叶松人工林树冠形状的模拟[J]. 林业科学, 40(5): 16-24. DOI:10.11707/j.1001-7488.20040503 |
李永慈, 唐守正. 2004. 用Mixed和Nlmixed过程建立混合生长模型[J]. 林业科学研究, 17(3): 279-283. |
宋喜芳, 李建平, 胡希远. 2009. 模型选择信息量准则AIC及其在方差分析中的应用[J]. 西北农林科技大学学报:自然科学版, 37(2): 88-92. |
唐守正, 李勇. 2002. 生物数学模型的统计学基础[M]. 北京: 科学出版社.
|
中华人民共和国国家质量监督检验检疫总局. 2004. 中华人民共和国国家标准——原木检验(GB/T 144—2003). 北京: 中国标准出版社.
|
Burnham K P, Anderson D R. 2002. Model selection and multimodel inference: a practical information-theoretic approach[M]. New York: Springer Science+Business Media, Inc.
|
DIN-EN. 1997. Deutsche Norm: Laub-rundholz-qualitats-sortierung-teil 1:eiche und buche. Deutsche Fassung EN 1316-1, Standard, German.
|
Hägg A, Weslien H. 1995. Samband mellan tallens(Pinus sylvestris) naturliga grenrensing och ståndorten. Department of Forest Production, Swedish University of Agricultural Science, Uppsala, Report, 246: 25.
|
Hein S, Lenk E, Kladtke J, et al. 2007a. Z-Baum orientierte Auslesedurchforstung in Buche (Fagus sylvatica L.):Auswirkungen auf Qualitat, Sortenstruk tur und Wertleistung[J]. Allgemeine Forstund Jagdzeitung, 178: 1-20. |
Hein S, Makinen H, Yue C, et al. 2007b. Modelling branch characteristics of Norway spruce from wide spacing in Germany[J]. Forest Ecology and Management, 242(2/3): 155-164. |
Hein S, Spiecker H. 2007c. Comparative analysis of occluded branch characteristics for Fraxinus excelsior and Acer pseudoplatanus with natural and artificial pruning[J]. Canadian Journal of Forest Science, 37(8): 1414-1426. DOI:10.1139/X06-308 |
Hein S. 2008. Knot attributes and occlusion of naturally pruned branches of Fagus sylvatica[J]. Forest Ecology and Management, 256(12): 2046-2057. DOI:10.1016/j.foreco.2008.07.033 |
Kantola A, Mäkinen H, Makelä A. 2007. Stem from and branches of Norway spruce as a sawn timber——predicted by a process based model[J]. Forest Ecology and Management, 241(1/3): 209-222. |
Lämsvä P, Kellomäki S, Vaisanen H. 1990. Nuorten mantyjen oksikkuuden riippuuvuus puuston rakenteesta ja kasvupaikan viljavuudesta[J]. Folia For, 746: 1-22. |
Li F R, Choi J K, Kim J H. 2001. The development of growth and yield models for the natural broadleaved-Korean pine forests in northeast China[J]. Journal of Korean Forestry Society, 90(5): 650-662. |
Moberg L. 2001. Models of internal knot properties for Picea abies[J]. Forest Ecology and Management, 147(2/3): 123-138. |
Mäkinen H, Song T. 2002. Evaluation of models for branch characteristics of Scots pine in Finland[J]. Forest Ecology and Management, 158(1/3): 25-39. |
Shigo A L. 1989. A New Tree Biology. Shigo and Trees Associates, Durham.
|
SAS 9.1. 2004. SAS/STAT 9.1 User's Guide. SAS Publishing, Cary.
|
Verbeke G, Molenberghs G. 2000. Linear Mixed Models for Longitudinal Data. Springer, New York/Berlin/Heidelberg, Geert Verbeke Geert Molenberghs.
|
Vestøl G I, Høibø O A. 2000. Prediction of knot diameter in Picea abies (L.) Karst[J]. European Journal of Wood and Wood Products, 59(1/2): 129-136. |