基于MODIS數(shù)據(jù)的額爾齊斯河上游河冰物候研究

胡婉嬪， 效存德， 謝愛紅， 沈永平， 雷華錦

（1.中國科學院西北生態(tài)環(huán)境資源研究院冰凍圈科學國家重點實驗室，甘肅蘭州730000；2.中國科學院大學，北京100049；3.北京師范大學地表過程與資源生態(tài)國家重點實驗室，北京100875）

0 引言

額爾齊斯河是發(fā)源于中國新疆阿爾泰山南坡的一條跨國河流，沿阿爾泰山南麓向西北，經(jīng)中國、哈薩克斯坦和俄羅斯，在俄羅斯境內(nèi)流入鄂畢河，后注入北冰洋，是中國唯一流入北冰洋的河流［1-6］。額爾齊斯河干流屬于荒原河流，水面寬闊，流速平緩，落差小，能行舟航運［7］。20世紀50至60年代，額爾齊斯河在中蘇之間運輸建筑材料、石油制品等貨物，貨運量從630萬噸增長至1 850萬噸，后因中蘇關系惡化，該河航運停止［8］。目前隨著全球氣候變暖，北極變暖的速度遠快于預期，北冰洋冰層融化速度加快，極可能在2050年前后出現(xiàn)夏季無冰狀態(tài)［9］，北極航道成為各國共同關注的熱點，隨著中國“一帶一路”倡議的穩(wěn)步推進，“冰上絲綢之路”逐漸走進大眾的視野。據(jù)報道，設濱海國際運輸走廊的邀請，建議雙方共同開發(fā)和利用海上通道特別是北極航道，打造“冰上絲綢之路”。2018年1月，中國政府發(fā)布《中國的北極政策》白皮書，提出中俄共建“冰上絲綢之路”，促進北極地區(qū)互聯(lián)互通和經(jīng)濟社會可持續(xù)發(fā)展，至此“冰上絲綢之路”得到完整闡釋［10-12］。有報道稱哈薩克斯坦和俄羅斯政府早在2004年已簽訂相關協(xié)議，清理并拓深哈國邊境至鄂木斯克市的河道，使額爾齊斯河實現(xiàn)哈俄之間的跨境通航［13］。若中國與哈薩克斯坦之間的額爾齊斯河上游能夠通航，中國貨物可以從新疆布爾津縣經(jīng)額爾齊斯河直接運往哈薩克斯坦后，到達北冰洋，運往歐美等國家。據(jù)哈薩克斯坦估計，若額爾齊斯河上游通航，不僅可以開發(fā)沿途的旅游業(yè)，而且水路貨運量全年可達120萬噸［14］。

北半球許多季節(jié)性封凍的河流和湖泊對航運和運輸業(yè)有很大影響，河湖中浮冰密度太大，冰塊太厚或嚴重封凍，將導致每年內(nèi)陸水道停運數(shù)月［15］。額爾齊斯河為季節(jié)性封凍河流，每年10月開始出現(xiàn)岸冰，12月封凍，至次年5月初完全開河，約有120天至160天封凍期。隨著“冰上絲綢之路”的建設，額爾齊斯河冰凍期及其變化研究顯得尤為必要。

當前，學者們基于遙感數(shù)據(jù)對河冰及其變化過程進行了大量研究，Chaouch等［16］結合一種基于閾值的決策樹圖像分類算法，利用MODIS圖像探測和繪制薩斯奎哈納河河冰，并確定河冰范圍；Chu等［17］使用MOD09GQ第二波段數(shù)據(jù)研究加拿大西北地區(qū)Slave River河沿岸的河流結冰過程；Unterschultz等［18］利 用RADARSAT-1合 成 孔 徑 雷 達（SAR）衛(wèi)星影像描述河流冰在冬季與崩解期間的特性的可行性。此外，在河湖冰信息提取方面，Sarah等［19］利用MODIS數(shù)據(jù)研究了河冰開河時間及開河的影響因素，Cai等［20］基于開放水域像素數(shù)的閾值提取青藏高原58個湖泊的凍結開始和破裂結束時間。

本文擬采用TIMESAT軟件，提取河冰物候信息，并與南灣水文站實測資料、布爾津站與哈巴河站氣象觀測資料進行對比分析，進而對河冰開封河時間、冰期持續(xù)時長及開封河速率等進行分析，以期為關于本河段通航所帶來的經(jīng)濟效益研究奠定科學基礎，為評估河段航運現(xiàn)實能力及未來潛力提供科學依據(jù)。由于資料有限，目前僅研究額爾齊斯河布爾津至齋桑泊河段河冰物候。

1 數(shù)據(jù)源與研究方法

1.1 研究區(qū)概況

額爾齊斯河全長4 248 km，在中國境內(nèi)546 km，流域面積5.7×104km2，年徑流量多達111×108m3。本文研究區(qū)為額爾齊斯河流域上游，范圍為從布爾津到齋桑泊的額爾齊斯河干流河段（47°40′25.88″～47°48′37.16″N，86°54′40.24″～84°31′12.73″E）（圖1），全長約194 km，其中中國境內(nèi)約115 km，哈薩克斯坦境內(nèi)約79 km，河流最寬處約350 m，最窄處約100 m。研究區(qū)在中國境內(nèi)屬于中國新疆阿勒泰平原地區(qū)，地處歐亞大陸中心腹地，緯度較高，屬于溫帶大陸性氣候，冬季漫長而寒冷，夏季短促、氣溫平和，最冷月通常在1月，平均氣溫為-16℃，最熱月在7月，平均氣溫21℃，地區(qū)降水量少，年均降水量131～223 mm，蒸發(fā)量大，晝夜溫差大，光照充足［21］。

圖1 額爾齊斯河研究區(qū)Fig.1 Irtysh river research area

1.2 數(shù)據(jù)集及預處理

MOD09GQ遙感數(shù)據(jù)源于美國國家航空航天局（NASA）網(wǎng) 站（https：//search.earthdata.nasa.gov/search），是MODIS的反射率產(chǎn)品，時間分辨率為1天，空間分辨率為250 m，提供了兩個波段的數(shù)據(jù)，分別是：第一波段紅光波段（620～670 nm）和第二波段近紅外波段（841～876 nm）（表1）［22］。研究表明近紅外波段對于區(qū)別水和雪/冰覆蓋最敏感，與紅光波段相比，近紅外波段更少受到大氣條件的影響［23］，故選取第二波段數(shù)據(jù)作為研究數(shù)據(jù)，其時間范圍為2001—2018年的每年9月1日至5月31日。

表1 MOD09GQ數(shù)據(jù)說明Table 1 Data specification of MOD09GQ

本文的氣象數(shù)據(jù)選用布爾津站（47°25′12″N，86°31′11″E）及哈巴河站（48°01′48″N，86°14′23″E）2001—2017年的逐日平均氣溫，用于分析封河期逐日平均氣溫的變化情況，及其與封河速率的關系。水文數(shù)據(jù)選用南灣水文站（48°00′N，85°41′E）2001—2018年的河冰開封河時間，用于對比驗證遙感數(shù)據(jù)所得結果。

研究區(qū)影像通過谷歌地球（Google Earth）和ArcGIS軟件處理得到，使用Google Earth手動矢量化研究區(qū)范圍（額爾齊斯河布爾津至齋桑泊段）得到研究區(qū)矢量圖層，然后對MOD09GQ數(shù)據(jù)進行裁剪，最后提取研究區(qū)河流的平均反射率值進行后續(xù)研究。

1.3 研究方法

1.3.1 研究路線

本研究將使用HANTS算法去除云噪聲，用TIMESAT軟件得到河冰物候，具體流程如圖3所示。

圖3 MODIS的MOD09GQ產(chǎn)品對河冰物候分析總流程Fig.3 General process for river ice phenology analysis of MODIS MOD09GQ

1.3.2 河冰物候定義

河冰物候研究涉及河冰周期性形成與消失，河冰物候受到氣候變化以及人類的影響，這與植被物候周期性變化相似，植被物候受季節(jié)氣候和年際變化的影響，因此植被物候的研究方法亦可用于研究河冰物候。當氣溫降至0℃時，河水開始凍結，并歷經(jīng)水內(nèi)冰、薄冰、岸冰、冰覆蓋和封凍等階段［24］。在河流開始結冰之前，MODIS影像的地表反射率有最低值，隨著河水開始凍結，反射率逐漸增大。由于光學圖像的反射率值取決于地表物體的表面屬性而不是內(nèi)部結構，同樣的，河冰反射率僅受河冰表面積雪等的影響，因此，本研究認為當反射率達到最大值時是冰與雪的最大覆蓋情況［17］。隨著早春氣溫升高，以及上游開河后流向下游的水流量增加，研究區(qū)河段開始開河，直至河流完全無冰，結束開河過程。本文定義通過MODIS探測到河冰出現(xiàn)時為開始封河，當河流無冰時為結束開河，從開始封河到河面最大冰/雪覆蓋所經(jīng)歷的時間為封河期，從河面有最大冰/雪覆蓋到結束開河所經(jīng)歷的時間為開河期（圖2）。

圖2 河冰物候定義Fig.2 Phenological definition of river ice

1.3.3 時間序列諧波分析法（HANTS）

MODIS數(shù)據(jù)具有時間分辨率高、空間覆蓋范圍廣等優(yōu)點，但其影像常被大量云覆蓋，尤其在冬季，部分像元長時間（如10天以上）連續(xù)被云遮擋，因此需要采用濾波技術降低噪聲，填補MODIS時間序列數(shù)據(jù)的空白。

時間序列的諧波分析法（HANTS算法）［25-26］是為了去除歸一化植被指數(shù)（NDVI）時間序列中的云污染值而開發(fā)的。HANTS算法的基本思想是計算一個傅立葉級數(shù)來對像素級觀測的時間序列建模，同時識別相對于時間序列模型的異常值。該算法消除了這些異常值，并用傅里葉級數(shù)給出的值將其替換。與快速傅里葉變換（FFT）算法相比，HANTS算法在頻率選擇和時間序列長度上更具靈活性，此外，由于樣本在時間上不要求等距，從時間序列中排除某些點更容易。通過去除時間序列中明顯的異常值，提取的諧波比直接的FFT算法可靠得多［27］。HANTS算法已經(jīng)成功的應用于重建無云的遙感植被數(shù)據(jù)集（例如歸一化植被指數(shù)（NDVI）、增強型植被指數(shù)（EVI）和葉面積指數(shù)（LAI））、重建地表溫度數(shù)據(jù)集以及重建河冰反射率時間序列［27-33，17］。因此，本文采用HANTS算法消除云及其他噪聲影響，填補MODIS時間序列數(shù)據(jù)的空白，重建MODIS時間序列。

1.3.4 TIMESAT

雖然遙感時間序列數(shù)據(jù)在監(jiān)測植被季節(jié)方面的價值已經(jīng)得到認可，但從這些序列中探索和提取季節(jié)參數(shù)的方法有限，因此，TIMESAT軟件是為了提取季節(jié)參數(shù)而開發(fā)［34］?，F(xiàn)已應用于諸多研究，例如表征物候學［35］、促進生態(tài)系統(tǒng)中的數(shù)據(jù)分類［36］以及監(jiān)測在火災季節(jié)時人類的足跡［37］等。

該軟件實現(xiàn)了三種基于最小二乘擬合的NDVI數(shù)據(jù)處理方法。第一種方法采用局部多項式函數(shù)進行擬合，該方法可分類為自適應Savitzky Golay濾波器（簡稱S-G濾波）。另外兩種方法是普通的最小二乘法，分別為非對稱高斯函數(shù)和雙邏輯函數(shù)，其數(shù)據(jù)適合于不同復雜度的模型函數(shù)。這三種處理方法都使用了季節(jié)性（單模態(tài)或雙模態(tài)）的初步定義，以及生長季節(jié)的近似時間［38-39］。本文使用河冰反射率數(shù)據(jù)代替NDVI數(shù)據(jù)，所得生長季節(jié)的近似時間即為開封河時間及最大冰/雪覆蓋時間。

2 結果與分析

2.1 濾波提取結果對比

S-G濾波器及高斯濾波器提取的河冰物候信息如圖4所示，通過比較分析兩個時間序列結果，高斯擬合所得時間序列更符合研究所需，最終選用其所得河冰物候信息。由于MODIS部分數(shù)據(jù)被云遮蓋（尤其是冬季影像），所得反射率值高于地物反射率，經(jīng)過HANTS算法去云以及TIMESAT軟件濾波之后，MOD09GQ數(shù)據(jù)中因為云及其他噪聲引起的部分高反射率值得到糾正（圖4）。

圖4 由MOD09GQ第二波段得到的額爾齊斯河布爾津至齋桑泊段河冰物候Fig.4 River ice phenology of Irtysh river from Burqin to Lake Zaysan section was obtained by MOD09GQ band 2

2.2 河冰物候信息

表2 列出了2001—2018年每年的開始封河時間、結束封河時間、最大冰/雪覆蓋時間、最大冰/雪覆蓋時所對應的MOD09GQ第二波段的反射率、冰期持續(xù)時間、封河速率以及開河速率，其中開封河速率計算方法是開河（或封河）時所對應的河冰反射率與最大冰/雪覆蓋時河冰反射率的差值，除以開河期（或封河期）時長。

2.3 結果驗證

根據(jù)表2得到冰/雪最大覆蓋日期，結合預覽高分1號影像上研究區(qū)的河冰變化情況，選取了2015年的2月2日、3月15日、4月8日以及4月12日的MOD09GQ第二波段的反射率數(shù)據(jù)進行對比分析。2月2日為2015年冰/雪最大覆蓋日期，預覽高分1號影像上顯示4月8日本文研究區(qū)上游有部分冰，而下游齋桑湖附近冰較多，4月12日研究區(qū)上游無冰，下游齋桑湖附近無影像。如圖5所示，2月2日的反射率值總體最大，之后3月15日、4月8日、4月12日的反射率值明顯降低。這表明，隨著河冰開河，MOD09GQ的第二波段反射率值逐漸降低，從而認為MOD09GQ的第二波段數(shù)據(jù)能有效監(jiān)測河冰變化。

圖5 2015年冰/雪最大覆蓋期到開河結束期間4天的MOD09GQ第二波段反射率值的變化情況Fig.5 The reflectivity of MOD09GQ second band is worth changing for 4 days from the maximum ice/snow cover period in 2015 to the end of break-up

表2 根據(jù)MOD09GQ第2波段數(shù)據(jù)（841～876 nm）得到的額爾齊斯河布爾津至齋桑泊段2001—2018年間河冰物候信息Table 2 According to MOD09GQ band 2 data（841～876 nm），the river ice phenology information of the Irtysh River from Burqin to Lake Zaysan section from 2001 to 2018 was obtained

由于研究區(qū)水文站較少，并且水文站監(jiān)測河流開封河時間定義與本文不同，因此本文選取額爾齊斯河上，位于布爾津與齋桑泊中間的南灣水文站數(shù)據(jù)進行對比驗證，分別分析由MOD09GQ數(shù)據(jù)及南灣站得到的2001—2018年間冰期持續(xù)時間的變化趨勢。由于本文定義MODIS能夠探測到河冰出現(xiàn)時為開始凍結，到河流無冰時為開河結束，而水文站定義在封凍期和開河期，當河冰達到一定百分比則認為開始封凍和結束開河，因此本研究所得冰凍期歷時大于南灣站數(shù)據(jù)所得冰凍期歷時。2001—2018年間，由MOD09GQ數(shù)據(jù)及南灣水文站所得冰期持續(xù)時間均略有變長趨勢［圖6（a）、6（b）］，并且通過線性分析得到其相關系數(shù)R=0.306，呈正相關［圖6（c）］，說明本文所得結果與南灣站實測數(shù)據(jù)結果的變化趨勢一致，均呈現(xiàn)冰期歷時延長的趨勢。并且水文站及本文研究結果中各年開封河時間變化趨勢具有一致性，如2013年南灣站、布爾津站及本文結果解凍日期同比2012年分別提早10天、21天及10天，封凍日期同比2012年分別晚30天、21天及31天，趨勢一致，天數(shù)相近。

圖6 根據(jù)MODIS數(shù)據(jù)所得冰期持續(xù)時間與南灣站數(shù)據(jù)所得冰期持續(xù)時間分別的變化趨勢及二者的相關關系Fig.6 based on MODIS data from the glacial ice duration and Nanwan station data duration，respectively，the change trend and the correlation of both［Hydrological station in the Nanwan station measured ice duration changes with time（a），by MOD09GQ data of glacial duration changes with time（b），Nanwan station with the duration of the glacial obtained in this paper（c）］

因此，本文認為通過MOD09GQ第二波段數(shù)據(jù)可以有效監(jiān)測額爾齊斯河布爾津至齋桑泊段的河冰物候信息。不僅如此，相比于水文站分布稀疏、只能監(jiān)測較短距離河段物候信息的缺點，遙感數(shù)據(jù)可以更加靈活地監(jiān)測不同河段河冰物候信息。然而遙感數(shù)據(jù)空間分辨率與時間分辨率不能兼得，因此想要監(jiān)測河冰物候信息，只能選取時間分辨率高而空間分辨率較低的數(shù)據(jù)，從而遙感數(shù)據(jù)更適宜于監(jiān)測大型河流的河冰物候情況。

2.4 開封河時間分析

表2 顯示，額爾齊斯河布爾津至齋桑泊段的開始封河時間集中于每年10月，冰/雪最大覆蓋時間主要集中于每年1月，結束開河時間主要集中于每年4月。最早開始封河年份為2016年，日期是9月29日，最晚開始封河年份為2013年，日期是11月13日最早封河比最晚封河時間提前45天；最早結束開河在2013年，其日期為3月28日，最晚結束開河為2002年，其日期為4月28日，最早開河比最晚開河時間提前31天，封河與開河時間均提前。分析發(fā)現(xiàn)，在2001—2018年這17年中，開始封河時間與結束開河時間總體均有提前趨勢［圖7（a）、7（b）］，但開始封河時間隨年份變化的斜率絕對值（0.623）大于結束開河時間隨年份變化斜率的絕對值（0.267）。因此可以認為，封河比開河提前的時間長，導致冰期持續(xù)時間將延長［圖6（b）］。

圖7 開始封河時間與結束開河時間變化趨勢（縱軸表示每一年的第幾天）Fig.7 The change trend of river freezing time and river breaking time（The vertical axis represents the day of each year），［The change of river break-up time（a），the change of river freeze-up time（b）］

2.5 封河與開河速率變化分析

開河速率與封河速率無明顯大小區(qū)別，但在速率變化上確有很大不同。通過線性分析，得到封河速率隨時間變化的斜率為-0.007，開河速率隨時間變化的斜率為0.021，表明2001—2018年，封河速率有變慢趨勢，開河速率有增快趨勢（如圖8）。若額爾齊斯河布爾津段至齋桑泊段通航，則每年秋冬季（即封河期），由于封河速度變慢，河冰形成緩慢，可通航時間將延后。而在每年春季（即開河期），由于開河速度加快，河冰將快速減少，河流達到通航條件時間將提前。因此，隨著秋冬季封凍期可通航時間延后，春季開河期可通航時間提前，每年的總通航時間將增加，可帶來更大經(jīng)濟效益。

圖8 開封河速率變化趨勢Fig.8 The change trend of river break-up and river freeze-up rates is that the river freeze-up rate decreases gradually（a）and the river break-up rate accelerates gradually（b）

2.6 開封河速率與溫度的關系

河冰開河機制復雜，受熱力、動力及河流形態(tài)（如多支流河流）等諸多因素影響，因此僅分析封河速率與溫度的相關性。本文使用布爾津和哈巴河兩個氣象觀測站的2001—2017年逐日平均氣溫，分析封河期平均溫度與封河速率的相關性，結果顯示，在封河期，布爾津站和哈巴河站觀測到的平均溫度隨時間變化的斜率分別為0.155、0.124，均為正值具有上升的趨勢［圖9（a）、9（b）］，這與高沈瞳等［4］得到的額爾齊斯河流域年均溫及夏、秋、冬溫將繼續(xù)保持升溫的趨勢的結果保持一致。圖9（c）、9（d）顯示，布爾津站封河期平均溫度與封河速率的相關系數(shù)R=-0.179；哈巴河站封河期平均溫度與封河速率的相關系數(shù)R=-0.299，封河期平均溫度均與封河速率呈負相關，因此，隨著封河期平均溫度升高，封河速率變慢，本文認為封河期平均溫度是封河速率變化的影響因子之一。

圖9 兩個氣象觀測站的封河期平均溫度變化情況及封河期平均溫度與封河速率的相關關系Fig.9 average temperature changes of two meteorological observation stations Burqin station（a）and Habahe station（b）during the river freeze-up period，and the correlation between the average temperature during the river freeze-up and the river freeze rate measured by Burqin station（c）and Habahe station（d）

3 結論與討論

本文基于MOD09GQ數(shù)據(jù)，使用HANTS算法及TIMESAT軟件，研究額爾齊斯河布爾津至齋桑泊河段的河冰物候，得到以下結論：

（1）通過HANTS算法和TIMESAT軟件的處理，有效去除了云噪聲對MODIS影像的影響，反演得到河冰物候信息。對比分析2015年2月2日（即最大雪/冰覆蓋日）、3月15日、4月8日及4月12的反射率可知，隨著河冰開河，河冰反射率明顯降低，并且與南灣水文站觀測數(shù)據(jù)的對比結果表明，MODIS數(shù)據(jù)所得與南灣水文站所得冰期持續(xù)時間的相關系數(shù)R=0.306呈正相關，本文研究結果與實測數(shù)據(jù)具有一致性，因此認為MOD09GQ的第二波段數(shù)據(jù)能有效監(jiān)測河冰冰情變化。

（2）在額爾齊斯河布爾津至齋桑泊河段，開始封河時間和結束開河時間均有提前的趨勢，最早開始封河時間為2016年9月29日，最晚開始封河時間為2013年11月13日，最早開始封河比最晚開始封河提前45天；最早結束開河時間為2013年3月28日，最晚結束開河時間為2002年4月28日，最早結束開河比最晚結束開河時間提前31天。封河與開河時間均提前，封河比開河提前的時間長，從而導致冰期持續(xù)時間有延長趨勢。

（3）該河段封河速率隨時間變化的斜率為-0.007，開河速率隨時間變化的斜率為0.021，呈現(xiàn)出封河速率逐漸下降、開河速率逐漸加快的趨勢。封河速率與布爾津站封河期平均溫度的相關系數(shù)R=-0.179，與哈巴河站封河期平均溫度的相關系數(shù)R=-0.299，封河速率與封河期平均溫度均呈負相關，因此，封河期平均溫度升高是引起封河期速率下降的一個因素。若未來該河段通航，則在封河期河流凍結速度變慢，可通航時間將延后，在開河期河冰加速消融，可通航時間將提前，每年總通航期延長，可帶來更大經(jīng)濟效益。

本研究所得結果雖與實測數(shù)據(jù)所得結果相符，但尚有諸多不足之處，MODIS影像時間分辨率高，可有效監(jiān)測河冰隨時間的變化情況，但其空間分辨率較低，使用MODIS數(shù)據(jù)研究存在河寬不足250 m的河流，不能監(jiān)測其具體結冰與破冰過程，并且可能存在混合像元的情況，為消除混合像元的影響，

使用整個河段的平均反射率進行研究，僅能描述所研究河段河冰整體隨時間變化的趨勢。此外，2015年4月16日高分1號影像缺失，并且直到4月24僅有一幅研究區(qū)上游影像，預覽該影像顯示研究區(qū)河段無河冰，未能觀測到4月16日是否研究區(qū)全河段開河。額爾齊斯河屬于跨國河流，可獲取水文資料有限，僅能通過水文站監(jiān)測所得冰期持續(xù)時間對本文結果進行驗證。下一步可使用微波遙感數(shù)據(jù)監(jiān)測河冰形成及破碎過程，結合使用無人機實際觀察河流封河及開河過程，了解河冰具體變化過程，并驗證MODIS數(shù)據(jù)所得結果。