基于MaxEnt模型的小黃花茶在中國的適生區域研究

代玉烜, 金 瑭, 徐海霞, 王 東, 王 麗

(四川大學生命科學學院， 成都 610065)

1 引 言

小黃花茶(CamelliaLuteofloraY. K. Li)，山茶科(Theaceae)山茶屬(Camellia)灌木或小喬木，是山茶屬小黃花茶組的唯一種.1981年首次在貴州赤水發現小黃花茶[1]，據統計，當地僅有1 700余株[2-3].2013年長寧環保局科考人員在四川長寧自然保護區發現了至少幾十株小黃花茶和至少兩處小黃花茶的生長地[4].2015年陳鋒課題組在四川古藺縣桂花鄉發現了小黃花茶的存在，該物種在四川的新發現為進一步的研究提供了基礎資料[5].經過研究者野外踏查及長寧林竹局涂為熙等科考人員統計，目前發現小黃花茶數量約為瀘州古藺1 200株、宜賓長寧800株.

1983年國家科委明令將小黃花茶作為我國禁止外流的特殊保護物種[2]；1988年小黃花茶被列為貴州省一級珍稀瀕危保護植物[3]；2013年在《中國生物多樣性紅色名錄一高等植物卷》中，其瀕危等級評定為易危；IUCN(名)收錄時，瀕危等級為VU，為我國特有種(環境保護部和中國科學院，2013)[6].小黃花茶在形態、結構上的特殊性決定了其在遺傳育種、系統發育、種質資源等方向的科研價值，以及在庭園綠化、園藝觀賞、生態建設等方面的經濟價值.但小黃花茶開花不多，結果較少，種子易受蟲害鳥食、自然生長較慢、移栽成活率低，加之近年來隨著人為砍伐和過度采食加劇了其群落的破碎化等，導致了小黃花茶目前處在瀕危之境.

在如今的氣候變化背景下，對珍稀瀕危植物的生物多樣性保護最重要的在于進行就地保護，利用原生生態環境使被保護的生物能夠更好地生存，以保證生物遺傳繁育的完整性.近年來物種分布模型在生物多樣性保護的多方面得到應用，如物種環境生態位的量化、未來環境氣候變化對物種分布的影響、自然保護區的選擇、傳染病的空間傳播、生物質能源開發、入侵物種擴散及其適生區的問題等研究[7-10].其中MaxEnt生態位模型只需要基于物種存在數據和整個研究區域的環境信息，則能夠分析統計后預測出物種的潛在生境[11-12]，樣本量少、精確度高、操作簡明是其優勢.而ArcGIS10.2不僅能簡便地進行參數地輸入，還能形象化地輸出結果，顯示規律[13].

小黃花茶被發現以來，引起了國內外科研人員和國際組織的高度重視.國內眾多學者對其進行了繁殖技術、種群結構、分布格局、生理生化特性等多方面的研究[14]，但對其適生區的預測研究、氣候變化的影響研究等尚未有報道.本文基于貴州省、四川省小黃花茶的野生生境調查研究、標本記錄和各環境變量，運用MaxEnt生態位模型和ArcGIS相結合預測其潛在的適生區，對該樹種的繁育保護有實際指導意義.

2 材料與方法

2.1 小黃花茶分布點數據

小黃花茶的分布點數據來源于資料搜集和野外調查.(一)通過查閱國內外公開發表的相關文獻、中國植物圖像庫(PPBC， http://ppbc.iplant.cn/)、中國國家標本資源平臺(NSII, http://nsii.org.cn/2017/home.php)等找出已報道的分布地點；(二)查閱有關報道以及政府網站的新聞消息；(三)同時利用野外調查數據補充分布點.目前小黃花茶已知分布點共計20個，通過Google Earth校正，識別并刪除重復點，從而保持1 km×1 km范圍內的一個點，以避免采樣偏差[15].將有效點按照MaxEnt軟件要求的格式錄入到Excel表中，添加為XY數據，保存為*.csv格式.

2.2 方 法

2.2.1 環境數據來源及預處理 (一)從WorldClim網站(http://www.worldclim.org/)下載30 s(約1 km)空間分辨率的生物氣候因子19個(Bio1～Bio19).(二)坡度和坡向、海拔數據從地理空間數據云(http://www.gscloud.cn/)下載的數字高程模型(DEM)通過ArcGIS軟件利用空間分析工具獲得.(三)將從國家基礎地理信息系統(http://nfgis.nsdi.gov.cn/)中下載的1∶400萬中國矢量地圖作為分析底圖，并進行相對應的環境圖層的剪裁提取[16].

考慮到物種分布模型的一個問題是由于變量之間的多重共線性而導致模型過度擬合[17]，會影響對于響應關系、貢獻率的評判[18].故在SPSS22.0軟件中先進行標準化處理，然后進行斯皮爾曼相關性分析，比較所使用的22種環境因子之間的相關性，選擇相關系數<0.75的變量.而生態極限往往和植物的實際分布相關.小黃花茶貴州分布地屬于中亞熱帶溫和濕潤氣候，海拔在700～1 060 m高處，其植株常生在坡度較陡的崖壁上、灌叢中及溝谷兩旁的密林內，多分布于雜木林和竹林中，喜溫暖潮濕[5].結合植物生存環境需求及其生理特性，在相關系數>0.75的變量中保留最重要的生態因子[17].篩選環境變量后通過ArcGIS將獲得的氣候、地形數據.tif格式，轉換為MaxEnt模型運行所需的*.asc格式[11]，建立面文件對環境數據進行掩膜處理，使模擬結果的實用性更高.

2.2.2 模型預測及數據處理 將小黃花茶的有效分布點數據和主要環境變量數據集同時導入MaxEnt模型中，將測試子集設置為分布點數據的25%，訓練子集為數據的75%.以線性二次抽樣、鉸鏈功能的方式進行10次重復運算，其余參數設置為默認值，算法在160次迭代(1 s)后收斂，運行計算物種的適生區.

采用刀切法(Jackknife)檢驗環境因子貢獻權重，設置受試者特征曲線(ROC，Receiver operating characteristic curve)，以假陽性率即1-特異性為橫坐標，代表傭金誤差；以真陽性率即靈敏度為縱坐標，表示沒有遺漏誤差.用篩選的所有主要變量建立模型(紅色條帶表示)；單獨用每一個變量建立模型(寶藍色條帶表示)；每次排除依次排除一個變量，用剩余的變量建立模型(青藍色條帶表示)，觀察條帶之間的變化.當藍色條帶越長代表增益值越高，說明該變量對物種分布貢獻越大；當青藍色條帶降低較多，則代表該變量可能含有較多信息，對物種分布較為重要；紅色條帶表示所有環境變量對建模的累積貢獻率[19].ROC曲線是通過繪制所有可能閾值的x軸上的1-特異性和y軸上的靈敏度而成的[12].其曲線下面積(AUC，Area under the ROC curve)的大小作為評價模型運算結果精準度的指標[20-21].AUC值在0.5～1，一般來說, 0.5≤AUC<0.6為差，0.6≤AUC<0.7為較差, 0.7≤AUC<0.8為一般，0.8≤AUC<0.9為好，0.9≤AUC<1.0為非常好[18,22].

運行ArcGIS軟件處理預測結果，將適生分布信息按自然簡單分級法顯示[23].將分析結果導入軟件中，依據Maxent模型最大約登指數及自然間斷點對其進行重分類，得到小黃花茶生境適宜性評價圖，在布局視圖中點擊插入對應的比例尺、圖例、指北針，完善分布圖，最后運用ArcMap軟件空間分析模塊工具提取各潛在適生區的面積[24].

3 結果與分析

3.1 MaxEnt模型精度評估和變量分析

經過校正，本研究最終運用于生態位模型運算的小黃花茶有效分布點共15個(圖1)，數據導入后形成表1.圖2顯示了22種環境因子之間的相關性，環境數據處理最終結果僅保留了8個因子(表2)，其中降水量因子占據3個，氣溫因子占據2個，這表明與氣候變化相關的氣候不穩定性和干旱壓力的增加將是當地物種損失的主要決定因素.在模型預測中，使用14個分布點記錄用于訓練，1個記錄用于測試.預測結果如圖3a，顯示小黃花茶適生區分布最大熵模型運算的10組訓練集(cha)AUC均值為0.947，測試集的均值為0.953，表示了實驗數據可信度高.圖3b顯示了本研究運算結果中ROC預測精度數值最大的1組(cha_7)，訓練子集AUC、預測子集AUC分別為0.947和0.997，非常接近1，根據Swets提出的模型評價標準[25-26]，表明該模型對其適生區的預測精度高，為優級. cha_7可作為分析對象.

刀切法(Do Jackknife)檢驗結果估算了環境變量對MaxEnt模型的相對貢獻(表2)，最熱季平均氣溫(Bio10)貢獻率最高(51.1%)，其次是最熱季降水量(Bio18)(23.6%)、最冷季降水量(Bio19)(16.1%)；而坡度、最濕月降水量(Bio13)、坡向、海拔(Dem)、最熱月最高氣溫(Bio5)的貢獻率分別為4.7%、1.9%、1.1%、0.9%、0.3%，累計僅達8.9%，表明了它們本身對小黃花茶分布沒有占據很大的影響.上述環境因子中降水因子占有3個，地形因子占3個，溫度因子占1個，這與已知研究的該樹種實際分布情況(喜分布于地勢較低、地形閉塞、溫濕度較高的特殊生境[27-29])是一致的.這8個環境因子在模型構建過程中起主要作用，也是影響小黃花茶潛在分布最重要的生態因子.

圖4是可變重要性的折刀測試的結果，顯示了8個環境變量對MaxEnt模型預測結果的影響.當單獨使用變量時，重要性最大的環境變量為最熱季降水量(Bio18)，訓練增益超過1.0；其次是最熱季平均氣溫(Bio10)、最冷季降水量(Bio19)以及最濕月降水量(Bio13)對模型影響較高，訓練增益超過0.6.結果顯示，環境溫濕度是影響小黃花茶的潛在適生區分布的重要因素，而夏季最高溫及坡度是小黃花茶發育發展的主要制約因子.這與目前所發現的小黃花茶生長環境相符合.

圖1 小黃花茶樣點分布圖Fig.1 Distribution points of Camellia Luteoflora

圖2 22個環境因子間的Spearman相關分析Fig.2 Spearman Correlation results among 22 environmental factors

表1 小黃花茶15個分布點的地理數據

表2 8個環境變量對MaxEnt模型的貢獻率和訓練集增益

圖3 模型檢驗AUC值

圖4 基于AUC的Jackknife法環境變量重要性測試圖

3.2 小黃花茶潛在適生區分析

借助ArcGIS處理，小黃花茶的適生區預測圖可按4個等級進行劃分，分級標準如下：<0.08為非適生區，0.08～0.27為低適生區，0.27～0.54為中適生區，>0.54為高適生區.從圖5可知，小黃花茶在當前氣候條件下預測的潛在適生區與現已知的實際生長地大致相符，可能存在尚未發現的生長區.圖中淺米色區域表示非適生區，綠色區域表示低適生區，黃色區域表示較適生區，紅色區域代表高適生區.小黃花茶適生區主要處于中國西南方103°E～109°E，25°N～30°N的地區，四川瀘州、宜賓，貴州遵義、仁懷、赤水，云南昭通、昆明以及重慶江津、綦江是其高度適生區，黔、川、渝交界處是小黃花茶的核心分布地.除此之外，預測顯示在貴州畢節、銅仁、黔東南苗族侗族自治州，四川自貢、廣安，云南曲靖、玉溪，重慶西部及北部有零星低適生區.

提取適生區面積結果顯示(表3)，小黃花茶在中國的總適生區面積約為21 371 km2，占我國陸地總面積的0.222%.其中高、中、低適生區的面積占比分別為0.028%、0.068%和0.126%，而非適生區面積占比高達99.778%.

表3 小黃花茶潛在適生區面積統計

3.3 適生地理分布與變量的關系

模型用響應曲線圖(圖6)展示了每個環境變量是如何影響MaxEnt預測的，該圖橫坐標為對應環境因子的數值，縱坐標即為分布概率.曲線反映了預測的分布適宜性對單一變量的依賴性，曲線變化可以為對應物種提供最佳生長所需環境閾值的有用信息[30-31].一般認為，當分布概率>0.50時，對應的環境因子值適合物種的生長.

當最熱月最高氣溫(Bio5)(圖6a)處于<4.5 ℃時，小黃花茶分布概率極低；隨著最高氣溫的上升，分布概率逐漸增大；達到27.5 ℃時，小黃花茶分布概率到達頂峰.即在最熱月時最高氣溫達到27.5 ℃時，小黃花茶的最適生長溫度到達.而最熱季平均氣溫(Bio10)(圖6b)的閾值范圍為11.1～18.0 ℃，最適閾值低于14.5 ℃，當氣溫高于25 ℃植物分布概率趨近于0.

最濕月降水量(Bio13)在290～360 mm時，小黃花茶的分布概率到達峰值(圖6c).以分布概率>0.50為適宜范圍，最濕月降水量為290～500 mm，此降雨區間小黃花茶生長較好.當最濕月降水量>700 mm后，小黃花茶分布概率幾乎為零.濕度過大可能更容易造成植物病蟲害，反而不利于小黃花茶的生長.最熱季降水量(Bio18)(圖6d)的閾值為4.4～6.6 mm，8.5 mm以上的最熱季降水量將會影響植物的生長.最冷季降水量(Bio19)(圖6e)適合目標物種生長的閾值為16.7～24.8 mm.

小黃花茶在海拔450～1 000 m的區域中分布概率較高(圖6f)，在此區間之后呈現出海拔越高分布概率越低的趨勢，且小黃花茶不分布于海拔>4 000 m的地方.海拔愈高，溫度愈低、濕度增加、光照增強、土壤含量變化，對小黃花茶的生長和分布產生了負影響.在坡向為-40°～190°范圍內小黃花茶的分布概率>0.50，說明小黃花茶喜愛分布于陰坡和半陰坡以及小部分的陽坡(圖6g).當坡度>20°時，小黃花茶具有較高的分布概率(>0.50)，隨著坡度的增加其分布概率也在增大，在坡度達到64°時，小黃花茶的分布概率到其極限值(圖6h)，可能是過于陡峭的斜坡土壤含水量流失較快，小黃花茶無法汲取充足水分生長，其選擇分布在更適宜生長的生態帶.

圖5 基于MaxEnt模型小黃花茶適生分布預測

4 討 論

4.1 小黃花茶潛在適生區分布

小黃花茶目前的分布面積較為狹小，其植株生長緩慢、病蟲害嚴重、種子自然萌發率低，加之人為干擾迫害加劇，小黃花茶數量呈現下降趨勢.小黃花茶在1983年被我國列為珍稀瀕危保護物種，世界自然基金會(WWF)在2003—2004年將貴州赤水小黃花茶的保護和發展納入基金會在中國的物種保護資助項目內，2013年小黃花茶瀕危等級評定為易危，目前國家和國際組織對小黃花茶的保護高度重視[7,32].物種適生區分布的精準預測對于物種的推廣和引種極其重要[33]，在明確其目前已知的野外生存現狀下，對其潛在的適生分布區進行預測十分重要.

在本次模擬中，雖然我們能夠測試的物種數量分布點少于20個,但一些研究同樣使用低樣本量成功模擬了物種潛在分布，證明我們在樣本量小的情況下使用Maxent是可行的.例如Deb等人對一種稀有蘭花(VandabicolorGriff.)僅使用了4個發生點來開發小樣本分布預測模型[15]；Ortega-Huerta等對樣本量只有2～20的11種鳥類和其地區環境條件的數據運用不同的模型預測了物種分布[34]；Pearson等對馬達加斯加特有的樣本大小為4～23個的13種壁虎(Uroplatusspp.)建立模型，并對其進行評估[35]；Adhikari等人對印度東北部一種極度瀕危的大葉冬青樹種進行再引種生境分布模擬，其樣本量為16個[36]；Groff等利用17個已驗證的位置開發了生態位模型，以評估俄勒岡州斑點蛙(Ranapretiosa)在加利福尼亞州的存在情況[37]；Hernandez等人的研究表明當樣本量小到5、10和25個時，四種建模方法(Bioclim、Domain、GARP和Maxent)中，Maxent是最有能力產生有用結果的[38].本次MaxEnt模型預測結果的高精確度證明了在樣本量小的情況下預測小黃花茶分布是可信的.

適生區預測結果顯示，小黃花茶主要適生于黔、川、渝及滇部分區域，這些區域主要位于長江流域一帶，氣候、土壤條件適宜，適合喜濕喜蔭的小黃花茶生長.這與目前所報道的小黃花茶發現地四川瀘州古藺、宜賓長寧竹海自然保護區、貴州赤水自然保護區相吻合，根據試驗者實地調查亦得知，瀘州敘永、貴州貴陽等地均有小黃花茶生長.

4.2 小黃花茶生態因子特征

刀切圖顯示了影響小黃花茶分布的環境因子中貢獻率前三分別為Bio10、Bio18、Bio19，重要性最大的環境變量為Bio18，其次是Bio10、Bio19、Bio13對模型影響較高.可以看出，最熱季的平均氣溫和降水量以及最冷季降水量是影響小黃花茶生長分布最重要的環境因子，這與野外調查的小黃花茶多生長于雜林中、溝谷內、崖邊、溪水旁的生長地環境相符合，其喜愛溫暖濕潤，充沛的降水和濕潤的環境利于其生長.刀切法的分析結果認為平均溫度和降水是影響小黃花茶分布的關鍵要素，最熱季平均溫度與小黃花茶分布密切相關，降水量過低或過高均不利于其分布，說明降水和溫度共同決定了其分布區氣候類型.

主要環境因子響應圖顯示了小黃花茶的生態適宜帶為最熱月最高氣溫為22.5～27.5 ℃，最熱季平均氣溫為11.1～18.0 ℃，最濕月降水量為290～500 mm，最熱季降水量為4.4～6.6 mm，最冷季降水量為16.7～24.8 mm，海拔為450～1 000 m，坡向在-40°～190°范圍內，20°<坡度<64°.植物對生態環境的選擇策略不同，并對其有一定的耐受度，不在度區間時植物的生存則會受到影響.因此在明確了小黃花茶分布的限制環境因子后，可以對其展開生態修復.

本研究基于MaxEnt模型對小黃花茶的潛在適生區進行了預測分析.MaxEnt模型對于物種適生區分布預測的準確度，主要是基于統計或理論推導的響應面，將野外觀測與環境預測變量聯系起來的經驗模型，其由提供給模型的物種實際分布數據和環境變量數據決定[32,39].

在實際分布數據方面，因小黃花茶分布的地域性以及環境特殊性，其樣本分布點數據量很少，預測結果可能存在一定誤差.根據前人的研究，當數據集過小時可以選擇兩種方法降低對預測結果的影響：(1)不劃分訓練和測試子集時，因為每個位置都可能提供唯一的、有價值的信息，即先從數據集中刪除每個局部性一次，并用剩余的局部性生成一個模型，然后在設定為10%的固定閾值下進行評估，最后使用P-Value Compute程序以測試模型預測是否優于排除地點的隨機分配[35]；(2)劃分訓練和測試子集時，對具有10%訓練存在的閾值規則的物種進行數次重復的模型運行，采用交叉驗證技術將樣本劃分為重復折疊，并將所有其他參數用作測試數據,構造出一個可行的模型[15,36].本研究采用了劃分訓練子集和測試子集的方式，尚未在設定為10%的固定閾值下進行評估，可能存在預測范圍有誤差.但是小樣本量開發的模型應該被解釋為識別與已知物種出現的環境條件相似的區域，而不是預測物種范圍的實際限制，因此本研究尚具有一定的區域參考價值.當然在后續的研究中，我們將進行多方面的建模，以獲得更高、更顯著的成功率.在環境變量方面，本文僅僅考慮了氣候及地形因子，而沒有考慮土壤類型、水分、植被、植物天敵等環境因素以及由自然或人為引起的外界擾動對其分布的影響，存在一定的局限性和模型預測結果的不確定性.如能綜合所有因素，預測結果將會更加準確.同時，在環境變量的篩選時，可以采用MaxEnt模型將生態環境因子數據多次計算，每次舍去貢獻率為0的生態因子，最終聯系刀切法選擇出貢獻率在5%以上的環境因子作為最終環境變量組來減小對植物生境適宜性預測的影響[40].

在研究區域方面，本研究以四川省、湖北省、重慶市、湖南省、云南省、貴州省、廣西壯族自治區等七省(自治區、直轄市)為研究對象，空間尺度適中、精度有待提高，并且僅在小黃花茶已知的現存分布區域內進行預測與評估，可以考慮在更小范圍區域展開更加精細的研究.從預測區域大小來看，在較小的區域內Maxent模型能夠給出顯著的預測結果, 而當小樣本數據用于更大的區域的預測時，預測模型即變得不穩定[36]，因此應該在較小面積上使用小樣本量來開發模型以降低對預測結果的影響.在下一步的研究中或可以收集更多與植物分布相關的自然、人為因素添加到模擬模型中，以提高預測精度，為小黃花茶的資源保護與開發利用提供更有效的參考信息，并為其他繁育研究、工作等提供科學依據.同時我們也將采集更多的物種分布信息以提高預測精確度.

在未來的研究中，我們也會預測其對于未來環境氣候條件變化的響應，以此更好地了解植物對生長環境因子的需求以及其對環境變化做出的應對機制[41].同時考慮從生態位角度探討珍稀瀕危植物的瀕危機制、種間競爭及群落共存特征，為針對性地制定珍稀瀕危物種保護措施提供有利信息[42].同時我們將開展地面調查工作，通過地面試驗評估模型的預測能力.

MaxEnt模型根據氣候、地形因子與小黃花茶之間的關系，找出其分布規律的最大熵，評估了影響植物生長分布的主導環境因子，預測了其潛在適生區，這對于指導植物的野生資源保護、繁育引種、開發利用等具有重要意義.小黃花茶作為中國特有植物，在我國分布范圍也相當狹窄，而了解其適生區范圍、適生環境閾值對發揮其高經濟價值、高科研價值具有重要意義.