基于MTVDI與DDI二元回歸模型對毛烏素沙地腹部土壤表層水分的研究

王思楠，李瑞平，韓 剛，田 鑫，王耀強(qiáng)，胡勇平

(內(nèi)蒙古農(nóng)業(yè)大學(xué)水利與土木建筑工程學(xué)院，內(nèi)蒙古 呼和浩特 010018)

土壤水分是監(jiān)測土地退化的重要指標(biāo)，同時也是衡量土壤干旱程度的重要指標(biāo)[1]。及時知曉土壤水分狀況，可以了解旱情程度[2]。遙感技術(shù)能夠有效、大面積、實(shí)時動態(tài)地獲取干旱地區(qū)旱情資料[3-5]，為各級政府與農(nóng)業(yè)生產(chǎn)部門提供決策依據(jù)。因此，利用遙感進(jìn)行旱情監(jiān)測是一個研究和應(yīng)用的熱點(diǎn)[6-8]。作為同時與歸一化植被指數(shù)(Normalized vegetation index，NDVI)和陸地表面溫度(Land surface temperature，LST)相關(guān)的溫度植被干旱指數(shù)(Temperature vegetation drought index，TVDI)可用于干旱監(jiān)測[9-12]，尤其是監(jiān)測特定年內(nèi)某一時期整個區(qū)域的相對干旱程度。

目前利用地表溫度-植被指數(shù)(Ts-NDVI)特征空間進(jìn)行旱情監(jiān)測已取得一定的進(jìn)展[13]，但是使用歸一化植被指數(shù)容易造成紅光飽和，背景的土壤噪聲也在一定程度上損害了NDVI的空間一致性。鑒于此，利用改進(jìn)的TVDI模型降低實(shí)際計算過程中的偏差，同時提高土壤水分反演計算的準(zhǔn)確性就顯得尤為重要。伍漫春等[14]在此基礎(chǔ)上采用土壤調(diào)節(jié)植被指數(shù)(Modified soil adjusted vegetation index，MSAVI)對其進(jìn)行改進(jìn)，證明地表溫度-土壤調(diào)節(jié)植被指數(shù)(Ts-MSAVI)能夠更好地反映區(qū)域土壤水分狀況，是一種更有效的土壤水分監(jiān)測方法[15]。

雖然TVDI 以及改進(jìn)的MTVDI都可以較好地監(jiān)測土壤水分，但是精度還不是很高。針對這一問題，本研究根據(jù)研究區(qū)的特點(diǎn)采用修正土壤調(diào)整植被指數(shù)MSAVI構(gòu)成的MTVDI指數(shù)，并利用野外同步實(shí)測的土壤水分?jǐn)?shù)據(jù)分析和比較不同土壤深度的監(jiān)測性能，探討荒漠化差值指數(shù)(Desertification difference index，DDI)對毛烏素沙地腹部土壤水分狀況的監(jiān)測精度的影響，能夠?qū)γ珵跛厣车剞r(nóng)牧業(yè)活動提供一定的指導(dǎo)，同時也為改善大尺度區(qū)域的氣候及改善毛烏素沙地生態(tài)環(huán)境打下基礎(chǔ)。

1 研究方法

1.1 研究區(qū)概況

研究區(qū)為毛烏素沙地腹地的烏審旗,它位于鄂爾多斯高原西南部，處于蒙、陜、寧經(jīng)濟(jì)發(fā)展的“金三角”地帶。東經(jīng)108°17′-109°40′,北緯37°38 ′-39°23′，面積11 645 km2。其所處位置為年平均氣溫6.8℃，多年平均降水量350～400 mm，多年平均蒸發(fā)量2 443 mm，海拔一般在1 300～1 400 m的干旱區(qū)，生態(tài)環(huán)境脆弱，沙生植被是本地區(qū)植被主體。

圖1 采樣點(diǎn)分布Fig.1 Sampling site distribution

圖2 子樣點(diǎn)采樣設(shè)計Fig.2 Sampling design at each site

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

1.2.1 遙感數(shù)據(jù) 研究用的2016年4月21日和9月28日兩期軌道號為128/33、128/34,覆蓋烏審旗的Landsat8 OLI遙感數(shù)據(jù)由地理空間數(shù)據(jù)云提供。對Landsat 影像數(shù)據(jù)進(jìn)行輻射定標(biāo),利用Modtran 4模型進(jìn)行大氣校正以及裁剪鑲嵌等一系列的處理。

1.2.2 土壤含水量數(shù)據(jù) 參照土地利用類型圖、荒漠化程度分布圖等資料，在2016年4月20日、21日與9月27日、28日借助手持GPS布設(shè)了23個樣區(qū)(圖1)，每個樣區(qū)中包含5個子樣點(diǎn)A1～A5共5組(1 km×1 km像元，如圖2)共計115個子樣點(diǎn)，重復(fù)0～10 cm、10～20 cm、20～30 cm的土層深度，逐點(diǎn)用土鉆采集樣品，并刮去樣本點(diǎn)土壤表層的浮土后迅速封裝。在實(shí)驗(yàn)室利用烘干法計算土壤含水量數(shù)據(jù)。

1.3 數(shù)據(jù)處理

1.3.1 構(gòu)建Albedo-NDVI特征空間 不同荒漠化信息相對應(yīng)的地表反照率Albedo與歸一化植被指數(shù)NDVI的特征空間中具有顯著的線性關(guān)系，計算公式如下：

(1)

A=0.356α2+0.130α4+0.373α5+0.085α6

+0.072α7-0.0018

(2)

式中,A為地表反照率，α2、α4、α5、α6、α7為經(jīng)過大氣校正后的第2、4、5、6、7波段的反射率。

為了進(jìn)一步確定公式(4)中的K值，在研究區(qū)選擇分布于不同沙漠化類型的300個點(diǎn)，并進(jìn)行歸一化處理，然后利用Albedo和NDVI兩組數(shù)據(jù)進(jìn)行回歸擬合(如圖3)，得到相應(yīng)的方程：

A=0.8438-0.5531×NDVI

(3)

1.3.2 荒漠化差值指數(shù)提取 Verstrate等[16]研究發(fā)現(xiàn)在代表荒漠化變化趨勢的垂直方向上劃分Albedo-NDVI特征空間，荒漠化差值指數(shù)DDI能夠有效地區(qū)分出不同程度的荒漠化土地，計算公式如下：

DDI=K×NDVI-Albedo

a×K=-1

(4)

式中,a為式(3)中垂線的斜率，由此，確定荒漠化差值指數(shù)DDI的最終表達(dá)式為：

DDI=1.808×NDVI-Albedo

(5)

1.3.3 植被指數(shù)-地表溫度特征空間 Sandholt等[17]在研究NDVI和Ts的散點(diǎn)圖時發(fā)現(xiàn)二者呈現(xiàn)三角形分布，在此基礎(chǔ)上，前人還發(fā)現(xiàn)NDVI和Ts構(gòu)成的特征空間可以間接表示土壤含水狀況，因此提出了溫度植被干旱指數(shù)的概念。本文使用TVDI由改進(jìn)型土壤調(diào)整植被指數(shù)和地表溫度[18]計算得到，其定義為公式：

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

(6)

圖3 Albedo與NDVI的回歸分析Fig.3 Albedo and NDVI regression analysis

Tsmax=MSAVIa+b
Tsmin=MSAVIc+d

(7)

式中,Ts為任意像元地表溫度；Tsmax為干邊上的地表溫度；Tsmin為濕邊上的地表溫度。a、b、c、d是干、濕邊通過線性擬合的模型參數(shù)。

2 結(jié)果分析

2.1 基于Landsat8數(shù)據(jù)的Ts和MSAVI構(gòu)建特征空間

通過Landsat8數(shù)據(jù)計算的植被指數(shù)對應(yīng)地表溫度的最大值和最小值，構(gòu)建毛烏素沙地腹部2期不同時相的Ts-MSAVI特征空間。根據(jù)公式(6)并擬合干、濕邊方程(表1)。結(jié)果表明，MSAVI對應(yīng)Ts的最大值和最小值呈近似線性關(guān)系，隨著植被指數(shù)的增大，地表溫度的最大值呈減小趨勢，地表溫度的最小值呈增大趨勢，二者差值呈減小趨勢。從圖4可以看出2期不同時相的散點(diǎn)圖的形狀相似[19]，其中9月的干、濕邊的R2都比4月的大。

2.2 Landsat8-MTVDI與土壤水分相關(guān)性分析

利用野外實(shí)測采樣點(diǎn)0～10 cm、10～20 cm、20～30 cm土層的土壤含水量值與MTVDI值進(jìn)行回歸分析用于反演精度驗(yàn)證。從圖5和圖6中可以看出不同遙感數(shù)據(jù)獲取的不同月份MTVDI和土壤各層含水量具有一定的負(fù)相關(guān)。即MTVDI越高，土壤含水量越低，滿足MTVDI指數(shù)值越大土壤水分越低的原理。4月份MTVDI指數(shù)與0～10 cm、10～20 cm、20～30 cm土壤含水量的R2值分別為0.656、0.646、0.637，整體高于9月份R2值0.457、0.436、0.431。不管是低植被覆蓋度的4月還是高植被覆蓋度的9月，在0～10 cm深度的R2值均大于10～20 cm、20～30 cm深度的R2值。說明MTVDI能夠較好地反映土壤表層的含水量狀況。

圖4 2016年4月(a)、9月(b)Landast8的Ts-MSAVI特征空間Fig.4 Ts-MSAVI feature space of Landast8 in April(a)and September(b)2016

表1 不同時相干濕邊擬合結(jié)果

圖5 2016年4月不同深度下土壤水分與MTVDI回歸分析Fig.5 The regression analysis of soil moisture and MTVDI at different depths in April 2016

圖6 2016年9月不同深度下土壤水分與MTVDI回歸分析Fig.6 The regression analysis of soil moisture and MTVDI at different depths in September 2016

2.3 基于MTVDI與DDI的二元線性回歸干旱監(jiān)測模型的分析

毛烏素沙地荒漠化程度的加劇，在土壤水分補(bǔ)給和散失的過程當(dāng)中又受到土地利用、地貌、地形等多種因素的影響。本研究考慮荒漠化因素的條件下，利用Landsat8遙感數(shù)據(jù)提取2016年4月、9月不同時相的MTVDI指數(shù)和DDI指數(shù)并與0～10 cm土層深度含水率建立二元線性監(jiān)測模型(SMC=a×MTVDI+b×DDI+c)，利用野外實(shí)測土壤0～10 cm土層含水量進(jìn)行驗(yàn)證。最后和單一使用MTVDI指數(shù)與0～10 cm土層深度含水量建立的一元線性監(jiān)測模型(SMC=a×MTVDI+c)進(jìn)行精度對比分析，結(jié)果如表2、3所示。

分析表2、3可知，2016年2個不同時相在0～10 cm土層深度的含水量的反演過程中，單獨(dú)利用MTVDI監(jiān)測含水量的相關(guān)系數(shù)分別為0.810、0.676，相對誤差分別為12.93%、14.35%，平均相對誤差為13.64%。當(dāng)引進(jìn)了DDI指數(shù)影響因子之后，二元線性監(jiān)測模型中相關(guān)系數(shù)分別為0.842、0.734，相對誤差分別為10.26%、11.64%，平均相對誤差為10.95%。二元線性監(jiān)測模型提高了反演土壤含水量的精度。說明在2個不同時相土壤表層0～10 cm深度的含水量受荒漠化信息的影響較大。二元線性回歸模型監(jiān)測0～10 cm土層深度的相對含水量空間分布，如圖7所示。

由圖7可知：上述2期土壤含水量反演圖均可反映烏審旗當(dāng)時的土壤含水量分布情況。分析表4可知，2016年2期數(shù)據(jù)中可得到該區(qū)域0～10 cm表層土壤含水量在4月份5%～10%之間的區(qū)域，占總面積的53.72 %以上，達(dá)到了6 256 km2，含水量偏低。0～10 cm表層土壤含水量在9月份10%～15%之間的區(qū)域，均占總面積的51.55 %以上，達(dá)到了6 003 km2，含水量偏高。與此同時根據(jù)研究區(qū)氣象局2016年氣象資料可知該年9月份降水偏高，與本文研究的結(jié)果基本一致。

表2 MTVDI模型對0～10 cm土層深度含水量的回歸分析

注：a,c為一元線性監(jiān)測模型的系數(shù)。

Note: “a” and “c” are the coefficients of a linear monitored model.

表3 MTVDI與DDI二元監(jiān)測模型對0～10 cm深度含水量的回歸分析

注：a,b,c為二元線性監(jiān)測模型的系數(shù)。

Note: “a”,“b” and “c” are the coefficients of binary linear monitored model.

圖7 2016年MTVDI與DDI二元監(jiān)測模型反演0～10 cm深度含水量空間分布Fig.7 Inversion of spatial distribution of water content in 0～10 cm soil depth in 2016 modeled by MTVDI and DDI

表4 土壤表層含水量面積分布

3 結(jié)論與討論

本研究建立了毛烏素沙地腹部MTVDI與DDI二元線性監(jiān)測模型與土壤表層0～10 cm含水率的關(guān)系模型，對該地區(qū)土壤表層0～10 cm含水率進(jìn)行反演與精度分析并統(tǒng)計表層土壤含水率的面積分布。

1)不同植被覆蓋度下的MTVDI均能反映不同深度的土壤含水量，且呈現(xiàn)不同程度的負(fù)相關(guān)，R2值4月份整體高于9月份，在0～10 cm深度的R2值均大于10～20 cm、20～30 cm深度的R2值，并且都高于0.4。說明MTVDI可以作為有效指示地表土壤水分狀況的指標(biāo)。

2)荒漠化程度加劇，地表覆蓋程度下降，地表能量與水分平衡發(fā)生變化，均可導(dǎo)致土壤水分發(fā)生改變。荒漠化指數(shù)DDI與MTVDI結(jié)合建立二元線性回歸模型監(jiān)測區(qū)域0～10 cm土層深度含水率，平均相對誤差值要比單獨(dú)TVDI小2.69 %。

3)通過二元線性回歸監(jiān)測模型反演的土壤表層含水量分布圖發(fā)現(xiàn)：該區(qū)域0～10 cm表層土壤含水量在4月份5%～10%之間的區(qū)域，占總面積的53.72%以上，達(dá)到了6256km2，含水量偏低；0～10 cm表層土壤含水量在9月份10%～15%之間的區(qū)域，達(dá)到了6 003 km2，均占總面積的51.55%以上，含水量偏高。

本研究在MTVDI指數(shù)的基礎(chǔ)上，考慮了毛烏素沙地腹部荒漠化信息，建立MTVDI指數(shù)和DDI指數(shù)的二元線性回歸監(jiān)測模型，在一定程度上提高了MTVDI指數(shù)指示土壤水分的合理性。從另一個角度看，MTVDI指數(shù)指示土壤水分的能力與DDI指數(shù)有很大關(guān)系，使用單一的MTVDI指數(shù)無法精確達(dá)到反映干旱特征的目的，這也促使干旱監(jiān)測由單因素向多因素綜合發(fā)展，綜合多指數(shù)的干旱監(jiān)測模型是研究復(fù)雜的干旱監(jiān)測問題的新途徑，在解決干旱監(jiān)測的復(fù)雜問題中有著較大的應(yīng)用潛力。本研究尚屬可行性研究，但單期的影像不能說明研究結(jié)果的普遍性，在后續(xù)的研究中將會進(jìn)一步考慮采用多期遙感數(shù)據(jù)做動態(tài)分析。