北京山區(qū)油松生長模型研究

北京山區(qū)油松生長模型研究

越來越多的森林工人關注林分生長和收獲模型的兼容性。Cluttter首先提出林分的生長量與收獲量應該是一致的。Sullivan等和Pienaar等利用系統(tǒng)模型解決林分生長與收獲模型的兼容性問題。在他們所建的模型中,利用限制性條件保證了林分生長預測從A1到A3與從A1到A2、然后再從A2到A3預測的結果一樣,即階段無偏性(step-invariance)。然而正是由于這些限制性條件,降低了系統(tǒng)模型的預測精度。Ochi等提出利用林分的年生長預測來克服這些限制性條件,雖然直觀上來看不是對林分生長與收獲模型兼容性的研究,但是保證了林分生長預測的階段無偏性。據(jù)查閱相關文獻,國內(nèi)對階段無偏性這方面的研究很少。而且,目前森林資源監(jiān)測的發(fā)展趨勢正向監(jiān)測周期的年度化發(fā)展,連續(xù)定期清查的數(shù)據(jù)很難準確的提供森林每年的生長量和枯損量以及森林資源消長動態(tài)。為此,筆者利用北京山區(qū)油松Pinustabulaeformis的定期調查數(shù)據(jù),建立油松林分的年生長預測模型,進而達到以下兩個目的:一、為我國林分生長預測階段無偏性的研究提供一種方法;二、提供林分每年的生長變化情況,為資源監(jiān)測周期年度化的發(fā)展提供重要的理論依據(jù)。1樣地的復測計算本研究采用的數(shù)據(jù)來源于北京市林業(yè)調查設計院,其中油松的固定樣地有63個,每個樣地面積為0.067hm2。樣地主要調查因子有:林木胸徑、方位角、林分年齡、林分優(yōu)勢平均高、郁閉度、水平距、坡向、坡位、坡度、海拔高度、土層厚度等因子。樣地每隔5年復測1次,本研究利用的數(shù)據(jù)是1986、1991、1996、2001年復測的一類清查數(shù)據(jù),利用這些復測數(shù)據(jù)可以組成間隔調查期5年的樣地有156塊,間隔10年的樣地有93塊,間隔15年的樣地有39塊的數(shù)據(jù)集。因此,建模樣地有110塊(間隔5年),檢驗樣地共有178塊(間隔5年的樣地有46塊,間隔10年的樣地有93塊,間隔15年的樣地有39塊)。根據(jù)Ochi等和Yue等在研究中的檢驗方法,本研究對3個間隔期進行分開檢驗。油松林分樣地分布情況及變量因子統(tǒng)計分別列于表1、表2。2林分健全at/at全林分生長方程主要包括林分優(yōu)勢高、林分株樹密度、林分平均胸徑、林分斷面積、林分蓄積等方程。本研究利用林分平方平均胸徑代替林分算術平均胸徑,因為在林分統(tǒng)計中常用林分平方平均胸徑,而不是算術平均胸徑。年生長模型在單木生長模型中很常見,研究者們通過利用定期調查數(shù)據(jù),來估計年生長模型的參數(shù)。Cao在對單木生長模型進行年生長預測時,提出了可變生長率法,該方法考慮了單木每年生長的變化,而不是按每年固定生長率進行的。張雄清等利用可變生長率法對北京地區(qū)油松一類清查數(shù)據(jù)作了研究。那么,林分每年的生長量也是有變化的,不應該固定不變。因此,利用可變生長率法,對林分優(yōu)勢高方程、林分株樹密度方程、林分平均胸徑方程、林分斷面積方程和林分蓄積方程的推導如下。(t+1)年時:Hi+1=exp[(At/At+1)ln(Ht)+(1-At/At+1)(α1+α2Ht)];(1)Hi+1=exp{[(At/At+1)ln(Nt)+(1-At/At+1)(β1+β2At+β3/ln(Nt)]};(2)Dt+1=exp{(At/At+1)ln(Dt)+(1-At/At+1)[χ1+χ2At/ln(Nt)+χ3Dt]};(3)Bt+1=Bt+exp{(At/At+1)ln(Bt)+(1-At/At+1)[δ1+δ2Ht+δ3At/ln(Nt)+δ4ln(Bt)]};(4)Mt+1=Mt+exp{(At/At+1)ln(Mt)+(1-At/At+1)[γ1+γ2Ht+γ3/ln(Nt)+γ4ln(Bt)]}。(5)(t+q)年時:Ht+q=exp[(At+q-1/At+q)ln(Ht+q-1)+(1-At+q-1/At+q)(α1+α2Ht+q-1)]};(6)Nt+q=exp{(At+q-1/At+q)ln(Nt+q-1)+(1-At+q-1/At+q)(β1+β2At+q-1+β3/ln(Nt+q-1)]};(7)Dt+q=exp{(At+q-1/At+q)ln(Dt+q-1)+(1-At+q-1/At+q)[χ1+χ2At+q-1/ln(Nt+q-1)+χ3Di+q-1]};(8)Bt+q=Bt+q-1+exp{(At+q-1/At+q)ln(Bt+q-1)+(1-At+q-1/At+q)[δ1+δ2Ht+q-1+δ3At+q-1/ln(Nt+q-1)+δ4ln(Bt+q-1]};(9)Mt+q=Mt+q-1+exp{(At+q-1/At+q)ln(Mt+q-1)+(1-At+q-1/At+q)[γ1+γ2Ht+q-1+γ3/ln(Nt+q-1)+γ4ln(Bt+q-1]}。(10)式(1)～(10)中:q為生長期(調查間隔期,5a);At為t年時林分的平均年齡(a);Ht為t年時林分的優(yōu)勢木平均高(m);Nt為t年時林分的公頃株樹(株/hm2);Dt為t年時林分平方平均胸徑(cm);Bt為t年時林分斷面積(m2/hm2);Mt為t年時林分蓄積量(m3/hm2);α1、α2,β1、β3,χ1～χ3,δ1～δ5,γ1～γ4為待估參數(shù)。3參數(shù)估計與模型評價3.1林分平均胸徑和密度相關性全林整體模型是一個大的模型體系,該系統(tǒng)中的各模型從不同的側面對生長過程進行了描述,這些模型之間是相互關聯(lián)的。如:林分蓄積與林分平均高生長、林分平均胸徑和林分密度息息相關(cov(εi,εj)≠0)。似乎不相關聯(lián)立方程組作為聯(lián)立方程組的一個特例,從系統(tǒng)方程上看,這些方程之間不同的因變量依賴的自變量不同或者部分相同,利用最小二乘法對這些模型進行獨立估計時,會產(chǎn)生系統(tǒng)偏差,而采用似乎不相關聯(lián)立估計則會解決系統(tǒng)偏差,保證參數(shù)估計的有效性。因此,本研究利用似乎不相關聯(lián)立估計全林分模型參數(shù),減少系統(tǒng)誤差所帶來的估計偏差,這樣就能夠提高參數(shù)估計的有效性和一致性。3.2立地因子及預測方法林分優(yōu)勢高模型、林分株樹密度模型、林分平均胸徑生長模型、林分斷面積生長模型、林分蓄積生長模型等可以通過統(tǒng)計量平均偏差(EMD)、平均絕對偏差(EMAD)、均方根誤差(ERMSE)和決定系數(shù)(R2)進行評價。它們的數(shù)學表達式分別為:EMD=∑(yi-y′i)/n。(11)EMAD=∑|yi-y′i|/n。(12)ERΜSE=√∑(yi-y′i)2n-p。(13)R2=1-∑(yi-y′i)2/∑(yi-ˉyi)2。(14)式(11)～(14)中:yi為實際值(林分優(yōu)勢高、林分株數(shù)密度、林分平均胸徑、林分斷面積和林分蓄積量等);y′i、ˉyi分別為它們的預測值和平均值;n為樣地個數(shù);p為參數(shù)個數(shù)。評價一個模型或者一種預測方法的優(yōu)劣,可以利用上述這幾個統(tǒng)計量來完成。平均偏差、平均絕對偏差和均方根誤差小,且決定系數(shù)大,則該模型或該預測方法為優(yōu)。本研究參數(shù)估計和模型評價均利用SAS軟件來完成。4林分樹密度模型的檢驗表3列出了用林分優(yōu)勢高模型、林分株樹密度模型、林分平均胸徑生長模型、林分斷面積生長模型和林分蓄積生長模型進行獨立估計時,生長模型殘差的相關系數(shù)矩陣。根據(jù)表中數(shù)據(jù)可知,其中4個生長模型的殘差在置信水平為0.01下相關性顯著,這也表明了林分平均胸徑、林分斷面積、林分蓄積和林分株樹密度這4個林分變量因子存在著一定的相關性,即:cov(εi,εj)≠0。因此,根據(jù)Borders的方法對這4個方程的參數(shù)進行似乎不相關聯(lián)立估計,這樣就能夠提高參數(shù)估計的有效性和一致性。各林分生長模型的參數(shù)估計及標準誤見表4。由表4可知,各參數(shù)估計值比較穩(wěn)定,估計有效。表5列出了建模樣地及檢驗樣地的模型統(tǒng)計量。從表5可知:在建模數(shù)據(jù)中,林分優(yōu)勢高模型的EMD為0.0977,EMAD為0.6099,ERMSE為0.8983,R2為0.8690;林分株樹密度模型的EMD為11.1066,EMAD為113.3270,ERMSE為185.0350,R2為0.8578;林分平均胸徑生長模型的EMD為0.0812,EMAD為0.4920,ERMSE為0.6577,R2為0.9355;林分斷面積生長模型的EMD為0.2836,EMAD為1.15485,ERMSE為1.6296,R2為0.9321;林分蓄積生長模型的EMD為1.2365,EMAD為4.5492,ERMSE為6.5330,R2為0.9397。經(jīng)過Kolmogorov-Smirnov檢驗,林分優(yōu)勢高模型的殘差、林分株樹密度模型的殘差、林分平均胸徑生長模型的殘差、林分斷面積生長模型的殘差和林分蓄積生長模型的殘差都服從正態(tài)分布,且分別對林分優(yōu)勢高預測值、林分株樹密度預測值、林分平均胸徑預測值、林分斷面積預測值和林分蓄積預測值都沒有明顯的估計偏差(見圖1)。表6列出了利用3類檢驗數(shù)據(jù)所得的各林分生長模型的評價統(tǒng)計量。q=5時,林分株樹密度模型的EMD為36.0396,EMAD為39.2258,ERMSE為60.3986,R2為0.9810;q=10時,林分株樹密度模型的EMD為52.1159,EMAD為197.1320,ERMSE為284.7610,R2為0.6492;q=15時,林分株樹密度模型的EMD為96.0569,EMAD為238.4860,ERMSE為330.0300,R2為0.4891。比較這3個間隔期林分株樹密度模型EMD、EMAD、ERMSE和R2值的大小,可知,EMD(36.0396,q=5)<EMD(52.1159,q=10)<EMD(96.0569,q=15),EMAD(39.2258,q=5)<EMAD(197.1320,q=10)<EMAD(238.4820,q=15),ERMSE(60.3986,q=5)<ERMSE(284.7610,q=10)<ERMSE(330.0300,q=15),R2(0.9810,q=5)>R2(0.6492,q=10)>R2(0.4891,q=15)。很明顯,隨著預測時間變長,預測偏差越大。而且,在本研究中,主要是針對林分的年生長預測,建立林分年生長預測模型,因此,預測間隔的時間越長,預測的偏差會更大。同樣,林分平均胸徑生長模型、林分斷面積生長模型和林分蓄積量生長模型也一樣。5標準ce—結論林分年