魏俊,朱坤,2,陳文德,2*,盧雨田,蔡艷琨,吳群,黃仲宣
(1.成都理工大學(xué)旅游與城鄉(xiāng)規(guī)劃學(xué)院,四川成都610059;2.四川省社會科學(xué)重點(diǎn)研究基地(擴(kuò)展)國家公園研究中心,四川成都610059)
香果樹(Emmenopterys henryiOliv.)隸屬茜草科(Rubiaceae)、香果樹屬(EmmenopterysOliv.)落葉大喬木,是古老孑遺植物。在《中華本草》和《新華本草綱要》中記載香果樹具有濕中和胃,降逆止嘔之功效。香果樹種子萌發(fā)力較低,主要靠根系或樹樁萌條更新,繁殖比較困難,尚未發(fā)現(xiàn)純林,只見其星散分布,野外資源量少,1991 年被認(rèn)定為近危保護(hù)植物納入《中國植物紅皮書》,1999 年被列為第一批中國特有二級珍稀瀕危保護(hù)植物,有研究認(rèn)為香果樹是香果樹屬中己知的唯一現(xiàn)存種[1],具有較高的經(jīng)濟(jì)和科研價(jià)值[2]。
目前,關(guān)于香果樹主要研究具體保護(hù)區(qū)等小區(qū)域的種群特征、群落結(jié)構(gòu)、生長動態(tài)、基因組特征、繁殖能力和病害等方面且以中東部省區(qū)為主[3-8],對西南地區(qū)香果樹的地理空間范圍分布情況的研究較少,其中張永華對全國香果樹的親緣地理學(xué)和景觀遺傳學(xué)及其生態(tài)位分化進(jìn)行了研究[1];周文娟等對貴州省境內(nèi)香果樹的分布現(xiàn)狀研究發(fā)現(xiàn)其在貴州屬廣布種,但多為零星分布[9]。
最大熵模型(Maximum entropy model,Maxent)是近年來用作預(yù)測物種在地理尺度上空間分布的重要方法之一。王娟娟等利用Maxent 模型對川貝母的潛在分布和適應(yīng)性進(jìn)行了研究[10];崔晉亮等利用該模型對藍(lán)莓的潛在分布區(qū)進(jìn)行了預(yù)測研究[11];張東方等采用最大熵模型研究了當(dāng)歸在全球的生態(tài)適宜區(qū)和生態(tài)特征[12];王國崢等采用4 種基于不同算法的模型(GARP、Bioclim、Domain 和Maxent)研究金錢松的潛在分布區(qū),結(jié)果表明Maxent 模型對植物潛在地理分布預(yù)測結(jié)果比其他模型更為精確[13]。因此,本研究采用最大熵模型,以對西南地區(qū)涼山山系樣地調(diào)研為基礎(chǔ),結(jié)合香果樹標(biāo)本數(shù)據(jù),對西南地區(qū)(川、云、貴、藏、渝)香果樹在歷史時(shí)期、當(dāng)前時(shí)期和未來時(shí)期適生區(qū)的空間分布情況進(jìn)行預(yù)測分析研究。
采用樣地調(diào)查法對涼山山系香果樹的分布情況進(jìn)行實(shí)地調(diào)研,獲得樣本數(shù)據(jù)點(diǎn)44 個(gè),在全球生物多樣性信息數(shù)據(jù)庫網(wǎng)絡(luò)(GBIF,http://www.gbif.org/)、中國數(shù)字植物標(biāo)本館(CVH,http://www.cvh.ac.cn/)和中國自然保護(hù)區(qū)資源平臺(http://www.papc.cn/)中共查詢到西南地區(qū)香果樹樣本點(diǎn)117 個(gè),其中貴州 19 個(gè)、云南 5 個(gè)、重慶 7 個(gè)、四川 86 個(gè),數(shù)字植物標(biāo)本館中顯示西藏地區(qū)有1 份香果樹標(biāo)本,但是沒有提供坐標(biāo)信息。加上實(shí)地調(diào)研的數(shù)據(jù)一共得到樣本點(diǎn)161 個(gè),統(tǒng)計(jì)發(fā)現(xiàn)有5 個(gè)重復(fù)樣本點(diǎn),去除后共有156個(gè)樣本點(diǎn)(圖1)。
圖1 西南地區(qū)香果樹分布情況Fig.1 Distribution of E.henryiin Southwest China
生態(tài)因子共36 個(gè),其中環(huán)境因子包含19 個(gè)生物氣候因子、7 類氣候因子和1 個(gè)地形因子,環(huán)境因子 數(shù) 據(jù) 來 源 于 Worldclim(http://www. worldclim.org/),此數(shù)據(jù)集為當(dāng)前條件(2000~2020 年)且分辨率為30S(約1km2)的生物氣候變量,這些變量主要反映溫度和降水的特點(diǎn)及其季節(jié)性變化特征,地形因子來源于地理空間數(shù)據(jù)云(http://www.gscloud.cn/),下載分辨率為90 m 的數(shù)字高程地圖(DEM),從中提取海拔,然后在ArcGIS10.2 中重采樣生成分辨率為30″的柵格數(shù)據(jù)圖層[14]。土壤數(shù)據(jù)來自聯(lián)合國糧食及農(nóng)業(yè)組織網(wǎng)站中的世界土壤數(shù)據(jù)庫(Har?monized world soil database v1.2|FAO SOILS POR?TAL),坐標(biāo)系為WGS-1984,柵格大小約為1km2,共選取了9個(gè)土壤因子。
1.3.1 數(shù)據(jù)預(yù)處理
表1 生態(tài)因子縮寫及含義Table 1 Abbreviations and meanings of ecological factors
從國家基礎(chǔ)地理信息系統(tǒng)(http:nfgis.nsdi.gov.cn/)下載獲得的1∶400 000 0 的西南地區(qū)(四川、重慶、云南、貴州、西藏)標(biāo)準(zhǔn)地圖作為分析底圖,再利用ArcGIS 10.2 軟件將西南地區(qū)的生物變量數(shù)據(jù)掩膜出來,將圖層的坐標(biāo)系定義為GCS-WGS-1984,并導(dǎo)出為ASCII格式。將香果樹樣本點(diǎn)數(shù)據(jù)文件保存為“.csv”格式。
1.3.2 生態(tài)因子篩選
由于環(huán)境變量之間具有多重共線性,會導(dǎo)致預(yù)測分布過度擬合[13],用 ArcGIS 10.2 中 ArcTool?box 下多元分析的波段集統(tǒng)計(jì)功能計(jì)算生態(tài)因子的協(xié)方差和相關(guān)矩陣進(jìn)行相關(guān)性檢驗(yàn),篩選出相關(guān)系數(shù)絕對值小于0.75 的生態(tài)因子,并依據(jù)香果樹的生長環(huán)境再次篩選,最終得到4 個(gè)土壤因子:表層土礫石含量、表層土粘土含量、表層土壤有機(jī)碳含量、表層土基本飽和度,13 個(gè)環(huán)境因子:海拔、1 月和 7 月平均降水量、1 月太陽輻射、7 月最高溫、1 月最低溫、年均溫、等溫性、溫度季節(jié)性變化的標(biāo)準(zhǔn)差、年均降水量、最干月降水量、降水量變異系數(shù)、最冷季度降水量,共篩選出17 個(gè)生態(tài)因子。
圖2 生物氣候因子相關(guān)系數(shù)圖Fig.2 Correlation coefficient of bioclimatic factors
1.3.3 最大熵模型
最大熵原理(The principle of maximum entropy)起源于信息論和統(tǒng)計(jì)力學(xué),是基于有限的已知信息對未知分布進(jìn)行無偏推斷的一種數(shù)學(xué)方法[15]。Phil?lips 團(tuán)隊(duì)基于生態(tài)位理論,考慮氣候、海拔、植被等環(huán)境因子,用最大熵原理作為統(tǒng)計(jì)推斷工具,構(gòu)建了物種地理尺度上空間分布的最大熵模型[16-17],編寫了一個(gè)可以免費(fèi)獲取的計(jì)算軟件(http://www.cs.princeton.edu/~schapire/maxent/)[15]。該模型是借助機(jī)器學(xué)習(xí)方法建立的對最大熵原理的應(yīng)用[16],Elith等[18]在2011年對其做了統(tǒng)計(jì)上學(xué)的解釋,使這個(gè)模型在生態(tài)學(xué)領(lǐng)域得到了更加廣泛的應(yīng)用。
最大熵模型軟件的具體操作為:將已知的物種分布數(shù)據(jù)的“.csv”格式文件和篩選出來的生態(tài)因子導(dǎo)入Maxent 軟件,定義結(jié)果輸出的位置以及生態(tài)因子圖層位置,勾選創(chuàng)建響應(yīng)曲線和刀切法復(fù)選框?qū)ψ兞恐匾赃M(jìn)行測試,設(shè)置采取25%的樣本點(diǎn)作為測試集,剩下75%的樣本點(diǎn)作為訓(xùn)練集,其余選項(xiàng)采用模型的默認(rèn)設(shè)定,重復(fù)運(yùn)行10次進(jìn)行建模。受試者工作特征曲線(Receiver operating characteristic?curve,ROC)分析法在物種潛在分布預(yù)測模型評價(jià)中應(yīng)用比較普遍,這一診斷試驗(yàn)評價(jià)指標(biāo)的認(rèn)可度較高[12,19]。計(jì)算 ROC 曲線下方的面積(Area under curve,AUC),其中AUC 的數(shù)值范圍為 0.6~0.7 時(shí)預(yù)測精度較低,0.7~0.8 時(shí)預(yù)測精度中等,0.8~0.9 時(shí)預(yù)測精度良好,大于0.9 時(shí)預(yù)測精度優(yōu)秀;并利用Jackknife刀切圖來評價(jià)各環(huán)境因子的權(quán)重[14]。
將篩選出來的17 個(gè)生態(tài)因子以及香果樹樣本點(diǎn)數(shù)據(jù)分別導(dǎo)入Maxent3.3.3軟件中,重復(fù)運(yùn)行10次上述方法建立最大熵模型,得到生態(tài)因子Maxent 模型平均訓(xùn)練集的AUC 值為0.956(圖3),表明模型預(yù)測結(jié)果為優(yōu)秀,滿足進(jìn)行香果樹地理分布和生態(tài)因子關(guān)系模擬研究的精度要求。
圖3 生態(tài)因子預(yù)測ROC曲線Fig.3 Receiver operating characteristic curve of ecological factor prediction
為了對各區(qū)域的適生程度加以區(qū)分,按照政府間氣候變化專門委員會(Intergovernmental panel on climate change,IPCC)關(guān)于可能性的劃分標(biāo)準(zhǔn)并結(jié)合文獻(xiàn)研究[10],將香果樹潛在分布預(yù)測區(qū)劃分為4個(gè)等級:P < 0.1 為不適宜區(qū);0.1 ≤ P < 0.4 低適宜區(qū);0.4 ≤ P < 0.6 為較適宜區(qū);P ≥ 0.6 為高適宜區(qū)。將最大熵模型模擬導(dǎo)出的“Oliv._avg.asc”文件導(dǎo)入ArcGIS 10.2,利用工具欄將“.asc”文件轉(zhuǎn)為柵格文件并進(jìn)行重分類,得到香果樹適生區(qū)分布圖(圖4)。
圖4 生態(tài)因子模擬香果樹潛在適生區(qū)Fig.4 Ecological factors simulating the potential adaptive area of E.henryi Oliv
利用ArcGIS10.2 統(tǒng)計(jì)分析工具,將模擬結(jié)果的柵格數(shù)據(jù)轉(zhuǎn)換為矢量數(shù)據(jù)并進(jìn)行面積計(jì)算,得到香果樹當(dāng)前潛在分布區(qū)及其適生區(qū)面積。生態(tài)因子模擬結(jié)果中香果樹分布的高適宜區(qū)面積約為3.27×104km2,中適宜區(qū)面積約為 16.78 × 104km2,低適宜區(qū)面積約為24.13 × 104km2,不適宜區(qū)面積約為198.82 × 104km2,其中適宜區(qū)總面積僅占西南地區(qū)總面積的18.19%。
根據(jù)參與模型建立的生態(tài)因子對最大熵模型的貢獻(xiàn)率(表2),可以判斷影響物種分布的主要生態(tài)因子。生態(tài)因子貢獻(xiàn)率結(jié)果顯示共有8個(gè)環(huán)境因子對模擬結(jié)果的貢獻(xiàn)率大于4%,由高到低依次為表層粘土含量(25.25%)、最干月降水量(14.15%)、海拔(13.51%)、表層土礫石含量(13.1%)、1 月太陽輻射強(qiáng)度(11.69%)、表層土基本飽和度(5%)、1 月最低溫(4.38%)、表層土壤有機(jī)碳含量(4.1%),累積貢獻(xiàn)率高達(dá)91.12%。
表2 生態(tài)因子貢獻(xiàn)率Table 2 The percent contribution of ecological factors
2.3.1 香果樹分布區(qū)環(huán)境特征
基于各氣候因子響應(yīng)曲線(圖5)可知,當(dāng)潛在分布概率大于0.5 時(shí),其對應(yīng)的生態(tài)因子值比較適合香果樹的生存[12],曲線中藍(lán)色部分表示10次重復(fù)運(yùn)算的波動范圍,紅線表示10次運(yùn)算的平均值。最干月降水量小于8 mm 時(shí)香果樹適生概率小于0.5;海拔在950~2 200 m 時(shí)其適生概率大于0.5;1 月太陽輻射在6 900~8 100 kJ/m2·day 時(shí)其適生概率大于0.5;1 月最低溫在-4~1℃時(shí)其適生概率大于0.5,香果樹在極端最低溫-15℃時(shí)都可能生存下來;7 月最高溫在21~32℃時(shí)其適生概率大于0.5,在35℃的高溫下仍具有9%的生存概率;降水量變異系數(shù)的范圍為60~95 時(shí),適生概率大于0.5;等溫性的適宜范圍在25~34。
圖5 香果樹的潛在分布概率與環(huán)境主導(dǎo)影響因子的關(guān)系Fig.5 Relationship between the potential distribution probability of E.henryi Oliva environmental dominant impact factors
2.3.2 香果樹分布區(qū)土壤特征
基于各土壤因子響應(yīng)曲線(圖6)可知適合香果樹生存的土壤因子范圍。表層土粘土含量在21%~40%時(shí)香果樹適生概率大于0.5;表層土礫石含量在4%的時(shí)候其適生概率達(dá)到了0.5,隨著表層土礫石含量升高,適生概率短暫下降,礫石含量達(dá)到5.5%之后其適生概率開始增加,當(dāng)?shù)[石含量增加到8%時(shí)適生概率最大,礫石含量適宜范圍為4%~10%;表層土基本飽和度在15~22之間時(shí)香果樹適生概率逐漸下降,之后香果樹存在概率逐漸增加,當(dāng)表層土基本飽和度達(dá)到33時(shí)香果樹存在概率最大,表層土基本飽和度適宜范圍為10~40;香果樹適生概率隨著表層土壤有機(jī)碳含量的增加逐漸增大,表層土壤有機(jī)碳含量達(dá)到1.5%的時(shí)候香果樹存在概率最大,表層土壤有機(jī)碳含量適宜范圍為0.5%~1.5%。
圖6 香果樹的潛在分布概率與土壤主導(dǎo)影響因子的關(guān)系Fig.6 Relationship between the potential distribution probability of E.henryi Oliva soil dominant impact factors
由于香果樹屬于落葉大喬木,其群落演替時(shí)間較長,因此基于歷史氣候數(shù)據(jù)和未來氣候數(shù)據(jù)對其近100 年的時(shí)空分布變化情況進(jìn)行預(yù)測。氣候數(shù)據(jù)來源于 Worldclim(http://www.worldclim.org/),其中歷史氣候數(shù)據(jù)是該網(wǎng)站發(fā)布的1 970~2 000 年世界氣候2.1 版氣候數(shù)據(jù);未來氣候數(shù)據(jù)(2030s、2050s、2090s)來自于該網(wǎng)站 BCC-CSM2-MR 模型ssp585 共享社會經(jīng)濟(jì)路徑,目前提供的未來氣候數(shù)據(jù)的空間分辨率最高為2.5 min,提供的氣候數(shù)據(jù)類型為月最低溫、月最高溫、月平均降水量和19 個(gè)生物氣候因子。篩選出1 月和7 月平均降水量、1月和7 月最低溫、1 月和7 月最高溫、年均溫、等溫性、溫度季節(jié)性變化的標(biāo)準(zhǔn)差、年均降水量、最干月降水量、降水量變異系數(shù)、最冷季度降水量,假設(shè)海拔因子沒有變化,沿用前面的DEM 數(shù)據(jù),共得到14項(xiàng)環(huán)境因子。
將4 個(gè)時(shí)期的數(shù)據(jù)按照前面的方式進(jìn)行處理,并分別用Maxent3.3.3 軟件進(jìn)行模擬預(yù)測,平均訓(xùn)練集的 AUC 值分別為 0.942、0.939、0.933、0.930,預(yù)測結(jié)果達(dá)到了優(yōu)秀水平。利用ArcGIS 10.2 進(jìn)行加工處理,得到 20 世紀(jì)末、2030s、2050s和2090s 各時(shí)期香果樹適生區(qū)分布情況(圖7),同樣利用統(tǒng)計(jì)分析工具計(jì)算4 個(gè)時(shí)期香果樹各等級適生區(qū)面積(圖8)、Maxent3.3.3 軟件輸出結(jié)果中包含環(huán)境變量對最大熵模型的相對貢獻(xiàn)的估計(jì)值(圖9)。
IPCC 第五次評估報(bào)告利用全球氣候耦合模式,模擬發(fā)布了在溫室氣體排放濃度最低情景下,到本世紀(jì)中葉全球平均地表溫度將上升0.4~1.6℃[20]。根據(jù)圖 8 可以看出,在未來的 3 個(gè)時(shí)期,隨著氣候的變化,香果樹在西南地區(qū)的不適宜分布區(qū)面積將會減少,低適宜區(qū)和中適宜區(qū)面積將會有所增加,高適宜區(qū)面積在2030s 和2050s 將會有所增加,但是在2090s 時(shí)期內(nèi)相比于前一時(shí)期將會有所減少。由圖7 可以看出,未來的三個(gè)時(shí)期內(nèi)低適宜區(qū)分布范圍在西藏的東南邊緣和云南的東北部將會有所增加,中適宜區(qū)分布范圍將主要在貴州省境內(nèi)有所增加,2030s 和2050s 高適宜區(qū)分布范圍將在四川的東北邊緣、重慶的東南部和貴州境內(nèi)有所增加,2090s 高適宜區(qū)分布范圍在貴州中西部地區(qū)將會有所減少。
圖7 不同時(shí)期香果樹適生區(qū)的預(yù)測分布圖Fig.7 Distribution map for predicting the suitable areas of E.henryi Oliv in different periods
圖8 不同時(shí)期預(yù)測適生區(qū)分布面積Fig.8 Distribution acreage of adaptive area predicted in different periods
根據(jù)圖9可以看出最干月降水量因子在各時(shí)期貢獻(xiàn)率最大,尤其是2030s、2050s 這兩個(gè)時(shí)期內(nèi)其貢獻(xiàn)率超過了35%;其次是海拔和等溫性因子,在四個(gè)時(shí)期內(nèi)貢獻(xiàn)率都超過了10%;1 月平均降水量、溫度季節(jié)性變化的標(biāo)準(zhǔn)差、7 月最高溫、7 月平均降水量這幾個(gè)因子的貢獻(xiàn)率最小,不超過1%。由此可知,影響香果樹適生區(qū)在西南地區(qū)分布范圍的主要環(huán)境因子是最干月降水量、海拔和等溫性,表明香果樹具有喜濕潤的特點(diǎn)。
圖9 環(huán)境因子在各時(shí)期貢獻(xiàn)率百分比Fig.9 The percent contribution of environmental factors in each period
將 ArcGIS 10.2 與 Maxent3.3.3 結(jié)合,通過 17 個(gè)生態(tài)因子對西南地區(qū)香果樹當(dāng)前的潛在分布區(qū)進(jìn)行預(yù)測,預(yù)測精度達(dá)到優(yōu)秀水平(AUC=0.956),利用ArcGIS 10.2 統(tǒng)計(jì)分析工具計(jì)算得到西南地區(qū)香果樹低適宜區(qū)面積約為21.54×104km2,占西南地區(qū)總面積的9.93%;中適宜區(qū)面積約為16.78 × 104km2,占西南地區(qū)總面積的6.91%;高適宜區(qū)面積約為3.28 × 104km2,占西南地區(qū)總面積的1.35%。香果樹分布中等適宜以上區(qū)域面積占比為8.26%。分析結(jié)果表明西南地區(qū)香果樹的適生區(qū)范圍很小,中高適宜區(qū)主要分布在四川南部、云南東北部、貴州省和重慶市全域范圍內(nèi),分布比較零散,西藏的東南邊緣的唐古拉山區(qū)分布有少量中低適宜區(qū)。
在影響香果樹分布的主導(dǎo)環(huán)境因子中,最干月降水量、海拔、1 月太陽輻射這三項(xiàng)的貢獻(xiàn)率都超過了10%。適宜區(qū)最干月降水量的閾值為5~43 mm,海拔的閾值為133~3 000 m,1 月太陽輻射的閾值為5 030~9 700 kJ/m2·day,1 月最低溫的閾值為-14℃~7℃。在影響香果樹分布的主導(dǎo)土壤因子中,表層粘土含量、表層土礫石含量的貢獻(xiàn)率超過了10%,其中表層土粘土含量的貢獻(xiàn)率達(dá)到了25.25%;適宜區(qū)表層粘土含量的閾值為15%~56%,表層土礫石含量的閾值為0~15%,表層土基本飽和度的閾值為10%~100%,高適宜區(qū)表層土壤有機(jī)碳含量的閾值為0.2%~2.5%。
對香果樹在西南地區(qū)近100 年的適生區(qū)的預(yù)測結(jié)果表明,在未來的幾十年內(nèi)香果樹的適生區(qū)面積整體會有小幅度增加,但主要是中低適宜區(qū)取代小部分不適宜區(qū),高適宜區(qū)面積將僅在2050s 有所增加,到2090s 則會減少至現(xiàn)有水平。其中最干月降水量、海拔和等溫性等環(huán)境因子是主要影響因子。
通過對西南地區(qū)香果樹的潛在分布及適生區(qū)分布的預(yù)測研究,為西南地區(qū)的中國二級珍稀瀕危保護(hù)物種香果樹的樹種分布調(diào)查和保護(hù)提供了較為詳細(xì)的依據(jù);對西南地區(qū)香果樹分布主導(dǎo)生態(tài)影響因子的分析,得出了香果樹生長的生境條件;近100 年的預(yù)測結(jié)果展示了只在氣候環(huán)境變化的影響下,西南地區(qū)香果樹適生區(qū)范圍的變化趨勢,為香果樹的植育保護(hù)提供了參考。