何 瀟，趙嘉誠(chéng)，曹 磊，李海奎*，徐期瑚

(1.中國(guó)林業(yè)科學(xué)研究院 資源信息研究所，北京 100091：2.廣東省林業(yè)調(diào)查規(guī)劃院，廣東 廣州 510520)

森林是陸地生態(tài)系統(tǒng)最大的碳庫(kù)，其固碳能力是評(píng)價(jià)全球大氣碳收支的重要參數(shù)。作為森林固碳能力的計(jì)算基礎(chǔ)和森林生態(tài)系統(tǒng)運(yùn)轉(zhuǎn)的物質(zhì)來源的森林生物量，在固碳增匯、促進(jìn)生態(tài)系統(tǒng)碳循環(huán)、維持全球碳平衡等方面發(fā)揮著重要作用[1]。近年來，隨著各國(guó)對(duì)全球氣候變化的重視程度不斷加深，森林生物量的測(cè)定和估算已成為當(dāng)代林業(yè)生產(chǎn)和科研的焦點(diǎn)[2]。我國(guó)對(duì)森林生物量的研究?jī)?nèi)容多為利用生物量模型、生物量估算參數(shù)、3S技術(shù)等方法在大尺度上對(duì)森林生物量及碳儲(chǔ)量進(jìn)行估算[3]來評(píng)價(jià)森林的固碳能力。

然而，對(duì)森林生物量及其固碳能力的評(píng)價(jià)不僅應(yīng)體現(xiàn)在空間尺度上，更應(yīng)該將森林不同時(shí)期的固碳效率考慮在內(nèi)——在對(duì)當(dāng)前階段森林固碳率有一個(gè)準(zhǔn)確的評(píng)價(jià)的基礎(chǔ)上，對(duì)未來幾年、幾十年內(nèi)，其生物量及固碳能力的變化也應(yīng)有一個(gè)合理的估計(jì)。這種動(dòng)態(tài)變化是森林經(jīng)營(yíng)方面的重要參考指標(biāo)[4]。但是，森林生長(zhǎng)周期長(zhǎng)、長(zhǎng)期定位監(jiān)測(cè)困難等特點(diǎn)，使得如何對(duì)森林生長(zhǎng)過程進(jìn)行評(píng)價(jià)與動(dòng)態(tài)預(yù)測(cè)成為林業(yè)工作中的重點(diǎn)和難點(diǎn)[5]。建立生長(zhǎng)模型，往往是解決這一問題的有效途徑。

生長(zhǎng)模型，描述了林木個(gè)體或者林分大小隨時(shí)間的變化[6]。單木生長(zhǎng)模型是以林分中各單株林木作為基本單位，它從競(jìng)爭(zhēng)機(jī)制出發(fā)，通常用來預(yù)測(cè)林木胸徑、斷面積或樹高等的生長(zhǎng)過程[7]。單木生長(zhǎng)模型可以直接判定各株單木的生長(zhǎng)狀況，以單木模型為基礎(chǔ)可以模擬林分的生長(zhǎng)，更加詳盡地描述林分的結(jié)構(gòu)和動(dòng)態(tài)[8-9]且有較好的適用性[10]。最早的生長(zhǎng)模型主要包括邏輯斯蒂(Logistic)、坎派茲(Gompertz)、單分子式(Mitscherlich)、理查茲(Richards)等在內(nèi)的基礎(chǔ)模型，后來逐步擴(kuò)展為考慮林木之間競(jìng)爭(zhēng)的，與距離有關(guān)的生長(zhǎng)模型(distance-dependence)和與距離無關(guān)的生長(zhǎng)模型(distance-independence)[11]。近年來，盡管我國(guó)學(xué)者已對(duì)單木生長(zhǎng)模型進(jìn)行了廣泛研究[12-13]，但林分年齡(尤其是天然林)的確定難度較大且不同時(shí)期林木生物量的測(cè)定成本較高，導(dǎo)致我國(guó)主要樹種生物量生長(zhǎng)模型一直極為缺乏，但有關(guān)研究表明，年齡對(duì)單木和林分生物量影響極顯著[14]。對(duì)于已經(jīng)建立好的可靠的生長(zhǎng)模型，盡管林木的年齡無法準(zhǔn)確確定，但兩個(gè)生長(zhǎng)期之間的時(shí)間間隔是可以確定的，在森林經(jīng)營(yíng)措施不改變的條件下，已知任意時(shí)間的森林生物量均可根據(jù)生長(zhǎng)方程直接推出若干年后的生物量。但由于不同的林木，受立地條件、氣候因子、人為活動(dòng)等多種因素的影響，同一樹種的單株生物量生長(zhǎng)過程往往不盡相同。因此，對(duì)于同一樹種不同生長(zhǎng)速度，一般需要分別建模，這樣不僅工作量較大，而且使模型不相容。而在模型中加入?yún)^(qū)分生長(zhǎng)環(huán)境特征的啞變量則可以解決以上問題[15]。

楓香(Liquidambarformosana)為金縷梅科楓香屬落葉喬木，是我國(guó)閩北林區(qū)重要的鄉(xiāng)土樹種，生長(zhǎng)迅速、適應(yīng)性強(qiáng)，有較強(qiáng)的耐火性和對(duì)毒性氣體的抗性，在觀賞、藥用、用材等方面都有著重要作用且頗具造林綠化的潛力。本文基于廣東省楓香解析木，利用啞變量技術(shù)，建立楓香單木干材生物量生長(zhǎng)模型，合理地分析該樹種不同生長(zhǎng)速度下生物量生長(zhǎng)的動(dòng)態(tài)變化，以期為指導(dǎo)森林經(jīng)營(yíng)工作及全面評(píng)價(jià)森林固碳能力提供科學(xué)依據(jù)。

1 材料與方法

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

本文根據(jù)楓香在廣東省的分布，從廣州、河源、惠州、茂名、梅州、清遠(yuǎn)、韶關(guān)、云浮、肇慶 9個(gè)地區(qū)中的一類固定樣地附近，使用典型抽樣的方法，抽取40株胸徑較大、樹高接近林分平均高、干型通直飽滿、無明顯病蟲害、梢頭完整的楓香樣木進(jìn)行樹干解析。伐倒之前記錄樹種、起源(其中25株屬于天然起源，15株為人工起源)，并準(zhǔn)確測(cè)量基部位置、胸徑、樹高、枝下高、東西南北4個(gè)方向的冠幅。楓香解析木的基本統(tǒng)計(jì)信息如表1所示。之后采用破壞性方法將樣木伐倒，分干材、樹皮、樹枝、樹葉稱鮮質(zhì)量，并分別抽取樣品帶回實(shí)驗(yàn)室，在85℃恒溫下烘干至恒質(zhì)量，根據(jù)樣品鮮質(zhì)量和干質(zhì)量分別推算樣木各部分干質(zhì)量并匯總得到地上部分干質(zhì)量。根據(jù)外業(yè)測(cè)量的每株樣木根頸處、0.5/10、1/10、2/10、3/10、4/10、5/10、6/10、7/10、8/10、9/10處的帶皮直徑與去皮直徑，計(jì)算出樹干的帶皮材積與去皮材積，繼而由式(1)求得每株樣木的干材生物量材積轉(zhuǎn)換系數(shù)R：

(1)

式中，Biomassi為第i株樣木實(shí)測(cè)干材生物量，Vi為第i株樣木去皮材積。

在轉(zhuǎn)換系數(shù)確定后，從1.3 m處開始以2 m為一個(gè)區(qū)分段，在每個(gè)區(qū)分段中央截取圓盤并包括根莖及梢底圓盤，將圓盤拋光后查數(shù)年輪，以5 a為一個(gè)齡級(jí)，分別量取東西和南北2個(gè)方向上各齡級(jí)的直徑，取其平均值作為該齡級(jí)的直徑，之后根據(jù)中央斷面區(qū)分求積式計(jì)算各個(gè)齡級(jí)的材積，假設(shè)不同齡級(jí)的材積與生物量的轉(zhuǎn)換系數(shù)相同，則可以通過式(2)求得每株樣木各齡級(jí)生物量Biomass：

Biomassit=Ri·Vit(i=1,2,3,…,40)

(2)

式中，Biomassit為第i株樣木t年齡干材生物量，Vit為第i株樣木t年齡去皮材積。

表1 楓香解析木基本統(tǒng)計(jì)量

1.2 方法

1.2.1 基礎(chǔ)模型 本文選取Logistic及其變形式、Gompertz和Richards這4種符合生物學(xué)意義且被國(guó)內(nèi)外學(xué)者廣泛使用[16-18]的生長(zhǎng)方程作基礎(chǔ)模型，其模型形式分別為：

Logistic方程：

(3)

Logistic方程變形式：

(4)

Gompertz方程：

B=A·e-me-rt+ε(A,m,r>0)

(5)

Richards方程：

B=A·(1-e-rt)c+ε(A,r,c>0)

(6)

式中，B為林木生物量；t為林木年齡；A為林木生長(zhǎng)的最大值參數(shù)，即A=Bmax；m是與初始值有關(guān)的參數(shù)；r是內(nèi)稟增長(zhǎng)率(最大生長(zhǎng)速率)參數(shù)；c，b1，b2為模型參數(shù)；ε為誤差項(xiàng)。

基礎(chǔ)模型中的參數(shù)A的確定，首先對(duì)所有樣本擬合生物量與胸徑、樹高的二元模型(B=a·DbHc)，然后帶入已知較大的胸徑以及中國(guó)植物志[19]記載的楓香最高樹高，得到理論最大生長(zhǎng)量A，最后估計(jì)其他參數(shù)。

1.2.2 啞變量模型 天然林與人工林的具有不同的立地條件、氣候條件，在人工林中林木還受到人為經(jīng)營(yíng)措施的影響，使得不同起源的生物量生長(zhǎng)情況有所不同。同一樹種即使在相同立地條件下，不同樹木之間生長(zhǎng)速度也會(huì)有所差異，這是受個(gè)體差異所影響的，林木受到氣候、土壤質(zhì)量、水分和溫度影響時(shí)，林木的生長(zhǎng)速度差異將更明顯。本文用調(diào)查時(shí)實(shí)際測(cè)量的各樣本樹高與徑階的關(guān)系，將所有樣木劃分為 “高”(n=11)、“中”(n=17)和“低”(n=12)3個(gè)立地等級(jí)，分等級(jí)建立生物量生長(zhǎng)模型。劃分具體步驟為：首先確定每株樣木所屬?gòu)诫A，然后用每株樣木的實(shí)測(cè)樹高除以該樣木的所屬?gòu)诫A，在同一徑階內(nèi)，把該比值的最大(小)值作為最大(小)區(qū)間的中間值，以此劃分3個(gè)區(qū)間，最后在每個(gè)徑階內(nèi)把樣木分為3個(gè)立地等級(jí)，得到3個(gè)立地等級(jí)。

從擬合出的4種基礎(chǔ)方程中選出1個(gè)最優(yōu)方程，在參數(shù)A上引入代表?xiàng)飨闫鹪吹膯∽兞?，在參?shù)R上引入?yún)^(qū)分不同生長(zhǎng)速度的啞變量，使得有不同起源和立地等級(jí)的楓香的生物量生長(zhǎng)用一個(gè)生物量生長(zhǎng)模型來表示，有利于減少工作量?；谧顑?yōu)基礎(chǔ)方程，分起源和立地等級(jí)建立啞變量模型。以Richards方程為例，啞變量模型形式如下：

B=(A+A0f0)·(1-e-(r+r1f1+r2f2)t)C+ε

(7)

式中，A0表示天然林與人工林下生物量生長(zhǎng)的最大值之差；ri(i=1,2)代表中、低立地等級(jí)與高立地等級(jí)下生物量的生長(zhǎng)速率之差；C為不同生長(zhǎng)速度下模型形狀參數(shù)的差異；f0為區(qū)分起源的啞變量，天然林下f0=0，人工林下f0=1；fi(i=1，2)為區(qū)分立地等級(jí)的啞變量，在高立地等級(jí)下f1=0，f2=0；中立地等級(jí)下f1=1，f2=0；低立地等級(jí)下f1=0，f2=1；ε為誤差項(xiàng)。其他變量含義同上。

1.2.3 模型評(píng)價(jià) 基于Forstat[20]和R軟件的nls和nlme[21]包進(jìn)行模型擬合并采用以下4個(gè)指標(biāo)進(jìn)行模型評(píng)價(jià)[22]，包括決定系數(shù)R2，均方根誤差RMSE，總相對(duì)誤差TRE和平均預(yù)估誤差MPE。計(jì)算公式分別為：

(8)

(9)

(10)

(11)

對(duì)于包含啞變量的參數(shù)的選擇及設(shè)置問題，一方面從統(tǒng)計(jì)學(xué)上對(duì)模型參數(shù)的t檢驗(yàn)是否顯著來考慮。首先將每個(gè)參數(shù)都設(shè)置啞變量進(jìn)行模型擬合，如若在顯著水平α=0.05下不顯著(P>0.05)則不在此參數(shù)上引入啞變量，調(diào)整后重新擬合，直到所有參數(shù)全部顯著為止。另一方面，用非線性額外平方和方法[23](non-linear extra sum of squares method)檢驗(yàn)啞變量的引入是否有意義。該檢驗(yàn)需要對(duì)全模型(full model)和刪減模型(reduced model)進(jìn)行擬合。全模型是指通過擴(kuò)充啞變量形式的參數(shù)，實(shí)現(xiàn)對(duì)不同生長(zhǎng)速度下都有各自對(duì)應(yīng)參數(shù)的完整模型。而刪減模型是全模型的一部分，它嵌套在全模型中，“刪減”掉一些沒有意義的啞變量，降低了模型參數(shù)個(gè)數(shù)，但幾乎不降低模型精度，使模型更佳簡(jiǎn)化[24]。這一過程通過對(duì)統(tǒng)計(jì)量F*進(jìn)行檢驗(yàn)來完成：

(12)

式中，SSER是刪減模型的殘差平方和，SSEF為全模型的殘差平方和，dfR和dfF分別為刪減模型和全模型的自由度。

2 結(jié)果與分析

2.1 最優(yōu)基礎(chǔ)模型

表2是4種基礎(chǔ)模型的擬合參數(shù)和指標(biāo)評(píng)價(jià)結(jié)果，最大生長(zhǎng)量A由胸徑、樹高的二元模型生物量模型(B=0.0214D1.78H1.12)得到。可以看出，所有模型各參數(shù)全部極顯著(P<0.01)。Logistic方程式(3)給出的最大生長(zhǎng)速率(0.04)大于Gompertz(式(5)，0.02)和Richards方程(式(6)，0.01)。Richards方程的R2達(dá)到0.6以上，其次是Gompertz方程，Logistic及其變形式R2的最小，僅0.5；RMSE也顯示出同樣的規(guī)律，Richards方程優(yōu)于Gompertz方程，Logistic及其變形式的RMSE最大，為57.33 kg； Richards和Gompertz方程的TRE<10%，而Logistic及其變形式的TRE最大，為11.80%，所有模型給出的TRE均<0，說明4種生長(zhǎng)曲線對(duì)樣木的生物量估計(jì)整體偏高；所有方程的MPE差別不大，均處在5%左右。從4個(gè)評(píng)價(jià)指標(biāo)上看，Richards方程在所有評(píng)價(jià)指標(biāo)上效果最優(yōu)，擬合結(jié)果較好，其次是Gompertz方程，而Logistic及其變形式的擬合效果較差，且二者在擬合效果上幾乎無任何差異，4種基礎(chǔ)生長(zhǎng)模型的擬合效果如圖1所示。綜合考慮，選擇Richards方程為最優(yōu)基礎(chǔ)模型。

表2 基礎(chǔ)生長(zhǎng)模型參數(shù)和評(píng)價(jià)結(jié)果

2.2 啞變量模型的評(píng)價(jià)

基于最優(yōu)基礎(chǔ)模型構(gòu)建啞變量模型，表3為啞變量模型各擬合參數(shù)和指標(biāo)評(píng)價(jià)結(jié)果，立地等級(jí)從高到低，帶有啞變量的Richards方程式(7)給出的生長(zhǎng)速率參數(shù)估計(jì)值也逐漸下降。與基礎(chǔ)模型相比，啞變量模型R2提高了11%，RMSE降低了7.95 kg，TRE也下降了1.85%，該值>0，說明啞變量模型對(duì)樣木的生物量估計(jì)偏低，MPE指標(biāo)進(jìn)一步減少。綜合來看，啞變量模型明顯優(yōu)于Richards基礎(chǔ)模型，說明將啞變量引入方程中，細(xì)化了模型的使用范圍，大幅提高了模型擬合效果。

從模型擬合預(yù)測(cè)值對(duì)實(shí)測(cè)值散點(diǎn)圖(圖1)中也可以看出，基礎(chǔ)模型(上半幅)的擬合值相對(duì)于理想標(biāo)準(zhǔn)線較為離散，年齡較小的林木干材生物量預(yù)測(cè)值幾乎全都處于標(biāo)準(zhǔn)線的上方，幾棵大樹的預(yù)測(cè)值距離標(biāo)準(zhǔn)線也比較遠(yuǎn)，整體擬合效果比較差；啞變量模型(下半幅)則較為緊密，年齡較小的林木均勻分布在標(biāo)準(zhǔn)線的兩側(cè)，擬合效果很好，僅有幾株年齡較大的林木擬合效果較差。啞變量模型和基礎(chǔ)模型都有一些偏離標(biāo)準(zhǔn)線預(yù)測(cè)值較大的值，但基礎(chǔ)模型的偏離程度遠(yuǎn)遠(yuǎn)大于啞變量模型。

2.3 不同起源與立地等級(jí)下的楓香干材生物量生長(zhǎng)差異

由表4可以看出，不同起源和立地等級(jí)的生長(zhǎng)特征有所差異，平均生長(zhǎng)量和連年生長(zhǎng)量的最大年齡均隨著立地等級(jí)的不斷降低而增加，數(shù)量成熟齡差距明顯；在立地等級(jí)相同時(shí)，連年生長(zhǎng)量與平均生長(zhǎng)量的最大值由于起源的不同也出現(xiàn)了差異，而這種差異在立地等級(jí)較高時(shí)表現(xiàn)的更為明顯；同時(shí)天然起源最大生長(zhǎng)量要大于人工起源的最大生長(zhǎng)量，這點(diǎn)由啞變量模型中A0的估計(jì)值<0也可以發(fā)現(xiàn)。圖3給出了天然林高立地等級(jí)的楓香干材生物量連年生長(zhǎng)量與平均生長(zhǎng)量的關(guān)系。在幼齡和中齡階段干材迅速生長(zhǎng)，連年增長(zhǎng)量與平均生長(zhǎng)量都隨年齡的增加而不斷增大，生長(zhǎng)率持續(xù)上升，但連年生長(zhǎng)量增加速度較快，最大連年生長(zhǎng)量(13.3 kg)在23 a處到達(dá)。連年生長(zhǎng)量與平均生長(zhǎng)量于38 a相交，此時(shí)平均生長(zhǎng)量(8.4 kg)達(dá)到最大值，之后生長(zhǎng)逐漸趨于緩慢，連年生長(zhǎng)量小于平均生長(zhǎng)量。

圖1 4種基礎(chǔ)生長(zhǎng)模型擬合效果

方程參數(shù)估計(jì)值P值R2RMSE/kgTRE/%MPE/%啞變量模型(式7)A455<0.01A0-166<0.01r0.071 5<0.01r10.019 3<0.01r20.025 4<0.01c5.148 1<0.010.7243.043.704.15

圖2 基礎(chǔ)模型和啞變量模型的干材生物量預(yù)測(cè)值對(duì)實(shí)測(cè)值散點(diǎn)

3 結(jié)論與討論

通過對(duì)不同起源、地區(qū)抽取的楓香樣木進(jìn)行樹干解析，建立了廣東省楓香干材生物量生長(zhǎng)模型。樣木之間生長(zhǎng)過程的差異主要是樣木所在林分的起源、立地質(zhì)量及其周圍環(huán)境等因素綜合影響的結(jié)果。人工林中由于存在各種撫育經(jīng)營(yíng)措施，使林木更早地進(jìn)入速生期，生長(zhǎng)速率應(yīng)高于同時(shí)期的天然林。增加起源和立地等級(jí)的啞變量模型的最大生長(zhǎng)速率參數(shù)差異顯著，與黃勇來[25]得到的楓香人工林較天然林生長(zhǎng)更為迅速結(jié)論一致。陳波濤[26]等通過對(duì)貴州省楓香樹的調(diào)查，發(fā)現(xiàn)林木由于受到氣候的影響，在土壤條件基本一致的情況下，低海拔地區(qū)下的生長(zhǎng)量大于高海拔地區(qū)，高積溫地區(qū)大于低積溫地區(qū)。這些結(jié)果也印證了劃分起源和立地等級(jí)來建模的合理性。我國(guó)南方楓香天然林一般于21～40a生長(zhǎng)迅速，51～70 a成熟；人工林于11～20 a生長(zhǎng)迅速，31～50 a成熟[27]。而本研究顯示，廣東省楓香干材生物量于23～36 a迅速生長(zhǎng)，在38～59 a成熟。楓香干材生物量在天然高立地等級(jí)的條件下，連年生長(zhǎng)量與平均生長(zhǎng)量于38 a相交，這與周火明[28]等通過對(duì)楓香天然次生林研究，得到楓香材積的連年生長(zhǎng)量在10～40 a時(shí)一直呈現(xiàn)較快的增長(zhǎng)趨勢(shì)，而在50 a時(shí)還未達(dá)到平均生長(zhǎng)量最大值的結(jié)果相近。

表4 不同起源和立地等級(jí)下的生物量最大生長(zhǎng)年齡和生長(zhǎng)量

圖3 天然林高立地質(zhì)量干材生物量連年生長(zhǎng)量與平均生長(zhǎng)量的關(guān)系

選用4種國(guó)內(nèi)外常用的理論生長(zhǎng)方程，擬合廣東省楓香干材生物量的生長(zhǎng)過程， Richards方程在所有基礎(chǔ)模型中效果最優(yōu)，擬合結(jié)果最好?；赗ichards方程的啞變量模型，啞變量模型則進(jìn)一步提升了精度，更加科學(xué)合理地描述了廣東省楓香干材生物量的生長(zhǎng)過程。使用啞變量模型進(jìn)一步分析不同起源和立地等級(jí)下楓香干材生物量的生長(zhǎng)過程：干材生物量生長(zhǎng)的最大年齡隨著立地等級(jí)的不斷降低而增加，數(shù)量成熟齡差距明顯；在立地等級(jí)相同而起源不同時(shí)，干材生物量的最大連年生長(zhǎng)量與平均生長(zhǎng)量不同，在高立地等級(jí)較時(shí)差異更大。以天然林高立地等級(jí)為例，其干材生物量的最大連年生長(zhǎng)量(13.3 kg)和最大平均生長(zhǎng)量(8.4 kg)分別在23 a和38 a處到達(dá)，此后生長(zhǎng)逐漸緩慢。

本模型的解析木建模樣本取自整個(gè)廣東省，適用范圍為整個(gè)廣東省，但由于生長(zhǎng)模型缺乏對(duì)立地條件的輸入，所以即使啞變量模型能夠提高模型精度，但72%的解釋力也不是很理想，此模型在某一林分下使用可能會(huì)產(chǎn)生較大誤差，若能增加建模樣本數(shù)量，獲得更多解析木，則模型的解釋力會(huì)隨之增強(qiáng)，預(yù)估精度將會(huì)有所提高，減小誤差。