李斌成，許業(yè)洲，袁 慧，侯義梅，李雙龍，杜超群

(1.湖北省林業(yè)科學(xué)研究院，湖北武漢430075；2.長江大學(xué)園藝園林學(xué)院，湖北荊州434025；3.建始縣林業(yè)局，湖北建始445300；4.恩施土家族苗族自治州林業(yè)科學(xué)研究所，湖北恩施445000)

立地指森林或其他植被類型生存的空間及與之相關(guān)的自然因子的綜合，在一定時(shí)間內(nèi)相對不變[1]。森林的立地質(zhì)量對森林生產(chǎn)力有直接影響，對森林立地質(zhì)量作出準(zhǔn)確評估，不僅可以方便管理，還能為改造森林提供科學(xué)依據(jù)。評價(jià)立地質(zhì)量的方法較多，但最有效、最客觀的方法是地位指數(shù)法[2]。對于立地指數(shù)模型的擬合，傳統(tǒng)研究大多使用固定基準(zhǔn)年齡的靜態(tài)模型，隨著統(tǒng)計(jì)學(xué)以及計(jì)算機(jī)科學(xué)的不斷發(fā)展，人們逐漸發(fā)展了基準(zhǔn)年齡可變的動(dòng)態(tài)方程來擬合立地指數(shù)模型，任何一個(gè)林齡都可以得出一條相同的立地指數(shù)曲線，動(dòng)態(tài)模型能夠明顯提高模型的擬合效果[3]。

近年來，國內(nèi)利用動(dòng)態(tài)模型來模擬優(yōu)勢高生長的相關(guān)研究較為成熟，段愛國等[4]利用代數(shù)差分法(algebraic difference approach，ADA)建立了杉木[Cunninghamialanceolata(Lamb.)Hook.]人工林多形地位指數(shù)模型，分別利用理論模型及其差分模型模擬了杉木樹高生長過程，擬合結(jié)果表明，差分方程精度較理論模型有所提高，且表現(xiàn)出更好的生物學(xué)解釋。唐誠[5]利用差分法建立了西南樺(BetulaalnoidesBuch.-Ham.ex D.Don)人工林生長模型，其擬合效果總體優(yōu)于理論模型。廣義代數(shù)差分法(generalized algebraic difference approach，GADA)也是動(dòng)態(tài)模型的一種，因其能夠建立具有可變漸進(jìn)值的多形地位指數(shù)曲線簇而在近些年廣泛被使用。曹元帥等[2]利用GADA研究了杉木人工林地位指數(shù)模型，得出應(yīng)用GADA構(gòu)建地位指數(shù)動(dòng)態(tài)模型有更好的模擬效果。一些學(xué)者的相關(guān)研究[6-11]也表明了GADA相較于ADA具有更好的擬合效果。

日本落葉松[Larixkaempferi(Lamb.)Carr]在我國栽培區(qū)域北起黑龍江省林口縣，南至湖南省城步縣和四川省雷波縣[3]，湖北省是日本落葉松最早引種區(qū)域之一，該樹種目前已發(fā)展成為鄂西亞高山區(qū)主要造林樹種和重要工業(yè)原料林樹種。但湖北省日本落葉松人工林的立地質(zhì)量評價(jià)和林分生產(chǎn)力預(yù)估等方面的研究較少，李忠國[12]利用導(dǎo)向曲線法編制了恩施土家族苗族自治州(恩施州)日本落葉松立地指數(shù)表，馬友平等[13]利用比例法編制了長嶺崗林場日本落葉松的立地指數(shù)表，這些研究范圍較窄，不足以代表湖北省日本落葉松人工林的立地狀況，而且研究方法較為傳統(tǒng)，研究結(jié)果精度與適應(yīng)性有待進(jìn)一步驗(yàn)證。以經(jīng)典理論生長[6-8]模型為基礎(chǔ)，分別采用ADA和GADA構(gòu)建動(dòng)態(tài)立地指數(shù)模型，比較兩種模型參數(shù)擬合結(jié)果及預(yù)測精度的差異性，選擇最優(yōu)模型編制湖北省日本落葉松人工林立地指數(shù)表，以期為湖北日本落葉松人工林營造、分類經(jīng)營及林分定向培育提供參考和支持。

1 材料與方法

1.1 試驗(yàn)地概況

湖北省地處中國中部， 轄區(qū)范圍地跨東經(jīng) 108°21′～116°07′， 北緯 29°01′～ 33°06′， 屬于亞熱帶季風(fēng)氣候區(qū)，年平均氣溫16℃，年平均降水量800～1 600 mm，氣候濕潤，雨熱同期。湖北省地貌以山地、丘陵為主，約占全省面積的80%，東、西、北三面環(huán)山而高高凸起，中間低平。日本落葉松作為外來引種樹種，主要分布于海拔1 300～2 100 m的鄂西亞高山區(qū)，取樣地海拔為1 075～2 007 m，多集中在1 650 m左右，坡度2°～37°，以緩坡(15°)為主；腐殖層厚度2～19 cm，土壤類型以黃棕壤為主。

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

在造林面積大于15 hm2的日本落葉松人工純林中按照試驗(yàn)設(shè)計(jì)要求，選取不同年齡、不同立地、代表性強(qiáng)、林相完整的林分，設(shè)置600 m2(20 m×30 m)的標(biāo)準(zhǔn)地。標(biāo)準(zhǔn)地取樣主要分布于恩施州，而十堰市較少，這是因?yàn)楹笔∪毡韭淙~松主要生長于恩施州境內(nèi)，十堰市造林?jǐn)?shù)量有限。標(biāo)準(zhǔn)地內(nèi)進(jìn)行每木檢尺，分別測量樹高、胸徑、冠幅、枝下高等生長量指標(biāo)，并記錄相應(yīng)的地理位置、立地環(huán)境及林分造林年度、造林密度等基本信息。每塊標(biāo)準(zhǔn)地中以生長量最大的5株優(yōu)勢木的樹高平均值做為優(yōu)勢高，并選取1株優(yōu)勢木伐倒，按2 m區(qū)分段截取樹干圓盤并進(jìn)行樹干解析。以86個(gè)標(biāo)準(zhǔn)地及其優(yōu)勢木解析資料為本研究基礎(chǔ)數(shù)據(jù)(表1)。

表1 標(biāo)準(zhǔn)地基本信息Table 1 Basic information of standard area

1.3 模型選擇

1.3.1 理論生長模型選擇 以Richards、Logistic、Korf、Mitscherlich以及Gompertz等5種經(jīng)典樹高生長模型為基礎(chǔ)模型，其表達(dá)式見表2。將解析數(shù)據(jù)整理成樹高-年齡成對值形式用于理論模型的擬合。

表2 基礎(chǔ)模型表達(dá)式Table 2 Basic model expressions

1.3.2 代數(shù)差分模型 利用R語言中的nlsLM函數(shù)進(jìn)行模型擬合，根據(jù)基礎(chǔ)模型的擬合結(jié)果，將解析數(shù)據(jù)整理成雙樹高-雙年齡的形式用于差分模型的擬合。選擇擬合效果最好的Richards模型，分別利用ADA和GADA構(gòu)建其動(dòng)態(tài)指數(shù)差分模型。模型拐點(diǎn)表達(dá)式見文獻(xiàn)[3-4]，模型表達(dá)式見表3。

表3 動(dòng)態(tài)指數(shù)差分模型表達(dá)式Table 3 Expression of the dynamic exponential difference model

1.3.3 模型統(tǒng)計(jì)與評估 選用決定系數(shù)(coefficient of determination，R2)、平均殘差(mean residual，MR)、均方根誤差(root mean square error，RMSE)、絕對殘差(absolute mean residual，AMR)、相對殘差(relative absolute residual，RAR)以及殘差平方和(residual sum of square，RSS)等統(tǒng)計(jì)指標(biāo)評估模型的擬合質(zhì)量(統(tǒng)計(jì)指標(biāo)計(jì)算公式見文獻(xiàn)[4])。

1.4 立地指數(shù)表的編制

1.4.1 標(biāo)準(zhǔn)年齡與指數(shù)級 根據(jù)樹高理論變動(dòng)系數(shù)[14]及前期對于日本落葉松的相關(guān)研究結(jié)果[12]確定標(biāo)準(zhǔn)年齡T＝20 a，并利用標(biāo)準(zhǔn)年齡時(shí)的樹高絕對變動(dòng)幅度及經(jīng)營水平確定級距為2 m，立地指數(shù)級為7個(gè)(12～24 m)。

1.4.2 立地指數(shù)表編制 選擇最優(yōu)的模型，設(shè)H＝H2，t＝t2，IS＝H1，T＝t1，其中IS、T分別表示立地指數(shù)和指數(shù)年齡，將各變量代入差分模型求出各年齡優(yōu)勢高，整理并編制立地指數(shù)表。

1.4.3 立地指數(shù)表檢驗(yàn) 分別將基于ADA、GADA編制的立地指數(shù)表中的數(shù)據(jù)整理繪制成立地指數(shù)曲線簇，并利用86個(gè)標(biāo)準(zhǔn)地優(yōu)勢木平均高進(jìn)行落點(diǎn)檢驗(yàn)。按照立地指數(shù)表查定林分立地質(zhì)量所在的指數(shù)級，采用χ2適合性檢驗(yàn)方法，統(tǒng)計(jì)合格次數(shù)(不跳級)與不合格次數(shù)(跳級)的值，并利用皮爾遜定理求出η值與χ2值并進(jìn)行比較，按照森林調(diào)查90%的精度，計(jì)算出合格與不合格次數(shù)的理論值[15]。

選擇30株年齡大于20 a的優(yōu)勢解析木，按基準(zhǔn)年齡時(shí)的樹高值以2 a為齡階進(jìn)行分組，利用公式(1)對各齡階平均值進(jìn)行換算，得到檢驗(yàn)樣本的立地指數(shù)。利用公式(2)～(3)計(jì)算不同年齡(齡組)和不同立地指數(shù)級的估計(jì)誤差值，進(jìn)行立地指數(shù)表的預(yù)測精度檢驗(yàn)。

式中:H0為基準(zhǔn)年齡時(shí)的樹高(m)；Hi為第i齡階平均實(shí)際樹高(m)；H0k為基準(zhǔn)年齡時(shí)的導(dǎo)向曲線樹高(m)；Hik為第i齡階導(dǎo)向曲線樹高(m)；H(i，z)為年齡為i，立地指數(shù)為z時(shí)的實(shí)際樹高(m)；H(i，z)為年齡為i，立地指數(shù)為z時(shí)的平均樹高(m)；Si為齡階為i的立地指數(shù)估測誤差(m)；Sz為立地指數(shù)為z的齡階估測誤差(m)；l為立地指數(shù)序列的組數(shù)；jz為第z組優(yōu)勢木觀測齡階數(shù)。

2 結(jié)果與分析

2.1 模型擬合結(jié)果與比較

理論模型的擬合結(jié)果見表4。除Korf(a＝150.739 0)、Mitscherlich(b＝-48.294 2)模型的參數(shù)值異常，不符合所表征的生物學(xué)涵義外，其他模型參數(shù)均合理，決定系數(shù)(R2)為0.926 4～0.949 3。其中，R2排序?yàn)?Korf＞Richards＞Mitscherlich＞Gompertz＞Logistic； AMR 排序?yàn)?Logistic＞Gompertz＞Richards＞Mitscherlich＞Korf； RMSE 排序?yàn)?Logistic＞Gompertz＞Richards＞Mitscherlich＞Korf。 在保證參數(shù)具有生物學(xué)意義的前提下，結(jié)合統(tǒng)計(jì)指標(biāo)的大小，可以最終R2(0.947 9)大且AMR(1.281 9)和RMSE(1.621 3)小的Richards模型為建?；A(chǔ)模型，并分別利用ADA和GADA構(gòu)建Richards模型的差分形式，結(jié)果見表5。

表4 理論模型擬合結(jié)果Table 4 Fitting results for the theoretical model

由表5可知，利用Richards基礎(chǔ)模型構(gòu)建的ADA和GADA模型的AMR、RAR的最大值分別為0.254 7、0.037 8，均小于基礎(chǔ)模型的AMR、RAR最小值1.260 2和0.172 9。也就是說，相較于基礎(chǔ)模型，ADA模型明顯提高了擬合精度。ADA模型之間(A1、A2、A3、A4)，各統(tǒng)計(jì)指標(biāo)值均無太大差異。以b或c為自由參數(shù)的模型(A1或A2)的評價(jià)結(jié)果要好于以a為自由參數(shù)的差分模型(A3)。GADA模型在擬合精度一樣的情況下，3參數(shù)模型(G1)的各項(xiàng)統(tǒng)計(jì)指標(biāo)均小于2參數(shù)模型(G2)。

表5 差分模型擬合結(jié)果Table 5 Fitting results for the difference model

選擇ADA模型中擬合效果靠前的A2和A4以及GADA模型的G1和G2進(jìn)行殘差分析(圖1)。A2、G2模型的殘差分布范圍在-1.5～1.5 m之間，A4模型的在-1.0～1.5 m之間，G1模型的在-1～1 m之間，但4個(gè)模型的殘差值都隨機(jī)分布在y＝0附近，沒有明顯的區(qū)別。

圖1 殘差分析Figure 1 Residual analysis diagrams

2.2 立地指數(shù)方程的選擇

利用ADA模型，以20 a為基準(zhǔn)年齡、2 m的指數(shù)級距生成立地指數(shù)曲線簇，立地指數(shù)范圍為12～24 m(圖2)。從圖2中可以看出，4個(gè)模型所構(gòu)建的模型均具有多條水平漸近線，并且能夠反映出立地質(zhì)量對優(yōu)勢高生長漸進(jìn)值和模型拐點(diǎn)的影響。但對比A2與A4模型，可以看出A4模型由于模型自身限制，在低齡階低立地指數(shù)時(shí)，會(huì)有部分?jǐn)?shù)據(jù)無法計(jì)算，出現(xiàn) “0”處不收斂現(xiàn)象，在低齡階時(shí)，A4模型呈現(xiàn)出離散狀態(tài)且在高立地指數(shù)時(shí)，優(yōu)勢高明顯偏大，并且A2與A4模型是以b和c為自由參數(shù)，所形成的立地指數(shù)曲線簇a值相同，表明所有曲線的漸進(jìn)線相同，這類模型無法反映立地質(zhì)量與優(yōu)勢高漸進(jìn)值成正比的關(guān)系。

圖2 模型立地指數(shù)曲線簇Figure 2 Site exponential curve cluster of the model

同樣，對比G1與G2模型可以看出，兩種模型均符合立地指數(shù)模型理想的生物學(xué)假設(shè)，但在齡階為2 a、立地指數(shù)為24 m時(shí)，G1模型的優(yōu)勢高為3.32 m，G2模型的為1.70 m，G2模型相較于G1更符合日本落葉松在湖北的實(shí)際生長，雖然在統(tǒng)計(jì)指標(biāo)的優(yōu)度上G2模型稍遜于G1，但從圖形表現(xiàn)以及參數(shù)數(shù)量來看，G2模型要優(yōu)于G1。綜上，最終選擇G2模型為最優(yōu)立地指數(shù)曲線方程，并進(jìn)行立地指數(shù)表的編制，方程表達(dá)式為:

式中:Ht代表年齡為t時(shí)的預(yù)測樹高(m)；Hk是已知的樹高(m)；T是基準(zhǔn)年齡(20 a)；X0和F0是引入的新參數(shù)。

2.3 立地指數(shù)表的編制與檢驗(yàn)

設(shè)定基準(zhǔn)年齡T＝20 a，基于GADA模型編制的湖北日本落葉松立地指數(shù)表見表6。

表6 基于GADA模型的日本落葉松立地指數(shù)表Table 6 Site index table of L.kaempferi based on GADA method

利用編表數(shù)據(jù)繪制立地指數(shù)曲線簇(圖3)，并利用全部86個(gè)標(biāo)準(zhǔn)地優(yōu)勢木平均高進(jìn)行落點(diǎn)檢驗(yàn)。由圖3可以看出，GADA模型曲線簇中僅有1個(gè)標(biāo)準(zhǔn)地落于曲線外，精度為98.8%。

隨機(jī)選擇30株(占總數(shù)34.9%)年齡大于基準(zhǔn)年齡20 a的日本落葉松人工林解析木數(shù)據(jù)進(jìn)行跳級檢驗(yàn)，檢驗(yàn)其各年齡樹高是否在同一指數(shù)級中，統(tǒng)計(jì)各解析木在相應(yīng)的指數(shù)級中各年齡樹高的合格數(shù)和不合格數(shù)。結(jié)果表明，有27對樹高-年齡數(shù)據(jù)存在跳級，其中25對跳1級，2對跳2級，289對不存在跳級現(xiàn)象，檢驗(yàn)合格(表7)。以α＝0.05，f＝1， 查表可得 χ2值， 得出兩種模型均存在η＜χ2，說明其對應(yīng)的立地指數(shù)表對于各優(yōu)勢解析木不同年齡的立地指數(shù)是一致的，對立地質(zhì)量的評判有一定的準(zhǔn)確性。但是根據(jù)合格數(shù)占總數(shù)的比例表明，GADA模型的合格數(shù)比例為91.46%。落點(diǎn)檢驗(yàn)和跳級檢驗(yàn)結(jié)果均能夠說明所編立地指數(shù)表預(yù)測精度高，能很好地反映湖北省日本落葉松生長實(shí)際情況。

表7 跳級檢驗(yàn)Table 7 Skip test list

將檢驗(yàn)樣本分組，計(jì)算各齡階樹高平均值，利用公式(1)進(jìn)行樹高值的轉(zhuǎn)換，再利用公式(2)～(3)計(jì)算不同齡階、不同立地指數(shù)的預(yù)測精度，結(jié)果見表8。結(jié)果表明，立地指數(shù)估測誤差介于0.22～1.39之間，且高立地指數(shù)時(shí)的估測誤差大于低立地指數(shù)；不同齡階的估測誤差介于0.25～1.20之間，總體呈現(xiàn)出中間齡階預(yù)估精度高。立地指數(shù)和不同齡階的預(yù)估精度絕大部分未超過1 m的誤差，達(dá)到了較好的預(yù)估效果。

表8 立地指數(shù)表預(yù)測精度檢驗(yàn)Table 8 Precision test for the prediction of the site index table

3 討論與結(jié)論

GADA構(gòu)建的立地指數(shù)動(dòng)態(tài)模型在滿足統(tǒng)計(jì)學(xué)和生物學(xué)要求的基礎(chǔ)上對樹高生長過程有很好的模擬效果[3，16]，國內(nèi)外學(xué)者對此均有相關(guān)研究。目前國內(nèi)對森林立地質(zhì)量評價(jià)多采用ADA代替?zhèn)鹘y(tǒng)導(dǎo)向曲線法，而GADA模型在避免ADA模型缺陷的基礎(chǔ)上，可以對基礎(chǔ)模型進(jìn)行擴(kuò)展，提高模擬精度，并且可獲得多重漸近線模型，因此使用范圍更廣[17]。通過對比ADA與GADA模型的擬合結(jié)果，可以得出GADA模型擬合精度更高、適用性更好，與曹元帥等[2]、趙磊等[18]得出的廣義差分型立地指數(shù)模型具有更高的預(yù)測性能，更適合作為立地指數(shù)模型推導(dǎo)方法的結(jié)論一致。差分模型的擬合結(jié)果表明，以b或c為自由參數(shù)的差分模型的擬合結(jié)果優(yōu)于以a為自由參數(shù)的差分模型，與相聰偉[3]得出的基礎(chǔ)模型相同時(shí)，以b或c為自由參數(shù)的差分模型擬合精度好于以a為自由參數(shù)的差分模型的結(jié)論相似，但是與段愛國等[4]得出的同一基礎(chǔ)模型，自由參數(shù)為b和a的差分模型模擬結(jié)果更好的結(jié)論不同，可能與樣本數(shù)據(jù)分布不均有關(guān)。GADA模型編制的立地指數(shù)表檢驗(yàn)結(jié)果預(yù)測精度高，落點(diǎn)檢驗(yàn)結(jié)果表明擬合精度為98.8%，略大于李忠國[12]關(guān)于恩施州日本落葉松的落點(diǎn)檢驗(yàn)精度(97.65%)，能很好地反映湖北省日本落葉松生長實(shí)際情況。

選擇立地指數(shù)表中齡階6～32 a、立地指數(shù)12～22 m部分與李忠國[12]編制的湖北省恩施州日本落葉松人工林立地指數(shù)表進(jìn)行比較，整體來看，各齡階誤差絕大部分未超過1 m，說明兩種立地指數(shù)表所反映的立地質(zhì)量十分接近。在齡階6～12 a和20～32 a時(shí)，編制的立地指數(shù)表中大部分立地質(zhì)量略高于后者，相差范圍在0.01～0.72 m；齡階14～18 a時(shí)，部分立地質(zhì)量略低于后者，但相差不大(0.01～0.06 m)。

對比馬友平等[13]有關(guān)長嶺崗日本落葉松立地指數(shù)表的研究，發(fā)現(xiàn)其研究結(jié)果呈現(xiàn)在低齡階時(shí)樹高偏高，高齡階時(shí)樹高偏低的現(xiàn)象。利用編制的立地指數(shù)表對湖北省已調(diào)查的日本落葉松人工林立地質(zhì)量進(jìn)行分析發(fā)現(xiàn)，立地質(zhì)量高(20～24 m)的林分占45.4%，立地質(zhì)量中等(16～20 m)的林分占48.9%，立地質(zhì)量較差(14～16 m)的林分占5.7%，相較于惠淑榮等[19]關(guān)于遼東地區(qū)日本落葉松立地質(zhì)量研究中得出的80%標(biāo)準(zhǔn)地適合或較適合發(fā)展日本落葉松，湖北省日本落葉松人工林立地質(zhì)量更優(yōu)，符合南方日本落葉松生長更快的研究結(jié)論。

只選擇Richards為基礎(chǔ)理論模型、1～2個(gè)自由參數(shù)的假設(shè)，利用其他基礎(chǔ)模型及自由參數(shù)假設(shè)構(gòu)建立地指數(shù)模型的效果還有待進(jìn)一步研究。

通過擬合比較，選擇Richards模型為優(yōu)勢高擬合的最優(yōu)基礎(chǔ)模型，基于最優(yōu)基礎(chǔ)模型分別構(gòu)建ADA和GADA的差分模型，通過評價(jià)和檢驗(yàn)，差分模型擬合效果優(yōu)于理論模型，差分模型之間統(tǒng)計(jì)指標(biāo)差異不大，殘差分析也沒有明顯區(qū)別。通過對兩種模型所形成的立地指數(shù)曲線簇的拐點(diǎn)和漸近線分析，最終選擇選擇參數(shù)a＝expX、c＝b1+1/X的GADA模型作為立地指數(shù)表編表方程，落點(diǎn)檢驗(yàn)精度和跳級檢驗(yàn)合格率分別為98.8%和91.5%。GADA模型編制的立地指數(shù)表精度高、適應(yīng)性好，適合用于評價(jià)湖北日本落葉松立地質(zhì)量。