基于MaxEnt模型的內蒙古自治區(qū)樟子松潛在分布研究

宋慧芳，黨曉宏，2*，高 永，蒙仲舉，孫艷麗

（1. 內蒙古農業(yè)大學沙漠治理學院，呼和浩特 010000；2. 內蒙古杭錦荒漠生態(tài)系統國家定位觀測研究站，內蒙古 鄂爾多斯 017000）

樟子松（Pinus sylvestris）生長慢但壽命長、抗逆性和耐寒性強，具有耐腐蝕、穩(wěn)定性強等特性[1-2]。1995 年科爾沁沙地首次成功引種樟子松[3]，并在三北防護林建設中得到廣泛應用，是我國北方水土保持和防風固沙工程中最重要的樹種之一[4]。但是由于環(huán)境和植物本身影響因素的干擾，北方樟子松人工林開始大面積退化[5]，特別是科爾沁沙地將近40%的樟子松防風固沙林出現生長不良癥狀[6]。目前樟子松的研究主要集中在培育技術[7]、土壤理化性質[8]和物種多樣性[9]等方面，對其潛在適生分布區(qū)的研究較少。

目前常見的用于物種潛在適生區(qū)預測模型有最大熵（MaxEnt）、生態(tài)位因子分析（ENFA）、遺傳算法（GARP）、生物氣候分析系統（BIOCLIM）和區(qū)域環(huán)境模型（DOMAIN）等[10]。對比上述這些模型，其中適用范圍較廣的是最大熵（MaxEnt）模型，即使物種樣本數據數量較少的情況下也能有較高的預測精度和準確的預測成果[11]。MaxEnt 模型根據物種已有地理分布位置和相應的多種環(huán)境因子，按照氣候相似原理，計算出特定生態(tài)位限制條件下物種分配定律的最佳狀態(tài)，即熵最大時該物種在研究區(qū)內的潛在分布范圍[12]。近年來，普遍運用在外來物種入侵預警[13-14]、疾病傳播途徑[15-16]及物種潛在地理分布預測[17-18]等領域。Phillips團隊運用Java語言編寫了MaxEnt模型軟件[19]，劉超等[20]基于MaxEnt模型結合ArcGIS 軟件對川梨（Pyrus pashia）自冰期以來不同時期的分布進行預測，認為隨歷史氣候的變遷，川梨由末次盛冰期經全新世中期進入當代，分布區(qū)面積呈現先下降后上升的趨勢。王浩等[21]分析野生防風（Saposhnikovia divaricata）分布信息及多個環(huán)境因子，預測防風生態(tài)適宜性區(qū)劃，為栽培防風的擴大生產提供基礎研究數據。秦委等[22]基于MaxEnt模型結合ArcGIS 軟件研究我國范圍內東南茜草（Rubia argyi）的潛在地理分布，預測結果可為東南茜草野生撫育與人工栽培的合理布局提供理論參考依據。

因此，為了探索樟子松潛在分布區(qū)域，并進行科學推廣種植，為樟子松的合理分區(qū)栽種提供基礎研究數據，本文采用實地調查與數據處理相結合的方法，以內蒙古自治區(qū)為例，將樟子松存在記錄分布點數據及其相關環(huán)境因子，如氣候、土壤和地形因子數據等利用ArcGIS軟件進行前期數據處理，利用MaxEnt 模型統一構建其潛在生境的空間適宜分布函數，最終得到樟子松潛在分布區(qū)域以及影響其分布的主要環(huán)境因子，對實現樟子松科學種植以及防風固沙林的穩(wěn)定發(fā)展具有重要意義。

1 材料和方法

1.1 研究區(qū)概況

內蒙古自治區(qū)(37.40°～53.38°N，97.27°～126.07°E)位于我國西北部，東西直線距離2 400 km，土地總面積118.3 萬km2，占中國土地總面積的12.3%。全區(qū)地勢較高，平均海拔高度在1 000 m以上，屬于蒙古高原型地貌，全年大風日數平均在10～40 d，年日照時數都在2 700 h 以上，屬于溫帶大陸性季風氣候。年降水量在100～500 mm 之間，受到氣候及地形地貌影響，內蒙古降水量在空間分布上呈現出東北高西南低的趨勢[23]，是我國典型的農牧帶過渡區(qū)，東部地區(qū)土壤肥沃，土地生產力較強適宜農業(yè)發(fā)展，中西部地區(qū)土壤水分結構相對較差，存在大量荒地和沙漠，生態(tài)系統以草原、森林、荒漠和裸地為主[24]。

1.2 數據來源

樟子松分布點數據來源于國家標本資源共享平臺（http://www.nsii.org.cn/）、中國數字植物標本館（http://www.cvh.ac.cn/）、全球生物多樣性信息服務網絡平臺（https://www.gbif.org）和相關文獻資料[25]和實地調查，年份范圍為1973—2020 年；氣候因子來源于世界氣候數據庫（http://www.worldclim.org），選擇1970—2000年生物氣候（Bioclimatic）變量圖層數據；地形因子來源于地理空間數據云平臺（http://www.gscloud.cn）；土壤因子來源于北京大學城市與環(huán)境學院地理數據平臺（https://geodata.pku.edu.cn/），年份范圍為1971—1981年。

1.3 數據處理

研究區(qū)內共獲取73條樟子松分布信息，對于有地名無坐標信息的分布點通過谷歌地圖補全其地理坐標，將以上坐標信息導入ArcGIS剔除超出研究區(qū)范圍、錯誤和重復的分布點，最終獲得45 條樟子松有效分布點，如圖1 所示，使用Excel 將點位信息統一為spiecis、latitude、longitude 3 列表并以.csv 格式保存；由于樟子松為深根性樹種，所以最終選擇共17條以S_開頭的下層土壤屬性（30～100 cm）土壤類型包括：土壤有效含水量、沙含量、黏土含量、土壤質地分類、土壤容重、有機碳含量、酸堿度、土壤的陽離子交換能力、基本飽和度、硫酸鹽含量和電導率等；結合ArcGIS處理內蒙古自治區(qū)海拔、坡度、坡向數據、19項氣候數據和17項土壤數據，共39項環(huán)境因子轉換為ASCⅡ格式用于MaxEnt 模型對內蒙古地區(qū)樟子松生境的模擬。

圖1 內蒙古地區(qū)樟子松存在記錄分布點位置圖Figure 1 Location map of P. sylvestris recorded in Inner Mongolia

1.4 MaxEnt模型參數設置

環(huán)境因子之間高度相關會使模型預測準確性降低，故用MaxEnt模型中的刀切法（Jackknife）來初步測定環(huán)境因子模擬貢獻率大小，為避免過度擬合要選擇保留貢獻率大的因子（貢獻率≥2%），在Arc-GIS 中利用多元分析工具，統計柵格圖層數據相關矩陣，比較相關系數，若|r|＜0.8，全部保留，若|r|≥0.8，對比MaxEnt模型中兩者貢獻率，剔除貢獻率相對較小的影響因子[26]，經過篩選本研究最終確定了6 個主要環(huán)境因子用于模型預測樟子松的潛在分布，如表1所示。

表1 環(huán)境變量對MaxEnt模型構建的貢獻率Table 1 Contribution rate of environmental variables to MaxEnt model construction

結合樟子松分布點數據和6 項環(huán)境數據構建MaxEnt模型，建立模型時隨機選取75%的樟子松存在記錄點作為訓練數據，剩余25%的數據則作為驗證數據[27]，將重復次數設置為10 次，以此來消除訓練數據（75%）和驗證數據（25%）產生的誤差[28]，最終得到預測結果。

MaxEnt 模型物種潛在分布預測運行結果中的物種工作特征曲線（receiver operating characteristic，ROC）能夠模擬預測進行自檢驗，采用刀切法（Jackknife）分別對每一個環(huán)境變量進行測定，且衡量每一項環(huán)境變量的貢獻率，可以分析不同環(huán)境變量對物種分布情況的作用。在ArcGIS 中將MaxEnt 模型運行結果的樟子松潛在適宜分布區(qū)進行重分類，適宜程度劃分為（適生指數）P＜0.05為非適生區(qū)、0.05≤P＜0.33 為中適生區(qū)和P≥0.33 為高適生區(qū)3 個等級[29]，進行深入系統的分析。

2 結果與分析

2.1 MaxEnt模型精度

MaxEnt 模型對物種潛在分布預測結果進行精度評價是利用受試者工作特征曲線來驗證的，ROC曲線以假陽性率為橫坐標，真陽性率為縱坐標，曲線與橫坐標所圍成的面積值為AUC（area under curve，AUC）。AUC值的范圍為0～1，值在0.5～0.6為不及格；0.6～0.7 為較差；0.7～0.8 為一般；0.8～0.9 為較好；0.9～1.0 為優(yōu)秀[30]，測試結果AUC 值越大效果越好，AUC值越大說明模型預測的精度越大。經過MaxEnt模型對樟子松分布信息及環(huán)境因子重復迭代運行10次且取其均值后，ROC曲線研究結果如圖2所示，AUC＞0.8，可以看出MaxEnt 模型預測精度較好，預測結果準確可靠，可以應用在樟子松潛在分布研究中。

圖2 樟子松潛在分布預測的ROC曲線Figure 2 ROC curve for predicting potential distribution of P. sylvestris

2.2 主導環(huán)境因子

從表1、圖3 看出，氣候因子最潮濕月份的降水量（BIO13）貢獻率最大，達到32.1%，最潮濕季節(jié)的平均溫度（BIO8）置換重要性最大，為47.2%，氣候因子累計貢獻率和置換重要性分別為88.7%、73.7%；而土壤因子（class）和地形因子（altitude）累計貢獻率、置換性為11.4%、26.2%，因此影響樟子松潛在分布的主要因子為氣候因子，其次是土壤、地形因子。根據刀切法測定的變量權重模擬時的訓練得分判定標準：“只含某變量”得分較高時，說明該因子重要，對物種分布貢獻較大。正則化訓練增益、測試增益和樟子松工作特征曲線下面積（AUC）值“只含某變量”得分最高前4 位的變量均為最潮濕月份的降水量（BIO13）、等溫性（BIO3）、最潮濕季節(jié)的平均溫度（BIO8）和最干燥季節(jié)的平均溫度（BIO9）。結合上述分析最終確定最潮濕月份的降水量（BIO13）、等溫性（BIO3）、最潮濕季節(jié)的平均溫度（BIO8）、最干燥季節(jié)的平均溫度（BIO9）為影響樟子松潛在分布的主導環(huán)境因子。

圖3 樟子松潛在分布刀切法預測圖Figure 3 Prediction of potential distribution of P.sylvestris by knife cutting method

分析圖4樟子松分布概率—環(huán)境因子的響應曲線表明，當樟子松分布概率大于0.5時，最潮濕月份降水量（BIO13）分布范圍在80～135 mm，其中最適宜范圍在105～115 mm之間；等溫性（BIO3）計算：[平均日較差（BIO2）/溫度年較差（BIO7）]×100，代表一年間溫度變化的穩(wěn)定狀態(tài)，可以看出等溫性值≤26 時為其適宜范圍，最適宜范圍低于23.0；最潮濕季節(jié)平均溫度（BIO8）范圍在12～19.5 ℃時為適宜變化范圍，17 ℃左右時樟子松分布概率出現最高值，為最適宜；最干燥季節(jié)的平均溫度（BIO9）在-30～-15 ℃和0 ℃以上時為適宜變化范圍，其中溫度回暖至0 ℃以上最適宜。當最潮濕季節(jié)的平均溫度過高，最干燥季節(jié)的平均溫度過低都會抑制樟子松的分布。

圖4 樟子松的分布概率—環(huán)境因子的響應曲線Figure 4 Distribution probability of P. Sylvestris -- Response curve of environmental factors

2.3 樟子松潛在分布范圍

結合MaxEnt模型適生指數P值的運行結果，參照IPCC[31]關于評估物種適宜生存“可能性”的劃分標準以及周炳江等[29]對云南榧樹（Torreya yunnanensis）、李昂等[32]對樟子松的潛在分布劃分方法，將樟子松在研究區(qū)域的潛在生境適宜程度劃分為：（適生指數）P＜0.05 為非適生區(qū)、0.05≤P＜0.33 為中適生區(qū)、P≥0.33為高適生區(qū)。根據模型預測結果如圖5、表2 所示，樟子松在內蒙古自治區(qū)適宜分布區(qū)域為呼倫貝爾、興安盟、錫林郭勒盟、通遼、赤峰、呼和浩特、巴彥淖爾和鄂爾多斯等地區(qū)，適宜分布面積約為86.997 8 萬km2，占總面積的72.73%。高適生分布區(qū)主要在呼倫貝爾市北部、興安盟北部、錫林郭勒盟南部、赤峰市北部和呼和浩特市，分布面積為30.811 7 萬km2，占總面積的25.76%，模型預測結果中阿拉善盟地區(qū)樟子松呈現零星破碎化分布，累計非適生區(qū)分布總面積共約32.623 6萬km2，占總面積的27.27%，從環(huán)境條件看阿拉善盟地區(qū)高溫缺水，目前不適宜引種人工樟子松。

表2 內蒙古地區(qū)樟子松潛在分布面積統計Table 2 Statistics of potential distribution area of Pinus Sylvestris in Inner Mongolia 萬km2

圖5 內蒙古地區(qū)樟子松潛在分布區(qū)劃Figure 5 Potential distribution regionalization of P.Sylvestris in Inner Mongolia

3 討論與結論

3.1 討論

目前樟子松林天然分布于我國大興安嶺和呼倫貝爾沙地草原[33]，分布最廣的位于內蒙古自治區(qū)呼倫貝爾市鄂溫克族自治旗境內紅花爾基樟子松林區(qū)。樟子松是東北大興安嶺山地地帶造林和固沙造林的先鋒樹種，能夠涵養(yǎng)水源，十分有利于防風固沙，因此新中國成立以后，樟子松的人工栽植區(qū)日益擴展，自此先后引入到內蒙古等地區(qū)[34]。同樣本研究預測結果表明樟子松主要位于內蒙古呼倫貝爾、興安盟、赤峰等東北部，以及呼和浩特市、包頭市和鄂爾多斯市等，這與張日升等[35]在三北防護林建設地區(qū)樟子松適宜分布區(qū)分析結果相一致。MaxEnt 模型預測結果與樟子松實際地理分布高度一致，表明MaxEnt模型對于樟子松的分布預測結果較為準確、可靠。研究區(qū)內高、中適生區(qū)面積分別為30.811 7 km2、56.186 1 km2，而阿拉善盟、鄂爾多斯西北部以及烏蘭察布北部樟子松分布面積較小。受環(huán)境制約，西部區(qū)域普遍干旱少雨，年均溫較高，降水較少且蒸發(fā)強烈，導致水分不足無霜期較長，地貌類型多為流動沙地、半流動沙地及戈壁，環(huán)境惡劣，導致樟子松分布較少。

本研究結果表明，影響樟子松潛在分布的主要因子是氣候因子，其次是土壤和地形因子。這正與李昂[32]對于北方樟子松適宜性分布研究結果高度一致。在氣候衍生變量中，最潮濕月份的降水量（BIO13）、等溫性（BIO3）、最潮濕季節(jié)的平均溫度（BIO8）和最干燥季節(jié)的平均溫度（BIO9）為影響樟子松分布的主要因子，其中溫度類型的主要氣候因子累計貢獻率高達69.4%，說明樟子松分布對于溫度變化響應更為敏感。而李蒙蒙等[36]研究也同樣表明，溫度是影響天然樟子松引種的主導氣候因子，其次是濕度、光照、海拔和風速。其中最潮濕月份的降水量（BIO13）最適宜范圍在105～115 mm之間，北方地區(qū)降水量月份較為集中，降水過多造成林分高溫高濕的環(huán)境條件林分生長勢降低[37]。最潮濕季節(jié)的平均溫度（BIO8）適宜值在17 ℃左右、過高時導致蒸發(fā)量增加，水分缺失，樟子松分布急劇降低。最干燥季節(jié)的平均溫度（BIO9）0 ℃以上最適宜，其干旱、低溫的環(huán)境條件導致樟子松生長速度減緩[36]。

由于森林生態(tài)系統對影響因子的響應較為復雜，預測結果與實際分布會存在一定的誤差。比如地形因素影響導致一些點位信息未被調查到，因此通過查閱盡可能多的網絡數據庫來減小偏差；今后可以從生物因素、種內和種間關系、人為因素等多方面綜合進行考量[38]；在篩選環(huán)境變量時通過綜合考量方式來降低模型誤差。未來還要將大量的研究數據進行模擬，情景進行分析，使得預測更精確，結果更可靠。

3.2 結論

①在影響樟子松分布的多重環(huán)境變量中，氣候因子對其生長響應更為敏感，其次土壤和地形因素。

②樟子松在內蒙古自治區(qū)適宜分布區(qū)域主要集中在中部和東北部，其中呼倫貝爾、興安盟、錫林郭勒盟、通遼、赤峰、呼和浩特、巴彥淖爾和鄂爾多斯等地分布最廣，阿拉善盟地區(qū)分布較少，其高溫缺水的環(huán)境條件抑制了樟子松的生長。

③按照本研究篩選出的主導環(huán)境因子，可以增加相對應的管護措施，其要注意最潮濕、干燥季節(jié)和最潮濕月份的平均溫度和降水量，可以在冬季來臨之前做埋土處理，在臨冬和臨春季節(jié)及時澆水，減小其蒸發(fā)，保證其生存環(huán)境內適宜的溫度及水分[39]，促進樟子松的生長發(fā)育繁殖。