基于最粗優勢木胸徑生長的湖南櫟類天然林立地質量評價模型*

何 靜 李新建 朱晉梅 朱光玉

(1.中南林業科技大學林學院 長沙 410004； 2.國家林業和草原局中南調查規劃院 長沙 410004； 3.湖南省林業局 長沙 410001)

立地生產力指一個特定林分生產地上林木的潛力，是林地和林分共同形成的生產單位(Skovsgaardetal.， 2008)。立地生產力評價是判斷或預測立地適宜性或潛在生產力的一項基礎性工作，在森林科學經營中具有重要作用(Berrilletal.， 2013)。森林立地生產力評價指標可歸納為林分高、林分蓄積或斷面積生長量、林分特征和立地屬性4類(朱光玉等， 2016)，基于4類評價指標產生了多種評價方法，其中以地位指數評價最為常用，即以某樹種林分基準年齡時的優勢平均高作為衡量森林立地生產潛力的指標(Fonwebanetal.， 1995； 孟憲宇， 2006)。

野外實地調查時，尤其在林區中，由于地形十分復雜、坡度變化較大、林木生長狀況和林分密度不同等因素，給樹高的準確測定帶來了很多困難，而林木胸徑與和樹高之間存在密切關系(Sharmaetal.， 2007； Leietal.， 2009； 王冬至等， 2016； 2019； 段光爽等， 2020)，且胸徑測量非常準確，因此，一些學者以優勢木胸徑為指標，通過構建模型或采用地位級法對林分生長與立地質量評價進行研究： 如應用多期胸徑-年齡數據編制胸徑立地指數表，用于天然林立地質量評價(劉洵等， 2019)； 應用樹木理論生長方程，以優勢木平均胸徑為指標，導出立地指數模型，從理論和技術上解決立地質量評價問題(李鐵華等， 1997)； 以最大胸徑生長率為指標建立的杉木(Cunninghamialanceolata)人工林立地質量評價數量化理論模型，不僅可以評估現有杉木人工林的生產力，而且能夠預測非林地的潛在生產力(楊海賓等， 2020)。可見，應用胸徑-年齡構建立地指數模型評價立地質量是可行的(李振芳等， 2016； 吳恒等， 2015； 馬煒等， 2013)。

混合效應模型是一種常用的處理局部變異性的建模方法，可在模型中反映其他因子的作用，針對不同立地類型，混合效應模型能夠根據混合效應參數構建優勢木胸徑生長曲線簇，從而反映立地類型對胸徑生長的影響程度(Fangetal.， 2001； Calamaetal.， 2004)。而且，混合效應模型既可以反映個體之間的差異，也可以反映總體的平均變化趨勢，利用混合效應構建多行立地指數曲線能夠提高立地指數估計精度(Calegarioetal.， 2005； 符利勇等， 2013； Zhuetal.， 2019)。

櫟類(Quercus)天然林是我國天然林的主要林分類型，作為一種闊葉用材林在湖南林業生產中發揮了重要作用(胡芳名等， 2000)。天然林經營管理相對復雜，天然混交林立地質量評價仍是當今世界林業研究的難點和焦點，對天然林生產力預估和評價研究非常必要。鑒于此，本研究以湖南櫟類天然林為研究對象，采用數量化方法Ⅰ篩選對優勢木胸徑生長影響顯著的立地因子，通過混合效應、k-means聚類構建含立地類型的優勢木胸徑生長模型，導出以最粗優勢木胸徑為指標的基于立地分級的湖南櫟類天然林立地質量評價模型，以期為天然林立地質量評價提供一種新方法。

1 研究區概況

湖南位于長江以南，屬長江中游地區，地處108°47′—114°15′E，24°38′—30°08′N，為大陸性亞熱帶季風濕潤氣候，年均氣溫15～18℃，年降水量1 200～1 700mm。土壤以黃壤和黃棕壤為主。植物資源豐富，包括青岡櫟(Cyclobalanopsisglauca)、栲(Castanopsisfargesii)、銀杏(Ginkgobiloba)、南方紅豆杉(Taxuschinensis)、水杉(Metasequoiaglyptostroboides)等珍貴樹種。

2 數據與方法

2.1 數據采集

在湖南五蓋山林場、八大公山自然保護區、平江縣蘆頭實驗林場等5個地點，設置20 m×30 m櫟類天然林固定樣地51塊，對樣地內胸徑大于5cm的活力木進行每木檢尺，調查因子包括胸徑、樹高、海拔、坡度、坡位、坡向、土壤類型等。對樣地數據進行預處理，選取每塊樣地的最粗優勢木，并利用生長錐鉆取年齡(AGE)。林分因子統計見表1。

表1 林分因子統計

2.2 研究方法

2.2.1 數量化方法 Ⅰ篩選因子 數量化方法Ⅰ指自變量中包含定性變量的回歸模型，表達式(周光亞等， 1979)如下：

(1)

式中：yi為變量y在i樣品中的觀測值；bjk為j項目k類目的得分；εi為隨機誤差。

運用Forstat 2.2“統計之林”軟件進行數量化方法Ⅰ分析，篩選可能影響櫟類天然林優勢木胸徑生長的因子。

2.2.2 基礎模型擬合與選擇 選取4種常見的樹木理論生長方程，分別對櫟類天然林最粗優勢木胸徑與年齡的相關關系進行擬合，候選基礎模型名稱、表達式和參數見表2，通過一系列評價指標，確定最符合櫟類天然林最粗優勢木胸徑生長的最優基礎模型。

2.2.3 立地分類 為方便劃分立地類型和建立模型，將海拔、坡度、坡向、坡位和土壤類型5種立地因子進行分級處理，分級標準見表3。

2.2.4k-means聚類 應用k-means聚類劃分立地類型。k-means算法實質上是一種將聚類視為密度估計問題的概率方法，假設樣本來自于如下形式的混合模型(雷小鋒等， 2008)：

p(x|θ)=∑j=1∶kP(Cj)pj(x|θj,Cj)。

(2)

式中：θ=(θ1,…,θk)為待估計的參數向量； 條件概率密度pj(x|θj,Cj)稱為分量密度，表示類別j的概率密度形式，且參數向量θj未知； 先驗概率P(Cj)稱為混合因子。

表2 候選基礎模型擬合結果①

①D為最粗優勢木胸徑(cm)，AGE為最粗優勢木年齡(年)，a、b、c為模型待求參數。Dis the diameter at breast height (cm) of the thickest dominant tree, AGE is the age (year) of the thickest dominant tree,a,bandcare the parameters to be calculated in the model.

表3 調查因子等級劃分

2.2.5 混合效應模型構建 混合效應模型(mixed effect models， MEMs)是根據回歸函數依賴于效應參數的非線性關系建立的數學模型，其中，效應參數包括固定效應參數和隨機效應參數，其一般形式(符利勇等， 2020)如下：

Yi=f(β,μi,Xi)+εi。

(3)

式中：Yi和Xi分別為第i個立地類型的因變量和自變量向量；εi為誤差向量；β、μi分別為固定效應參數和隨機效應參數。

2.2.6 模型評價與檢驗 采用赤池信息量(Akaike information criterion, AIC)、貝葉斯信息量(Bayesian information criterion, BIC)、確定系數(determination coefficient,R2)、平均絕對誤差(mean absolute error, MAE)和均方根誤差(root mean square error, RMSE)評價模型，指標計算公式如下:

AIC=-2lnl+2p;

(4)

BIC=-2lnl+lnn×p;

(5)

(6)

(7)

(8)

3 結果與分析

3.1 立地因子篩選

林分最粗優勢木胸徑生長受多種立地因子影響，根據最粗優勢木胸徑與不同立地因子的關系，采用數量化方法Ⅰ對其進行顯著性分析，結果(表4)表明，對優勢木胸徑生長影響顯著的立地因子包括海拔、坡度、坡向和坡位，顯著性順序為海拔>坡度>坡向>坡位。

3.2 基礎模型擬合與選擇

基于櫟類天然林51塊樣地實測數據，分別以最粗優勢木胸徑和年齡為因變量和自變量，構建林分最粗優勢木胸徑生長基礎模型，采用表2中的基本理論生長方程通過R語言的nls函數擬合最粗優勢木胸徑和年齡的相關關系。由表5可知，4種理論生長方程的確定系數(R2)均在0.7左右，其中Richards模型的擬合效果最好，R2為0.7318，平均絕對誤差(MAE)和均方根誤差(RMSE)分別為5.4426和6.8791，故本研究選擇擬合效果最好的Richards模型為最粗優勢木胸徑生長基礎模型(M4)。

表4 立地因子顯著性檢驗①

表5 基礎模型擬合結果

3.3 基于立地類型效應的非線性混合效應模型構建與分析

不同立地條件下，林分最粗優勢木胸徑生長可能存在差異，故將海拔、坡度、坡向和坡位按標準(表3)分級、組合構成720種初始立地類型，具體劃分情況見表6。

表6 初始立地類型與立地類型組合劃分

本研究只用到初始立地類型中的38種，將組合后的初始立地類型作為隨機效應加入基礎模型構建混合效應模型，模型擬合結果(表7)顯示，只在參數a上收斂，含立地類型隨機效應的櫟類天然林立地指數模型表達式為：

D=[a+μi)×(1-exp(-c×(AGE)]^b+ε。

(M4.1)

式中：a、b、c為模型參數；μi為隨機效應參數；D為林分最粗優勢木胸徑； AGE為林分最粗優勢木年齡；ε為誤差項。

混合效應模型的R2相比基礎模型有所提高，由原來的0.7318升至0.9016，含立地類型的混合效應模型明顯優于基礎模型。

殘差分布是模型異質性檢驗最直觀的工具之一，對基礎模型M4和混合效應模型M4.1的擬合結果進行殘差分析(圖1)可知，模型殘差均呈隨機分布，且M4.1殘差分布更為集中、均勻，可見混合效應能夠提高模型精度。

表7 混合效應模型擬合結果

圖1 基礎模型M4和混合效應模型M4.1殘差分布

3.4 基于立地類型聚類的非線性混合效應模型構建與分析

不同初始立地類型下，最粗優勢木胸徑生長狀況不同，生長過程存在差異。本研究51塊樣地共有38種初始立地類型，為簡化初始立地類型數量，將模型M4.1擬合得到的a參數值應用k-means聚類，在滿足生長規律的基礎上，以聚類精度≥95%為標準，將影響效果相近或相同的初始立地類型聚類成立地類型組。由表8可知，當聚類成6個立地類型組時，精度達96.4%，滿足精度要求。

表8 立地類型聚類結果

將聚類后的立地類型組作為隨機效應加入最優基礎模型構建含立地類型組的混合效應模型(M4.2),由表9可知，與基礎模型M4相比，M4.2的AIC、BIC均有所減小，Loglik有所上升；R2從0.7318升至0.9269，擬合精度提高26.7%，說明聚類后的立地類型組對櫟類天然林最粗優勢木胸徑生長影響顯著。模型M4.2的表達式為：

Dj=aj×[1-exp(-c×AGE)]^b+ε。

(M4.2)

式中：Dj為第j個立地類型組的最粗優勢木胸徑；aj為第j個立地類型組的參數值；b、c為模型參數；ε為誤差。

表9 模型M4.2的精度評價與參數估計

對模型M4.2的擬合結果進行殘差分析(圖2)可知，模型殘差在橫坐標上下浮動，與圖1對比，分布更加集中、均勻，說明立地類型組作為隨機效應加入模型能夠有效改進殘差，提高模型精度。

圖2 混合效應模型M4.2殘差分布

利用模型M4.2繪制基于立地類型組的胸徑生長曲線簇，由圖3可知，曲線擬合效果較好，進一步證實模型的精度和適用性。

圖3 基于立地類型組的胸徑生長曲線簇

3.5 基于混合效應的立地質量評價模型

模型M4.2兩邊同時取對數：

lnDj=b×lnaj×[1-exp(-c×AGE)]，

(9)

得

(10)

因此，

SQEIM-DBH=aj×[1-exp(-c×AGE0)]^b+ε。

(11)

式中： SQEIM-DBH為基于最粗優勢木胸徑生長的立地質量評價指標；AGE0為基準年齡。

整理式(11)導出以最粗優勢木胸徑為指標的基于立地分級的立地質量評價模型為：

SQEIM-DBH=aj×[1-exp(-c×AGE0)]^

(12)

代入參數值得到最終模型為：

SQEIM-DBH=aj×[1-exp(-0.03×AGE0)]^

(13)

根據模型(13)，當給定基準年齡(AGE0)、優勢木年齡(AGE)和胸徑(D)時，在確定立地類型組下得到a參數值，即可得到立地質量評價指標SQEIM-DBH。

3.6 立地質量評價模型適用性檢驗

林分斷面積是反映立地質量的常用指標，以林分每公頃斷面積為因變量，林分密度指數、年齡、立地質量評價指標為自變量，采用方差分析驗證立地質量評價模型的可行性，結果見表10。

表10 顯著性檢驗

由表10可知，立地質量評價指標與林分斷面積顯著相關，以最粗優勢木胸徑為指標評價櫟類天然林立地質量、預估林地生產力在理論上具有可行性。

4 討論

立地質量是森林生產力的基礎，對森林更新、樹種選擇、地力維持和經營管理至關重要。在異齡林或混交林中，林分密度不影響優勢木的樹高-胸徑關系(李清順等， 2021； Vanclayetal.， 1988； Huangetal.， 1993)，通常利用優勢木樹高-年齡(立地指數法)、樹高-胸徑(立地形法)評價立地質量，其中立地指數法應用最為廣泛。但野外實地調查時，胸徑測量相比樹高更準確，且本研究通過相關分析也得到林分最粗優勢木胸徑生長與株數密度無顯著相關性，而與立地因子顯著相關，因此可以根據優勢木胸徑生長評價林地立地狀況。本研究表明，影響湖南櫟類天然林優勢木胸徑生長的因子包括海拔、坡度、坡位和坡向，其中海拔的影響最顯著，坡度、坡位、坡向的影響次之。海拔主要影響氣溫，坡位、坡度主要影響土壤水分和養分再分配，坡向主要影響光照條件，再結合降雨、氣溫等氣象條件，不同立地條件下的樹木生長會存在差異(鄧秀秀等， 2016； 高若楠等， 2019)。

混合效應模型廣泛用于林分生長建模中，不僅可以提高模型精度，還能夠提升模型適用性。齊戰濤等(2021)基于混合效應模型構建含氣候效應的杉木斷面積模型，為森林經營者在氣候影響下采取合理的經營措施提供了理論參考。本研究構建含立地類型的湖南櫟類天然林最粗優勢木胸徑生長模型，通過分析隨機效應參數發現， 立地類型影響最粗優勢木胸徑生長最大值。由殘差分布可知，將初始立地類型作為隨機效應加入模型構建混合效應模型的殘差分布明顯更加集中、均勻，說明混合效應模型精度更優，能夠更好擬合湖南櫟類天然林最粗優勢木胸徑生長。不同立地類型對優勢木生長的影響差異可通過立地類型混合效應參數值體現出來。為簡化立地類型數量，應用k-means聚類將影響效果相近或相同的初始立地類型聚類成6個立地類型組，并作為隨機效應加入最優基礎模型構建非線性混合效應模型，R2為0.9269，相比基礎模型提高26.7%，相比含初始立地類型的混合模型提高2.8%，說明含立地類型組的非線性混合效應模型能夠更加精確擬合湖南櫟類天然林最粗優勢木胸徑生長狀況。6個立地類型組對應6組參數值，根據不同參數值a可得到不同立地類型組下優勢木胸徑生長狀況，就本研究結果來看，立地類型組1的最粗優勢木胸徑生長最佳，立地類型組6的最粗優勢木胸徑生長最差，因此，基于立地類型的湖南櫟類天然林最粗優勢木胸徑生長模型可為經營者在不同立地條件下的經營管理提供參考。

5 結論