史阿榮，李書恒，侯 麗，陳 蘭，白紅英

(1.西北大學城市與環(huán)境學院，陜西 西安 710127；2.西北大學陜西省地表系統(tǒng)與環(huán)境承載力重點實驗室，陜西 西安 710127)

樹木的生長輪，即樹輪資料具有定年準確、連續(xù)性強、分辨率高以及易于獲取復(fù)本等特點[1]，采用樹木年輪寬度值作為代用指標進行樹木生長的研究，有利于探究樹木生長量變化趨勢[2]。由于海拔、坡向和坡度等因子的不同，樹木的生長也隨之不同，隨著樹木年輪多樣化和深入化的研究，樹木年輪資料對研究樹木生長規(guī)律、分析地形因子對樹木生長的影響有著十分重要的作用。太白紅杉(Larixchinensis)是秦嶺高山林線森林群落的建群種[3]，研究表明太白紅杉的生長發(fā)育受環(huán)境的影響較顯著，如苗莉云對太白紅杉群落與環(huán)境因子進行逐步多元回歸和典范相關(guān)分析顯示，海拔和坡度決定太白紅杉群落物種多樣性的主要環(huán)境因子[4]。郭華用主成分分析研究太白山、光頭山樣地以及總體生境的環(huán)境因子對太白紅杉種群格局的影響，結(jié)果顯示太白紅杉種群格局特征主要影響因子有海拔高度、坡度和坡向[5]等。以上研究表明環(huán)境制約因子中的海拔、坡度和坡向?qū)淠旧L影響較大，因此對樹木生長規(guī)律的研究必須考慮外部的生長因素，即地形因子的制約。外部的生長因素影響生長變化，內(nèi)部的生長動態(tài)反映樹木的生長趨勢，關(guān)于太白紅杉生長趨勢的相關(guān)研究中，建立生長方程和生長模型是當前的主流研究方向，王孝安利用樹輪資料構(gòu)建4種競爭指數(shù)的太白紅杉單木生長的最優(yōu)模型[6]，段仁燕采用樹輪資料和多元回歸的方法構(gòu)建的太白紅杉的生長模型[7]，瞿亞明利用杉木樹干材積資料建立杉木連年生長量模型[8]；肖東耀根據(jù)海拔范圍和樹輪資料，利用5種常用的生長模型擬合太白紅杉平均胸徑生長過程[9]。以上研究考慮單木生長量、連年生長量結(jié)合回歸分析時，需要添加輔助信息并建立生長曲線，個別參數(shù)獲取相對困難。近年來，在樹木生長量研究中出現(xiàn)了新的研究方向，即根據(jù)樹種斷面積增長量(BAI，Basal Areas Increment)研究樹木生長變化。研究表明該模型可較準確的模擬BAI的生長趨勢，如王寶欽用ARMA(0,0)模型發(fā)現(xiàn)34個BAI序列中有30個櫟樹斷面積生長量可以得到最好的模擬，而另4個用ARMA(1,0)模型可得到最好的模擬，所有黃杉和挪威云杉BAI序列都能用ARMA(1,0)模型來模擬[10]。張荷觀采用ARIMA模型研究白樺、馬尾松、云杉、白楊和紙皮樺的生長量并進行模擬預(yù)測，結(jié)果顯示ARIMA模型預(yù)測效果更好[11]，目前利用ARIMA模型研究太白紅杉BAI生長規(guī)律尚未見報道。鑒于以上研究結(jié)果的優(yōu)良性，證明該模型可根據(jù)單個變量的取值，并對其自身的變動進行預(yù)測，無須添加任何的輔助信息。本文基于太白紅杉年輪寬度采用時間序列ARIMA模型模擬太白紅杉BAI序列生長過程，并探討B(tài)AI與立地因子的相關(guān)關(guān)系。希望能為保護和研究太白紅杉提供新的方法和思路，為進一步探知太白紅杉的生長規(guī)律奠定實踐基礎(chǔ)。

1 數(shù)據(jù)研究區(qū)

太白紅杉也叫太白落葉松，是松科落葉松屬紅杉組的植物[12]。形成于古冰川遺跡之上，是我國陜西唯一天然分布的落葉松，是太白山林線上限建群樹種[13]。太白紅杉集中分布在太白山(N33°49′～34°10′，E107°19′～107°58′)地區(qū),另外在長安光頭山、牛背梁自然保護區(qū)、長青自然保護區(qū)等山頭也有少量分布[14]。采集對象木樣本位于樹種適生區(qū)[3]3 000 m～3 300 m(圖1)。氣象站數(shù)據(jù)參考太白縣氣象站(N34°02′，E107°19′，海拔 1 543.6 m)。

共采集60棵樹芯，將采集好的樹芯樣本按照國際通用基本處理程序[15]，結(jié)合計算機程序COFECHA[16]。采樣點信息匯總(見表1)。

2 研究方法

2.1 計算方法

利用實測輪寬數(shù)據(jù)與胸徑序列得出杉木每年的胸徑面積增長量(BAI)數(shù)學公式[10]：

表1太白山采樣點信息

Tab. 1 Sampling information of the Taibai Mountains

圖1 采樣點及氣象站點的分布Fig.1 Distribution of sample sites and meteorological stations

(1)

2.2 模擬步驟

ARIMA模型的又稱為自回歸移動平均模型，由自回歸過程AR(p)、差分次數(shù)(d)、移動平均過程MA(q)組成。是將預(yù)測對象隨時間推移而形成的數(shù)據(jù)序列視為1個隨機序列，用一定的數(shù)學模型來近似描述這個序列，步驟[11]如下：(1)根據(jù)時間序列的散點圖、自相關(guān)函數(shù)和偏自相關(guān)函數(shù)圖以單位根檢驗其方差、趨勢及其季節(jié)性變化規(guī)律，對序列的平穩(wěn)性進行識別。(2)若數(shù)據(jù)隨著時間的變化存在一定的增長或下降趨勢，則需對數(shù)據(jù)進行差分，若差分后數(shù)據(jù)存在異方差，則需采取對因變量取log值等處理。(3)根據(jù)平穩(wěn)序列的自相關(guān)和偏自相關(guān)圖拖尾和截尾，判斷p、q的值并將其帶入模型進行參數(shù)檢驗，檢驗?zāi)Ｐ褪欠窬哂薪y(tǒng)計意義。(4)進行白噪聲檢驗，診斷殘差序列是否為白噪聲。(5)最后利用已通過檢驗的模型進行預(yù)測分析。

3 結(jié)果與討論

3.1 海拔、坡向和坡度因子與太白紅杉BAI的空間分布

依照南北坡向為劃分基礎(chǔ)，分別繪制兩坡向內(nèi)不同海拔與坡度范圍的太白紅杉分布圖，如圖2所示，在海拔 3 000 m～3 300 m與坡度10°～30°中，南坡樣本數(shù)量明顯高于北坡，從BAI生長幅度來看，南北坡海拔范圍相近，坡度差異明顯。

對比發(fā)現(xiàn)：南坡5°～25°范圍太白紅杉BAImax=53.91 cm2，BAImin=4.64 cm2，BAImean=23.64 cm2；北坡25°～30°范圍太白紅杉BAImax=48.03 cm2，BAImin=8.59 cm2，BAImean=21.04 cm2。由此可得，太白紅杉較適宜生長區(qū)應(yīng)位于海拔 3 080 m～3 200 m和坡度15°～25°的地區(qū)。

3.2 海拔、坡度和坡向因子與太白紅杉BAI的關(guān)系

為研究海拔、坡度與坡向3種地形因子與BAI的相關(guān)關(guān)系，采用雙變量相關(guān)性分析BAI、海拔、坡向與坡度數(shù)據(jù)(見表2)，顯著性檢驗采用雙側(cè)檢驗(T)。由表2可知，平均BAI與海拔、坡度負相關(guān)性較顯著，南坡BAI與海拔、坡度呈負相關(guān)(r=-0.640**，P<0.01),另外平均BAI與南北坡向均為正相關(guān)性關(guān)系，與BAI南坡向呈顯著正相關(guān)性(r=0.496*，P<0.05)，與北坡向正相關(guān)性較弱(r=0.371，P<0.05)。這表明海拔越高、坡度越大則BAI的生長速度就越慢，這與太白紅杉的生長習性相似。

3.3 ARIMA模型模擬分析太白紅杉BAI的生長趨勢

樹木的BAI隨著生長年份的增加而增加， 采用時間序列ARIMA模型分析太白紅杉BAI生長過程，為了消除數(shù)據(jù)的不平穩(wěn)和減弱BAI波動幅度，對BAI進行一階差分且BAI取log值，后對log(BAI)進行一階差分，數(shù)據(jù)趨于平穩(wěn)(見圖3)。

表2太白紅杉BAI與地形因子的相關(guān)分析

Tab. 2 Correlation statistics of BAI and slope, slope direction and elevation of Larix chinensis

注: * * P < 0. 01; * P < 0. 05

圖2 不同海拔、坡度和坡向范圍的太白紅杉分布圖Fig. 2 Distribution of Larix chinensis at different altitudes, slopes and slope aspects

得到平穩(wěn)的時間序列之后，對其建立ARIMA模型。由上可知，在本文中，在d=2時，數(shù)據(jù)達到平穩(wěn)。由圖4自相關(guān)與偏自相關(guān)圖可知，為確保因變量之間、誤差與因變量沒有相關(guān)性而得到平穩(wěn)的時間序列，自相關(guān)與偏自相關(guān)系數(shù)的取值縮小至0左右。因而建立ARIMA模型，選取BIC值最小為最終模型， 模型BIC(BIC=kln(n)-2ln(L)) 值越小越好，k為模型參數(shù)個數(shù)，n為樣本數(shù)量，L為似然函數(shù)。kln(n)在維數(shù)過大且訓練樣本數(shù)據(jù)相對較少的情況下，可以有效避免出現(xiàn)維度錯亂現(xiàn)象，同時也避免了選擇模型是偏向于維數(shù)大的模型。L為似然函數(shù)，依據(jù)使似然函數(shù)最大選擇模型，則BIC要最小)。

圖3 太白紅杉log(BAI)一階差分后與生長年份的序列圖Fig. 3 Sequence diagram of first order difference and growth year of Larix chinensis log (BAI)

由表3和表4得，ARIMA(1，1，2)模型BIC值最小且各參數(shù)的sig.值均小于0.05，表明各參數(shù)均通過顯著性檢驗。根據(jù)上述擬合參數(shù)，組合的模型公式表示為：

Y=0.011+1.562yt-1-0.562yt-2+εt-0.73εt-1-0.269εt-2

(2)

式中Y為 log(BAI)，且該模型的殘差為白噪聲。

表3 ARIMA模型統(tǒng)計表

表4 ARIMA模型參數(shù)統(tǒng)計量

Tab. 4 Parameter statistics of ARIMA model

根據(jù)ARIMA(1，1，2)模型繪制太白紅杉BAI生長模擬圖(見圖5)。從圖5中可以看出，置信區(qū)間長度較小，預(yù)測精度較高。在1840年左右出現(xiàn)歷史最低值，1997年左右出現(xiàn)歷史最高值。造成特殊值出現(xiàn)原因有兩方面：一方面，特殊年份的最高最低生長值出現(xiàn)，是由于1840年～1890年屬于氣候寒冷期，降水冬多夏少，而1997年我國大部分地區(qū)平均氣溫偏高，為北方近30 a降水最少年份[18]，夏季北熱南涼，Rolland[19]在研究中提到，夏季溫度升高對樹木徑向生長的抑制作用更為明顯，高海拔樹木在強烈的太陽輻射下，其代謝會受到影響，進而導致樹木自身缺水，使正常生長受到抑制?？梢?，除了中間競爭和地形因子影響以外，還受到當年氣候的制約而導致太白紅杉BAI出現(xiàn)歷史最低值和最高值。另一方面，太白紅杉隨著生長年限增加生長從快速增長到之后的增長速度減緩，這是受K選擇對策的影響，太白紅杉的幼齡樹皆具有較強的生存力，外加成林的保護關(guān)系使幼齡樹的BAI增長速率持續(xù)加快，但隨著年齡的增長種內(nèi)的保護關(guān)系演變成競爭關(guān)系導致BAI的生長速率開始減慢。從大的時間尺度來看，太白紅杉BAI生長高峰約10 a1次，BAI年生長超過3.0 cm2約有12 a以上，低于1.0 cm2約有9 a左右，大幅度生長動蕩期約兩個階段：第1個階段為中齡80 a～120 a，持續(xù)40 a；第2個階段為老齡177 a～200 a，持續(xù)23 a，兩階段間隔57 a，后55 a有小幅度降低階段，持續(xù)20 a，之后趨于平穩(wěn)。

圖5 ARIMA模型模擬太白紅杉280年BAI生長趨勢圖Fig. 5 ARIMA simulation of BAI growth trend of Larix chinensis in 280 years

利用時間序列ARIMA模型和相關(guān)性分析法研究了太白紅杉BAI與地形因子(海拔、坡度、坡向)的關(guān)系和太白紅杉BAI的生長過程，由于太白紅杉BAI與地形因子海拔、坡度負相關(guān)顯著，認為太白紅杉具有生長的范圍限制，海拔與坡度越大則BAI的生長速度就越慢，較適宜生長區(qū)應(yīng)位于海拔 3 080 m～3 200 m和坡度15°～25°的地區(qū)；ARIMA模型顯示研究區(qū)年齡最大的太白紅杉的生長趨勢約10 a左右經(jīng)歷1次生長高峰，幼齡太白紅杉BAI處快速增長的階段，中齡后太白紅杉BAI增速減緩，但整體生長趨勢呈現(xiàn)較平穩(wěn)上升的態(tài)勢。

參考文獻：

[1] 吳祥定.樹木年輪分析在環(huán)境變化研究中的應(yīng)用[J].第四紀研究.1990,(2):188～196.

[2] 吳祥定,邵雪梅.采用樹輪寬度資料分析氣候變化對樹木生長量影響的嘗試[J].地理學報,1996(s1): 92～101.

[3] 康永祥,劉婧輝,代拴發(fā)，等.太白山不同海拔太白紅杉年輪生長對氣候變化的響應(yīng)[J].西北農(nóng)林科技大學學報:自然科學版,2010,38(12):141～147.

[4] 苗莉云.太白紅杉群落交錯帶生態(tài)學特性的研究[D].陜西師范大學,2005.

[5] 郭華.秦嶺太白紅杉種群空間分布格局研究[D].陜西師范大學,2004.

[6] 王孝安，段仁燕，王明利．太白紅杉單木胸徑生長模型的研究[J]．武漢植物學研究2005,23(2) 157～162．

[7] 段仁燕，王孝安，黃敏毅.太白紅杉徑向生的預(yù)測模型[J].廣西植物.2009.3.29(2):212 ～216.

[8] 翟亞明,楊倩.長嶺杉木樹干材積連年生長量的模擬模型研究[J].江西農(nóng)業(yè)學報,2015(6):27～29.

[9] 肖東耀,廖超英,楊曉娟,等.太白山北坡太白紅杉(Larixchinensis)胸徑生長規(guī)律研究[J].西北林學院學報, 2012,27 (4): 194～198.

[10] 王寶欽.用年輪數(shù)據(jù)估計生長趨勢[J].林業(yè)資源管理,1992(1):71～75.

[11] 張荷觀.時間序列的ARIMA模型在預(yù)測樹木生長量中的應(yīng)用[J].林業(yè)科學,1986,22(1):94～100.

[12] 王戰(zhàn).中國落葉松林[M].中國林業(yè)出版社,1992.

[13] 張仰渠.陜西森林[M],北京:中國林業(yè)出版社,153～156.

[14] 孫毓,王麗麗,陳津,等.中國落葉松屬樹木年輪生長特性及其對氣候變化的響應(yīng)[J].中國科學:地球科學, 2010,40 (5): 645～653.

[15] Stokes M A, Smiley T L. An introduction to tree-ring dating[M].University of Arizona Press,1968.

[16] Holmes R L.Computer-Assisted Quality Control in Tree-Ring Dating and Measurement[J]. Tree-Ring Bull,1983,44(3):69～75.

[17] 祁建,馬克明,張育新.北京東靈山不同坡位遼東櫟(Quercusliaotungensis)葉屬性的比較[J].生態(tài)學報,2008,28(1):122～128.

[18] 竺可楨.中國近五千年來氣候變遷的初步研究[J].中國科學,1973(2):168～189.

[19] Rolland C.Tree-ring and climate relationships for Abies Alba in the internal Alps[J].1993,53.

[20] 吳祥定.樹木年輪與氣候變化[M].北京:氣象出版社,1990.