王桂榮，王彥珺，王尚濤，孫志蓉

北京中醫(yī)藥大學(xué)中藥學(xué)院，北京 102488

當(dāng)歸為傘形科植物當(dāng)歸Angelica sinensis(Oliv.)Diels的干燥根[1]，其用途廣泛[2]，除藥用外，還廣泛應(yīng)用于保健品、食品領(lǐng)域，也常作為日化產(chǎn)品的原料。目前人類活動(dòng)導(dǎo)致大氣中二氧化碳濃度顯著升高，溫室效應(yīng)加劇，極端氣候頻發(fā)[3]，且我國已成為全球氣候變化敏感區(qū)[4]。研究表明，道地藥材產(chǎn)區(qū)的變遷與氣候變化有關(guān)[5]，但盲目引種導(dǎo)致藥材良莠不齊，制約當(dāng)歸產(chǎn)量和質(zhì)量的提升[6]。目前當(dāng)歸區(qū)劃研究主要集中在當(dāng)前潛在適生區(qū)預(yù)測，鮮有對(duì)未來不同氣候情景下適生區(qū)的預(yù)測。本研究利用最大熵（MaxEnt）模型和地理信息系統(tǒng)（GIS）對(duì)未來不同氣候條件下當(dāng)歸適生區(qū)進(jìn)行預(yù)測，分析影響當(dāng)歸分布的主要環(huán)境因子，預(yù)測適生區(qū)遷移方向，為當(dāng)歸生產(chǎn)規(guī)劃提供參考。

1 數(shù)據(jù)收集及處理

1.1 當(dāng)歸分布信息

檢索全球生物多樣性信息平臺(tái)（https://www.gbif.org/）、中國數(shù)字植物標(biāo)本館（http://www.cvh.ac.cn/）、國家標(biāo)本平臺(tái)（http://www.nsii.org.cn），查閱相關(guān)文獻(xiàn)[7-29]，收集當(dāng)歸在我國的地理分布信息。對(duì)于無經(jīng)緯度的樣本，根據(jù)文獻(xiàn)和網(wǎng)站記錄的位置信息，結(jié)合satelite map（https://map.jiqrxx.com/jingweidu/）確定經(jīng)緯度；對(duì)于經(jīng)緯度單位不一致的，結(jié)合相關(guān)網(wǎng)站（https://www.osgeo.cn/app/s2703）進(jìn)行轉(zhuǎn)化。共收集到921條當(dāng)歸地理分布信息，去除經(jīng)緯度重復(fù)和錯(cuò)誤的信息，采用緩沖區(qū)分析法對(duì)獲得的分布數(shù)據(jù)進(jìn)行篩選（緩沖區(qū)半徑2.5 km，當(dāng)分布點(diǎn)距離＜5 km時(shí)保留其中一個(gè)點(diǎn)），最終得到當(dāng)歸樣點(diǎn)信息218 條，并轉(zhuǎn)化成MaxEnt模型使用的.csv格式。本文所用中國行政區(qū)劃圖下載于國家基礎(chǔ)地理信息系統(tǒng)（http://nfgis.nsdi.gov.cn/），分辨率為1∶400萬。

1.2 環(huán)境因子

當(dāng)前（1970－2000 年）和未來（2041－2060 年、2061－2080年、2081－2100年）氣候數(shù)據(jù)來自世界氣候數(shù)據(jù)庫（http://www.worldclim.org/），空間分辨率2.5 min，包括22個(gè)環(huán)境因子（19個(gè)氣候變量及海拔、坡度和坡向）數(shù)據(jù)。選擇共享社會(huì)經(jīng)濟(jì)路徑SSP126和SSP585，其中SSP126代表低強(qiáng)迫、溫室氣體低排放情景，SSP585為高強(qiáng)迫、溫室氣體大量排放情景。采用ArcGIS10.2軟件對(duì)環(huán)境數(shù)據(jù)進(jìn)行提取、轉(zhuǎn)換，保存為MaxEnt模型運(yùn)行所需的.asc格式。

2 方法

2.1 環(huán)境因子篩選

將218 條樣點(diǎn)信息和22 個(gè)環(huán)境因子數(shù)據(jù)導(dǎo)入MaxEnt3.4.1軟件進(jìn)行迭代運(yùn)算，重復(fù)運(yùn)行10次，逐步舍去貢獻(xiàn)率為0的因子，直到所有因子貢獻(xiàn)率均大于0。由于不同因子間存在一定相關(guān)性，可能導(dǎo)致過擬合而影響判斷，因此將分布點(diǎn)對(duì)應(yīng)的環(huán)境因子（貢獻(xiàn)率＞0）數(shù)據(jù)提取后進(jìn)行Pearson分析，當(dāng)2個(gè)環(huán)境因子的相關(guān)系數(shù)|r|＞0.85時(shí)，剔除貢獻(xiàn)率較小者，得到影響當(dāng)歸潛在分布的環(huán)境因子。

2.2 模型構(gòu)建及評(píng)價(jià)

將218條當(dāng)歸樣點(diǎn)信息和篩選的環(huán)境因子數(shù)據(jù)導(dǎo)入MaxEnt3.4.1軟件，隨機(jī)選取75%分布點(diǎn)作為訓(xùn)練集建立預(yù)測模型，剩余25%分布點(diǎn)作為測試集驗(yàn)證模型，刀切法檢驗(yàn)各變量權(quán)重，選擇創(chuàng)建環(huán)境變量響應(yīng)曲線，其余參數(shù)為默認(rèn)值，重復(fù)10次，結(jié)果輸出形式為logistic，柵格數(shù)值為生存概率（P值），輸出為.asc格式文件。

通過MaxEnt3.4.1軟件建模，以（1-特異度）為橫坐標(biāo)、靈敏度為縱坐標(biāo)繪制受試者工作特征（ROC）曲線，用曲線下面積（AUC）評(píng)價(jià)模型準(zhǔn)確性。AUC 取值范圍0.5～1，越接近1表明模型預(yù)測精度越高。具體評(píng)價(jià)標(biāo)準(zhǔn)：0.5≤AUC＜0.6為擬合失敗，0.6≤AUC＜0.7為擬合效果較差，0.7≤AUC＜0.8為擬合效果一般，0.8≤AUC＜0.9為擬合效果好，0.9≤AUC＜1為擬合效果極好。

2.3 適生區(qū)等級(jí)劃分

將“2.2”項(xiàng)下輸出的.asc 格式文件導(dǎo)入ArcGIS10.2軟件，按P值進(jìn)行重分類，采用人工分級(jí)法分為5個(gè)等級(jí)：0～0.10為非適生區(qū)，0.11～0.30為低適生區(qū)，0.31～0.50為中適生區(qū)，0.51～0.70為高適生區(qū)，0.71～1.00為極高適生區(qū)[30-31]。利用ArcGIS柵格計(jì)算每個(gè)適生區(qū)的比例權(quán)重，再根據(jù)國土實(shí)際面積換算為每個(gè)適生區(qū)的面積。

2.4 適生區(qū)面積變化及質(zhì)心遷移

將“2.2”項(xiàng)下輸出的.asc 格式文件導(dǎo)入ArcGIS10.2軟件，使用SDM工具對(duì)模型運(yùn)行結(jié)果進(jìn)行二進(jìn)制分類，將適生區(qū)分割后進(jìn)行比較。為方便查看面積動(dòng)態(tài)變化，使用SDM工具下的MaxEnt工具進(jìn)行路徑轉(zhuǎn)化，將生成的.asc格式文件導(dǎo)入ArcGIS10.2中添加唯一值（-1、0、1、2）：-1代表增加的適生區(qū)面積，0代表未改變的非適生區(qū)面積，1代表未改變的適生區(qū)面積，2 代表減少的適生區(qū)面積。同時(shí)，運(yùn)用SDM工具計(jì)算當(dāng)前和未來時(shí)期2個(gè)路徑下各適生區(qū)的質(zhì)心坐標(biāo)，并繪制遷移趨勢圖，利用在線工具（https://www.box3.cn/tools/lbs.html）計(jì)算直線距離。

3 結(jié)果與分析

3.1 環(huán)境因子

影響當(dāng)歸分布的9個(gè)環(huán)境因子及其貢獻(xiàn)率見表1。貢獻(xiàn)率為海拔＞年均降水量＞最冷季度平均溫度＞溫度季節(jié)性變化的標(biāo)準(zhǔn)差＞等溫性＞晝夜溫差月均值＞坡向＞坡度＞降水量變異系數(shù)，其中海拔、年均降水量、最冷季度平均溫度和溫度季節(jié)性變化的標(biāo)準(zhǔn)差累計(jì)貢獻(xiàn)率達(dá)92.8%。由刀切法檢驗(yàn)結(jié)果（見圖1）可知，僅使用單一變量時(shí)，影響正規(guī)化訓(xùn)練增益的環(huán)境因子為海拔＞降水變異系數(shù)＞溫度季節(jié)性變化標(biāo)準(zhǔn)差＞最冷季度平均溫度。

圖1 影響當(dāng)歸分布環(huán)境因子刀切圖

表1 影響當(dāng)歸分布環(huán)境因子貢獻(xiàn)率

總體來看，影響當(dāng)歸分布的環(huán)境因子包括地形、溫度和降水，其中地形因子（海拔、坡度、坡向）貢獻(xiàn)率達(dá)50.8%，溫度因子（最冷季度平均溫度、溫度季節(jié)性變化的標(biāo)準(zhǔn)差、等溫性、晝夜溫差月均值）貢獻(xiàn)率為24.7%，降水因子（年均降水量、降水量變異系數(shù)）貢獻(xiàn)率為24.6%。

3.2 模型精度評(píng)價(jià)

當(dāng)前氣候條件下MaxEnt 預(yù)測模型AUC=0.950，未來氣候條件下預(yù)測模型AUC=0.949，表明模型擬合效果極好，結(jié)果可信度較高。見圖2、圖3。

圖2 當(dāng)前氣候條件下當(dāng)歸預(yù)測模型ROC曲線

圖3 未來氣候條件下當(dāng)歸預(yù)測模型ROC曲線

3.3 當(dāng)歸適生區(qū)分布

當(dāng)前氣候條件下，當(dāng)歸總適生區(qū)面積約為148.13×104km2。其中極高適生區(qū)面積4.8×104km2，主要分布在甘肅南部，青海東部、西藏東部、寧夏南部、四川中部、湖北西部、重慶北部、貴州西部和陜西中南部也有零星分布；高適生區(qū)面積約17.74×104km2，主要分布在甘肅南部及中部、青海東部、寧夏南部、云南北部和四川盆地西部，陜西中南部、重慶北部、湖北東部，西藏東部和貴州西部也有零星分布；中適生區(qū)面積約35.38×104km2，主要分布在甘肅南部及中部、四川北部及南部、云南中部，青海東部、西藏東部、寧夏陜西交界處和陜西中南部、重慶北部及中部、湖北西部、貴州西部也有零星分布；低適生區(qū)面積約90.21×104km2，主要分布在甘肅東部、青海東部、陜西中部及南部、寧夏中部、山西大部、河南西部、貴州大部、四川盆地周邊和湖北重慶交界處，福建、浙江、江西、湖南、廣東和廣西也有零星分布。見圖4。

圖4 當(dāng)前氣候條件下當(dāng)歸適生區(qū)分布

由表2可知，SSP126氣候情景下，當(dāng)歸總適生區(qū)面積呈先增加后減少再增加的變化趨勢，2041－2060年總適生區(qū)面積達(dá)到最大值，之后有所降低，到2081－2100年，不同等級(jí)適生區(qū)面積又有所增加，其中低適生區(qū)面積較當(dāng)前增加0.76%，中適生區(qū)面積增加0.68%，高適生區(qū)面積增加2.33%，極高適生區(qū)面積增加16.67%。SSP585氣候情景下，當(dāng)歸總適生區(qū)面積先增加后減少，2041－2060年總適生區(qū)面積達(dá)到最大值，主要增加的是低適生區(qū)面積，到2081－2100年，除低適生區(qū)面積增加1.35%外，其他等級(jí)適生區(qū)面積均有不同程度減少，其中中適生區(qū)面積較當(dāng)前減少5.17%，高適生區(qū)面積減少9.86%，極高適生區(qū)面積減少9.17%，總適生區(qū)面積減少1.89%。不同氣候情景下當(dāng)歸適生區(qū)分布預(yù)測見圖5、圖6。

圖5 未來SSP126氣候情景下當(dāng)歸適生區(qū)分布預(yù)測

圖6 未來SSP585氣候情景下當(dāng)歸適生區(qū)分布預(yù)測

表2 當(dāng)前和未來氣候情景下當(dāng)歸適生區(qū)面積（×104 km2）

3.4 當(dāng)歸適生區(qū)時(shí)空變化

不同氣候情景下當(dāng)歸適生區(qū)面積變化見表3。擴(kuò)張表示與當(dāng)前氣候環(huán)境相比適生區(qū)新增，不變表示與當(dāng)前氣候環(huán)境相比適生區(qū)未變化，收縮表示與當(dāng)前氣候環(huán)境件相比適生區(qū)減少。SSP126 氣候情景下，除2061－2080年外，擴(kuò)張面積均大于收縮面積，在此氣候情景下當(dāng)歸未來適生區(qū)面積總體呈增長趨勢。SSP585氣候情景下，2041－2060年擴(kuò)張面積大于收縮面積，在后續(xù)2個(gè)時(shí)期，收縮面積均大于擴(kuò)張面積，在此氣候情景下未來當(dāng)歸適生區(qū)面積逐漸縮小。由圖7可知，SSP126氣候情景下，甘肅中部和南部、寧夏南部、陜西中部和南部、青海東部、四川中部、云南北部、貴州大部、重慶湖北交界處適生區(qū)面積基本不變；山西北部、青海東部、福建浙江交界處有新增適生區(qū)；陜西中部、寧夏中部、山西中部及南部、四川西部、云南南部、貴州東部、福建南部適生區(qū)有縮小情況，且隨著時(shí)間推移，南方收縮面積較北方更多。由圖8可知，與SSP126氣候情景相比，SSP585氣候情景下青海東部擴(kuò)張面積更大，適生區(qū)不斷向高海拔遷移，且收縮面積遠(yuǎn)大于SSP126氣候情景，同時(shí)隨著時(shí)間推移，北方收縮面積較南方更大。

圖7 未來SSP126情景下當(dāng)歸適生區(qū)變化

圖8 未來SSP585情景下當(dāng)歸適生區(qū)變化

表3 不同氣候情景下當(dāng)歸適生區(qū)面積變化（×104 km2）

3.5 當(dāng)歸適生區(qū)質(zhì)心遷移

當(dāng)前當(dāng)歸適生區(qū)質(zhì)心位于四川省德陽市中江縣，坐標(biāo)30.986 189 46°N、104.633 186 1°E。SSP126氣候情景下，2041－2060年質(zhì)心從當(dāng)前點(diǎn)向北移動(dòng)22.76 km，仍位于四川省德陽市中江縣，坐標(biāo)31.190 946 08°N、104.628 298°E；2061－2080年質(zhì)心向西北移動(dòng)55.79 km，位于四川省德陽市綿竹市，坐標(biāo)31.387 679 03°N、104.087 755 4°E；2081－2100年質(zhì)心繼續(xù)向西北方向移動(dòng)，距離23.22 km，位于四川省德陽市綿竹市，坐標(biāo)31.580 069 78°N、103.992 185 6°E。SSP585氣候情景下，當(dāng)歸適生區(qū)質(zhì)心總體向西北方向遷移，2041－2060年位于四川省阿壩藏族羌族自治州茂縣，距當(dāng)前質(zhì)心98.23 km，坐標(biāo)31.590 599 79°N、103.878 383 3°E；2061－2080年質(zhì)心向南遷移，但距離較小，為7.73 km，位于四川省德陽市旌陽區(qū)，坐標(biāo)31.521 036 51°N、103.877 561 4°E；2081－2100 年質(zhì)心繼續(xù)向南遷移，距離21.62 km，位于四川省德陽市旌陽區(qū)，坐標(biāo)31.332 438 44°N、103.933 367 6°E。總體來看，當(dāng)歸適生區(qū)向高緯度地區(qū)遷移。見圖9。

圖9 未來氣候條件下當(dāng)歸適生區(qū)質(zhì)心遷移趨勢

4 討論

本研究建立的MaxEnt模型AUC值均大于0.9，說明研究結(jié)果準(zhǔn)確可靠。海拔、溫度、降水共同影響當(dāng)歸適生區(qū)分布，其中海拔對(duì)當(dāng)歸適生區(qū)分布影響最大，貢獻(xiàn)率為49.1%。研究發(fā)現(xiàn)，當(dāng)歸主要分布在第二階梯的云貴高原和黃土高原，海拔較低的四川盆地和第一階梯幾乎沒有分布，而海拔更高的青藏高原也無分布，但在新疆天山山脈和東南丘陵地帶有零星低適生區(qū)分布，可見低海拔或海拔過高均不適合當(dāng)歸生長。

未來不同氣候情景下，當(dāng)歸適生區(qū)面積有不同程度增加或減少，均在2041－2060年達(dá)到最大值，其中高適生區(qū)和極高適生區(qū)受氣候影響最為明顯，在SSP126氣候情景下2個(gè)適生區(qū)面積均有所增加，但在SSP585氣候情景下2個(gè)適生區(qū)面積均有所減少。當(dāng)歸適宜生長在高海拔、寒冷陰濕的地方，在其生長過程中存在低溫春化現(xiàn)象[32]，在SSP585氣候情景下，溫室效應(yīng)加劇，溫度升高，部分地區(qū)當(dāng)歸抽薹無法成藥的問題可能會(huì)愈發(fā)嚴(yán)重。目前當(dāng)歸主要分布在我國西北地區(qū)，降水過多或不足都會(huì)影響藥材的產(chǎn)量和質(zhì)量。在SSP126氣候情景下，當(dāng)歸適生區(qū)面積波動(dòng)范圍較小，山西北部、青海東部、福建浙江交界處有新增適生區(qū)，隨著時(shí)間推移，適生區(qū)面積減少的部分主要集中在我國南部，如云南東部、貴州東部及南部等地；在SSP585氣候情景下，當(dāng)歸適生區(qū)面積波動(dòng)范圍較大，新增適生區(qū)較SSP126氣候情景多，且主要集中在青海東部的高海拔地區(qū)，適生區(qū)面積減少部分主要集中在我國北部，如寧夏中部、甘肅中部及東部、陜西西部等地，且逐漸向甘肅省中高適生區(qū)靠近。當(dāng)歸適生區(qū)質(zhì)心整體北移，向高緯度方向移動(dòng)，高排放造成的溫室效應(yīng)對(duì)當(dāng)歸的適生區(qū)面積及分布造成了一定影響，因此，環(huán)境保護(hù)對(duì)氣候穩(wěn)定及當(dāng)歸的可持續(xù)開發(fā)利用至關(guān)重要。

根據(jù)本研究結(jié)果，當(dāng)歸高適生區(qū)和極高適生區(qū)對(duì)氣候變化更為敏感，甘肅南部面積較為穩(wěn)定，可繼續(xù)進(jìn)行當(dāng)歸大規(guī)模栽培生產(chǎn)，青海東部適生區(qū)增加，可適當(dāng)引種當(dāng)歸。當(dāng)歸適生區(qū)分布受氣候、土壤和地形等環(huán)境因子的影響，本研究僅考慮了氣候、地形條件，未對(duì)土壤、植被及人為因素進(jìn)行考察，因此預(yù)測結(jié)果可能存在一定偏差，后續(xù)將進(jìn)一步完善。