內(nèi)蒙古地表凍融指數(shù)動態(tài)變化與驅(qū)動因素分析

張昊琛， 薩楚拉， 孟凡浩， 羅 敏， 王牧蘭，高紅豆， ADIYA Saruulzaya

（1.內(nèi)蒙古師范大學(xué)地理科學(xué)學(xué)院，內(nèi)蒙古 呼和浩特 010022；2.內(nèi)蒙古自治區(qū)遙感與地理信息系統(tǒng)重點實驗室，內(nèi)蒙古 呼和浩特 010022；3.蒙古科學(xué)院地理與地球生態(tài)研究所，蒙古 烏蘭巴托 15170）

凍土一般是指溫度在0 ℃或0 ℃以下，并含有冰的各種巖土和土壤，包括永久凍土(PF)和季節(jié)性凍土(SFG)。我國是世界凍土大區(qū)，多年凍土和季節(jié)凍土分別占國土面積的22%和53.5%［1］。IPCC 第5次評估顯示，全球地表溫度在1880—2012年期間升高了0.85 ℃［2］。全球變暖在高緯度和高海拔區(qū)域加劇，對凍土的形成有較大的影響，尤其是在我國北方高緯度，溫度升高幅度大于全國平均升溫水平［3-5］。隨著氣候變化，我國凍土大致體現(xiàn)在：最大凍土深度下降、凍結(jié)日期推遲、消融日期提前、凍土下限上升［6-7］。20 世紀(jì)70 年代至21 世紀(jì)初，我國凍土面積減少約18.6%，面積由2.15×106km2減少到1.75×106km2［8］。按照對未來的天氣情況預(yù)報，未來冰凍圈各組成部分面積將持續(xù)減少，而多年凍土的活動層厚度將會進一步增加，季節(jié)凍土的凍結(jié)深度也將進一步下降［9-10］。

凍融指數(shù)是研究凍土特性的重要參數(shù)，也是判定其存在和長期變化的一個重要指標(biāo)［11］。凍土研究中凍融指數(shù)有助于繪制和預(yù)測永久凍土分布，估算土壤凍結(jié)深度，并計算活躍層厚度，凍結(jié)指數(shù)的下降和融化指數(shù)的上升在一定程度上說明凍土環(huán)境的退化［12］。利用季節(jié)凍土區(qū)的凍結(jié)指數(shù)，可以預(yù)測季節(jié)凍融深度［13］；利用多年凍土區(qū)的凍融指數(shù)，不僅能夠繪制出凍土的分布圖，還能對多年凍土的活動層厚度進行預(yù)測［14-15］；Peng 等［12］的研究結(jié)果顯示，北半球高緯度和高海拔地區(qū)的凍結(jié)指數(shù)下降幅度最大，凍結(jié)指數(shù)的下降速率要比融化指數(shù)的上升速率快得多；Shi 等［16］研究表明，從1901—2015 年，環(huán)北極地區(qū)的凍結(jié)和融化指數(shù)分別呈現(xiàn)降低、上升趨勢，在1988年之后，其變化趨勢更為明顯；Liu等［17］基于地面溫度觀測資料和再分析產(chǎn)品計算了雅魯藏布江上游凍融指數(shù)，并基于凍融指數(shù)模擬了凍土面積的變化，效果良好；Wu等［18］對青藏高原凍融指數(shù)的變化規(guī)律進行了分析，發(fā)現(xiàn)多年凍土區(qū)的凍結(jié)指數(shù)比季節(jié)凍土區(qū)的凍結(jié)指數(shù)降低速度要快，而融化指數(shù)上升速率更慢；Wang 等［19］研究認(rèn)為，黃河源區(qū)海拔高度對凍融指數(shù)的空間分布有一定的影響，地表凍融指數(shù)的變化比空氣凍融指數(shù)的變化要大。

內(nèi)蒙古地區(qū)的氣象觀測站數(shù)量較少，有些站點建設(shè)時間較晚，總體分布稀疏，個別站點存在數(shù)據(jù)缺失問題，氣候變化條件下的凍融指數(shù)變化研究較為欠缺。僅有的研究主要聚焦在典型草原或多年凍土區(qū)，且研究對象為空氣凍融指數(shù)［20-21］，不能完全代表內(nèi)蒙古全區(qū)的地表凍融特征，并且研究區(qū)地表凍融指數(shù)變化的驅(qū)動因素有待探究。鑒于此，本文基于內(nèi)蒙古氣象站1980—2019 年日均地表溫度數(shù)據(jù)、結(jié)合中國第一代全球陸面再分析產(chǎn)品（CRA）數(shù)據(jù)以及NDVI 數(shù)據(jù)，對內(nèi)蒙古近40 a 地表凍融指數(shù)時空變化特征及驅(qū)動因素進行分析，該研究能夠為下一步估算內(nèi)蒙古季節(jié)凍土凍結(jié)深度、多年凍土活動層厚度、模擬凍土面積變化等提供重要參數(shù)數(shù)據(jù)，以及為內(nèi)蒙古農(nóng)牧業(yè)生產(chǎn)布局、草場生態(tài)評價及凍害防治等方面提供一定的科學(xué)依據(jù)。

1 數(shù)據(jù)與方法

1.1 研究區(qū)概況

內(nèi)蒙古位于126°04′～97°12′E，37°24′～53°23′N之間（圖1），地處北部邊疆，總面積118.3×104km2，占全國總面積的12.3%。地貌以高原為主，兼有平原、沙地、丘陵與山地，大部分地區(qū)海拔在1000 m以上［22］。氣候?qū)贉貛Т箨懶约撅L(fēng)氣候，處于干旱與半干旱氣候帶，年平均氣溫0～8 ℃，且由東北向西南遞增；全區(qū)降水量少且分布不均，年總降水量50～450 mm，由東北向西南遞減［23］。氣候分異明顯使得研究區(qū)域從東北向西南呈現(xiàn)出寒溫帶針葉林、落葉闊葉林、草甸草原、典型草原、荒漠草原等植被分布［22］。研究區(qū)位于環(huán)北極地區(qū)多年凍土區(qū)歐亞大陸多年凍土南界。冬季嚴(yán)寒，夏季高溫，且局地日夜溫差大，春冬季土壤凍融現(xiàn)象突出，夏秋季土壤干濕交替明顯［24］。多年凍土主要分布于呼倫貝爾地區(qū)，其余地區(qū)廣泛分布季節(jié)凍土，淺層表層土壤冬季結(jié)冰，夏季完全解凍。

圖1 內(nèi)蒙古氣象站點/凍土區(qū)分布及土地覆蓋類型Fig.1 Distribution of meteorological stations/permafrost regions and land cover types in Inner Mongolia

1.2 數(shù)據(jù)來源與處理

使用的45 個氣象臺站1980—2019 年日平均0 cm 地表溫度觀測數(shù)據(jù)來自于國家氣象信息中心的中國地面氣候資料日值數(shù)據(jù)集（http://data.cma.cn/）。為更好地體現(xiàn)研究區(qū)的凍融特征，盡量采用過去40 a來不同觀測站的共同時間數(shù)據(jù)來求取地表凍融指數(shù)。對缺測數(shù)據(jù)采用如下方式內(nèi)插［25］：缺測1 d，對其前后各1 d數(shù)據(jù)的平均值內(nèi)插；缺測2 d，在缺測當(dāng)天選擇該日前2 d 的數(shù)值取平均，缺測第2 d選取該日后2 d的數(shù)值取平均；如果有連續(xù)2 d以上的缺測，可以利用對相鄰站同一月份的溫度平均值做線性相關(guān)性檢查，再選取溫度相關(guān)性最佳的站點用線性回歸的溫度估計值內(nèi)插。采用反距離插值法(inverse distance weighted)對內(nèi)蒙古45個氣象站點的地表凍融指數(shù)進行空間插值，插值分辨率為34 km。

凍土類型分區(qū)依據(jù)青藏高原科學(xué)數(shù)據(jù)中心(http://data.tpdc.ac.cn)提供的1:1000凍土分區(qū)和類型地圖資料進行劃分，通過對凍土類型整合，整體劃分為多年凍土區(qū)和季節(jié)凍土區(qū)，用來對地表凍融指數(shù)和影響因素分區(qū)討論。0.1 m和0.4 m表層土壤含水量、降水量、雪深數(shù)據(jù)均源自我國國家氣象信息中心的中國第一代全球大氣和陸面再分析產(chǎn)品CRA/Land-逐月陸面產(chǎn)品（http://data.cma.cn/），空間分辨率為34 km，用來分析內(nèi)蒙古地表凍融指數(shù)影響因素時空變化。

GIMMS NDVI3g 數(shù)據(jù)由美國馬里蘭大學(xué)GLCF(Global Land Cover Facility)研究組生產(chǎn)的15 d 合成產(chǎn) 品 數(shù) 據(jù) (https://ecocast.arc.nasa.gov/data/pub/gimms/)，MODIS NDVI 數(shù)據(jù)為美國國家航天航空局(NASA)提供的MODIS 16 d合成產(chǎn)品數(shù)據(jù)(https://mo?dis.gsfc.nasa.gov/data/)。在GIMMS NDVI 和MODIS NDVI 的基礎(chǔ)上，采用線性回歸分析，對2 種數(shù)據(jù)進行時間序列的重建，從而使數(shù)據(jù)重采樣到同一個時空尺度上。使用重合時間的GIMMS NDVI 和MO?DIS NDVI 進行建模，根據(jù)建模結(jié)果對GIMMS 數(shù)據(jù)作延長，獲取1982—2019年NDVI數(shù)據(jù)。對于NDVI數(shù)據(jù)的預(yù)處理，主要分為坐標(biāo)變換、數(shù)值裁剪、數(shù)據(jù)轉(zhuǎn)換、異常值刪除、月最大值獲取等過程，使參與計算的所有NDVI 值在0～1 范圍內(nèi)，根據(jù)研究區(qū)土地覆蓋分類圖（圖1），結(jié)合研究區(qū)實際情況，認(rèn)為ND?VI＜0.15 的地區(qū)無植被覆蓋。考慮到數(shù)據(jù)空間分辨率的一致性，對NDVI數(shù)據(jù)重采樣為34 km。

1.3 研究方法

（1）地表凍融指數(shù)計算?；贙lene 等［26］和Frauenfeld等［15］對凍結(jié)期和融化期的限制，界定凍結(jié)期介于每年7 月1 日至次年6 月30 日，融化期介于每年1月1日至12月31日，以最大可能統(tǒng)計一年中寒冷期（溫度＜0 ℃）和溫暖期（溫度＞0 ℃）的溫度累加值，并確保負(fù)溫和正溫計算分別在整個寒冷期和溫暖期內(nèi)。年凍結(jié)指數(shù)計算如下：

式中：FI和TI分別為凍結(jié)指數(shù)和融化指數(shù)(℃·d)；Ti和Tj分別為負(fù)值逐日地表溫度和正值逐日地表溫度，單位為℃；Nf和NT為年內(nèi)地表溫度小于和大于0 ℃的天數(shù)，單位為d。

（2）趨勢分析。利用一元線性回歸方法計算凍融指數(shù)和影響因素的年際變化趨勢，并通過最小二乘法擬合的斜率反映，當(dāng)斜率大于0時，說明有增大的傾向，當(dāng)斜率小于0時，說明有下降的傾向。其變化趨勢的計算公式為：

式中：Slope表示地表凍融指數(shù)或影響因素的變化趨勢；n為時間序列長度；為研究對象X第i年的平均值。

（3）相關(guān)分析法。通過計算地表凍融指數(shù)與對應(yīng)的影響因素相關(guān)系數(shù)，表示其與地表凍融指數(shù)的相關(guān)性。計算公式如下：

式中：n表示研究時間段的累計年數(shù)；rxy表示x和y的相關(guān)系數(shù)；xi和yi分別表示第i年x和y的值；和表示2 個要素樣本值的平均值；rxy在-1～1 之間，隨著絕對值的增大，地表凍融指數(shù)與影響因素的變化關(guān)系也相應(yīng)較好；反之越差，當(dāng)rxy＞0 表示相關(guān)性為正，rxy＜0，表示相關(guān)性為負(fù)［27］。最后采用T檢驗法進行顯著性檢驗，將地表凍融指數(shù)與影響因素的相關(guān)程度劃分為顯著相關(guān)（0.01≤P＜0.05）和不顯著相關(guān)（P≥0.05）。

（4）灰色關(guān)聯(lián)分析。采用灰色關(guān)聯(lián)分析法，研究影響地表凍融指數(shù)的主要因子及重要性?；疑P(guān)聯(lián)度可以測量2 個系統(tǒng)之間的關(guān)系，并確定相關(guān)因子中的主控因子［28］。地表凍融指數(shù)的時間序列作為目標(biāo)變量，將0.1 m土壤含水量、0.4 m土壤含水量、降水量等影響因素作為影響變量?；疑P(guān)聯(lián)系數(shù)表達為：

2 結(jié)果與分析

2.1 內(nèi)蒙古地表凍融指數(shù)時空變化

2.1.1 地表凍融指數(shù)的時間變化特征 通過對內(nèi)蒙古地區(qū)1980—2019 年的45 個觀測站SFI（地表凍結(jié)指數(shù)）和STI（地表融化指數(shù)）進行分析，得到研究區(qū)40 a 整體和不同凍土區(qū)SFI 與STI 的時間變化特征。如圖2a和圖2b所示，整體上SFI近40 a呈現(xiàn)顯著下降趨勢，多年凍土區(qū)SFI 下降速率高于季節(jié)凍土區(qū)，內(nèi)蒙古SFI 多年變化范圍介于956.1～1848.3 ℃·d，最低值和最高值分別發(fā)生于2016年和1984年，年均值為1457.0 ℃·d，變化率為-156.4 ℃·d·（10a）-1(P＜0.01)，表明土壤的季節(jié)凍結(jié)潛力有連續(xù)降低的趨勢［30］。進一步分析各站點SFI 變化特征，近40 a 所有站點的SFI 均呈下降趨勢，但各站SFI變化特征差異也較大，變化最小值和最大值分別出現(xiàn)于阿拉善盟阿拉善右旗站點［-20.7 ℃·d·（10a）-1]和興安盟阿爾山站點［-711.1 ℃·d·（10a）-1]。如圖2c 和圖2d 所示，整體上STI 近40 a 呈現(xiàn)顯著上升趨勢，多年凍土區(qū)STI上升速率高于季節(jié)凍土區(qū)，內(nèi)蒙古STI 多年變化范圍介于3717.6～4442.3 ℃·d，最低值和最高值分別發(fā)生于1984 年和2007 年，年均值為4040.9 ℃·d，上升速率為152.4 ℃·d·（10a）-1(P＜0.01)，說明土壤融化期間最大融化深度增加［30］。1980—2019年間，所有站點的STI均呈上升趨勢，各站STI 變化特征存在差異，變化特征最小值和最大值分別出現(xiàn)于赤峰市赤峰站點［43.1 ℃·d·（10a）-1］和赤峰市巴林左旗站點［221.1 ℃·d·（10a）-1］。

圖2 1980—2019年內(nèi)蒙古整體和不同凍土區(qū)地表凍融指數(shù)年際變化趨勢Fig.2 Interannual variation trend of surface freezing-thawing indices in Inner Mongolia and different permafrost regions from 1980 to 2019

2.1.2 地表凍融指數(shù)空間分布與變化 內(nèi)蒙古SFI多年平均值的空間分布特征如圖3a 所示，受經(jīng)度、緯度和海拔的綜合影響，SFI 均值總體呈現(xiàn)西低東高、由西南向東北逐漸增加的特征，其中SFI受緯度的影響最大（表1），相關(guān)系數(shù)達到0.92（P＜0.01），每上升1°，SFI 增加165.2 ℃·d，該分布特征同環(huán)北極凍融指數(shù)分布變化特征一致［12,16］。SFI 高值區(qū)位于呼倫貝爾北部地區(qū)，SFI 低值區(qū)位于阿拉善西部地區(qū)，與凍土區(qū)分布狀況類似，多年凍土區(qū)SFI 普遍大于季節(jié)凍土區(qū)。內(nèi)蒙古SFI 年際變化率空間分布特征見圖3b，SFI 年際變化率的變化范圍介于-10.3～-710.6 ℃·d·（10a）-1，SFI 變化率在空間上自西南向東北呈增加趨勢，表明內(nèi)蒙古東北地區(qū)SFI 下降速率較西南地區(qū)更快。

內(nèi)蒙古STI 多年平均值的空間分布特征見圖3c，STI 均值的空間分布特征呈現(xiàn)西部高，東部低，由西南向東北逐漸減小的趨勢，多年凍土區(qū)的STI普遍小于季節(jié)凍土區(qū)，在呼和浩特、赤峰等地出現(xiàn)了點狀高值區(qū)，這與近幾年農(nóng)業(yè)開發(fā)、城市擴張、工業(yè)發(fā)展等社會經(jīng)濟活動密切相關(guān)。同時，上述地區(qū)處于生態(tài)環(huán)境較為脆弱，人為活動對環(huán)境的影響更加敏感。STI同緯度的相關(guān)系數(shù)最高（表1），為-0.73（P＜0.01），每上升1°STI 減少141.1 ℃·d，STI 高值區(qū)位于阿拉善西北部地區(qū)，STI 低值區(qū)位于呼倫貝爾中部地區(qū)。內(nèi)蒙古STI年際變化率空間分布特征見圖3d，STI上升速率較高的地區(qū)（Slope＞150 ℃·d·（10a）-1）主要分布在研究區(qū)中部、中東部和東北部分地區(qū)，可能是因為STI 上升速率較高的地區(qū)受人類活動等因素的影響，其次，呼倫貝爾高原大部分地區(qū)的凍融指數(shù)變化率最大，這也證實了該區(qū)域（原島狀凍土區(qū)）目前處于季節(jié)性凍土的狀態(tài)［31］。

表1 地表凍融指數(shù)與地理因子的相關(guān)系數(shù)Tab.1 Correlation coefficient between surface freezingthawing index and geographical factors

圖3 1980—2019年內(nèi)蒙古地表凍融指數(shù)均值和變化趨勢Fig.3 Mean value and variation trend of surface freezing-thawing index in Inner Mongolia from 1980 to 2019

2.2 內(nèi)蒙古地表凍融指數(shù)影響因素時空變化

2.2.1 內(nèi)蒙古地表凍融指數(shù)影響因素時間變化 如圖4 所示，近40 a 研究區(qū)0.1 m 和0.4 m 土壤含水量年均值分別介于0.1401～0.1673 m3·m-3和0.1388～0.1779 m3·m-3之間，0.4 m 土壤含水量略大于0.1 m土壤含水量，不同深度土壤含水量的總體變化趨勢均顯示出向上態(tài)勢，平均年際變化率為0.001 m3·m-3·（10a）-1，而0.1 m 和0.4 m 土壤含水量并不是呈現(xiàn)一致性的向上態(tài)勢，而是隨氣候波動經(jīng)歷了“上升-下降-上升”3 個波動性變化階段，且0.1 m 和0.4 m 深度的土壤含水量年際變化趨勢大體相同：1980—1998年為波動上升階段，1998—2007年呈波動下降趨勢，2007—2019 年又呈波動上升趨勢；降水量以10.0 mm·（10a）-1的速率呈不顯著增加趨勢(P＞0.05)，年均值介于240.1～415.9 mm 之間，降水量與0.1 m和0.4 m土壤含水量的年際波動變化趨勢極為相似，同樣經(jīng)歷了“上升-下降-上升”3 個波動變化階段，1980—1998 年、2007—2019 年屬于降水偏豐的年代；1998—2007 年屬于降水偏枯的年代；研究區(qū)雪深以11.4～32.8 mm 的起伏趨勢為主，多年積雪厚度為22.1 mm，年平均雪深以-0.82 mm·（10a）-1的速率呈微弱減少趨勢(P＞0.05)；基于1982—2019年內(nèi)蒙古NDVI 數(shù)據(jù)，發(fā)現(xiàn)內(nèi)蒙古NDVI 總體以0.007·（10a）-1的速率呈顯著增加趨勢(P＜0.05)，年均值介于0.41～0.47之間。

圖4 1980—2019年影響因素年際變化趨勢Fig.4 Interannual variation trend of impact factors from 1980 to 2019

2.2.2 內(nèi)蒙古地表凍融指數(shù)影響因素空間分布與變化 如圖5 所示，研究區(qū)內(nèi)各影響因素空間分布整體上存在明顯的緯度地帶性和區(qū)域性規(guī)律，且均呈現(xiàn)出由西南向東北遞減的發(fā)展趨勢，與SFI 空間分布特征一致，與STI 空間分布特征相反。但上述要素在部分地區(qū)空間分布和空間變化存在一定差異：0.1 m和0.4 m土壤含水量呈現(xiàn)出東北部和中西部土壤水分相對較高、西部相對較低的空間格局，可能受到不同土地覆蓋類型、積雪分布及人為活動等多重因素的共同作用；由于大興安嶺、陰山、賀蘭山等地勢因素的阻礙，東南地區(qū)的水汽抬升，容易造成降雨，西北山區(qū)的水汽含量偏低，因此東南地區(qū)的降水量比西北地區(qū)大；雪深與降水量分布特征基本一致，但與溫度相反，研究區(qū)東北部降水量較多、年均地表溫度較低且森林覆蓋率高，有利于積雪的形成和積累，中部和西部地區(qū)降水較少，平均地表溫度較高，植被覆蓋率較低，不利于積雪的形成與累積，使得雪深較淺且分布均勻；內(nèi)蒙古NDVI由西南向東北呈增加趨勢，部分地區(qū)NDVI 存在區(qū)域性特征和空間異質(zhì)性較明顯，可能與灌溉和開墾農(nóng)田等人類活動有關(guān)。

圖5 1980—2019年地表凍融指數(shù)影響因素多年平均Fig.5 Annual average of surface freezing-thawing index influencing factors from 1980 to 2019

如圖6所示，近40 a來內(nèi)蒙古大部分地區(qū)0.1 m和0.4 m 土壤含水量年際變化呈上升態(tài)勢，0.1 m 和0.4 m 土壤含水量呈減少態(tài)勢的區(qū)域分別僅占研究區(qū)整體的19.9%和36.1%，且大部分集中在研究區(qū)西北部和東北部的多年凍土區(qū)；內(nèi)蒙古中部和西部的年降水量有所增加，特別是中偏南地區(qū)降水量增加較多，東部降水量則有所減少，特別是在呼倫貝爾西北部、赤峰中部以及通遼中部地區(qū)降水量減少情況較為突出；研究區(qū)38.5%的區(qū)域雪深呈下降態(tài)勢，且10.8%的區(qū)域雪深平均下降速度在0.5 mm·（10a）-1以上，主要分布在呼倫貝爾東北部，并位于多年凍土區(qū)；61.5%的區(qū)域雪深呈現(xiàn)增加趨勢，主要分布于西南部的沙漠及中部和東部的部分草原，多位于季節(jié)凍土區(qū)，且以上區(qū)域的雪深平均增加速率僅為0.4 mm·（10a）-1，這一變化特征可能與研究區(qū)近年來溫度和降水變化有關(guān)；研究時段內(nèi)蒙古有75.2%區(qū)域的NDVI 呈現(xiàn)增長的趨勢，僅24.8%的區(qū)域呈現(xiàn)下降的趨勢，下降區(qū)域以錫林郭勒盟大部、赤峰市北部、通遼市北部、興安盟南部、呼倫貝爾市林區(qū)為主，這一變化特征可能與上述地區(qū)近年來溫度升高、降水減少及人類活動有關(guān)。

圖6 1980—2019年地表凍融指數(shù)影響因素年際變化Fig.6 Interannual variation of surface freezing-thawing index influencing factors from 1980 to 2019

2.3 地表凍融指數(shù)與影響因素相關(guān)分析

地表凍融指數(shù)的變化受多種因素的影響，具有較強的隨機性，適于進行灰色關(guān)聯(lián)度分析。如表2所示，內(nèi)蒙古地表凍融指數(shù)年均值與0.1 m 土壤含水量、0.4 m土壤含水量、降水量、雪深、NDVI的灰色關(guān)聯(lián)度系數(shù)皆大于0.65，說明地表凍融指數(shù)變化與各影響因素存在較大的相關(guān)性。SFI與影響因素的灰色關(guān)聯(lián)序為：0.4 m 土壤含水量＞0.1 m 土壤含水量＞NDVI＞降水量＞雪深，說明了0.4 m 土壤含水量對SFI變化的影響較大，是影響內(nèi)蒙古整體SFI變化的主控因素，0.1 m 土壤含水量、NDVI、降水量次之，最后是雪深；STI 與影響因素的灰色關(guān)聯(lián)序為：NDVI＞0.1 m 土壤含水量＞0.4 m 土壤含水量＞降水量＞雪深，表明NDVI 對STI 變化的影響較大，是影響內(nèi)蒙古整體STI變化的主控因素，0.1 m土壤含水量、0.4 m土壤含水量、降水量次之，最后是雪深。

表2 地表凍融指數(shù)與影響因素灰色關(guān)聯(lián)系分析Tab.2 Grey correlation analysis between surface freezing-thawing index and influencing factors

為探討地表凍融指數(shù)對影響因素變化的動態(tài)響應(yīng)，基于像元尺度逐像元計算了內(nèi)蒙古地表凍融指數(shù)與影響因素的相關(guān)關(guān)系。SFI 與各影響因素（除雪深外）相關(guān)關(guān)系有超過65%的地區(qū)呈負(fù)相關(guān)（表3），如圖7所示，從空間分布來看，在多年凍土區(qū)大部分地區(qū)呈現(xiàn)正相關(guān)關(guān)系，但在季節(jié)凍土區(qū)大部分地區(qū)卻呈現(xiàn)負(fù)相關(guān)關(guān)系，造成不同凍土區(qū)SFI 與影響因素相關(guān)性差異的原因可能與影響因素的時空變化差異有關(guān)。其中，0.1 m土壤含水量與SFI顯著負(fù)相關(guān)(P＜0.05)的區(qū)域大部分在季節(jié)凍土區(qū)的中部和東南部，顯著正相關(guān)(P＜0.05)區(qū)域集中分布在多年凍土區(qū)北部小部分地區(qū)，0.4 m 土壤含水量與SFI 顯著負(fù)相關(guān)(P＜0.05)的地區(qū)在季節(jié)凍土區(qū)的南部和東南部均有分布，顯著正相關(guān)(P＜0.05)的地區(qū)則在多年凍土區(qū)北部呈現(xiàn)片狀分布；降水量與SFI顯著負(fù)相關(guān)(P＜0.05)的地區(qū)主要在鄂爾多斯市以東、呼和浩特市和錫林郭勒盟以南，呈正相關(guān)關(guān)系的區(qū)域主要在多年凍土區(qū)北部和阿拉善中西部；NDVI 與SFI顯著負(fù)相關(guān)(P＜0.05)的區(qū)域主要在鄂爾多斯東部、呼和浩特市南部和多年凍土區(qū)部分林地，顯著正相關(guān)(P＜0.05)區(qū)域主要在多年凍土區(qū)部分林地和濕地地區(qū)；雪深與SFI顯著負(fù)相關(guān)(P＜0.05)的區(qū)域主要在阿拉善西南部和季節(jié)凍土區(qū)中部地區(qū)零星分布，顯著正相關(guān)(P＜0.05)區(qū)域集中在多年凍土區(qū)北部和阿拉善西北部小片地區(qū)。

圖7 SFI與影響因素相關(guān)性T檢驗Fig.7 Correlation T test between SFI and influencing factors

表3 SFI與影響因素相關(guān)系數(shù)面積占比Tab.3 Area ratio of correlation coefficient between SFI and influencing factors

STI 與影響因素的相關(guān)關(guān)系以負(fù)相關(guān)為主（表4）。如圖8所示，雖然STI與各影響因素在多年凍土區(qū)及大部分地區(qū)呈現(xiàn)負(fù)相關(guān)關(guān)系，但在季節(jié)凍土區(qū)部分地區(qū)卻呈現(xiàn)正相關(guān)關(guān)系，造成不同凍土區(qū)STI與影響因素相關(guān)性差異的原因可能與影響因素的時空變化差異有關(guān)。其中，0.1 m 土壤含水量與STI顯著負(fù)相關(guān)(P＜0.05)區(qū)域集中在多年凍土區(qū)東部和多年凍土區(qū)周圍的盟市，顯著正相關(guān)(P＜0.05)的區(qū)域在黃河流域內(nèi)蒙古地區(qū)有所分布；0.4 m土壤含水量與STI 顯著負(fù)相關(guān)(P＜0.05)區(qū)域主要在多年凍土區(qū)北部和阿拉善西部零星地區(qū)，顯著正相關(guān)(P＜0.05)的區(qū)域主要在阿拉善東部、巴彥淖爾市大部和鄂爾多斯市周邊；降水量與STI顯著負(fù)相關(guān)(P＜0.05)區(qū)域主要在錫林郭勒盟東部、興安盟東部及多年凍土區(qū)東北部，呈正相關(guān)關(guān)系的區(qū)域主要在季節(jié)凍土區(qū)中部和西部；雪深與STI顯著負(fù)相關(guān)(P＜0.05)區(qū)域主要在錫林郭勒盟西北部、赤峰市北部及多年凍土區(qū)東北部，顯著正相關(guān)(P＜0.05)的區(qū)域集中在巴彥淖爾市中部。NDVI 與SFI 顯著負(fù)相關(guān)(P＜0.05)的區(qū)域主要在赤峰市北部、興安盟、多年凍土區(qū)部分林地和錫林郭勒盟荒漠草原地區(qū)，顯著正相關(guān)(P＜0.05)的區(qū)域主要在黃河流域內(nèi)蒙古地區(qū)、通遼市和多年凍土區(qū)北部。

圖8 STI與影響因素相關(guān)性T檢驗Fig.8 Correlation T test between STI and influencing factors

表4 STI與影響因素相關(guān)系數(shù)面積占比Tab.4 Area ratio of correlation coefficient between STI and influencing factors

3 討論

氣候變暖背景下，內(nèi)蒙古1980—2019年間地表凍結(jié)和融化指數(shù)分別顯著上升和下降，多年凍土區(qū)地表凍融指數(shù)變化速率較季節(jié)凍土區(qū)更快。研究結(jié)果基本與Peng等［12］、Wu等［18］、李昀赟等［32］分別在北半球、青藏高原、東北地區(qū)的結(jié)論相吻合，說明在內(nèi)蒙古整體凍土環(huán)境退化明顯，多年凍土區(qū)冬季溫度變暖趨勢更為顯著，其凍土退化更嚴(yán)重，該區(qū)地表凍融指數(shù)的顯著變化可能導(dǎo)致多年凍土區(qū)面積縮減，凍結(jié)深度減少，活動層厚度增加等一系列凍土退化現(xiàn)象的發(fā)生［12］；季節(jié)凍土區(qū)的地表凍融指數(shù)變化速率同樣不容忽視，該區(qū)是內(nèi)蒙古農(nóng)牧業(yè)生產(chǎn)主要布局地區(qū)，凍融指數(shù)的變化同農(nóng)作物的種植，畜牧業(yè)的發(fā)展息息相關(guān)。氣候變暖雖然是造成該區(qū)多年凍土區(qū)凍土退化顯著的根本性因素，但是多種人為經(jīng)濟活動對研究區(qū)多年凍土退化起到了加速促進作用［24］。近40 a研究區(qū)表層土壤含水量、降水量的增加趨勢表明，內(nèi)蒙古地區(qū)的凍土環(huán)境正向暖濕化發(fā)展，與馬梓策等［33］得出的研究結(jié)論相一致。但是基于本文研究發(fā)現(xiàn)，不同凍土區(qū)的發(fā)展情況存在較大差異，多年凍土區(qū)正向暖干化發(fā)展，與研究區(qū)整體和季節(jié)凍土區(qū)大部分地區(qū)變化趨勢相反。地面溫度升高、降水減少、潛在蒸發(fā)量增加是導(dǎo)致土壤含水量下降的主要因素［34］。地表凍融指數(shù)是根據(jù)地表溫度計算所得，與氣溫和地表溫度有密切的聯(lián)系，但表層土壤含水量、降水量、雪深和NDVI 同樣對凍融指數(shù)的變化起著至關(guān)重要的作用，總體上，表層土壤含水量和NDVI對地表凍融指數(shù)的變化起主導(dǎo)作用。此外，凍土環(huán)境退化對生物多樣性的作用是負(fù)面的。溫度增高、雨量減少的情況對草原植物生長發(fā)育不利，對森林、濕地植物有促進作用［35］；季節(jié)凍土地帶冬季增溫明顯，農(nóng)作物凍害損失減少，從而促進了糧食作物單產(chǎn)的增加和質(zhì)量提高，利于家畜健康過冬度春，同時也會使病蟲害及蟲卵易于過冬，容易造成病蟲害的防治困難［36］。

20世紀(jì)80年代前研究區(qū)氣象數(shù)據(jù)缺測較多，本研究只采用了1980—2019年的氣象資料，下一步還需依據(jù)更長時間序列和更多的站點數(shù)據(jù)開展研究。在對地表凍融指數(shù)處理后，需要根據(jù)凍融指數(shù)，進一步確定研究區(qū)多年凍土空間分布范圍和活動層厚度及季節(jié)凍土凍結(jié)深度的時空變化特征，對內(nèi)蒙古地區(qū)凍土對氣候變化的響應(yīng)作出較為全面、精確的預(yù)測。目前，由于中國全球陸面再分析產(chǎn)品(CRA/Land)數(shù)據(jù)在凍土研究中仍具有空間分辨率略小的缺陷，在今后凍土變化研究中，擬將與空間分辨率更高的再分析資料及實際觀測資料相結(jié)合，以盡可能降低研究結(jié)論的不確定性。

4 結(jié)論

本文基于內(nèi)蒙古氣象站1980—2019 年日均地表溫度數(shù)據(jù)、結(jié)合中國第一代全球陸面再分析產(chǎn)品（CRA）數(shù)據(jù)以及NDVI 數(shù)據(jù)，借助GIS 空間分析功能，通過RStudio編程處理分析了內(nèi)蒙古近40 a地表凍融指數(shù)時空變化格局及關(guān)鍵影響要素，主要結(jié)論如下：

（1）1980—2019年間，內(nèi)蒙古SFI與STI分別呈現(xiàn)明顯減少與增加的態(tài)勢，緯度是制約內(nèi)蒙古地表凍融指數(shù)空間分布的重要因素，與SFI 呈現(xiàn)極明顯正相關(guān)的關(guān)系，則與STI 相反。多年凍土區(qū)地表凍融指數(shù)變化速率較季節(jié)凍土區(qū)更為劇烈，對氣候變暖的響應(yīng)更為敏感。

（2）近40 a 研究區(qū)0.1 m 和0.4 m 土壤含水量、降水量、NDVI 均呈增加趨勢，僅雪深呈減少趨勢。其中不同深度土壤含水量平均年際變化率為0.001 m3·m-3·（10a）-1；降水量以10.0 mm·（10a）-1的速率呈不顯著增加趨勢；NDVI總體以0.007·（10a）-1的速率呈顯著增加趨勢；年平均雪深以-0.82 mm·（10a）-1的速率呈微弱減少趨勢。

（3）多年凍土區(qū)0.1 m 和0.4 m 土壤含水量、降水量、NDVI 及雪深的多年均值普遍大于季節(jié)凍土區(qū)，同樣影響因素在年際變化上也表現(xiàn)出不同的空間差異性，多年凍土區(qū)呈暖干化發(fā)展趨勢，季節(jié)凍土區(qū)呈暖濕化發(fā)展趨勢。

（4）地表凍融指數(shù)與影響因素以負(fù)相關(guān)關(guān)系為主，SFI 與影響因素在多年凍土區(qū)大部呈正相關(guān)關(guān)系，在季節(jié)凍土區(qū)大部呈負(fù)相關(guān)關(guān)系，STI反之。內(nèi)蒙古地表凍融指數(shù)變化受影響因素共同驅(qū)動，而0.4 m土壤含水量是影響SFI變化的主導(dǎo)因素，NDVI是影響STI變化的主導(dǎo)因素。