臧 顥，黃錦程，劉洪生，歐陽勛志，寧金魁*

（1.江西農(nóng)業(yè)大學 林學院，江西 南昌 330045；2.江西農(nóng)業(yè)大學 鄱陽湖流域森林生態(tài)系統(tǒng)保護與修復國家林業(yè)和草原局重點實驗室，江西 南昌 330045；3.江西省崇義縣林業(yè)局，江西 贛州 341300）

【研究意義】杉木（Cunninghamia lanceolata（Lamb.）Hook.）是我國最重要的用材樹種之一。據(jù)第9 次全國森林資源清查結(jié)果[1]，杉木人工林總面積達990.20 萬hm2、蓄積為7.55 億m3，占全國人工林總面積的17.33%和總蓄積的22.30%。林分斷面積是林業(yè)調(diào)查中的主要因子以及林分生長和收獲預估體系中的重要變量[2-3]，因而，研究杉木林分斷面積生長對預估杉木林的收獲有重要意義?！厩叭搜芯窟M展】為了構(gòu)建林分斷面積模型，一些學者采用多種統(tǒng)計方法進行了研究?；诮魇?1 塊樣地數(shù)據(jù)，杜紀山等[4]采用Richards 和Schumacher 方程作為基礎(chǔ)模型，以年齡、林分密度指標和立地質(zhì)量指標為自變量構(gòu)建了杉木林分斷面積生長模型，并采用最小二乘法估計了模型參數(shù)。倪成才等[5]基于美國佐治亞州148 個火炬松人工林固定樣地的調(diào)查數(shù)據(jù)，采用代數(shù)差分方程構(gòu)建了林分斷面積模型，并采用最大似然和限制性最大似然估計了模型的參數(shù)。Zhao 等[6]收集了江西省、廣東省和廣西省共583 個杉木人工林固定樣地數(shù)據(jù)，采用最小二乘法和3 水平混合效應模型對含優(yōu)勢高和每公頃株數(shù)的林分斷面積生長模型進行了參數(shù)估計。Zhang 等[7]以江西省分宜縣杉木人工林為研究對象，采用貝葉斯模型平均法將林分模型、單木模型和直徑分布模型聯(lián)合起來，并對林分斷面積進行了預測。張雄清等[8]采用廣義代數(shù)差分法構(gòu)建杉木林分斷面積生長模型，并采用最小二乘法和貝葉斯法分別估計了模型參數(shù)。【本研究切入點】盡管最小二乘法、混合效應模型、貝葉斯等方法被廣泛用于構(gòu)建林分斷面積生長模型，但這些統(tǒng)計方法在應用時都存在著一定的假設(shè)前提，諸如參數(shù)的特定分布形式、自變量的獨立性等，且需要預先給定因變量和自變量之間的關(guān)系式[9-11]。而森林的生長是一個復雜的過程，它受到立地條件[12-14]、氣候變化[15-16]、林分結(jié)構(gòu)[17-18]及其交互作用的綜合影響，因而上述假設(shè)通常難以滿足[19]。由于能處理復雜的變量關(guān)系，且不需要預設(shè)模型的形式和輸入數(shù)據(jù)的分布等優(yōu)點[20-21]，近些年來，機器學習被逐漸用于模擬森林生長模型中，而增強回歸樹作為機器學習的一種，其預估精度較高且估計結(jié)果比較穩(wěn)定[22-24]，因而在林業(yè)上有著較多的應用?！緮M解決的關(guān)鍵問題】本研究以江西省崇義縣杉木人工林為研究對象，采用增強回歸樹構(gòu)建杉木人工林林分斷面積生長模型，探討了增強回歸樹算法中不同參數(shù)的變化對模型的影響，并分析了不同林分密度指標和不同立地指標及其交互作用對林分斷面積的模擬影響。

1 研究區(qū)概況

研究區(qū)位于江西省贛州市崇義縣，屬于南嶺的北端。崇義縣的地理坐標為113°55′E～114°38′E，25°24′N～25°55′N。全縣屬亞熱帶季風氣候，年平均氣溫17.9 ℃，極端最高39.9 ℃，最低-8 ℃，年均降水量1 603.1 mm，無霜期307 d。全縣總面積2 206.27 km2，主要地貌為低山、丘陵，土壤類型有黃壤、紅壤等。崇義縣屬常綠闊葉林生物氣候帶，主要樹種有絲栗栲（Castanopsis fargesiiFranch.）、米櫧（Castanopsis carlesii（Hemsl.）Hayata）、杉木（Cunninghamia lanceolata（Lamb.）Hook.）、馬尾松（Pinus massonianaLamb.）、苦楝（Melia azedarachL.）、南酸棗（Choerospondias axillaris（Roxb.）Burtt.et Hill）等。

2 材料與方法

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

數(shù)據(jù)來源于2016 年江西省贛州市崇義縣的森林資源二類調(diào)查數(shù)據(jù)，本研究選擇了起源為人工且優(yōu)勢樹種為杉木的喬木林小班進行分析，共計8 780個小班，主要因子統(tǒng)計信息見表1。

表1 主要因子統(tǒng)計Tab.1 Summary statistics of primary factors

2.2 研究方法

2.2.1 確定自變量 基于前人的研究結(jié)果[25-27]，立地、密度和年齡是林分斷面積的主要影響因子。因此，本研究的自變量也按這3方面進行選擇，其中，立地指標選取了地位指數(shù)和地位級指數(shù)，密度指標選取了每公頃株數(shù)和林分密度指數(shù)，年齡為平均年齡。地位級指數(shù)、地位指數(shù)和林分密度指數(shù)的計算如下：

（1）地位指數(shù)。依據(jù)蔡學林等[28]對贛南山地杉木人工林構(gòu)建的數(shù)量化地位指數(shù)模型和小班的地形因子信息，計算各小班的地位指數(shù)，基準年齡為20 a。

（2）地位級指數(shù)。地位級指數(shù)是基準年齡下的林分平均高，可采用差分方程進行計算[29]。計算步驟為：1）采用理論生長方程構(gòu)建林分平均高模型，如式①；2）假定立地條件的高低僅影響林分平均高所能達到的最大值，即方程中反映理論最大值的參數(shù)為自由參數(shù)，此時，年齡取基準年齡則方程計算得到的就是地位級指數(shù)，如式②所示，進而可采用差分方程對該參數(shù)進行差分，最終形式如式③所示。

式中，Hi和ti分別表示第i個小班的林分平均高和年齡，SCIi表示第i個小班的地位級指數(shù)，t0為基準年齡，本研究中取20年，β0為反映理論最大值的參數(shù)，βj為其他待估參數(shù)，β0i表示隨小班立地條件變化的自由參數(shù)，f（）表示理論方程，g（）表示對應的差分方程。

本研究中選取了常見的理論生長方程[30]，各方程形式如下：

式中，β1、β2為待估參數(shù)，其他符號如上所述。

（3）林分密度指數(shù)。依據(jù)Reineke[31]給出的公式計算各小班的林分密度指數(shù)，計算公式如下：

式中：SDIi、Ni和Dgi為第i個小班的林分密度指數(shù)、每公頃株數(shù)和林分平均胸徑，D0表示標準平均直徑，本研究中取25 cm[32]。

Pretzsch 等[33]通過對歐洲中部的4 種森林類型進行分析，發(fā)現(xiàn)樹木的徑向生長會受到立地和密度的交互作用影響，這一影響隨樹種的變化而變化，因此，本研究在構(gòu)建林分斷面積模型時考慮了不同立地指標和不同密度指標的交互作用，并于不考慮交互作用的模型進行對比，以判斷立地和密度的交互作用對杉木人工林林分斷面積的影響。通過8種自變量組合（如表2所示），分別模擬林分斷面積，構(gòu)成8種模型，以模型的模擬效果對立地指標、密度指標及交互作用的影響。

表2 不同的自變量組合Tab.2 Combinations of different independent variables

2.2.2 增強回歸樹構(gòu)建林分斷面積模型 增強回歸樹是一種基于集成學習的決策樹模型，通過對訓練樣本進行多次不放回隨機抽樣，生成一系列訓練子樣本，并用回歸樹對每一個子樣本進行擬合，擬合過程按順序進行，將當前的模型添加到之前的模型中，該過程一直持續(xù)到指定的迭代次數(shù)為止，最終的預測結(jié)果則依據(jù)各回歸樹擬合效果進行加權(quán)[11,24]。

增強回歸樹建模有5 個輸入?yún)?shù)需要設(shè)定，分別為損失函數(shù)（distribution）、交互深度（interaction.depth）、收縮參數(shù)（shrinkage）、回歸樹數(shù)量（n.trees）、子取樣比例（bag.fraction）。針對回歸問題，損失函數(shù)一般設(shè)為平方誤差。結(jié)合前人的研究結(jié)果[24,34]，子取樣比例取0.5。而交互深度、收縮參數(shù)和回歸樹數(shù)量的確定，不同的研究中無統(tǒng)一的結(jié)果[11,21,24]。為探討交互深度、收縮參數(shù)和回歸樹的數(shù)量的設(shè)置對擬合結(jié)果的影響，本研究中對這3個參數(shù)設(shè)置不同的數(shù)值，通過對比不同數(shù)值組合的模型評價指標，對模型進行優(yōu)化：交互深度取1、2、3、4、5、6、7、8、9、10，收縮參數(shù)取0.1、0.05、0.01、0.005、0.001，回歸樹的數(shù)量取100、200、300、400、500、600、700、800、900、1000。

2.2.3 模型評價 本研究采用5折交叉驗證進行模型評價：1）將數(shù)據(jù)隨機分成5個子樣本；2）每次用其中4個子樣本作為建模樣本，最后1個子樣本作為檢驗樣本，依據(jù)檢驗樣本計算評價指標；3）重復上一步，直到5個子樣本都進行了模型檢驗；4）計算5次模型檢驗得到的評價指標的平均值，并以此判斷模型優(yōu)劣。

采用3種指標評價模型的預估精度：MAB，平均絕對誤差；RMA，平均相對誤差絕對值；RMSE，均方根誤差，計算公式如下[21,35-36]：

式中，Gi和分別表示第i個小班的林分斷面積實測值和預測值，n為樣本數(shù)。

本研究中的所有計算均在R 軟件中完成，其中，增強回歸樹調(diào)用“gbm”包中的“gbm”函數(shù)實現(xiàn)的，增強回歸樹的參數(shù)優(yōu)化是在“caret”包中的“train”函數(shù)實現(xiàn)的。

3 結(jié)果與分析

3.1 地位級指數(shù)的計算

由表3可知，7種理論方程擬合的效果差異不大，其中，式⑩擬合的林分平均高效果最好。相較其他6 種理論方程，式⑩的MAB減少了0.06%～7.22%，RMSE減少了0.01%～3.91%，因此選擇式⑩基礎(chǔ)模型構(gòu)建差分方程，并計算地位級指數(shù)。計算地位級指數(shù)的差分形式如下：

式中各符號如前所述。

3.2 增強回歸樹參數(shù)優(yōu)化

圖1 顯示了8 種林分斷面積模型的參數(shù)優(yōu)化情況?？偟膩砜?，隨著回歸樹數(shù)量的增加，8 種模型的MAB、RMA和RMSE均呈現(xiàn)減小的趨勢。這一趨勢在不同交互深度下沒有明顯差異，但隨著收縮參數(shù)增加，這一趨勢逐漸減弱，當收縮參數(shù)取0.05 和0.1 時，隨著回歸樹數(shù)量的增加，MAB、RMA和RMSE的變化無明顯趨勢。

綜合考慮MAB、RMA和RMSE，8 種模型的最優(yōu)參數(shù)值及其交叉驗證的評價指標見表4。可知：含有林分密度指數(shù)的模型明顯優(yōu)于含有每公頃株數(shù)的模型，含有林分密度指數(shù)的模型的MAB減少了80.05%～84.15%，RMA減少了80.42%～84.66%，RMSE減少了77.94%～82.27%；含有地位級指數(shù)的模型優(yōu)于含有地位指數(shù)的模型，含有地位級指數(shù)的模型的MAB減少了11.32%～25.84%，RMA減少了11.68%～22.96%，RMSE 減少了6.95%～25.24%；含有立地和密度交互作用的模型略優(yōu)于不含交互作用的模型，含交互作用的模型的MAB減少了4.56%～12.48%，RMA減少了5.88%～15.48%，RMSE減少了1.07%～9.34%，僅M3（地位指數(shù)與每公頃株數(shù)的無交互模型）優(yōu)于M7（地位指數(shù)與每公頃株數(shù)的有交互模型）。M6（考慮了地位級指數(shù)、林分密度指數(shù)及其交互作用）的各項評價指標均為最小，可用于預估崇義縣杉木人工林林分斷面積。

表3 基于交叉驗證的林分平均高模型評價Tab.3 The evaluation of stand mean height model using cross-validation

圖1 8種林分斷面積模型的增強回歸樹參數(shù)優(yōu)化Fig.1 Tuning boosted regression tree for 8 stand basal area models

表4 8種林分斷面積模型的參數(shù)優(yōu)化結(jié)果Tab.4 The tuning results for 8 stand basal area models

4 結(jié)論與討論

本研究以江西省崇義縣的杉木人工林為對象，采用增強回歸樹構(gòu)建了林分斷面積模型，通過交叉驗證分析了不同立地指標、不同密度指標及其交互作用對模擬的影響。研究發(fā)現(xiàn)，經(jīng)過參數(shù)優(yōu)化的增強回歸樹能較好的擬合林分斷面積，最優(yōu)模型的MAB為0.504 7 m2/hm2，RMA為2.29%，RMSE為0.900 8 m2/hm2。

地位指數(shù)是基準年齡時的優(yōu)勢高[37]，被廣泛應用于林分斷面積模型中以反映立地質(zhì)量對林分斷面積的影響[4,9,27]。因為連續(xù)調(diào)查固定樣地數(shù)據(jù)通常缺乏林分優(yōu)勢高數(shù)據(jù)，只有林分平均高數(shù)據(jù)[38-39]，因此有學者使用基于林分平均高計算的地位級指數(shù)評價立地質(zhì)量[29]。同時，也有學者研究發(fā)現(xiàn)林分平均高和優(yōu)勢高之間存在著相當穩(wěn)定的數(shù)量關(guān)系[40-41]，因而使用基于林分平均高的地位級指數(shù)反映立地質(zhì)量是一種可行的方法。本研究對比了地位指數(shù)和地位級指數(shù)對林分斷面積模擬的效果，發(fā)現(xiàn)地位級指數(shù)對林分斷面積的模擬效果良好，這說明用基于林分平均高的地位級指數(shù)能有效的反映立地質(zhì)量對林分斷面積的影響。此外，盡管本研究中，地位級指數(shù)對林分斷面積的模擬優(yōu)于地位指數(shù)，這可能是因為本研究是采用前人構(gòu)建的模型推算得到的地位指數(shù)，即存在著一定的誤差，進而影響了地位指數(shù)對林分斷面積的模擬效果。

每公頃株數(shù)和林分密度指數(shù)是被廣泛應用于森林生長收獲模型中[17,42-43]，本研究對比了每公頃株數(shù)和林分密度指數(shù)對杉木人工林林分斷面積的模擬效果，結(jié)果顯示林分密度指數(shù)對林分斷面積的模擬效果明顯優(yōu)于每公頃株數(shù)，這與一些學者的研究結(jié)果相同[9,27]。盡管每公頃株數(shù)描述了單位面積上的林木株數(shù)，但卻沒有包含樹木大小信息，隨著樹木的徑向生長，可能存在每公頃株數(shù)較少但由于樹木胸徑較大，導致樹木個體之間競爭更強的情況，這在一定程度上也會影響到林分斷面積的生長。林分密度指數(shù)則是將當前的株數(shù)轉(zhuǎn)換到標準胸徑下的株數(shù)，即考慮了林分平均胸徑的影響，這可能是林分密度指數(shù)的模擬效果優(yōu)于每公頃株數(shù)的原因。

為了探討立地和競爭是否存在對林分斷面積的交互作用，本研究對比了含交互作用項與不含交互作用項的林分斷面積模型的模擬效果，結(jié)果顯示，除含地位指數(shù)和每公頃株數(shù)的交互作用項的模型外，其余含交互作用項的林分斷面積模型擬合精度均略高于不含交互作用項的模型，這說明了立地和林分密度的交互作用對杉木人工林林分斷面積有一定影響。這與Pretzsch 等[33]的研究發(fā)現(xiàn)類似。因此，如果需要精確模擬森林的徑向生長時，有必要考慮立地和林分密度的交互作用。

由于機器學習算法對模型的形式和輸入數(shù)據(jù)的分布沒有假設(shè)前提，且能處理復雜的變量關(guān)系、擬合效果好[11,20]，所以逐漸應用于森林生長收獲預估中。但也有學者認為機器學習存在“黑箱”問題，無法給出因變量和自變量的關(guān)系式，但其實每種機器學習算法存在各自的工作原理，而因變量和自變量的關(guān)系也可以通過模擬的方式給出，偏依賴圖也可以給出各自變量對因變量的重要性[21]。相對傳統(tǒng)的統(tǒng)計方法更關(guān)注于模型的基本形式（線性、非線性等）和變量的關(guān)系（是否存在共線性、自變量參數(shù)是否顯著）等問題，采用機器學習方法構(gòu)建森林生長模型時需要更多的考慮參數(shù)值的確定和過擬合問題，參數(shù)值的確定可以通過對比不同參數(shù)取值組合的擬合效果來確定，但過于追求擬合效果也可能會導致過擬合問題，即對建模樣本較好，但在檢驗樣本表現(xiàn)較差的情況[11]，因此建議使用機器學習構(gòu)建模型時需要進行模型檢驗。

綜上所述，增強回歸樹作為一種機器學習算法，可應用于林分斷面積模型的構(gòu)建，不同的交互深度、收縮參數(shù)及回歸樹數(shù)量對擬合精度影響較大，有必要進行參數(shù)優(yōu)化，同時不同的自變量組合對擬合精度也有較大的影響，基于林分平均高的地位級指數(shù)和林分密度指數(shù)及其交互作用可作為模型自變量以提高模型的預測效果。