基于分位數回歸的針闊混交林樹高與胸徑的關系

張冬燕，王冬至，李 曉，高雨珊，李天宇，陳 靜

（1.河北農業大學 林學院，河北 保定 071000；2.河北農業大學 經濟管理學院，河北 保定 071000）

樹高和胸徑不僅是用來預測林分蓄積量[1]、生物量[2]、立地生產力[2-3]及林分結構[4]的重要變量，而且是森林資源調查及經營效果評價的重要因子。在標準地調查過程中，樹高測量難度較大且觀測成本高，其觀測誤差也相對較大，給精準林業質量提升帶來了一定困難[2,4]，而胸徑觀測方便且精度較高。因此，根據標準地調查數據，建立樹高與胸徑關系預測模型，可降低調查成本，提高預測精度[5]，這對于森林質量精準提升具有重要意義。當前林業研究多采用線性或非線性樹高與胸徑關系模型來模擬預測兩者之間的關系，其參數估計方法多采用最小二乘法來模擬，然而最小二乘法是基于均值回歸，利用變量均值來擬合模型參數[6]，該方法要求調查數據需滿足獨立正態同分布等條件。在林業調查中，樹高與胸徑觀測數據不能滿足該要求，而分位數回歸對調查數據沒有嚴格要求[7]，利用變量條件分位數來建模[6]，對具有尖峰、厚尾、異方差顯著的數據擬合效果更加穩健[8-10]。分位數回歸理論框架是由KOENKER[11]提出，已在醫學[12]、經濟學[13]、教育與政策[14-15]及自然資源管理等領域進行了研究與應用。在林業相關研究中，分位數回歸被應用于模擬林分自疏邊界線[16]、直徑分布規律[17]、林分密度指數[18]及森林病蟲害[19]等方面研究。?Z?ELIK等[20]基于分位數回歸建立了樹高與胸徑關系模型，高慧淋等[21]采用此方法建立了長白落葉松Larix olgensis人工林最大林分密度線模型，提高了模型預測精度及適用性。然而在華北暖溫帶針闊混交林中，如何基于一個分位數回歸模型，預測不同樹種樹高與胸徑關系是亟待解決的科學問題。在混交林中為了描述樹種結構對樹木生長影響，部分學者[22-25]采用啞變量方法構建了不同間伐方式、不同地域樹高曲線及生長量預測模型。然而基于包含啞變量的非線性分位數回歸方法來構建不同樹種樹高與胸徑關系模型的研究較少。因此，本研究以河北省塞罕壩華北落葉松Larix principisrupprechtii與白樺Betula platyphylla針闊混交林為研究對象，基于啞變量方法和分位數回歸相結合方法，構建混交林不同樹種分位數回歸模型，為精確描述樹高與胸徑的關系提供依據。

1 研究區概況

河北省塞罕壩機械林場(41°22′～42°58′N，116°53′～118°31′E)位于河北省最北部，地勢北高南低，屬華北暖溫帶立地類型區，林場總面積約9.2×104hm2，總蓄積約8.1×106m3。土壤類型以褐色森林土、棕色森林土及風沙土等為主；成土母質主要為坡積物、殘積物及洪積物等；極端最高氣溫為33.4 ℃，最低氣溫-43.3 ℃，年均氣溫-1.3 ℃，年均無霜期64 d，年均降水量460.3 mm，是典型的半干旱半濕潤寒溫性大陸季風氣候。研究區植被類型豐富，主要喬木樹種有華北落葉松、白樺、樟子松Pinus sylvestris、云杉Picea asperata等，主要灌木樹種有山刺玫Rosa daverica、胡枝子Lespedeza bicolor、沙棘Hippophae rhamnoides等，主要草本植物有地榆Sanguisorba offcinalis、唐松草Thalictrum aquilegifolium、曼陀羅Datura stramonium等。

2 研究方法

2.1 數據來源

在北曼店、大喚起、陰河、千層板和第三鄉等5個林場設立了83塊標準地(30 m×30 m)，對標準地內各林分因子(林分密度、平均高、平均胸徑、樹種斷面積、林分總斷面積、優勢高等)和立地因子(海拔、坡度、坡向、坡位、土層厚度等)進行調查，共調查立木10 104株(華北落葉松5 258株，白樺4 846株)，林分年齡分布為24～45 a，不同標準地混交度分布為0.39～0.62。研究過程中，分樹種將觀測數據分別按3∶1分為建模數據(62塊標準地)和檢驗數據(21塊標準地)，基本信息如表1和表2所示。

表1 模型建立數據Table 1 Data of establishment model

表2 模型檢驗數據Table 2 Data of test model

2.2 基礎模型

在描述樹木生長及樹高與胸徑關系的近百種不同模型中，Richard方程不但具有可解釋的生物學意義，而且具有易收斂且靈活性高等特性。部分研究基于Richard方程構建了不同林分類型樹高與胸徑關系的預測模型，均取得了較好的預測結果[18,20,26-29]。因此，本研究以Richard方程作為構建華北落葉松與白樺針闊混交林樹高與胸徑關系基礎模型，模型表達如式(1)所示。

式(1)中：Hij為第i個樣地第j株樹的樹高(m)；dij為第i個樣地第j株樹的胸徑(cm)；a、b、c為基礎模型的參數；εij為誤差項。

2.3 啞變量

為了解決模型預測精度的影響，可以在模型中加入啞變量[23,30-32]。包含啞變量的樹高與胸徑關系預測模型，不僅可以實現模型對不同樹種相容性，而且在一定程度上可以提供模型預測精度及適用性，包含啞變量的樹高與胸徑關系預測模型表達如式(2)所示。

式(2)中：Mi為啞變量，當M1=1、M2=0時為華北落葉松，當M1=0、M2=1時為白樺；ai、bi、ci為模型參數；εij為誤差項。

2.4 分位數回歸

由于分位數回歸對模型誤差不需要嚴格假設條件，因此本研究基于Richard方程，選取5個分位點(τ=0.1、0.3、0.5、0.7、0.9)構建不同樹種的樹高與胸徑關系預測模型，利用加權最小一乘法可以得到不同分位點參數估計值，具體見式(3)。

式(3)中：S為不同分位點估計值；τ、Hij分別為第i個樣地第j株樹在不同分位點τ樹高預測值與樹高值(m)；dij第i個樣地第j株樹胸徑(cm)；τ為分位點。

2.5 模型檢驗與評價

統計分析均基于SPSS 24.0和SAS 9.4中的PROC NLIN和PROC NLP完成，基于模型確定系數(R2)、平均差(MD)、平均絕對誤差(MAD)對模型擬合精度及適用性進行評價與比較。

式(4)～(6)中：Hij、和分別為樹高觀測值、預測值和平均值；m為標準地數量；n為標準地株數。

3 結果與分析

3.1 林分特征分析

圖1為不同樹種胸徑與樹高的關系。華北落葉松和白樺的樹高分別為4～18和6～16 m，胸徑分別為6～32和6～28 cm。

圖1 不同樹種胸徑與樹高分布Figure 1 Height-diameter distribution of different tree species

3.2 基礎模型和分位數模型擬合與評價

從表3可見：在華北落葉松與白樺針闊混交林中，基于分位數回歸的不同樹種不同分位點確定系數均大于傳統回歸方法，平均差及平均絕對誤差均小于傳統回歸方法。在確定的5個分位點中，當分位點τ=0.7時，華北落葉松與白樺的樹高與胸徑關系預測模型精度最高。基于不同分位點預測模型建立了不同樹種在各分位點殘差分布圖(圖2)，確定當分位點位為τ=0.7時，華北落葉松與白樺樹高與胸徑關系模型能夠更好描述兩者之間的關系。

表3 基礎模型和分位數回歸模型擬合與評價Table 3 Fitting and evaluation of basic model and quantile regression model

圖2 不同樹種各分位點殘差分布Figure 2 Residual distribution of each quantile of different tree species

3.3 不同樹種樹高模擬

基于不同分位點預測模型，分別對華北落葉松和白樺的樹高與胸徑關系進行了模擬(圖3)。不同樹種在不同分位點樹高與胸徑關系預測趨勢及范圍基本一致，表明包含啞變量的分位數回歸模型預測效果較好。

圖3 不同分位點樹高與胸徑的關系Figure 3 Relationship between tree height and DBH of different quantiles

4 討論

在描述樹高與胸徑關系的多種線性和非線性預測模型中，通常采用具有生物學意義且靈活性較高的Richard方程來研究不同林分類型樹高與胸徑的關系[26,33]，因此，本研究將Richard方程作為構建華北落葉松與白樺混交林樹高與胸徑關系基礎模型。對紅松Pinus koraiensis人工林、土耳其松Pinus brutia和黎巴嫩雪松Cedrus libani混交林的研究[20,26-28]表明：Richard方程是描述其樹高與胸徑關系的最優模型。

在構建樹木生長及生物量預測模型中，包含啞變量模型具有更高的預測精度[22,25,31-32]。本研究在華北落葉松與白樺針闊混交林中，基于非線性分位數回歸構建的包含啞變量樹高與胸徑的關系模型，其精度高于傳統回歸預測模型，這與人工林最大密度線確定[21]及樹高與胸徑關系模型[20]的研究結論相一致，表明非線性分位數回歸較傳統回歸方法更加穩定，可用于人工林和混交林立地潛在生產力的評價。

5 結論

本研究以塞罕壩華北落葉松與白樺針闊混交林調查數據為基礎，確定Richard模型為描述不同樹種樹高與胸徑關系的基礎模型，在基礎模型中構造一個表示樹種的啞變量，并利用分位數回歸在一個模型中同時估計不同樹種及不同分位點的樹高與胸徑關系模型參數，經過檢驗不同樹種分位數回歸模型均能較好反映樹高與胸徑的關系，當分位點τ=0.7時，分位數回歸模型預測精度最高，擬合效果最好。