王帥玲，龍時(shí)勝，曾思齊，楊盛揚(yáng)

（中南林業(yè)科技大學(xué) 林學(xué)院，湖南 長沙 410004）

櫟類Quercus sp.在我國森林資源中占有重要地位，也是我國天然林的主要樹種，分布范圍廣，種類豐富。據(jù)第八次全國森林資源清查結(jié)果顯示，櫟林面積占全國森林總面積的10.15%、占天然林總面積的13.7%[1]。櫟林為人類提供了涵養(yǎng)水源、保持水土等生態(tài)價(jià)值和優(yōu)質(zhì)用材、培養(yǎng)食用菌等經(jīng)濟(jì)價(jià)值，且在湖南省分布廣、面積大、資源豐富。因大多數(shù)櫟林處于自生自滅的狀態(tài)，存在林分經(jīng)營不合理等問題，需要建立林分?jǐn)嗝娣e生長模型，對(duì)在自然條件下的林分生長動(dòng)態(tài)進(jìn)行模擬，探究櫟林生長規(guī)律，為其合理、科學(xué)經(jīng)營提供理論依據(jù)。

林分?jǐn)嗝娣e生長模型之所以在研究森林生長規(guī)律中起著重要的基礎(chǔ)作用，是因?yàn)榱址謹(jǐn)嗝娣e是主要的林分測(cè)樹因子之一，具有易測(cè)定性、較高的穩(wěn)定性及預(yù)估性，常被用于林分生長收獲預(yù)估及森林資源動(dòng)態(tài)估測(cè)[2]。國內(nèi)外對(duì)斷面積生長模型進(jìn)行了大量研究，主要可分為單木斷面積生長模型和林分?jǐn)嗝娣e生長模型兩大類[3-5]，林分?jǐn)嗝娣e生長模型常用Richards 和Schumacher 兩種模型形式，構(gòu)建含有年齡、立地質(zhì)量和密度3 個(gè)變量的非線性回歸模型或非線性混合效應(yīng)模型[6-8]。由于傳統(tǒng)模型未考慮立地條件等因子的差異對(duì)林木生長的影響，致使模型預(yù)估精度不高[9]。因此，一些學(xué)者通過引入啞變量或隨機(jī)效應(yīng)構(gòu)建混合效應(yīng)模型對(duì)傳統(tǒng)模型進(jìn)行改進(jìn)，研究表明，啞變量模型和混合效應(yīng)模型的預(yù)估精度優(yōu)于傳統(tǒng)模型[10-12]。目前，對(duì)櫟類林分?jǐn)嗝娣e生長模型的研究多以立地類型、林分類型為隨機(jī)效應(yīng)或啞變量[3,13]，還未出現(xiàn)以林分密度為隨機(jī)變量的研究。林分密度對(duì)林分生長量具有重要影響，因此本研究以林分密度為隨機(jī)效應(yīng)構(gòu)建林分?jǐn)嗝娣e生長混合效應(yīng)模型，并與以林分密度為連續(xù)變量構(gòu)建的再參數(shù)化模型進(jìn)行比較，選出預(yù)估精度最高的模型，以期為櫟類次生林的經(jīng)營措施和提質(zhì)增量提供理論依據(jù)。

1 研究區(qū)概況

湖南省的經(jīng)緯度位置為108°47′～114°15′E，24°38′～30°08′N，屬于長江中游地區(qū)。全境東西寬667 km，南北長774 km，是云貴高原向江南丘陵及南嶺山地向江漢平原過渡地帶，海拔在24～2 099 m 之間。境內(nèi)地形復(fù)雜多樣，主要為丘陵和山地地貌，土壤以紅壤和黃壤為主。屬于大陸性亞熱帶季風(fēng)濕潤氣候，境內(nèi)光熱水資源充足，全年日照時(shí)數(shù)為1 300～1 800 h，年平均氣溫16～19℃，無霜期天數(shù)在253～311 d 之間，年均降水量為1 200～1 700 mm，氣候年際、年內(nèi)變化較大。研究區(qū)內(nèi)植物資源豐富，主要喬木樹種有馬尾松Pinus massoniana Lamb、青岡櫟Cyclobalanopsis glauca、甜櫧Castanopsis eyrei、樟樹Cimamomum camphora等樹種，也有銀杏Ginkgo biloba L、珙桐Davidia involucrata Baill等珍貴樹種；灌木和草本植物有山茶Camellia japonica、鹿角杜鵑Rhododendron latoucheae、麥冬Ophiopogon japonicus、蘭花Cymbidiu-m等。

2 材料與方法

2.1 數(shù)據(jù)來源與整理

本研究所用數(shù)據(jù)源于6 期湖南省森林資源連續(xù)清查數(shù)據(jù)（1989—2014年），篩選數(shù)據(jù)遵從以下原則：選擇櫟類為優(yōu)勢(shì)樹種或主要樹種的樣地；相鄰兩期株數(shù)變動(dòng)控制在20%以內(nèi)；只考慮林分?jǐn)嗝娣e正增長的數(shù)據(jù)。由于國家森林資源連續(xù)清查數(shù)據(jù)中包括林木胸徑、樹高、株數(shù)等測(cè)樹因子。因此需要通過對(duì)選用的數(shù)據(jù)進(jìn)行處理計(jì)算得到林分年齡和林分?jǐn)嗝娣e，因林木年齡和林木胸徑具有相關(guān)性，首先利用固定樣地分樹種構(gòu)建胸徑與年齡的生長方程，之后利用各樹種平均胸徑預(yù)估各樹種的平均年齡，最后用斷面積加權(quán)法算出林分平均年齡[14]。將林分年齡劃分齡階，用3 倍標(biāo)準(zhǔn)差法剔除出林分?jǐn)嗝娣e異常數(shù)據(jù)。最終在125塊樣地中隨機(jī)選出80 塊樣地130 組數(shù)據(jù)作為建模樣本構(gòu)建林分?jǐn)嗝娣e生長模型，剩余45 塊樣地65組數(shù)據(jù)作為檢驗(yàn)樣本。相關(guān)統(tǒng)計(jì)量見表1。

表1 建模數(shù)據(jù)與檢驗(yàn)數(shù)據(jù)統(tǒng)計(jì)Table 1 Statistics of modeling data and inspection data

2.2 基礎(chǔ)模型的選取

生長方程有經(jīng)驗(yàn)方程和理論方程兩類，由于理論生長方程是基于生物學(xué)原理做出假設(shè)模擬生物體大小或者生長的微分方程，具有邏輯性強(qiáng)、參數(shù)可做生物學(xué)解釋的特點(diǎn)。因此，本研究選取常用的7 種方程為備選基礎(chǔ)模型[15]，以林分?jǐn)嗝娣e（G）為因變量，林分平均年齡（A）為自變量，模型表達(dá)式見表2。比較基礎(chǔ)模型擬合結(jié)果，選出最優(yōu)基礎(chǔ)模型，構(gòu)建含有林分年齡、林分密度及立地質(zhì)量3 個(gè)因子為變量的林分?jǐn)嗝娣e生長模型。

表2 基礎(chǔ)模型表達(dá)式?Table 2 The expression of basic model

2.3 混合效應(yīng)模型與參數(shù)化模型構(gòu)建

林分?jǐn)嗝娣e的生長不僅受到林分年齡的影響，還受到立地質(zhì)量和林分密度的影響[10]。在模型中，立地質(zhì)量常用地位指數(shù)或林分優(yōu)勢(shì)木平均高來表示[3,6]，因林分優(yōu)勢(shì)木平均高數(shù)據(jù)不全，所以將地位指數(shù)（SI）[16]作為反映立地質(zhì)量的因子加入到基礎(chǔ)模型。在立地條件相同時(shí)，林分密度不同，林分?jǐn)嗝娣e生長狀態(tài)不同[8]。林分調(diào)查因子中包含郁閉度和株數(shù)兩種直觀反映林分密度的指標(biāo)，因郁閉度數(shù)據(jù)不完整，故用株數(shù)作為隨機(jī)變量加入到基礎(chǔ)模型中構(gòu)建混合效應(yīng)模型?；旌闲?yīng)模型是根據(jù)回歸函數(shù)基于固定效應(yīng)參數(shù)和隨機(jī)效應(yīng)參數(shù)的非線性關(guān)系建立的[17]，所以構(gòu)建混合效應(yīng)模型要先確定隨機(jī)效應(yīng)與固定效應(yīng)的組合形式及數(shù)量，在保證預(yù)估精度的同時(shí)，還應(yīng)避免參數(shù)過多化，以免造成模型不收斂[18]。將株數(shù)進(jìn)行定性化處理，利用定性代碼表示，并將定性數(shù)據(jù)轉(zhuǎn)化為（0，1），根據(jù)樣地林分特征，將株數(shù)劃為6 個(gè)密度等級(jí)，相鄰的兩個(gè)等級(jí)間隔為400 株（表3），定性因子和混合效應(yīng)模型表達(dá)式為：

表3 林分密度等級(jí)劃分Table 3 The classification of stand density

式（1）中，x表示不同的林分密度等級(jí)，i=1,2,3,4,5,6。

式（2）中：Yij、Xij分別為i樣地第j次觀測(cè)的因變量和自變量，εij表示殘差向量，f為擬合的基礎(chǔ)模型，β和ui分別為固定和隨機(jī)參數(shù)向量，aij、bij分別為固定和隨機(jī)效應(yīng)的設(shè)計(jì)矩陣設(shè)計(jì)矩陣，m為樣地?cái)?shù)，ni表示第i塊樣地的觀測(cè)次數(shù)，σ2、ψ分別為方差協(xié)方差矩陣和隨機(jī)效應(yīng)協(xié)方差矩陣。

根據(jù)以林分株數(shù)為隨機(jī)變量的混合效應(yīng)模型形式，構(gòu)建以林分株數(shù)為連續(xù)變量的林分?jǐn)嗝娣e生長參數(shù)化模型，對(duì)混合效應(yīng)模型和參數(shù)化模型進(jìn)行適用性比較，選用預(yù)估精度最高的模型作為櫟類次生林林分?jǐn)嗝娣e生長模型。

2.4 模型精度比較及適用性檢驗(yàn)方法

對(duì)備選模型擬合結(jié)果進(jìn)行評(píng)價(jià)后選取最佳基礎(chǔ)模型，用確定系數(shù)R2和預(yù)估精度P值及均方根誤差（RMSE）3個(gè)評(píng)價(jià)指標(biāo)對(duì)模型擬合結(jié)果進(jìn)行評(píng)價(jià)，R2和P越大且RMSE 越小越好；對(duì)具有相同參數(shù)個(gè)數(shù)的混合效應(yīng)模型選用赤池信息準(zhǔn)則（AIC）、貝葉斯信息準(zhǔn)則（BIC）及對(duì)數(shù)似然值（LL）進(jìn)行評(píng)價(jià)，其中AIC 和BIC 的值越小、LL 的值越大，模型的擬合效果越好，具有不同參數(shù)個(gè)數(shù)的模型還應(yīng)進(jìn)行似然比（LRT）檢驗(yàn)（ρ＜0.000 1，差異顯著）。用未參與建模的45 塊樣地65 組數(shù)據(jù)對(duì)所構(gòu)建的3 種模型進(jìn)行適用性檢驗(yàn)，驗(yàn)證預(yù)測(cè)值和實(shí)際值是否接近，模型是否有實(shí)際使用效果，檢驗(yàn)指標(biāo)選用確定系數(shù)R2、預(yù)估精度P值、平均誤差（ME）、均方根誤差（RMSE）。表達(dá)式如下：

式（3）～（9）中：n為樣本數(shù)，yi為斷面積實(shí)測(cè)值；為斷面積預(yù)測(cè)值；為斷面積平均預(yù)估值，n為樣本數(shù)，ρ為模型中參數(shù)個(gè)數(shù)，SSE 為殘差平方和，ln 為自然對(duì)數(shù)，置信水平α=0.05。

3 結(jié)果與分析

3.1 基礎(chǔ)模型構(gòu)建

以篩選出的80 塊樣地130 組樣本為基礎(chǔ)，分別用備選生長理論方程建立基礎(chǔ)模型，模型擬合結(jié)果見表4。從表中可以看出7 個(gè)模型的確定系數(shù)R2在0.297～0.314 之間，預(yù)估精度P值在92.76%～92.89%，RMSE 的范圍為6.211～6.286。模型M3 的確定系數(shù)R2（0.314）和預(yù)估精度P值（92.89%）均最大，均方根誤差（RMSE=6.211）最小，因此將Richards 作為斷面積生長基礎(chǔ)模型。

表4 基礎(chǔ)模型擬合結(jié)果Table 4 The fitting result of basic model

3.2 混合效應(yīng)模型與參數(shù)化模型構(gòu)建

在構(gòu)建混合效應(yīng)模型前，應(yīng)先確定地位指數(shù)（SI）的位置，其中參數(shù)a代表斷面積生長的最大值，b為生長速率參數(shù)，c為形狀參數(shù)，參數(shù)a和c與立地質(zhì)量有關(guān)[3,19]，因此將SI 分別放到模型M3 中不同參數(shù)上，從收斂的組合中選擇擬合效果最好的模型，從表5中可以看出將SI 加到參數(shù)a上要優(yōu)于加到參數(shù)c上，最終將地位指數(shù)加到參數(shù)a上。其表達(dá)式為：

表5 加立地指數(shù)擬合結(jié)果Table 5 The fitting result of add site index

式（10）中，a，b，c，d為參數(shù)；SI 為地位指數(shù)；A為年齡；G為斷面積。

對(duì)加入SI 的模型引入林分株數(shù)（隨機(jī)變量）構(gòu)建混合效應(yīng)模型，將林分株數(shù)分別加入到不同參數(shù)上及其組合上，有10 種組合形式收斂（表6）。對(duì)具有相同參數(shù)個(gè)數(shù)的模型，通過比較AIC、BIC和LL 選出最優(yōu)模型，對(duì)具有不同參數(shù)個(gè)數(shù)的模型還需用似然比（LRT）檢驗(yàn)兩個(gè)模型之間差異是否顯著（ρ＜0.000 1 時(shí)，差異顯著），從表6中可以看出，含有1 個(gè)隨機(jī)參數(shù)時(shí)，模型2 的AIC、BIC 最小，LL 值最大，擬合效果最好；含有2 個(gè)隨機(jī)參數(shù)時(shí)，模型7 的AIC、BIC 最小，LL 值最大，擬合效果最好；為了避免參數(shù)過多導(dǎo)致模型復(fù)雜化，對(duì)具有不同參數(shù)個(gè)數(shù)的模型進(jìn)行似然比檢驗(yàn)，模型2 與基礎(chǔ)模型和模型7 差異均顯著，因此選擇參數(shù)較少的模型2 作為湖南櫟類次生林林分?jǐn)嗝娣e生長混合效應(yīng)模型，模型參數(shù)值見表7。表達(dá)式為：

表6 含密度等級(jí)混合效應(yīng)模型擬合結(jié)果Table 6 The fitting results of mixed effect model with density grades

表7 混合效應(yīng)模型參數(shù)值Table 7 The parameter value of mixed effect model

式（11）中，a，b，c，d為固定參數(shù)；β為隨機(jī)參數(shù)；SI 為地位指數(shù)；A為年齡；G為斷面積。

根據(jù)混合效應(yīng)模型的表達(dá)式可知，SI 對(duì)斷面積生長的最大值影響顯著，林分株數(shù)則主要影響斷面積的生長速度，因此以林分株數(shù)為連續(xù)變量構(gòu)建參數(shù)化模型，其擬合精度較高，確定系數(shù)R2為0.964，預(yù)估精度P值為98.36%，均方根誤差（RMSE）為1.427，3 項(xiàng)指標(biāo)均優(yōu)于混合效應(yīng)模型，因此將林分株數(shù)為連續(xù)變量的參數(shù)化模型作為湖南櫟類次生林林分?jǐn)嗝娣e生長模型，表達(dá)式為：

3.3 模型的檢驗(yàn)與評(píng)價(jià)

用未參加建模的樣本，對(duì)模型進(jìn)行檢驗(yàn)，看其是否適用，并對(duì)3 種模型進(jìn)行比較（表8）。從表中可以看出：相較于基礎(chǔ)模型，混合效應(yīng)模型及參數(shù)化模型預(yù)估精度P值從92.92%提升到98.36%（建模數(shù)據(jù)）、90.91%提升到97.74%（檢驗(yàn)數(shù)據(jù)），確定系數(shù)R2從0.314 提升到0.964（建模數(shù)據(jù)）、0.400 提升到0.964（檢驗(yàn)數(shù)據(jù)），平均誤差值（ME）、均方根誤差（RMSE）都有所減?。粎?shù)化模型預(yù)估精度優(yōu)于混合效應(yīng)模型，但各項(xiàng)評(píng)價(jià)指標(biāo)差異較小。

表8 3 種模型的評(píng)價(jià)指標(biāo)Table 8 The evaluation index of three models

用檢驗(yàn)數(shù)據(jù)對(duì)3 種模型分別作林分?jǐn)嗝娣e的殘差分布圖（圖1），可以直觀清晰地看出殘差分布范圍：參數(shù)化模型＜混合效應(yīng)模型＜基礎(chǔ)模型，且殘差分布基本均勻分布在橫軸兩側(cè)，說明混合效應(yīng)模型和參數(shù)化模型均優(yōu)于基礎(chǔ)模型，且參數(shù)化模型擬合效果最佳。

圖1 3 種模型殘差分布Fig.1 The residual distribution of three kinds of model model

3.4 林分?jǐn)嗝娣e生長曲線

根據(jù)構(gòu)建的湖南櫟類次生林林分?jǐn)嗝娣e生長參數(shù)化模型，將不同株數(shù)及地位指數(shù)分別代入到參數(shù)化模型中，便可得到各立地指數(shù)在不同林分密度下林分?jǐn)嗝娣e生長曲線（圖2）。從圖2中可以看出：林分年齡和林分密度相同時(shí)，地位指數(shù)級(jí)越大林分?jǐn)嗝娣e總生長量越大；當(dāng)?shù)匚恢笖?shù)和林分年齡相同時(shí)，密度越大林分?jǐn)嗝娣e總生長量越大；9、12 指數(shù)級(jí)在6 種密度條件下林分?jǐn)嗝娣e總生長量都隨著林分年齡的增加而增加，且密度越大林分?jǐn)嗝娣e呈“S”型生長規(guī)律越明顯，相鄰兩個(gè)密度間斷面積總生長量的差距越小。

圖2 不同林分密度下林分?jǐn)嗝娣e生長曲線Fig.2 The growth curve of stand basal area with different stand density

3.5 林分?jǐn)嗝娣e定期生長量

通過構(gòu)建的櫟類林分?jǐn)嗝娣e生長參數(shù)化模型，預(yù)測(cè)間隔期為5 a 的林分?jǐn)嗝娣e生長量，圖3為各地位指數(shù)在不同林分密度下，林分?jǐn)嗝娣e5 a 間的生長量曲線，其符合林木斷面積定期生長量的規(guī)律：定期生長量隨著林分年齡的增加而增加，到定期生長峰值后逐漸下降。從圖3中可以看出：在9、12 指數(shù)級(jí)中，林分年齡相同時(shí)，林分密度越大斷面積5 a 間生長量越大，到最大生長量后密度越大生長量下降速度越快；隨著林分年齡的增加，不同的初始密度條件下5 a 間生長量的差距由小變大至峰值后減小，直至密度對(duì)5 a 定期生長量無顯著影響；9 和12 指數(shù)級(jí)中5 a 間斷面積生長量最大值出現(xiàn)的時(shí)間均為：密度I 未達(dá)到最大值，密度II 為45 a，密度III 為35 a，密度IV、V 為25 a，密度VI 為20 a，由此可以看出，林分密度越大，5 a 的生長量最大值出現(xiàn)的越早（20 a 時(shí)）。

圖3 不同林分密度下5年間林分?jǐn)嗝娣e生長量曲線Fig.3 The growth curve of stand basal area in 5 years under different stand density

4 結(jié)論與討論

4.1 結(jié) 論

將湖南省6 期國家連續(xù)清查數(shù)據(jù)（1989—2014年）中以櫟類為優(yōu)勢(shì)樹種或主要樹種的樣地作為研究對(duì)象，以Richards生長模型作為基礎(chǔ)模型，通過構(gòu)建加入林分株數(shù)、地位指數(shù)的混合效應(yīng)模型及參數(shù)化模型，得出以下結(jié)論：

1）相較于基礎(chǔ)模型，混合效應(yīng)模型和參數(shù)化模型對(duì)櫟類林分?jǐn)嗝娣e生長的預(yù)估精度顯著提高，預(yù)估效果較好。參數(shù)化模型的預(yù)估精度比混合效應(yīng)模型稍高，最終將參數(shù)化模型作為湖南櫟類次生林林分?jǐn)嗝娣e生長模型。

2）9、12 地位指數(shù)級(jí)均表現(xiàn)為在各密度下林分?jǐn)嗝娣e總生長量隨著林分年齡的增加而增加；密度越大林分?jǐn)嗝娣e呈“S”型生長規(guī)律越明顯，相鄰兩個(gè)密度間斷面積總生長量的差距越小；在林分年齡和地位指數(shù)級(jí)相同時(shí)，隨著林分密度的增加，林分?jǐn)嗝娣e總生長量增加；在相同林分年齡和密度下，地位指數(shù)越大林分?jǐn)嗝娣e總生長量越大。

3）在9、12 地位指數(shù)級(jí)中均表現(xiàn)為林分年齡相同時(shí)，林分密度越大林分?jǐn)嗝娣e5 a 間生長量越大，到最大生長量后密度越大生長量下降速度越快；林分密度太低則沒有明顯的定期生長量峰值，密度越大，峰值出現(xiàn)的林分年齡越小，在林分年齡為20～45 a 之間；在相同林分年齡和密度下，地位指數(shù)越大林分5 a 的生長量越大。

4.2 討 論

本研究構(gòu)建以林分密度為隨機(jī)變量的湖南櫟類次生林林分?jǐn)嗝娣e混合效應(yīng)模型的主要目的，是為了分析不同密度下櫟類林分?jǐn)嗝娣e生長狀況，為制定營林措施提供科學(xué)依據(jù)。備選模型中Richards 生長模型擬合效果最佳，這與彭其龍等[4]研究成果相符?；旌闲?yīng)模型和參數(shù)化模型的預(yù)估精度優(yōu)于基礎(chǔ)模型這一結(jié)論與符利勇等[11]構(gòu)建的混合效應(yīng)模型和冉啟香等[20]構(gòu)建的林分?jǐn)嗝娣e生長模型相符。

本研究所構(gòu)建的以林分株數(shù)為連續(xù)變量的參數(shù)化模型，其總生長量曲線呈“S”型，定期生長量符合“隨著年齡增長而逐漸增加，到一定年齡后逐漸下降”的規(guī)律[15]。對(duì)斷面積生長量的研究較少，本文關(guān)于林分?jǐn)嗝娣e總生長量和5 a 的生長量的結(jié)論與孫洪剛[21]的結(jié)論相符。

朱光玉等[3]構(gòu)建的林分?jǐn)嗝娣e生長模型中，用林分密度指數(shù)擬合效果優(yōu)于林分株數(shù)，本文沒有對(duì)這兩個(gè)表示密度的指標(biāo)進(jìn)行比較。本文建立的林分?jǐn)嗝娣e生長參數(shù)化模型具有較高的預(yù)估精度，可以為櫟類次生林制定合理的營林措施提供依據(jù)，但沒有考慮林分樹種組成以及對(duì)林分株數(shù)和林分密度指數(shù)的對(duì)比，在下一步研究中可以加入林分樹種組成和對(duì)林分株數(shù)、林分密度指數(shù)分別建模的比較對(duì)模型進(jìn)行完善。