洪李斌，卿蘊(yùn)賢，田佳赫，康潔敏，盧偉

(河北農(nóng)業(yè)大學(xué) 林學(xué)院，河北 保定 071000)

胸徑作為林業(yè)中最基本的測(cè)樹(shù)因子之一，是林分生長(zhǎng)狀況和生長(zhǎng)潛力判定的重要依據(jù)，也是估算森林蓄積量、生物量、碳儲(chǔ)量等指標(biāo)的重要因子，為預(yù)估森林碳匯潛力，提高森林經(jīng)營(yíng)管理水平，實(shí)現(xiàn)我國(guó)碳達(dá)峰碳中和的林業(yè)發(fā)展目標(biāo)提供重要參考[1-7]。塞罕壩機(jī)械林場(chǎng)擁有全球最大的人工林基地，是世界生態(tài)文明建設(shè)史上的典范，2018年塞罕壩林場(chǎng)實(shí)現(xiàn)了第一筆造林碳匯交易，標(biāo)志著塞罕壩林場(chǎng)在碳匯項(xiàng)目上邁出實(shí)質(zhì)性的步伐，有力詮釋了“綠水青山就是金山銀山”的理念，實(shí)現(xiàn)了森林生態(tài)效益和經(jīng)濟(jì)效益雙贏。華北落葉松(Larixprincipis-rupprechtii)是塞罕壩地區(qū)主要造林樹(shù)種，占林場(chǎng)有林地總面積的70.6%[8]，對(duì)塞罕壩地區(qū)水源涵養(yǎng)、防風(fēng)固沙和維護(hù)該地區(qū)生態(tài)系統(tǒng)健康具有至關(guān)重要的作用[9]。因此，構(gòu)建塞罕壩地區(qū)華北落葉松胸徑預(yù)測(cè)模型對(duì)估算本地區(qū)森林碳匯儲(chǔ)量，預(yù)測(cè)森林碳匯潛力具有重要意義。

針對(duì)單木胸徑預(yù)測(cè)模型已有大量研究，于士濤基于塞罕壩林內(nèi)160株華北落葉松解析木數(shù)據(jù)構(gòu)建了華北落葉松估測(cè)胸徑傳統(tǒng)經(jīng)驗(yàn)?zāi)Ｐ蚚10]；張冬燕在塞罕壩地區(qū)利用固定樣地調(diào)查數(shù)據(jù)基于啞變量方法構(gòu)建包含不同齡組的混合效應(yīng)胸徑生長(zhǎng)模型[11]；周澤宇基于塞罕壩地區(qū)的標(biāo)準(zhǔn)地調(diào)查數(shù)據(jù)構(gòu)建了華北落葉松人工林直徑分布線性混合效應(yīng)模型[12]。樹(shù)皮生長(zhǎng)量與木材生長(zhǎng)量的速度并不相同，利用胸徑尺測(cè)量數(shù)據(jù)構(gòu)建胸徑預(yù)測(cè)模型會(huì)引入樹(shù)皮因子的誤差，F(xiàn)inch研究表明如果不考慮樹(shù)皮生長(zhǎng)量，那么樹(shù)木的帶皮胸徑生長(zhǎng)量和材積的估算就會(huì)產(chǎn)生0～6%的誤差[13]。傳統(tǒng)經(jīng)驗(yàn)?zāi)Ｐ投嗉僭O(shè)立地條件不變，而樹(shù)木之間的競(jìng)爭(zhēng)也會(huì)使同一年齡的林木胸徑產(chǎn)生一定的差異，競(jìng)爭(zhēng)能力與其對(duì)周圍生存資源的占有及利用存在聯(lián)系[14-15]。以往研究還表明樹(shù)木生長(zhǎng)和周圍樹(shù)木的相對(duì)空間位置和相互之間的距離顯著相關(guān)，即林木的胸徑生長(zhǎng)受到鄰木大小和密度制約[16-17]。混合效應(yīng)模型是一種同時(shí)考慮時(shí)間變化和區(qū)域分布效應(yīng)的有效模型，現(xiàn)已在單木以及林分等研究中取得了廣泛的應(yīng)用[18-23]。因此開(kāi)展對(duì)華北落葉松單木去皮胸徑預(yù)測(cè)模型的研究可以為林業(yè)調(diào)查提供技術(shù)支持，對(duì)華北落葉松的經(jīng)營(yíng)管理十分必要。

本研究以河北省塞罕壩機(jī)械林場(chǎng)的華北落葉松人工林為研究對(duì)象，于2018年7月在37塊20 m×20 m的固定樣地中鉆取木芯數(shù)據(jù)，與傳統(tǒng)胸徑尺測(cè)量胸徑相比，這種方法能夠避免樹(shù)皮因子的影響，可以得到確切的年齡與去皮胸徑數(shù)據(jù)。通過(guò)引入樣地競(jìng)爭(zhēng)因子林分密度和單木競(jìng)爭(zhēng)因子相對(duì)直徑來(lái)構(gòu)建塞罕壩地區(qū)華北落葉松單木去皮胸徑預(yù)測(cè)的混合效應(yīng)模型，以期能為塞罕壩地區(qū)華北落葉松的生長(zhǎng)和經(jīng)營(yíng)提供參考依據(jù)。

1 數(shù)據(jù)與方法

1.1 研究區(qū)概況和數(shù)據(jù)來(lái)源

河北省塞罕壩機(jī)械林場(chǎng)(E 116°53′～117°39′，N 41°92′～42°36′)地處陰山山脈、大興安嶺余脈向西南延伸與燕山余脈的交合點(diǎn)，地形上分南部壩下和北部壩上2部分。該區(qū)域?qū)俸疁匦源箨懠撅L(fēng)氣候，海拔1 010～1 939.9 m。年平均氣溫-1.3 ℃，年均降水量460 mm，蒸發(fā)量1 339.2 mm，年均無(wú)霜期64 d，積雪時(shí)間長(zhǎng)達(dá)7個(gè)月。土壤主為灰色森林土，該區(qū)域的主要喬木樹(shù)種有華北落葉松、樟子松、白樺、云杉等[24]。

在塞罕壩機(jī)械林場(chǎng)華北落葉松純林中設(shè)置37塊固定樣地，樣地大小為20 m×20 m。于2018年7月在每個(gè)樣地中選擇長(zhǎng)勢(shì)良好的華北落葉松，使用口徑5.5 mm的生長(zhǎng)錐在距離地面1.3 m高度處，從2個(gè)不同的方向鉆取2個(gè)樣芯，并原則上保證每個(gè)樣芯通過(guò)髓心以及保留樹(shù)皮。將取好樹(shù)芯粘到木質(zhì)凹槽中，經(jīng)晾干、打磨后進(jìn)行交叉定年。測(cè)量并記錄每個(gè)年齡的去皮胸徑。建模數(shù)據(jù)概況見(jiàn)表1。

表1 建模數(shù)據(jù)統(tǒng)計(jì)表Table 1 Statistics table of modeling data

1.2 研究方法

1.2.1 變量選擇及基礎(chǔ)模型的構(gòu)建 林木胸徑估測(cè)模型通常是基于林齡、樹(shù)高、冠幅等自變量與胸徑關(guān)系進(jìn)行構(gòu)建。然而，在實(shí)際測(cè)量中，樹(shù)木高度和冠幅的數(shù)據(jù)采集成本相對(duì)較高且精度不高。相比之下，林齡是一個(gè)較為容易獲取的數(shù)據(jù)，尤其在人工林中。因此，本文首先選擇了14個(gè)以年齡為唯一自變量的年齡-去皮胸徑模型來(lái)估測(cè)去皮胸徑，從中選擇合適的模型進(jìn)行混合效應(yīng)模型的構(gòu)建。具體模型形式見(jiàn)表2。

表2 14個(gè)基于林齡-去皮胸徑的胸徑預(yù)測(cè)模型Table 2 The 14 models based on the age-inside bark diameter at breast height for the prediction of DBH

1.2.2 混合效應(yīng)模型的構(gòu)建 利用方差膨脹系數(shù)(Variance inflation factor，VIF) 檢驗(yàn)?zāi)Ｐ投嘀毓簿€性，依據(jù)VIF值<10原則，挑選出3個(gè)變量，分別為林分密度(D)、林木直徑與林分平均直徑之比(RH)和立地指數(shù)(SI)，由于變量的增加使混合效應(yīng)模型過(guò)復(fù)雜，故依據(jù)Pearson相關(guān)系數(shù)刪去對(duì)方程貢獻(xiàn)較小的立地指數(shù)(SI)。選取林分密度(D)以及林木直徑與林分平均直徑之比(RH)作為隨機(jī)效應(yīng)構(gòu)建華北落葉松去皮胸徑估測(cè)模型。

單水平線性混合效應(yīng)模型：

DBH=a0+a0m0i+(a1+m1i)A+εik(εik)～N(0，Ri) (m0im1i)T～N(0，Dm)

(1)

DBH=a0+a0u0j+(a1+u1j)A+εjk(εjk)～N(0，Rj) (u0ju1j)T～N(0，Du)

(2)

嵌套兩水平線性混合效應(yīng)模型：

圖書(shū)館閱讀推廣平臺(tái)后臺(tái)設(shè)計(jì)包括系統(tǒng)管理、用戶管理、閱讀推廣服務(wù)、資源鏈接、讀者積分管理、日志管理及其他管理模塊。具體的設(shè)計(jì)如圖2所示。

DBH=a0+a0m0i+a0u0ij+(a1+m1i+u1ij)A+εijk(εijk)～N(0，Rij) (m0im1i)T～N(0，Dm) (u0iju1ij)T～N(0，Du)

(3)

式中：A為林木年齡，i為密度水平編號(hào)，j為相對(duì)直徑編號(hào)，k為樣木編號(hào)，m0i、m1i為密度水平隨機(jī)效應(yīng)參數(shù)，u0i、u1i為相對(duì)直徑水平隨機(jī)效應(yīng)參數(shù)；εik、εjk、εijk為殘差向量，Dm、Du隨機(jī)參數(shù)的方差-協(xié)方差矩陣。在混合效應(yīng)模型構(gòu)建的過(guò)程中需要明確3個(gè)問(wèn)題。

(1)模型參數(shù)的確定

(4)

(5)

(6)

(7)

(2)誤差方差-協(xié)方差結(jié)構(gòu)的確定。

林業(yè)中常用于描述方差-協(xié)方差結(jié)構(gòu)表達(dá)式為：

(8)

本研究采用廣義正定矩陣(General positive-definite matrix)來(lái)作為隨機(jī)效應(yīng)的方差-協(xié)方差結(jié)構(gòu)。其結(jié)構(gòu)如下：

(9)

以上計(jì)算過(guò)程均在R軟件(版本4.0.3)中完成，華北落葉松去皮胸徑估測(cè)基礎(chǔ)模型由R中的“l(fā)m”和“nls”函數(shù)完成；混合效應(yīng)模型其建模和預(yù)測(cè)有R包的“NLME”完成。

2 結(jié)果與分析

2.1 基礎(chǔ)模型的擬合

在R軟件中將建模數(shù)據(jù)代入到14個(gè)林齡-去皮胸徑估測(cè)模型中，對(duì)各個(gè)模型進(jìn)行擬合，擬合結(jié)果見(jiàn)表3。

表3 基礎(chǔ)模型的建模結(jié)果Table 3 Modeling results of base models

由表3可知，14個(gè)基礎(chǔ)模型中除了Model.3以外的都達(dá)到了0.8以上，表明各模型對(duì)去皮胸徑的擬合效果較好。其中，Model.1、Model.4、Model.5、Model.7、Model.12、Model.13、Model.14決定系數(shù)均達(dá)到了0.85，說(shuō)明相對(duì)于其他7個(gè)模型擬合效果更好。其中，Model.1的赤池信息準(zhǔn)則(AIC)值低于其余模型，同時(shí)考慮到實(shí)際工作中適用性，本研究選擇Model.1作為基礎(chǔ)模型。

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

考慮不同隨機(jī)參數(shù)的組合，對(duì)基礎(chǔ)模型進(jìn)行擬合，根據(jù)表4可見(jiàn)，基于競(jìng)爭(zhēng)效應(yīng)的華北落葉松去皮胸徑估測(cè)混合效應(yīng)模型共有8種收斂的情況，MEM.2、MEM.4、MEM.6、MEM.7以及MEM.8模擬效果優(yōu)于傳統(tǒng)模型，其中MEM.8模擬精度最高。具體MEM.8參數(shù)估計(jì)見(jiàn)表5。

表4 去皮胸徑估測(cè)線性混合模型收斂結(jié)果Table 4 Simulated results of estimation inside bark diameter at breast height (DBH) on linear mixed model

表5 MEM.8去皮胸徑估測(cè)線性混合模型的參數(shù)估計(jì)Table 5 Parameter estimation of linear mixed model of MEM.8 inside bark diameter at breast height (DBH) estimation

2.3 傳統(tǒng)模型和混合效應(yīng)模型的比較

本研究通過(guò)對(duì)傳統(tǒng)模型以及不同水平的混合效應(yīng)模型進(jìn)行比較分析，探尋能在塞罕壩地區(qū)更精確預(yù)測(cè)華北落葉松去皮胸徑的模型形式,詳見(jiàn)表6和圖1。

表6 不同模型擬合統(tǒng)計(jì)量Table 6 Fitting statistics of different models

(a)

(b)圖1 傳統(tǒng)模型(a)與混合效應(yīng)模型(b)的的殘差圖Figure 1 Residual plots of the traditional model (a) and the mixed effects model (b)

由表6可知，相對(duì)直徑水平混合效應(yīng)模型比密度水平混合效應(yīng)模型更能體現(xiàn)競(jìng)爭(zhēng)對(duì)林木生長(zhǎng)的影響；與傳統(tǒng)模型和單水平混合效應(yīng)模型相比，嵌套兩水平的混合模型擬合精度顯著提高。根據(jù)圖1可以看出，混合效應(yīng)模型的殘差離散程度較小，說(shuō)明混合效應(yīng)模型的擬合效果較好。

3 討論

林木生長(zhǎng)是一個(gè)復(fù)雜的過(guò)程，個(gè)體之的差異往往受到立地條件、光照、溫度以及林木相互之間競(jìng)爭(zhēng)等多方面的的影響，故僅通過(guò)年齡推斷去皮胸徑是不準(zhǔn)確的。本研究在僅考慮年齡的情況下，其均方根誤差(RMSE)達(dá)到2 cm左右，誤差較大，這與周潤(rùn)凱等人的研究結(jié)果一致[26]。因此在構(gòu)建去皮胸徑預(yù)測(cè)模型時(shí)必須要考慮林木的生存環(huán)境，而劉泰瑞在關(guān)帝山的研究表明，華北落葉松受不同競(jìng)爭(zhēng)壓力后胸徑生長(zhǎng)量差異顯著，由此得知林木胸徑對(duì)競(jìng)爭(zhēng)壓力的響應(yīng)往往超過(guò)其對(duì)立地條件、光照、溫度等其他因素的響應(yīng)[27-28]?？紤]到所要構(gòu)建模型能兼顧實(shí)用性、復(fù)雜性和精確性的要求，因此除年齡外，競(jìng)爭(zhēng)是必須的建模因子。在引入競(jìng)爭(zhēng)因子構(gòu)建華北落葉松去皮胸徑混合效應(yīng)模型后，預(yù)測(cè)精度被大幅提高，尤其是嵌套兩水平的混合效應(yīng)模型，其均方根誤差(RMSE)達(dá)到了0.29 cm，均方根誤差降低了85%，本結(jié)論與王玉濤、楊惠濱的研究結(jié)果相一致[29-30]?；旌闲?yīng)模型的擬合精度優(yōu)于傳統(tǒng)經(jīng)驗(yàn)?zāi)Ｐ?，其原因主要是混合效?yīng)模型既可以反映總體的平均變化趨勢(shì)，又可以提高數(shù)據(jù)方差、協(xié)方差等多種信息來(lái)反映個(gè)體之間的差異[31]。李春明、符利勇等人運(yùn)用混合效應(yīng)模型顯著提高了單木胸徑生長(zhǎng)量的擬合精度[32-33]。

為進(jìn)一步探討競(jìng)爭(zhēng)因子對(duì)華北落葉松去皮胸徑-年齡關(guān)系的影響，分別建立密度水平、相對(duì)直徑水平以及嵌套兩水平的混合效應(yīng)模型。結(jié)果表明，單水平混合模型中，相對(duì)直徑水平擬合精度優(yōu)于密度水平，而嵌套兩水平的混合效應(yīng)模型擬合效果達(dá)到了最優(yōu)，明顯提高了對(duì)去皮胸徑的預(yù)測(cè)精度。說(shuō)明在預(yù)測(cè)單木去皮胸徑時(shí)，僅考慮林分尺度的競(jìng)爭(zhēng)因子是不夠的，更需要考慮單木之間競(jìng)爭(zhēng)的差異。嵌套兩水平混合效應(yīng)模型擬合精度進(jìn)一步提高，體現(xiàn)對(duì)林分尺度競(jìng)爭(zhēng)因子分組后考慮單木競(jìng)爭(zhēng)效應(yīng)對(duì)去皮胸徑生長(zhǎng)影響更為準(zhǔn)確，且嵌套兩水平混合效應(yīng)模型還可以分析不同隨機(jī)效應(yīng)對(duì)去皮胸徑生長(zhǎng)影響程度。本研究數(shù)據(jù)來(lái)源于塞罕壩機(jī)械林場(chǎng)，在實(shí)際應(yīng)用中需要考慮區(qū)域性差異以及實(shí)際需要，選擇不同層次進(jìn)行模擬和預(yù)測(cè)。

本研究數(shù)據(jù)中缺乏大于42年的華北落葉松木芯數(shù)據(jù)，具有一定的局限性，后續(xù)應(yīng)當(dāng)對(duì)數(shù)據(jù)進(jìn)行補(bǔ)全，完善模型。盡管在引入競(jìng)爭(zhēng)因子后，模型達(dá)到了很高的精度，但桑杰通過(guò)建立氣候影響下的青岡櫟胸徑生長(zhǎng)預(yù)測(cè)混合效應(yīng)模型，模擬不同氣候條件下青岡櫟的生長(zhǎng)趨勢(shì)，為華北落葉松胸徑預(yù)測(cè)模型的研究提供了新的思路[34]。以后的研究需考慮將競(jìng)爭(zhēng)指數(shù)作為參數(shù)直接用于構(gòu)建胸徑模型，再引入氣候等因子進(jìn)一步從機(jī)理上解釋胸徑的增長(zhǎng)，探討光照、溫度等氣候因子對(duì)預(yù)測(cè)模型的敏感性和貢獻(xiàn)性。最后，通過(guò)研究華北落葉松去皮胸徑隨年齡變化的規(guī)律，可進(jìn)一步研究華北落葉松帶皮胸徑與去皮胸徑的關(guān)系，探尋樹(shù)皮生長(zhǎng)量隨林齡增加的變化規(guī)律，達(dá)到改進(jìn)華北落葉松胸徑生長(zhǎng)模型的目的。

4 結(jié)論