基于Landsat8 OLI遙感影像的延河流域土壤水分反演研究

姚 靜， 薛超玉， 焦 峰,3*

(1. 西北農(nóng)林科技大學(xué)水土保持研究所， 陜西 楊凌712100； 2. 中科宇圖科技股份有限公司， 北京100000；3. 中國科學(xué)院水利部水土保持研究所， 陜西 楊凌712100)

土壤水分是大氣和陸地之間物質(zhì)能量交換中的關(guān)鍵因子，是氣候、植被、地形及土壤因素等自然條件的綜合反映，同時也是水文學(xué)、氣象學(xué)、土壤學(xué)、生態(tài)學(xué)以及農(nóng)業(yè)科學(xué)等領(lǐng)域衡量土地退化和土壤干旱程度的重要指標(biāo)[1-5]。因此對于區(qū)域或全球范圍的土壤水分進行實時、有效地監(jiān)測和反演研究具有很重要意義，這也一直是國內(nèi)外研究的熱點問題[1,4,6]。

傳統(tǒng)的土壤水分監(jiān)測方法存在采樣速度慢、樣點較稀疏，而且需花費大量人力物力、代表范圍有限、數(shù)據(jù)收集時效性差等不足，使其難以滿足大范圍土壤水分分析研究的需要[7]。隨著遙感技術(shù)的不斷發(fā)展，特別是“3S”集成技術(shù)的日益成熟，不僅使大面積土壤水分實時動態(tài)監(jiān)測成為可能，而且提高了土壤水分與旱情遙感監(jiān)測的可行性和實用精度[8-9]。國內(nèi)外用遙感技術(shù)監(jiān)測土壤水分反演方法主要有熱慣量法[10-11]、植被指數(shù)法(Ts-NDVI)[12]和微波遙感法[13-14]。其中，溫度植被干旱指數(shù)(TVDI)是目前可見光和熱紅外波段遙感監(jiān)測土壤水分理論相對成熟，應(yīng)用比較廣泛的一種方法[4]。TVDI耦合了地表溫度(Ts)和植被指數(shù)(NDVI)信息，通過Ts-NDVI特征空間的變化特征分析土壤水分狀況以及時空動態(tài)演變規(guī)律。一方面，國內(nèi)相關(guān)領(lǐng)域的學(xué)者利用TVDI指數(shù)進行了全國級、區(qū)域級、省級、市級等不同空間尺度上的土壤水分反演、旱情監(jiān)測、作物產(chǎn)量估計[1-3,15-21]等相關(guān)研究，并取得一定成果，但關(guān)于黃土高原地區(qū)的土壤水分反演研究較少；另一方面，隨著Landsat OLI/TIRS系列數(shù)據(jù)的開放，為區(qū)域土壤水分的遙感監(jiān)測和反演提供了新的數(shù)據(jù)源，但目前以Landsat 8為數(shù)據(jù)源進行的土壤水分反演的相關(guān)研究較少。

黃土高原地區(qū)由于長期的水力、重力等自然作用的侵蝕以及人類不合理的開墾，水土流失嚴(yán)重，土壤貧瘠，有效養(yǎng)分缺乏，基礎(chǔ)肥力低，導(dǎo)致植被和生態(tài)系統(tǒng)嚴(yán)重退化[22-23]。延河流域是黃土高原典型丘陵溝壑區(qū)，一直以來也是黃土高原水土流失重點治理地區(qū)之一。土壤水分既是制約該地區(qū)植被恢復(fù)和重建的主要限制因子[24],也是決定土壤生產(chǎn)力的重要因素之一[25]，充分認識延河流域的土壤水分的空間分布狀況對該地區(qū)的生態(tài)系統(tǒng)的發(fā)展、植被恢復(fù)和重建都具有重要的意義。故本研究選擇黃土丘陵區(qū)延河流域作為研究區(qū)，以2013年的延河流域Landsat8 OLI遙感影像為數(shù)據(jù)源，利用TDVI對研究區(qū)進行土壤水分狀況反演，分析探究該流域各土層的土壤水分分布狀況，為該流域的土壤水分監(jiān)測與生態(tài)恢復(fù)提供依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

延河流域(36°23'～37°17' N，108°45'～110°28' E)發(fā)源于陜西省靖邊縣天賜灣鄉(xiāng)的周山，屬黃河中游河口鎮(zhèn)-龍門區(qū)間的一級支流。流域總面積7 687 km2，干流總長度286.9 km，河網(wǎng)密度約4.7 km·km-2，海拔為495～1 795 m。該區(qū)氣候?qū)倥瘻貛О敫珊荡箨懠撅L(fēng)氣候，年平均氣溫9℃，多年平均降水量為495.6 mm，降雨多集中在7—9月份。延河流域?qū)儆邳S土丘陵溝壑區(qū)第二副區(qū)，地形破碎，主要土壤類型為黃綿土，土質(zhì)比較疏松，抗蝕抗沖性差，土壤侵蝕嚴(yán)重，生態(tài)環(huán)境脆弱[26-28]。該地區(qū)植被類型隨環(huán)境梯度變化明顯，從南向北主要物種有刺槐(Robiniapseudoacacia)、遼東櫟(Quercuswutaishanica)、檸條(CaraganaKorshinskii)、白羊草(Bothriochloaischaemum)、百里香(Thymusvulgaris)、長芒草(StipaBungeana)等。

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

1.2.1數(shù)據(jù)源與預(yù)處理 本研究采用2013年9月6日和2013年9月13日的Landsat8 OLI影像數(shù)據(jù)，均為雨季后影像，分辨率30×30 m，且云量≤2%。利用ENVI5.3軟件對影像進行預(yù)處理，主要包括輻射定標(biāo)、FLAASH模型大氣校正(能見度參數(shù)來自https://www.wunderground.com/)、幾何校正、鑲嵌和裁剪等，其中幾何校正誤差在0.5個像元之內(nèi)，精度符合研究標(biāo)準(zhǔn)。本研究基于降雨、溫度、植被分帶等因子在延河流域均勻選取10個野外實測區(qū)，包括張家灘、化子坪、譚家營、安塞、招安、縣南溝、高橋、萬花山、甘谷驛、張家灘。2013年6月雨季后開始進行野外測定，采用烘干法來測定土壤水分。土壤水分采樣點為178個，選取148個作為土壤水分反演模型回歸點，30個作為驗證點。

圖1 土壤水分野外驗證點分布圖Fig.1 Distribution of soil moisture verification point

1.2.2植被覆蓋度 (1)計算NDVI的值

歸一化植被指數(shù)(NDVI)是地表植被覆蓋狀況的遙感指標(biāo)，可以將它定義為近紅外波段的反射值與紅光波段的反射值之差比上兩者之和，歸一化植被指數(shù)(NDVI)的公式如下：

NVDI=(ρNIR-ρRED)/(ρNIR+ρRED)

(1)

式中，為band5的反射率即近紅外波段，為band4的反射率即紅光波段。

(2)基于像元二分模型提取植被覆蓋度，公式如下：

FC=(NDVI-NDVIsoil)/(NDVIveg-NDVIsoil)

(2)

式中，NDVIveg為純植被像元的NDVI值，NDVIsoil為純裸土像元的NDVI值。

1.2.3地表溫度 采用輻射傳導(dǎo)方程計算得到同溫度下的黑體輻射亮度后，根據(jù)普朗克定律反函數(shù)，計算得到真實的地表溫度。公式如下：

(3)

式中，Ts為真實地表溫度，K1和K2均為常量校訂系數(shù)，K1為774.89W/(m2·μm·sr)，K2為1321.08 K。

1.2.4植被覆蓋干旱指數(shù) Nemani和Moran等[29-30]發(fā)現(xiàn)Ts和NDVI的散點圖為梯形，即Ts-NDVI特征空間。Sandholt等[12]在Moran等提出的Ts-NDVI梯形特征空間的基礎(chǔ)上將特征空間從梯形簡化為三角形，提出了溫度植被干旱指數(shù)TVDI(Temperature-Vegetation Dryness Index)的概念，公式為：

TVDI=(Ts-Tsmin)/(Tsmax-Tsmin)

(4)

Tsmin=a1+b1*NDVI

(5)

Tsmax=a2+b2*NDVI

(6)

式中，Ts為地表溫度，Tsmin(濕邊)為相同NDVI條件下對應(yīng)的最小地表溫度(℃)，Tsmax(干邊)為相同NDVI條件下對應(yīng)的最大地表溫度(℃)。將公式干邊和濕邊帶入TVDI的計算模型中，得到TVDI的計算公式：

(7)

式中，a1、b1、a2、b2為Ts-NDVI特征空間干濕邊擬合方程的斜率與截距。

TVDI值的范圍在[0，1]之間，干邊所對應(yīng)TVDI值為1，濕邊所對應(yīng)TVDI值為0，TVDI的值隨著地表溫度的升高而增大，土壤水分隨TVDI值增大而減小，因而土壤水分與溫度、TVDI值都呈負相關(guān)關(guān)系。

2 結(jié)果與分析

2.1 Ts-NDVI特征空間及干濕邊方程

本文選取IDL編程語言，在ENVI5.3-IDL環(huán)境下利用ZONALMAX、ZONALMIN函數(shù)計算TVDI特征空間的干濕邊擬合方程，NDVI作為橫坐標(biāo)(自變量X)，Ts作為縱坐標(biāo)(因變量Y)。設(shè)置NDVI的步長值為0.01，將NDVI(0～1)等分100份，利用NDVI單波段影像查找固定范圍值的索引，分別提取相同NDVI值對應(yīng)的地表溫度的最大值和最小值。在NDVI大于0時，隨著NDVI的增大，地表溫度的最大值呈減小趨勢；同時隨著植被覆蓋度增加，地表溫度的最小值呈增大趨勢，并且地表溫度的最大值、最小值與NDVI呈線性關(guān)系。

根據(jù)上述干、濕邊的計算模型，利用2013年影像數(shù)據(jù)Ts-NDVI干、濕邊特征空間的趨勢擬合方程，得到2013年影像數(shù)據(jù)的干邊擬合回歸方程為y=-25.9840x+303.542(R2=0.90)；濕邊擬合回歸方程為y=13.9418x+268.711(R2=0.75)。其中，濕邊的斜率大于0，干邊的斜率小于0，干、濕邊斜率的變化可以歸因于蒸散、冠層傳導(dǎo)度以及土壤水分含量。從線性擬合的效果來看，干邊線性擬合的效果更好。

2.2 延河流域土壤水分的反演

2.2.1土壤表層TDVI指數(shù) 延河流域溫度植被干旱指數(shù)反演結(jié)果表明(圖2)：微旱和干旱的分布范圍廣且多分布于延河流域上游杏子河和延河主河道上游，王堯水庫和化子坪西北部，TVDI的數(shù)值范圍主要集中在0.6～0.8；濕潤地區(qū)的分布范圍較小且多分布于延安以南的南川河流域等地，呈現(xiàn)邊緣式分布方式，TVDI的數(shù)值范圍為0～0.5。

不同的地表覆蓋類型(圖3)，土壤水分存在很大的差異，各地表覆蓋類型的TVDI呈現(xiàn)出不同的分布格局，總體表現(xiàn)趨勢為：低覆蓋度>中覆蓋度>

高覆蓋度。延河流域西北地區(qū)的植被覆蓋度低，水土流失嚴(yán)重，土壤的蓄水性能較差，土壤比較干旱；在城市建成區(qū)及流域的沿岸地帶，多為城市建筑用地土壤透水性差，呈干旱狀態(tài)；在植被覆蓋度良好的區(qū)域，土壤蓄水性能好，土壤水分狀況正常。

2.2.2不同深度土壤水分反演回歸模型 將5個土層深度(0～20 cm，20～40 cm，40～60 cm，60～80 cm，80～100 cm)各回歸樣點與對應(yīng)的TVDI作線性回歸分析，得各個土層深度的線性回歸方程(表1)。

圖2 延河流域TVDI空間分布圖Fig.2 Spatial distribution of TVDI

圖3 延河流域植被覆蓋度Fig.3 Yanhe River basin Vegetation Coverage

表1 各土層深度TVDI回歸方程Table 1 TVDI regression equation for different soil depth

各土層深度的實測土壤水分與TVDI的反演值具有較好的負相關(guān)性，這表明TVDI指數(shù)適用于延河流域的土壤水分信息反演。其中0～20 cm，20～40 cm，40～60 cm土層的土壤水分與TVDI反演值的相關(guān)系數(shù)分別為0.65、0.78和0.64，其中20～40 cm土層的土壤水分實測值與TVDI反演值的相關(guān)性較強；60～80 cm，80～100 cm土層的相關(guān)系數(shù)分別為0.52和0.37，其中60～80 cm土層實測土壤水分與TVDI反演值相關(guān)性較弱，80～100 cm土層的相關(guān)性最差，標(biāo)準(zhǔn)差最高。

2.2.3不同深度土壤水分反演結(jié)果 0～20 cm深度的土壤含水量最低，20～40 cm深度的土壤含水量最高，40～60 cm深度的土壤含水量相對較高，60～80 cm，80～100 cm深度的土壤含水量較低。其中：0～20 cm土層(圖4A)小于5%的土壤含水量在西北和東南部的邊緣地區(qū)呈現(xiàn)較小的片狀分布；5%～7%和7%～13%的土壤含水量在延河流域西北部、中部、東南部均有大面積分布；大于18%的土壤含水量主要分布在水體和南部的林地區(qū)域。20～40 cm土層(圖4B)土壤含水量小于5%只在西北和東南部的邊角地帶有所分布；5%～7%的土壤含水量分布區(qū)域大致與0～20 cm中5%～7%的土壤含水量的分布區(qū)域一致，但其分布范圍大大縮小；7%～13%和13%～18%土壤含水量的分布范圍擴大，其主要分布在延河流域西北部、中部、東南部；大于18%的土壤含水量在中部和東南部的分布范圍明顯擴大。40～60 cm土層(圖4C)小于5%的土壤含水量主要分布在西北部的邊角地區(qū)；5%～7%的土壤含水量在西北和東南部地區(qū)呈零星片狀分布；7%～13%的土壤含水量在整個延河流域均有分布；大于18%的土壤含水量主要分布在水體和南部的林地區(qū)域，中部和東南部。60～80 cm(圖4D)和80～100 cm(圖4E)土層的土壤含水量分布圖相對一致，小于5%和5%～7%的土壤含水量主要分布在西北部、南部和東南部的邊緣地區(qū)，分布范圍較廣；7%～13%和13%～18%的土壤含水量分布范圍縮?。淮笥?8%的土壤含水量的分布范圍明顯縮小，其集中分布在南部邊緣地區(qū)和流域水體區(qū)域。

2.3 精度檢驗

利用未參與回歸模型建立的30個野外采樣點數(shù)據(jù)對提取的各深度土層的土壤含水量反演值進行精度驗證，選取的驗證點均勻的分布在周家山小流域、陳家坬小流域、毛堡則小流域，精度驗證(表2)如下：

圖4 延河流域不同深度土壤含水量分布圖Fig.4 Distribution of soil water content at different depths in the Yanhe River Basin

表2 各層土壤含水量反演值精度驗證Table 2 Accuracy verification of inversion value of soil moisture in each layer

注：表中Sc代表土壤含水量實測值，F(xiàn)y代表土壤含水量反演值，Jd代表精度

Note：Sc represents the measured value of soil water content,Fy represents the inversion value of soil water content,and Jd represents the accuracy

0～20 cm，20～40 cm，40～60 cm，60～80 cm和80～100 cm土層深度的土壤水分反演值的平均相對精度分別為70.2%，80.6%，72.6%，67.3%和55.4%，其中20～40 cm土層深度的土壤水分反演值的相對精度最高；60～80 cm的土壤水分反演值的整體精度最低，且波動性最大。

結(jié)果表明TVDI法能夠?qū)\層土層的土壤水分進行有效的監(jiān)測，特別是能夠更穩(wěn)定反映和指示20～40 cm土層的土壤水分狀況，但是不能實現(xiàn)對更深土層土壤水分信息提取，這與夏燕秋等[1]、王純枝等[31]和姚春生等[32]的研究結(jié)果一致。

3 討論

土壤水分是黃土高原關(guān)鍵帶中的關(guān)鍵要素，同時又是黃土高原植被恢復(fù)與生態(tài)環(huán)境重建的決定因子[33-36]。故本研究基于Landsat8 OLI數(shù)據(jù)，選用溫度植被干旱指數(shù)法對延河流域的土壤水分進行反演，結(jié)果表明延河流域微旱和干旱的分布范圍廣且多分布于西北地區(qū)，分析陜西省水利廳網(wǎng)站雨情簡報中的降雨數(shù)據(jù)可以發(fā)現(xiàn)[37]，2013年暴雨分布和降雨空間分布表現(xiàn)出由南到北逐漸降低的趨勢，這與延河流域的地表土壤水分的空間分布規(guī)律有很大的一致性，已有研究也表明了降水格局是影響延安-安塞地區(qū)土壤水分空間格局主導(dǎo)因素[38]。

利用植被干旱指數(shù)法進行土壤水分反演，其反演精度在20～40 cm最高，這與夏燕秋等[1]、王純枝等[31]和姚春生等[32]的研究結(jié)果一致。該地區(qū)近地表的表土層是降水、土壤水、熱量、生物利用的一個交換面(作用面)，除受上述因素影響外，植被、土壤、地形等對土壤水分動態(tài)變化也有著重要影響，故淺層土壤水分不穩(wěn)定，變動較大[22,39]。由于土層較深，土壤水分對光學(xué)影像熱紅外波段的敏感度下降，與TVDI指數(shù)的相關(guān)性較弱，導(dǎo)致深層土壤的土壤水分反演精度也較低[1]。

目前，利用遙感影像進行土壤水分的遙感監(jiān)測和反演已經(jīng)取得一定的成果，但單期的遙感影像的研究結(jié)果僅能反映特定時期的土壤水分狀況，而未能反映土壤水分的動態(tài)變化。因此，今后對延河流域土壤水分進行長時間序列的動態(tài)監(jiān)測和反演時，應(yīng)考慮多種傳感器及多時相的遙感影像。

4 結(jié)論

延河流域微旱和干旱的分布范圍廣且多分布于西北地區(qū)，不同的地表覆蓋類型，土壤水分的總體趨勢表現(xiàn)為：低覆蓋度>中覆蓋度>高覆蓋度；各土層深度的實測土壤水分與TVDI的反演值具有較好的負相關(guān)性，表明TVDI指數(shù)法適用于延河流域的土壤水分反演；土壤水分反演值的精度檢驗結(jié)果中，20～40 cm的土壤水分反演值為80.57%，整體精度最高，表明TVDI更能穩(wěn)定反映和指示地表20～40 cm土壤水分狀況，但對于提取更深土層深度的土壤水分信息不能滿足精度要求。