王漢文，原喜忠，雷勝友，郭長志，褚志成，陳虹宇，陳鵬輝

(長安大學(xué) 公路學(xué)院，陜西 西安 710064)

0 引 言

在地球系統(tǒng)中，土壤濕度狀態(tài)在陸氣能量交換中起到了重要的控制作用。傳統(tǒng)依靠地面站點統(tǒng)計所計算得到的干旱濕度指數(shù)雖然可以較準(zhǔn)確地表達(dá)附近地區(qū)的干濕情況，但是獲取數(shù)據(jù)成本較高，并且往往數(shù)據(jù)的數(shù)量較少，分布存在較大的差異性，在時間序列方面，很難表達(dá)土壤濕度的動態(tài)變化。隨著遙感衛(wèi)星技術(shù)不斷發(fā)展，并伴隨諸多相關(guān)研究都表明其反演得到的濕度情況可以較為準(zhǔn)確地反映實際情況，遙感以其宏觀、動態(tài)、高效、受地面條件限制少等優(yōu)點為開展區(qū)域尺度下大面積土壤濕度監(jiān)測提供了有效手段[1-2]。在之前的研究中，大多數(shù)學(xué)者是通過建立植被指數(shù)(NDVI)與地表溫度(Ts)的綜合模型進(jìn)行土壤濕度監(jiān)測，比如溫度植被指數(shù)(TVI)、植被供水指數(shù)(VSWI)、水分虧缺指數(shù)(WDI)、條件植被溫度指數(shù)(VTCI)和溫度植被干旱指數(shù)(TVDI)等[3]，其中VTCI與TVDI運用較為廣泛，它們均與表層土壤濕度有較好相關(guān)性[4-5]，并且反演速度較快，因此本文采用其中之一的TVDI進(jìn)行相關(guān)研究。最新發(fā)射的Landsat 8衛(wèi)星比以前的Landsat系列衛(wèi)星在掃描方式、波段設(shè)置、輻射分辨性能等方面都有較大改進(jìn)[6]，并且由于其所測數(shù)據(jù)多為開源資料，所以它所提供的影像已成為現(xiàn)階段多數(shù)土壤濕度研究的數(shù)據(jù)源。運用遙感影像進(jìn)行濕度反演以及旱情動態(tài)監(jiān)測的研究，之前的學(xué)者已經(jīng)對諸多地區(qū)進(jìn)行了探索，效果較為理想，但是針對具體地區(qū)不同土地類型一年內(nèi)的濕度變化趨勢，卻鮮有人分析。為此，本文使用Landsat 8 OLI/TIRS產(chǎn)品數(shù)據(jù)，采用溫度植被干旱指數(shù)(TVDI)對研究區(qū)一年內(nèi)的各個月的土壤濕度進(jìn)行反演，結(jié)合不同土地類型分析林地、草地等土地類型的年濕度變化趨勢。

1 資料與方法

1.1 研究區(qū)概況

論文研究的矩形區(qū)域位于廣東省北部的韶關(guān)市內(nèi)113°36′36″～114°04′12″E，24°30′10″～24°58′12″N，含武江區(qū)與曲江區(qū)兩大部分，研究區(qū)域總面積約2 927.8 km2，其影像位置見圖1。

圖1 研究區(qū)Landsat 8遙感影像位置

該區(qū)域?qū)儆谥衼啛釒駶櫺图撅L(fēng)氣候，春季陰雨連綿，秋季降水偏少，冬季寒冷，夏季偏熱。年均降雨1 400～2 400 mm，年平均氣溫為18.8～21.6 ℃。研究區(qū)地形以丘陵為主，河谷盆地分布其中，海拔變化較大，為65～1 200 m。土地覆被主要有林地、草地以及耕地和市區(qū)的建筑等。韶關(guān)是全國重點林區(qū)，而且還擁有十分豐富的耕地資源，對其進(jìn)行土壤濕度監(jiān)測反演以及不同土地類型的濕度變化趨勢研究，對了解區(qū)域內(nèi)部不同植被的分布發(fā)育情況和保護(hù)森林、耕地資源，發(fā)展綠色生態(tài)有著極其重要的意義。

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

1.2.1 數(shù)據(jù)源與預(yù)處理

研究所使用的Landsat 8影像數(shù)據(jù)來源于地理空間數(shù)據(jù)平臺，其詳細(xì)信息可以參見表1。如果2017年內(nèi)某月沒有可使用的影像數(shù)據(jù)，也應(yīng)選取接近年份對應(yīng)月份的影像。將影像加載到軟件ENVI 5.3中進(jìn)行裁剪，幾何校正、輻射定標(biāo)以及FLAASH大氣校正等預(yù)處理，相關(guān)精度均符合研究標(biāo)準(zhǔn)。將處理后的影像進(jìn)行監(jiān)督分類，分類后處理以及精度檢驗之后得到分辨率為30 m的土地類型分類圖(圖2)，5種土地類型為草地、耕地、林地、水體、建筑用地及裸地，再對這5種用地類型進(jìn)行掩膜提取處理，并對TVDI適用的用地類型進(jìn)行分區(qū)統(tǒng)計分析，計算出相應(yīng)的平均值來表示這一地類的TVDI值。

表1 遙感數(shù)據(jù)信息

圖2 土地類型分類圖

1.2.2 歸一化植被指數(shù)

歸一化植被指數(shù)(NDVI)是反映土地覆蓋植被狀況的一種遙感指標(biāo)。該指數(shù)在一定程度上能消除地形和群落結(jié)構(gòu)的陰影、輻射干擾及太陽高度角和大氣所帶來的噪聲[6]。其被定義為近紅外波段的反射值與紅光波段的反射值之差除以兩者之和，公式為

(1)

式中：ρNIR為近紅外波段反射率；ρRED為紅光波段反射率。

NDVI的取值為-1～1，負(fù)值表示對可見光高反射的物質(zhì)，接近0表示地物為巖石、裸地等[7]，正值表示存在植被覆蓋，且值越大，覆蓋度越大。

1.2.3 地表溫度

地表溫度是反映地球表面能量狀態(tài)的一個重要指標(biāo)，在研究區(qū)干濕情況反演中，其作用是要與歸一化植被指數(shù)一同構(gòu)建Ts-NDVI特征空間[8]。本研究中使用輻射傳導(dǎo)方程法(大氣校正法)來計算地表溫度，其原理是利用衛(wèi)星拍攝影像的熱輻射強(qiáng)度減去大氣對地表的熱輻射強(qiáng)度，最終得到地表的熱輻射強(qiáng)度，由于不同物體表面的比輻射率不同，通過確定不同地表的比輻射率，之后根據(jù)普朗克公式即可反算出地表的溫度[9]。輻射傳輸方程的表達(dá)式為

Lλ=[ε·B(Ts)+(1-ε)L↓]·τ+L↑， (2)

式中：Lλ為Landsat 8的熱紅外輻射亮度值；ε為地表比輻射率；B(Ts)為黑體熱輻射強(qiáng)度，是地表溫度Ts的函數(shù)；τ為大氣在熱紅外波段的透過率；L↑,L↓分別為大氣上行輻射強(qiáng)度與大氣下行熱輻射強(qiáng)度。τ,L↑,L↓這3個大氣剖面參數(shù)可通過NASA提供的網(wǎng)站查詢得到。公式變換得到B(Ts)表達(dá)式，即：

B(Ts)=[Lλ-L↑-τ·(1-ε)L↓]/(τε)， (3)

將其代入普朗克公式，得到地表溫度值

(4)

式中：K1，K2為熱紅外波段的定標(biāo)常數(shù)，K1=774.89 W/(m2·μm·sr)，K2=1 321.08 K。

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

溫度植被干旱指數(shù)是利用Ts-NDVI特征空間提取的水分脅迫指標(biāo)，它是估算陸面表層土壤濕度的一種重要方法[3，10]。對于某一區(qū)域而言，隨著植被覆蓋變化，植物蒸騰所轉(zhuǎn)化為的潛熱與顯熱能量就處于動態(tài)平衡之中，進(jìn)而表現(xiàn)為地面溫度的動態(tài)變化，又由于土壤的濕度狀態(tài)與地面溫度存在關(guān)系，于是就能間接反演出土壤的濕度狀態(tài)。在之前的研究中，S.M.Moran等[11]發(fā)現(xiàn)Ts-NDVI特征空間呈梯形，而I Sandholt等[10]則在此基礎(chǔ)上將其簡化為三角形，并提出了溫度植被干旱指數(shù)TVDI(Temperature Vegetation Dryness Index)[12]的概念，計算式為

(5)

式中：TS為影像中任意像元對應(yīng)的地表溫度；TSmax，TSmin分別為影像中同一NDVI值對應(yīng)的最高與最低地表溫度。

在簡化的特征空間中，對干邊(TSmax)、濕邊(TSmin)進(jìn)行線性擬合，即

TSmin=a1+b1·NDVI,

(6)

TSmax=a2+b2·NDVI,

(7)

式中：a1，b1為TVDI濕邊方程的系數(shù)；a2，b2為TVDI干邊方程的系數(shù)。將式(6)～(7)代入式(5)得

(8)

TVDI取值為[0，1]，其中濕邊的TVDI值為0，干邊的TVDI值為1，取值越接近0表示土壤濕度接近正常濕潤土壤，越接近1，表示土壤越干旱缺水[13-15]，故土壤濕度與TVDI取值呈負(fù)相關(guān)關(guān)系。以TVDI作為土壤濕度分級標(biāo)準(zhǔn)，將土壤濕度分為5個層次[16]，如表2所示。

表2 TVDI土壤濕度分級標(biāo)準(zhǔn)

2 結(jié)果與分析

2.1 NDVI結(jié)果分析

通過逐月計算得到研究區(qū)每月的NDVI空間分布圖，共12幅。均采用密度分級得到4個等級[14]，以10月份的NDVI空間分布圖為例(圖3)。經(jīng)與土地類型分類圖比對得到：NDVI值小于0的區(qū)域多為水體區(qū)域，經(jīng)像素統(tǒng)計后得到水域面積約占總面積的1.34%；NDVI值在0～0.3之間多為建筑用地、裸地以及稀疏草地等低植被覆蓋的區(qū)域，面積約占總面積的4.11%，其中0～0.2基本為建筑用地及裸地[5]；NDVI值在0.3～0.6之間主要是長勢較差的耕地、草地等較高植被覆蓋的區(qū)域，面積約占總面積的12.90%；對于NDVI值大于0.6的高植被覆蓋區(qū)域，主要包括林地以及長勢較好的耕地、草地，其約占總面積的81.65%，其中大于0.9的部分只占總區(qū)域面積的1.56%。NDVI空間分布效果基本可以較好地反映土地類型的基本情況。

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

2.2.1 特征空間中NDVI范圍的確定

NDVI<0對應(yīng)的地表類型為水體(可以認(rèn)定濕度為100%)，在特征空間構(gòu)建時不予討論[4]。根據(jù)相關(guān)研究[5，17-18]發(fā)現(xiàn)，當(dāng)植物覆蓋度較小時，NDVI值很難指示該地區(qū)的植物生物量，當(dāng)植物覆蓋度太大時，由于植物覆蓋密集，NDVI值增加延緩而呈現(xiàn)飽和狀態(tài)，對植物覆被檢測的靈敏度將降低，無法準(zhǔn)確顯示植被覆蓋情況，故根據(jù)NDVI像元分布直方圖及多次擬合嘗試，對于0≤NDVI<0.2及NDVI>0.9的值予以刪除。綜上所述，構(gòu)建特征空間的NDVI取值為0.2～0.9。該區(qū)間的選擇也與前人的研究[5，18]接近。

2.2.2 特征空間構(gòu)建及干濕邊方程擬合

將反演得到的TS以及NDVI數(shù)據(jù)導(dǎo)入ArcGIS軟件中，并提取出選定區(qū)間內(nèi)某一NDVI及其所對應(yīng)的最大和最小地表溫度TSmax，TSmin。以一定步長的NDVI為橫坐標(biāo)(自變量x)，反演得到的地表溫度為縱坐標(biāo)(應(yīng)變量y)，構(gòu)建了一年內(nèi)各個月份的Ts-NDVI特征空間，并擬合出了各個時相的干濕邊方程，如圖4所示。

圖3 10月份NDVI等級分布

從圖4可以看出，不同時相的Ts-NDVI散點圖基本都具有相似的形狀，呈三角形或梯形，與前人研究[19-21]結(jié)論相同。總體而言，地表溫度最大值TSmax隨NDVI增大呈減小的趨勢，地表溫度最小值TSmin，除幾個月因刪除了NDVI大于0.9后的散點圖導(dǎo)致有減小的趨勢外，其他多數(shù)均有增大的趨勢，并且地表溫度的最大與最小值的差值逐漸縮小。

根據(jù)Ts-NDVI特征空間，擬合得到各個時相的干濕邊方程，見表3，其中干邊的斜率均小于0，濕邊斜率大多大于0，這也間接反映了地表溫度隨NDVI值變化的趨勢情況。表3中干邊方程的比例系數(shù)的絕對值均大于濕邊方程的比例系數(shù)的絕對值，說明對該地區(qū)而言，干邊所擬合的地表溫度對植被指數(shù)的變化更加敏感。從擬合效果而言，干邊的擬合效果較濕邊整體較好，R2的平均值達(dá)到了0.837。

2.3 TVDI等級圖分布

根據(jù)表3中的干濕邊方程，計算得到不同時相的TVDI值，以TVDI值作為土壤濕度分級的標(biāo)準(zhǔn)，結(jié)合表2將其分為5個層次，由于在構(gòu)建特征空間時未考慮NDVI<0.2及NDVI>0.9的部分，為減少土壤濕度分布誤差，故再疊加水體、建筑用地及裸地并將NDVI>0.9的區(qū)域歸為無數(shù)據(jù)區(qū)，進(jìn)而得到本研究區(qū)年內(nèi)各個月的土壤濕度分布(圖5)。

2.3.1 土壤濕度空間分布特征

對于研究區(qū)域的土壤濕度狀況，整體呈現(xiàn)出西側(cè)較為干旱，東側(cè)較為濕潤。結(jié)合土地類型分類圖可以看出，東側(cè)大多為林地或者草地，地表的黏性土層也相對西側(cè)較厚，故儲存水分的能力相對較高，相應(yīng)地，植被覆蓋指數(shù)圖也能反映這一方面；另一方面，東側(cè)的地勢較高且多為山地丘陵，反演得到的地表溫度較低，植物的蒸發(fā)蒸騰能力較低，進(jìn)而地表的土壤濕度相應(yīng)較高。局部區(qū)域可以看出，建筑用地及裸地周圍的土壤濕度明顯較低，林地密集的區(qū)域土壤濕度則相對較高，草地和耕地之間的土壤濕度比較與時間有一定的關(guān)系，主要是存在人為因素導(dǎo)致耕地土壤濕度在某些特定時間段內(nèi)發(fā)生較大變化。而對于山地丘陵區(qū)域，通過提取陰陽坡位置處的土壤濕度發(fā)現(xiàn)，相同高程陰坡的TVDI值要比陽坡低0.2左右，即濕度大約相差一個層次。

2.3.2 土壤濕度時間分布特征

土壤濕度變化在一年內(nèi)大體呈現(xiàn)出先降低再升高的循環(huán)，濕度最低即干旱程度最嚴(yán)重(極干旱)主要集中在4—6月，通過ArcGIS相關(guān)的分類統(tǒng)計工具，計算得到了4—6月極干旱區(qū)的面積比例分別為43.67%，30.62%，14.65%，其他月份極干旱區(qū)則均小于10%。研究區(qū)西側(cè)大片區(qū)域表現(xiàn)為極干旱的類型，結(jié)合土地類型、作物生長、氣候變化等因素可知，西側(cè)多為耕地及草地。農(nóng)作物一般4月播種，生長期多為4，5兩月，而草地的快速生長期一般也在4月，對水分需求都較大，加之該月份內(nèi)區(qū)域降水較少并且地表溫度也處于上升趨勢，所以導(dǎo)致了大面積極干旱的情況。為了直觀了解一年內(nèi)研究區(qū)的土壤濕度整體情況，按照計算的TVDI值，將表2中的極濕潤(0

圖4 年內(nèi)各月的TS-NDVI特征空間

表3 年內(nèi)各月Ts-NDVI特征空間干濕邊擬合方程

圖5 年內(nèi)各月份的TVDI土壤濕度等級分布

圖6 一年內(nèi)濕度變化整體趨勢Fig.6 Overall trend of humidity change of the year

2.4 不同用地類型的土壤濕度趨勢變化

結(jié)合之前的土地類型分類圖2以及一年內(nèi)逐月的TVDI土壤濕度等級分布圖5掩膜提取出各個用地類型，并對林地、耕地、草地每月的濕度情況進(jìn)行分區(qū)統(tǒng)計分析，計算得到TVDI的月平均值，進(jìn)而繪制出這3種不同用地類型的土壤濕度趨勢如圖7所示。

從圖7可以看出，幾種用地類型的土壤濕度情況在一年的變化整體表現(xiàn)為：3月，TVDI值迅速增大，土壤濕度減??；4—8月濕度基本維持不變；9月，TVDI值又減小，相應(yīng)土壤濕度則呈現(xiàn)增大趨勢；之后的10月至來年2月各個用地類型的TVDI值又基本保持恒定。這一規(guī)律也基本符合該地區(qū)植被在生長周期內(nèi)與土壤水分的相互作用關(guān)系。局部區(qū)域，各個用地類型土壤濕度的情況又可以分為兩大階段：第一階段是3—10月，土壤濕度情況表現(xiàn)為林地>草地>耕地；第二階段是一年內(nèi)的其余時間，表現(xiàn)為林地>耕地>草地。兩個階段中濕度變化主要發(fā)生在耕地與草地之間。對于耕地草地之間的變化，經(jīng)查閱資料得知，與該地區(qū)農(nóng)作物播種后因生長需要吸收較多土壤水分以及在該時間段內(nèi)降雨量大大減小有密切關(guān)系。對于林地，由于其植被覆蓋相對于其他類型一直較高，地表土壤多為黏性土且性質(zhì)穩(wěn)定，保水能力較強(qiáng)，故TVDI值相對較低，濕度一直比其他用地類型高。通過計算一年內(nèi)林地、草地、耕地TVDI值的變化幅度值(TVDImax-TVDImin)分別為0.226，0.241和0.324，得到3種用地類型的土壤濕度變化幅度情況為耕地>草地>林地。其中林地土壤濕度變化幅度最小，土壤水分的保持能力最強(qiáng)，耕地由于受到較大的人為因素影響，加之保水能力較差，土壤濕度變化幅度最大。

圖7 不同用地類型的土壤濕度趨勢

綜上，幾種用地類型在一年內(nèi)土壤濕度情況為林地>耕地或草地；土壤濕度變化幅度情況為耕地>草地>林地。幾種用地類型在同月中TVDI值差別不是很大，變化趨勢基本相同，且均符合大區(qū)域范圍下的土壤濕度變化情況，因此具有一定的參考性。

3 討 論

本研究討論了一年中每月的濕度情況，相對于其他研究而言，研究時期較長，更具有代表性和動態(tài)特性，但由于研究區(qū)內(nèi)Landsat 8一個月中可用的影像數(shù)據(jù)十分有限，用某期的影像進(jìn)行土壤濕度反演，進(jìn)而用于代表該月的土壤濕度情況，其中可能存在一些偶然性因素導(dǎo)致較大的土壤濕度誤差。

由于NDVI在低值區(qū)易受到土壤背景因素的影響，在高值區(qū)對高覆蓋植被反應(yīng)不敏感，故會導(dǎo)致TVDI的計算存在一定誤差[18]。雖然本文中采取排除的方法提高準(zhǔn)確性，但處理過程過于麻煩，故可以考慮用EVI[22-23]及SAVI[24-25]取代NDVI來構(gòu)建特征空間。

在使用TVDI進(jìn)行不同土地類型土壤濕度研究時受到了諸多條件限制，今后可以考慮使用其他相關(guān)性較好的濕度指數(shù)，例如VTCI等，用以驗證TVDI的準(zhǔn)確性并對比研究濕度指數(shù)對不同地類的適用性。

研究中對于用地類型的分類數(shù)量不多，且精度不高。例如未對林地中的二級分類進(jìn)行詳細(xì)區(qū)分，今后可以考慮選取地類豐富的區(qū)域，進(jìn)行更加詳細(xì)的分類與深入研究。

4 結(jié) 論

本文利用Landsat 8相關(guān)數(shù)據(jù)對廣東省北部韶關(guān)市部分區(qū)域進(jìn)行了土地類型分類及土壤濕度監(jiān)測反演研究，得到了一年內(nèi)各個月具有代表性的歸一化植被指數(shù)NDVI、地表溫度TS和溫度植被干旱指數(shù)TVDI等數(shù)值，并結(jié)合土地類型分類，得到了以下結(jié)論。

(1)對于研究區(qū)的土壤濕度空間分布，整體呈現(xiàn)西側(cè)干旱，東側(cè)濕潤的情況，局部區(qū)域與人為因素、氣候條件以及地形地貌存在一定聯(lián)系。

(2)對于研究區(qū)土壤濕度時間分布，在一年內(nèi)大體呈現(xiàn)先降低至基本恒定，再升高并保持穩(wěn)定的循環(huán)。濕度最低主要集中在4—6月。結(jié)合TVDI將濕度情況分為潮濕、干燥兩大類，4—9月，干燥區(qū)面積大于潮濕區(qū)面積，一年中其他月份，潮濕區(qū)面積大于干燥區(qū)面積。

(3)基于TVDI，幾種用地類型在一年內(nèi)土壤濕度情況，3—10月，林地>草地>耕地，其余時間為林地>耕地>草地，其中耕地與草地的濕度大小情況與研究區(qū)內(nèi)人為因素及氣候關(guān)系密切。通過計算TVDImax-TVDImin，得到一年內(nèi)幾種用地類型的土壤濕度變化幅度情況為耕地>草地>林地。