常曉敏,余新曉,*,賈國(guó)棟,孫立博,劉自強(qiáng),楊照國(guó),閆騰飛

1 北京林業(yè)大學(xué),水土保持與荒漠化防治教育部重點(diǎn)實(shí)驗(yàn)室, 北京 100083 2 豐寧縣大灘林場(chǎng), 石家莊 068357 3 張北縣林業(yè)局, 張家口 076450

近年來全球氣候變化背景下大范圍樹木死亡逐漸引起人們的廣泛關(guān)注[1-2]。地區(qū)水平上大面積的樹木死亡不僅影響碳平衡,還影響森林生態(tài)系統(tǒng)服務(wù)功能[3]。我國(guó)的三北地區(qū)氣候也正經(jīng)歷著暖干化變化[4]。該地區(qū)三北防護(hù)林經(jīng)過40多年的建設(shè),為治理風(fēng)沙、改善京津冀環(huán)境質(zhì)量發(fā)揮了重要作用[5]。但是近十年來,楊樹防護(hù)林出現(xiàn)了大面積的枯梢和死亡現(xiàn)象?？萆伊址置娣e接近80%,枯死和瀕死的約占總面積的1/3[6],防護(hù)功能明顯下降。水分是該區(qū)域植物生長(zhǎng)的主要限制因子[7],同時(shí)粗放經(jīng)營(yíng)也是導(dǎo)致防護(hù)林生長(zhǎng)發(fā)育不良甚至衰退的重要原因[4]。研究楊樹林的死亡原因和影響因素,對(duì)全球氣候變化背景下三北防護(hù)林的林分經(jīng)營(yíng)管理是十分必要的。

國(guó)內(nèi)外很多學(xué)者從造林方式[8]、林分結(jié)構(gòu)[9-11]、立地條件[12-13]、氣候環(huán)境條件[14-16]等多方面對(duì)樹木死亡的原因開展了大量研究。Qiu等[9]建立了林木相對(duì)斷面積、胸徑、密度、土壤理化性質(zhì)和氣候因子與樹木死亡關(guān)系模型,認(rèn)為林分密度過大,或者林木相對(duì)斷面積增長(zhǎng)量達(dá)不到一定值時(shí),樹木死亡率增加。陳清等[17]研究了林齡和林分?jǐn)嗝娣e對(duì)樹木死亡的影響。Fan等[18]研究了胸徑、密度、坡度、坡向、海拔等因子與樹木死亡的關(guān)系,發(fā)現(xiàn)胸徑是影響樹木死亡最重要的因子。De Toledo等[11]總結(jié)了影響樹木死亡的土壤和地質(zhì)因子。以上研究表明影響樹木死亡的原因有很多,包括胸徑、密度、林木相對(duì)斷面積、土壤理化性質(zhì)等。然而不同徑級(jí)的林木,導(dǎo)致其死亡的原因是否相同尚不清楚。

本文以張北壩上地區(qū)的小葉楊(PopulussimoniiCarr)為研究對(duì)象,以林木胸徑、林木相對(duì)斷面積、林齡、林分密度、不同土層土壤含水量和土壤容重等11個(gè)因子作為輸入變量,運(yùn)用隨機(jī)森林建立不同徑級(jí)樹木死亡模型,對(duì)影響樹木死亡的變量重要性進(jìn)行評(píng)估,揭示不同徑級(jí)楊樹死亡的主要原因,這對(duì)于今后壩上地區(qū)防護(hù)林的恢復(fù)重建和經(jīng)營(yíng)管護(hù)具有重要的意義。

1 材料與方法

1.1 研究區(qū)概況

試驗(yàn)地點(diǎn)位于河北省壩上地區(qū)的張北縣,屬溫帶大陸性氣候,多年平均氣溫3.2℃,極端最高氣溫30.8℃,極端最低氣溫-30.3℃,年降水量不足400mm。每年約有50—70d風(fēng)力達(dá)到六級(jí)以上,以偏西風(fēng)和西北風(fēng)為主,土壤風(fēng)蝕嚴(yán)重。地帶性土壤為栗鈣土和沙壤土。該地區(qū)楊樹防護(hù)林自20世紀(jì)60年代開始營(yíng)造,小葉楊為主要樹種。近年來,楊樹出現(xiàn)大面積的枯梢枯死現(xiàn)象。喬木樹種還有榆樹、落葉松等,主要灌木有沙棘、檸條等。

1.2 數(shù)據(jù)來源與處理

1.2.1數(shù)據(jù)來源

張北的防護(hù)林以楊樹純林為主。2016年7—8月和2017年7—8月,在對(duì)張北全縣18個(gè)鄉(xiāng)鎮(zhèn)進(jìn)行踏查的基礎(chǔ)上,目測(cè)立地類型、胸徑、樹高、林齡、密度等主要調(diào)查因子后,以此為輪廓設(shè)置了100塊20m×20m的標(biāo)準(zhǔn)地進(jìn)行調(diào)查,基本上涵蓋了張北防護(hù)林不同立地類型、林分密度、林齡、胸徑分布范圍等。當(dāng)調(diào)查對(duì)象為帶狀防護(hù)林時(shí),為減少邊緣效應(yīng),設(shè)置矩形標(biāo)準(zhǔn)地,并根據(jù)實(shí)際情況改變標(biāo)準(zhǔn)地的大小,使得標(biāo)準(zhǔn)地內(nèi)的林木特征和立地條件盡量一致。

記錄標(biāo)準(zhǔn)地的地理位置、海拔、地形等因子。由于張北的防護(hù)林都在平坦區(qū)域,因此沒有記錄坡度、坡向、坡位等信息。對(duì)標(biāo)準(zhǔn)地內(nèi)的楊樹進(jìn)行每木檢尺,記錄胸徑、樹高、林木枯梢情況、林木死亡情況等。楊樹防護(hù)林大多為同齡林,因此每塊標(biāo)準(zhǔn)地選取3株標(biāo)準(zhǔn)木,應(yīng)用生長(zhǎng)錐測(cè)定林木的年齡,取其平均值代表該林分的年齡。在每個(gè)標(biāo)準(zhǔn)地內(nèi)按照一定方向選取3個(gè)有代表性的位置挖掘土壤剖面,分別在每一剖面按照0—10、10—20、20—60cm采集土壤樣本,并把同一層次土壤樣本進(jìn)行混合。利用烘干法測(cè)定土壤含水量,環(huán)刀法測(cè)定土壤容重、孔隙度等土壤物理性質(zhì)。由于土壤含水量是變化值,因此對(duì)樣地進(jìn)行了多次調(diào)查并取其平均值。該地區(qū)降水頻率和強(qiáng)度都較小,為研究一般情況下土壤含水量對(duì)樹木死亡的影響,調(diào)查時(shí)間選擇在長(zhǎng)期未降水的時(shí)期,避免降水對(duì)土壤含水量造成的影響。

1.2.2數(shù)據(jù)處理

通過實(shí)地調(diào)查發(fā)現(xiàn),張北的楊樹林以成熟林和過熟林居多,即大徑級(jí)林木較多。為了使不同徑級(jí)的林木樣本量相近,去除掉部分大徑級(jí)的林木和部分生長(zhǎng)健康無枯死木的林分,只保留部分生長(zhǎng)健康的林分,最終選取26塊樣地,基本涵蓋了張北地區(qū)楊樹林的徑級(jí)范圍,并考慮到不同的立地類型和林分密度等因子。

研究顯示與樹木死亡有顯著相關(guān)性的因子有林木徑級(jí)大小、林齡、密度、土壤理化性質(zhì)、立地等[9,18]。張北的楊樹防護(hù)林分布于平坦區(qū)域,海拔的分布范圍較窄。研究表明張北防護(hù)林的退化死亡情況與土壤養(yǎng)分沒有明顯的相關(guān)關(guān)系[6]。因此選擇表征林齡、林木徑級(jí)大小、林木疏密程度、土壤物理性質(zhì)的因子。表征林木徑級(jí)大小的因子有胸徑和林木相對(duì)斷面積,表征林木疏密程度的因子有林分密度和林分?jǐn)嗝娣e。土壤物理性質(zhì)包括0—20cm土壤含水量、20—60cm土壤含水量、土壤平均含水量、0—20cm土壤容重、20—60cm土壤容重、土壤平均容重。林木相對(duì)斷面積即林木斷面積與林分平均斷面積的比。

應(yīng)用主成分分析(PCA)的方法對(duì)11個(gè)解釋變量降維,經(jīng)過降維后得到4個(gè)主成分,累計(jì)方差貢獻(xiàn)率為87%。PC1主要由土壤含水量和土壤容重決定,說明土壤物理性質(zhì)；PC2主要由林木相對(duì)斷面積決定；PC3主要由密度決定；PC4說明林齡。將林木分為5個(gè)徑級(jí)(<10cm、10—15cm、15—20cm、20—25cm、>25cm),利用邏輯斯蒂?gòu)V義線性混合效應(yīng)模型分別不同徑級(jí)建立樹木死亡與4個(gè)主成分因子的回歸模型。因變量為二元分類變量,枯死記為1,無枯死記為0。由于空間自相關(guān)性的影響,距離相近的林木其生長(zhǎng)與死亡可能互相影響,因此將林木所屬的樣地號(hào)作為模型的隨機(jī)效應(yīng)。固定效應(yīng)系數(shù)b表示模型中固定效應(yīng)對(duì)因變量影響大小的參數(shù),將固定效應(yīng)系數(shù)和模型顯著性相近的劃分為一類,以此作為劃分徑級(jí)大小的標(biāo)準(zhǔn)。

應(yīng)用隨機(jī)森林研究解釋變量對(duì)樹木死亡的相對(duì)重要性,建立不同徑級(jí)樹木死亡預(yù)測(cè)模型。

所有的統(tǒng)計(jì)分析應(yīng)用R軟件實(shí)現(xiàn)。

1.3 研究方法

1.3.1邏輯斯蒂?gòu)V義線性混合模型

設(shè)樹木死亡的概率為P(二項(xiàng)分類因變量Y=1),樹木無枯死的概率為1-P(二項(xiàng)分類因變量Y=0)。對(duì)P進(jìn)行Logit變換,

(1)

概率P與解釋變量的關(guān)系為

(2)

1.3.2隨機(jī)森林算法

隨機(jī)森林算法是由Breiman于2001提出的一種基于分類回歸樹的機(jī)器學(xué)習(xí)算法。如果因變量有n個(gè)觀測(cè)值,有k個(gè)解釋變量與之相關(guān),它利用Bootstrap重抽樣的方法從n個(gè)原數(shù)據(jù)中隨機(jī)抽取n個(gè)觀測(cè)值,選擇部分變量作為分類樹節(jié)點(diǎn),從而構(gòu)建幾百甚至幾千個(gè)分類樹組合成隨機(jī)森林,將重復(fù)程度最高的樹作為最終結(jié)果。利用隨機(jī)森林算法可以對(duì)解釋變量的相對(duì)重要性進(jìn)行評(píng)價(jià),對(duì)解釋變量進(jìn)行預(yù)測(cè)。平均準(zhǔn)確率降低度是衡量將一個(gè)變量的取值變?yōu)殡S機(jī)數(shù)后隨機(jī)森林預(yù)測(cè)準(zhǔn)確性降低程度的指標(biāo),該值越大,說明該變量的重要性越大。

2結(jié)果與分析2.1不同徑級(jí)樹木死亡與環(huán)境因子關(guān)系

根據(jù)環(huán)境因子對(duì)各徑級(jí)樹木死亡影響的相似性,可將樹木徑級(jí)劃分為<10cm、10—25cm、>25cm(表1)。各徑級(jí)林木死亡率與土壤含水量、土壤容重(PC1)和林齡(PC4)的關(guān)系較弱。徑級(jí)為10—25cm的林木,其死亡率與林木相對(duì)斷面積(PC2)具有極顯著的正相關(guān)關(guān)系(P<0.001),與密度(PC3)具有顯著的負(fù)相關(guān)關(guān)系(P<0.05),固定效應(yīng)系數(shù)分別為1.011和-0.499。徑級(jí)<10cm和>25cm的林木,死亡率與解釋變量之間無顯著的相關(guān)關(guān)系。徑級(jí)為10—25cm時(shí),固定效應(yīng)解釋的變異最大,R2為0.297。不同徑級(jí)的回歸模型中,固定效應(yīng)解釋的變異小于隨機(jī)效應(yīng)解釋的變異。

表1 不同徑級(jí)樹木死亡對(duì)影響因子的響應(yīng)

Table 1 Results of logistic generalized linear mixed-effects model(GLMM)relating gradients of four axes to tree mortality in five DBH size classes

徑級(jí)Classes of DBH /cm固定效應(yīng)系數(shù)Fixed effects coefficients (b)PC1PC2PC3PC4固定效應(yīng)R2Marginal R2隨機(jī)效應(yīng)R2Random effect R2<10-0.0250.355-0.2240.3660.1050.5210—25-0.1271.011???-0.499?0.4520.2970.322>25-0.4861.997-1.9820.1710.1730.766

*P<0.05；***P<0.001;PC1主要由土壤含水量和土壤容重決定；PC2主要由林木相對(duì)斷面積決定；PC3主要由密度決定；PC4主要由林齡決定;DBH: 胸徑Diameter at breast height

2.2 隨機(jī)森林模型預(yù)測(cè)精度

運(yùn)用隨機(jī)森林分別不同徑級(jí)建立樹木死亡預(yù)測(cè)模型。表2的混淆矩陣顯示了隨機(jī)森林模型的預(yù)測(cè)結(jié)果。不同徑級(jí)模型預(yù)測(cè)的總體誤判率為10.7%—31.25%,解釋了樹木死亡70%—90%的變異。其中,模型對(duì)于楊樹防護(hù)林無枯死的誤判率為4.8%—13.1%,對(duì)于枯死的誤判率為30.2%—44.8%。預(yù)測(cè)精度較為滿意。

表2 隨機(jī)森林模型混淆矩陣

2.3 不同徑級(jí)樹木死亡影響因子重要性評(píng)估

分別不同徑級(jí),運(yùn)用隨機(jī)森林評(píng)估11個(gè)影響因子對(duì)樹木死亡的相對(duì)重要性,并比較哪些變量的影響是相似的。由圖1可知,胸徑<10cm時(shí),林木相對(duì)斷面積對(duì)死亡的影響遠(yuǎn)遠(yuǎn)大于其他因子,其次為胸徑和0—20cm土壤水分,密度和林分?jǐn)嗝娣e的影響相對(duì)較小。當(dāng)林木相對(duì)斷面積<0.22或>0.45時(shí),樹木死亡率較大(圖2)。

胸徑為10—25cm時(shí),密度對(duì)樹木死亡的影響最大,其次20—60cm土壤水分、林木相對(duì)斷面積、林分?jǐn)嗝娣e對(duì)樹木死亡的影響相似。當(dāng)密度>600株/hm2時(shí),樹木死亡率隨著密度的增加明顯上升(圖2)。

胸徑>25cm時(shí),20—60cm土壤水分和林分?jǐn)嗝娣e對(duì)樹木死亡的影響最大且相近,其次為林分密度、林木相對(duì)斷面積、0—20cm土壤容重的影響相近,胸徑對(duì)樹木死亡的影響最小。當(dāng)20—60cm土壤質(zhì)量含水量>5%時(shí),樹木死亡率明顯下降。

圖1 影響樹木死亡的自變量的重要性排序Fig.1 Variable importance measure 1：密度Density；2：林木相對(duì)斷面積Relative basal area；3：林分?jǐn)嗝娣eStand basal area；4：胸徑Diameter at breast height；5：林齡Stand age；6: 0—20cm土壤容重，0—20 cm soil bulk density；7: 20—60cm土壤容重，20—60 cm soil bulk density；8：土壤平均容重Average soil bulk density；9: 0—20cm土壤水分，0—20 cm soil moisture content；10: 20—60cm土壤水分，20—60 cm soil moisture content；11：平均土壤水分Average soil moisture content。

圖2 不同徑級(jí)相對(duì)重要變量對(duì)樹木死亡的影響Fig.2 Partial effect of variables on tree mortality in different tree size classes

3 討論

3.1 不同徑級(jí)樹木死亡對(duì)環(huán)境因子的響應(yīng)

本文的研究表明林木相對(duì)斷面積、密度和土壤水分對(duì)該區(qū)域樹木死亡的影響較大,可能與該區(qū)域林木可獲得的資源及受到的大風(fēng)、凍害等惡劣的自然環(huán)境有關(guān)。林分中每株樹木獲得的資源一方面由其獲得資源的能力即競(jìng)爭(zhēng)力決定,一方面由環(huán)境可提供的資源多少?zèng)Q定。胸徑和林木相對(duì)斷面積較大時(shí),林木競(jìng)爭(zhēng)資源的能力較大[19]。林分密度和林分?jǐn)嗝娣e較大時(shí),環(huán)境可提供給每株樹木的資源較少[20]。本文的研究顯示不同徑級(jí)的林木,林木相對(duì)斷面積對(duì)死亡的影響大于胸徑對(duì)死亡的影響,說明林木相對(duì)斷面積可以更好的描述林木之間的競(jìng)爭(zhēng)力,相對(duì)越大的樹木競(jìng)爭(zhēng)力越強(qiáng),這與Luo等[21]的研究結(jié)果相近。

不同徑級(jí)的林木,各影響因子對(duì)樹木死亡的影響權(quán)重不同。胸徑<10cm時(shí),林木相對(duì)斷面積對(duì)樹木死亡的影響最大,其次為胸徑和0—20cm土壤水分,密度和林分?jǐn)嗝娣e的影響最小。說明該階段樹木大小是影響樹木死亡最主要的因子。Zhu等[12]的研究也發(fā)現(xiàn)樹木大小對(duì)小徑級(jí)樹木的死亡有顯著影響。本研究中當(dāng)0.23<林木相對(duì)斷面積<0.43時(shí),樹木死亡率最低(圖2)。當(dāng)林木相對(duì)斷面積<0.23時(shí),樹木死亡率較高。壩上地區(qū)屬半干旱地區(qū),水資源嚴(yán)重缺乏,當(dāng)林木相對(duì)較小時(shí),一方面其競(jìng)爭(zhēng)資源的能力較弱,包括水資源、土壤資源、光照資源等。楊樹屬于速生物種,木材密度較低[22],若得不到充足的水分、營(yíng)養(yǎng)物質(zhì)和光照,將難以維持碳平衡,根系和主干的碳儲(chǔ)存減少,樹干密度降低,抵御病蟲害的能力隨之降低,死亡率增加。另一方面,周圍高大的林木會(huì)造成陰影,不利于矮小的林木獲得光照。同時(shí),有研究表明當(dāng)樹木年均斷面積增長(zhǎng)量達(dá)不到一定值時(shí),死亡率會(huì)增加[9]。本研究中當(dāng)林木相對(duì)斷面積>0.43時(shí),林木死亡率增大。原因可能是壩上地區(qū)大風(fēng)天氣頻繁,相對(duì)較高的林木更易受大風(fēng)凍害的影響,大風(fēng)造成落葉增加,葉損失會(huì)減少樹木光合作用產(chǎn)物,增加非結(jié)構(gòu)性碳消耗,減少樹木體內(nèi)的碳儲(chǔ)備,當(dāng)非結(jié)構(gòu)性碳消耗到一定程度時(shí)會(huì)導(dǎo)致樹木因碳饑餓而死亡[23]。這與Coomes[24]和De Toledo[11]的研究結(jié)果相似,大徑級(jí)林木更易受大風(fēng)、凍害等外界環(huán)境干擾,死亡率相對(duì)較高。

胸徑為10—25cm時(shí),密度對(duì)樹木死亡的影響最大,其次為20—60cm土壤水分、林木相對(duì)斷面積、林分?jǐn)嗝娣e對(duì)樹木死亡的影響相近。說明林分所能提供的資源多少成為該階段樹木生長(zhǎng)的主要限制因子。由于根系的生長(zhǎng),這一階段林木主要利用20—60cm土壤水分。苗博等人[25]的研究中,在生長(zhǎng)季節(jié)胸徑在該范圍的楊樹主要利用30—80cm土壤水分,與本文的研究結(jié)果相似。隨機(jī)森林的結(jié)果顯示,在這一階段密度越大,林木的死亡率越高。因?yàn)槊芏仍酱?來自周圍林木的資源競(jìng)爭(zhēng)越強(qiáng),林分所能提供給單株林木的資源越少,同時(shí)林木也會(huì)缺乏生長(zhǎng)所需的空間。說明這一階段的楊樹都已成林,需要進(jìn)行疏伐。然而邏輯斯蒂?gòu)V義線性模型的擬合結(jié)果顯示,這一階段密度與樹木死亡呈現(xiàn)負(fù)相關(guān)關(guān)系?？赡苁怯捎诒狙芯渴菍?duì)樹木死亡和4個(gè)主成分之間進(jìn)行了擬合回歸,PC3不能完全代表密度。同時(shí)也可能是由于該區(qū)域在冬春季節(jié)大風(fēng)、凍害頻繁,當(dāng)林分的密度較大時(shí),有利于林分抵御不利的自然條件。因此該徑級(jí)范圍內(nèi),密度與樹木死亡的關(guān)系還需進(jìn)一步的探究,從而得到樹木生長(zhǎng)的最適密度。

胸徑>25cm時(shí),20—60cm土壤水分和林分?jǐn)嗝娣e對(duì)樹木死亡的影響最大且相近。說明該階段林木主要利用20—60cm土壤水分,且林分所能提供的水資源多少成為導(dǎo)致樹木死亡的主要因子。在該階段,密度是導(dǎo)致樹木死亡的第三個(gè)主要因子,其影響小于林分?jǐn)嗝娣e,因?yàn)榱址謹(jǐn)嗝娣e是一個(gè)隨著林分生長(zhǎng)而變化的密度指標(biāo),不僅能表示林木的擁擠程度,還能反應(yīng)林分對(duì)資源的需求。

在本研究中林齡對(duì)不同徑級(jí)樹木死亡的影響都較小,影響樹木生長(zhǎng)的主要原因是水資源缺少,大風(fēng)、凍害等惡劣的自然環(huán)境,以及人類不合理的經(jīng)營(yíng)。這與部分學(xué)者認(rèn)為的壩上地區(qū)樹木死亡是樹齡老化的觀點(diǎn)不一致[26],在以后需要對(duì)楊樹是否達(dá)到老齡化這個(gè)問題進(jìn)行進(jìn)一步的研究。

3.2 模型預(yù)測(cè)精度

本研究運(yùn)用隨機(jī)森林算法建立的樹木死亡預(yù)測(cè)模型可以解釋樹木死亡70%—90%的變異,Botter等[27]、De Toledo等[11]、Luo等[21]建立的簡(jiǎn)單線性模型或回歸樹模型可以解釋樹木死亡20%—80%的變異,說明本研究中建立的樹木死亡模型擬合效果較好。邏輯斯蒂回歸對(duì)自變量的多元共線性非常敏感,要求自變量之間相互獨(dú)立,因此本研究在運(yùn)用模型之前將11個(gè)影響因子做了主成分分析,而隨機(jī)森林算法可以評(píng)估所有變量的重要性,不需要顧慮一般回歸問題面臨的多元共線性問題。

本研究只選取了表示林分結(jié)構(gòu)和土壤性質(zhì)的幾個(gè)因子作為解釋變量,在以后的研究中可以加入氣候因子。另外本文只研究了0—60cm土壤水分與樹木死亡的關(guān)系,以后可以增加更深層次土壤水分對(duì)樹木死亡的影響。

4 結(jié)論

本研究基于隨機(jī)森林,分析了壩上張北地區(qū)導(dǎo)致不同徑級(jí)樹木死亡的主要影響因子,并得到了關(guān)鍵因子的影響閾值。結(jié)果表明,林木相對(duì)斷面積、密度、林分?jǐn)嗝娣e、土壤水分是影響該區(qū)域樹木死亡的主要因子。林木胸徑<10cm時(shí),當(dāng)林木相對(duì)斷面積在0.23—0.43之間時(shí),樹木死亡率較低；胸徑為10—25 cm時(shí),密度>600株/hm2,樹木死亡率明顯上升；胸徑>25cm,當(dāng)20—60cm土壤質(zhì)量含水量>5%時(shí),樹木死亡率明顯下降。不同徑級(jí)樹木死亡的主要影響因子不同,但都表明水資源是影響該區(qū)域樹木生長(zhǎng)的關(guān)鍵限制因子。不同階段的林木,有其不同的適宜密度?；谒Y源承載力的林分密度動(dòng)態(tài)調(diào)控是該區(qū)域防護(hù)林經(jīng)營(yíng)的關(guān)鍵。要加強(qiáng)森林的經(jīng)營(yíng)管護(hù),且要遵循自然規(guī)律,不同生長(zhǎng)發(fā)育階段的林木經(jīng)營(yíng)措施不同,才能獲得森林的可持續(xù)發(fā)展。