基于無人機-衛星遙感升尺度的土壤鹽漬化監測方法

陳俊英 王新濤 張智韜 韓 佳 姚志華 魏廣飛

(1.西北農林科技大學旱區農業水土工程教育部重點實驗室， 陜西楊凌 712100；2.西北農林科技大學水利與建筑工程學院， 陜西楊凌 712100)

0 引言

土壤鹽漬化是一個全球生態環境問題[1]，不僅會對作物生長產生負面影響，還會引起土地退化[2]。因此，快速準確獲取土壤鹽漬化信息對灌區的可持續發展具有重大意義。

衛星遙感具有快速和范圍廣等優點，是大范圍監測土壤鹽漬化的有效途徑，目前已得到了廣泛的應用[3-4]。如LOBELL等[5]通過MODIS影像計算EVI和NDVI，進行區域尺度的鹽分監測；姚遠等[6]利用Landsat TM7影像計算植被指數，進行北疆地區土壤鹽分反演。但衛星遙感數據存在高時間分辨率和高空間分辨率矛盾的問題[7]，為此，國內外學者開展了大量尺度轉換的研究[8-9]。LIANG[10]通過尺度上推分析了從30 m到1 km反照率和LAI的變化規律，發現其隨空間尺度的變化呈線性變化；張萬昌等[11]基于統計分析的方法實現了30 m ETM+影像反演LAI向900 m的轉換，并對LAI產品進行了校正，效果較好。上述方法多是基于衛星空間尺度之間的轉換，在大范圍監測地表動態研究中，常用的衛星遙感影像空間分辨率相對較低[12],對于灌區鹽漬化的精準監測存在一定的局限性。無人機作為新型遙感平臺，具有高時效、高空間分辨率、云下低空飛行、高機動性等優點[13]，可以迅速、準確地完成一定區域內的鹽漬化監測任務。無人機在大范圍監測土壤鹽漬化方面存在缺陷，將無人機遙感的高空間分辨率與衛星遙感的大范圍監測相結合，則可實現高精度大尺度的土壤鹽漬化監測。因此，開展無人機到衛星尺度轉換研究對土壤鹽漬化監測具有非常重要的理論意義，而目前針對這一方面的研究還鮮有報道。

本文以河套灌區沙壕渠灌域裸土期不同鹽漬化程度的表層土壤為研究對象，采用MLR模型、SR模型和RR模型，分別構建基于無人機數據和GF-1衛星遙感數據的土壤鹽漬化監測模型，并通過改進的TsHARP尺度轉換方法，將無人機多光譜數據(空間分辨率為6.5 cm)建立的趨勢面應用到GF-1衛星(空間分辨率為16 m)中，經過轉換殘差校正后，對升尺度結果進行定性和定量評價，以期為提高衛星遙感監測土壤鹽漬化精度提供理論支持。

1 材料與方法

1.1 研究區概況

河套灌區解放閘灌域沙壕渠灌域為本研究的試驗區，如圖1所示，其面積約為5 333.36 hm2，位于西北高原地帶，平均海拔1 000 m以上，屬于典型的干旱半干旱大陸性氣候，土壤類型為粉質粘壤土，平均干容重1.466 t/m3。地下水埋深1.5 m左右，年降雨量為66.3～200 mm，年蒸發量約2 000 mm，蒸降比在10以上。降水量少蒸發量大[14]，再加上氣候條件、土質因素及不合理灌溉方式導致了此地區相對嚴重的土壤鹽漬化問題。

圖1 研究區示意圖Fig.1 Map of research area

1.2 數據來源與處理

1.2.1無人機多光譜遙感圖像數據獲取

本次試驗所使用的無人機為大疆公司生產的M600型六旋翼無人機，其攜帶的多光譜相機傳感器為Mirco-MCA (簡稱u-MCA)，它具有遠程觸發的特點，包括6個波段的光譜采集通道，分別是藍波段(中心波長490 nm)、綠波段(中心波長550 nm)、紅波段(中心波長680 nm)、紅邊波段(中心波長720 nm)、近紅外波段1(中心波長800 nm)、近紅外波段2(中心波長900 nm)。試驗選取不同鹽漬化程度的5個區域，分別對其進行了無人機飛行試驗，如圖1所示。無人機飛行模式按照提前規劃的航線飛行，拍照模式為等時間間隔，主航線間圖像重疊率設置為80%，多光譜影像獲取時間為2018年5月1—5日，選擇天氣晴朗的日期進行飛行，無人機飛行高度為120 m，多光譜相機鏡頭垂直向下，此時影像所對應的地面分辨率為6.5 cm。在研究區內布設了參考板，以便對遙感影像進行輻射標定。

使用與多光譜相機配套的處理軟件(PixelWrench2)對獲取的影像進行配準與合成，導出為6波段多光譜影像的tif格式。將采樣點的經緯度導入ENVI5.3.1軟件，提取tif圖像中對應采樣點的灰度，利用參考板進一步計算出對應的反射率。

1.2.2GF-1多光譜遙感影像的獲取及預處理

以國產高分一號衛星影像(GF-1 WFV相機)為數據源。衛星影像的成像時間與實測土壤含鹽量日期同步，高分一號衛星數據的重訪周期為4 d，空間分辨率為16 m，包括4個波段，分別為藍波段(450～520 nm)、綠波段(520～590 nm)、紅波段(630～690 nm)和近紅外波段(770～890 nm)。

在ENVI5.3.1軟件中對下載的影像進行幾何精校正、輻射定標、大氣校正等預處理，經剪裁后得到研究區域衛星影像圖，將野外實測采樣點的經緯度導入到GF-1衛星遙感影像中，導出各采樣點對應像元的反射率，用于構建土壤鹽分含量的反演模型。

1.2.3鹽分數據的獲取

本文實地采樣時間與無人機影像獲取時間相同，在5個采樣區域均勻設置80個表層土壤采樣點(圖2只列舉了其中1個區域的采樣點)。采樣點大多數為無植被覆蓋的表層土，采樣方法為五點法，將采集的土樣干燥研磨處理后，配置土水比為1∶5的土壤溶液，經攪拌、靜置、沉淀、過濾后，采用電導率儀(DDS-307型，上海佑科儀器公司)測定土壤溶液電導率，對每個采樣點的5個土樣電導率取平均值作為該樣點處的電導率EC1∶5，并通過經驗公式計算土壤含鹽量：S=(0.288 2EC1∶5+0.018 3)×100%[15]。剔除采樣點中的8個異常值，剩余72個樣本進行本次試驗土壤含鹽量的反演，將樣本從大到小進行排序，建模集和驗證集的劃分按2∶1的比例進行等間隔取樣，可保證建模樣本和驗證樣本范圍一致且分布均勻。統計情況見表1。

圖2 采樣點分布圖(其中1個區域)Fig.2 Soil sampling locations(one of the areas)

表1 土壤鹽分數據統計結果Tab.1 Summary statistics of soil salinity

1.3 改進光譜指數的計算

研究表明，僅用個別波段或多個單波段數據提取地物信息有明顯的局限性，進行光譜指數分析運算可有效提高地表觀測數據對反演模型參數的敏感性[16]。因此許多學者建立了基于光譜指數特征空間的土壤鹽漬化監測模型[17-18]，根據王飛等[19]的研究，發現通過NDVI與SI構建的SDI(Salinization detection index)指數可以較好地反映鹽漬化情況。基于此，本文在傳統的鹽分指數的計算基礎上，構建了NDVI-S1光譜指數，用于本次土壤鹽漬化反演。計算公式為

(1)

其中

INDVI=(ρNIR-ρRED)/(ρNIR+ρRED)

S1=ρBLUE/ρRED

式中INDVI——歸一化植被指數

S1——鹽分指數

INDVI-S1——NDVI-S1光譜指數

ρBLUE、ρRED、ρNIR——無人機和GF-1衛星對應的藍波段、紅波段和近紅外波段反射率

1.4 空間升尺度原理

TsHARP方法常用于遙感地表溫度降尺度轉換中[7]，此方法假設地表溫度與NDVI的關系在各個尺度上不變，通過引入NDVI構造趨勢面，實現對地表溫度的尺度轉換。本文通過改進這種方法來進行土壤鹽漬化的升尺度研究：首先建立無人機尺度上土壤鹽分與趨勢面因子間的關系式

S0.065=F0.065(B0.065)

(2)

式中S0.065——無人機尺度上采用趨勢面因子反演的土壤含鹽量

B0.065——無人機尺度上的趨勢面因子，即光譜指數

F0.065——趨勢面反演函數，同樣適用于升尺度到GF-1衛星16 m空間分辨率上土壤鹽分與趨勢面因子間的反演

考慮到趨勢面可能受到土壤含水率等因素的影響，趨勢面因子很難完全反映土壤鹽分的分布情況，表現在無人機高分辨率0.065 m尺度上的轉換殘差ΔS，在本文中，也認為等于衛星尺度上的轉換殘差ΔS1，計算式為

ΔS=S-S0.065=ΔS1

(3)

尺度轉換后的土壤含鹽量，應該由無人機尺度上建立的趨勢面應用到GF-1數據上計算的土壤含鹽量和衛星尺度上的轉換殘差構成，計算式為

S16=F16(B16)+ΔS1=F0.065(B16)+ΔS1

(4)

式中S16——尺度轉換后的土壤含鹽量

B16——16 m衛星數據上的變量因子

1.5 升尺度精度評價指標

本文采用決定系數R2(Coefficient of determination)、均方根誤差(Root mean square error，RMSE)和平均絕對誤差(Mean absolute error，MAE)[20]這3個指標來綜合評價升尺度后模型精度。R2能夠說明轉換前后數據擬合的精度高低，RMSE和MAE能夠反映升尺度結果與實測數據之間的偏差程度。R2越接近1，同時RMSE和MAE越小，說明模型效果越好。

2 結果分析

2.1 光譜波段和光譜指數與實測鹽分數據相關性分析

將無人機和GF-1衛星對應的光譜波段分別與實測土壤鹽分數據進行相關性分析，如表2所示。

由表2可以看出，無人機遙感影像中，除B3和B4波段外，其他4個波段與土壤鹽分均表現出較高的相關性，其中B1和B6波段在0.01水平上顯著相關，相關系數均大于0.43，B5波段與土壤含鹽量的相關性也相對較好，其相關系數為0.367， B2波段的相關系數最低，僅為0.233；對于GF-1衛星遙感影像， B1和B5波段與土壤鹽分的相關性在0.01水平上顯著相關，且相關系數均達到了0.35以上。綜上所述，無人機和衛星遙感影像的B1、B5波段與土壤含鹽量的相關性最好。

表2 光譜波段與表層土壤鹽分相關系數Tab.2 Correlation coefficient between each band data and surface soil salinity

注：*表示相關性顯著(P<0.05)，** 表示相關性極顯著(P<0.01)，下同。

將實測土壤鹽分數據與常用的光譜指數相關性進行分析，其結果如表3所示。

表3 光譜指數與表層土壤鹽分相關系數Tab.3 Correlation coefficient between salt index and surface soil salinity

由表3可以看出，通過無人機遙感影像計算的鹽分指數SI和S5與土壤鹽分的相關性較好，相關系數分別為0.356和0.441。通過GF-1衛星遙感影像計算的鹽分指數SI和S5與土壤鹽分的相關性也相對較高，相關系數分別為0.306和0.315。對于改進光譜指數(NDVI-S1)，在無人機遙感數據和GF-1衛星遙感數據上均呈現良好的相關性，相關系數分別為0.476和0.547。

綜上所述，本文選取在0.01水平上顯著相關的2個敏感波段(B1、B5)和3個敏感光譜指數(SI、S5和NDVI-S1)作為模型的輸入變量，用于建立遙感影像與土壤含鹽量的定量關系。

2.2 多重共線性分析

多元回歸分析方法在數理統計中有著廣泛的應用，但變量因子間的近線性關系會影響回歸方程的穩定性[26]，所以有必要對變量因子間進行共線性診斷。本文選用相關系數矩陣和方差膨脹因子(VIF)兩種方法來綜合評估5個變量因子之間線性關系的強弱。

2.2.1相關系數矩陣

通過統計分析軟件SPSS 23.0對無人機數據和GF-1數據的變量因子分別進行相關系數矩陣分析，結果見表4和表5。

表4 無人機數據相關系數矩陣Tab.4 Correlation coefficient matrix of UAV data

表5 GF-1數據相關系數矩陣Tab.5 Correlation coefficient matrix of GF-1 data

當兩個變量因子間的相關系數越接近于1時，認為其線性關系越強。其中，相關系數超過0.9時，認為存在共線性問題，在0.8以上時可能會有問題[27]。由表4可以直觀地看出， B5變量因子與其他幾個變量因子的相關系數小于0.8，共線性關系較弱；其余變量因子的相關系數均大于0.8，表示可能存在共線性問題。從表5可以看出，5個變量因子之間，存在相關系數大于0.8的情況。綜上所述，無人機數據變量因子和GF-1數據變量因子間都存在共線性問題。

2.2.2方差膨脹因子

為進一步說明變量因子間的共線性問題，采用方差膨脹因子VIF來表達其嚴重性。利用SPSS 23.0軟件對無人機遙感數據和GF-1遙感數據的變量因子分別進行方差膨脹因子VIF分析，結果見表6。

表6 方差膨脹因子VIF統計Tab.6 Spectral variance expansion factor VIF statistics

一般認為某個變量因子的VIF大于10時，則表明該變量與其他自變量間具有較強的共線性問題[28]，由表6可以看出，對于無人機數據，只有B5波段的VIF小于10，其他幾個變量因子的VIF均大于10。對于GF-1數據，除B5波段和改進光譜指數，其他3個變量因子的VIF均大于10，其中改進光譜指數的VIF最低，為4.079，比一般鹽分指數的VIF有明顯的降低。充分說明了無人機數據和GF-1數據各變量因子間存在嚴重的共線性問題，這種問題會導致模型的穩定性相對較差，預測的精度也會隨之降低。本研究通過使用MLR模型，以及可以有效解決共線性問題的SR模型和RR模型進行表層土壤的鹽分反演。

2.3 不同數據源含鹽量回歸模型分析

采用統計分析軟件SPSS，以樣本的2/3進行建模，1/3進行驗證。對GF-1衛星數據的變量因子進行建模和驗證分析，結果如表7所示，Y為土壤含鹽量；X1為GF-1衛星波段1反射率；X2為GF-1衛星波段5反射率；X3為GF-1衛星計算的SI光譜指數；X4為GF-1衛星計算的S5光譜指數；X5為GF-1衛星計算的NDVI-S1光譜指數；MLRS、SRS、RRS下標S表示衛星數據反演。

表7 GF-1衛星變量因子回歸模型Tab.7 Monitoring models of GF-1 satellite variable factor

由表8(MLRU、SRU、RRU的下標U表示無人機數據反演)可以看出，3種反演土壤表層鹽分的模型都具有統計學意義(P<0.001)，且均表現為極顯著，進一步表明了高分辨率多光譜信息可以進行表層土壤鹽分的預測。但是模型間也存在一定的差異，其中RRU模型的效果相對較差，建模集和驗證集R2最小，分別為0.403和0.369，RMSE最大，為0.237%。MLRU模型和SRU模型的建模和驗證效果相對較好，其中建模集和驗證集R2都在0.45以上，

表8 土壤鹽分與無人機變量因子回歸模型Tab.8 Monitoring models of soil salt and UAV variable factor

RMSE相差較小，分別為0.195%和0.202%。但是SRU模型的F明顯高于MLRU模型，且SRU模型縮減變量因子，可在一定程度上減弱共線性問題、減少了計算量，可簡潔高效地得到表層土壤含鹽量情況，故SRU模型為本次無人機遙感反演土壤含鹽量的最優模型。

結合表7、8可得，基于無人機遙感數據建立的裸土期表層土壤鹽分反演模型精度均優于GF-1衛星數據建立的表層土壤鹽分反演模型。利用無人機遙感數據建立的土壤鹽分反演最優模型為SRU模型，利用GF-1衛星遙感數據建立的土壤鹽分反演最優模型為MLRS模型。

2.4 轉換殘差分析

由于無人機尺度上的轉換殘差等于衛星尺度的轉換殘差，因此將無人機數據在3種回歸模型下生成的轉換殘差與衛星數據變量因子進行分析。經過對變量因子的篩選，選取NDVI-S1和B1進行轉換殘差的擬合，如圖3所示。由圖可以看出，RRU轉換殘差與(NDVI-S1)-B1的擬合效果總體相對較好，決定系數R2為0.586 3，其他兩種轉換殘差下，決定系數R2分別為0.456 7和0.523 3。通過圖3擬合的多項式生成衛星尺度上轉換殘差的面圖像，加上在無人機尺度上建立的趨勢面應用在GF-1衛星尺度生成的面圖像，設為升尺度后的圖像。

圖3 轉換殘差與光譜指數擬合圖Fig.3 Fitting graphs of conversion residual

2.5 升尺度分析

2.5.1升尺度定性分析

升尺度后的土壤鹽分影像圖與表7所示的MLRS模型反演的土壤鹽分圖進行密度散點圖擬合，如圖4所示，可見二者具有一定的擬合度，R2為0.204。為了進一步定性評價升尺度后的土壤鹽分影像圖，本文通過直方圖來對比升尺度后的土壤鹽分影像圖與MLRS模型反演的鹽漬化影像圖之間的關系，如圖5所示。

圖4 轉換前后土壤鹽分密度散點圖Fig.4 Soil salt density scatter plot before and after scaling

由圖5可以看出，MLRS模型反演的鹽漬化影像圖，含鹽量在0～1.75%之間，峰值出現在0.4%附近，且像元數量為2 000左右；而升尺度后的鹽漬化影像圖，含鹽量在0～2.0%之間，峰值位于0.25%附近，像元數量為6 500左右。可見相比MLRS反演的鹽漬化影像圖，升尺度后的土壤含鹽量值相對偏小，降低了試驗區的鹽漬化情況。

圖5 轉換前后土壤鹽分直方圖Fig.5 Histograms of soil salt before and after scaling

2.5.2升尺度定量分析

為了對升尺度結果進行定量評價，分別對模型計算值與實測值進行了對比分析，如表9所示。

由表9可以看出，模型1和模型2相差不大，建模集和驗證集R2都約為0.72，且RMSE和MAE都相近，可見兩者具有相同的建模和預測效果。模型3相對最優，雖然建模集和驗證集R2與其他兩種模型較為接近，但RMSE和MAE在建模集和驗證集下都最小。

表9 預測值與實測值對比分析Tab.9 Comparison of predicted and true values

注：ΔS1、ΔS2、ΔS3分別為不同模型下的土壤鹽分轉換殘差。

2.6 土壤鹽分預測

根據升尺度后最優的模型3進行試驗區表層土壤的含鹽量估測，通過鹽漬土劃分標準，將試驗區域劃分為4類鹽漬化類型，即非鹽漬土(0～0.2%)、輕度鹽漬土(0.2%～0.5%)、重度鹽漬土(0.5%～1.0%)、鹽漬土(>1.0%)[15]，如圖6所示。

圖6 升尺度16 m鹽漬化分級圖Fig.6 Classification map of salinization by ascending scale of 16 m

3 討論

土壤鹽漬化問題嚴重制約當地農業經濟的發展，對其進行監測十分重要，遙感技術能夠較好地解決這一問題。本文對GF-1衛星數據和無人機數據的各波段和光譜指數與實測鹽分數據進行了相關性分析，發現兩者的藍波段、近紅外波段、鹽分指數SI、鹽分指數S5與土壤鹽分的相關性較好，這與馮娟等[29-30]研究成果一致。新構建的光譜指數NDVI-S1與土壤鹽分的相關性可以提高0.03左右。通過3種回歸模型對試驗區進行了土壤含鹽量的遙感反演，發現無人機數據監測模型精度均優于GF-1衛星數據，原因是GF-1衛星數據分辨率低且受外界因素影響更大，無人機數據在獲取時，可以人為選擇較好的外界環境來進行飛行試驗，得到成像品質較好的影像數據。

本文通過改進的TsHARP尺度轉換方法，提高GF-1衛星遙感對土壤鹽漬化的監測精度，具有較好的效果。這是因為TsHARP尺度轉換方法有效結合衛星數據和無人機數據，通過無人機尺度上生成的轉換殘差對GF-1衛星反演模型進行校正，獲得像元尺度土壤鹽分的“相對真值”，提高了衛星遙感數據對定量遙感反演模型的可靠性。與傳統基于重采樣的尺度上推方法相比，引入高分辨率下的光譜信息構造趨勢面實現升尺度轉換，能較好地保持原始遙感數據的信息量及空間結構特征，可對今后的尺度上推研究提供一定的參考。

先前的研究表明，地表異質性和模型非線性均是尺度效應產生的原因[31]，本文沒有考慮尺度效應問題，只是通過尺度轉換方法實現GF-1衛星監測土壤鹽分精度的提高，這是本文的不足之處。另外，試驗區水熱耦合情況復雜且多變，制約了土壤鹽分的有效監測，因此本文所得的結果也僅限于本次試驗，在其他月份和其他地區的適用性還有待進一步探索。

4 結論

(1)通過土壤鹽分與無人機數據和GF-1數據的光譜波段、鹽分指數進行相關性分析發現，兩者的藍波段B1、近紅外波段B5、鹽分指數SI、鹽分指數S5和新構建的光譜指數NDVI-S1與土壤鹽分具有較好的相關性，相關系數均在0.3以上，其中NDVI-S1可達到0.47以上。

(2)對比2種數據源在3種回歸模型下反演表層土壤鹽分時發現，無人機數據建立的模型精度均優于GF-1衛星數據；在無人機高分辨率數據下，逐步回歸模型相對最優，在GF-1衛星低分辨率數據下，多元線性回歸模型相對最優。