昆明市4個主要針葉樹種林分斷面積生長模型研究

吳 恒 朱麗艷 李 華 羅春林 吳雪瓊

(國家林業(yè)局昆明勘察設計院，云南 昆明 650216)

林分立地質(zhì)量的準確評價是森林科學經(jīng)營的前提，地位級法作為立地質(zhì)量評價的可靠手段在林業(yè)生產(chǎn)實踐中被廣泛的運用[1-2]。林分斷面積生長模型以林分總體特征為基礎預估整個林分的生長和收獲量。 Richards和Schumacher是可變密度生長模型的2種主要形式，Richards模型適應性強、擬合優(yōu)度大[3]，Schumacher模型形式簡單、實用性強[4]。對不同生長模型的性能比較將有助于模型形式的不斷改進。

基于生物學理論或假設而推導出來的理論生長模型，其參數(shù)有一定的生物學意義并且適用性較強[2-5]，林分斷面積生長模型參數(shù)求解常將非線性方程線性化，再利用最小二乘法或極大似然估計法等求解[6-7]。最小二乘法和極大似然估計法都是建立在具有連續(xù)可導的光滑搜索空間的假設基礎上，且是在梯度方向上尋優(yōu)的局部搜索技術[8]，可能陷入局部最優(yōu)解，參數(shù)求解時搜索速度隨反演參數(shù)增多呈級數(shù)減慢，誤差傳遞導致不收斂等缺點[9-12]，參數(shù)求解誤差影響著林分斷面積生長模型的預測精度。昆明市地處云貴高原中部，分布有云貴高原典型森林植被類型，但長期以來缺少對林區(qū)主要針葉樹種立地質(zhì)量的準確評價和生長過程的模型研究。采用優(yōu)化算法進行林分斷面積生長模型參數(shù)估計，對提高模型預測精度和確保參數(shù)穩(wěn)定性具有積極意義。本研究旨在建立昆明市4種主要針葉樹種地位級表和林分斷面積生長模型，基于不同優(yōu)化算法建立針葉樹種林分斷面積生長模型，比較不同優(yōu)化算法效率，分析模型形式和參數(shù)求解穩(wěn)定性，為優(yōu)化算法在生長模型參數(shù)估計中的運用提供參考。

1 材料與方法

1.1 數(shù)據(jù)來源與處理

基于昆明市森林資源二類調(diào)查數(shù)據(jù)(2015—2016)，選取昆明市主要針葉樹種華山松 (Pinusarmandii)、云南松 (Pinusyunnanensis)、油杉 (Keteleeriafortunei) 和柏木 (Cupressusfunebris) 林分為研究對象，篩選樹種組成系數(shù)大于6成的林分作為 “相對純林” 小班[13]，計算林分平均年齡、平均樹高、平均胸徑、每公頃斷面積平均值，以 ±3倍標準差為上下限剔除異常值[1]，各樹種滿足建模要求的數(shù)據(jù)樣本量分別為華山松19 275組、云南松35 821組、油杉2 847組、柏木3 799組，樣本基本信息見表1。

表1 樣本基本信息Table 1 Basic information of sample

1.2 立地質(zhì)量評價方法

采用洪玲霞等[14]在建立蒙古櫟 (Quercusmongolica) 林全林分生長模型時使用的地位級指數(shù) (SCI) 作為立地質(zhì)量評價指標，建立各針葉樹種地位級模型；采用Richards、Schumacher和Compertz模型擬合林分平均年齡和平均樹高間的關系，模型形式見公式 (1)～(3)，參數(shù)a1表示林分平均樹高生長的上限水平，參數(shù)a2與平均樹高生長速率有關，參數(shù)a3決定了樹高生長曲線形狀。依據(jù)林分平均年齡與樹高關系建立地位級曲線簇，編制各針葉樹種地位級表，并進行落點檢驗，將檢驗樣本作散點圖，繪制到地位級曲線簇中，算出散點落在曲線簇內(nèi)的概率，即為地位級表能夠解釋樹高生長的概率值。

HT=a1(1-e-a2A)a3

(1)

HT=a1e-a2e-a3A

(3)

式中：HT為林分平均高；A為林分平均年齡；a1、a2和a3為參數(shù)。

1.3 建模方法

森林的生長和收獲取決于年齡、立地、林分密度等3個主要因子，本研究結(jié)合數(shù)據(jù)特征利用林分密度指數(shù) (SDI) 作為競爭指標[15]，引入到生長模型中擬合不同立地條件下林分的生長過程，其計算方法見公式 (4)。

備選模型有Schumacher、Richards和Korf 3種，模型形式見公式 (5)～(7)。模型的參數(shù)估計為在最小二乘法意義下極小化離差平方和，采用麥夸特算法 (LM)[16]、差分進化算法 (DE)[12]、遺傳算法 (GA)[17]、模擬退火算法 (SA)[11]和粒子群算法 (PSO)[8-10]，運用MatlabR2014b編寫各算法代碼進行最優(yōu)值求解，得到不同優(yōu)化算法下林分斷面積生長模型參數(shù)，根據(jù)擬合優(yōu)度 (R2)、均方根誤差 (RMSE) 綜合考慮擬合參數(shù)穩(wěn)定性，篩選一組最優(yōu)參數(shù)估計值作為模型擬合結(jié)果。

式中：G為林分每公頃斷面積；A為林分平均年齡；SCI為地位級指數(shù)；SDI為林分密度指數(shù)；a1、a2、a3、a4和a5為參數(shù)。

1.4 模型擬合結(jié)果評價方法

采用殘差對模型擬合結(jié)果進行分析，比較模型形式對林分生長預測的偏差。采用不同算法求解參數(shù)迭代次數(shù)比較算法效能，并用各算法擬合不同生長模型參數(shù)間歐式距離分析求解參數(shù)穩(wěn)定性。不同擬合參數(shù)間歐式距離的計算方法見公式 (8)。

(8)

2 結(jié)果與分析

2.1 林分年齡與樹高的模型擬合

立地質(zhì)量的準確評價是林分斷面積生長模型構建的基礎，也是林分科學經(jīng)營的前提。根據(jù)林分平均年齡和平均樹高曲線擬合結(jié)果 (表2) 可知，Richards模型和Compertz模型擬合優(yōu)度均高于Schumacher模型，Richards模型擬合華山松、云南松和柏木的生長上限水平、生長速率和曲線形狀相近，符合樹高生長的生物學規(guī)律，因此4個樹種均采用Richards模型。

根據(jù)林分平均年齡和平均樹高建立的地位級曲線簇，落點檢驗見圖1。華山松落點檢驗值為97.9%、云南松為98.3%、油杉為98.1%、柏木為98.9%。建立的地位級曲線簇通過落點檢驗，精度滿足立地質(zhì)量評價適用性檢驗，可編制各樹種地位級表，用于指導林業(yè)生產(chǎn)。

2.2 林分斷面積生長模型參數(shù)擬合

采用不同優(yōu)化算法擬合林分斷面積生長模型參數(shù)，擬合結(jié)果見表3。各生長模型R2介于0.87～0.92、平均值為0.90，RMSE介于1.48～1.99、平均值為1.80，符合模型擬合精度要求。4個針葉樹種模型形式均采用Richards模型，參數(shù)a1和a2與立地質(zhì)量的生長上限有關，參數(shù)a4和a5與林分密度的控制的生長速率有關。本研究中華山松和云南松林分生長上限參數(shù)相似，且均高于油杉和柏木，云南松的生長速率參數(shù)與其他針葉樹種存在差異。華山松、云南松、油杉和柏木模型預測值和觀測值的相關系數(shù)分別為0.95、0.94、0.93和0.96。華山松、云南松、油杉和柏木模型的t檢驗顯著性p值分別為0.93、0.88、0.97和0.58，表明4個針葉樹種模型預測值和觀測值間差異均不顯著。

表2 林分年齡-樹高曲線模型參數(shù)Table 2 Parameters of age and height growth model

圖1 各樹種地位級曲線簇及落點檢驗Fig.1 Site class curves for different tree species and falling point test results

表3 各針葉樹種林分斷面積生長模型參數(shù)擬合結(jié)果Table 3 Fitting results of stand basal area growth models for different coniferous species

2.3 模型擬合結(jié)果評價

2.3.1參數(shù)穩(wěn)定性分析

Richards模型擬合4個針葉樹種R2平均值為0.90，Schumacher模型擬合4個針葉樹種R2平均值為0.79，Korf模型擬合4個針葉樹種擬合優(yōu)度平均值為0.88。不同模型形式殘差分布見圖2。

圖2 不同模型形式殘差分布Fig.2 Distribution of residuals of different models

根據(jù)不同模型形式的殘差分布圖 (圖2) 可知，Richards模型擬合林分每公頃斷面積生長預測值與觀測值間殘差分布合理，隨林分每公頃斷面積的增大，殘差變動幅度逐漸擴大，未表現(xiàn)出系統(tǒng)性偏差；Schumacher模型擬合林分每公頃斷面積生長預測值與觀測值間殘差分布不合理，林分每公頃斷面積介于0～5 m2和20～40 m2時，預測值系統(tǒng)性低于觀測值，林分每公頃斷面積介于5～20 m2時，隨林分每公頃斷面積的增大，殘差變動幅度逐漸減??；Korf模型擬合林分每公頃斷面積生長預測值與觀測值間殘差分布較合理，殘差變動幅度大于Richards模型，林分每公頃斷面積介于0～5 m2時，預測值系統(tǒng)性高于觀測值。因此，針葉樹種更適宜采用Richards模型。

不同優(yōu)化算法求解各針葉樹種Richards模型各參數(shù)擬合結(jié)果間歐式距離為1 988.44，其中華山松斷面積生長采用不同優(yōu)化算法求解參數(shù)間歐式距離最大為7 912.49；Schumacher模型各參數(shù)擬合結(jié)果間歐式距離為4.08；Korf模型各參數(shù)擬合結(jié)果間歐式距離為8.44。采用Schumacher模型和Korf模型擬合參數(shù)穩(wěn)定性優(yōu)于Richards模型，但Richards模型的擬合優(yōu)度大于Schumacher模型和Korf模型。

2.3.2算法效能分析

不同優(yōu)化算法求解斷面積生長模型參數(shù)迭代次數(shù)見表4。

表4 不同優(yōu)化算法求解斷面積生長模型參數(shù)迭代次數(shù)Table 4 Iterations of parameter solving for basal area growth models by different algorithms

由表4可知，LM、DE、GA、SA和PSO迭代次數(shù)平均值分別為33、1 042、3 137、5 000和2 003。LM迭代次數(shù)平均迭代次數(shù)最小，SA迭代次數(shù)均達到了算法設置的最大迭代次數(shù) (5 000次)。不同模型形式對優(yōu)化算法求解參數(shù)效能也存在差異，Richards模型擬合迭代次數(shù)普遍高于Korf模型和Schumacher模型，Schumacher模型擬合迭代次數(shù)普遍最小。PSO迭代次數(shù)介于中間，但求解參數(shù)的擬合優(yōu)度較差。SA、DE對不同模型的參數(shù)擬合優(yōu)度較高。

3 結(jié)論與討論

采用不同優(yōu)化算法對昆明市4個主要針葉樹種進行林分斷面生長模型構建，比較不同優(yōu)化算法求解參數(shù)效率和穩(wěn)定性，得出以下結(jié)論：

1) Richards模型擬合華山松、云南松和柏木的生長上限水平、生長速率和曲線形狀相近，符合樹高生長的生物學規(guī)律；建立的地位級曲線簇通過落點檢驗，精度滿足立地質(zhì)量評價適用性檢驗。

2) 建立了華山松、云南松、油杉和柏木可變密度斷面積生長模型，模型R2平均值為0.90，均方根誤差平均值為1.80，模型精度符合要求。

3) Richards模型對針葉樹種生長擬合優(yōu)度大于Schumacher模型和Korf模型，Schumacher模型采用不同算法求解參數(shù)穩(wěn)定性大于Richards模型和Korf模型。

4) 各優(yōu)化算法求解生長模型參數(shù)效率由高到低順序為LM、DE、PSO、GA、SA，但PSO求解參數(shù)的擬合優(yōu)度較差。

林分斷面積生長模型中，Richards模型[18]、Schumacher模型[19-20]和Korf模型[4]具有較好的數(shù)學和生物學意義，得到普遍認可，但本研究表明，Richards模型的擬合優(yōu)度和解釋性更好，但參數(shù)求解較困難且不穩(wěn)定；Schumacher模型擬合參數(shù)更易求解且比較穩(wěn)定，但存在系統(tǒng)性偏差，在實際運用中需要對模型形式進行修正；Korf模型作為Schumacher模型的廣義形式R2和參數(shù)穩(wěn)定性均介于中間水平，但仍然存在特定階段的系統(tǒng)性偏差，相較于Schumacher模型有較好的改善。利用再參數(shù)化方法建立生長模型能夠兼顧Richards模型和Schumacher模型各自的優(yōu)點[4]。

林分斷面積生長模型是復雜的非線性回歸模型[4]，常用的非線性模型參數(shù)估計方法是基于估計目標函數(shù)對優(yōu)化變量梯度信息的優(yōu)化，要求優(yōu)化模型具有連續(xù)性和可導性[8]，而生長模型的系統(tǒng)非線性特點增加了梯度信息局部搜索方法尋找全局最優(yōu)解的難度和誤差傳遞的不收斂，從而影響生長模型精度[9-12]。運用單一的優(yōu)化算法對非線性模型參數(shù)進行估計限制了不同優(yōu)化算法間效能比較和參數(shù)穩(wěn)定性分析[21]，降低了各優(yōu)化算法在實際運用中的可選性。對比分析不同優(yōu)化算法擬合林分斷面積生長模型的性能，為林分斷面積生長模型參數(shù)求解方法提供了依據(jù)。

林分立地質(zhì)量評價和生長模型是森林經(jīng)營的重要工具。全林分生長模型能較好的模擬人工純林的生長[22]，采用針葉樹種組成系數(shù)大于等于6成的小班作為建模樣本單元，符合 “相對純林” 的概念[1,13]，能夠用于全林分生長模型的構建。建立的全林分生長模型擬合優(yōu)度介于0.87～0.92，具有較好的預測效能[4]。采用建立的全林分生長模型提取林分小班因子中年齡、平均樹高、公頃株數(shù)實現(xiàn)對小班樹高、胸徑、每公頃斷面積、蓄積等因子的更新，與采用林分蓄積生長率預測林分蓄積相比，傳統(tǒng)生長率更新林分因子的方法不能反映小班間生長差異性，不能實現(xiàn)樹高、胸徑和公頃株數(shù)等林分因子的更新變化。建立全林分生長模型，對探索無干擾林分動態(tài)監(jiān)測和林地 “一張圖” 蓄積量等因子更新具有重要作用。地位級指數(shù)作為林分立地質(zhì)量評價的指標在缺乏優(yōu)勢樹高數(shù)據(jù)的情況下用于林分立地質(zhì)量的評價具有可行性，本研究建立的華山松、云南松、油杉和柏木的地位級指數(shù)模型落點檢驗值均高于97%，具有較好的普適性和外推效能[13]，能夠用于昆明市針葉樹種造林立地質(zhì)量評價與劃分。