氣候變化情景下基于最大熵模型的青海云杉潛在分布格局模擬

曹雪萍,王婧如,魯松松,張曉瑋,*

1 甘肅農(nóng)業(yè)大學(xué)林學(xué)院, 蘭州 730070 2 蘭州大學(xué)生命科學(xué)學(xué)院, 蘭州 730000

青海云杉(Piceacrassifolia)為松科云杉屬多年生常綠針葉喬木,是我國(guó)青藏高原東北緣特有樹種,主要分布于中國(guó)甘肅、寧夏、內(nèi)蒙古、青海等海拔1600—3800米地帶[1]。青海云杉因其樹形高大通直,材質(zhì)輕軟、紋理直、有彈性,是重要的工業(yè)用材,同時(shí)其抗旱性較強(qiáng),生長(zhǎng)迅速,是青海東部、甘肅北部山區(qū)和祁連山區(qū)優(yōu)良的造林樹種[2-3]。賀蘭山區(qū)是青海云杉的遺傳多樣性中心,而其分布中心主要位于祁連山北坡[4]。在祁連山北坡,青海云杉多呈帶狀、塊狀與草原鑲嵌形成森林-草原復(fù)合生態(tài)系統(tǒng)。上述山區(qū)是我國(guó)西北地區(qū)重要的生態(tài)安全屏障,作為其森林的優(yōu)勢(shì)樹種,青海云杉在維系我國(guó)西北地區(qū)生態(tài)平衡、水土保持、水源涵養(yǎng)和生物多樣性等方面發(fā)揮著重要作用[5- 8]。目前的研究主要集中在青海云杉生理學(xué)、生態(tài)學(xué)、育苗造林、病蟲害防治,以及青海云杉林的林學(xué)特征、水源涵養(yǎng)和生態(tài)效益等領(lǐng)域[9- 11],但對(duì)該物種地理分布格局的研究主要集中于祁連山地區(qū)[8,12-13],而未見有關(guān)該物種整體分布區(qū)的相關(guān)模擬研究。且有關(guān)研究表明,祁連山中部正處于相對(duì)干旱和相對(duì)高溫時(shí)期,并呈現(xiàn)出向暖干化發(fā)展的趨勢(shì)[14],而降水減少和氣候變暖均會(huì)嚴(yán)重影響青海云杉的生長(zhǎng)[15-17]。

植物地理分布是由氣候、水文、土壤、人類活動(dòng)等多種因素長(zhǎng)期共同作用的結(jié)果,其中,氣候是決定植物在區(qū)域尺度上地理分布的最主要因素[18-20]。通過研究氣候與物種分布的關(guān)系,可找出影響物種分布的主導(dǎo)因子,可探索物種分布形成的原因、確定物種的潛在分布區(qū)以及分析未來氣候變化下物種分布區(qū)的變化情況,因此,物種的空間分布模擬一直是生物地理學(xué)的重要研究?jī)?nèi)容之一[12,21]。目前,因相關(guān)研究的側(cè)重點(diǎn)不同,已建立了多種物種分布模型[21-24],其中以最大熵模型(Maxent)的預(yù)測(cè)能力較高,被廣泛應(yīng)用于瀕危物種、入侵物種和重要經(jīng)濟(jì)物種的適生區(qū)預(yù)測(cè)和物種保護(hù)方面[25-26]。例如,獨(dú)葉草(Kingdoniaunioflora)、毛紅椿(Toonaciliatevar.pubescens)、紅豆杉(Taxuschinensis)和水葫蘆(Eichhorniacrassipes)等植物的研究[27- 30]。

本研究以青海云杉為研究對(duì)象,基于氣候相似性原理,利用最大熵模型篩選影響青海云杉地理分布的主導(dǎo)因子,給出其潛在的地理分布范圍,并預(yù)測(cè)其在未來氣候變化背景下潛在分布區(qū)的變化,從而深化未來氣候變暖對(duì)其潛在分布區(qū)變化影響的認(rèn)識(shí),以期為青海云杉的經(jīng)營(yíng)管理及應(yīng)對(duì)氣候變化提供依據(jù)。

1 材料方法

1.1 收集物種現(xiàn)有分布數(shù)據(jù)

青海云杉的現(xiàn)有分布樣點(diǎn)通過在線查閱中國(guó)數(shù)字植物標(biāo)本館(http://www.cvh.ac.cn/)中已有記載的標(biāo)本樣點(diǎn)和收集實(shí)驗(yàn)室已有野外采集樣點(diǎn)而得,并參考《中國(guó)植物志》中所述該物種的分布范圍,對(duì)分布點(diǎn)進(jìn)一步篩選,去除不確定和重復(fù)樣點(diǎn),最終確定69個(gè)樣點(diǎn),用于后續(xù)分析(見表1)。

表1 青海云杉樣點(diǎn)地理信息

1.2 當(dāng)前及未來潛在分布區(qū)預(yù)測(cè)

1.2.1獲取氣候數(shù)據(jù)

當(dāng)前氣候圖層下載自WorldClim 網(wǎng)站(http://www.worldclim.org/),該網(wǎng)站提供了1970—2000年30間與降雨和氣溫有關(guān)的19個(gè)氣候變量：年均溫(bio1),月均溫度變幅(bio2),等溫性(bio3),溫度季節(jié)性變化(bio4),最暖月最高溫(bio5),最冷月最低溫(bio6),溫度年變幅(bio7),最濕季均溫(bio8),最干季均溫(bio9),最暖季均溫(bio10),最冷季均溫(bio11),年降雨(bio12),最濕月降雨(bio13),最干月降雨(bio14),降雨季節(jié)性變化(bio15),最濕季降雨(bio16),最干季降雨(bio17),最暖季降雨(bio18),最冷季降雨(bio19)。

未來20世紀(jì)50年代和80年代(2050s和2080s)的氣候圖層下載自CCCFS網(wǎng)站(http://www.ccafs-climate.org),氣候增溫情景則參照了Ren等人的方法[31],從CMIP5計(jì)劃中選取了3種大氣環(huán)流模型(General circulation models, GCMs: BCC-CSM1- 1、 CCCma_CanESM2和CSIRO-Mk3.6.0)和3種代表性濃度路徑情景(Representative concentration pathways, RCPs: RCP2.6、RCP4.5和RCP 8.5))。這3種情景代表了不同的溫室氣體排放策略：RCP2.6下通過提高能源使用率和使用可再生新能源替代化石燃料等方法降低溫室氣體排放量,到2100年控制中國(guó)平均氣溫增幅低于2℃(緩和情景)；RCP4.5和RCP8.5則分別為中度穩(wěn)定和高度溫室氣體排放情景,兩種情景下,到2100年中國(guó)平均氣溫分別增加2.88℃和5.51℃[32-33]。以上數(shù)據(jù)圖層均采用2.5′的空間分辨率。

1.2.2潛在分布區(qū)預(yù)測(cè)

圖1 基于Maxent 模型預(yù)測(cè)青海云杉當(dāng)前潛在分布的受試者工作特征曲線(ROC)曲線Fig.1 The receiver operating characteristic curve (ROC) curve predicting the current potential distribution of Picea crassifolia based on Maxent圖中AUC表示受試者工作特征曲線下面積Area under the receiver operating characteristic curve

物種潛在分布區(qū)的預(yù)測(cè)基于最大熵模型(Maxent)來完成。Maxent通過物種現(xiàn)有分布數(shù)據(jù)可構(gòu)建物種當(dāng)前潛在分布模型,并能預(yù)測(cè)物種在氣候變化下未來的潛在分布,是目前應(yīng)用最廣泛的物種分布模型[34-35]。將以上收集的青海云杉現(xiàn)有分布數(shù)據(jù)和當(dāng)前及未來2050s和2080s氣候圖層導(dǎo)入軟件Maxent,使用軟件默認(rèn)參數(shù),計(jì)算得到青海云杉當(dāng)前和未來潛在分布模型。在此過程中,從數(shù)據(jù)基中抽取75%作為訓(xùn)練數(shù)據(jù),剩余25%作為交互數(shù)據(jù)以評(píng)估模型的有效性。模擬結(jié)果的精確度通過受試者工作特征曲線(Receiver operating characteristic curve,ROC)進(jìn)行評(píng)估,用曲線下面積(Area under the curve, AUC)表示,AUC值介于0—1.0,其值為0.5—0.7時(shí)模型預(yù)測(cè)精度較差,0.7—0.9時(shí)為一般,值高于0.9時(shí)預(yù)測(cè)精度較好。各氣候因子對(duì)青海云杉物種潛在分布的貢獻(xiàn)率通過Jackknife模塊進(jìn)行評(píng)估。

此外,結(jié)合植物志描述的實(shí)際分布區(qū)和野外調(diào)查結(jié)果,本研究將青海云杉當(dāng)前和未來潛在分布區(qū)域劃為4個(gè)等級(jí)：<0.1為非適生區(qū)；0.1—0.3為低度適生區(qū)；0.3—0.5為中度適生區(qū)；0.5—0.8為高度適生區(qū)。通過分別計(jì)算各等級(jí)對(duì)應(yīng)分布區(qū)面積,以明確氣候變化對(duì)青海云杉實(shí)際分布的影響。另外為對(duì)比青海云杉未來潛在分布區(qū)總面積與當(dāng)前的差異,本研究在進(jìn)行模型運(yùn)算時(shí)采用了第10分位訓(xùn)練存在閾值(10thpercentile training presence threshold),以獲取物種分布的二元分布圖像(僅指示適生與非適生兩種狀態(tài))。以上分布區(qū)面積的計(jì)算和圖形繪制皆在ArcGIS 10.3中完成。

2 結(jié)果

2.1 影響我國(guó)青海云杉地理分布的主要?dú)夂蛞蜃?/h3>

基于Maxent 模型預(yù)測(cè)青海云杉當(dāng)前潛在分布的ROC曲線顯示,其訓(xùn)練集和測(cè)試集的曲線下面積(AUC)均高于0.99,表明該模型對(duì)青海云杉潛在分布預(yù)測(cè)的準(zhǔn)確度很高(圖1)。通過Jackknife模塊對(duì)各氣候因子的貢獻(xiàn)率進(jìn)行評(píng)估可知,等溫性(bio3)和最冷季降雨(bio19)的訓(xùn)練得分值為最高,其次為最干季降雨(bio17)、最干月降雨(bio14)、降雨季節(jié)性變化(bio15)、最冷月最低溫(bio6)、最干季均溫(bio9)和最冷季均溫(bio11)(圖2)。繪制以上8個(gè)氣候因子對(duì)青海云杉適生度的響應(yīng)曲線可以發(fā)現(xiàn),這些曲線均存在顯著峰值,且最干月、最干季和最冷季降雨的適生范圍顯著低于其他因子,表明一年中的最低降雨量是影響青海云杉分布的主導(dǎo)因子(圖3)。

圖2 各氣候因子對(duì)青海云杉當(dāng)前潛在分布模型的Jackknife檢驗(yàn)得分Fig.2 Jackknife test for climate factor contributions in the current potential distribution of Picea crassifolia圖中All variables為全部環(huán)境變量; bio1: 年均溫Annual mean temperature; bio2: 月均溫度變幅Mean diurnal range; bio3: 等溫性Isothermality; bio4: 溫度季節(jié)性變化Temperature seasonality; bio5: 最暖月最高溫Max temperature of warmest month; bio6: 最冷月最低溫Min temperature of coldest month; bio7: 溫度年變幅Temperature annual range; bio8: 最濕季均溫Mean temperature of wettest quarter; bio9: 最干季均溫Mean temperature of driest quarter; bio10: 最暖季均溫Mean temperature of warmest quarter; bio11: 最冷季均溫Mean temperature of coldest quarter; bio12: 年降雨Annual precipitation; bio13: 最濕月降雨P(guān)recipitation of wettest month; bio14: 最干月降雨P(guān)recipitation of driest month; bio15: 降雨季節(jié)性變化Precipitation seasonality; bio16: 最濕季降雨P(guān)recipitation of wettest quarter; bio17: 最干季降雨P(guān)recipitation of driest quarter; bio18: 最暖季降雨P(guān)recipitation of warmest quarter; bio19: 最冷季降雨P(guān)recipitation of coldest quarter

圖3 主要?dú)夂蛞蜃拥捻憫?yīng)曲線Fig.3 Response curves of main climate factors

2.2 青海云杉未來潛在分布區(qū)變化

在3種情景、3種大氣環(huán)流模型(GCMs)和兩個(gè)時(shí)間段(2050s和2080s)下,通過模擬得到18個(gè)青海云杉的未來潛在分布模型,且所有模型的平均AUC測(cè)試值均高于0.99,表明模型的預(yù)測(cè)精度極高。對(duì)比潛在分布區(qū)面積變化發(fā)現(xiàn)(表2),盡管在3種增溫情景下,青海云杉在2050s和2080s的總分布區(qū)面積與當(dāng)前相比變化并不大(約1%),但不同適生等級(jí)的潛在分布面積變化顯著。其中,0.1—0.3和0.3—0.5(低度和中度適生區(qū))兩個(gè)等級(jí)的分布受氣候增溫影響顯著,在低度和高度增溫情景(RCP2.6和RCP8.5)下這兩個(gè)適生等級(jí)在2080s的分布區(qū)面積大幅下降(>30%),而中度增溫情景(RCP4.5)下,低度和中度適生區(qū)面積則顯著增加(>10%)。與此相對(duì),青海云杉高度適生區(qū)的面積在所有情景下,與當(dāng)前相比均呈減少趨勢(shì)(表2)。從分布圖上看,3種增溫情景下,與當(dāng)前相比,在2050s青海云杉在青藏高原中部及南部的潛在分布區(qū)(低度適生區(qū))大幅減少,而增加了其在中國(guó)東部的潛在分布(低度和中度適生區(qū)),呈現(xiàn)出潛在分布區(qū)向東遷移的趨勢(shì)；在2080s,青海云杉在青藏高原的潛在分布(低度適生區(qū))進(jìn)一步減少,而東部邊緣地區(qū)(低度和中度適生區(qū))依然呈擴(kuò)張趨勢(shì)(圖4)。盡管如此,青海云杉的核心分布區(qū)域(高度適生區(qū))仍然以青海東部、甘肅北部為主,無明顯變遷趨勢(shì)(圖4)。

3 討論

植物在區(qū)域尺度上的地理分布主要受氣候因子的制約,其中水熱條件起主導(dǎo)作用。本研究通過Maxent模型的Jackknife模塊對(duì)影響青海云杉分布的氣候因子進(jìn)行評(píng)估后,證明了年最低降雨量是影響青海云杉分布的主要?dú)夂蛞蜃印Ｅc此不同,Xu等人利用最大熵模型模擬青海云杉在祁連山的分布后得出影響其分布的主要?dú)夂蛞蜃邮亲钆伦罡邷睾妥顫窦揪鶞豙36]。這種結(jié)果的出現(xiàn)可能與預(yù)測(cè)分布尺度的差異有關(guān)：在本研究中,對(duì)青海云杉潛在分布的預(yù)測(cè)是基于全球尺度,遠(yuǎn)大于祁連山尺度。前人研究表明,青海云杉對(duì)氣候因子的響應(yīng)受海拔因素的影響,在森林上限(海拔3300 m)溫度是影響青海云杉徑向生長(zhǎng)的主要因素,而在森林下線(海拔2770 m)其徑向生長(zhǎng)主要受降雨量的影響[37]。因此,我們推測(cè)在區(qū)域水平上,青海云杉的分布受到溫度的影響,但其總體分布仍然受限于最低年降雨量。

Maxent模型預(yù)測(cè)結(jié)果顯示,當(dāng)前氣候條件下,青海云杉的潛在分布區(qū)主要集中于：青海東部、甘肅東南部、寧夏大部分地區(qū)、西藏東部、四川西部山區(qū)以及陜西、新疆和內(nèi)蒙古部分地區(qū)。此預(yù)測(cè)分布范圍大于青海云杉的實(shí)際分布范圍。原因可能來自于多方面：1)隨著人類活動(dòng)的加劇導(dǎo)致以青海云杉為優(yōu)勢(shì)樹種的森林生態(tài)系統(tǒng)遭受了嚴(yán)重的破壞[38],僅1989年到1998年期間祁連山保護(hù)區(qū)的森林植被就減少了33.5%[39]。別強(qiáng)等[8]通過分布1960年至2000年期間祁連山保護(hù)區(qū)青海云杉林分布的動(dòng)態(tài)變化發(fā)現(xiàn),自20世紀(jì)60年代至本世紀(jì)初青海云杉面積顯著減少,減少率最大時(shí)達(dá)到11.7%；2)青海云杉作為青藏高原東北邊緣特有樹種,主要集中分布于海拔2500—3200 m的范圍[38],而本研究重點(diǎn)探討氣候因子與物種分布的關(guān)系,沒有考慮分布海拔、坡向等地理空間因素,可能會(huì)造成預(yù)測(cè)的潛在分布區(qū)域大于實(shí)際已知生長(zhǎng)區(qū)域的結(jié)果；3)Maxent模型基于氣候相似性進(jìn)行模擬,但并未考慮到物種現(xiàn)實(shí)分布時(shí)受到的其他因素,如物種的生長(zhǎng)特性、擴(kuò)散和遷移能力、自然更新能力、生物間的相互作用等[40-41]。盡管如此,利用Maxent模型得出的青海云杉潛在分布區(qū)域均具有相似的氣候特點(diǎn),說明青海云杉林具有很高的恢復(fù)潛力[42],可為上述地區(qū)進(jìn)行森林生態(tài)系統(tǒng)修復(fù)與恢復(fù)提供一定的樹種參考價(jià)值。

氣候變化背景下,盡管青海云杉的總潛在分布面積波動(dòng)較小(約1%),但其在不同適生等級(jí)的潛在分布面積及分布范圍均發(fā)生了改變。這種差異的產(chǎn)生可能與本研究中對(duì)潛在分布的等級(jí)劃分有關(guān),總適生面積是通過Maxent第10分位訓(xùn)練存在閾值生成的二元圖計(jì)算而來,而適生等級(jí)則是人為劃分,與前者的適生范圍存在一定差異。在不同適生等級(jí),雖然中度及以下的適生區(qū)在未來3種增溫情景下兩個(gè)時(shí)間段(2050s和2080s)的面積變化不一致,但高度適生區(qū)的面積均呈下降趨勢(shì),這表明未來氣候增溫仍然影響了青海云杉的分布。前人對(duì)祁連山青海云杉過去近30年對(duì)氣候響應(yīng)的調(diào)查分析發(fā)現(xiàn),青海云杉在全球變暖背景下增強(qiáng)了對(duì)溫度的敏感性,這種敏感性與區(qū)域內(nèi)水分含量密切相關(guān),干旱會(huì)增加青海云杉與溫度負(fù)相關(guān)的敏感性[43]。且有研究表明青海云杉的分布地未來氣候?qū)⑾驕嘏稍镛D(zhuǎn)化[14],這也間接支持了與本研究的結(jié)果,即青海云杉的核心分布范圍呈下降趨勢(shì)。梁和延對(duì)我國(guó)21世紀(jì)末期氣候變化預(yù)測(cè)發(fā)現(xiàn),在未來氣候顯著變暖情景(RCP8.5)下,中西部地區(qū)的夏季降雨量明顯下降[33]。這也解釋了本研究中青海云杉在該情境下本世紀(jì)末期(2080s)在中國(guó)中西部的分布范圍顯著低于其他情景及時(shí)間段的現(xiàn)象(圖4)。此外,在不同增溫情境下,不同適生等級(jí)的青海云杉適生區(qū)在2050s和2080s的潛在分布區(qū)面積變化趨勢(shì)并不一致,這可能與大氣環(huán)流模型存在地區(qū)氣候變幅差異有關(guān)[44-45]。

表2 青海云杉不同適生等級(jí)未來潛在分布與當(dāng)前相比面積變化百分比/%

Table 2 The percentage change in potential environmental area forPiceacrassifoliain the future when compared with the current at different level of fitness

適生等級(jí)Level of fitnessRCP 2.6RCP 4.5RCP 8.52050s2080s2050s2080s2050s2080s<0.1-8.10+26.62-19.23-9.46-9.24+43.560.1—0.3+5.83-13.81+2.20+10.19-1.83-22.500.3—0.5+20.54-32.14+13.22+20.764.74-37.890.5—0.8-9.48-16.98-21.73-3.72-12.80-6.25總面積變化 Percentage change in total area-0.87-0.71-0.53-0.77-1.13-0.34

“+”和“-”分別指示“增加”和“減少”

圖4 在3種增溫情景下,青海云杉當(dāng)前及未來不同適生等級(jí)的潛在分布區(qū)變化Fig.4 The change of potential distribution for Picea crassifolia in the future when compared with the current at different fitness under three greenhouse gas emission scenarios

4 結(jié)論

本研究基于最大熵模型證明了年最低降雨量是影響青海云杉分布的主要?dú)夂蛞蜃樱辉谖磥須夂蜃兓榫跋?青海云杉分布范圍雖然有減小趨勢(shì),但其高度適生區(qū)仍以祁連山、賀蘭山區(qū)為主。因此在未來西部高山地區(qū),特別是作為我國(guó)重要生態(tài)屏障的祁連山、賀蘭山區(qū),青海云杉依然具有重要的經(jīng)濟(jì)價(jià)值并將持續(xù)其生態(tài)服務(wù)功能。