基于Sentinel-1/2的土壤鹽分含量反演研究

馬 馳,劉曉波

(1.遼寧省交通高等?？茖W校，遼寧 沈陽 110122；2.中國電建集團北京勘測設計研究院有限公司，北京 100024)

土壤鹽堿化是土壤環(huán)境災害的重要表現(xiàn)，一般發(fā)生在地勢低平且水位較高、排水不暢的干旱與半干旱地區(qū)，已成為世界性有待解決的生態(tài)環(huán)境難題[1-2]。因此，實時監(jiān)測、評估土壤鹽堿化嚴重地區(qū)的土壤含鹽量對于區(qū)域土壤環(huán)境的恢復與保護具有重要意義。傳統(tǒng)的監(jiān)測區(qū)域土壤含鹽量的方法多為實地調(diào)查土壤采樣，并在實驗室化驗土樣的含鹽量。此方法雖然精度較高，但費時費力，只能獲得一定數(shù)量的點狀數(shù)據(jù)。遙感技術(shù)具有數(shù)據(jù)獲取方便迅捷、影像覆蓋范圍廣泛、數(shù)據(jù)成本低廉等優(yōu)點，已廣泛應用于土地鹽堿化的監(jiān)測與評估中[3-5]。

合成孔徑雷達(Synthetic aperture radar，SAR)具有全天時、全天候等特點，對植被、云霧、雨雪具有一定的穿透能力，能夠更好地探測地表信息[6-7]。然而，于合成孔徑雷達屬于主動微波遙感，SAR影像的后向散射系數(shù)容易受到地形、極化方式、入射角等因素的影響，降低了對地表地物識別的精度[8-10]。光學遙感影像能夠反映地表地物的光譜信息，在一定程度上可以彌補SAR數(shù)據(jù)在地物識別中的缺陷。研究表明，將光學遙感影像與SAR影像進行融合應用，可以增加融合影像中的有用信息，減少地物識別過程中的差異性與不確定性，從而增強對地物的識別能力。

本文試驗利用Gram-Schmidt(GS)融合算法將Sentinel-1A雙極化雷達數(shù)據(jù)與Sentinel-2A多光譜影像進行融合，以我國典型蘇打鹽堿土區(qū)——吉林省白城市為研究對象，結(jié)合野外調(diào)查的土壤采樣實驗室化驗數(shù)據(jù)，建立白城市土壤含鹽量的反演模型，實現(xiàn)對研究區(qū)土壤含鹽量的監(jiān)測與制圖，為區(qū)域土壤環(huán)境監(jiān)測提供數(shù)據(jù)支持，為Sentinel-1A雷達數(shù)據(jù)與Sentinel-2A多光譜數(shù)據(jù)融合及應用提供參考。

1 材料與方法

1.1 研究區(qū)概況

白城市地處吉林省西部，是世界三大蘇打鹽堿土區(qū)域之一，面積約2.6萬km2，氣候?qū)儆谥袦貛Т箨懶约撅L氣候。白城市西部、北部地區(qū)為山地地貌，海拔較高，中部、東部和南部地區(qū)為平原地貌，地勢低平。多年來，受自然條件、土壤質(zhì)地、人為因素等影響，白城市土壤鹽堿化嚴重，實時監(jiān)測土壤含鹽量對于精確改良土壤質(zhì)量意義重大。

1.2 技術(shù)路線

將Sentinel-1A SAR影像與Sentinel-2A多光譜影像融合，反演白城市土壤含鹽量的技術(shù)路線如圖1。

圖1 技術(shù)路線圖Fig.1 Technology roadmap

1.3 研究區(qū)土壤采樣及化驗

根據(jù)白城市氣候及地表植被情況，選擇2019年春季5月5日至11日對研究區(qū)進行土壤采樣，此時研究區(qū)地表無冰雪且綠色植被稀少。首先，參考研究區(qū)遙感影像設計采樣路線，采樣點選擇在視野開闊的裸土區(qū)域，且均勻地分布于研究區(qū)內(nèi)，采樣點及采樣路線如圖2(見257頁)所示。采樣過程采用多點采樣法，即：在采樣位置20×20 m范圍內(nèi)采集4～5點表層(0～15 cm)土壤樣品，混合后收取約1 kg土壤裝入采集袋，同時利用GPS接收機測定采樣位置坐標，共采得土壤樣品54個。在實驗室化驗土壤樣品的含鹽量，為研究區(qū)表層土壤含鹽量的反演提供數(shù)據(jù)。本試驗將土壤樣品隨機分組：建模樣品38個、檢驗樣品16個，其含鹽量的描述統(tǒng)計見表1。

圖2 采樣路線圖Fig.2 Sampling roadmap

表1 土壤樣品分組及含鹽量統(tǒng)計

1.4 遙感數(shù)據(jù)選擇與處理

為了真實反映土壤采樣時刻研究區(qū)地表信息，選擇與采樣時間相近的遙感影像數(shù)據(jù)。其中，選取Sentinel-1A SAR影像2景，成像時間為2019年5月8日；選取Sentinel-2A多光譜遙感影像共4景，成像時間為2019年5月7日。Sentinel-1 SAR影像預處理采用歐空局發(fā)布的SNAP軟件，主要完成SAR影像的輻射定標、噪聲去除、Lee濾波、地理編碼等，生成具有WGS 84/UTM坐標系的后向散射系數(shù)影像。利用歐空局發(fā)布的Sen2cor和SNAP兩種軟件讀取Sentinel-2A多光譜遙感影像并對其進行影像預處理工作，主要包括輻射定標和大氣校正，生成具有WGS 84/UTM坐標系的地表反射率影像。Sentinel-2A多光譜影像的幾何校正和裁剪等工作在ENVI5.3軟件中完成。

1.5 遙感影像融合

本文采用Gram-Schmidt變換將Sentinel-1A SAR影像與Sentinel-2A多光譜遙感影像進行融合處理。GS變換能夠保證影像融合前后波譜信息的一致性，是一種高保真的遙感影像融合方法[11-12]。GS波段融合表示如下：

(1)

式中，B表示原始波段，T表示已轉(zhuǎn)換波段的編號，μT表示T波段的均值，φ(BT,GSl)表示高分辨率全色波段的協(xié)方差。

(2)

式中，C為圖像總列數(shù)，R為圖像的總行數(shù)。

(3)

計算得到平均值μT、方差δT以及協(xié)方差φ后，按以下公式調(diào)整高分圖像：

(4)

最后，利用Gram-Schmidt逆變換生成融合圖像。融合后影像的反射率表示為：

(5)

融合工作主要包括：將預處理后的SAR波段作為模擬的全色波段，將預處理后的Sentinel-2A影像作為多光譜波段；對全色和多光譜波段進行GS變換，其中模擬的全色波段為第一波段，并將GS變換后的第一波段替代全色波段；最后，應用GS反變換生成融合波段，實現(xiàn)SAR與多光譜影像的融合。

2 結(jié)果與分析

2.1 影像的后向散射系數(shù)或反射率與土壤含鹽量的相關(guān)性分析

2.1.1 SAR影像與土壤含鹽量的相關(guān)性 參考前人的研究[10,13-14]，將Sentinel-1A SAR影像后向散射系數(shù)及其波段組合與研究區(qū)土壤含鹽量進行相關(guān)性分析，本文將Sentinel-1A SAR影像的VH、VV波段后向散射系數(shù)進行函數(shù)變換及波段組合。波段的函數(shù)變換包括VH、VV兩個波段的倒數(shù)、對數(shù)、冪、指數(shù)、二次方等，波段組合包括VH和VV兩個波段求商、兩個波段的平方求商、兩個波段的和與兩個波段的差求商、兩個波段的平方和與平方差求商等。相關(guān)性分析結(jié)果顯示(表2)，Sentinel-1A SAR影像的VV波段后向散射系數(shù)與研究區(qū)土壤含鹽量的相關(guān)性較好，相關(guān)系數(shù)為r=0.628；將兩個波段進行函數(shù)變換后，VV波段后向散射系數(shù)的二次方與土壤含鹽量相關(guān)性有較明顯的改善，相關(guān)系數(shù)達到0.641；將波段進行組合并對波段組合進行函數(shù)變換，并未提高與土壤含鹽量的相關(guān)性。

表2 影像的后向散射系數(shù)或反射率與土壤含鹽量的相關(guān)性分析

2.1.2 Sentinel-2A影像反射率與土壤含鹽量相關(guān)性分析 由表2可見，Sentinel-2A遙感影像反射率與研究區(qū)土壤含鹽量呈正相關(guān)性，其中第5波段相關(guān)性最好，達到了r=0.696；反射率經(jīng)過適當?shù)暮瘮?shù)變換可以改善與土壤含鹽量的相關(guān)性，其中，反射率的冪函數(shù)變換與土壤含鹽量的相關(guān)性在第5波段達到峰值，相關(guān)系數(shù)為r=0.768，也是與土壤含鹽量相關(guān)性最好的函數(shù)變換形式；第11、12兩個波段的反射率及其函數(shù)變換與研究區(qū)土壤含鹽量的相關(guān)性較低，均未達到顯著相關(guān)性標準；本文參考前人研究成果引入多種鹽分指數(shù)，并分析鹽分指數(shù)與研究區(qū)土壤含鹽量的相關(guān)性。結(jié)果顯示，鹽分指數(shù)(B2B4)/B3與土壤含鹽量相關(guān)性較高，達到r=0.822。

2.2 反演模型的建立

在54個土壤樣本中隨機選取38個作為建模樣本，在SPSS軟件中以建模樣本土壤含鹽量為因變量、遙感影像的反射率(反射率的變換)為自變量，利用回歸分析的方法建立研究區(qū)土壤含鹽量的反演模型(表3)。建模結(jié)果顯示，以Sentinel-1A SAR影像VV波段的后向散射系數(shù)建立的研究區(qū)土壤含鹽量反演模型Y=13.937X+4.255，模型決定系數(shù)為0.405，均方根誤差RMSE為2.626 g·kg-1；以VV波段建立的曲線模型Y=19.763X-6.224X2+3.285，模型精度有所改善，模型決定系數(shù)為0.418；以Sentinel-2A影像建立的研究區(qū)土壤含鹽量反演模型中，以鹽分指數(shù)(B2B4)/B3建立的冪函數(shù)模型精度最好，模型的決定系數(shù)R2達到0.675，均方根誤差為2.286 g·kg-1；SAR影像與多光譜影像經(jīng)過GS融合后建立的土壤含鹽量反演模型，模型精度被顯著改善。其中，將SAR影像的VV波段分別與Sentinel-2A影像的第2、第3、第4波段融合后，以鹽分指數(shù)(D2D4)/D3建立的曲線模型Y=86.260X-66.206X2-5.312精度最好，模型決定系數(shù)R2達到0.791，均方根誤差為1.884 g·kg-1。

表3 土壤含鹽量反演模型

2.3 模型精度檢驗

將除去建模樣本而剩余的16個土壤樣本作為模型檢驗樣本，參考建模結(jié)果，選擇以Sentinel-1A SAR影像VV波段與Sentinel-2A第2、3、4波段融合建立的鹽分指數(shù)(D2D4)/D3為因變量，建立的土壤含鹽量反演模型Y=86.260X-66.206X2-5.312，反演檢驗樣本土壤含鹽量，以式(6)計算檢驗樣本土壤含鹽量的相對誤差，相對誤差如圖3。

圖3 相對誤差散點圖Fig.3 Relative error scatter diagram

δ=(L測-L反)/L測

(6)

式中，δ為相對誤差，L測為檢驗樣本土壤含鹽量的實驗室測試值，L反為檢驗樣本土壤含鹽量的反演值。

檢驗樣本土壤含鹽量的相對誤差統(tǒng)計結(jié)果顯示，16個檢驗樣本含鹽量實驗室測試值與反演值最大相對誤差為-36.49%，最小相對誤差為4.76%，其中，共有12個檢驗樣本的相對誤差在-20%～20%區(qū)間，占總檢驗樣本的75%，表明土壤含鹽量反演模型Y=86.260X-66.206X2-5.312具有較高的反演精度。

2.4 研究區(qū)土壤含鹽量制圖

將Sentinel-1A的VV極化波段與Sentinel-2A第2、3、4波段融合并引入鹽分指數(shù)(D2D4)/D3建立的白城市土壤含鹽量反演模型Y=86.260X-66.206X2-5.312，反演研究區(qū)的土壤含鹽量，見圖4。制圖結(jié)果顯示，研究區(qū)西南部、南部及北部地區(qū)土壤鹽分含量較高，局部地區(qū)土壤含鹽量甚至高于5 g·kg-1，研究區(qū)中部、東部以及西北部地區(qū)土壤含鹽量普遍低于1 g·kg-1。

3 討 論

本文以Sentinel-1A、Sentinel-2A主動與被動遙感影像為數(shù)據(jù)源，以我國典型蘇打鹽堿土區(qū)——白城市為研究對象，結(jié)合研究區(qū)調(diào)查、土壤采樣化驗數(shù)據(jù)反演了白城市土壤鹽分含量。劉全明等[13-14]研究表明，土壤中鹽分含量的升高會影響土壤的介電常數(shù)以及表層土壤的后向散射系數(shù)，為利用雷達遙感監(jiān)測土壤含鹽量變化提供了可能。本文將Sentinel-1A遙感影像VH、VV波段后向散射系數(shù)與研究區(qū)土壤含鹽量進行相關(guān)性分析，分析結(jié)果顯示VH、VV波段的后向散射系數(shù)與土壤含鹽量呈正相關(guān)，且VV波段的相關(guān)性較好，相關(guān)系數(shù)達到0.628，其顯著性P<0.001，表明在超過99%的置信區(qū)間內(nèi)兩者顯著相關(guān)。將SAR影像與多波段影像融合后，融合波段與土壤含鹽量的相關(guān)性較SAR影像及多波段影像均有顯著提升。本文研究結(jié)果表明，Sentinel-2A各波段反射率與土壤含鹽量相關(guān)系數(shù)在第5波段達到峰值，為r=0.696，而SAR影像與Sentinel-2A的第2波段融合后與土壤含鹽量相關(guān)系數(shù)達到0.753。將SAR影像后向散射系數(shù)及多波段反射率進行適當?shù)臄?shù)學變換后與土壤含鹽量的相關(guān)性有不同程度的改善，其中，將VV波段后向散射系數(shù)進行二次方變換后與土壤含鹽量相關(guān)系數(shù)由0.628提升至0.641；將多波段影像第5波段進行冪函數(shù)變換后，與土壤含鹽量的相關(guān)系數(shù)由0.696提升至0.768；將融合波段進行二次方變換后與土壤含鹽量的相關(guān)系數(shù)由0.735提升至0.820。前人研究發(fā)現(xiàn)[15-17]，遙感影像的可見光與近紅外波段對含鹽土壤較敏感，將影像反射率與土壤含鹽量相關(guān)聯(lián)構(gòu)成適當?shù)柠}分指數(shù)，可以有效地監(jiān)測表層裸土含鹽量。本文引入(B4-B8)/(B4+B8)、(B2B4)/B3等多種鹽分指數(shù)，與土壤含鹽量進行相關(guān)性分析，結(jié)果顯示：利用Sentinel-2A第2、3、4波段建立的鹽分指數(shù)(B2B4)/B3與研究區(qū)土壤含鹽量相關(guān)系數(shù)達到r=0.822，SAR影像VV波段與Sentinel-2A第2、3、4波段融合建立的(D2D4)/D3鹽分指數(shù)與土壤含鹽量相關(guān)系數(shù)達到r=0.889；龐治國等[18-19]對覆蓋白城市的遙感影像進行人機交互式解譯提取土地鹽堿化信息，并對鹽堿化土地面積進行統(tǒng)計分析時發(fā)現(xiàn)，鹽堿化土壤主要分布于白城市南部的通榆縣、北部的鎮(zhèn)賚縣以及東部的大安市，在對研究區(qū)進行調(diào)研、土壤采樣時發(fā)現(xiàn)，白城市東部、南部及北部地區(qū)中小型湖泊眾多，受氣候及人為影響，湖泊周圍土地鹽堿化嚴重。

本文以融合波段鹽分指數(shù)(D2D4)/D3建立的土壤含鹽量反演模型Y=86.260X-66.206X2-5.312，反演白城市表層裸土含鹽量，獲得了較好的反演效果，究其原因：(1)本文選取的遙感影像與土壤采樣時間同步，能夠真實反映地表地物信息；(2)在遙感影像獲取時刻及實地調(diào)查時刻、土壤采樣時刻，研究區(qū)地表無植被及冰雪覆蓋；(3)采用Gram-Schmidt(GS)變換將SAR影像與多光譜影像進行融合，避免了傳統(tǒng)融合方法造成個別波段信息過量集中及光譜響應范圍不一致等問題，較好保持了影像的紋理和光譜信息，提升了融合影像光譜信息的豐富度；(4)引入合適的鹽分指數(shù)改善了遙感影像對土壤中鹽分的敏感度。

然而，本試驗未考慮研究區(qū)土壤中有機質(zhì)及含水量對遙感影像光譜信息的影響，在今后的實踐中將予以研究。

4 結(jié) 論

本文利用Sentinel-1A的SAR影像、Sentinel-2A多光譜影像及兩者的融合影像反演白城市裸土含鹽量，獲得了以下結(jié)論：

(1)Sentinel-1A遙感影像的VH、VV波段后向散射系數(shù)與白城市表層土壤含鹽量呈顯著正相關(guān)性，表明Sentinel-1A影像可作為研究區(qū)土壤鹽堿化監(jiān)測的數(shù)據(jù)源。

(2)Gram-Schmidt(GS)影像融合方法適用于Sentinel-1A雷達影像與Sentinel-2A多光譜影像融合，融合后的遙感影像與研究區(qū)土壤含鹽量的相關(guān)性較融合前得到明顯改善。

(3)引入合適的鹽分指數(shù)可以有效提升Sentinel-2A多光譜影像及Sentinel-1/2A 融合影像與研究區(qū)土壤含鹽量的相關(guān)性，利用融合影像的第2、3、4波段建立的鹽分指數(shù)(D2D4)/D3與研究土壤含鹽量相關(guān)系數(shù)達到0.889，表明利用Sentinel-1A SAR影像與Sentinel-2A多光譜影像進行融合是一種監(jiān)測土壤含鹽量的有效方法。

(4)反演結(jié)果顯示，白城市西部、南部的中小型湖泊周圍土壤鹽堿化較嚴重，與實地調(diào)研結(jié)果相吻合。