基于啞變量的濕地松樹皮厚度模型研究

吳宏煒，黃光燦，張偉志，嚴銘海，莊崇洋，江希鈿

(福建農(nóng)林大學(xué)林學(xué)院，福建 福州 350002)

濕地松(PinuselliottiiEngelm)原分布于古巴、美國東南部等區(qū)域[1]，因其適應(yīng)性強，干型通直，木材質(zhì)量好[2-3]，抗侵染能力優(yōu)于馬尾松[4-5]，現(xiàn)已成為我國南方低丘陵地區(qū)的重要造林樹種之一[6-8]。濕地松的樹皮用途很廣，有研究表明其樹皮提取物具有一定的藥用價值，同時也是制作樹脂膠的原料[9]。樹皮厚度可用于計算出材率、木材材積和樹皮材積[10]。國內(nèi)外許多學(xué)者對樹皮厚度進行了研究，涉及的樹種有西南樺、銀杉、云杉、白樺和厚樸等[11-15]，有研究發(fā)現(xiàn)不同樹種的樹皮生長規(guī)律存在差異[13]，通常針葉樹的樹皮厚度比闊葉樹的大[10]。對于樹皮厚度模型的研究，擬合的因子主要有胸徑、樹高和年齡[12]，但較少有考慮立地質(zhì)量對樹皮厚度的影響。如果考慮立地質(zhì)量對樹皮厚度的影響，則需對生長在不同立地質(zhì)量等級的林木分別建立樹皮厚度模型，在實際工作中這可能會耗費更多人力物力，且可能會出現(xiàn)模型不兼容的問題，因此采用啞變量模型可為兼顧模型擬合效果和成本提供途徑[16]。

已有研究將立地質(zhì)量等級、競爭狀態(tài)、經(jīng)營管理措施等定性因子作為啞變量引入到模型中，使模型的擬合效果得到提升。朱光玉等[17]在研究櫟類天然林林分斷面積生長模型時引入了立地質(zhì)量等級為啞變量，華偉平等[18]在黃山松地位級指數(shù)模型的研制中引入立地質(zhì)量等級為啞變量，曹夢等[19]在單木胸徑和樹高生長模型研究中引入不同競爭狀態(tài)為啞變量，ZENG[20]在生物量的研究中引入森林起源為啞變量，王金池等[16]在云南松林分蓄積量生長模型的研究中引入間伐與未間伐指標為啞變量。因此本研究通過引入立地質(zhì)量等級作為啞變量，構(gòu)建濕地松樹皮厚度啞變量模型，以達到提升模型精度的目的。

本文通過參考前人研究所得的樹皮厚度模型，利用濕地松樣地和樣木數(shù)據(jù)，得到了濕地松帶皮直徑、去皮直徑和相關(guān)的樹皮厚度數(shù)據(jù)，利用含熵權(quán)值的TOPSIS綜合分析法選出最優(yōu)的樹皮厚度模型，在此基礎(chǔ)上，引入立地質(zhì)量作為啞變量建模，構(gòu)建了基于啞變量的濕地松樹皮厚度模型，為測算濕地松木材材積、出材率以及樹皮蓄積量提供依據(jù)。

1 研究區(qū)概況

寧德市地處福建省東北部(東經(jīng)118°32′～120°44′，北緯26°18′～27°04′)，土地面積13 400 km2，以丘陵山地為主兼沿海小平原,屬于亞熱帶海洋性季風(fēng)氣候，溫暖濕潤，降水充足。年平均氣溫13.4～20.2 ℃，年平均降雨量1 250～2 350 mm,無霜期270.4 d，日照時間1 637.7 h。

2 研究方法

2.1 材料來源

在寧德市濕地松人工林中設(shè)置樣地共86塊，包括臨時樣地與固定樣地，樣地形狀為方形與矩形，面積為0.04～0.67 hm2，設(shè)置于不同立地質(zhì)量等級的林分中。結(jié)合伐區(qū)設(shè)計數(shù)據(jù)，從中選取涵蓋不同年齡、胸徑與立地質(zhì)量的樣木258株，所在樣地基本信息見表1。伐倒后用皮尺(精度0.01 m)測量樹高與冠長，打枝后以2 m為區(qū)分段，用鋼尺(精度為0.1 cm)于樹干上0.3、1.0、1.3、3.0、5.0 m處依次測量帶皮直徑(outside bark diameter，dob)和去皮直徑(inside bark diameter，dib)。隨機抽取樣木193株(約占75%)用于構(gòu)建模型，剩余65株用于模型檢驗(約占25%)。樣木數(shù)據(jù)的描述性統(tǒng)計結(jié)果見表2，模型變量包括林分年齡(tree age，t)、樹高(tree height，H)、胸徑(diameter at breast height，DBH)、胸高處樹皮厚度(bark thickness at breast height, BBT)、任意樹高處樹皮率(該處兩倍樹皮厚度與帶皮直徑的比值，Bark ratio at any height，z)、為任意樹高處樹皮厚度(bark thickness at any height，BT)、相對樹皮厚度(任意樹高處樹皮厚度與胸徑處樹皮厚度的比值，relative bark thickness，RBT)、樹冠長(crown length，CL)和樹冠率(冠長與樹高的比值，crown rate，CR)等。

表1 樣木所在標準地信息Table 1 Information of sample plots

2.2 啞變量定義

啞變量(虛擬變量)是定性變量，通常取值為0或1，這種方法稱為定性因子(0,1)化展開，即將變量δ(x，i)表示為:δ(x，i)=0或1，當x是第i等級時為1，否者為0,因此稱變量δ(x，i)為啞變量。本研究引入立地質(zhì)量等級作為啞變量，根據(jù)福建省森林資源清查中的相關(guān)技術(shù)標準和實際林業(yè)生產(chǎn)情況將立地質(zhì)量分為Ⅰ肥沃、Ⅱ較肥沃、Ⅲ中等肥沃、Ⅳ貧瘠4個等級。將立地質(zhì)量等級用定性代碼0或1表示，第i種立地質(zhì)量等級表示為Si，定性數(shù)據(jù)Si轉(zhuǎn)化為(0，1)形式:Si為0或1，當x時第i等級時為1，否者為0。其中i=1、2、3、4；S1、S2、S3、S4分別是Ⅰ肥沃、Ⅱ較肥沃、Ⅲ中等肥沃、Ⅳ貧瘠的定性代碼。

表2 濕地松人工林樹皮厚度建模和檢驗數(shù)據(jù)Table 2 Fitting and valiating data of bark thickness models for P.elliottii plantations

2.3 樹皮厚度模型選擇

通過查閱相關(guān)資料，參考唐誠等[21]、嚴銘海等[22]建立的17個模型，并用于本研究，模型表達式見表3。其中模型3、模型4和模型5屬于對數(shù)回歸模型，需要計算校正因子進行修正，本文采用目前應(yīng)用較多的校正因子(correction factor，CF)[23]，其計算公式如下：

CF=exp (s2/2)

(1)

式中：CF為校正因子；s2為回歸方程樣本方差。

表3 選用模型表達式Table 3 The model forms

2.4 模型評價指標

采用最小二乘法對17個模型進行擬合，對模型參數(shù)進行顯著性分析，選取所有參數(shù)均達到顯著的模型，采用赤池信息準則(Akaike information criterion，AIC)、貝葉斯準則(Bayesian information criterion，BIC)、偏差(deviation，B)、絕對偏差(absolute deviation，AB)和決定系數(shù)(decision coefficient，R2)5個指標做模型質(zhì)量綜合評價，計算公式如下：

(2)

(3)

(4)

BIC=-2ln+lnnp

(5)

AIC=-2ln+2p

(6)

2.5 模型優(yōu)選方法

本研究采用含熵權(quán)值的TOPSIS法，在計算模型評價指標權(quán)重時運用熵值法取代一般的主觀權(quán)重法，避免了人為主觀性的影響[24-26]，其計算步驟如下。

(1)根據(jù)模型指標求解結(jié)果構(gòu)建m個評價對象、n個評價指標的判斷矩陣R，公式如下：

R=(xij)m×n(i=1,2,…，m；j=1,2,…，n)

(7)

(2)將判斷矩陣歸一化處理，得到無量綱化后的矩陣B，求解出的指標值范圍為0～1，其值越大代表指標越優(yōu)，最優(yōu)解的值為1，最差解值為0。

高優(yōu)指標(效益型指標)公式如下：

(8)

低優(yōu)指標(成本型指標)公式如下：

(9)

(3)計算熵權(quán)值Wj

(10)

(11)

(12)

式中：xij為模型評價指標值；xj，max為指標最大值；xj，min為指標最小值；rij為評價指標歸一化值；n為評價指標個數(shù)；fij為一致性指標；Hj為指標信息熵；Wj為熵權(quán)值。

3 結(jié)果與分析

3.1 樹皮厚度模型擬合結(jié)果

3.2 模型質(zhì)量評價

模型擬合效果與預(yù)估精度直接相關(guān)，選用的5個評價指標從不同方面體現(xiàn)模型的優(yōu)度，分別對4類模型的評價指標進行綜合分析，模型指標計算結(jié)果及含熵權(quán)值TOPSIS法結(jié)果見表5，各評價指標熵權(quán)值見表6。

表4 模型參數(shù)及顯著性檢驗Table 4 Model parameters and its significance test

表5 模型評價指標及TOPSIS法Table 5 Model statistical index and TOPSIS

表6 各評價指標熵權(quán)值Table 6 Entropy weight value of the model statistical index

在模型評價指標無量綱化后，指標取值范圍為0～1，其值越接近1代表指標優(yōu)度越高，其值越接近0代表指標優(yōu)度越低。模型3評價指標的B值在胸高處樹皮厚度模型中最小為0.000 002，該指標無量綱化后取值為1，代表該模型的評價指標B優(yōu)度在胸高處樹皮厚度模型中最高。模型6評價指標的AIC值在任意高度處樹皮厚度模型中最大為6 121.3，該指標無量綱化后取值為0，代表該模型的評價指標AIC優(yōu)度在任意高度處樹皮厚度模型中最低。

綜合本文所選的各評價指標無量綱化值及各評價指標的熵權(quán)值得到的結(jié)果表明，胸高處樹皮斷面積模型相對最優(yōu)解距離分別為0.194 3>0.034 3>0.000 9，表明模型質(zhì)量模型3>模型4>模型1，評價指標中AB熵權(quán)值最高為0.211 0，B熵權(quán)值最低為0.188 9。任意高度處樹皮厚度模型相對最優(yōu)解距離分別為0.115 3>0.086 5>0.000 4>0.000 0表明模型質(zhì)量模型9>模型10>模型8>模型6，評價指標中AB熵權(quán)值最高為0.206 3，R2熵權(quán)值最低為0.194 2。相對樹皮厚度模型相對最優(yōu)解距離分別0.200 4>0.042 2>0.000 4表明模型質(zhì)量模型12>模型13>模型14,評價指標中AB熵權(quán)值最高為0.215 1，B熵權(quán)值最低為0.190 6。去皮直徑模型相對最優(yōu)解距離分別0.001 5>0.000 4>0.000 0表明模型質(zhì)量模型16>模型15>模型17,評價指標中AB熵權(quán)值最高為0.205 0，R2熵權(quán)值最低為0.185 0。由此，選出模型3、模型9、模型12、模型16這4個相對最優(yōu)解距離最小的最優(yōu)模型。

3.3 啞變量模型及其擬合結(jié)果

在這4個最優(yōu)模型引入啞變量，經(jīng)過嘗試發(fā)現(xiàn)，在不同參數(shù)中加入啞變量時模型的決定系數(shù)有差異，模型3、模型9、模型12、模型16分別在參數(shù)α、α、β2、α中加入啞變量后所得模型不存在未顯著參數(shù)且擬合效果最佳，因此，啞變量模型的形式見表7。其中bi為啞變量Si的參數(shù)(i=1，2，3，4)，當立地質(zhì)量為Ⅰ肥沃時S1=1、S2=0、S3=0、S4=0，Ⅱ較肥沃時S1=0、S2=1、S3=0、S4=0，Ⅲ中等肥沃時S1=0、S2=0、S3=1、S4=0，Ⅳ貧瘠時S1=0、S2=0、S3=0、S4=1。以建模數(shù)據(jù)擬合以上4個啞變量模型，其參數(shù)求解值與各評價指標值見表8和表9。

表7 啞變量模型表達式Table 7 The form of dummy variable models

表8 啞變量模型參數(shù)及顯著性檢驗Table 8 Dummy variable model parameters and its significance test

由表9可知，4類啞變量模型的參數(shù)均顯著，4類啞變量模型較對應(yīng)的最優(yōu)基礎(chǔ)模型決定系數(shù)R2均得到提升，其中模型18的R2為0.592，模型19的R2為0.862，模型20的R2為0.868，模型21的R2為0.997。啞變量模型較最優(yōu)基礎(chǔ)模型AIC、BIC等低優(yōu)指標均有降低，其中模型18為 341.3、356.7；模型19為5 045.7、5 102.2；模型20為919.7、955.4；模型21為13 011.1、13 053.3，說明所建的啞變量樹皮厚度模型的擬合效果在篩選出的最優(yōu)模型的基礎(chǔ)上得到了提升。

3.4 模型擬合效果

根據(jù)用于檢驗的65株樣木實測數(shù)據(jù)，來進一步檢驗啞變量模型對生長于4種立地質(zhì)量等級的濕地松樹皮厚度預(yù)估效果。采取配對t檢驗的方法檢驗?zāi)Ｐ蛿M合效果，假設(shè)模型預(yù)估值與實測值之間存在顯著差異。若P<0.05 ，則假設(shè)成立，表示模型的擬合效果不好；若P>0.05，則假設(shè)不成立，表示模型擬合效果優(yōu)良(表10)。

結(jié)果表明，4類濕地松樹皮厚度啞變量模型的P值均大于0.05，其值分別為0.811、0.955、0.805、1.000。配對t檢驗結(jié)果表明本文所建立的4類濕地松樹皮厚度啞變量模型的預(yù)估值與實測值之間無顯著差異，表明所建模型可以做為計算生長于Ⅰ肥沃、Ⅱ較肥沃、Ⅲ中等肥沃、Ⅳ貧瘠4種不同立地質(zhì)量下濕地松樹皮厚度的參考。

表9 模型評價指標計算Table 9 Calculation of model statistical indexes

表10 模型配對t檢驗Table 10 Paired t test for the models

4 討論與結(jié)論

本研究對擬合出的模型進行顯著性檢驗，剔除含不顯著參數(shù)模型后運用5個評價指標進行綜合評價。依據(jù)含熵權(quán)值TOPSIS法評價結(jié)果篩選出最優(yōu)模型3、模型9、模型12、模型16，在這4個模型的基礎(chǔ)上引入啞變量，4個啞變量模型較基礎(chǔ)模型擬合效果均得到提升。啞變量模型經(jīng)配對t檢驗后，其P值分別為0.811、0.955、0.805、1.000，檢驗結(jié)果表明4類啞變量模型對檢驗樣本的預(yù)估值與其實測值不存在明顯差異，可以做為計算不同立地質(zhì)量下濕地松樹皮厚度的參考。對4類樹皮厚度啞變量模型的Si參數(shù)值進行分析，發(fā)現(xiàn)模型18的啞變量的Si參數(shù)值b4Ⅱ較肥沃>Ⅲ中等肥沃>Ⅳ貧瘠；模型19的啞變量的Si參數(shù)值b4Ⅱ較肥沃>Ⅲ中等肥沃>Ⅳ貧瘠；模型20的啞變量的Si參數(shù)值b2Ⅰ肥沃>Ⅳ貧瘠>Ⅱ較肥沃；模型21的啞變量的Si參數(shù)值b4Ⅱ較肥沃>Ⅲ中等肥沃>Ⅳ貧瘠。其中模型18、模型19、模型21均是在常數(shù)項加入啞變量的擬合效果最佳，依此可見立地質(zhì)量對濕地松胸高處樹皮厚度、任意樹高處樹皮厚度和去皮直徑可能有比較直接的影響而且是加法關(guān)系。模型20在常數(shù)項加入啞變量后出現(xiàn)了參數(shù)不顯著的情況，在指數(shù)項引入啞變量后模型的擬合效果最好。

基礎(chǔ)模型選優(yōu)采用含熵權(quán)值的TOPSIS法，選用的各項評價指標的熵權(quán)值由其本身的值決定，該方法具有較好的客觀性[24]，評價過程避免了人為主觀的確定各個指標權(quán)重而影響綜合評價的結(jié)果。本文中4種類型的模型各項評價指標熵權(quán)值均有不同，其中不同指標在同一類型模型中的重要性不一致且同一指標在不同類型模型中的重要性也存在差異，因此有必要采用改進的含熵權(quán)值的TOPSIS法對模型進行綜合評價。