夏星旎, 趙 林, 汪凌霄, 周華云, 劉世博,李成燁, 趙建婷, 李智斌
(1. 南京信息工程大學(xué) 地理科學(xué)學(xué)院,江蘇 南京 210044; 2. 常熟市王淦昌高級中學(xué),江蘇 蘇州 215531; 3. 中國科學(xué)院 西北生態(tài)環(huán)境資源研究院 冰凍圈科學(xué)國家重點實驗室 藏北高原冰凍圈特殊環(huán)境與災(zāi)害國家野外科學(xué)觀測研究站,甘肅 蘭州 730000)
多年凍土,是指至少連續(xù)存在兩年且溫度在0 ℃及以下的巖土層,是巖石圈與大氣圈等圈層綜合作用的產(chǎn)物,大多分布在高海拔和高緯度嚴(yán)寒地帶[1]。我國的青藏高原多年凍土區(qū)是世界上中低緯度區(qū)域面積最大的高海拔多年凍土區(qū),面積達(dá)1.06×106km2,占青藏高原總面積的40%[2]。近年來青藏高原變暖趨勢明顯,多年凍土退化呈現(xiàn)出加速趨勢。多年凍土加速退化引發(fā)高寒生態(tài)系統(tǒng)退化[3]、生態(tài)系統(tǒng)碳釋放量增加[4],影響水循環(huán)[5-6],并給凍土區(qū)工程穩(wěn)定性帶來危害。
活動層位于多年凍土層之上,在青藏高原地區(qū)活動層厚度通常為2~3 m[7-8],夏季融化冬季凍結(jié)。活動層作為多年凍土區(qū)地氣間水熱交換的重要土層,在季節(jié)融化和凍結(jié)過程中,冰水之間發(fā)生的相變會導(dǎo)致地表出現(xiàn)季節(jié)性的凍脹和融沉變形。在氣候變暖的背景下,活動層厚度加深,多年凍土上限的地下冰逐漸融化,地下冰融化釋水,引起地表發(fā)生長期的沉降形變,并有部分水分參與到活動層凍融的年循環(huán)過程中[9-10]。多年凍土區(qū)的地面形變過程蘊含著多年凍土活動層內(nèi)部的變化狀況及地下冰變化信息,可將地面形變信息作為監(jiān)測多年凍土退化的“窗口”。
傳統(tǒng)的地表變形測量技術(shù)有水準(zhǔn)實測、埋設(shè)變形儀器測量和GPS 測量等,能夠?qū)θ舾呻x散點或小范圍區(qū)域進(jìn)行監(jiān)測,監(jiān)測精度高,但青藏高原平均海拔在4 000 m 以上,氣候條件惡劣,環(huán)境復(fù)雜,由于觀測條件限制,僅僅依靠人工和站點觀測完成大時空尺度的地表形變精確監(jiān)測極其困難[11]。隨著遙感技術(shù)的興起,利用合成孔徑雷達(dá)干涉測量技術(shù)(InSAR)對多年凍土區(qū)進(jìn)行全天候大范圍監(jiān)測已成為可能[12-13]。參與干涉的SAR 數(shù)據(jù)在形變發(fā)生前后分別獲取,提取干涉相位,利用DEM 消除地形相位后,將相位差轉(zhuǎn)為地表形變信息[14-18]。1999 年Wang 等[19-20]首次使用D-InSAR 技術(shù)應(yīng)用C 波段ERS-1 SAR數(shù)據(jù)在阿拉斯加Brooks Range地區(qū)監(jiān)測到了冬季地面凍脹和夏季融沉現(xiàn)象,開啟了使用D-InSAR 監(jiān)測多年凍土的地面形變研究。但是,傳統(tǒng)D-InSAR 測量會容易受到時空失相關(guān)和大氣延遲的影響,隨著哨兵1號(Sentinel-1)發(fā)射成功,SAR數(shù)據(jù)獲取數(shù)目增多,發(fā)展了時序InSAR(MT-In-SAR)技術(shù),通過統(tǒng)計模型和時空濾波的方法進(jìn)一步降低軌道誤差、DEM 誤差、大氣噪聲等的影響,主要為永久散射體干涉測量[19](PS-InSAR)和小基線集干涉測量[21](SBAS-InSAR)兩方面。在青藏高原自然地貌監(jiān)測中,SBAS-InSAR 相對于永久散射體干涉測量PS-InSAR 更有優(yōu)勢,因為PS-InSAR 依靠分布密集的永久散射體(如建筑物,青藏鐵路/公路的路基等),而SBAS-InSAR通過采用多幅主影像組合的方式,選擇多幅相干性高的干涉影像對,可以最大程度消除時間失相關(guān)與體散射失相關(guān)的影響[21]。不少學(xué)者利用Sentinel-1 數(shù)據(jù),采用SBAS-InSAR 技術(shù)對青藏高原腹地[22]、唐古拉南北坡[23]、兩道河地區(qū)[24]、五道梁地區(qū)[25]、鹽湖地區(qū)[26]、青藏鐵路沿線北麓河至風(fēng)火山南部至沱沱河[27]進(jìn)行地表形變監(jiān)測,并對地表形變特征進(jìn)行分析。隨著研究深入,已經(jīng)初步表明地表形變長期趨勢主要受地下冰變化影響[28],季節(jié)性形變主要受活動層內(nèi)水分含量控制[29-31]。
目前在青藏高原的凍土區(qū)地表形變監(jiān)測研究多集中在高原腹地沿青藏高原工程走廊沿線,對于多年凍土邊界地區(qū)研究很少。西大灘位于青藏高原多年凍土區(qū)北部邊界,是多年凍土和季節(jié)凍土交錯存在的過渡區(qū)。伴隨著高原氣候的暖濕化[3],該地區(qū)凍土的退化特征明顯,活動層增厚、多年凍土溫度升高,局部的島狀多年凍土甚至出現(xiàn)消失的跡象。因此,在青藏高原多年凍土區(qū)北界通過遙感手段獲取地表形變信息,綜合分析形變的空間分布特征,不僅可以探明氣候變暖條件下多年凍土的退化狀態(tài),也可為凍土分類及制圖提供支持。本文采用SBAS-InSAR 技術(shù)獲取2014—2020 年西大灘—昆侖埡口多年凍土區(qū)的地表形變時間序列結(jié)果,計算并分析了該地區(qū)凍土地表的長期和季節(jié)特征,并對西大灘多年凍土區(qū)北部邊界地區(qū)的形變特征做了細(xì)致分析。本研究首次嘗試了地表形變信息在高原多年凍土分布邊界判定中的應(yīng)用,并分析了西大灘地區(qū)與昆侖埡口地區(qū)不同地貌類型下多年凍土地表形變特征的差異,表明凍土地表形變的季節(jié)性形變量和長期趨勢可為多年凍土退化狀態(tài)及分布邊界判定提供新的視角和方法。
西大灘—昆侖山埡口地區(qū)位于青藏工程走廊沿線多年凍土北界,西大灘地處東昆侖山區(qū),是一個近東西走向的斷陷盆地,地勢自東向西升高,平均海拔4 480 m,谷底兩側(cè)的高山山谷中發(fā)育有山谷冰川、冰斗冰川和懸冰川[32-33],青藏鐵路自北向南經(jīng)過該區(qū)(圖1)。由于高空西風(fēng)帶的控制,氣候嚴(yán)寒干燥,年平均氣溫在-4.0 ℃以下,年降水量在400 mm左右,且集中在5—9月。區(qū)域植被條件較差,多年凍土和季節(jié)凍土并存[34]。多年凍土分布的下界海拔為4 360 m 左右[7],多年凍土溫度較高,年均地溫為-0.5~0.2 ℃,多年凍土厚度通常小于20 m[8],活動層厚度為1.5~4.5 m。隨著海拔向西逐漸升高,多年凍土溫度和活動層厚度呈下降趨勢,而多年凍土下限和多年凍土厚度呈增加趨勢[35]。研究區(qū)西南部為楚瑪爾河高平原,其熱融現(xiàn)象顯著,發(fā)育有大量熱融湖塘。
圖1 研究區(qū)示意圖(底圖為Sentinel-2影像真彩色合成數(shù)據(jù),多年凍土分布圖引自文獻(xiàn)[2])Fig. 1 Schematic diagram of the study area (The base image is Sentinel-2 true-color-composited image.The permafrost distribution is derived from Reference [2])
研究使用C 波段Sentinel-1 SAR 數(shù)據(jù),選取2014 年11 月3 日至2020 年12 月25 日第77 軌道共137 景Sentinel-1 降軌數(shù)據(jù),其中6 個為Sentinel-1B數(shù)據(jù),其余均為Sentinel-1A 數(shù)據(jù),以獲取地表形變信息。圖2 所示為技術(shù)路線圖,SBAS-InSAR 形變處理具體流程如下:
圖2 SBAS-InSAR多年凍土地表形變監(jiān)測技術(shù)路線圖Fig. 2 SBAS-InSAR workflow for surface deformation monitoring over permafrost terrain
(1)生成差分干涉相對。每景SAR 影像與其時間序列中的兩幅后續(xù)SAR 影像構(gòu)建干涉像對,例如時序上第n個影像與第n+1 和n+2 個影像分別生成干涉圖,對干涉圖進(jìn)行方位向:距離向3∶13 多視處理,以形成正方形像元同時降低噪聲。研究區(qū)入射角約為40°,多視后方位向格網(wǎng)大小約為42 m,地距向格網(wǎng)大小約為47 m。所有干涉圖的時間基線小于48天,空間基線小于300 m。去除SRTM DEM 模擬的地形相位后,得到差分干涉圖,對差分相位采用最小費用流算法(SNAPHU)進(jìn)行解纏,得到解纏后的差分相位,此步驟使用ISCE 軟件完成(https://github.com/isce-framework/isce2)。
(2)檢查干涉集,選擇相干性高與解纏錯誤少的差分干涉結(jié)果組成小基線干涉集,對于干涉效果差的像對予以刪除。
(3)選擇穩(wěn)定的(不受凍脹融沉影響的)裸露基巖作為空間上形變信息參考點,并利用閉合相位信息評估相位解纏結(jié)果并對解纏錯誤進(jìn)行修正。
(4)小基線集形變時間序列求解,使用解纏后的干涉圖相位方差倒數(shù)作為權(quán)重系數(shù),對干涉圖進(jìn)行加權(quán),給予高質(zhì)量(低方差)的干涉圖比低質(zhì)量(高方差)的干涉圖更大的權(quán)重,使用加權(quán)最小二乘法解算相位時間序列,從而更準(zhǔn)確地估計形變時間序列[27]。
(5)對使用大氣ERA5 再分析數(shù)據(jù)模擬干涉像對兩個SAR 影像由于對流層大氣情況不同造成的大氣延遲相位差別,進(jìn)行對流層大氣相位延遲改正[36]。
(6)大范圍的相位趨勢由電離層延遲和軌道誤差造成,通過二階多項式對相位趨勢進(jìn)行擬合,移除該相位梯度。
(7)地形殘差校正,消除由DEM 不準(zhǔn)確帶來的誤差,最終得到形變過程;以上小基線集形變映射采用時序InSAR 處理軟件MintPy 完成[27](https://github.com/insarlab/MintPy)。
InSAR 所獲取形變信息為相對于參考點形變,圖1 中白色圓點為選取參考點位置,經(jīng)緯度坐標(biāo):93°58′22″ E,35°38′24″ N。利用線性模型表示長期形變趨勢[式(1)],根據(jù)青藏高原凍土區(qū)凍融過程分析表明,凍脹值通常在1—2 月達(dá)到最大,融沉通常在8—10 月達(dá)到最大[29,37]。在減去長期形變趨勢后[式(2)],取1—2 月形變最大值及8—10 月形變最小值,計算兩者差值得到當(dāng)年季節(jié)性形變量,再通過研究時段內(nèi)多年度平均得到平均季節(jié)性形變量(振幅)。
式中:D(t)為對應(yīng)t時刻(相對起始日期的累積天數(shù))的累積形變值;v為長期形變速率;D'(t)為減去長期形變趨勢后的形變時間序列;max(D'(Jan,F(xiàn)eb))和min(D'(Aug,Sep,Oct))表示每年1—2 月地形高點,8—10月地形低點;S為減去長期趨勢后的季節(jié)性形變量,在空間網(wǎng)格中逐像素計算。
DEM數(shù)據(jù)、熱融湖塘數(shù)據(jù)、中國1∶400萬數(shù)字地貌數(shù)據(jù)和高分辨率國家土壤信息格網(wǎng)屬性數(shù)據(jù)用于分析形變的空間分布特征。30 m SRTM DEM 數(shù)據(jù)從美國USGS 官網(wǎng)(https://lpdaac. usgs. gov)獲取。青藏高原熱融湖塘數(shù)據(jù)[38-40]來源于國家青藏高原科學(xué)數(shù)據(jù)中心(http://data.tpdc.ac.cn)。中國1∶400萬數(shù)字地貌數(shù)據(jù)集[41]來源于時空三極環(huán)境大數(shù)據(jù)平臺(http://poles.tpdc.ac.cn)。中國高分辨率國家土壤信息格網(wǎng)基本屬性數(shù)據(jù)集(2010—2018)數(shù)據(jù)來源于國家科技基礎(chǔ)條件平臺—國家地球系統(tǒng)科學(xué)數(shù)據(jù)中心—土壤分中心(http://soil.geodata.cn)[42]。
考慮到不同土層深度的土壤質(zhì)地存在差異,無法代表整體的土壤質(zhì)地情況,因此選取250 m 分辨率的砂粒、粉粒、黏粒數(shù)據(jù),將其0~5 cm、5~15 cm、15~30 cm、30~60 cm 和60~100 cm 的柵格數(shù)據(jù)分別乘上對應(yīng)的土層深度作為權(quán)重,計算得到0~100 cm深度內(nèi)的土壤質(zhì)地信息。
圖3 展示了2014 年11 月3 日至2020 年12 月25日研究期間內(nèi)24 個時段地表累積形變(以2014 年11 月3 日第一幅影像形變?yōu)?)的空間分布,涵蓋了不同年份凍結(jié)季與融化季的形變結(jié)果。圖3中地表抬升用冷色調(diào)正值表示,地表沉降用暖色調(diào)負(fù)值表示,季節(jié)性凍脹融沉形變及長期沉降趨勢明顯,并且呈現(xiàn)出很強的空間異質(zhì)性。例如在2015-02-07、2015-05-02、2016-02-14 與2016-04-26 等日期,能觀察到明顯的地表凍脹現(xiàn)象,而伴隨著地表長期沉降,在2018年之后,冬季地表的累積形變量為負(fù)值,例如在研究時段末尾日期2020-12-25相對于起始時段沉降明顯,其中南部楚瑪爾河高平原地表沉降最為明顯,累積沉降量在70 mm 以上。圖4 為在形變時間序列基礎(chǔ)上提取的長期形變趨勢以及季節(jié)性形變量,圖3~4中白色顯示失相干地區(qū)。
圖3 研究區(qū)時間序列形變(參考2014年11月3日)Fig. 3 Deformation maps in the study area in different acquisition dates, all with reference to November 3, 2014
圖4 研究區(qū)地表形變空間分布圖Fig. 4 Spatial distribution of surface deformation in the study area: spatial distribution of long-term deformation velocities (a); spatial distribution of seasonal deformation (b)
圖4東南和西南部(楚瑪爾河高平原)地表沉降最為明顯,長期沉降速率值大于15 mm·a-1,在沉降明顯的地區(qū)季節(jié)性形變也相應(yīng)較大,季節(jié)性形變量大于12 mm。圖4中部(昆侖埡口)西側(cè)(圖1),地表長期沉降明顯,但是季節(jié)性形變量不大。在北部格爾木河河谷,地表略微抬升,可能沉積作用也可能是土壤水分變化對InSAR 結(jié)果造成影響。西大灘谷地長期形變速率多位于-5~2.5 mm·a-1之間,西大灘東部位于河流下游,長期形變速率在0~5 mm·a-1之間。昆侖山埡口的地表主要表現(xiàn)為沉降趨勢,形變速率主要在-15~0 mm·a-1之間,位于埡口地區(qū)的青藏鐵路附近沉降速率大于5 mm·a-1。2015—2020年研究區(qū)季節(jié)性形變量基本位于0~40 mm 之間,空間分布差異較大[圖4(b)]。在昆侖山脈北部和西大灘溝谷地區(qū),季節(jié)性形變量多在0~10 mm 之間。在昆侖山埡口地區(qū)和楚瑪爾河高平原以北地區(qū),活動層含水量更高,季節(jié)性形變也更顯著,形變量大于10 mm。
西大灘谷地位于多年凍土北界,多年凍土與季節(jié)凍土共同存在,多年凍土呈島狀分布且厚度較薄。多年凍土特性隨海拔變化表現(xiàn)出明顯的差異,且所處退化歷史階段不同。為探究多年凍土區(qū)北界西大灘谷地地表形變的空間差異,沿青藏鐵路在西大灘溝谷自西向東繪制橫向剖線HH′,并根據(jù)西大灘多年凍土分布邊界探測剖面數(shù)據(jù),自北向南繪制縱向剖線AA′、BB′、CC′、DD′、EE′,剖線位置如圖5 所示。HH′剖線形變結(jié)果如圖5(e)~(f)所示,南北向剖線處長期形變速率如圖6(b)~(f)所示。為了便于對形變結(jié)果進(jìn)行解譯,圖5(b)和圖6(a)分別展示了西大灘谷地模型模擬的凍土類型和10 m深度處2015—2019年期間平均地溫,數(shù)據(jù)來自Zhao等[35]近期發(fā)布的模擬結(jié)果。西大灘北坡和南坡連續(xù)多年凍土下限模擬分別為約4 525 m 和4 732 m。多年凍土面積約占西大灘地區(qū)的80%(其中,33.93%為連續(xù)多年凍土,46.07%為不連續(xù)多年凍土)。圖5(b)中季節(jié)凍土I 表示該區(qū)1970 年代已為季節(jié)凍土,季節(jié)凍土II 表示該區(qū)1970 年代為不連續(xù)多年凍土,但是2009 年已退化為季節(jié)凍土,其地溫要比季節(jié)凍土I低。
圖5 西大灘溝谷地區(qū)東西向與南北向剖線位置及地表形變情況Fig. 5 The position of the east-west and north-south profile lines and surface deformation in the Xidatan valley area: Sentinel-2 image (a); type of frozen soil (b); long-term deformation velocity (c); seasonal deformation (d); variation of long-term deformation rate with altitude (e); variation of seasonal deformation with altitude (f)
圖6 西大灘溝谷地區(qū)南北向剖線的形變速率Fig. 6 The deformation velocities of the south-north profiles in the Xidatan valley area: average ground temperature at 10 m depth of regional permafrost during 2015—2020 (a); AA′ profile (b); BB′ profile (c); CC′ profile (d); DD′ profile (e); EE′ profile(f); (The starting point of the profile corresponding to A, B, C, D, E is displayed on the left side of the x-axis, and the endpoint of the profile corresponding to A′, B′, C′, D′, E′ is displayed on the right side of the x-axis)
圖5展示東西向HH′剖線長期形變速率和季節(jié)性形變量隨海拔變化情況。HH′剖線海拔自西向東從4 580 m 遞減至4 100 m,坡度均小于5.6°,根據(jù)長期形變速率和季節(jié)性形變速率,可將剖線大致分為三段,在圖5中以紅色線段示意。
第①段多年凍土區(qū):在HH′剖線海拔高于4 520 m 處(HH′剖線西側(cè)),地表長期形變速率多為0~1 mm·a-1,季節(jié)性形變量小于6 mm。從圖5(a)影像圖中可以看出該段為河漫灘,在風(fēng)化作用下,山頂基巖風(fēng)化形成的碎屑物質(zhì)在坡腳處堆積可能使得地表輕微上升,流水的沉積作用也會導(dǎo)致地表輕微抬升,使得該地區(qū)在剖線的三段中沉降量最小,且有微弱的抬升信號。
第②段不連續(xù)多年凍土區(qū):HH′剖線海拔位于4 370~4 520 m 的地區(qū),地表長期形變速率在-4 mm·a-1附近,表現(xiàn)出明顯的沉降趨勢,季節(jié)性形變量多高于6 mm,且波動較大,空間差異明顯。
第③段季節(jié)凍土區(qū):HH′剖線海拔低于4 370 m的地區(qū),地表長期形變速率在-3~-1 mm·a-1之間,季節(jié)性形變量多為6.5 mm左右,波動明顯小于第二段。
紅色實線為根據(jù)形變信息估計的多年凍土和季節(jié)凍土的分界線,紅線左側(cè)(第①段和第②段)為多年凍土分布區(qū),多年凍土存在及退化使得該區(qū)段具有更大的季節(jié)性形變量同時表現(xiàn)出更高的沉降速率。根據(jù)長期形變速率和季節(jié)性形變速率將剖線大致分為三段,與區(qū)域模擬的多年凍土分布吻合。其中第①段對應(yīng)圖5(b)中所示連續(xù)多年凍土,第②段對應(yīng)非連續(xù)多年凍土,第③段對應(yīng)季節(jié)凍土區(qū)。
南北向剖線AA′(海拔4 536~4 617 m),BB′(海拔4 498~4 553 m),CC′(海拔4 427~4 545 m),DD′(海拔4 379~4 498m)和EE′(海拔4 317~4 434 m)的形變細(xì)節(jié)如圖7 所示,A、B、C、D、E位于山體南坡,A′、B′、C′、D′、E′位于北坡,山體南坡溫度普遍高于山體北坡溫度。根據(jù)該區(qū)域凍土地溫的模擬結(jié)果,AA′、BB′、CC′、DD′和EE′剖線2015—2020期間10 m 深度處的平均地溫分別為-1.1~-0.4 ℃、-1.4~0.6 ℃、-2.9~0.2 ℃、-1.3~0.8 ℃、-1.3~2.7 ℃。圖6(a)剖線AA′位于非連續(xù)多年凍土至連續(xù)多年凍土區(qū),是5條剖線中溫度最低的一條,形變速率大多位于-3~1 mm·a-1區(qū)間。圖6(b)剖線BB′位于季節(jié)凍土區(qū)-非連續(xù)多年凍土-連續(xù)多年凍土區(qū),長期形變速率由季節(jié)凍土區(qū)的0 mm·a-1附近增大至多年凍土區(qū)的-7 mm·a-1。圖6(c)剖線CC′位于非連續(xù)多年凍土-連續(xù)多年凍土區(qū),形變速率在-3 mm·a-1左右,由于地溫差異不大,其形變速率差異也相應(yīng)較小,其中南端C′處位于山坡處,形變速率出現(xiàn)正值由于雷達(dá)視線對地體斜坡方向位移不敏感出現(xiàn)誤差。圖6(d)剖線DD′和EE′都位于季節(jié)凍土區(qū)-非連續(xù)多年凍土區(qū),形變速率在-7~-2 mm·a-1之間。BB′、DD′和EE′三條剖線由季節(jié)凍土區(qū)向非連續(xù)多年凍土區(qū)過渡,長期沉降速率也隨之增大。AA′、CC′兩條剖線海拔較高、多年凍土地溫較低,該處多年凍土退化以升溫為主,地表沉降不明顯,對應(yīng)地下冰融化不顯著[43-44]。BB′、DD′和EE′剖線海拔較低,剖線南端屬于高溫多年凍土區(qū),地表沉降速率明顯高于海拔較高的AA′、CC′剖線,對應(yīng)較為明顯的地下冰的融化[28,43]。
圖7 昆侖山埡口—楚瑪爾河高平原地區(qū)地貌土壤類型及地表形變情況Fig. 7 Types of geomorphic soil and surface deformation in Kunlun Pass-Chumar River High Plain: the location and natural surface condition of points P1, P2, P3, and P4 (a); morphological and geomorphological types in the study area (b);sand content, 1, 2, 3, 4 in the figure indicates the selected sectors (c); silt content (d);long-term deformation velocity (e); seasonal deformation (f)
基于東西向HH′與南北向AA′、BB′、CC′、DD′和EE′共6 條剖線的分析結(jié)果表明,在多年凍土邊界地區(qū),地表相近情況下,多年凍土區(qū)的長期沉降速率和季節(jié)性的形變量比季節(jié)凍土區(qū)要高,在不連續(xù)多年凍土區(qū)形變空間差異較大,此外,高溫多年凍土地表沉降相比低溫多年凍土更為顯著。由此分析,地表形變特征差異不僅可以被用來評估多年凍土退化狀態(tài),也可作為一種數(shù)據(jù)來協(xié)同輔助判定多年凍土分布邊界。
圖7展示了昆侖山埡口—楚瑪爾河高平原地區(qū)不同地貌土壤類型下的形變空間分布特征。與西大灘谷地相比,該區(qū)域長期形變速率與季節(jié)性形變量都明顯增大。結(jié)合圖7(b)~(d)表明不同地貌單元與土壤質(zhì)地下地表長期形變速率和季節(jié)性形變量存在明顯的差異。對不同地貌單元下的形變特征進(jìn)行統(tǒng)計如圖8 所示。統(tǒng)計表明,類型中1-高海拔丘陵、5-高海拔沖洪積湖積臺地、9-高海拔洪積平原、10-高海拔冰水沉積平原表現(xiàn)出顯著的地表沉降趨勢和高季節(jié)性形變量,平均沉降速率高于5 mm·a-1,平均季節(jié)性形變量高于10 mm,其中類型5-高海拔沖洪積湖積臺地的季節(jié)性形變量最大,平均值高達(dá)20 mm。高海拔丘陵在所有類型中長期沉降量最大,沉降速率平均值接近10 mm·a-1,而平均季節(jié)性形變量約10 mm,小于平原地區(qū)。類型2與3即海拔中起伏山地與小起伏山地的長期形變量和季節(jié)性形變量都偏小,長期形變速率趨近于0,季節(jié)性形變量通常小于6 mm。
圖8 不同地貌單元下的形變特征箱線圖Fig. 8 Box plots of deformation characteristics under different geomorphic units: box plots of long-term deformation velocity under different geomorphic units (a); box plots of seasonal deformation under different geomorphic units (b)
從圖7(a)和圖1 中的Sentinel-2 影像上可知該區(qū)域東南部楚瑪爾河高平原上有大量熱融湖塘分布。在高含冰量的多年凍土區(qū),在平坦低洼處(坡度通常小于3°),富含冰的多年凍土層發(fā)生退化或局部地下冰融化,水分在地勢低洼處匯聚,形成積水洼地并且不斷積聚和擴(kuò)張,地表有可能形成熱融湖塘和熱融洼地[45-47]。隨著原有地表水熱狀況的失衡,多年凍土不斷融化所產(chǎn)生的水分匯集起來,進(jìn)而加速熱融湖塘的擴(kuò)張[48-49]。等到湖塘底部的多年凍土徹底融化,凍結(jié)層內(nèi)部水相互貫通,形成貫通融區(qū)[9],非貫穿型湖塘轉(zhuǎn)變?yōu)樨灤┬秃?,新的多年凍土可能會在湖盆底部形成?0]。此地區(qū)內(nèi)眾多熱融湖的存在是多年凍土退化的直接證據(jù),此前研究也表明有熱融湖塘發(fā)育的地區(qū)通常能觀察到明顯的地表沉降[51]。為進(jìn)一步探明熱融湖塘地貌與形變特征的聯(lián)系,在區(qū)域內(nèi)選取了4 個斑塊[圖7(c)]進(jìn)行分析,這4 個斑塊在土壤質(zhì)地圖和形變圖中都有明顯的邊界,且土壤質(zhì)地、熱融湖塘特征、形變特征方面各有差異,其中圖斑2 和圖斑3 土壤質(zhì)地相近,圖斑1和圖斑4土壤質(zhì)地相近,但圖斑1和圖斑2中的熱融湖塘面積明顯小于圖斑3和圖斑4,且在特性特征中也各有差異。4 個圖斑的湖塘面積、分布、土壤質(zhì)地及形變特征信息總結(jié)見表1。圖7 可以看出,所選4 個圖斑的季節(jié)性形變量和長期沉降量都高于區(qū)域整體,即熱融湖塘發(fā)育使得區(qū)域季節(jié)性形變量和長期沉降速率增大。其中圖斑1 和圖斑2 內(nèi)的湖塘更小、分布密度更高,圖斑3 和圖斑4 內(nèi)的湖塘面積明顯增大,與此同時圖斑3和圖斑4的長期沉降趨勢比圖斑1 和圖斑2 偏大,而季節(jié)性形變量偏小。圖斑2與圖斑3的土壤質(zhì)地相近,圖斑3中的湖塘明顯比圖斑2更大,圖斑3的季節(jié)性形變量比圖斑2偏小且沉降量更大,平均沉降速率大于10 mm·a-1。
表1 斑塊內(nèi)湖塘、土壤質(zhì)地與形變信息Table 1 Information statistics of ponds, soil texture and deformation in patch
圖9 展示P1、P2、P3、P4 四點形變時間序列。P1、P2 點季節(jié)性形變量約為15 mm,長期沉降速率2~3 mm·a-1。P3 點與P1、P2 相比,沉降趨勢更加顯著,所處地理位置靠近楚瑪爾河,活動層含水量相對更高,因此季節(jié)性幅度更大,在2017 年之后季節(jié)性形變量顯著增大,達(dá)到25 mm 以上。P4點同樣靠近楚瑪爾河,長期沉降速率比P3 更大,其季節(jié)性形變量在2015—2017 年間增大,但2017 年之后卻輕微減小。
圖9 P1 (a),P2 (b),P3 (c),P4 (d)四點形變時間序列(e)Fig. 9 Deformation time series of P1 (a), P2 (b), P3 (c), P4 (d) (from November 3, 2014 to December 25, 2020) (e)
基于以上分析推斷,在熱融湖塘發(fā)展的早期湖塘面積較小,但地下冰融化使得區(qū)域季節(jié)性形變量增大,隨著熱融湖塘擴(kuò)張,使得區(qū)域長期沉降速率加劇,熱融湖塘進(jìn)一步發(fā)展后,伴隨著非貫穿型湖塘轉(zhuǎn)變?yōu)樨灤┬秃梁螅瑓^(qū)域季節(jié)性形變量可能降低。
研究采用SBAS-InSAR 技術(shù)對西大灘—昆侖山埡口地區(qū)進(jìn)行長期連續(xù)的形變監(jiān)測,得到了研究區(qū)2014—2020 年形變序列,并且計算獲得多年凍土長期形變速率和其季節(jié)性形變量結(jié)果,發(fā)現(xiàn)研究區(qū)長期形變速率基本位于-30~10 mm·a-1之間,形變速率大小的空間分布差異較大,地表總體呈現(xiàn)出下降的趨勢。
在多年凍土北界西大灘溝谷地區(qū),地表相近情況下,多年凍土區(qū)的長期沉降速率和季節(jié)性的形變量比季節(jié)凍土區(qū)要高,在不連續(xù)多年凍土區(qū)形變空間差異較大。此外,高溫多年凍土地表沉降相比低溫多年凍土更為顯著。由此分析,地表形變特征差異不僅可以被用來評估多年凍土退化狀態(tài),也可作為一種數(shù)據(jù)協(xié)同輔助判定多年凍土分布邊界。
與西大灘谷地相比,昆侖山埡口地區(qū)和楚瑪爾河高平原區(qū)域的長期形變速率與季節(jié)性形變量都明顯增大,且區(qū)域內(nèi)平原地區(qū)的沉降速率和季節(jié)性形變量通常高于起伏山地。結(jié)合土壤地貌信息和熱融湖塘數(shù)據(jù)進(jìn)一步研究形變規(guī)律表明,熱融湖塘發(fā)展的早期面積較小,但地下冰融化使得區(qū)域季節(jié)性形變量增大,隨著熱融湖塘擴(kuò)張,使得區(qū)域長期沉降速率加劇,熱融湖塘進(jìn)一步發(fā)展后,伴隨著非貫穿型湖塘轉(zhuǎn)變?yōu)樨灤┬秃梁?,區(qū)域季節(jié)性形變量可能降低。