廣西桉樹地徑與胸徑多元混合效應(yīng)模型研建

曾春陽，劉 峰，徐慶玲，李 恒，莫祝坤，黃正深

（1.廣西壯族自治區(qū)林業(yè)勘測設(shè)計院，廣西南寧 530011；2.凌云縣林業(yè)局，廣西凌云 533199；3.自然資源部北部灣經(jīng)濟(jì)區(qū)自然資源監(jiān)測評價工程技術(shù)創(chuàng)新中心，廣西南寧 530011）

桉樹（Eucalyptusspp.）為世界3 大速生樹種之一，具有生長快、產(chǎn)量高等優(yōu)點(diǎn)，在我國南方地區(qū)被廣泛種植[1]。作為廣西主要商品林樹種之一，截至2021年，廣西桉樹人工林種植面積約303萬公頃，蓄積量約1.86 億立方米（《廣西壯族自治區(qū)2021 年林業(yè)生態(tài)資源狀況報告》），為我國的木材供應(yīng)提供了強(qiáng)大支撐。在廣西主要桉樹種植區(qū)，違法采伐林木及風(fēng)雪和火災(zāi)等自然災(zāi)害是桉樹人工林發(fā)展的限制因素[2]。在對林木損失進(jìn)行評估的過程中，如何根據(jù)遺留地徑估算林木材積是目前林業(yè)司法鑒定面臨的一項技術(shù)難題。

在實(shí)際生產(chǎn)中，常利用一元材積表獲得林木材積；林木采伐后無法測得胸徑，通過建立地徑與胸徑的回歸關(guān)系模型估測胸徑是解決被伐木蓄積量測算問題的一種有效方法[3]。李輝[4]采用11 種常見的一元數(shù)學(xué)模型擬合馬尾松（Pinusmassoniana）地徑與胸徑的相關(guān)關(guān)系，結(jié)果表明冪函數(shù)模型為最優(yōu)模型。高德祥等[5]對比線性方程、對數(shù)方程、二次方程、三次方程和冪函數(shù)5個模型的擬合效果，結(jié)果顯示三次方程為擬合滇西南地區(qū)巨尾桉（E.grandis×E.urophylla）根徑與胸徑相關(guān)關(guān)系的最優(yōu)數(shù)學(xué)模型。李海龍[6]以呂梁山南段天然油松（Pinustabuliformis）林調(diào)查數(shù)據(jù)為樣本，運(yùn)用模型優(yōu)選法，確定一元線性模型為擬合天然油松胸徑與根徑相關(guān)關(guān)系的最佳回歸模型。趙浩彥等[7]運(yùn)用多模型對比選優(yōu)法，確定一元線性模型在擬合南京地區(qū)馬尾松根徑與胸徑相關(guān)關(guān)系方面優(yōu)度最高；并以林齡、林分密度等樣地因子為混合效應(yīng)因子構(gòu)建多元混合效應(yīng)模型，結(jié)果表明多元混合效應(yīng)模型具有更高的擬合優(yōu)度和精度。

本研究以廣西桉樹人工林為研究對象，組織建模樣本和驗(yàn)證樣本，分別擬合一元模型和加入樣地因子的多元混合效應(yīng)模型，建立廣西桉樹地徑與胸徑相關(guān)關(guān)系模型，可為廣西被伐桉樹材積計算提供理論支撐。

1 材料與方法

1.1 樣本采集區(qū)

根據(jù)桉樹人工林的地理位置、氣候、地形地貌和土壤等條件[8]，將廣西桉樹人工林劃分為6個樣本采集區(qū)，分別為桂南區(qū)、桂西區(qū)、桂中區(qū)、桂東區(qū)、桂北區(qū)和沿海區(qū)[9]，共70個樣地（表1）。桂南區(qū)水熱條件良好（年均降水量2 100 mm，年均氣溫22.5 ℃），多為丘陵和低山地貌，土壤厚度為厚至中等，土壤養(yǎng)分含量較高，為廣西桉樹主產(chǎn)區(qū)。桂西區(qū)海拔較高（平均海拔1 120 m），雨量較少（年均降水量1 115 mm），年均氣溫22.1 ℃，為廣西桉樹一般產(chǎn)區(qū)。桂中區(qū)多為丘陵和低山地貌，土壤厚度為厚至中等，土壤養(yǎng)分含量較高，年均降水量1 600 mm，年均氣溫21.5 ℃，為廣西桉樹主產(chǎn)區(qū)。桂東區(qū)土壤養(yǎng)分含量較高，氣候適宜（年均降水量1 465 mm，年均氣溫21.1 ℃），為廣西桉樹主產(chǎn)區(qū)。桂北區(qū)年均降水量1 949 mm，年均氣溫18.9 ℃，為廣西桉樹邊緣產(chǎn)區(qū)。沿海區(qū)雨量充沛，年均降水量2 200 mm，年均氣溫23.0 ℃，主要為沿海臺地地貌，地勢平緩，土壤較深厚、疏松，但土壤養(yǎng)分含量較低，需通過施肥促進(jìn)桉樹生長[10]，為廣西桉樹主產(chǎn)區(qū)。

表1 樣本采集區(qū)分布Tab.1 Distributions of sample collecting areas

1.2 樣地設(shè)置與調(diào)查

目前，廣西種植的桉樹主要有尾葉桉（E.urophylla）和巨尾桉，鄧恩桉（E.dunnii）、大花序桉（E.cloeziana）等種植面積較少。本研究選擇尾葉桉、巨尾桉為建模樹種。在各個區(qū)域采用典型抽樣法選取有代表性的林分設(shè)置臨時樣地[11]。樣地設(shè)置面積不小于667 m2；在密度較小的林分，樣地面積擴(kuò)大至1 000 m2，確保樣木株數(shù)不少于50 株。選擇與林緣距離10 m 以上、林分內(nèi)林木生長正常且有代表性的地段設(shè)置樣地。作業(yè)過程中，采用羅盤儀定向，采用測繩引線等確定樣地邊界，確保邊界線閉合差不超過1/200[12]。樣地林齡為1 ～17 年，林分密度為800 ～2 700 株/公頃，海拔為15 ～800 m，包括上坡、中坡、下坡、全坡和平地5種坡位類型。

1.3 樣木采集

共采集樣木929株。在70塊樣地中采集538株樣木，其中桂南區(qū)102 株，桂西區(qū)92 株，桂東區(qū)91株，桂北區(qū)83株，桂中區(qū)85株，沿海區(qū)85株，用于模型擬合；選取5 塊樣地進(jìn)行每木檢尺，共391 株，用于模型適用性檢驗(yàn)。為降低地徑變異性對模型精度的影響，將林木主干離地面5 cm 高處直徑作為地徑值[13]。

1.4 數(shù)據(jù)預(yù)處理

結(jié)合經(jīng)驗(yàn)判斷法[14]、繪圖分析法[15]和3倍標(biāo)準(zhǔn)差法[13]，剔除胸徑、地徑異常數(shù)據(jù)。

1.5 模型選擇

1.5.1 一元模型

比較線性、對數(shù)、逆函數(shù)、雙曲線、三次曲線、復(fù)合函數(shù)、冪函數(shù)、S 曲線、生長曲線、指數(shù)函數(shù)和Logistic 函數(shù)11 種常用一元模型的擬合效果，選出適用于擬合廣西桉樹地徑與胸徑相關(guān)關(guān)系的模型。不同模型公式為[7]：

式中，x為自變量，即林木地徑（cm）；y為因變量，即林木胸徑（cm）；c0、c1、c2和c3為模型參數(shù)。

1.5.2 多元混合效應(yīng)模型

為探索樣地因子對地徑與胸徑相關(guān)關(guān)系模型參數(shù)的影響，采用混合效應(yīng)模型建模方法，構(gòu)建以胸徑為自變量、地徑和樣地因子為自變量的多元混合效應(yīng)模型。

多元混合效應(yīng)模型包括微觀方程和宏觀方程。微觀方程用于分析單株樣木的變量關(guān)系，宏觀方程用于分析整個樣地水平因子與變量的關(guān)系。結(jié)合兩種方程可更好地體現(xiàn)變量間的關(guān)系[7]。

宏觀方程：

式中，x為樣木地徑（cm）；y為樣木胸徑（cm）；β0和β1分別為微觀方程的截距和斜率；ε 為微觀方程的殘差項，且ε ～N（0，σ2）；γ00、γ01、γ02、γ03和γ0n為模型參數(shù)；x1、x2、x3和xn為樣地因子；μ為宏觀方程的殘差項，且μ～N（0，σμ2）。采用限制性最大似然法估計兩種方程殘差的方差；采用Wald統(tǒng)計量判斷微觀方程參數(shù)在宏觀方程中是否具有隨機(jī)效應(yīng)。

1.6 模型評價指標(biāo)

采用建模數(shù)據(jù)擬合模型后，需對模型進(jìn)行評估。關(guān)鍵的評價指標(biāo)包括確定指數(shù)（R2）、赤池信息量準(zhǔn)則（Akaike Information Criterion，AIC）、估計值的標(biāo)準(zhǔn)差（Standard Error of Estimate，SEE）、總相對誤差（Total Relative Error，TRE，%）、平均系統(tǒng)誤差（Mean Systematic Error, MSE，%）、平均預(yù)估誤差（Mean Prediction Error,MPE，%）和平均百分標(biāo)準(zhǔn)誤差（Mean Percentage Standard Error, MPSE，%）7 項[16-17]。通過指標(biāo)分析，可以評估模型的解釋能力、估計精度和預(yù)測精度等，從而選出最優(yōu)模型。各指標(biāo)計算公式為[16-17]：

式中，yi為實(shí)際調(diào)查值；i為模型預(yù)估值；為樣本平均值；l為模型函數(shù)最大似然函數(shù)值；p為參數(shù)個數(shù)；n為樣本單元數(shù)；tα為置信水平α 時的t值。R2可衡量模型解釋因變量變化的能力，其值越接近1.000，模型解釋能力越強(qiáng)。AIC 可同時評估模型的復(fù)雜度和精確度，其值越小，模型擬合程度越高。SEE 可衡量估計值的精度，其值越小，估計精度越高。TRE、MSE、MPE 和MPSE 可衡量預(yù)測精度，其值越接近0，預(yù)測精度越高。一般要求R2在0.800 以上，MPE 小于5%，TRE 和MSE 在±3%范圍內(nèi)，MPSE不超過10%。

使用驗(yàn)證樣本數(shù)據(jù)檢驗(yàn)?zāi)Ｐ偷倪m用性。檢驗(yàn)指標(biāo)包括絕對殘差（Bias）和標(biāo)準(zhǔn)偏差（RMSE），其值越小，模型精度越高，計算公式為[7]：

式中，yi為實(shí)際調(diào)查值；i為模型預(yù)估值；n為樣本單元數(shù)。

1.7 數(shù)據(jù)處理

采用SPSS軟件進(jìn)行一元回歸分析，擬合不同林木地徑與胸徑的相關(guān)關(guān)系；采用R 語言中的函數(shù)包（nlme）建立混合效應(yīng)模型。

2 結(jié)果與分析

2.1 一元模型

2.1.1 模型擬合結(jié)果

除逆函數(shù)模型外，其他模型的確定指數(shù)均大于0.800，表明地徑與胸徑的相關(guān)性較強(qiáng)（表2）。其中，冪函數(shù)、線性、雙曲線、三次曲線和S 曲線模型的擬合效果較好，確定指數(shù)均大于0.900。

表2 一元模型擬合結(jié)果Tab.2 Fitting results of one-variable models

對比相關(guān)指數(shù)，選擇線性、雙曲線、三次曲線、冪函數(shù)和S 曲線模型進(jìn)行擬合指標(biāo)計算。5 個模型的總相對誤差和平均系統(tǒng)誤差均小于3%，說明這些模型不存在趨勢性系統(tǒng)誤差；平均預(yù)估誤差均小于5%，平均百分標(biāo)準(zhǔn)誤差均小于10%，說明這些模型可以較準(zhǔn)確地估算林木胸徑（表3）。綜合考慮相關(guān)指數(shù)和赤池信息量準(zhǔn)則，一元線性模型更適合作為擬合廣西桉樹地徑與胸徑相關(guān)關(guān)系的模型，其回歸方程為y= 0.016 + 0.747x。

表3 模型評價結(jié)果Tab.3 Evaluation results of models

2.1.2 模型檢驗(yàn)

5 個模型的絕對殘差均小于1.000，標(biāo)準(zhǔn)偏差均小于2.000，均可用于估測廣西桉樹的胸徑（表4）。其中，一元線性模型的絕對殘差最?。?0.001），標(biāo)準(zhǔn)偏差最小（0.976），精度較高，適用性強(qiáng)。

表4 模型檢驗(yàn)Tab.4 Model test

2.2 多元混合效應(yīng)模型

2.2.1 模型擬合結(jié)果

將林齡、林分密度、海拔和坡位4個樣地因子作為隨機(jī)參數(shù)進(jìn)行非線性混合擬合。結(jié)果表明，廣西桉樹的胸徑與林齡、海拔、坡位和地徑均呈正相關(guān)關(guān)系，與林分密度呈負(fù)相關(guān)關(guān)系（表5 ～6）。Wald統(tǒng)計量檢驗(yàn)結(jié)果表明，微觀方程的截距在宏觀水平上存在顯著隨機(jī)效應(yīng)，斜率在宏觀水平上的隨機(jī)效應(yīng)不顯著。多元混合效應(yīng)模型擬合和評價結(jié)果顯示，其確定指數(shù)（R2=0.966）、赤池信息量準(zhǔn)則（AIC=434.7）、估計值的標(biāo)準(zhǔn)差（SEE = 0.946）、總相對誤差（TRE=0.004%）、平均系統(tǒng)誤差（MSE = 0.000 05%）、平均預(yù)估誤差（MPE = 1.067%）和平均百分標(biāo)準(zhǔn)誤差（MPSE = 5.156%）均比一元線性模型表現(xiàn)好，說明多元混合效應(yīng)模型擬合效果更佳，精度更高，能更準(zhǔn)確地估算林木胸徑。

表5 多元混合效應(yīng)模型擬合結(jié)果Tab.5 Fitting results of multivariate mixed effect model

表6 多元混合效應(yīng)模型評價結(jié)果Tab.6 Evaluation results of multivariate mixed effect model

γ00、γ01、γ02、γ03和γ04為樣地隨機(jī)因子參數(shù)；β1為地徑參數(shù)；x1、x2、x3和x4分別為林齡、林分密度、海拔和坡位。γ00,γ01,γ02,γ03andγ04are random factor parameters of sample plots;β1is ground diameter parameter;x1,x2,x3andx4are stand age,stand density,elevation and slope position,respectively.

2.2.2 模型檢驗(yàn)

多元混合效應(yīng)模型的絕對殘差為0.001，標(biāo)準(zhǔn)偏差為0.930，均小于1.000，不存在趨勢性系統(tǒng)誤差，可用于估測桉樹胸徑。與一元線性模型相比，多元混合效應(yīng)模型的檢驗(yàn)效果更優(yōu)，模型適用性更強(qiáng)。

對模型殘差進(jìn)行Shapiro-Wilk 正態(tài)性檢驗(yàn)，一元線性模型P= 0.004，多元混合效應(yīng)模型P=0.059；從模型殘差分布圖可以看出，多元混合效應(yīng)模型殘差分布范圍小于一元線性模型，表明在預(yù)估胸徑時，多元混合效應(yīng)模型精度更高（圖1）。

圖1 一元線性模型(a)和多元混合效應(yīng)模型(b)殘差分布Fig.1 Residual distributions of one-dimensional linear model(a)and multivariate mixed effect model(b)

3 討論與結(jié)論

研究結(jié)果顯示，一元線性模型能較好地反映廣西桉樹地徑與胸徑的相關(guān)關(guān)系，該結(jié)論與多項研究結(jié)果一致[5，18]。這一發(fā)現(xiàn)對于預(yù)測桉樹生長和發(fā)展趨勢有實(shí)際意義，可為林業(yè)經(jīng)營提供參考。地徑與胸徑的相關(guān)關(guān)系可反映桉樹生長速度和形態(tài)特征，也可為林木選擇和林分結(jié)構(gòu)調(diào)整提供參考。

不同立地條件下，桉樹人工林的生長情況存在差異。在建立多元混合效應(yīng)模型時，選擇林分密度、林齡、海拔和坡位4個樣地因子為隨機(jī)因子。林分密度影響林木間的資源分配和生長競爭，進(jìn)而影響地徑與胸徑生長。隨林齡增長，林木胸徑增大，林木間的競爭進(jìn)一步加劇。海拔和坡位通過影響水熱條件使林木表現(xiàn)出不同的生長規(guī)律[19]。結(jié)果表明，多元混合效應(yīng)模型的擬合精度優(yōu)于一元線性模型。這一結(jié)果驗(yàn)證了樣地因子對桉樹生長的重要性[20]，可為合理的林分規(guī)劃和管理提供參考。

桉樹地徑與胸徑的相關(guān)關(guān)系還受多種因素影響。本研究未考慮氣溫、降水等自然影響因素，及施肥、撫育等人為影響因素[21]。未來可將這些因素納入模型，評估它們對桉樹生長的影響，更全面地了解桉樹生態(tài)特征。

本研究對廣西桉樹地徑與胸徑關(guān)系進(jìn)行深入研究，在一元模型的基礎(chǔ)上，考慮林分密度、林齡、海拔和坡位樣地因子對變量相關(guān)關(guān)系的影響，構(gòu)建地徑-胸徑多元混合效應(yīng)模型。該模型可更好地預(yù)測桉樹胸徑，有助于廣西桉樹人工林管理和經(jīng)營，實(shí)現(xiàn)可持續(xù)發(fā)展目標(biāo)。

利益沖突：所有作者聲明無利益沖突。

作者貢獻(xiàn)聲明：曾春陽負(fù)責(zé)試驗(yàn)調(diào)查、數(shù)據(jù)收集與分析、論文撰寫和文獻(xiàn)檢索；劉峰負(fù)責(zé)項目實(shí)施及研究計劃和試驗(yàn)設(shè)計的制定與執(zhí)行；徐慶玲負(fù)責(zé)桂南、桂中區(qū)的數(shù)據(jù)收集與分析；李恒負(fù)責(zé)桂北、桂東區(qū)的數(shù)據(jù)收集與分析；莫祝坤負(fù)責(zé)沿海區(qū)的數(shù)據(jù)收集與分析；黃正深負(fù)責(zé)桂西區(qū)的數(shù)據(jù)收集與分析。