賀 鵬，張會(huì)儒，雷相東，李曉冬

（中國(guó)林業(yè)科學(xué)研究院 資源信息研究所，北京 100091）

基于泛克里格方法的樹(shù)種多樣性空間估計(jì)研究

賀 鵬，張會(huì)儒，雷相東，李曉冬

（中國(guó)林業(yè)科學(xué)研究院 資源信息研究所，北京 100091）

生物多樣性保護(hù)已成為森林可持續(xù)經(jīng)營(yíng)的一個(gè)重要目標(biāo)。生物多樣性的空間分布估計(jì)對(duì)于森林經(jīng)營(yíng)者了解生物多樣性的空間格局和變化趨勢(shì)具有重要意義。本研究主要基于汪清林業(yè)局二類調(diào)查數(shù)據(jù)和局級(jí)固定樣地?cái)?shù)據(jù)，以樹(shù)種多樣性為例，利用泛克里格方法對(duì)全局樹(shù)種多樣性進(jìn)行空間分布估計(jì)。結(jié)果表明，泛克里格方法的預(yù)測(cè)精度較高（R2=0.697 7），樹(shù)種多樣性在固定樣地的抽樣尺度下，存在一定的空間自相關(guān)性，預(yù)測(cè)的樹(shù)種多樣性最高值為1.82，最低值為0.221，空間分布具有明顯的空間異質(zhì)性，呈現(xiàn)出北部低，東部，南部和西部高的空間分布格局。研究結(jié)果為區(qū)域尺度基于固定樣地和地統(tǒng)計(jì)學(xué)的樹(shù)種多樣性空間分布估計(jì)提供了方法和參照

樹(shù)種多樣性；空間分布格局；地統(tǒng)計(jì)學(xué)；泛克里格方法

生物多樣性保護(hù)已成為森林可持續(xù)經(jīng)營(yíng)的一個(gè)重要目標(biāo)。森林經(jīng)營(yíng)措施與生物多樣性密切相關(guān)，所以森林經(jīng)營(yíng)者需要了解森林經(jīng)營(yíng)活動(dòng)后生物多樣性的變化， 確保物種和種群不會(huì)因森林采伐和更新而處于危險(xiǎn)境地， 從而對(duì)經(jīng)營(yíng)和保護(hù)做出正確決策[1]。因此準(zhǔn)確地估計(jì)生物多樣性和分析生物多樣性的空間格局對(duì)于森林經(jīng)營(yíng)者從整體上了解生物多樣性的分布格局及變化趨勢(shì)和制定適應(yīng)性經(jīng)營(yíng)策略具有重要意義。

地統(tǒng)計(jì)學(xué)的理論和方法可以定量地描述和解釋空間異質(zhì)性或者空間相關(guān)性；也可以建立各種空間預(yù)測(cè)模型，并進(jìn)行空間數(shù)據(jù)插值和估計(jì)；同時(shí)能夠?qū)臻g格局的尺度、幾何形狀、變異方向進(jìn)行定量地分析和有效的估計(jì)；有利于我們深刻地了解生命有機(jī)體（個(gè)體、種群和群落）的空間分布情況和空間異質(zhì)的機(jī)制[2]。因此地統(tǒng)計(jì)學(xué)在林業(yè)中得到了比較廣泛地應(yīng)用，許多學(xué)者經(jīng)常運(yùn)用其對(duì)森林結(jié)構(gòu)進(jìn)行分析，例如胸徑地理變異規(guī)律分析[3]，胸高斷面積空間插值和空間格局分析[4-6]，蓄積空間插值和空間格局分析[7-8]，和生物量空間插值和空間格局分析[9-10]等。

生物多樣性受到各種生物因素和非生物因素的影響，所以生物多樣性不是純隨機(jī)變量，而是區(qū)域化變量，具有隨機(jī)性和結(jié)構(gòu)性雙重屬性。因此本研究以樹(shù)種多樣性為例，基于吉林省汪清林業(yè)局局級(jí)固定樣地和二類調(diào)查數(shù)據(jù)，利用地統(tǒng)計(jì)學(xué)的原理和方法對(duì)整個(gè)林業(yè)局13個(gè)林場(chǎng)的樹(shù)種多樣性進(jìn)行空間估計(jì)。實(shí)現(xiàn)由樣地“點(diǎn)”上樹(shù)種多樣性信息擴(kuò)展到整個(gè)區(qū)域“面”上樹(shù)種多樣性信息。同時(shí)利用插值得到整個(gè)林場(chǎng)的樹(shù)種多樣性空間格局分布圖，為在區(qū)域尺度內(nèi)，基于固定樣地和地統(tǒng)計(jì)學(xué)的樹(shù)種多樣性空間分布估計(jì)提供方法和參照。

1 研究區(qū)概況

汪清林業(yè)局位于吉林省延邊自治州的東部，所處的地理坐標(biāo)為東經(jīng) 123°56′～ 131°04′，北緯43°05′～ 43°40′，東西跨度約為 85 km， 南北跨度約為60 km。該局下設(shè)13個(gè)林場(chǎng)，全局經(jīng)營(yíng)面積為304 173 hm2。汪清林業(yè)局屬于長(zhǎng)白山系的中低丘陵地區(qū)， 海拔為360～1 477 m，氣候?qū)儆跍貛Т箨懶约撅L(fēng)氣候，平均氣溫為3.9 ℃， 最低溫度為-37.5 ℃，平均降水量為547 mm，且主要集中在夏季，占全年總降水量的80%。該區(qū)森林覆蓋率極高且植被種類繁多，主要以天然林為主，占有林地面積的93.7%。而天然林中以闊葉混交林、蒙古櫟林和針闊混交林為主，這三林分占全局有林地面積的73.8%、蓄積的77.6%。主要針葉樹(shù)種有長(zhǎng)白落葉松Larix olgensis、冷杉Abies nephrolepis、紅松Pinus koraiensis、云杉Picea jezoensis等。主要闊葉樹(shù)種蒙古櫟Quercus mongolica、白樺Betula platyphylla、楊樹(shù)Populus ussuriensis、楓樺Betula costata、色木Acer mono、水曲柳Fraxinus mandshurica、黃菠蘿Phellodendron amurense等。

2 數(shù)據(jù)與方法

2.1 數(shù)據(jù)來(lái)源

研究數(shù)據(jù)來(lái)源于吉林省汪清林業(yè)局2007年二類調(diào)查數(shù)據(jù)，包括局級(jí)固定樣地?cái)?shù)據(jù)。其中固定樣地共1 391個(gè)（剔除固定樣地?cái)?shù)據(jù)中蓄積為0的樣地后），為系統(tǒng)布設(shè)樣地，間距為1 km×2 km，樣地面積為0.06 hm2，樣地空間分布圖如圖1a（圖中樣地顏色越深表示樹(shù)種多樣性越高）。固定樣地主要調(diào)查因子為樹(shù)種、胸徑、坡向、坡度、海拔等，起測(cè)徑階為5 cm。

汪清林業(yè)局林場(chǎng)圖件：林相圖及相關(guān)的資源數(shù)據(jù)庫(kù)。圖1b表示全局的純林和混交林的空間分布圖。

圖1 樣地分布（a）和純林及混交林分布（b）Fig. 1 Spatial distribution of plots (a) and distribution of forest styles (b)

2.2 樣地樹(shù)種多樣性指標(biāo)計(jì)算及統(tǒng)計(jì)

根據(jù)固定樣地每木檢尺的數(shù)據(jù)，選取Shanonwiener指數(shù)（公式1）計(jì)算每個(gè)樣地的樹(shù)種多樣性指數(shù)。計(jì)算后樣地樹(shù)種多樣性指數(shù)的統(tǒng)計(jì)描述如表1，同時(shí)隨機(jī)預(yù)留出10%的樣地用來(lái)檢驗(yàn)方法的預(yù)測(cè)精度。

表1 固定樣地樹(shù)種多樣性指數(shù)統(tǒng)計(jì)Table 1 Statistical characteristics of tree species diversity indexes for studied area

式（1）中HSpecies表示樹(shù)種多樣性指數(shù)，pi表示各樹(shù)種胸高斷面積占林分胸高斷面積的比例。

2.3 泛克里格方法

在泛克里格方法（Universal Kriging）中， 預(yù)測(cè)值z(mì)(s0)由兩部分組成。一部分是m(s0)估計(jì)值，通過(guò)不同環(huán)境解釋變量對(duì)主變量的回歸估計(jì)，s0表示待估測(cè)點(diǎn)位置。第二部分是通過(guò)克里格插值殘差e(s0)。泛克里格估計(jì)的公式如下[11]：

式（3）中，βk表示趨勢(shì)模型的估計(jì)參數(shù)，本研究中選擇普通最小二乘法進(jìn)行估計(jì)。wi表示不同采樣點(diǎn)位置殘差的權(quán)重，此權(quán)重由殘差半變異函數(shù)確定。e表示回歸模型的殘差。本文中所有回歸估計(jì)和泛克里格插值均在R統(tǒng)計(jì)軟件gstat包中進(jìn)行[12]。

地統(tǒng)計(jì)建模中的第一步是對(duì)樣本計(jì)算試驗(yàn)半變異函數(shù)值。采用以下公式4：

式（4）中：r(h)為回歸殘差的變異函數(shù)；N(h)是在某一空間方位上相距距離為h的固定樣地對(duì)數(shù)；z(xi)和z(xi+h)分別是固定樣地xi和與xi相隔距離h處的回歸殘差值。采用不同的模型，如指數(shù)模型，球狀模型等對(duì)試驗(yàn)半變異函數(shù)值進(jìn)行擬合，得到理論半變異函數(shù)模型。一般半變異函數(shù)r(h)是一個(gè)單調(diào)遞增函數(shù)，隨距離h增加而增加，當(dāng)超出一定距離（變程a）時(shí)，r(h)不再增加，此時(shí)的值稱為“基臺(tái)值（C）”。而理論變異函數(shù)在原點(diǎn)的取值稱為“塊金值（C0）”，塊金值表示存在小于抽樣尺度的變異或者存在測(cè)量誤差?；_(tái)值和塊金值之差稱為“偏基臺(tái)值（C1）”[13]。

本研究中對(duì)樹(shù)種多樣性趨勢(shì)模型采用的解釋變量為每個(gè)樣地的坐標(biāo)值和森林類型，森林類型變量為一個(gè)分類變量，取0表示純林，取1表示混交林。趨勢(shì)模型的具體表達(dá)式如下：

m(s)=β0+β1x2(s)+ β2y2(s)+ β3xy(s)+β4x(s)+β5y(s)+ β6f(s) （5）

式（5）中：β表示各個(gè)解釋變量的系數(shù)，f(s)的值取決于s點(diǎn)的森林類型。

2.4 模型評(píng)價(jià)和檢驗(yàn)

采用交叉檢驗(yàn)（cross-validation）來(lái)評(píng)價(jià)變異函數(shù)模型擬合效果。交叉檢驗(yàn)分析重復(fù)從已知數(shù)據(jù)集中刪除一個(gè)采樣點(diǎn)，用剩余的采樣點(diǎn)估測(cè)刪除點(diǎn)的數(shù)值[14-16]。對(duì)于泛克里格方法預(yù)測(cè)的精度， 采用隨機(jī)預(yù)留的140個(gè)檢驗(yàn)樣本進(jìn)行檢驗(yàn)。本研究中通過(guò)計(jì)算平均誤差（ME）， 平均絕對(duì)誤差（MAE），均方根誤差（RMSE）和決定系數(shù)（R2）來(lái)對(duì)模型進(jìn)行檢驗(yàn)和預(yù)測(cè)精度進(jìn)行評(píng)估。

各式中Xf和X0分別表示樹(shù)種多樣性的預(yù)測(cè)值和實(shí)測(cè)值，表示實(shí)測(cè)值的平均值。

3 結(jié)果與分析

3.1 回歸模型和殘差半變異函數(shù)分析

表2 回歸模型各變量參數(shù)的估計(jì)值Table 2 Estimation of parameters of regression models

表2中的結(jié)果為通過(guò)普通最小二乘法估計(jì)回歸模型式（4）中各個(gè)變量的系數(shù)。β0截距包含了變量f取0（純林）時(shí)的信息。表明樹(shù)種多樣性與純林存在著負(fù)相關(guān)關(guān)系，β6為變量f取1（混交林）的系數(shù)，表明樹(shù)種多樣性與混交林存在正相關(guān)關(guān)系。

通過(guò)趨勢(shì)回歸模型對(duì)每個(gè)樣地的樹(shù)種多樣性成分進(jìn)行回歸估計(jì)，從而得到每個(gè)樣地的樹(shù)種多樣性殘差，在此基礎(chǔ)上對(duì)殘差進(jìn)行半變異函數(shù)結(jié)構(gòu)分析。本研究采用了不同類型了半變異函數(shù)理論模型對(duì)殘差的試驗(yàn)半變異函數(shù)進(jìn)行擬合，包括指數(shù)模型，高斯模型，球狀模型，圓形模型。其中指數(shù)模型擬合效果最好，塊金值（C0）為0.070 31，偏基臺(tái)值（C1）為 0.022 87， 變程（a）達(dá) 3 000 m。這表明在變程范圍內(nèi)，當(dāng)前固定樣地抽樣尺度下各個(gè)樣地之間存在著一定的空間自相關(guān)性。表3為理論半變異模型的參數(shù)及交叉檢驗(yàn)統(tǒng)計(jì)量。圖3表示通過(guò)試驗(yàn)半變異函數(shù)值擬合的理論半變異函數(shù)模型。 從表3和圖2中均可以得出趨勢(shì)回歸模型殘差存在著空間自相關(guān)性。并且交叉檢驗(yàn)表明理論半變異函數(shù)擬合效果比較理想。

表3 理論半變異函數(shù)參數(shù)和交叉檢驗(yàn)統(tǒng)計(jì)量Table 3 Parameters of theoretical semivariance model and statistical quantity of cross-validation

圖2 趨勢(shì)回歸模型殘差的半變異函數(shù)Fig.2 Semivariances of residuals of regression model

3.2 檢驗(yàn)

泛克里格空間預(yù)測(cè)的精度，采用隨機(jī)預(yù)留的140個(gè)樣地的實(shí)測(cè)值與預(yù)測(cè)進(jìn)行統(tǒng)計(jì)比較，結(jié)果如圖 3， y=0.408 78+0.699 72x， R2=0.697 7，表明泛克里格方法預(yù)測(cè)精度較高，是一種可行的方法。

圖3 檢驗(yàn)樣本的實(shí)測(cè)值與預(yù)測(cè)值統(tǒng)計(jì)比較Fig. 3 Statistical comparison between measured and predicted diversity in validation dataset

3.3 空間分布估計(jì)及分析

通過(guò)泛克里格方法空間插值得到整個(gè)全局的樹(shù)種多樣性空間分布格局圖（圖4a）和標(biāo)準(zhǔn)誤差預(yù)測(cè)圖（圖4b）。從圖4a可看出，樹(shù)種多樣性空間分布具有明顯的異質(zhì)性，預(yù)測(cè)的空間分布圖中樹(shù)種多樣性最高值為1.82，最低值為0.221。在北部，西北部和南部邊緣地區(qū)，樹(shù)種多樣性較低。這些的地區(qū)主要以蒙古櫟純林為主，這是導(dǎo)致這些地區(qū)物種多樣性普遍偏低的主要原因。而且主要道路和耕地兩邊的樹(shù)種多樣性也普遍偏低，這些地區(qū)可能受到人為干擾較大。而東部，西部，中部和南部的樹(shù)種多樣比較高，主要因?yàn)檫@些地區(qū)的森林類型以闊葉混交林，針闊混交林和針葉混交林為主。從整體上看，汪清林業(yè)局地區(qū)森林的樹(shù)種多樣性比較高，這與林業(yè)局主要以天然林為主，占有林地面積的93.7%，天然林保護(hù)類型主要以限伐區(qū)和禁伐區(qū)為主， 因此受到人為干擾比較小。從圖4b可以看出，在林場(chǎng)邊緣地帶預(yù)測(cè)標(biāo)準(zhǔn)誤差較大，這與樣地?cái)?shù)量較少有關(guān)，影響插值精度。同時(shí)標(biāo)準(zhǔn)誤差預(yù)測(cè)圖，也給出來(lái)對(duì)于每個(gè)點(diǎn)的預(yù)測(cè)的可信度。標(biāo)準(zhǔn)誤差越高的地方表明需要加密樣地，來(lái)提高預(yù)測(cè)精度。

圖4 樹(shù)種生物多樣性空間分布（a）和標(biāo)準(zhǔn)預(yù)測(cè)誤差分布（b）Fig. 4 Prediction spatial distribution (a) and standard error distribution (b)

4 結(jié)論與討論

樹(shù)種多樣性的空間分布估計(jì)可以為森林經(jīng)營(yíng)者了解森林經(jīng)營(yíng)活動(dòng)對(duì)森林生態(tài)系統(tǒng)生物多樣性的影響和制定合理的可持續(xù)經(jīng)營(yíng)的方案提供依據(jù)。本文基于汪清林業(yè)局二類調(diào)查數(shù)據(jù)和局級(jí)固定樣地?cái)?shù)據(jù)，采用泛克里格方法對(duì)整個(gè)汪清林業(yè)局13個(gè)林場(chǎng)的樹(shù)種多樣性指數(shù)進(jìn)行估計(jì)，并得到樹(shù)種多樣性的空間格局分布圖。研究結(jié)果表明泛克里格方法空間插值方法是一種可行的方法，它可以解決尺度推移的問(wèn)題，根據(jù)二類調(diào)查固定樣地的數(shù)據(jù)，就可以得到整個(gè)林業(yè)局甚至是更大區(qū)域水平上的樹(shù)種多樣性估計(jì)，現(xiàn)實(shí)了由“點(diǎn)”上的樹(shù)種多樣性信息擴(kuò)展到“面”上樹(shù)種多樣性信息。

從空間估計(jì)的結(jié)果來(lái)看，汪清林業(yè)局地區(qū)森林的樹(shù)種多樣性具有明顯的空間格局和空間異質(zhì)性，這與當(dāng)?shù)氐膶?shí)際情況相符合，蒙古櫟林中樹(shù)種多樣性較小，而闊葉混交林，針葉混交林和針闊混交林中樹(shù)種多樣性較高。

采用泛克里格方法對(duì)樹(shù)種多樣性進(jìn)行估計(jì)，首先對(duì)趨勢(shì)進(jìn)行回歸分析，然后對(duì)殘差進(jìn)行半變異函數(shù)結(jié)構(gòu)分析， 進(jìn)而對(duì)樣地實(shí)驗(yàn)半變異函數(shù)值擬合理論半變異函數(shù)模型，最后進(jìn)行局部空間估計(jì)。這種方法不僅考慮了解釋變量（空間位置和森林類型）對(duì)樹(shù)種多樣性的影響，同時(shí)也考慮了樹(shù)種多樣性的空間自相性。今后還可以結(jié)合不同的環(huán)境變量（如海拔，氣候因子等）和遙感信息進(jìn)行研究，可以有效地提高預(yù)測(cè)精度。

[1] 雷相東, 唐守正. 林分結(jié)構(gòu)多樣性指標(biāo)研究綜述[J]. 林業(yè)科學(xué)，2002,38(3)∶140-146.

[2] 馮益明. 空間統(tǒng)計(jì)學(xué)理論及其在林業(yè)中的應(yīng)用[M]. 北京∶ 中國(guó)林業(yè)出版社, 2008，1-123.

[3] 洪 偉, 吳承幀. 杉木種源胸徑生長(zhǎng)地理變異規(guī)律的研究[J]. 植物生態(tài)學(xué)報(bào) , 1998，22(2)∶ 186-192.

[4]Malhi Y, Wood D, Baker T, et al. The regional variation of aboveground live biomass in old-growth Amazonian forests [J]. Global Change Biology, 2006,12∶1107-1138.

[5]Meng Q, Cieszewski C, Madden M. Large area forest inventory using Landsat ETM+∶ A geostatistical approach[J]. ΙSPRS Journal of Photogrammetry and Remote Sensing,2009, 64∶27-36.

[6]Akhavan R, Kia-Daliri H. Spatial variability and estimation of tree attributes in a plantation forest in the Caspian region of Ιran using geostatistical analysis[J]. Caspian Journal of Environmental Sciences,2010, 8(2)∶ 163-172.

[7]李明陽(yáng) , 劉 敏 , 劉米蘭 . 基于 GΙS 的森林調(diào)查因子地統(tǒng)計(jì)學(xué)分析 [J]. 南京林業(yè)大學(xué)學(xué)報(bào)：自然科學(xué)版 , 2010, (6)∶ 66-70.

[8]Sales M H, Souza Jr C M, Kyriakidis P C. Ιmproving spatial distribution estimation of forest biomass with geostatistics∶ A case study for Rond?nia, Brazil[J]. Ecological Modelling, 2007,205(1/2)∶ 221-230.

[9]Du H Q, Zhou G M, Fan W Y, et al. Spatial heterogeneity and carbon contribution of aboveground biomass of moso bamboo by using geostatistical theory[J]. Plant Ecol, 2010,207∶131-139.

[10] 劉曉梅 , 布仁倉(cāng) , 鄧華衛(wèi) , 等 . 基于地統(tǒng)計(jì)學(xué)豐林自然保護(hù)區(qū)森林生物量估測(cè)及空間格局分析[J]. 生態(tài)學(xué)報(bào), 2011，31(16)∶ 4783-4790.

[11] Cressie N A C. Statistic for spatial data[M]. John Wiley&Sons, New York,1993. 151-172.

[12] Pebesma E J. Multivariable Geostatistics in S∶ the Gstat Package[J]. Computer & Geoderma, 2004,67∶ 215-226.

[13] 王政權(quán). 地統(tǒng)計(jì)學(xué)及在生態(tài)學(xué)中的應(yīng)用[M]. 北京∶ 科學(xué)出版社 , 1999.

[14] 湯國(guó)安 , 楊 昕 . ArcGΙS 地理信息系統(tǒng)空間分析實(shí)驗(yàn)教程[M]. 北京：科學(xué)出版社 , 2006.

[15] 羅 梅 , 鄭小賢 . 八達(dá)嶺遼東櫟 -油松混交林空間結(jié)構(gòu)及其多樣性 [J]. 中南林業(yè)科技大學(xué)學(xué)報(bào) , 2012, 32(9)∶ 55-58.

[16] 李際平 , 趙春燕 , 袁曉紅 , 等 . 杉闊間物種多樣性影響因子灰色關(guān)聯(lián)分析 [J]. 中南林業(yè)科技大學(xué)學(xué)報(bào) , 2011, 31(10)∶ 10-14.

Estimation of spatial distribution of tree species diversity based on Universal Krige Model

HE Peng, ZHANG Hui-ru, LEΙ Xiang-dong, LΙ Xiao-dong
(Ιnstitute of Forest Resource Ιnformation Techniques, Chinese Academy of Forestry, Beijing 100091, China)

Biodiversity conservation is one of major goals of forest sustainable management. Estimation of spatial distribution is extremely signif i cant for forest manager to get better sense of the spatial pattern and dynamical changes of biodiversity. By taking the tree species diversity as an example and using Universal Krige (UK) , the spatial distribution of tree species diversity based on bureaulevel permanent plot data of forest inventory in Wangqing forest enterprise, Jilin Province were estimated. The results demonstrate that the UK for the studied area was a feasible method with a higher prediction accuracy (R2=0.6977); Spatially, tree species diversity varied obviously in this region from 1.82 to 0.221; Ιn addition, the spatial distribution map shows that there was an obvious trend in the studied area∶ lower tree species diversity in north and higher in east, west and south. Meanwhile, the results provide a method and reference for the region-scale tree species diversity estimation based on the permanent plot and geo-statistics.

tree species diversity; spatial distribution pattern; geo-statistics; universal kriging

S718

A

1673-923X(2013)12-0067-05

2013-03-25

林業(yè)公益性行業(yè)專項(xiàng)“我國(guó)典型森林類型健康經(jīng)營(yíng)關(guān)鍵技術(shù)研究”（201004002）

賀 鵬（1988-），男， 碩士，主要研究方向：森林資源經(jīng)營(yíng)管理與決策

張會(huì)儒，博士，研究員，主要研究方向：森林可持續(xù)經(jīng)營(yíng)研究； E-mail：huiru@caf.ac.cn

[本文編校：吳 彬 ]