陳明星, 張玉虎
(首都師范大學 資源環境與旅游學院, 北京 100048)
土壤濕度是水文學、氣象學以及農業科學研究領域的一個重要參數。大型農灌區土壤濕度狀況的監測是農業用水管理以及農作物旱情預報的一個重要內容[1],一直以來受到學者們的廣泛關注[2-4]。遙感手段相對于傳統的土壤濕度監測手段,能夠快速獲取大面積區域的土壤濕度狀況,對于農業旱情監測具有較大的實用價值[5]。國內外學者對土壤濕度的遙感監測方法做了許多研究,提出了不同的模型和方法,總體上分為光學遙感和微波遙感兩種方式,主要有:表觀熱慣量法[6-7]、溫度植被干旱指數[8-9]、植被供水指數[10-11]、微波極化差異指數[12]、積分方程模型[13]、Dobson模型[14]等。微波遙感方式僅適合小尺度土壤濕度反演[15],光學遙感手段中部分方法在實際應用中也具有一定的局限性,如表觀熱慣量方法僅適合于裸土地區土壤濕度的估算[16]。在眾多的模型與方法中,通過LST-NDVI 特征空間計算TVDI反演土壤濕度的方法,由于較強的適用性在國內外得到了廣泛的研究與關注。齊述華等[17]、曹雷等[18]研究表明TVDI可有效反演區域土壤水分,且精度較高,J Chen[19]和Li Z等[20]用TVDI方法分別研究了我國黃淮海平原和陜北半干旱黃土高原的土壤濕度狀況,均得到了較好的效果;N.T.Sona等[21]探討了利用TVDI反演的土壤濕度結果對農業干旱監測的適用性,驗證結果顯示該方法效果較好。NR Patel等[22]評估了TVDI監測土壤濕度狀況的潛力,發現TVDI在植被覆蓋稀疏時能夠較好的監測土壤濕度的時間變化。當前大多數研究在TVDI模型計算中僅采用歸一化植被指數,但是歸一化植被指數對于不同時期植被覆蓋差異的敏感性不同,影響了TVDI監測土壤濕度的準確性。
本文考慮對植被覆蓋度敏感性不同的4種植被指數,以三江平原為研究區,計算不同植被指數下的TVDI,并用土壤相對濕度和降水量數據對TVDI監測土壤濕度狀況的效果進行驗證。對比不同植被指數下TVDI模型對土壤濕度狀況的監測效果,確定反演土壤濕度效果最佳的植被指數,為三江平原地區土壤濕度的監測提供參考依據。
三江平原位于中國黑龍江省東北部,介于43°50′02″—48°24′41″N,129°11′49″—134°46′37″E,地理位置如圖1所示。是由松花江、黑龍江和烏蘇里江沖積形成的低平原,土地面積為1.09×105km2,地勢總體特征是西南高東北低。本區屬于溫帶季風氣候,夏季炎熱多雨,冬季寒冷干燥。1月份平均氣溫低于-18℃,7月平均氣溫21~22℃,年降水量500~650 mm,且集中在夏季,屬于濕潤、半濕潤氣候區。行政區域包括佳木斯、雞西、鶴崗、雙鴨山、七臺河等地級市以及牡丹江市所屬的穆棱縣和哈爾濱所屬的依蘭縣,共計23個縣(市、區)。三江平原是我國重要的糧食產區,主要農作物為水稻、玉米、大豆等。

圖1 研究區位置
溫度植被干旱指數(Temperature Vegetation Dryness Index,TVDI)最早是由Sandholt等[23]提出,TVDI主要考慮了歸一化植被指數(NDVI)和地表溫度(LST)這兩個描述土壤表層特征的重要參數[24],以此構建NDVI--LST的特征空間,Sandholt等認為NDVI-LST的特征空間中存在一系列反映土壤濕度的等值線,據此提出了溫度植被干旱指數(TVDI)的概念(圖2)。其公式[25-26]為:
(1)
式中:LSTmin為NDVI相同值對應的最低地表溫度,為NDVI-LST特征空間的濕邊,LSTmax為NDVI相同值對應的最高地表溫度,為NDVI-LST特征空間的干邊。根據影像像元構造的特征空間,通過線性擬合獲得特征空間中的干邊和濕邊方程,其公式[25]為:
LSTmax=α1+b1×NDVI
(2)
LSTmin=α2+b2×NDVI
(3)
式中:α1,b1,α2,b2分別是特征空間中擬合的干邊和濕邊方程的系數。TVDI值在0~1之間,在干邊上TVDI為1,在濕邊上TVDI為0。TVDI越趨向于0,表示土壤濕度越高;TVDI越趨向于1,表示土壤濕度越低。

圖2 LST-NDVI特征空間
本研究采用的遙感數據來自于美國USGS網站提供的MODIS Terra產品數據(https:∥lpdaac.usgs.gov/dataset_discovery/modis),包括每8 d合成的空間分辨率為1 000 m的地表溫度產品(MOD11A2)和空間分辨率為500 m的地表反射率數據(MOD09A1)。數據時相覆蓋了2013年5—9月不同作物生長期的植被覆蓋情況(東北地區植被覆蓋度在一年內呈先上升后下降趨勢,最高的時間段為7—8月份[27])。本次研究使用的氣象數據和土壤濕度數據來源于中國氣象數據網(http:∥data.cma.cn),時間為2013年5月至9月。研究區共有12個土壤墑情站和22個氣象站(圖1),選取12個與土壤墑情站對應的氣象站。土壤墑情站提供了每旬的20 cm深度的土壤相對濕度數據,氣象站提供了逐日降水量數據。遙感影像的選取盡可能保證與土壤相對數據時間上的對應,并覆蓋5—9月的每個月份。
將MOD09A1,MOD11A2影像進行預處理后得到500 m分辨率的地表反射率和地表溫度數據,反射率數據包括:紅光波段(Red)和近紅外(NIR)波段和藍光波段(Blue)波段。使用ArcGIS軟件分別計算得到歸一化植被指數(NDVI)、增強型植被指數(EVI)、改進的修正土壤調整植被指數(MSAVI)以及比值植被指數(RVI),在ENVI中采用基于IDL語言開發的TVDI計算工具,輸入植被指數和地表溫度數據獲取構建特征空間的植被指數和地表溫度值及計算出的TVDI數字影像。不同植被指數[25,28]的計算公式見表1。在預處理過程中,由于8月份的遙感影像云量較大,經去云處理后影像缺失范圍較大,故不做分析。

表1 不同植被指數計算公式
注:ρNIR為近紅光波段反射率,ρRED為紅光波段反射率,ρBLUE為藍光波段反射率。
TVDI在計算過程中一般不考慮水體存在的區域,利用植被指數值將這部分區域剔除。為了進一步降低誤差,根據像元統計及目視解譯,結合不同時期的植被指數值,選取不同的閾值(5月為0.04、6月為0.05、7月為0.1、9月為0.1)對臨湖、臨河或河中沙洲等距離水體較近,但由于植被指數計算誤差存在少量水體或土壤水分飽和的裸地所在區域的少量像元進行剔除。在Excel軟件中,以植被指數為橫坐標,LST值為縱坐標,構建不同植被指數的LST-Ⅵ的特征空間,結果如圖3。
根據圖3,特征空間形狀總體相似,地表溫度的最大值均表現為隨著植被指數值的增大,先有輕微上升之后逐漸下降,但總體上隨植被指數增加呈減小趨勢,上升區間主要處于植被指數小于0.2區間,這部分主要為裸地或有稀疏植被分布的地區,由于植被覆蓋度較低,很難對區域內的地表溫度產生調控作用。地表溫度的最小值總體上呈現隨植被指數的增加而增大的趨勢,但是上升速率較小,地表溫度最小值所在區域可能代表了濕地、徑流量不大的水系等[26],這些區域的地表溫度在自然條件下處于比較穩定的水平。

圖3 2013年5-9月LST-Ⅵ特征空間
同一時期,不同植被指數構造的特征空間存在一定區別,LST-NDVI構造的特征空間在不同時期相較于其他植被指數較差,除5月外,特征空間的形狀與理論的三角形均存在差異。LST-RVI構造的特征空間在4種植被指數中效果最好,LST-EVI和LST-MSAVI構造的特征空間區別不大,效果均較好。特征空間總體上符合理論的三角形關系,與前人的相關研究結果一致[17-18],構造的特征空間可以用來獲取TVDI。
在理想狀態下,當地表覆蓋情況滿足從裸土到植被完全覆蓋均勻變化,土壤表層含水量從凋萎系數到田間持水量均勻變化時,特征空間的干、濕邊能擬合成直線狀[29]。根據特征空間的構造結果,擬合不同植被指數的特征空間對應的干邊和濕邊方程,情況見表2。
從表2可以看出,特征空間的干邊方程的擬合效果總體較好,R2大部分都在0.8以上;濕邊方程的擬合效果總體較差,R2大部分低于0.5,三江平原存在大面積的濕地及森林,可能是由于此類覆蓋物對地表溫度的調控作用所導致的。植被指數對于干、濕邊方程的擬合效果存在一定影響,總體上不同時期LST-RVI的干、濕邊方程擬合效果最好,其次為LST-EVI和LST-MSAVI,LST-NDVI的擬合效果相對其他3種植被指數最差,干、濕邊方程的R2均低于其他植被指數??傮w而言,干、濕邊方程的擬合結果能夠滿足TVDI計算的需要。

表2 不同植被指數的干邊和濕邊方程
2.3.1 TVDI與土壤相對濕度的回歸分析 根據三江平原地區12個土壤墑情站點測量的土壤相對濕度數據,以及由站點坐標原位提取的TVDI值,在SPSS中利用最小二乘法進行回歸分析??紤]到設置的觀測站點能夠代表相應地區的氣候特點,所以利用站點觀測數據進行驗證是可行的[17]?;貧w分析的結果見圖4。結果表明,4個時期不同植被指數計算的TVDI與土壤相對濕度均表現隨著TVDI值的上升,土壤濕度均呈下降趨勢,即二者之間為負相關關系,TVDI值越高,土壤濕度越低。
根據回歸分析結果,在植被覆蓋度較低的5月份,4種植被指數計算的TVDI與土壤相對濕度回歸分析的決定系數R2均在0.52左右,總體上相差不大,且均達到0.01的顯著水平。6月份植被覆蓋度有所增加,該時期R2均在0.2左右,植被指數的差異對于回歸分析影響較小。7月份植被覆蓋度較高,該時期計算的TVDI_EVI和TVDI_MSAVI與土壤相對濕度的回歸分析結果相較于TVDI_NDVI和TVDI_RVI有顯著優勢,R2明顯高于TVDI_NDVI和TVDI_RVI,且達到0.05的顯著水平,TVDI_NDVI和TVDI_RVI的回歸分析結果未達顯著水平。相關研究也表明EVI和MSAVI 對高植被覆蓋區域更為敏感,由此計算的TVDI比其他植被指數更能反映高植被覆蓋時的土壤濕度狀況[26,30]。9月植被覆蓋度處于較低水平,決定系數R2均在0.55左右,植被指數的差異對于回歸分析結果影響較小,該時間段TVDI反演的土壤濕度和實測的土壤相對濕度相關性均較高,且達到0.01的顯著水平。
對比4種植被指數計算的TVDI反演的土壤濕度和實測的土壤相對濕度的回歸分析結果。在低植被覆蓋的5月和9月,決定系數R2均顯著高于植被覆蓋度較高的6月和7月的分析結果,低植被覆蓋時的反演效果好于高植被覆蓋時。不同植被指數在5月、6月和9月份時R2相差較小,反演效果差異較??;在7月份,增強型植被指數(EVI)和修正土壤調節植被指數(MSAVI)的R2高于另外兩種植被指數,反演效果有所提升。
2.3.2 降水量與TVDI的相關分析 土壤濕度在一定程度上受降水量影響[31],研究中采用降水量探討其與TVDI和土壤相對濕度的關系,進一步確定TVDI方法監測土壤濕度狀況的準確性。主要是將降水量與土壤相對濕度和TVDI做相關分析,結果見表3。
根據相關分析結果,7月份的降水量與土壤相對濕度呈顯著正相關,即降水量增加,土壤相對濕度會增大,其他時間段均為不顯著正相關。綜合表明降水量對土壤相對濕度存在一定影響。
不同植被指數計算的TVDI和降水量的相關分析結果中,僅5月份的NDVI和EVI計算的TVDI通過了顯著性檢驗,其余均未通過顯著性檢驗。在5月、6月和9月期間,植被指數的不同對相關分析的結果并未造成明顯的差異,同一時期不同植被指數的相關系數值均較為接近;在7月份,基于EVI和MSAVI計算的TVDI值與降水量的相關系數明顯高于另外兩種植被指數??傮w情況與前文TVDI和土壤相對濕度的回歸分析結果一致,均表現為除7月外,其他月份不同植被指數計算TVDI驗證效果差異較小。降水量與不同植被指數計算的TVDI值的相關系數均為負值,即降雨量越大,TVDI值越小;結合降雨量與土壤濕度的相關分析結果(降雨量越大,土壤濕度越高),表明TVDI值越小,土壤濕度越大,與前面分析結果一致。

表3 2013年5-9月降水量分別與土壤相對濕度及TVDI相關系數
注:*表示在0.05水平上顯著相關。

圖4 2013年5-9月土壤相對濕度與不同植被指數TVDI的擬合效果比較
根據以上分析結果,MSAVI和EVI計算的TVDI在分別在5月和7月回歸分析R2最高,且特征空間擬合效果較好;NDVI計算的TVDI在6月和9月回歸分析R2最高,但特征空間擬合效果較差,EVI和RVI的R2分別僅次于NDVI且特征空間擬合效果較好。綜合考慮,在不同的時間選取相應的植被指數計算TVDI(5月為MSAVI,6月為EVI,7月為EVI,9月為RVI)并繪制土壤濕度等級圖,分析不同時間三江平原的土壤濕度狀況。根據TVDI值將研究區的土壤濕度狀況劃分為7個等級[32],分別為:極度濕潤(0~0.1)、中度濕潤(0.1~0.2)、輕度濕潤(0.2~0.4)、正常(0.4~0.6)、輕度干旱(0.6~0.8)、中度干旱(0.8~0.9)、極度干旱(0.9~1)。三江平原土壤濕度分布情況見圖5。

圖5 不同時期土壤濕度分布
由TVDI反映的三江平原地區的土壤濕度分布狀況可知,5月全區受旱范圍比較廣泛,主要分布在三江平原的西部,受旱地區多為輕度干旱,中部的集賢縣、西南部雞東縣以及密山市西部地區出現了中度干旱,同時南部地區存在小范圍的極度干旱。東部地區多為輕度濕潤或正常情況,土壤濕度狀況分區明顯,呈現顯著的從東北到西南由濕潤向干旱變化的特征。6月全區大部分區域土壤濕度狀況正常,中西部及西南部區域存在輕度干旱,中部存在零散的中度干旱,范圍較小,主要集中在富錦市、寶清縣和集賢縣。7月全區大部分地區為正常,中部及南部地區存在零散輕度干旱,中部的富錦市部分地區存在中度干旱,在東部饒河縣的西部地區發生了中度干旱,范圍較小,全區總體上土壤濕度狀況良好。9月份全區受旱范圍廣泛,除三江平原東北部以及西北部地區外,全區其他地方均發生不同程度的干旱,干旱等級主要為輕度干旱,西部的勃利縣以及中部的寶清縣部分區域為中度干旱,研究區的土壤濕度狀況整體呈現偏干旱的情勢。
對比4個時期的土壤濕度分布圖,在三江平原東部及東北部,土壤濕度一直保持為正常或輕度濕潤。相關研究[33]表明三江平原東部及東北部地物類型主要為水田或沼澤等濕地及森林,可能是森林具有的對地表溫度的調控作用[34],間接影響了土壤濕度。濕地區域由于土壤水分長期處于過飽和狀態,遙感手段反演的土壤濕度狀況顯示為正常或輕度濕潤。
針對單一植被指數對不同植被覆蓋狀況下敏感性不同的問題,本文對比分析了不同植被指數下TVDI對土壤濕度狀況的監測效果,得到了合理的結論,TVDI能夠適用于三江平原的土壤濕度狀況監測,與前人研究[9,18,24,29]結論一致。在不同時期選用何時的植被指數計算TVDI能夠更好的反演土壤濕度狀況。但是在研究中仍然存在不足之處有待提高:(1) 本研究在不同植被覆蓋狀況下選取不同的植被指數計算TVDI,雖然能夠在一定程度上提高土壤濕度狀況反演的準確性,但是土壤濕度狀況受多種因素的影響,下一步可以考慮建立綜合多種影響因素的土壤濕度狀況監測模型,進一步提高土壤濕度狀況監測的準確性。(2) 考慮到遙感影像和土壤濕度數據的時間對應及云量問題,本研究選用了四期遙感影像進行研究,對于不同植被覆蓋狀況為定性研究,下一步可以考慮實測土壤濕度,選用多期遙感影像進行分析,對于植被覆蓋狀況進行定量化,從而提高TVDI模型在時間上的適用性。(3) 基于MODIS遙感數據能夠監測大尺度土壤濕度狀況,本研究中也取得了較好的結果,但是MODIS影響容易受云和雨雪天氣影響,且分辨率較低,后期可以考慮利用雷達影像在小區域做更為細致的研究。
(1) TVDI方法能夠用來反演三江平原地區土壤濕度狀況,基于不同植被指數的TVDI與土壤濕度基本上呈負相關關系,即TVDI值越大,土壤濕度越低。(2) 不同植被指數計算的TVDI在5月、6月、9月這3個時期,回歸分析R2數值相近(5月均在0.52左右、6月均在0.20左右、7月均在0.55左右),即不同植被指數的TVDI在這3個時期對土壤濕度的反演效果差異較小,4種植被指數均適合用來反演這3個時間段的土壤濕度;7月份的EVI和MSAVI計算的TVDI的R2(均在0.36左右)明顯高于植被指數NDVI和RVI的R2(均在0.15左右),表明EVI和MSAVI計算的TVDI在7月份對土壤濕度的反演效果好于其他兩種植被指數, EVI和MSAVI更適合用來反演該時期的土壤濕度。(3) 三江平原2013年5—9月大部分地區濕度狀況屬于正?;蜉p度干旱,干旱現象主要出現在中部和西南部;東部及東北部地區由于森林和濕地的覆蓋,在不同時期均保持在正常和輕度濕潤狀況。