基于MaxEnt優化模型的閩楠潛在適宜分布預測

葉興狀,張明珠,賴文峰,楊淼淼,范輝華,張國防,陳世品,劉 寶,*

1 福建農林大學林學院,福州 350002

2 福建師范大學地理科學學院,福州 350007

3 福建省林業科學研究院特色資源研究所,福州 350012

氣候變暖影響全球水熱分布和土壤微氣候,而植物分布格局主要受水份、光照、土壤、海拔、坡向和CO2濃度等因子的制約[1]。近年來,氣候變暖對植物分布的影響是植物學、生態學及其交叉學科的研究熱點之一。第四次工業革命以來,全球氣候變暖趨勢加劇,導致其中部分物種的適生區面積大大減小和生境的破碎化,致使許多物種瀕危甚至滅絕,嚴重影響生態安全和生物多樣性[2- 3]。據統計,目前世界上約六分之一的物種受不同程度的威脅,有研究表明,隨著未來全球氣溫的升高,瀕危物種的滅絕風險將大大增加[4]。因此,在全球氣候變暖的背景下,探索瀕危物種的分布格局如何受未來氣候變化的影響,將有助于制定瀕危物種的保護策略及種質資源的保存計劃,以期規避滅絕風險,維持人類賴以生存的生態環境安全。

物種分布模型是研究物種適生區域受氣候環境變化影響的重要工具[5]。目前,國內外學者采用不同的算法開發了如DOMAIN(領域模型)、GARP(遺傳算法模型)和MaxEnt(最大熵值模型)等物種分布模型,其中,MaxEnt模型用于處理群體性溫度因子及降水因子的失真度小、穩定性好,是一種理想的預測模型,也是目前應用最普遍的物種分布模型之一[6],已被廣泛應用于瀕危植物適生區變化預測,如紅豆樹[7](Ormosiahosiei)、半楓荷[8](Semiliquidambarcathayensis)和四合木[9](Tetraenamongolica)等。因此,采用MaxEnt模型研究氣候變化對瀕危植物適生區的影響,進而制定相應的保護策略具有重要參考價值。

閩楠(Phoebebournei)屬樟科(Lauraceae)楠屬(Phoebe)植物,是我國特有二級瀕危珍稀樹種,其樹干高大通直,材質致密,是不可多得的珍貴用材樹種,同時還是實現碳中和的理想樹種,具有很高的經濟和生態價值[10-11]。天然閩楠林因歷年遭受嚴重的濫砍濫伐,導致其野外生境破碎化現象異常嚴重,現有天然林資源幾近枯竭[10- 12]。與歷史記載相比,其分布范圍嚴重收縮,除福建沙縣羅卜巖自然保護區和湖南平江等少數地點的天然林外,野外已難見集中成片的閩楠天然林[13- 14],因此,加大對天然閩楠資源的保護力度已迫在眉睫。此外,我國南方林區多為杉木或馬尾松人工純林,林分結構單一、病蟲害日趨嚴重、森林生態系統功能穩定性差,我國戰略儲備林規劃指出,應對戰略儲備樹種采用培育復層異齡混交林等措施,以期提高我國森林資源質量,而閩楠是理想的闊葉造林樹種之一[15]。然而,目前閩楠的研究主要集中在育苗造林、群落生態、良種選育和遺傳多樣性等方面[15- 19],目前還尚未有關于閩楠適生區預測的相關報道,而適生區模擬預測將有利于揭示閩楠未來分布變化規律及影響其分布的關鍵環境因子,對閩楠核心種質資源庫完善、資源保育和引種栽培具有重要意義。

本研究在野外實地調查和前人對閩楠資源分布調查的基礎上,應用基于Enmeval數據包輔助選擇MaxEnt最優參數設置來模擬推測閩楠的地理分布,同時分析影響閩楠分布的關鍵環境變量,并進一步分析其適生區的環境特征,以期為閩楠的資源保護、天然種群恢復和引種栽培提供科學依據。

1 材料與方法

1.1 研究區域

閩楠是我國特有珍貴樹種,其天然分布北界自湖北沿神農架、大別山直至江西北部山丘,南界為南嶺山脈的廣東、廣西丘陵至花坪國家自然保護區,東界為浙江雁蕩山、烏巖嶺至福建東北部丘陵,西界自重慶西部山區、貴州中西部山區至廣西西北部山區[17]。閩楠天然分布區的地理坐標范圍為23°—30°N,105°—122°E[16],考慮到要為科學引種提供參考,因此以100°—130°E,10°—35°N為預測區域。

1.2 樣本數據的搜集與篩選

課題組在1998—2019年陸續對福建、浙江、江西、湖南、廣西和廣東等6省部分閩楠天然種群進行了實地調查,搜集分布記錄133條；查閱已出版的文獻資料,搜集閩楠分布記錄46條[11- 19]；檢索全球物種多樣性信息庫(http://www.gbif.org)、國家標本平臺(http://www.nsii.org.cn)、教學標本資源共享平臺(http://mnh.scu.edu.cn/main.aspx)、中國數字植物標本館(http://www.cvh.org.cn)分別搜集閩楠記錄170條、635條、113條和653條,合計1750條。參照李垚[20]的方法對上述全部記錄進行篩選,刪掉非天然分布記錄及重復點,并刪除2000年之前的記錄點,為降低群集效應帶來的誤差,每個網格(10 km×10 km)只保留一個分布點,最終得到有效樣本186個。

1.3 環境變量篩選和數據處理

氣候數據選擇來自于WorldClim數據庫(http://worldclim.org)的第六次國際耦合模式比較計劃(CMIP6)模式的共享經濟路徑。為提高預測準確性,選取現代及未來兩個時期(2050s和2090s)各3個氣候情景(SSP126、SSP245和SSP585)進行等模擬,空間分辨率為2.5 arc-minutes。本研究氣候數據采用第二代國家氣候中心中等分辨率氣候系統模式(BCC-CSM2-MR)生成,該模式對中國降水與氣溫的模擬較為可靠[21]。采用Spearman Correlation Analysis和VIF方差膨脹因子,選擇相關系數小于0.8的因子。生態極限和植物的實際分布密切相關,野生閩楠多見于山地溝谷闊葉林中,海拔在40—1 220 m,且多分布于海拔600 m以下丘陵及低山,結合植物生存環境需求及其生理特性,在相關系數大于0.8的變量中保留最重要的生態因子[17]。最終,保留9個變量用于建模,包括年均溫(bio1)、晝夜溫差日均值(bio2)、等溫性(bio3)、最暖月最高溫(bio5)、最冷月最低溫(bio6)、最濕季度平均溫度(bio8)、最干月降水量(bio14)、降水量變異系數(bio15)和最暖季度降水量(bio18)。

1.4 模型建立、優化和評估

為使閩楠出現的概率接近正態分布,選擇75%的數據用于模型訓練,剩下的數據用于模型測試,其他值為默認值。使用Rv3.61中的Emneval數據包優化MaxEnt模型[22],將調控倍頻設置為0.5—5,每次間隔0.5,一共10種調控倍頻；MaxEnt模型提供了5種特征,分別是線性特征、二次型特征、片段化特征、 乘積型特征和閾值性特征。采用8個特征組合。Enmeval 數據包將上述 80種參數組合進行測試,最終采用最小信息準則AICc(The minimum information criterion AICc value, delta.AICc)模型,采用訓練AUC和測試AUC值之差(Difference between the AUC values, avg.AUCDIFF)檢驗模型的擬合度與復雜度[23]。

優化完成后,采用優化參數進行模擬和預測閩楠不同時期適生區。利用受試者曲線下的面積評估MaxEnt預測的準確性。AUC值的取值范圍為0—1,值越大意味著模擬可信度越高,AUC=0.5表示該模型判斷完全失敗,小于0.8表示可信度低,0.9—1.0表示極準確[7- 9]。采用平均邏輯值將適宜度分為(Ⅰ：0—20%；Ⅱ：20%—40%；Ⅲ：40%—60%；Ⅳ：60%—100%)4級,分別表示非適生區、一般適生區、中度適生區和高度適生區。

1.5 物種適生區空間格局變化

本研究將邏輯值Q≥40%作為潛在適生區,即適生區包括了中度適生區和高度適生區,用“1”表示；將邏輯值Q<40%作為非適生區,用“0”表示。本研究定義了4種物種適生區變化的類型：新增適生區、喪失適生區、保留適生區和非適生區。未來面積變化均基于當前閩楠適生區面積來比較計算。矩陣值0→1為新增適生區,1→0為喪失適生區,1→1為保留適生區,0→0為非適生區。具體方法參考張殷波[24]的方法,最后將矩陣變化值加載入Arc GIS10.4中,實現閩楠適生區空間格局變化的可視化表達。

1.6 多元環境相似度面和最不相似變量分析

參考李垚[19]的方法,裁剪閩楠現代潛在適生區的環境變量作為參考圖層,采用多元環境相似度面和最不相似變量分析探討未來6個溫室氣體排放情境下的氣候異常區域及引起潛在適生區改變的關鍵氣候因子。多元相似度能夠反映特定時期某點的氣候狀況與參考圖層氣候狀況的相似程度,負值表明該點至少有一個環境變量的數值超出了參考圖層對應值的范圍,稱之為氣候異常點,最大值100代表該點的氣候完全正常。某點的最不相似變量指的是相似度最低(即異常程度最高)的變量,異常程度最高的環境變量很可能是引起地理分布變遷的關鍵因素[24]。

2 結果與分析

2.1 模型優化及準確性評價

本研究基于閩楠186個分布點和9個生態因子,運用Enmeval 數據包調用MaxEnt對閩楠現代及未來的潛在分布區進行預測。當FC(特征組合)=PT(乘積型特征),RM(調控倍頻)=0.5時,delta.AICc=0,該參數下模為最優模型(表1)。優化后的最小信息準則AICc(the minimum information criterion AICc value, delta.AICc)和訓練集AUC與測試集AUC差值(Difference between the AUC values, avg.diff.AUC)均低于默認設置,表明優化后的參數有利于減低模型的擬合度和復雜度,故本研究選取FC=PT,RM=0.5作為建模的參數設置。在該參數下,MaxEnt預測結果顯示受試者工作曲線AUC值為0.9846±0.0037,表示預測結果極準確。

表1 MaxEnt模型Enmeval優化評價指標Table 1 Evaluation metrics of MaxEnt model generated by Enmeval

2.2 環境變量的重要性

從表2可知,貢獻率最重要的6個環境因子累計占比達86%,依次分別為最干月降水量(bio14,37.4%),適宜范圍為大于17.93 mm；降水量變異系數(bio15,15.8%),適宜范圍為大于46.75—68.35；最冷月最低溫(bio6,12%),適宜范圍為-3.32—6.31℃；最暖季降水量(bio18,11.1%),適宜范圍為大于72.32 mm；濕季度平均溫度(bio8,10.2%),適宜范圍為40.69—81.86；晝夜溫差日均值(bio2,9.6%),適宜范圍為小于8.89 ℃。置換重要值位居前6位的最冷月最低溫(bio6)、最干月降水量(bio14)、最濕季度平均溫度(bio8)、海拔(elev)、等溫性(bio3)和晝夜溫差日均值(bio2)。

表2 閩楠主要環境因子參數Table 2 Environmental variables and their main parameters of P. bournei

2.3 現代和未來潛在適生區

186條分布記錄的平均邏輯值為0.49,最高值為0.76(湖南藍山),最低值為0.11(廣西田林),分別屬于高度適宜區、中度適宜區、一般適宜區和不適宜區的分布點比例分別為21.76%、57.41%、16.66%和4.17%。從表3可知,當前潛在適生區的面積為54.32×104km2,其中高度適生區面積和中度適生區面積分別為3.37×104km2和21.60×104km2。如圖1,當前氣候條件下閩楠高度適生區主要位于廣西東北部、湖南南部、廣東北部和福建中北部局部。中度適生區主要分布于福建中東部、江西中南部、貴州東南部、湖南西部及南部、廣東北部和廣西東北部。

表3 不同時期閩楠適生區面積變化/(×104 km2)Table 3 Changes of suitable area of P. bournei in different periods

圖1 閩楠現代適生區Fig.1 Modern suitable habitats area of P. bournei

2.4 未來潛在適生區空間變換格局

將6個未來不同氣候情景下閩楠適生區的空間格局與現代閩楠適生區進行對比分析(表4、圖2和圖3),結果表明：新增適生區在2050s-SSP585排放情景下分布面積最大,新增率達到81.32%,新增面積主要集中在湖南中西部、貴州東部；中高度適生區約為42.17×104km2,為現代的1.69倍,保留面積21.78×104km2,保留率達87.21%；潛在適生區喪失面積最少,喪失3.01×104km2,喪失率僅為12.06%,主要喪失區域為廣西北部、貴州中部和江西中部。

表4 不同氣候情景下閩楠適生區空間變化Table 4 Changes in suitable distribution area of P. bournei different climate scenarios

圖2 MaxEnt模型預測的不同氣候情景下閩楠潛在適生區Fig.2 MaxEnt model predicted potential suitable growth areas of P. bournei under different climate scenarios

圖3 不同氣候情景下閩楠適生區空間變換格局 Fig.3 Suitable distribution of P. bournei under different climate change

2090s-SSP126排放情景下,新增面積最小,新增率為29.32%,新增面積主要集中于湖南南部、江西西部以及福建中西部,中高度適生區為22.42×104km2；該時期潛在喪失面積最大,達到10.07×104km2,喪失率為40.32%,主要喪失區域包括貴州中南部、廣西東北部、湖南西部、江西南部和浙江南部。

2090s-SSP245和2090s-SSP585情景下,保留適生區面積分別為第二及第三大,保留面積分別為17.22、16.73×104km2,保留率分別為68.96%、66.99%。兩個情景的保留區均主要位于閩東北、桂東北和湘南。兩個情景的新增適生區均主要位于黔東北及湘西,喪失區均主要位于黔東南、桂東北及浙東南。前者適生區較現代增加3%,而后者較現代減少3.1%。

2050s-SSP126和2050s-SSP245情景下,保留適生區面積幾乎相同,均較小,保留面積分別到16.45×104、16.57×104km2,保留率分別為65.89%、66.35%。這兩個時期喪失區較大,前者喪失區主要位于閩浙交界、贛中南部及湘西,后者喪失區主要位于黔東、贛中南部及湘西。前者新增適生區主要位于閩西、湘西和黔東,而后者新增適生區為閩西和湘西,前者適生區比現代略有增加,后者適生區比現代增長6%。

2.5 多元環境相似度面和最不相似變量分析

如圖4所示,在2050s-SSP126、2050s-SSP245、2050s-SSP585、2090s-SSP126、2090s-SSP245和2090s-SSP585等6個排放條件下,閩楠186個現代分布點的平均多元相似度分別為2.53、19.17、10.03、15.46、7.46和-2.06,多元相似度為負值的點比例分別為18.98%、5.56%、18.98%、6.48%、26.85%和61.11%,表明2090s-SSP585排放條件下,氣候異常程度最高,2050s-SSP126情景次之,2050s-SSP245情景下氣候異常程度最低。在2050s-SSP126情景下,現代適生區內的氣候異常區域主要出現在中北部,最不相似變量為bio3和bio5。在2090s-SSP126情景下,現代適生區內的氣候異常區主要出現在中北部,最不相似變量為bio5、bio13和bio15。在2050s-SSP245情景下,現代適生區內的氣候異常區主要出現在中北部,最不相似變量為bio2、bio5和bio15。在2090s-SSP245情景下,現代適生區內的氣候異常區主要出現在中北部和南部,最不相似變量為bio5和bio2。在2050s-SSP585情景下,現代適生區內的氣候異常區主要出現在中部、南部和北部,最不相似變量為bio5和bio15。在2090s-SSP585情景下,現代適生區內除西部和東部少數區域外,其他均為氣候異常區,最不相似變量為bio5(圖5)。

圖4 MaxEnt模型預測的不同排放情景下閩楠多元相似面(Mess)Fig.4 Multiple similarity surfaces (Mess) of P. bournei under different emission scenarios predicted by MaxEnt model

圖5 MaxEnt模型預測的不同排放情景下閩楠最不相似變量(Mod)Fig.5 The most dissimilar variable (Mod) predicted by MaxEnt model in different emission scenarios of P. bournei

3 討論

3.1 模擬結果的可靠性

生態位模型模擬結果的可靠性主要取決于三個方面：首先是模型的選擇,不同模型預測精度差異較大,經驗證,MaxEnt模型預測可靠性較其他模型更高[26- 28]；其次是樣本采集范圍及樣本覆蓋均勻度,采樣點的選擇越接近物種真實分布情況,預測精度越高；第三是環境因子的來源、種類及其精度,模擬時選取因子種類和數量越多,預測精度就越高。調用Enmeval數據包對MaxEnt進行優化可降低過擬合程度和抽樣偏差、從而提高預測精度[29]。本研究經優化后,AUC =0.9846±0.0037,大于0.9,表明預測精度很高。此外,本研究預測的現代閩楠潛在適生區范圍與葛永金得到的結論相近[23],因此模擬結果是可靠的。

3.2 制約閩楠分布的主要環境因子

本研究的8種主要評估途徑中,降水因子重要性列入第一位6次、第二位3次,溫度因子4次列入第一位、第二位6次,因此,降水和溫度共同制約著閩楠現代分布格局。雖然單從重要值來看,分布格局主導因子為最干月降雨量(bio14),但閩楠適宜的溫暖指數為141.7—197.1℃·月[30],且本研究表明溫度因子的置換重要值(Permutation importance)累計占比達86.7%,說明溫度也是制約閩楠分布的重要因子。本研究在進行潛在適生面積預測時,雖然未將城鎮建設和工農業用地排除,但這并不影響本次模擬的意義,對發掘潛在閩楠天然分布區和指導閩楠人工林營建均具有一定的參考價值。相關研究證明,MaxEnt模擬預測存在潛在適生區偏大的現象[7- 8],導致這一現象是由于多種因素綜合作用的結果,物種的空間分布不僅與降水和溫度有關,還受紫外線、生物、人為活動和物種擴散等因素制約,未來的研究中,增加上述變量可能會提高預測的準確性。

貢獻率、置換重要值、單獨使用該因子的正則化訓練增益和單獨使用該因子的測試增益均表明,最干月降水量是影響閩楠潛在地理分布的首要因素,這與葛永金的研究結論相一致[31]。吳顯坤等[32]研究發現影響浙江楠(Phoebechekiangensis)地理分布的首要因素是最干季降水量,而閩楠和浙江楠同為楠屬瀕危植物,親緣關系較近,而且最干月降水量月與最干季降水量兩因子相關性很高。因此本研究得出最干月降水量(bio14)是制約閩楠分布格局的關鍵因子的結論是可靠合理的。此外,本研究發現最冷月最低溫、最冷季降水量、年平均氣溫、降水季節性和晝夜溫差月均值也是制約閩楠分布格局的主要環境因子。

氣候因子并非相互獨立地影響閩楠的生長與分布。多元環境相似性面和最不相似變量分析表明,6個情景閩楠分布區北側氣候異常程度均比其他地區強烈,北側最不相似變量主要是最暖月最高氣溫和最濕季平均氣溫,除2050s-SSP585情景外,與此對應的閩楠未來其他5個情景新增潛在適生區均與喪失區面積無顯著性差異,溫度因子(主要指最冷月最低溫和最熱月最高溫)是制約閩楠適生區南擴的重要因子。然而,閩楠分布變遷是由多個因子共同作用的結果,這一現象尤其在2050s-SSP245和2090s-SSP126情景下表現尤為突出,溫度和降水因子共同制約其地理分布的變遷。

3.3 閩楠潛在適生區的變化格局

在未來6種溫室氣體排放情景下,MaxEnt預測的閩楠潛在分布區都將北移, 其中在2050s-SSP585情景下,其潛在適生區面積擴大較明顯,這可能是由于閩楠為喜暖濕樹種,其現代適生區北部毗鄰區氣溫和降水的增加,可能為其提供了更適宜的生長環境。相關研究也表明楠木[33](Phoebezhennan)在2050s的潛在適生區面積比現代更大；桫欏[34](Alsophilaspinulosa)潛在適生區在2090s將逐漸向西北遷移,且4個情景面積也比現代更大。然而,預測植物潛在分布還需結合土壤、地形、植被、人為干擾和種間競爭等因素,綜合生物與非生物因子,才能更有效地揭示閩楠的未來潛在適生區。未來各情景下,閩楠生境破碎化有加劇的趨勢,特別是2090s-SSP126情景下,閩楠潛在適生區喪失面積占現代適生區面積的40.32%,適生區明顯北移,生境破碎化最嚴重,說明氣候變暖對閩楠分布構成較大威脅。本研究表明,6種未來溫室氣體排放情景下,現代閩楠南端潛在適生區均出現不同程度喪失的現象,且閩楠適生區均往北遷移,黔桂兩省(區)閩楠喪失區比較大,而黔桂又是閩楠古樹分布的集中區,因此,開展閩楠的抗逆性研究已刻不容緩。此外,本研究還表明,雖然閩楠潛在喪失區與新增適生區無顯著差異(除2050s-SSP585情景外),但閩楠未來生境將變的更加破碎化,這一結果與浙江楠、丁香葉忍冬(Loniceraoblata)和藤棗(Eleutharrhenamacrocarpa)等瀕危植物的研究結論相一致[32,35- 36],說明全球氣候變暖將加劇閩楠資源的喪失,亟待從基因組和蛋白組等分子層面提出有效的閩楠種質資源保護方案,以減緩日益變暖的氣候狀況對閩楠生境的影響。此外,種源傳播在不受阻礙的情況下,閩楠有可能向重慶東南部、湖南、湖北、臺灣中部及安徽南部等北部潛在邊緣分布區擴展,因此潛在適生區可將閩楠作為改造針葉人工純林重要樹種之一。

3.4 閩楠的資源保護

前人研究發現[31,37],閩楠致瀕的可能原因有：首先,種子(幼苗)在高濕度生境中受到土壤病原菌的侵染而易變質,而且容易遭受動物捕食,加之閩楠種子休眠特性,從而導致較低的野外閩楠種子發芽率和幼苗存活率低；其次,閩楠是中等耐陰樹種,在高郁閉度的天然林環境下,可能成為制約其幼苗生長的主要因素,導致其幼苗的競爭力低于群落其他耐陰樹種幼苗；再次,由于閩楠經濟價值高,歷史上長期遭受大規模濫砍盜伐。野外調查發現,大多數閩楠古樹并未分布在已建立的自然保護區內,而是以風水林的形式破碎化地分布在村鎮和寺廟旁,加上閩楠種子主要靠重力散布,導致基因流較低[33-34,37],而本研究模擬結果表明,隨著氣候變暖,未來閩楠分布破碎化現象越來越嚴重,閩楠資源保護面臨嚴峻挑戰。綜合閩楠現代適生區圖,未來6個溫室氣體排放情景下適生區圖和空間變換格局,發現閩楠潛在適生區的穩定區為閩中北部、桂東北部、贛中南部、湘贛交界、湘南和黔東等山地溝谷闊葉林,尤其閩中北部、湘南和桂東北山區適生度最高,建議將上述地區列入優先保護區域。劉丹[19]利用EST-SSR和g-SSR兩種分子標記,均得出閩楠具有較高的遺傳多樣性,然而其天然分布破碎化程度較高,說明造成其瀕危的原因可能是由于濫砍濫伐、生境破壞和其特殊的自身繁育特性導致的[18, 38- 39]。因此,建議應盡快對閩楠天然種群采取遷地保護和就地保護相結合的立體保護模式,如閩楠盜伐現象應當引起重視[40-41],同時應建立閩楠保護小區,并定期巡查。本研究團隊已在福建南平來舟林場和順昌洋口林場建立閩楠種子園及種質資源庫,下一步將加快貴州、廣西和湖南等地閩楠種質資源的收集與保存工作,以期為更好地保護和利用閩楠種質資源奠定堅實基礎。