基于溫度植被指數(shù)TVDI的拉薩地區(qū)土壤濕度特征分析

鐘 偉，盧宏瑋，管延龍，仇理化

（1.華北電力大學(xué)水利與水電工程學(xué)院，北京 102206；2.中國科學(xué)院地理科學(xué)與資源研究所陸地水循環(huán)及地表過程重點(diǎn)實(shí)驗(yàn)室，北京 100101；3.南方科技大學(xué)環(huán)境科學(xué)與工程學(xué)院，廣東深圳 518055；4.華北電力大學(xué)新能源學(xué)院，北京 102206）

0 引 言

快速發(fā)展的遙感技術(shù)為土壤濕度的獲取提供了重要的途徑，相較于傳統(tǒng)監(jiān)測方法，遙感技術(shù)具有數(shù)據(jù)覆蓋范圍廣、多時向、多光譜、高分辨率、易獲取等優(yōu)點(diǎn)，被廣泛應(yīng)用于區(qū)域旱情監(jiān)測研究中［6-10］。目前，基于遙感技術(shù)的土壤濕度估計方法主要有表觀熱慣量法、蒸散量計算和溫度植被干旱指數(shù)等［11-14］，其中溫度植被干旱指數(shù)TVDI有著明確具體的物理含義，受土壤背景和數(shù)據(jù)分辨率的影響較?。?5］，被廣泛應(yīng)用于土壤濕度估計。TVDI指數(shù)在國內(nèi)外已有大量相關(guān)研究，楊曦等［16］研究發(fā)現(xiàn)，TVDI指數(shù)能夠較好地反映華北平原表層土壤的干濕狀況，趙杰鵬等［17］改進(jìn)了Sandholt［14］提出的TVDI模型，并發(fā)現(xiàn)改進(jìn)后的模型較為明顯地提高了新疆地區(qū)的土壤濕度反演精度。Patel 和Sun［18，19］等基于MODIS 遙感數(shù)據(jù)分析發(fā)現(xiàn)TVDI指數(shù)與表層土壤濕度具有顯著的負(fù)相關(guān)關(guān)系。部分研究發(fā)現(xiàn)［20-24］，植被指數(shù)對于TVDI指數(shù)模擬精度影響較大，不同的植被指數(shù)各有優(yōu)劣。由于不同季節(jié)和不同下墊面條件，植被指數(shù)對于植被覆蓋的敏感性具有顯著差異，導(dǎo)致TVDI不同植被覆蓋下的監(jiān)測能力存在較大的差異性。

早期TVDI的研究主要集中在平原農(nóng)業(yè)區(qū)、旱區(qū)［16，17，23］等，而針對高寒、地形條件復(fù)雜區(qū)域的探索仍較為缺乏。楊玲和楊艷昭［25］發(fā)現(xiàn)西遼河流域干旱區(qū)域主要位于流域的上中游地區(qū)，濕潤區(qū)域主要沿河岸分布，且不同土地覆被類型的土壤濕度差異較大。Jian 等［26］研究指出，在不同土壤深度下，黃淮海平原TVDI與實(shí)測土壤濕度之間表現(xiàn)為顯著負(fù)相關(guān)，且在10～20 cm之間的關(guān)系最為密切。李躍鵬等［27］指出烏魯木齊地區(qū)不同下墊面條件下，地表溫度與植被指數(shù)均呈顯著負(fù)相關(guān)關(guān)系，TVDI與ERA-interim Level 2層體積含水量數(shù)據(jù)擬合最好。拉珍等［28］發(fā)現(xiàn)西藏一江兩河地區(qū)TVDI空間分布特征與降水空間分布一致性較好，實(shí)測地表溫度與TVDI兩者具有較好的線性正相關(guān)性。這些研究表明，不同的地理條件可能導(dǎo)致TVDI的區(qū)域適用性存在顯著差異。

青藏高原是地球上獨(dú)特的地質(zhì)、地理和生態(tài)單元，獨(dú)特的地形和高寒氣候條件使得青藏高原成為全球氣候變化的敏感區(qū)。土壤濕度作為陸-氣系統(tǒng)相互作用的重要物理參數(shù)，在青藏高原等生態(tài)脆弱區(qū)的水-碳循環(huán)中扮演著至關(guān)重要的作用。盡管已有針對高原土壤濕度相關(guān)研究［10，24，28］，但缺乏針對下墊面條件復(fù)雜以及海拔差異大等情況的深入探討。拉薩地區(qū)作為青藏高原高寒氣候的代表性區(qū)域，地表類型復(fù)雜，地面監(jiān)測站點(diǎn)稀少，土壤濕度特征及影響因素尚待揭示。因此本文基于MODIS 高分辨率地表溫度和植被指數(shù)產(chǎn)品等，結(jié)合第二次青藏科學(xué)考察采樣實(shí)測數(shù)據(jù)，以及區(qū)域DEM 數(shù)據(jù)、氣象數(shù)據(jù)和土地利用數(shù)據(jù)，引入溫度植被指數(shù)TVDI對拉薩地區(qū)土壤濕度狀況進(jìn)行分析，進(jìn)一步探究土地利用類型、海拔以及氣象要素等因素對該地區(qū)土壤濕度的影響。

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

1.1 研究區(qū)概況

拉薩地處青藏高原南部，西藏自治區(qū)東南部，雅魯藏布江支流拉薩河北岸，位于東經(jīng)91°06′，北緯29°36′。全年多晴朗天氣，年日照時間3 000 h 以上，年降雨量500 mm 左右，屬于典型高原溫帶半干旱季風(fēng)氣候。區(qū)域內(nèi)海拔跨度較大，地勢由北向東南傾斜，高山河谷相間分布［圖1（a）］，地表類型復(fù)雜［圖1（b）］。區(qū)域內(nèi)生態(tài)系統(tǒng)單一脆弱，農(nóng)業(yè)主要集中在地勢平坦，灌溉條件較好的河谷地區(qū)。土地利用類型包括了耕地、林地、草地、濕地、建設(shè)用地以及未利用土地［表1，圖1（b）］，草地是最主要的土地利用類型，呈較均勻地分散狀分布，面積占比為56.62%，其次為未開發(fā)土地，占21.61%，耕地、林地和濕地面積占比較小，濕地主要位于北部海拔較低區(qū)域，林地主要分布在東南部山谷盆地。面積占比最小的為建設(shè)用地，僅為0.39%。

表1 拉薩地區(qū)主要土地利用類型面積及占比Tab.1 Area and proportion of main land use types in Lhasa

1.2 數(shù)據(jù)來源

（1）遙感數(shù)據(jù)和氣象數(shù)據(jù)。本文使用遙感數(shù)據(jù)來源于NASA 數(shù)據(jù)中心（https：//earthdata.nasa.gov），包括植被指數(shù)產(chǎn)品（MOD13A2）、地表溫度產(chǎn)品（MOD11A2）和短波地表反照率產(chǎn)品（MCD43C3）。數(shù)據(jù)選取時段為2018年7-9月，其中增強(qiáng)型植被指數(shù)EVI共6 期，時間分辨率為16 d，空間分辨率為1 km。地表溫度LST 和地表反照率Albedo 共12 期，時間分辨率均為8 d，空間分辨率分別為1 km和0.05°。遙感影像采用MRT等軟件進(jìn)行預(yù)處理后用于本文研究。研究區(qū)2018年土地利用遙感監(jiān)測數(shù)據(jù)來源于中科院資源環(huán)境科學(xué)與數(shù)據(jù)中心（http：//www.resdc.cn），拉薩地區(qū)2018年夏季日均氣溫和降水?dāng)?shù)據(jù)來源于中國氣象數(shù)據(jù)網(wǎng)（http：//data.cma.cn/），經(jīng)處理后得到16 d 平均氣溫和累積降水量用于本文研究。

增強(qiáng)型植被指數(shù)EVI的紅光和近紅外探測波段范圍設(shè)置更窄，改進(jìn)了針對稀疏植被的探測效果，同時，引入了藍(lán)光波段對大氣氣溶膠的散射和土壤背景進(jìn)行了矯正，減少了大氣和土壤背景的影響［20，22，29，30］，因此本研究采用EVI代替歸一化植被指數(shù)NDVI?？紤]到植被指數(shù)為負(fù)表示湖泊積雪等，本研究在統(tǒng)計計算中暫不納入EVI值為負(fù)的像元。利用短波地表反照率產(chǎn)品（albedo）提供的黑空以及白空地表反照率數(shù)據(jù)，通過線性加權(quán)平均法得到短波波段的地表真實(shí)反照率?；贛RT（MODIS Reprojection Tool）將LST和EVI數(shù)據(jù)統(tǒng)一時間尺度為16 d，計算各時期TVDI并構(gòu)建LST-EVT二維特征空間。

（2）采樣數(shù)據(jù)。繼20 世紀(jì)70年代進(jìn)行首次青藏高原綜合科學(xué)考察研究后，第二次青藏高原綜合科考從2017年8月全面啟動，在收集大量數(shù)據(jù)資料的基礎(chǔ)上，對青藏高原的水、生態(tài)、人類活動等環(huán)境問題進(jìn)行考察研究，為青藏高原生態(tài)保護(hù)和社會經(jīng)濟(jì)發(fā)展提供科學(xué)依據(jù)。本文所用實(shí)際土壤濕度數(shù)據(jù)來源于2018年7月考察支隊在西藏拉薩河流域進(jìn)行農(nóng)牧業(yè)發(fā)展以及環(huán)境問題考察時的樣點(diǎn)實(shí)測數(shù)據(jù)，選取分布在研究區(qū)內(nèi)的17個樣點(diǎn)，具體位置以及分布情況見圖1（b）和表2，覆蓋了耕地，草地，林地以及濕地等主要土地利用類型。本文擬合分析方法采用最小二乘法擬合，顯著性檢驗(yàn)采用t檢驗(yàn)法。

表2 采樣點(diǎn)位置以及土地覆被類型Tab.2 Sampling point location and land cover type

1.3 TVDI計算

Sandhol等［14］基于地表溫度與植被指數(shù)之間的對應(yīng)關(guān)系，構(gòu)建了LST-NDVI特征空間模型，指出溫度植被干旱指數(shù)TVDI能夠有效地反演區(qū)域陸表干濕變化。其中，特征空間的上邊界被定義為干邊，表示研究區(qū)土壤干旱特征，將下邊界定義為濕邊，反映土壤的濕潤狀況。具體計算公式如下：

高校思想政治教育中融入中國傳統(tǒng)文化能夠有效向大學(xué)生傳遞正確價值觀念，從而促進(jìn)大學(xué)生在日常生活中加強(qiáng)精神文明建設(shè)，提升大學(xué)生愛國思想。通過中國傳統(tǒng)文化的融入，能夠使大學(xué)生在接受思想政治教育時更深刻理解國家對自身的重要性，從心底產(chǎn)生愛國情感，并積極主動投入社會傳統(tǒng)文明建設(shè)中，從而加強(qiáng)我國社會傳統(tǒng)文明建設(shè)，滿足我國現(xiàn)代社會發(fā)展需求。

式中：TVDI為任意像元的溫度植被干旱指數(shù)；TS和NDVI為分別為任意像元的地表溫度和植被指數(shù)，TSmax、TSmin分別為某一NDVI對應(yīng)的最大最小地表溫度，即分別代表了干邊和濕邊；a1、a2和b1、b2分別為干邊和濕邊擬合方程的系數(shù)。

其中，干邊對應(yīng)TVDI值為1，濕邊對應(yīng)TVDI值為0，計算得到任一像元TVDI值介于0 和1 之間，TVDI值越大，對應(yīng)的像元土壤濕度越低；TVDI值越小，對應(yīng)的像元土壤濕度越高［31］。

2 結(jié)果與分析

2.1 TVDI與實(shí)測土壤濕度關(guān)系

圖2所示為實(shí)測土壤體積含水量與TVDI之間的線性擬合。具體而言，研究區(qū)內(nèi)17 個樣點(diǎn)實(shí)測表層土壤體積含水率與TVDI呈顯著負(fù)相關(guān)關(guān)系（p＜0.05），隨著土壤濕度的增加，TVDI呈現(xiàn)逐漸下降的趨勢，擬合系數(shù)r2達(dá)到0.71，擬合效果良好。所選取的樣點(diǎn)基本覆蓋了研究區(qū)主要的土地利用類型，表明在不同的下墊面條件下，總體上TVDI值越高土壤濕度越低。擬合結(jié)果有效反映了TVDI在反演表層土壤干濕狀況方面具有良好的表現(xiàn)，可以作為研究區(qū)地表土壤濕度監(jiān)測的有效手段和補(bǔ)充。

2.2 地表溫度與植被指數(shù)關(guān)系

圖3所示為不同土地利用類型的地表溫度LST與植被指數(shù)NDVI的擬合結(jié)果。對于耕地，隨著EVI的上升，地表溫度逐漸下降，整體上呈顯著的負(fù)相關(guān)（p＜0.05）。對于林地和草地，隨著EVI的增加，LST呈逐漸上升趨勢，都通過了0.05 的顯著檢驗(yàn)。對于濕地，地表溫度LST和植被指數(shù)EVI同樣呈顯著正相關(guān)關(guān)系（p＜0.05），上升趨勢高于林地和草地。不同土地利用類型的擬合結(jié)果表明，拉薩地區(qū)植被變化與地表溫度變化之間存在明顯的響應(yīng)過程，空間異質(zhì)性強(qiáng)烈。

圖4所示為不同土地利用類型之間地表溫度與海拔變化的擬合關(guān)系?？傮w上，不同土地利用類型下LST均與海拔均存在顯著的負(fù)相關(guān)關(guān)系（p＜0.05），隨著海拔高度的升高，地表溫度逐漸下降。具體而言，耕地主要集中分布于海拔3 500～4 500 m范圍內(nèi)，少部分位于4 500 m 以上，跨度相對較小，地表溫度表現(xiàn)為隨植被指數(shù)增加而顯著降低。林地和草地基本上均勻分布于海拔3 500～5 500 m 之間以及3 500～6 000 m 之間，垂直分異明顯，表明林地和草地多分布于地形坡度較大的山區(qū)，地表溫度受海拔梯度變化影響顯著，擬合系數(shù)明顯高于耕地。對于濕地而言，整體面積所占比例較小，多分布于海拔3 500～3 800 m和海拔5 000～6 000 m處之間，空間分布異質(zhì)性較為明顯。結(jié)合研究區(qū)地形條件和土地利用類型分布（圖1），研究區(qū)地表溫度LST與植被指數(shù)EVI呈顯著的負(fù)相關(guān)，空間異質(zhì)性較為明顯。

表3 不同土地利用LST與海拔擬合方程Tab.3 Fitting equation of LST and altitude over different land use types

圖5和表4所示為研究區(qū)2018年夏季7-9月每16 d的LSTEVI特征空間及其各時段干濕邊擬合方程及擬合系數(shù)。除07-28-08-12 時段外，各時段干邊斜率均顯著小于0，總體上LSTmax與EVI呈顯著負(fù)相關(guān)關(guān)系（p＜0.05）。各時段濕邊斜率均顯著大于0，總體上LSTmin與EVI呈顯著正相關(guān)關(guān)系（p＜0.05）。大部分時段的最高和最低地表溫度與植被指數(shù)擬合效果比較好（表4），r2基本大于0.4，最大達(dá)到0.83，能較好地擬合出干邊和濕邊方程。在EVI值較低（0＜EVI＜0.1）時，最高地溫表現(xiàn)為隨著植被指數(shù)增加而逐漸升高，這可能是在植被覆蓋處于較低水平時，植被指數(shù)描述植被覆蓋的動態(tài)變化不夠準(zhǔn)確，從而引起該區(qū)間地表溫度對植被指數(shù)的變化不夠敏感［18］。

表4 7-9月各時期LST-EVI特征空間干濕邊擬合方程及擬合系數(shù)Tab.4 Fitting equations and fitting coefficients for dry and wet edges of LST-EVI feature space from July to September

2.3 土壤濕度時空格局及影響因素

圖6所示為研究區(qū)2018年7-9月每16 d的TVDI空間分布。在7月和8月，大部分地區(qū)TVDI處于0～0.2 之間［圖6（a）～（d）］，地表整體較為濕潤。進(jìn)入9月，TVDI開始升高至0.2～0.6 區(qū)間，土壤濕度降低，少部分地區(qū)TVDI大于0.6，處于較為干燥的狀態(tài)。整體上看，研究區(qū)夏季土壤濕度處于較為良好的水平，并隨時間推移呈現(xiàn)逐漸變干的趨勢。

參考已有研究方法［27］的相關(guān)閾值將研究區(qū)土壤干濕狀況劃分為極濕潤至極干旱五個等級（表5），基于7-9月6期的平均TVDI，得到夏季地表干濕情況空間分布如圖7所示。從空間上看，大部分地區(qū)處于極濕潤和濕潤的狀態(tài)，分別占研究區(qū)面積的57.01%和22.27%。土壤濕度正常的區(qū)域占比13.14%，而干旱和極干旱區(qū)域主要分布在北部，西南以及中南部分地區(qū)，面積占比分別為7.06%和0.51%。

表5 基于TVDI閾值干旱等級劃分標(biāo)準(zhǔn)Tab.5 The classification of TVDI and drought

圖8和表6所示為研究區(qū)07-28-08-12 時段所計算的TVDI，考慮不同下墊面的TVDI與海拔擬合結(jié)果。隨著海拔的不斷升高，幾種主要土地利用類型TVDI均與海拔呈顯著負(fù)相關(guān)關(guān)系（p＜0.05），耕地、林地、草地和濕地的TVDI下降率分別達(dá)到0.151、0.235、0.281 和0.211/km，表明研究區(qū)土壤濕度隨海拔的升高而顯著增加。

表6 07-28-08-12時段不同土地利用類型TVDI與海拔擬合方程Tab.6 Fitting equation of TVDI and altitude over different land use types during 07-28-08-12

研究區(qū)7-9月平均植被指數(shù)和地表反照率空間分布如圖9所示，結(jié)合地形和下墊面情況（圖1），西部和北部高山地區(qū)以裸土為主，植被覆蓋較低，對應(yīng)地表反照率最高達(dá)0.365，造成土壤涵水能力不足，即使有降水補(bǔ)充但由于水分散失較快，導(dǎo)致土壤濕度處于較低水平。而對于面積占比較大的林地和草地，EVI指數(shù)處于0.15～0.5 之間，植被覆蓋相對完整，對應(yīng)地表反照率較低，雖然植被蒸散消耗水分，但相比裸土直接蒸發(fā)更為緩和，土壤濕度能保持較好的水平。

圖10 所示為夏季7-9月每16 d 的平均氣溫和累積降水量以及平均TVDI。整體來看，TVDI與氣溫正相關(guān)變化，與降水呈負(fù)相關(guān)變化。當(dāng)氣溫較低，降水充足的時候［圖10（08-13-08-28）］，TVDI下降至最低值，表明土壤含水率較高，較為濕潤。而當(dāng)氣溫升高，降水減少時［圖10（08-29-09-13）］，TVDI上升至最高點(diǎn)，土壤濕度迅速降低。TVDI與氣溫和降水的相關(guān)系數(shù)r（皮爾遜線性相關(guān)）分別為0.41和-0.67，表明降水對于研究區(qū)土壤濕度的影響較氣溫更為明顯。

3 結(jié) 論

基于MODIS 地表溫度（LST）和增強(qiáng)型植被指數(shù)產(chǎn)品（EVI）以及土地利用類型數(shù)據(jù)等，構(gòu)建了拉薩地區(qū)LST-EVI特征空間和溫度植被干旱指數(shù)TVDI，驗(yàn)證了TVDI在拉薩地區(qū)的適用性，分析了土壤濕度時空格局，進(jìn)一步探究了土地利用類型、海拔以及氣溫降水等因素對TVDI的影響，主要結(jié)論如下。

（1）TVDI與實(shí)測表層土壤體積含水率擬合效果良好，呈顯著負(fù)相關(guān)（r2=0.71，p＜0.05），表明TVDI能較好地反映研究區(qū)土壤濕度狀況。

（2）不同土地利用類型地表溫度與植被指數(shù)的關(guān)系具有明顯的差異性。耕地地表溫度隨著植被指數(shù)增加而顯著減小，林地、草地和濕地受海拔影響，地表溫度與植被指數(shù)呈顯著正相關(guān)關(guān)系。研究區(qū)夏季6 期LST-EVI特征空間，其散點(diǎn)圖呈三角形狀，均能夠較好地擬合出干邊和濕邊（r2基本大于0.4），均通過了顯著性水平為0.05的顯著性檢驗(yàn)。

（3）拉薩地區(qū)夏季整體上土壤濕度處于較為良好的水平，并隨時間推移呈現(xiàn)逐漸變干的趨勢?？臻g上，大部分地區(qū)處于極濕潤（0～0.2）和濕潤（0.2～0.4）的狀態(tài)，干旱（0.6～0.8）和極干旱（0.8～1.0）區(qū)域主要分布在北部，西南以及中南部分地區(qū)，面積占比僅分別為7.06%和0.51%。不同土地利用下土壤濕度隨海拔的升高均顯著增加（p＜0.05），TVDI基本與氣溫呈同向變化，與降水呈反向變化，降水對于土壤濕度影響更為明顯。

青藏高原地表類型復(fù)雜，受地形大氣環(huán)流等因素影響，植被覆蓋水平和垂直分異明顯，屬于典型高原氣候區(qū)。區(qū)內(nèi)監(jiān)測站點(diǎn)稀少，依靠傳統(tǒng)方式難以全面動態(tài)地掌握土壤濕度變化情況。本文以拉薩作為代表性區(qū)域，借助遙感技術(shù)反演了研究區(qū)土壤濕度，并進(jìn)一步考慮了其時空格局以及影響因素，以期為地表類型復(fù)雜的高寒地區(qū)土壤濕度研究提供參考。由于實(shí)測數(shù)據(jù)缺乏等客觀原因限制，本研究對研究區(qū)表層土壤濕度進(jìn)行了初步分析，未深入考慮不同深度層土壤濕度變異情況，且對于氣候要素的影響分析探討略淺，植被的動態(tài)響應(yīng)也尚未充分考慮，這是本研究的不足之處，也將是下一階段的工作重點(diǎn)。□