吳 見(jiàn),劉民士,李偉濤
(滁州學(xué)院 地理信息與旅游學(xué)院,安徽 滁州239000)
土壤是農(nóng)業(yè)生產(chǎn)的基礎(chǔ),土壤的水分含量、質(zhì)地、有機(jī)質(zhì)含量[1-2]等性質(zhì)會(huì)影響作物生長(zhǎng),其中,土壤水分不僅是植物生長(zhǎng)和發(fā)育的必要條件,也是精準(zhǔn)農(nóng)業(yè)中重要的信息。因此,土壤含水量監(jiān)測(cè)一直是學(xué)者們研究的重點(diǎn)[3-4]。傳統(tǒng)的土壤含水量測(cè)定方法,如張力計(jì)法、烘干法等均以點(diǎn)測(cè)量為基礎(chǔ),準(zhǔn)確性較高,但周期長(zhǎng)、工作量大。遙感技術(shù)主要以可見(jiàn)光—近紅外、熱紅外和微波波段對(duì)土壤含水量進(jìn)行反演,熱紅外波段反映的是地表很薄的一層信息,但也可利用熱慣量和輻射平衡估測(cè)表層以下的土壤水分,高光譜遙感具有精細(xì)的光譜分辨率,能根據(jù)可見(jiàn)光和近紅外波段的光譜信息,分析影響土壤含水量變化的光譜特征[5-6],建立光譜特征與 土壤含水量的 關(guān)系模型[7-8],進(jìn)而定量反演土壤含水量。
目前,土壤含水量遙感監(jiān)測(cè)研究大多是針對(duì)裸土進(jìn)行的,對(duì)于衛(wèi)星遙感來(lái)說(shuō),實(shí)際情況往往是植被與土壤同時(shí)存在于像元中,土壤含水量的監(jiān)測(cè)不可避免地會(huì)受到植被光譜信息的影響。實(shí)際研究發(fā)現(xiàn)土壤含水量變化引起的光譜差異與植被含水量變化引起的光譜差異并不相同[9],因此剔除植被光譜信息對(duì)土壤水分的干擾,是監(jiān)測(cè)土壤含水量的有效途徑。本研究對(duì)高光譜影像進(jìn)行光譜分解,剔除植被光譜的干擾,同時(shí)對(duì)土壤貢獻(xiàn)的光譜進(jìn)一步處理,建立土壤含水量反演模型。
懷柔區(qū)地處北京東北部,總面積達(dá)2 128.7km2,是華北平原、內(nèi)蒙古高原以及燕山山脈的過(guò)渡地帶。該區(qū)北部環(huán)山,南部是草原,地形包括深山、淺山、平原3種類型,從南向北綿延128km,且山地面積占整個(gè)懷柔區(qū)的88.7%,宜林山場(chǎng)森林面積達(dá)41%。地勢(shì)呈北高南低,海拔在34~1 661m,變化較大。該區(qū)氣候冬天受西伯利亞冷空氣影響,寒冷干燥,夏天受海洋氣團(tuán)影響,溫?zé)釢駶?rùn),春秋季節(jié)短且干旱多風(fēng)。年平均氣溫在6~12℃,常年平均降水量在470~850mm。土壤主要包括草墊土、褐土、棕壤、風(fēng)砂土等類型。
本研究選取北京懷柔部分地區(qū)2001年5月19日的EO-1Hyperion高光譜數(shù)據(jù),該數(shù)據(jù)共有242波段,其中1~70波段是可見(jiàn)光近紅外波段(VNIR),71~242波段是短波紅外波段(SWIR),光譜分辨率為10nm,空間分辨率為30m。
首先刪除2個(gè)重復(fù),20個(gè)受水汽影響嚴(yán)重及44個(gè)未定標(biāo)波段,剩余176個(gè)波段;然后對(duì)剩余波段進(jìn)行處理,包括壞線修復(fù)、條紋去除以及smile效應(yīng)去除;最后,對(duì)處理后的圖像進(jìn)行檢驗(yàn),繼續(xù)刪除質(zhì)量差的波段7個(gè),剩余169個(gè)波段。被刪除的具體波段為1~7,58~78,121~129,166~180,185~186,224~242。利用FLAASH軟件對(duì)剩余的169個(gè)波段進(jìn)行大氣糾正,得到反射率圖像。大氣糾正后,采用1∶50 000地形圖對(duì)影像進(jìn)行幾何糾正,總誤差是0.35個(gè)像元。
目前,全受限的線性光譜混合分解模型是最為常用的方法,該模型中像元的光譜表達(dá)成各基本組分光譜及其面積比率的線性組合[10],即:
式中:Rb——b波段的光譜信息;N——基本組分?jǐn)?shù)目;fi——基本組分i占所在像元的面積比例;Ri,b——基本組分i在b波段的光譜信息;eb——b波段的誤差。應(yīng)用該模型的關(guān)鍵是盡可能減少模型中每個(gè)像元的誤差,即使均方根誤差RMSE最小化[10]:
式中:M——波段數(shù)。解該模型最常用的方法是最小二乘法,即求得式(2)RMSE最小的解。
將影像或光譜庫(kù)中選取的基本組分光譜(r1,r2,…,rn)根據(jù)已設(shè)定的比例(f1,f2,…,fn)組合,獲取若干測(cè)試光譜S1=f1r1+f2r2+…+fnrn,然后將待分解的像元光譜,即目標(biāo)光譜S2與測(cè)試光譜S1進(jìn)行匹配[10],把匹配效果最佳的測(cè)試光譜的基本組分比例作為影像分解結(jié)果。因此,可將混合像元分解視為帶約束的非線性最優(yōu)問(wèn)題[10]:
式中:G——目標(biāo)函數(shù),描述目標(biāo)與測(cè)試光譜差異的匹配模型。由最小二乘法原理可定義G為測(cè)試與目標(biāo)光譜間的歐式距離[10]:
因S1b-S2b=eb,根據(jù)式(2)可知GLS= ■M ×RMSE,說(shuō)明最小二乘法也屬于光譜匹配的非線性最優(yōu)化算法,但GLS采用歐式距離描述光譜匹配程度時(shí),易受局部奇異值等因素影響。基于此,本研究引用楊偉等[10]提出的基于相關(guān)系數(shù)匹配的分解算法,設(shè)目標(biāo)函數(shù)GSCM為:
式中:SCM——測(cè)試與目標(biāo)光譜的相關(guān)系數(shù)。
采用相關(guān)系數(shù)SCM描述目標(biāo)與測(cè)試光譜之間的匹配程度,以SCM最大值的測(cè)試光譜作為混合像元分解的結(jié)果。
高光譜混合像元分解后,得到土壤和植被比例,以及裸土和純植被基本組分的光譜信息,因此可進(jìn)一步剔除植被光譜的干擾,計(jì)算各像元中土壤貢獻(xiàn)的反射率信息。獲取土壤貢獻(xiàn)的反射率信息后,為突出土壤光譜的吸收特征,本研究采用光譜包絡(luò)線去除[11]和一階微分[12]兩種方法對(duì)光譜曲線做進(jìn)一步處理。
2.2.1 反射率光譜信息 由于不同像元對(duì)應(yīng)的土壤含水量不同,在不同含水量條件下裸土基本組分光譜是有差異的,因此不能采用土壤比例信息與裸土基本組分光譜直接計(jì)算不同比例土壤貢獻(xiàn)的反射率信息。由于研究區(qū)內(nèi)植被類型較為單一,因此設(shè)定研究區(qū)域純植被基本組分光譜不變,從原始影像混合光譜中將植被光譜剔除,進(jìn)而獲取土壤貢獻(xiàn)的反射率光譜信息,由像元二分模型[13]:
可知,剔除植被光譜后的土壤光譜信息為:
式中:Rsb,Rv——土壤和植被在波段b貢獻(xiàn)的光譜信息;fs,fv——土壤和植被在像元中的比例。
2.2.2 光譜一階微分 一階微分法是常用的光譜處理技術(shù),不僅能夠識(shí)別重疊光譜,增強(qiáng)曲線坡度上的光譜細(xì)微變化,使光譜吸收峰參數(shù)便于提取,而且可以有效消除系統(tǒng)誤差,降低大氣吸收、散射、輻射等噪聲的影響[12]。一階微分公式為[12]:
式中:λb——波段b;FDRλb——波段b和b+1之間的一階微分值;Rλb,Rλb+1——波段b和b+1的光譜反射率;Δλ——波段b和b+1的波長(zhǎng)之差。
2.2.3 包絡(luò)線去除法 包絡(luò)線去除法能夠增強(qiáng)光譜曲線的反射、吸收特征,且反射率歸一化在0~1.0之間,是一種有效突出目標(biāo)吸收特征的光譜分析法,同時(shí)能夠?qū)⒐庾V吸收特征歸一化在同一背景上,方便不同光譜曲線之間的比較[11]。“包絡(luò)線”即以直線逐點(diǎn)連接光譜曲線上的峰值點(diǎn),且峰值點(diǎn)上的折線外角大于180°。光譜去包絡(luò)是將原始曲線上的值除以包絡(luò)線上對(duì)應(yīng)的光譜值,計(jì)算公式為[11]:
式中:λb——波段b;Rcb——波段b的包絡(luò)線去除值;Rb——波段b的反射率值;Rstart,Rend——光譜曲線起始、末端點(diǎn)的反射率值;λstart,λend——光譜曲線起始、末端點(diǎn)的波長(zhǎng);K——光譜曲線起始、末端波段間的斜率。
2.3.1 土壤含水量與反射率光譜相關(guān)性分析 從圖1可知,反射率光譜與土壤含水量在450~1 050nm的相關(guān)系數(shù)為正值,相關(guān)系數(shù)的峰值分別出現(xiàn)在波段456~465,496~505,786~795,1 004~1 013nm處,相關(guān)系數(shù)值分別為0.47,0.48,0.37,0.29;在1 050~2 380nm的相關(guān)系數(shù)為負(fù)值,相關(guān)系數(shù)的峰值分別出現(xiàn)在波段1 134~1 143,1 324~1 333,1 467~1 476,2 012~2 021,2 302~2 311nm處,相關(guān)系數(shù)值分別為-0.42,-0.38,-0.51,-0.61,-0.56。為避免“過(guò)度擬合”現(xiàn)象及共線性問(wèn)題,并非全部選取相關(guān)系數(shù)值最大的波段,而是根據(jù)光譜曲線特征劃分區(qū)域,選取各區(qū)域中具有代表性或相關(guān)系數(shù)絕對(duì)值最大的波段作為敏感波段,即496~505,786~795,1 134~1 143,1 467~1 476,2 012~2 021nm處。
圖1 土壤含水量與反射率光譜相關(guān)系數(shù)
2.3.2 土壤含水量與一階微分光譜相關(guān)性分析 一階微分光譜與土壤含水量的相關(guān)系數(shù)在部分波段比反射率光譜明顯增加,部分波段相對(duì)減小,且一階微分光譜與土壤含水量之間相關(guān)性時(shí)正時(shí)負(fù),不呈單一相關(guān)(如圖2所示)。相關(guān)系數(shù)的正峰值分別出現(xiàn)在波段1 004~1 013,1 174~1 183,1 194~1 203,1 467~1 476,1 497~1 506,1 787~1 796,1 797~1 806nm處,相關(guān)系數(shù)值分別為0.37,0.52,0.39,0.34,0.23,0.76,0.69;相關(guān)系數(shù)的負(fù)峰值分別出現(xiàn)在波段1 064~1 073,1 264~1 273,1 324~1 333,1 527~1 536,1 952~1 961,1 972~1 981nm處,相關(guān)系數(shù)值分別為-0.25,-0.48,-0.40,-0.25,-0.58,-0.47。選取各區(qū)域中具有代表性或相關(guān)系數(shù)絕對(duì)值最大的波段作為敏感波段,即1 004~1 013,1 174~1 183,1 324~1 333,1 787~1 796,1 952~1 961nm處。
2.3.3 土壤含水量與包絡(luò)線去除光譜相關(guān)性分析從圖3可以看出,包絡(luò)線去除光譜與土壤含水量在1 560~1 805nm為負(fù)相關(guān),其余波長(zhǎng)范圍內(nèi)為正相關(guān)。相關(guān)系數(shù)的峰值分別出現(xiàn)在波段656~665,736~745,1 014~1 023,1 274~1 283,1 467~1 476,1 517~1 526,2 062~2 071,2 242~2 251nm處,相關(guān)系數(shù)值分別為0.66,0.67,0.66,0.62,0.26,0.26,0.54,0.49。選取各區(qū)域中具有代表性或相關(guān)系數(shù)絕對(duì)值最大的波段作為敏感波段,即656~665,1 014~1 023,1 517~1 526,2 062~2 071,2 242~2 251nm波段。
圖2 土壤含水量與一階微分光譜相關(guān)系數(shù)
圖3 土壤含水量與包絡(luò)線去除光譜相關(guān)系數(shù)
選取反射率、一階微分、包絡(luò)線去除光譜與土壤含水量敏感的波段進(jìn)行逐步回歸分析。通過(guò)比較,建立土壤含水量反演多項(xiàng)式(表1—3),選擇R2值最大的反演模型,作為土壤光譜與土壤含水量的關(guān)系模型。
表1 反射率光譜與土壤含水量逐步回歸分析模型
表2 一階微分光譜與土壤含水量逐步回歸分析模型
表3 包絡(luò)線去除光譜與土壤含水量逐步回歸分析模型
由表1—3可知,預(yù)測(cè)土壤含水量的最佳模型是以波段X661,X1019和X2067的土壤包絡(luò)線去除光譜為自變量建立的回歸方程:
其預(yù)測(cè)R2值為0.85。采用實(shí)測(cè)含水量數(shù)據(jù)與式(12)得到的預(yù)測(cè)值進(jìn)行擬合,檢驗(yàn)結(jié)果的R2值為0.80(如圖4所示)。
圖4 土壤含水量預(yù)測(cè)與真實(shí)值比較
若不對(duì)原始影像進(jìn)行分解,直接采用未剔除植被光譜的像元反射率、一階微分、包絡(luò)線去除光譜建立土壤含水量估測(cè)模型,得到的最佳模型是以波段X541,X979和X1632的一階微分光譜為自變量建立的回歸方程:
其預(yù)測(cè)R2值為0.36。采用實(shí)測(cè)含水量數(shù)據(jù)與式(13)得到的估測(cè)值進(jìn)行擬合,檢驗(yàn)結(jié)果的R2值僅為0.06。因此,剔除植被光譜干擾能夠更有效地反演土壤含水量。
(1)將混合像元分解問(wèn)題歸結(jié)為一個(gè)基于光譜匹配的非線性最優(yōu)化問(wèn)題,并針對(duì)最小二乘法的不足,引用了一種基于相關(guān)系數(shù)匹配的混合像元分解技術(shù)。高光譜數(shù)據(jù)本身提供了豐富的光譜信息,可以預(yù)見(jiàn)在高光譜數(shù)據(jù)中利用基于相關(guān)系數(shù)匹配的混合像元分解算法能夠得到更為精確的結(jié)果。
(2)剔除植被光譜干擾后,土壤反射率、一階微分、包絡(luò)線去除光譜的部分波段能夠敏感地反映土壤含水量的變化,最佳模型的估測(cè)R2值為0.85;若直接采用原始影像光譜建立估測(cè)模型,估測(cè)R2值僅為0.36。
(3)通過(guò)高光譜影像分解剔除植被光譜干擾估測(cè)土壤含水量的方法是可行的,可為今后遙感估測(cè)土壤含水量的研究提供參考。
[1] 周萍,王潤(rùn)生,閻柏琨,等.高光譜遙感土壤有機(jī)質(zhì)信息提取研究[J].地理科學(xué)進(jìn)展,2008,27(5):27-34.
[2] 陳榮毅,張?jiān)?,潘伯榮,等.古爾班通古特沙漠土壤養(yǎng)分空間分異與干擾的關(guān)系[J].中國(guó)沙漠,2007,27(2):257-265.
[3] 宋韜,鮑一丹,何勇.利用光譜數(shù)據(jù)快速檢測(cè)土壤含水量的方法研究[J].光譜與光譜學(xué)分析,2009,29(3):675-677.
[4] Bedidi A,Cervelle B,Madeira J,et al.Moisture effects on visible spectral characteristics of lateritic soils[J].Soil Science,1992(153):129-141.
[5] 孫建英,李民贊,唐寧,等.東北黑土的光譜特性及其與土壤參數(shù)的相關(guān)性分析[J].光譜學(xué)與光譜分析,2007,27(8):1502-1505.
[6] 劉煥軍,張柏,宋開(kāi)山,等.黑土土壤水分光譜響應(yīng)特征與模型[J].中國(guó)科學(xué)院研究生院學(xué)報(bào),2008,25(4):503-509.
[7] 王昌佐,王紀(jì)華,王錦地,等.裸土表層含水量高光譜遙感的最佳波段選擇[J].遙感信息,2003(4):33-36.
[8] 王靜,何挺,李玉環(huán).基于高光譜遙感技術(shù)的土地質(zhì)量信息挖掘研究[J].遙感學(xué)報(bào),2005,9(4):438-445.
[9] 吳代暉,范聞捷,崔要奎,等.高光譜遙感監(jiān)測(cè)土壤含水量研究進(jìn)展[J].光譜學(xué)與光譜分析,2010,30(11):3067-3071.
[10] 楊偉,陳晉,松下文經(jīng),等.基于相關(guān)系數(shù)匹配的混合像元分解算法[J].遙感學(xué)報(bào),2008,12(3):454-461.
[11] 丁麗霞,王志輝,葛宏立.基于包絡(luò)線法的不同樹種葉片高光譜特征分析[J].浙江林學(xué)院學(xué)報(bào),2010,27(6):809-814.
[12] 尹業(yè)彪,李霞,趙釗,等.沙質(zhì)土壤含水率高光譜預(yù)測(cè)模型建立及分析[J].遙感技術(shù)與應(yīng)用,2011,26(3):355-359.
[13] 牛寶茹,劉俊蓉,王政偉.干旱半干旱地區(qū)植被覆蓋度遙感信息提取研究[J].武漢大學(xué)學(xué)報(bào):信息科學(xué)版,2005,30(1):27-30.