劉 馨，宋小寧?，冷 佩，夏 龍

(1 中國科學院大學資源與環(huán)境學院， 北京 100049； 2 中國農業(yè)科學院農業(yè)資源與農業(yè)區(qū)劃研究所， 北京 100081)

干旱是嚴重的自然災害之一，持續(xù)的干旱災害會直接影響農業(yè)生產和經濟發(fā)展，長期干旱甚至會造成植被減少、河水斷流、沙漠化現(xiàn)象[1]。有證據(jù)表明，自20世紀70年代以來，干旱等極端氣候現(xiàn)象頻次不斷增加[2]。因此，加強干旱監(jiān)測，對于農業(yè)管理、旱災防治、水資源管理及保護都有著極強的現(xiàn)實意義[3]。

黃河源是黃河重要的水源補給區(qū)，近30年來由于人類活動及氣候變化的影響，黃河源區(qū)氣溫上升，下游頻繁斷流，極大地影響著中下游地區(qū)的生態(tài)環(huán)境和生產生活。因此對黃河源區(qū)進行區(qū)域尺度上的干濕狀況時空分析對于水資源管理和生態(tài)環(huán)境保護具有重要意義。

干旱作為一種復雜的現(xiàn)象難以直接觀測，因此通常采用干旱指標進行評估?；谟^測手段的不同，干旱監(jiān)測評估主要可分為兩類：基于站點的觀測指標和基于遙感的觀測指標[2]。相比傳統(tǒng)基于站點的觀測，遙感技術能夠實時動態(tài)地獲取地表參數(shù)的時空分布信息，從而表征地表干濕狀況，具有覆蓋范圍廣，多時相等優(yōu)點，是目前干旱監(jiān)測的重要手段之一[4]。

干旱過程與土壤干濕狀況以及作物水分虧缺有著緊密聯(lián)系。20世紀90年代開始，基于可見光-近紅外波段反演的植被指數(shù)和基于熱紅外反演的地表溫度形成的三角形或梯形的特征空間逐漸發(fā)展用于旱情遙感監(jiān)測。由于相同光照條件下，植被溫度在空間上沒有變化，三角形內的溫度變化只反映土壤表面的干燥程度。隨著植被覆蓋增加，溫度變化范圍逐漸縮小，因此溫度/植被特征空間內的散點隨著植被增加而傾斜形成三角形或梯形。三角形邊界反映物理邊界限制，即裸露土壤或完全的植被覆蓋，以及完全干燥或完全濕潤的土壤狀況[5]。三角形概念最初由Price[6]提出，后由Carlson等[7-8]詳細闡述。Moran等[9]提出用三角形空間內作物水分脅迫指數(shù)等值線評估水分虧缺。隨后，學者們陸續(xù)利用溫度/植被特征空間研究土壤水分反演[10-12]。Sandholt等[12]在水分虧缺指數(shù)WDI (water deficit index) 基礎上發(fā)展了溫度植被干旱指數(shù)TVDI (temperature vegetation dryness index)。該指數(shù)假設土壤水分是溫度變化的主要來源，指數(shù)的計算無需地面的輔助數(shù)據(jù)，通過光學和熱紅外遙感數(shù)據(jù)的反演就可直接計算。該指數(shù)已被廣泛用于土壤水分反演和干旱狀況的監(jiān)測[12-16]，其結果表明，TVDI表征的干旱程度及范圍與實地監(jiān)測氣象數(shù)據(jù)和土壤水分數(shù)據(jù)吻合[17-18]。因此，本研究基于LST-NDVI特征空間通過計算TVDI來研究黃河源區(qū)土壤相對干濕狀況，進而分析干旱時空變化，為有效對黃河源地區(qū)進行水資源管理和生態(tài)環(huán)境保護提供科學依據(jù)。

1 研究區(qū)與數(shù)據(jù)源

1.1 研究區(qū)

黃河源區(qū)(95°50′E～103°30′E，32°30′N～36°10′N)為河源至唐乃亥之間的匯水區(qū)域(圖1)，流經青海、甘肅、和四川3省[19]。該區(qū)域集水面積為1.22×105km2，占黃河流域面積的15.3%，多年平均徑流量為2.051×1011m3，占全流域的34.5%，是黃河流域重要的水源地和產流區(qū)，在黃河流域有著不可或缺的重要作用[18]。源區(qū)平均海拔4 000 m以上，為低山谷和湖盆組合地貌。植被類型以高寒草甸、高寒沼澤草甸和高寒草原為主，土壤類型包括高山草甸土、沼澤化草甸土、高山荒漠土、高山草原土等[20]。

圖1 黃河源地區(qū)示意圖Fig.1 Location of the source area of Yellow River(SAYR)

黃河源區(qū)氣候屬典型高原大陸性高寒氣候，空氣稀薄、輻射強烈且日照時間長。受西南季風影響比較明顯，氣溫和降水從東南向西北呈遞減趨勢。源區(qū)年均氣溫在5 ℃左右，源區(qū)多年平均降水量為426 mm，年內變化較大且分布極度不均勻，主要集中在6—9月，約占全年降水量的75%～90%甚至以上[20]。

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

本研究使用遙感數(shù)據(jù)為2007—2016年7—9月(時相第193～257天)黃河源區(qū)的地表溫度和NDVI的MODIS影像數(shù)據(jù)，分別為16天合成的1 km分辨率的MOD13A2的NDVI產品數(shù)據(jù)和8天合成的1 km分辨率的MOD11A2的LST產品數(shù)據(jù)。數(shù)據(jù)處理過程中使用MRT對上述數(shù)據(jù)進行幾何校正、投影變換的批量處理，使用IDL將數(shù)據(jù)批量轉換成實際地表溫度和NDVI的數(shù)值，同時對NDVI數(shù)據(jù)水體進行掩膜處理，以消除水體對特征空間可能造成的影響?？紤]到植被指數(shù)的值對植被生長狀況的敏感性問題，只選取NDVI處于0.1～0.9之間的數(shù)據(jù)進行干濕邊反演計算。此外，利用Arcgis裁剪出產品數(shù)據(jù)中所需的研究區(qū)域部分，剔除異常值。同時，為保持時相一致性，將8天合成的MOD11A2的LST產品數(shù)據(jù)合成為每16天的LST數(shù)據(jù)。

為了對構建的TVDI作為旱情指標進行可信度檢驗，本研究利用黃河源瑪曲地區(qū)(瑪曲站點連續(xù)原位土壤水分與土壤溫度的測量——CEOP-AEGIS項目)站點提供的表層土壤水分觀測數(shù)據(jù)與TVDI進行相關性分析。圖1(b)為這20個站點的分布位置(具體信息見表1)，站點監(jiān)測頻率為每15 min一次。光學-熱紅外遙感所在波段穿透性弱，因此研究選擇表層土壤水分進行TVDI相關性分析。眾多研究選取的表層土壤水分深度不盡相同，大多數(shù)研究選擇0～10 cm深度范圍，而有的研究選取0～20 cm深度范圍內的土壤水分數(shù)據(jù)，均與TVDI有較好相關性[18，21-22]。本研究參照黃河源區(qū)的土壤水分相關研究，同時考慮到實測數(shù)據(jù)的可獲取性，選擇5 cm深度的土壤水分數(shù)據(jù)。由于衛(wèi)星過境的時間在11時左右，為保持實測土壤水分與MODIS數(shù)據(jù)時相一致，選取對應的每16天11時的土壤水分測量值計算平均土壤水分。

2 研究方法

2.1 溫度植被干旱指數(shù)

LST和NDVI 特征空間可以表征地表植被和干濕狀況，從而實現(xiàn)旱情監(jiān)測，在國內外有了廣泛的研究[8]。基于熱紅外遙感的地表溫度和基于可見光-近紅外波段的植被指數(shù)改善了單一使用表征地表干濕狀況的不足，可以有效地判斷作物長勢和干旱程度，同時解決了植物在水分虧缺時表征現(xiàn)象的滯后性。2002年Sandholt等[12]提出TVDI(溫度-植被干旱指數(shù))估測地表土壤水分狀況，由地表溫度和NDVI計算而得(計算示意圖如圖2)，其形成的三角形或梯形的空間特征表征地表的干濕狀況。計算公式如下

(1)

表1 實測站點信息Table 1 In-situ site information

式中：Ts為像元溫度；Tsmin為相同NDVI值的最小地表溫度，對應特征空間的濕邊；a+bNDVI對應特征空間的干邊，為相同NDVI值的最大地表溫度，a和b為干邊擬合線性關系的系數(shù)。干邊對應的TVDI為1，表示極度干旱；濕邊對應的TVDI為0，表示極度濕潤。TVDI范圍在0～1.0之間，地表干旱程度與TVDI數(shù)值呈正比[12,23]。

圖2 TVDI計算示意圖[13]Fig.2 TVDI triangle algorithm[13]

2.2 LST/NDVI特征空間構建

本研究使用批量重投影、拼接、裁剪等方法處理LST和NDVI產品數(shù)據(jù)，并利用處理好的產品數(shù)據(jù)構建出LST/NDVI特征空間。計算TVDI時，根據(jù)像元的NDVI值確定對應的溫度變化，從而確定該像元對應的干邊和濕邊，再根據(jù)像元的地表溫度值以及確定的干濕邊計算出TVDI。干濕邊的確定參照Garcia等[23]的方法，干邊對應公式(1)中(a+bNDVI)，濕邊與X軸平行，對應NDVI所在像元最低溫度的均值。如前所述，研究區(qū)植被覆蓋和土壤水分范圍涵括所有極端情況且像元足夠多時，LST/NDVI特征空間呈三角形，但實際情況中到達模擬干邊的溫度時，植被為減少水分流失會主動關閉氣孔，因此土壤水分濕度不會因為高溫而迅速減少。在NDVI高植區(qū)域的干邊像元數(shù)值也不可能為0，特征空間內的散點多呈梯形分布特征。

圖3為基于MODIS的地表溫度和NDVI產品數(shù)據(jù)擬合的TVDI梯形特征空間，選取2007—2016年間的部分擬合結果為例。由圖3可知，干濕邊的擬合效果較好，R2均為0.9以上，最高的甚至達到0.99。算法假定當特征空間內植被覆蓋度及土壤水分含量變化范圍較大時，特征空間呈梯形?？臻g內干濕邊為直線。LST/NDVI梯形特征空間對應的4個頂點在圖3中分別為：1)水分飽和的裸土；2)干燥裸土；3)水分充足的全覆蓋植被；4)水分虧缺的全覆蓋植被。

2.3 干旱等級劃分標準

本研究獲取TVDI以后，根據(jù)相關文獻提出的基于TVDI的干旱等級劃分標準[16, 19, 24-26]，將干旱程度做了等級劃分，便于后續(xù)對黃河源區(qū)域的干濕狀況進行時空分析。具體劃分指標如表2。

圖3 LST/NDVI特征空間Fig.3 LST/NDVI space

等級類型TVDI1極濕潤0～0.22濕潤0.2～0.43正常0.4～0.64干旱(輕度)0.6～0.85重旱0.8～1.0

3 結果與分析

3.1 TVDI作為旱情指標的可信度評價

由TVDI的物理意義可知，TVDI與土壤水分具有負相關性；眾多研究也表明TVDI與土壤水分存在負相關的線性關系[12-15, 27-28]。本研究利用TVDI作為干旱指標對黃河源區(qū)進行10年的干濕狀況時空分析，因此在分析前利用站點所測的表層土壤水分(5 cm)與TVDI進行相關性分析,有利于檢驗TVDI作為干旱指標表征黃河源區(qū)域尺度土壤干濕狀況變化的有效性。

由于云、水等因素影響導致部分LST和NDVI產品數(shù)據(jù)缺失，加上極端天氣、設備條件等因素影響站點的實測數(shù)據(jù)，研究選取時間尺度上橫向和縱向兩組數(shù)據(jù)：2010年D193～D257的5個時相的土壤水分數(shù)據(jù)，以及2008—2010年D209的3年土壤水分數(shù)據(jù)與對應的TVDI進行相關性分析，增強了驗證的可信度。由圖4可知，TVDI與土壤水分呈顯著的負相關線性關系，其中圖4 (a)為TVDI與2008—2010年期間D209時相的土壤水分的線性關系，圖4 (b)為TVDI與2010年D193～D257的5個時相的土壤水分的線性關系。圖4 (a)中線性擬合的R為0.71，圖4 (b)線性擬合的R為0.70，說明TVDI與土壤水分具有較好的相關關系。圖4 (a)、4(b)兩組線性關系略有差異，可能是由于時相選取不同，不同時相具有的不同天氣狀況(如云層、降雨等)會對遙感反演結果造成一定的影響；且氣候條件和設備條件也會影響實地站點的數(shù)據(jù)精度，這些因素都使得TVDI作為干旱指標可能會與實際地表干濕狀況有一定的偏差。但是總體來看，圖4兩組線性關系良好，說明TVDI很夠有效表征土壤表層干濕狀況。

3.2 黃河源區(qū)土壤干濕狀況時空格局分析

3.2.1 土壤濕度時間格局分析

圖5(e)為10年來TVDI自D193至D257五個時相的變化規(guī)律圖。逐時相來看，除2008年D257外，其余時相的趨勢大致一致。10年間的TVDI均在0.4～0.6，屬干旱分級中的正常狀況。各時相旱情大體上波動趨勢一致，說明在旱情嚴重的年份，不同時相旱情都有一定程度的加重，旱情緩解的年份，不同時相的旱情也有一定程度的緩解。各月旱情在10年期間的波動變化各有區(qū)別，也是由于不同時相當?shù)氐奶鞖庾兓?，如降雨、云層的影響，造成對旱情變化的短期影響?/p>

圖4 TVDI與土壤水分線性關系Fig.4 Linear relationship between TVDI and soil moisture

圖5 2007—2016年TVDI時間格局分析Fig.5 TVDI analysis of temporal patterns from 2007 to 2016

此外，根據(jù)中國科學院資源環(huán)境科學數(shù)據(jù)中心發(fā)布的2010年土地利用分類圖，本研究選取高、中、低覆蓋度草地區(qū)域，逐時相分析TVDI 在不同植被狀況下對土壤干濕狀況的指示。為控制降雨量、地形、地表覆蓋類型等變量因素，選取東部地區(qū)相離較近的實測站點所在區(qū)域。東部地區(qū)地勢平坦且植被覆蓋度較廣，選取的3個實測站點(CST_01, NST_01, NST_07)地形均為河谷區(qū)域，地表覆蓋類型為草地。TVDI變化規(guī)律如圖5(a)～5(d)所示。圖5(a)、(b)、(c)分別為中、低、高覆蓋度草地的TVDI趨勢圖，圖5 (d)為10年內D241時相上TVDI在3種不同覆蓋度草地的變化趨勢圖。東北部地區(qū)整體來看處于干旱-重旱。2007年—2016年間，在中低草地覆蓋區(qū)，7月上旬旱情最為嚴重，中覆蓋區(qū)有5年時間TVDI數(shù)值都在0.8及以上，為重旱程度，10年期間TVDI均在0.6以上，為干旱程度；低覆蓋區(qū)有7年時間TVDI達0.8以上，為重旱程度，10年期間TVDI均在0.6以上，為干旱程度。其次是9月上旬，中覆蓋區(qū)有4年時間TVDI數(shù)值在0.8以上，為重旱程度，10年期間TVDI均在0.6以上，為干旱程度；低覆蓋區(qū)有6年時間達0.8及以上，為重旱程度；10年期間TVDI均在0.6以上，為干旱程度。中低植被覆蓋區(qū)10年來旱情均較為嚴重，TVDI數(shù)值大致均在0.6～1間浮動，其中低植被覆蓋區(qū)所有時相的TVDI在10年內數(shù)值均在0.6及以上，旱情持續(xù)。各時相旱情在中、低覆蓋度區(qū)域TVDI年際趨勢相似，各月旱情變化基本具有一致性，而浮動差異可能是由于短期天氣變化造成的一定影響。圖5 (d)可以看出，NST_01與CST_01站點所在中、低植被覆蓋區(qū)域與NST_07所在高植被覆蓋區(qū)域10年間TVDI趨勢有所差異。中、低植被覆蓋區(qū)域在10年間TVDI變化大體相似，而高植被覆蓋區(qū)域TVDI波動較大，在2008、2012、2016年與中、低植被覆蓋區(qū)域的TVDI變化呈現(xiàn)出相反的趨勢。高植被覆蓋區(qū)域7月上旬有5年旱情嚴重，TVDI達0.8以上，其余時相旱情較輕，均僅有2～3年TVDI達0.8以上。盡管嚴重干旱的情況較中低覆蓋區(qū)域緩解，但是干旱依然在10年間持續(xù)發(fā)生。D193和D257在10年間TVDI均在0.6～1之間小幅波動，D209僅有1年TVDI略低于0.6；D225和D241分別有3年和2年的TVDI數(shù)值低于0.6，但也均大于0.5。干旱與植被覆蓋度互為因果，植被具有水分涵養(yǎng)的作用，植被覆蓋度低的地方易水土流失，引發(fā)旱情，而干旱持續(xù)發(fā)生時亦會造成植被覆蓋減少的情況[1]。

3.2.2 土壤濕度空間格局分析

本研究采用空間統(tǒng)計的分析方法，計算各等級旱情(參照表2)區(qū)域占黃河源區(qū)總面積的比例，從而分析旱情的時空分布特征。圖6(a)為2007—2016年間D193～D257時相的旱情占比頻譜圖，橫坐標為2007—2016年的D193～D257順序各時相，縱坐標為各旱情面積占比；圖6(b) 為2007—2016年間D193～D257時相的旱情占比直方圖，橫坐標為2007—2016年的D193～D257順序各時相，縱坐標表示各旱情面積占比及其總和。圖6(a)、6(b)結合分析可得，10年來7—9月正常干濕狀況的區(qū)域所占面積最多，最高占比達48%；干旱區(qū)域和濕潤區(qū)域面積交替占比第二，在10%～40%之間浮動；而極端情況即極濕潤和重旱區(qū)域占比均較少，在0～17%之間浮動。各旱情面積占比波動規(guī)律較為明顯，而干旱與濕潤等級的區(qū)域面積波動有著明顯的相反趨勢。

本文中以2016年7月上旬—9月上旬(D193～D257時相)的TVDI圖像為例進行說明。旱情的空間變化特征顯著：每年7月上旬西部地區(qū)旱情較輕，東北部地區(qū)和東南部受旱區(qū)域多為嚴重干旱。7月下旬，干旱區(qū)域逐漸向西部遷移，東南部旱情緩解，東北部仍然遭受嚴重干旱。8月上旬東北部旱情緩解，東南、西北兩部旱情有加重趨勢。8月下旬開始，西北部旱情緩解，且逐漸濕潤，旱情逐漸向東南部轉移。9月上旬西北部旱情繼續(xù)緩解，基本處于濕潤-正常的地表狀態(tài)，而東北、東南區(qū)域保持較嚴重的干旱狀態(tài)。中部區(qū)域在7—9月一直處于濕潤-正常的穩(wěn)定狀態(tài)。分布趨勢基本吻合中西部水分整體水分較充足，干旱比較嚴重的區(qū)域主要在東北部和東南部的分布規(guī)律。

4 結論與討論

本研究以黃河源為研究區(qū)，利用2007—2016年D193～D257的MODIS的1 km分辨率的地表溫度(LST)和植被指數(shù)(NDVI)的產品數(shù)據(jù)反演10年溫度-植被干旱指數(shù)(TVDI)。分析結果表明：

1)本研究提取NDVI對應的最大、最小地表溫度構建LST/NDVI特征空間，10年間D193～D257時相的LST/NDVI散點圖均構成梯形特征空間，符合前人的特征空間研究成果。

圖6 2016年7—9月(D193～D257)黃河源區(qū)旱情變化圖Fig.6 Variation of drought situation from July to September (D193-D257) in 2016

2)通過研究區(qū)站點監(jiān)測的實地土壤水分與反演的TVDI建立線性關系進行驗證，結果表明實地土壤水分與TVDI具有良好的線性負相關關系。2008—2010年實測土壤水分與TVDI線性關系的相關系數(shù)達到0.7，說明TVDI能作為干旱指標有效地指示黃河源區(qū)的干濕狀況。

3)10年間D193～D257的TVDI均值屬干旱分級中的正常狀況；各時相TVDI大體上波動趨勢一致，說明在旱情嚴重的年份，不同時相上旱情都有一定程度的加重，旱情緩解的年份，不同時相的旱情也有一定程度的緩解。東北部地區(qū)整體來看處于干旱-重旱，年年有旱情發(fā)生，7月上旬(D193)的干旱情況均最為嚴重。中、低植被覆蓋區(qū)域的干旱趨勢較為相似，10年來旱情均較為嚴重，其中低植被覆蓋區(qū)所有時相在十年間旱情持續(xù)。高植被覆蓋區(qū)域嚴重干旱的情況較中低覆蓋區(qū)域緩解，但是干旱依然在10年間持續(xù)發(fā)生。

4)旱情的空間變化特征顯著，7月上旬至9月上旬期間(D193～D257)，西部尤其是中部區(qū)域旱情較輕，嚴重旱情多集中于東北部和東南部區(qū)域，分布趨勢基本吻合中西部整體上土壤水分較充足，干旱比較嚴重的區(qū)域主要在東北部和東南部的分布規(guī)律。

根據(jù)中國科學院資源環(huán)境科學數(shù)據(jù)中心(http:∥www.resdc.cn/data.aspx?DATAID=125)提供的中國生態(tài)地理分區(qū)圖可知，黃河源區(qū)并非處于干旱區(qū)，而TVDI的旱情分級也只是提供黃河源區(qū)的相對干濕狀況信息，并非絕對的旱情。本研究以TVDI為指標，研究黃河源區(qū)土壤相對干濕狀況，為有效對黃河源地區(qū)進行水資源管理和生態(tài)環(huán)境保護提供科學依據(jù)。TVDI結合基于熱紅外遙感的地表溫度和基于可見光-近紅外波段的植被指數(shù)兩種指示因子，改善了單一使用表征地表干濕狀況的不足，同時解決了植物在水分虧缺時表征現(xiàn)象的滯后性。由于本研究使用NDVI、LST以及監(jiān)測站點土壤水分觀測值的時相不同，最終將所有數(shù)據(jù)統(tǒng)一對應NDVI時相，3種數(shù)據(jù)在時間對應上存在偏差，從而讓反演結果產生誤差。實地站點僅有2008—2010年3年的實地土壤水分數(shù)據(jù)，且由于極端天氣，站點設備等原因，許多實測值不具有對比驗證意義，可能也對遙感反演的數(shù)據(jù)驗證存在一定的影響。后續(xù)將對TVDI反演結果進行進一步的研究，如參照Hoffmann等[29]提出的方法，考慮改進干濕邊的擬合，為提高遙感影像反映地區(qū)干濕狀況進行更深一步的工作。