黃興召 豆玉萍 黃慶豐 陳東升

(1. 安徽農(nóng)業(yè)大學(xué)林學(xué)與園林學(xué)院，合肥 230036；2. 中國(guó)林業(yè)科學(xué)研究院林業(yè)研究所，北京 100091)

立木地上生物量是森林生物量的重要組成部分，如何準(zhǔn)確地估算立木地上生物量是研究森林生態(tài)系統(tǒng)物質(zhì)循環(huán)的基礎(chǔ)[1?2]。異速生物量模型作為預(yù)估生物量的常用方法，已被廣泛應(yīng)用于各類研究中[3?6]。異速生物量方程參數(shù)a和b隨著研究區(qū)域和樹種的變化而變化，通過(guò)收獲法，利用異速生物量方程預(yù)估生物量時(shí)，只有單獨(dú)建模才能提高預(yù)估精度[5?8]。貝葉斯方法作為其中的一種擬合方法，通過(guò)收集已建立方程的參數(shù)信息，將該信息轉(zhuǎn)化為參數(shù)的概率分布，以參數(shù)的均值和協(xié)方差矩陣為先驗(yàn)信息進(jìn)行擬合，越來(lái)越受到人們的關(guān)注[9?12]。但是，利用先驗(yàn)信息對(duì)貝葉斯方法進(jìn)行擬合時(shí)，研究者收集參數(shù)a和b的先驗(yàn)信息具有隨機(jī)性和異質(zhì)性，由此對(duì)模型擬合和預(yù)估效果產(chǎn)生的影響尚未得到驗(yàn)證，也沒(méi)有闡述如何選取參數(shù)a和b的先驗(yàn)信息來(lái)提升方程的擬合和預(yù)估效果。

收集已發(fā)表文獻(xiàn)異速生物量方程參數(shù)a和b的信息，是分析先驗(yàn)信息異質(zhì)性對(duì)立木生物量預(yù)估影響的基礎(chǔ)。Ter-Mikaelian和Korzukhin[13]建立了北美65個(gè)樹種的生物量預(yù)估方程、Zianis和Mencuccini[3]建立了279個(gè)異速生物量預(yù)估方程、Henry等[14]收集了850個(gè)立木及各組分的生物量預(yù)估方程、Zapata-Cuartas等[8]收集了134個(gè)的立木地上異速生物量方程、Rojas-García等[15]收集了487個(gè)立木及各組分的生物量預(yù)估方程。這些方程都是通過(guò)常規(guī)方法擬合得到參數(shù)a和b的值，為貝葉斯方法提供了豐富的先驗(yàn)信息。分析發(fā)現(xiàn)，立木地上部分生物量預(yù)估方程參數(shù)b與緯度呈顯著相關(guān)，且在不同屬之間存在顯著差異[8,13]。研究者在利用貝葉斯方法預(yù)估立木地上部分生物量的過(guò)程中，討論了樣本大小對(duì)方程擬合效果和預(yù)估效果的影響，缺乏反映樹種特征的先驗(yàn)信息的擬合比較[8,10,16]。樹木生長(zhǎng)與生物量的積累受環(huán)境[17,18]和生物學(xué)特性[19?20]的影響，在方程擬合時(shí)先驗(yàn)信息也會(huì)受到影響[21]。因此，檢驗(yàn)參數(shù)a和b先驗(yàn)信息的氣候區(qū)和生物學(xué)特性對(duì)立木生物量預(yù)估的影響顯得尤為重要。

通過(guò)分析已發(fā)表文獻(xiàn)關(guān)于立木地上部分異速生物量方程的參數(shù)a和b信息，將參數(shù)的先驗(yàn)信息按主要?dú)夂騾^(qū)和樹種進(jìn)行分類，以日本落葉松（Larix kaempferi）立木地上部分生物量為擬合數(shù)據(jù)，利用貝葉斯方法擬合異速生物量方程，比較先驗(yàn)信息的氣候差異和樹種差異對(duì)立木生物量擬合和預(yù)估的影響，探討先驗(yàn)信息對(duì)立木生物量預(yù)估的作用規(guī)律，為利用貝葉斯方法擬合異速生物量方程時(shí)選擇先驗(yàn)信息提供了理論和實(shí)踐的依據(jù)。

1 研究方法

1.1 異速生物量方程

異速生物量方程作為預(yù)估立木生物量的常用模型，其方程表達(dá)式為：

其中，Mi為第i株樣本的地上部分生物量；Di為第i株樣本的胸徑；a和b為方程參數(shù)；ei為誤差項(xiàng)。

貝葉斯方法將公式（1）中的參數(shù)a、b定義為θ，需要知道θ的概率密度函數(shù)。接著，根據(jù)先驗(yàn)信息確定參數(shù)a和b的先驗(yàn)分布 π(θ)，再結(jié)合樣本信息可得出后驗(yàn)分布 π(θ|D)。 π(θ|D)的計(jì)算公式為：

1.2 先驗(yàn)信息

利用貝葉斯方法擬合異速生物量方程時(shí)，以參數(shù)a、b的均值和協(xié)方差矩陣為先驗(yàn)信息。為研究不同先驗(yàn)信息對(duì)異速生物量方程擬合效果的影響，本文通過(guò)分析[8,13]立木地上部分異速生物量方程參數(shù)a和b的數(shù)據(jù)集，按照氣候和屬進(jìn)行分類（表1）。以主要?dú)夂騾^(qū)為單位，將參數(shù)a和b分為寒帶、溫帶和熱帶三種類型，分別計(jì)算它們的均值和協(xié)方差矩陣。同時(shí)，以屬為單位，將參數(shù)a和b分為櫟屬（Quercus）、樺木屬（Betula）、楊屬（Populus）、槭屬（Acer）、桉屬（Eucalyptus）和松屬（Pinus）等六個(gè)屬，計(jì)算其均值和協(xié)方差矩陣。以氣候區(qū)和屬為單位進(jìn)行分類，是進(jìn)一步研究先驗(yàn)信息異質(zhì)性對(duì)異速生物量方程擬合效果影響的關(guān)鍵。

表1 異速生物量方程先驗(yàn)信息分類Table1Priorinformationclassificationofallometric biomassequation

2 方程擬合

2.1 擬合數(shù)據(jù)

以遼寧省清源縣大孤家林場(chǎng)（42°16?～42°23? N，124°48?～124°55? E）日本落葉松（Larixkaempeferi）立木地上部分生物量數(shù)據(jù)為擬合數(shù)據(jù)[22]，該區(qū)屬于典型的溫帶季風(fēng)氣候?；跇拥卣{(diào)查，在標(biāo)準(zhǔn)地附近選取標(biāo)準(zhǔn)木，利用收獲法測(cè)定樣木的樹干、樹皮、樹枝和樹葉各組分的生物量，求出立木地上部分生物量。本數(shù)據(jù)由幼齡林、中齡林、近熟林以及成熟林4個(gè)林齡段共45株立木的胸徑和地上部分生物量組成（表2）。本文在建模過(guò)程中，利用無(wú)重復(fù)抽樣，從45株總樣木中抽取30株樣木，進(jìn)行1 000次無(wú)重復(fù)擬合，每次擬合未被抽取的15株作為驗(yàn)證數(shù)據(jù)來(lái)檢驗(yàn)預(yù)估效果。

表2 日本落葉松樣地信息Table2PlotsinformationofLarixkaempferiplantations

2.2 預(yù)估效果與預(yù)估精度檢驗(yàn)

使用不同的先驗(yàn)信息，利用貝葉斯方法擬合異速生物量方程時(shí)，采用決定系數(shù)（R2）、平均偏差（Mean Bias , MB）和平均均方根差（Mean Root Mean Square Error, MRMSE）來(lái)檢驗(yàn)方程的擬合結(jié)果和預(yù)估效果。R2、MB和MRMSE的表達(dá)式分別為：

其中，yi為實(shí)際觀測(cè)值；yˉ為樣本平均值；y?i為模型預(yù)測(cè)值；n為樣本數(shù)。數(shù)據(jù)統(tǒng)計(jì)分析和繪圖使用R（Development Core Team, 2009）和Excel軟件。以不同區(qū)域和不同樹種參數(shù)a、b的均值以及協(xié)方差矩陣為先驗(yàn)信息，應(yīng)用R的MCMCglmm包[23]對(duì)異速生物量方程分別進(jìn)行擬合。

3 結(jié)果與分析

3.1 先驗(yàn)信息氣候區(qū)差異對(duì)異速生物量方程擬合的影響

為檢驗(yàn)先驗(yàn)信息區(qū)域差異對(duì)異速生物量方程擬合的影響，分別以寒帶、溫帶、熱帶以及總數(shù)據(jù)集參數(shù)a、b的均值及協(xié)方差矩陣為先驗(yàn)信息，利用無(wú)重復(fù)抽樣，從45株總樣木中抽取30株樣木，進(jìn)行1 000次無(wú)重復(fù)擬合（圖1）。寒帶和熱帶先驗(yàn)信息擬合后參數(shù)a、b值的變化范圍相似；溫帶和總數(shù)據(jù)集先驗(yàn)信息擬合后參數(shù)a、b值的變化范圍相似。以熱帶和寒帶為先驗(yàn)信息進(jìn)行1 000次無(wú)重復(fù)擬合，參數(shù)a的平均值均為?2.794，95%的置信區(qū)間分別為（?3.185, ?2.386）和（?3.197,?2.368）；參數(shù)b的平均值均分別為2.640和2.638，95%的置信區(qū)間分別為（2.484, 2.788）和（2.476, 2.792）。以溫帶和總數(shù)據(jù)集為先驗(yàn)信息進(jìn)行1 000次無(wú)重復(fù)擬合，參數(shù)a的平均值分別為?2.836和?2.826，95%的置信區(qū)間分別為（?3.240,?2.418）和（?3.228, ?2.407）；參數(shù)b的平均值均分別為2.655和2.651，95%的置信區(qū)間分別為（2.496, 2.809）和（2.491, 2.805）（表3）。這說(shuō)明使用不同區(qū)域的先驗(yàn)信息擬合異速生物量方程會(huì)影響參數(shù)a、b值的穩(wěn)定性。

圖1 不同氣候區(qū)先驗(yàn)信息預(yù)估參數(shù)a和b的擬合結(jié)果Fig.1Simulationresultsfortheestimationofthe parametersaandbusingdifferentprior informationofregion

決定系數(shù)R2反映方程的擬合效果（表3），寒帶、溫帶、熱帶以及總數(shù)據(jù)集四種先驗(yàn)信息1 000次無(wú)重復(fù)擬合，方程的平均決定系數(shù)均為0.934。平均偏差MB和平均均方根差MRMSE反映了方程的預(yù)估效果，四種先驗(yàn)信息1 000次無(wú)重復(fù)擬合，平均偏差分別為?1.392、?0.597、?0.929、?0.715，平均均方根差分別為0.347、0.344、0.347、0.345。由此可知，區(qū)域先驗(yàn)信息主要影響方程的預(yù)估效果。

表3 不同氣候區(qū)先驗(yàn)信息擬合異速生物量方程的結(jié)果對(duì)比Table3Resultcomparisonallometricbiomassequationsindifferentpriorinformationofregions

為了進(jìn)一步比較寒帶、溫帶、熱帶以及總數(shù)據(jù)集四種先驗(yàn)信息對(duì)方程預(yù)估效果的影響，對(duì)BE和RMSE分別進(jìn)行方差分析（圖2）；總數(shù)據(jù)集先驗(yàn)信息的MB和MRMSE與溫帶的MB和MRMSE無(wú)顯著差異；寒帶與熱帶的MB和MRMSE也無(wú)顯著差異；其他先驗(yàn)信息的MB和MRMSE之間存在顯著差異（P<0.01）?？倲?shù)據(jù)集和溫帶先驗(yàn)信息的擬合效果優(yōu)于寒帶和熱帶，驗(yàn)證數(shù)據(jù)集所在區(qū)域?yàn)闇貛?，因此，利用貝葉斯方法擬合異速生物量方程時(shí)，應(yīng)充分考慮先驗(yàn)信息的氣候區(qū)特征。

圖2 不同氣候區(qū)先驗(yàn)信息擬合的平均偏差和平均均方根差比較Fig.2Simulationresultsforthemeanbiasandmean rootmeansquareerrorusingdifferent priorinformationofregion

3.2 先驗(yàn)信息屬差異對(duì)異速生物量方程擬合的影響

以櫟屬、樺木屬、楊屬、槭屬、桉屬、松屬和總數(shù)據(jù)集為先驗(yàn)信息，在樣木數(shù)量為30的情況下，也進(jìn)行1 000次無(wú)重復(fù)擬合（圖3）。桉屬、櫟屬、松屬以及總數(shù)據(jù)集的參數(shù)a的變化范圍相似，其平均值分別為?2.859、?2.866、?2.855、?2.826；樺木屬和楊屬的參數(shù)a的變化范圍相似，其平均值分別為?2.728和?2.712；槭屬的參數(shù)變化范圍尤為突出，其平均值為?2.492。參數(shù)b的變化范圍也表現(xiàn)出相同的規(guī)律，桉屬、櫟屬、松屬以及總數(shù)據(jù)集參數(shù)b的平均值分別為2.663、2.666、2.662、2.651；樺木屬和楊屬參數(shù)b的平均值分別為2.616和2.605；槭屬參數(shù)b的平均值為2.528。

圖3 不同屬先驗(yàn)信息預(yù)估參數(shù)a和b的擬合結(jié)果Fig.3Simulationresultsfortheestimationofthe parametersaandbusingdifferent priorinformationofgenus

以櫟屬、樺木屬、楊屬、槭屬、桉屬、松屬和總數(shù)據(jù)集為先驗(yàn)信息，進(jìn)行異速生物量方程擬合（表4）表4。槭屬的平均決定系數(shù)為0.935，其他屬的平均決定系數(shù)均為0.934。櫟屬、樺木屬、楊屬、槭屬、桉屬、松屬和總數(shù)據(jù)集的平均偏差分別為?0.205、?1.287、?2.926、?3.726、?0.367、?0.494、?0.715；平均均方根差分別為0.343、0.352、0.354、0.391、0.344、0.343、0.345。由此可知，樹種先驗(yàn)信息主要影響方程的預(yù)估效果。

表4 不同屬先驗(yàn)信息擬合異速生物量方程的結(jié)果對(duì)比Table4Resultcomparisonallometricbiomassequationsindifferentpriorinformationofgenus

為了進(jìn)一步比較六個(gè)屬以及總數(shù)據(jù)集七種先驗(yàn)信息對(duì)方程預(yù)估效果的影響，對(duì)MB和MRMSE分別進(jìn)行方差分析（圖4）。櫟屬和桉屬、松屬和桉屬之間的MB無(wú)顯著差異，其他屬之間的MB均存在顯著差異（P<0.01），即櫟屬和桉屬為先驗(yàn)信息的MB>松屬>總數(shù)>樺木屬>楊屬>槭屬。櫟屬和桉屬、松屬和桉屬、總數(shù)據(jù)集和桉屬之間的MRMSE無(wú)顯著差異，其他屬之間的MRMSE存在顯著差異（P<0.01），即松屬和櫟屬為先驗(yàn)信息的MRMSE>總數(shù)據(jù)集和桉屬>樺木屬>楊屬>槭屬。因此，樹種的生物學(xué)特性也是影響生物量預(yù)估的重要因子。

圖4 不同屬先驗(yàn)信息擬合的平均偏差和平均均方根差比較Fig.4Simulationresultsforthemeanbiasandmean rootmeansquareerrorusingdifferent priorinformationofgenus

4 討論

4.1 區(qū)域差異對(duì)異速生物量方程擬合的影響

利用貝葉斯方法擬合生物量方程，先驗(yàn)信息尤為重要[8?9,16]，通過(guò)收集立木地上部分異速生物量方程參數(shù)信息，分析表明參數(shù)a和b符合二元正態(tài)分布。擬合異速生物量方程過(guò)程中，貝葉斯原理使用先驗(yàn)信息，而傳統(tǒng)方法只能使用立木數(shù)據(jù)。因此，貝葉斯方法可以實(shí)現(xiàn)少量的樣木高精度預(yù)估生物量，黃興召等[10]和Zell等[16]均證明，當(dāng)樣木的數(shù)量小于50株時(shí)，貝葉斯方法的預(yù)估效果明顯優(yōu)于傳統(tǒng)方法。立木地上部分異速生物量方程參數(shù)a和b的值有明顯的區(qū)域特征，參數(shù)b與緯度呈線性關(guān)系[10]。本文將參數(shù)a和b劃分為寒帶、溫帶和熱帶三個(gè)主要?dú)夂騾^(qū)，基于信息量度的相似概率分布[22,24?25]，以三個(gè)氣候區(qū)和總數(shù)據(jù)集參數(shù)a、b的均值和協(xié)方差矩陣為先驗(yàn)信息，利用無(wú)重復(fù)抽樣，從總樣木中抽取30株立木生物量數(shù)據(jù)，進(jìn)行1 000次無(wú)重復(fù)擬合，得出每次擬合結(jié)果中參數(shù)a、參數(shù)b、R2、MB和MRMSE的值。

將寒帶、溫帶和熱帶三個(gè)主要?dú)夂騾^(qū)先驗(yàn)信息的擬合結(jié)果與總數(shù)據(jù)集先驗(yàn)信息的擬合結(jié)果進(jìn)行比較，發(fā)現(xiàn)R2均為0.934無(wú)顯著差異，說(shuō)明先驗(yàn)信息的氣候差異對(duì)異速生物量方程的擬合效果無(wú)影響。MB和MRMSE差異分析表明，先驗(yàn)信息的氣候差異對(duì)異速生物量方程的預(yù)估效果有顯著影響。本文中擬合數(shù)據(jù)的研究地點(diǎn)屬于典型的溫帶季風(fēng)氣候，溫帶的先驗(yàn)信息擬合效果（顯著）優(yōu)于寒帶和熱帶，與總數(shù)據(jù)集無(wú)差異，驗(yàn)證了先驗(yàn)信息的氣候差異對(duì)立木生物量預(yù)估的影響，使用立木生長(zhǎng)主要?dú)夂騾^(qū)的先驗(yàn)信息，有助于提升立木生物量的預(yù)估精度。

通過(guò)模型擬合，證實(shí)了先驗(yàn)信息的氣候區(qū)差異性對(duì)立木生物量預(yù)估的影響。氣候是環(huán)境因子最直接的外在表現(xiàn)，其綜合了降雨、光照、地形地貌等環(huán)境因子，決定了森林生物量和生產(chǎn)力的大小[26?27]。利用異速生物量方程進(jìn)行預(yù)估時(shí)，相同胸徑不同區(qū)域的立木生物量相差達(dá)到500倍[3,28]，相同區(qū)域內(nèi)的立木生物量相差也達(dá)到99倍[14,29]。Deans等[30]和Henry等[31]分別利用異速生物量方程估算熱帶森林生物量，發(fā)現(xiàn)同一氣候區(qū)內(nèi)立木生物量預(yù)估方程也存在差異，如何通過(guò)先驗(yàn)信息解釋預(yù)估方程的差異是下一步的研究重點(diǎn)。另外，氣候區(qū)之間的過(guò)渡帶森林的生物量和生產(chǎn)力均有增加的趨勢(shì)，進(jìn)一步驗(yàn)證先驗(yàn)信息區(qū)域差異對(duì)過(guò)渡區(qū)生物量預(yù)估的影響，對(duì)闡述氣候?qū)ι稚锪亢蜕a(chǎn)力的影響有重要的作用。

4.2 樹種生物學(xué)特性對(duì)異速生物量方程擬合的影響

大量研究從針葉林和闊葉林[3,13]、天然林和人工林[13,15]、常綠林和落葉林[14,32]等森林類型來(lái)比較立木生物量和預(yù)估方程的差異。樹種之間生物量的差異多數(shù)從木材密度來(lái)分析對(duì)比不同演替過(guò)程中木材密度的差異對(duì)立木生物量的影響[14,33]。本文利用櫟屬、樺木屬、楊屬、槭屬、桉屬和松屬六個(gè)屬的先驗(yàn)信息擬合異速生物量方程，從先驗(yàn)信息的生物學(xué)特性分析方程的預(yù)估效果。結(jié)果表明，槭屬平均決定系數(shù)R2為0.935，其他的平均決定系數(shù)R2均為0.934，樹種先驗(yàn)信息的差異對(duì)方程R2無(wú)顯著影響。平均偏差MB分別為? 0.205、?1.287、?2.926、?3.726、?0.367、?0.494；平均均方根差MRMSE分別為0.343、0.352、0.354、0.391、0.344、0.343，樹種先驗(yàn)信息的差異對(duì)MB和MRMSE影響顯著（P<0.01），與先驗(yàn)信息的氣候區(qū)差異的結(jié)果一致。

黃興召等[22]研究結(jié)果表明，a和b的值在針葉樹和闊葉樹之間有明顯差異，本研究結(jié)果顯示以松屬和櫟屬為先驗(yàn)信息的擬合效果>總數(shù)據(jù)集和桉屬>樺木屬>楊屬>槭屬。櫟屬、樺木屬、楊屬、槭屬、桉屬和松屬分別屬于殼斗科（Fagaceae）、樺木科（Betulaceae）、楊柳科（Salicaceae）、槭樹科（Aceraceae）、桃金娘科（Myrtaceae）和松科（Pinaceae）。但是，落葉松屬與松屬的生物學(xué)關(guān)系>落葉松屬與桉屬>落葉松屬與槭屬>落葉松屬與楊屬>落葉松屬與櫟屬>落葉松屬與樺木屬。由此推斷，樹種先驗(yàn)信息影響異速生物量的預(yù)估效果，與樹種之間的生物學(xué)關(guān)系不一致，其作用規(guī)律有待進(jìn)一步研究。本文選擇先驗(yàn)信息的氣候區(qū)差異和樹種差異對(duì)立木生物量預(yù)估的影響，驗(yàn)證了氣候因子和生物學(xué)特性對(duì)先驗(yàn)信息的影響和方程預(yù)估的作用，但先驗(yàn)信息差異性對(duì)交錯(cuò)區(qū)立木生物量預(yù)估的影響及與樹種間生物學(xué)特性的關(guān)系，還有待進(jìn)一步研究。

5 結(jié)論

先驗(yàn)信息的區(qū)域差異和樹種差異對(duì)異速生物量方程的擬合效果R2影響不顯著，對(duì)方程的擬合效果MB和MRMSE影響顯著（P<0.01），先驗(yàn)信息的差異性影響立木地上部分生物量的預(yù)估效果。

先驗(yàn)信息區(qū)域差異對(duì)立木生物量預(yù)估的影響為總數(shù)據(jù)集和溫帶的MB和MRMSE優(yōu)于寒帶和熱帶的MB和MRMSE。擬合的立木生物量研究地點(diǎn)位于遼寧，屬于典型的溫帶季風(fēng)氣候，先驗(yàn)信息的區(qū)域作用效果與立木所在區(qū)域一致。

先驗(yàn)信息區(qū)域差異對(duì)立木生物量預(yù)估效果的影響程度為：松屬和櫟屬>總數(shù)據(jù)集和桉屬>樺木屬>楊屬>槭屬。但是，落葉松屬與松屬的生物學(xué)關(guān)系>落葉松屬與桉屬>落葉松屬與槭屬>落葉松屬與楊屬>落葉松屬與櫟屬>落葉松屬與樺木屬，先驗(yàn)信息生物學(xué)特性的預(yù)估效果表現(xiàn)規(guī)律與樹種的生物學(xué)關(guān)系不一致。