朱念福，張懷清*，崔澤宇，楊廷棟，李永亮，劉 華

(1.中國林業(yè)科學(xué)研究院資源信息研究所，北京 100091；2.國家林業(yè)和草原局林業(yè)遙感與信息技術(shù)重點實驗室，北京 100091)

枝下高是指樹冠第一個活枝距離地面的高度，是反映樹木生長活力和冠形特征的重要指標(biāo)[1-2]。目前林木生長可視化模擬主要包括結(jié)構(gòu)功能模型對不同發(fā)育階段的單木細(xì)節(jié)形態(tài)模擬[3-5]，以及基于林木生長模型和三維渲染模擬等形式[6-12]，僅通過樹木模型整體的伸縮和替換來實現(xiàn)，均未考慮枝下高自身的變化規(guī)律，無法很好體現(xiàn)林木形態(tài)的差異性，也無法表現(xiàn)林木生長過程中枝下高的動態(tài)變化。

目前，國內(nèi)外學(xué)者對枝下高模型研究中主要考慮單木因子(樹高、胸徑等)、林分統(tǒng)計因子(林分?jǐn)嗝娣e、冠幅競爭因子等)和立地指數(shù)對枝下高的影響[13-19]。如Ritchie等[15]以冠幅競爭因子、林分?jǐn)嗝娣e和立地指數(shù)等預(yù)測變量，構(gòu)建了美國太平洋西北地區(qū)14個樹種的枝下高模型；Rijal等[17]引入樹木尺度因子和競爭因子，構(gòu)建了北美阿卡迪地區(qū)13個樹種的枝下高模型；馬載陽等[20]構(gòu)建了湖南省黃豐橋國有林場的杉木枝下高模型；盧康寧[21]分析了杉木分枝生物量與分枝結(jié)構(gòu)的關(guān)系，構(gòu)建與分枝輪數(shù)有關(guān)的杉木枝下高模型，這些模型均需要林木分枝的生長輪數(shù)和年齡信息支持，未考慮林木空間結(jié)構(gòu)對枝下高的影響。枝下高變化主要通過自然整枝實現(xiàn)，除了自身因素，還受周圍環(huán)境競爭影響，而對象木與其鄰近木組成的空間結(jié)構(gòu)對林木生長和自然整枝具有較大影響[22-25]。

本研究以杉木(Cunninghamialanceolata)人工林為例，考慮林木自身屬性(胸徑、樹高等)和林木所處空間結(jié)構(gòu)單元對枝下高的影響，構(gòu)建出基于空間結(jié)構(gòu)的枝下高模型。應(yīng)用林分三維模型實時生成方法，構(gòu)建出枝干可控的杉木三維模型。根據(jù)單木初始胸徑，通過單木胸徑連年生長量模型、樹高曲線模型、冠幅面積估計模型和考慮空間結(jié)構(gòu)的枝下高模型擬合林木生長狀態(tài)，結(jié)合Unity3D渲染引擎，最終實現(xiàn)杉木枝下高的可視化模擬。

1 材料與方法

1.1 試驗區(qū)概況

試驗區(qū)位于湖南省攸縣黃豐橋國有林場(113°04′～113°43′E，26°43′～27°06′ N)，平均氣溫17.8 ℃，年均降水量為1 410.8 mm，無霜期292 d，年平均日照時間1 612 h，屬于亞熱帶季風(fēng)濕潤氣候。地貌以丘陵、山地為主，境內(nèi)最高海拔1 270 m，最低海拔115 m，坡度一般為20°～35°。森林覆蓋率90%以上，樹種以杉木、松類為主。

1.2 數(shù)據(jù)來源

采用常規(guī)的測樹方法，使用卷尺測量胸徑(DBH)，使用測高桿測量樹高和枝下高(HCB)，使用塔尺測量冠幅(CW)，利用全站儀測量樣地內(nèi)杉木的相對位置坐標(biāo)(X,Y,Z)。根據(jù)2012年調(diào)查數(shù)據(jù)，選取地位指數(shù)[25]相同的6塊杉木人工林臨時樣地，共計2 498組。調(diào)查數(shù)據(jù)按照3∶1的比例隨機劃分?jǐn)M合數(shù)據(jù)(1 874組)和檢驗數(shù)據(jù)(624組)。樣地核心區(qū)(距邊界5 m)屬性分布如表1所示。

表1 樣地數(shù)據(jù)統(tǒng)計量

1.3 空間結(jié)構(gòu)因子構(gòu)建

根據(jù)惠剛盈等[24]提出的最佳空間結(jié)構(gòu)單元，選擇最近4株鄰木參與空間結(jié)構(gòu)因子的構(gòu)建，本研究構(gòu)建了3個空間結(jié)構(gòu)因子：水平空間結(jié)構(gòu)參數(shù)、垂直空間結(jié)構(gòu)參數(shù)和空間結(jié)構(gòu)單元平均距離。

1.3.1 水平空間結(jié)構(gòu)參數(shù)

在水平空間維度，對象木主要受到鄰近木的冠幅擠壓[26-27]。鄰近木的冠幅擠壓抑制了對象木的樹冠生長,促進(jìn)對象木自然整枝。而且同樣程度的冠幅擠壓，對于不同大小林木的影響是有差異的[28]，可以用鄰近木與對象木胸徑的比值來反映這種差異。據(jù)此建立表示對象木受最近4株鄰木冠幅擠壓的水平空間結(jié)構(gòu)參數(shù)PHi(圖1A)，表達(dá)式如下：

(1)

式中：Cwi為對象木i的冠幅；Cwj為鄰近木j的冠幅；Dj為鄰近木j的胸徑；Di為對象木i的胸徑；dij為對象木i與鄰近木j的水平距離，n為鄰近木株數(shù)，n=4。

1.3.2 垂直空間結(jié)構(gòu)參數(shù)

在垂直空間維度，由于鄰近木的上方遮蓋，抑制了對象木的生長，促進(jìn)對象木自然整枝[23]。據(jù)此建立表示對象木受最近4株鄰木上方遮蓋的垂直空間結(jié)構(gòu)參數(shù)PVi(圖1B)，表達(dá)式如下：

(2)

式中：Hj為鄰近木j的樹高；Hi為對象木i的樹高；ΔZ為鄰近木j與對象木i的空間位置高度差。

1.3.3 空間結(jié)構(gòu)單元平均距離

已有研究表明林分密度對枝下高有較為顯著的影響[29]，因此對象木與鄰近木之間的距離顯然對對象木的枝下高有影響。據(jù)此建立表示對象木與最近4株鄰木的空間結(jié)構(gòu)單元平均距離(dDISi)(圖1C)，表達(dá)式如下：

(3)

1.4 模型選擇

本研究選擇5個常用枝下高(HCB)與樹高(H)、其他變量(X)的關(guān)系模型作為候選基礎(chǔ)模型，各模型表達(dá)式見表2。

表2 枝下高模型

樹高和胸徑是預(yù)測枝下高最重要的兩個變量[13,19]。因此，本研究首先擬合只含樹高和胸徑變量的枝下高模型，選出最優(yōu)的枝下高基礎(chǔ)模型。本研究構(gòu)建的空間結(jié)構(gòu)參數(shù)表示了林木受到鄰近木的競爭影響，一定程度上可以作為競爭因子。考慮競爭因子在常用枝下高模型中的表達(dá)形式，將空間結(jié)構(gòu)因子線性加入所選擇的枝下高基礎(chǔ)模型中，選取綜合指標(biāo)較好的枝下高模型。考慮空間結(jié)構(gòu)的枝下高模型表示如下：

(4)

式中：HCB為枝下高，DBH為胸徑，f(H,X)為枝下高模型表達(dá)形式，H為樹高，X為除樹高外其他變量，Y為空間結(jié)構(gòu)因子及組合，a、b、c、c1、c2、c3為方程參數(shù)。然后通過調(diào)查數(shù)據(jù)進(jìn)行擬合和檢驗，最終選出綜合指標(biāo)較好的因子構(gòu)建考慮空間結(jié)構(gòu)的枝下高模型。

1.5 模型評價與檢驗

本研究采用決定系數(shù)(R2)、估計值的標(biāo)準(zhǔn)差[SEE，式中記為σ(SEE)]進(jìn)行模型評價；采用平均誤差[ME,式中記為σ(ME)]，平均絕對誤差[MAE，式中記為σ(MAE)]、總相對誤差[TRE，式中記為σ(TRE)]和平均系統(tǒng)誤差[MSE，式中記為σ(MSE)]進(jìn)行模型檢驗[30]。其中，R2和SEE是回歸模型的最常用指標(biāo)，反映了模型的擬合優(yōu)度。TRE和MSE是反映擬合效果的重要指標(biāo)，二者都應(yīng)該控制在一定范圍內(nèi)(如±3%或±5%)。6個指標(biāo)的具體計算公式如下：

(5)

(6)

(7)

(8)

(9)

(10)

1.6 枝下高模擬

為了表示杉木枝下高的變化，除了本研究需要構(gòu)建的枝下高模型，還使用了其他研究在同一區(qū)域擬合得到的3個相關(guān)杉木生長模型：單木胸徑連年生長量模型、樹高曲線模型、冠幅面積估計模型。單木胸徑連年生長量模型[31]表示如下：

(11)

樹高曲線模型[32]表示如下：

(12)

冠幅面積估計模型[33]表示如下：

CWa=DBH/(2.861 7-0.034 98DBH)。

(13)

式中：t為時間；RS為相對植距；RD為相對優(yōu)勢度；SI為地位指數(shù)；H為樹高，DBH為胸徑，CWa為冠幅面積。

李永亮等[34]基于林分三維模型實時生成方法，構(gòu)建了林分中形態(tài)各異的杉木三維模型，且三維模型的干、枝和葉是分離的，較好地表現(xiàn)了林分中杉木的形態(tài)結(jié)構(gòu)，但無法表現(xiàn)杉木枝下高的變化。

本研究僅構(gòu)建樹干和樹枝分離的杉木三維模型，使用上述方法實現(xiàn)林分三維模型構(gòu)建，實現(xiàn)枝下高動態(tài)調(diào)整功能。使用下述方法模擬樹干和枝干的變化。

1)樹木胸徑：

DTrunk2=DTrunk1+fD(DTrunk1)。

(14)

式中：DTrunk2為下一時刻的樹木胸徑，DTrunk1為當(dāng)前時刻的樹木胸徑，fD(DTrunk1)為單木胸徑連年生長量模型擬合的胸徑變化量。

2)樹木樹高：

HTrunk2=fH(DTrunk2)。

(15)

式中：HTrunk2為下一時刻的樹木樹高，fH(DTrunk2)為樹高曲線模型擬合的杉木樹高。

3)枝干高度：

HCB2=fHCB(HTrunk2,DTrunk2);

(16)

HBranch2i=HCB2+(HBranch1i-HCB1)·(HTrunk2/HTrunk1)。

(17)

式中：HCB2為下一時刻枝下高，HCB1為當(dāng)前時刻枝下高，fHCB(HTrunk2,DTrunk2)為枝下高模型擬合的枝下高，HBranch2i為下一時刻樹枝i的高度，HBranch1i為當(dāng)前時刻樹枝i的高度，HTrunk1為當(dāng)前時刻樹木樹高。

4)枝干長度：

CW2=fCW(DTrunk2);

(18)

LBranch2i=LBranch1i·(CW2/CW1)。

(19)

式中：CW2為下一時刻冠幅，CW1為當(dāng)前時刻冠幅，fCW(DTrunk2)為冠幅面積估計模型擬合的冠幅，LBranchi2i為下一時刻樹枝i的長度，LBranch1i為當(dāng)前時刻樹枝i的長度。

5)枝干基徑：

DBranch2i=DBranch1i·(CW2/CW1)·(HTrunk2i/HBranch2i)。

(20)

式中：DBranch2i為下一時刻樹枝i的基徑，DBranch1i為當(dāng)前時刻樹枝i的基徑。

2 結(jié)果與分析

2.1 枝下高模型選取和檢驗

最初擬合只含H和DBH為變量的枝下高模型，選擇最優(yōu)的枝下高基礎(chǔ)模型。枝下高基礎(chǔ)模型的參數(shù)估計值、擬合優(yōu)度和檢驗結(jié)果見表3。

表3 枝下高基礎(chǔ)模型的擬合優(yōu)度及檢驗結(jié)果

由表3可知，5個模型的參數(shù)估計值標(biāo)準(zhǔn)誤均較小，說明模型各參數(shù)的穩(wěn)定性較好。5個模型的相關(guān)系數(shù)(R2)相差不大，均在0.717左右，估計值標(biāo)準(zhǔn)差SEE為1.40～1.41 m。5個模型檢驗結(jié)果區(qū)別不明顯，所有模型的平均誤差(ME)為-0.29～-0.30 m，平均相對誤差(MAE)均小于1.02 m，總相對誤差(TRE)為-2.90%～-3.00%，平均系統(tǒng)誤差(MSE)為-2.90%～-3.20%。其中，模型A的總相對誤差(TRE)和平均系統(tǒng)誤差(MSE)均最小，分別為-2.906%和-2.931%，平均誤差(ME)也較小，為-0.293 m，且模型A的表達(dá)形式相對簡單，符合生物學(xué)意義，因此選用模型A作為加入空間結(jié)構(gòu)參數(shù)的枝下高基礎(chǔ)模型。

將空間結(jié)構(gòu)因子加入模型A中，空間結(jié)構(gòu)因子及其組合對應(yīng)的擬合優(yōu)度和模型檢驗結(jié)果見表4。由表4可知，加入空間結(jié)構(gòu)因子后，枝下高模型的擬合優(yōu)度得到了提升，各模型檢驗指標(biāo)減小。其中PV、PH和dDIS對枝下高模型擬合效果的提升依次減小。加入的空間結(jié)構(gòu)參數(shù)越多，對枝下高模型的擬合優(yōu)度提升越大，同時，枝下高模型變得越復(fù)雜。

表4 加入空間結(jié)構(gòu)參數(shù)的枝下高模型擬合優(yōu)度及檢驗結(jié)果

實際上，構(gòu)建的3個空間結(jié)構(gòu)參數(shù)是從不同角度分析對象木受鄰近木的競爭影響，本質(zhì)上都是反映周圍林木與對象木對陽光和生長空間的競爭?？紤]到模型中的變量越少越便于應(yīng)用，因此，選擇影響效果較為明顯的垂直空間結(jié)構(gòu)參數(shù)PV，枝下高模型的R2從0.717提升到了0.741，σ(SEE)從1.407 m減小到了1.321 m,σ(TRE)和σ(MSE)均在±3%范圍內(nèi)。最終建立的枝下高模型如下所示：

HCB=H/[1+exp(0.046 42-0.000 21DBH-0.237 22PV)]。

(21)

根據(jù)枝下高模型的預(yù)測值，繪制殘差分布圖(圖2)，由圖2可知，殘差的散點隨機分布，模型擬合效果較好。

2.2 杉木枝下高模擬

選擇調(diào)查數(shù)據(jù)，利用胸徑連年生長量模型獲取每年的胸徑數(shù)據(jù)，利用樹高曲線模型獲取樹高數(shù)據(jù)，冠幅面積估計模型獲取冠幅數(shù)據(jù)，利用枝下高模型獲取枝下高數(shù)據(jù)。

根據(jù)林木當(dāng)前狀態(tài)，擬合10年后、20年后林木的生長狀態(tài)。為了方便觀察，顯示基于前5條數(shù)據(jù)構(gòu)建的杉木枝下高數(shù)據(jù)，并選取其中編號1～5號的數(shù)據(jù)(表5)。由表5可知，枝下高與樹高相關(guān)性較高，樹高較高的杉木，其枝下高也較高，且杉木受到垂直方向上的遮蓋越嚴(yán)重，其枝下高也越高。

表5 枝下高預(yù)測數(shù)據(jù)

3 討 論

林木除了受自身因素調(diào)控，還受外部環(huán)境影響，而林木屬性(樹高、胸徑、冠幅、枝下高等)則是這些影響因素綜合作用的結(jié)果。樹高和胸徑是描述林木健康狀況的兩個重要屬性，與枝下高息息相關(guān)[35]，一定程度上可以反映枝下高與林木自身因素的關(guān)系。枝下高與外部環(huán)境的關(guān)系則可以通過立地指數(shù)、林分密度等統(tǒng)計因子、林分空間結(jié)構(gòu)參數(shù)等來表征，但本質(zhì)上是反映林木間競爭的優(yōu)劣程度[26]。與類似研究相比[1-5,13,19]，本研究構(gòu)建的枝下高模型參數(shù)少，擬合精度稍低：一是因為樣地較少，不能提供完整的林分統(tǒng)計信息；二是因為在數(shù)據(jù)預(yù)分析階段，發(fā)現(xiàn)與距離無關(guān)的競爭因子等林分因子對枝下高模型并無明顯改善，且林分因子在生物學(xué)上的解釋性較差。實際上，林木主要受周圍幾株鄰近木的競爭影響，考慮林分空間結(jié)構(gòu)單元對林木枝下高的影響更為直接。

本研究的枝下高模型僅與林木自身屬性和林木所處空間結(jié)構(gòu)單元有關(guān)，無需調(diào)查較大范圍內(nèi)的所有林木信息，便于在樣地狀況復(fù)雜、無法全面調(diào)查林木數(shù)據(jù)時，能夠基于較小范圍內(nèi)的空間結(jié)構(gòu)單元，獲得受競爭影響的枝下高數(shù)據(jù)。

在枝下高可視化模擬中，本研究基于一期數(shù)據(jù)預(yù)測了枝下高未來10年和20年后的狀態(tài)，由于沒有真實數(shù)據(jù)進(jìn)行對比，因此無法對模擬結(jié)果的準(zhǔn)確性進(jìn)行檢驗。但枝下高模型反映了杉木林分中枝下高與樹高、胸徑和垂直空間結(jié)構(gòu)參數(shù)的關(guān)系，其中垂直空間結(jié)構(gòu)參數(shù)是基于對象木和相鄰木的樹高得來的，因此在胸徑、樹高模型擬合效果較好的基礎(chǔ)上，枝下高模型的預(yù)測準(zhǔn)確性在可接受的范圍內(nèi)。結(jié)合他人研究的杉木生長模型，實現(xiàn)了杉木枝下高可視化模擬，效果較為逼真，更符合現(xiàn)實。

為了提高杉木枝下高模型的預(yù)測精度，在今后的外業(yè)調(diào)查中可以考慮利用大區(qū)域的多個樣地進(jìn)行研究；由于林木胸徑和樹高具有極為顯著的相關(guān)性，因此聯(lián)立樹高模型和枝下高模型，形成僅與胸徑和林木距離相關(guān)的枝下高模型是今后研究的方向。