楊麗萍 任 杰 王 宇 張 靜 王 彤 李凱旋
(長安大學地質工程與測繪學院,西安 710054)
土壤鹽漬化(Soil salinization)是破壞土地資源、嚴重威脅生態系統安全與穩定的全球性問題[1]。中國西部地區鹽漬土面積占全國鹽漬化土地總面積的69.03%[2],土壤鹽漬化問題突出。黑河流域的鹽漬化問題在中國西部干旱區具有典型性和代表性,其下游的額濟納盆地尤為嚴重,快速準確地獲取大范圍的土壤鹽分信息對鹽漬化監測與防治具有重要意義。
遙感技術具有數據獲取速度快、覆蓋范圍廣、成本低等特點,已被廣泛用于土壤鹽漬化監測研究。基于光譜波段或植被指數等[3]單一參數的鹽漬化信息提取精度有限,隨著遙感數據種類的日益豐富,多源、多參數聯合的反演方法得到了廣泛應用。
光學遙感通過地物反射輻射特征的變化,模擬各種遙感指數、地表溫度等參數與土壤鹽分之間的關系,實現土壤鹽分監測。文獻[4-6]基于光學影像,通過多參數聯合的方式提升了鹽分信息的提取精度,但未考慮各類參數的鹽分反演效率和相互比較優勢。
雷達遙感具有全天時、全天候成像的特點,且微波具有一定的穿透能力,可以與光學遙感優勢互補,在土壤鹽漬化監測方面發揮著重要作用[7-8]。近年來,全極化雷達遙感技術快速發展,全極化雷達影像由于包含了散射目標的幾何特征、后向散射特征及豐富的極化特征信息,在土壤鹽分反演中備受關注。賈殿紀[9]基于Radarsat-2全極化雷達影像,通過構建土壤鹽分與極化特征參數間的偏最小二乘模型,實現了青海湖流域的土壤鹽分反演,效果良好。然而現階段利用全極化雷達影像開展土壤鹽分定量反演的研究雖然取得了一定的進展,但工作較為有限。
光學遙感不僅可以通過光譜信息反映土壤鹽漬化狀況,同時也有助于微波遙感信息的解譯和定量反演,而全極化雷達遙感能夠提供豐富的地物特征信息,因此,集成光學與全極化雷達遙感技術的土壤鹽分協同反演具有十分重要的科學意義與應用價值[10]。在利用多源遙感數據的鹽分反演中,前人采用了眾多的鹽分指示環境變量,如波段反射率[6]、植被指數、鹽分指數[4]、后向散射系數[7-8]、極化特征參數[9,11]、地表溫度[5-6]以及地形因子[6,12]等,這些變量在鹽分反演效率和相互比較優勢方面,具有一定程度的不確定性和易混淆性。因此,有必要對眾多鹽分指示環境變量進行全面的綜合評價,以確定鹽分監測的優勢敏感變量,進而建立鹽分監測的有效模型。在利用多源遙感參數反演土壤鹽分的方法上,RF和SVM等機器學習算法因具有能夠處理復雜的非線性問題和不受制于輸入環境變量的類型與個數的特點而被廣泛應用,且取得了較好效果[4,6],為土壤鹽分的定量反演研究提供了切實可行的思路。
本文以黑河下游額濟納旗東南的居延澤地區為研究區,基于Sentinel-2、Radarsat-2、Landsat-8和SRTM DEM數據提取波段反射率、植被指數、鹽分指數、極化雷達參數、地表溫度和地形因子共6類82個變量,采用變量優選策略篩選各類變量及其組合的最優變量,完成鹽分敏感變量的優選,進而構建RF與SVM鹽分預測模型,以期實現居延澤地區土壤鹽分的定量反演。
居延澤位于內蒙古自治區阿拉善盟額濟納旗的東南部,北鄰阿爾泰山脈,南與巴丹吉林沙漠相接,由東、西居延澤2個子盆地構成(圖1),西北角存有殘留湖沼——天鵝湖。該區域地勢東北高,西南低,平均海拔900 m,由于深居歐亞大陸腹地,氣候干燥,降雨稀少,為典型的大陸型氣候。

圖1 研究區示意圖
根據研究區的土壤類型和土地覆蓋類型,考慮交通可達性,于2017年8月14—17日共采集42個鹽分樣點(圖1)。在各采樣點人工開挖長×寬×深約50 cm×50 cm×50 cm的探坑,采集0~10 cm的土壤樣品。土壤樣品經風干、研磨,稱取20 g樣品(已過1 mm篩),加100 mL超純水5倍稀釋,封口后上搖床3~5 min,過濾,并使用儀器測定過濾后原液的濃度,依據原液質量濃度是否大于500 mg/L考慮配置稀釋液,最后用sensION5型電導率儀測定濾液或稀釋液的電導率,并轉換為土壤含鹽量(SSC)[13]。
本文采用光學與雷達影像。光學影像采用Sentinel-2和Landsat-8影像,Sentinel-2影像成像時間為2017年8月24日,數據來源于哥白尼數據中心(https:∥scihub.copernicus.eu/),數據等級為L1C級,空間分辨率為10、20、60 m。利用配套處理插件Sen2Cor將經過輻射校正與幾何校正的L1C級數據轉為L2A級數據,完成大氣校正等預處理。Landsat-8影像成像時間為2017年9月18日,數據來源于地理空間數據云(http:∥www.gscloud.cn),空間分辨率為30 m。訂購了野外期間一景精細全極化模式Radarsat-2雷達影像,成像時間為2017年8月17日,幅寬25 km×25 km,空間分辨率8 m。使用ENVI軟件的SARScape模塊和PolSARpro軟件完成原始SLC影像的多視、濾波、地理編碼、輻射定標和正射校正等預處理。DEM數據采用SRTM DEM,空間分辨率為30 m,數據來源于美國地質調查局(http:∥www.usgs.gov/)。為便于后續分析,以上影像均重采樣至10 m。
前人研究表明,植被指數、鹽分指數以及地形因子可以為土壤鹽漬化監測提供有效信息[6,14]。地表溫度影響水分中可溶性鹽的析出積聚,是干旱半干旱區土壤鹽漬化監測的重要參數[15]。極化目標分解是從極化SAR數據中提取目標散射特征的方法,根據不同的分解機制,前人提出了An&Yang分解[16]、Cloude-Pottier分解[17]、Freeman-Durden分解[18]、H-A-α分解[19]、Krogager分解[20]和Yamaguchi分解[21]等方法。分解后的特征物理含義明確,在一定程度上可體現不同目標之間散射機理的差異,能夠提供豐富的地物特征信息。此外,土壤鹽分的積累會導致土壤介電常數虛部發生變化,進而引起雷達后向散射系數的變化[22]。因此,極化雷達的后向散射系數及極化分解后的極化特征參數也是土壤鹽漬化監測的重要參數。鑒于此,本文選取波段反射率、植被指數、鹽分指數、極化雷達參數、地表溫度和地形因子共6類82個變量參與土壤鹽分預測建模。其中,波段反射率、植被指數、鹽分指數由Sentinel-2數據獲取,極化雷達參數由Radarsat-2數據獲取,地表溫度由Landsat-8數據獲取,地形因子利用DEM數據獲取。具體變量及計算方法見表1。

表1 環境變量
在土壤鹽分反演中,并不是所有環境變量均可以提供有用信息而參與模型構建,模型應該避免不相關變量的影響[32]。研究表明,基于迭代刪除潛在不相關的環境變量,有利于減少不確定性,提高預測精度[33]。本文參與模型構建的變量以各類變量及其組合進行分組,每組均會根據迭代刪除產生最優變量集合。變量優選策略參照文獻[32-33],具體為:
(1)將每組變量輸入RF算法,基于平均精度減少量(Mean decrease accuracy,MDA),對環境變量進行重要性排序。
(2)刪除重要性排序的最后一個環境變量,將剩余的變量輸入RF算法重新訓練并排序,隨后再刪除最后一位變量。重復刪除最不重要的變量,直至剩余2個環境變量,循環結束。
(3)每次算法運行均會計算出決定系數(Determination coefficient,R2)和均方根誤差(Root mean square error,RMSE),最終以RMSE為主、R2為輔綜合判斷最優變量集合。
(4)重復步驟(1)~(3),直到遍歷所有分組,得到每組的最優變量集合。
RF是2001年BREIMAN[34]提出的一種機器學習算法,以CART決策樹為弱學習器,通過自舉法(Bootstrap)重采樣技術,從原始訓練樣本集N中有放回地重復隨機抽取n個樣本生成新的訓練樣本集,重復采樣生成多個新的訓練樣本集,用新的樣本集訓練決策樹組成RF,并將所有弱學習器得到的結果進行算術平均作為最終的預測結果。其中,在每個弱學習器構建Bagging集成、生成決策樹時,所選特征均為隨機選出的少數特征,從而使RF不需要額外剪枝即可取得較好的泛化能力和抗過擬合能力。該算法僅需定義2個參數[33],即終端節點樹和每個節點隨機選擇的變量個數。通過反復計算與比較,設置終端節點樹ntree為1 500,選擇的變量個數m為自變量數目的1/3。
SVM是以VC維(Vapnik-Chervonenkis dimension)和結構風險最小化為理論基礎的機器學習算法,該算法基于核函數映射將低維空間非線性問題轉換為高維空間的線性問題,利用一個超平面根據最大化間隔值,將輸入數據劃分到n維的特征空間,從而實現分類、回歸任務。通過引入不敏感函數,SVM在曲線擬合中得以應用并發展為支持向量回歸。作為模型的重要參數,核函數與懲罰系數c、參數g的選取對模型精度有很大影響[35]。本文選擇徑向基函數作為核函數,c、g通過交叉驗證來計算,并通過調用基于libsvm網格劃分(Grid search)的參數尋優函數來尋找最優解。
交叉驗證對機器學習有很強的指導與驗證意義,K折交叉驗證[36](K-fold cross validation)是廣為使用的交叉驗證方法之一,尤其在數據集較小的情況下,能夠最大程度地提高數據利用率,使模型性能更為優越。本文采用10折交叉驗證法,將樣本劃分為10個互斥子集,每次將9個子集作為訓練集以擬合模型,剩余1個子集作為驗證集用于評估所建模型的預測能力,驗證過程重復10次,每次將不同的子集作為驗證集,故每個模型會得到10次模型驗證結果,取其平均值,作為該模型的最終驗證結果。本文采用R2和RMSE為指標評價預測模型的效果與性能。當R2趨近1,RMSE趨近0,表示模型效果最佳。
研究區實測土壤含鹽量在0.201~119 g/kg之間,土壤含鹽量變化較大,平均值為23.898 g/kg,變異系數為1.175,屬于強變異性。根據土壤鹽分分級標準[37],鹽漬土可分為非鹽漬化(含鹽量小于1 g/kg)、輕度鹽漬化(含鹽量1~2 g/kg)、中度鹽漬化(含鹽量2~4 g/kg)、重度鹽漬化(含鹽量4~10 g/kg)和鹽土(含鹽量大于10 g/kg)5個等級。參考以上標準,研究區實測鹽分等級分布如表2和圖2所示。由表2可知,該地區非鹽土樣本與鹽土樣本約各占50%,而非鹽土樣本又以中度和重度鹽漬化土居多。在圖2中,小提琴寬度表示數據分布密度,越寬表示數據分布密度越大,高度表示土壤鹽分分布范圍。

表2 實測鹽分等級分布

圖2 實測鹽分分布圖
利用Random forests軟件包對變量進行了重要性分析。就所有變量(圖3,圖中以平均精度減少(MDA)衡量環境變量重要性)而言,地表溫度的重要性最高,為12.416%,且遠遠高于其他變量;Kkh的重要性次之,為4.121%;地形深度指數TWI的重要性居第三,為3.984%;波段反射率、植被指數及鹽分指數的重要性靠后且均小于2.3%。綜合考慮,6類變量的重要性從大到小依次為地表溫度、地形因子、極化雷達參數、鹽分指數、植被指數、波段反射率。通過變量優選策略進一步選擇各組最優變量,結果如表3所示。

圖3 變量重要性
由表3可知,優選的波段反射率重要變量有B6、B8a、B11和B12等,揭示了紅邊和短波紅外波段在土壤鹽分監測中的重要性,與TAGHADOSI等[6]的研究結果相一致。在植被指數中,CRSI、ERVI、ENDVI和SAVI與土壤鹽分存在密切關系,與周曉紅等[38]的研究中增強型植被指數與鹽分相關性最強的結果相一致,同時也與王飛等[39]在不同綠洲區得到的CRSI和增強型植被指數(ENDVI、EVI、EEVI)在干旱區具有普適性的結論相吻合。在鹽分指數中,優選的重要變量有SI2re3、S5re1、S1re2和S2re3等,而這些變量均有紅邊波段參與計算,再次證明了紅邊波段在土壤鹽分監測中的重要性,與馬國林等[12]的研究中RESI63、RESI31、RESI12和RESI16紅邊光譜指數重要性居前的結果一致。優選的極化雷達參數重要變量有An4Odd、An4Vol、FOdd、YOdd、Kkh、H和Kks等,其中,FOdd與Kks同樣在依力亞斯江·努爾麥麥提等[11]的研究中作為優選變量。各變量有明確的物理含義,可以在一定程度上體現不同目標之間散射機理的差異,尤其是重度與中-輕度鹽漬地之間的差異[40-41],對土壤鹽分預測有重要作用。地形因子優選的重要變量有CA、VD、CNBL和TWI等,以水文相關的地形因子為主,這與居延澤地區為古湖盆區有關,而CA和TWI也在馬國林等[12]的研究中作為重要性居前的優選地形指數。

表3 優選變量
以波段反射率、植被指數、鹽分指數、極化雷達參數、地表溫度和地形因子6類變量及其組合方案的優選變量作為RF和SVM模型輸入的自變量,土壤含鹽量作為目標變量,建立土壤鹽分RF和SVM預測模型,不同變量模型的鹽分預測精度如表4所示。

表4 土壤鹽分反演精度
由表4可知,在單一變量反演方案中,極化雷達參數與地形因子對土壤鹽分預測有較大貢獻,極化雷達參數(方案4)建立的RF模型驗證集R2為0.701,RMSE為21.837%,SVM模型驗證集R2為0.504,RMSE為21.771%。地形因子(方案6)建立的RF模型驗證集R2為0.704,RMSE為20.257%,SVM模型驗證集R2為0.519,RMSE為21.065%。其次地表溫度(方案5)對土壤鹽分反演也較為重要,RF驗證集R2稍低,為0.609,RMSE為21.302%,SVM驗證集R2為0.528,RMSE為23.538%。鹽分指數、植被指數與波段反射率(方案3、2、1)建立的RF模型精度略低且相差不大,R2分別為0.621、0.593、0.599,RMSE分別為28.453%、25.866%、26.326%,SVM模型趨勢與RF模型相似,R2分別為0.451、0.490、0.428,RMSE分別為29.055%、22.740%、22.296%,與陳俊英等[4]研究優選鹽分與植被指數組以及張智韜等[42]研究敏感波段組與光譜指數組對鹽分反演精度影響的結果相吻合。極化雷達參數包含了豐富的極化信息,不同地物的散射機理不同,反映的極化信息也不盡相同,尤其是重度與中-輕度鹽漬地[40],因而可以較好地預測土壤鹽分。居延澤地區為古湖區,湖水在匯聚、流動和退縮過程中會受到地形影響,并最終在低洼地帶匯集。湖水干涸后,地下水毛細上升使得地表鹽分大量析出。因此,土壤鹽分的分布與地形有密切關系。地表溫度同樣影響鹽分分布,土壤中的可溶性鹽隨湖水流動,在溫度與地形的影響下,鹽分朝著某一方向匯聚,隨著高溫條件下水分的快速蒸發,鹽分便會積聚在土壤表面。由于植物鹽脅迫的影響,植被指數在較高植被覆蓋下對土壤鹽分變化敏感[30],而本研究區極度干旱,植被稀疏,導致植被指數對土壤鹽分預測精度較低。鹽分指數與植被指數恰恰相反,由對鹽分敏感的波段組合計算得到,主要依據土壤的光譜響應來反映土壤鹽分變化,地表植被越少,對光譜響應影響越小,鹽分指數對土壤鹽分的監測效果越好,因而鹽分指數對土壤鹽分的預測精度高于植被指數。
在多變量反演方案中,波段反射率、植被指數、鹽分指數這3類變量聯合反演的精度幾乎沒有提升,波段反射率與植被指數聯合反演(方案7)的RF驗證集R2為0.619,RMSE為26.446%,SVM驗證集R2為0.468,RMSE為22.197%,而加入鹽分指數后(方案8),RF驗證集R2降低0.018,RMSE升高0.587個百分點,SVM驗證集R2降低0.229,RMSE升高6.775個百分點,或許是因為植被指數與鹽分指數都是由對鹽分敏感的波段計算的光譜指數,這3類變量本質上都是波段信息組合,相互之間存在信息冗余與重疊。而對植被指數與極化雷達參數聯合反演的方案13而言,RF驗證集R2比兩者單獨反演分別提升0.202、0.095,RMSE比前者單獨反演降低3.352個百分點,比后者單獨反演升高0.667個百分點,SVM驗證集R2比兩者單獨反演分別提升0.061、0.047,RMSE比兩者單獨反演分別升高1.719、2.688個百分點,說明光譜指數和極化雷達參數之間不存在信息冗余,體現出光學遙感與微波遙感信息互補的優勢。當光譜指數加入極化雷達參數、地表溫度和地形因子聯合反演時(方案9~11),RF驗證集R2分別提升0.048、0.087、0.202,RMSE分別降低1.584、6.354、8.745個百分點,SVM驗證集R2分別提升0.009、0.160、0.186,方案9的RMSE升高1.726個百分點,方案10、11分別降低0.747、3.674個百分點,表明極化雷達參數、地表溫度及地形因子包含豐富的信息,是土壤鹽分反演的重要參數,這與前文單一變量鹽分反演精度較高的結果是一致的,同時與6類變量的重要性排序結果吻合。當地表溫度與地形因子聯合反演時(方案12),RF驗證集R2為0.804,RMSE為16.748%,SVM驗證集R2為0.629,RMSE為18.271%,但在加入植被指數或鹽分指數后(方案14、15),RF驗證集R2分別下降0.102、0.139,RMSE分別升高2.034、2.823個百分點,SVM驗證集R2分別下降0.067、0.154,RMSE分別升高2.833、3.278個百分點,與極化雷達參數、地表溫度和地形因子聯合反演(方案16)加入植被指數(方案17)的變化趨勢一致,這或許是因為在聯合反演時,各變量間互相影響、互相牽制,植被指數和鹽分指數的加入對地表溫度及地形因子的聯合反演有一定的影響。
從不同方案的優選變量來看,各變量之間相互影響、相互制約,聯合反演并不是將重要性最高的變量組合在一起。從不同方案優選變量的個數結合預測精度也可以說明,反演時變量并不是越多越好,引入過多的自變量會產生信息冗余,造成過擬合,使模型預測精度下降[42]。由上文可知,RF模型的R2均高于SVM模型,RF模型的RMSE總體也均低于SVM模型,且兩者在各方案模型上的變化趨勢一致,表明RF模型比SVM模型對土壤鹽分的預測效果更好。張智韜等[42]在內蒙古河套灌區耕地的土壤鹽分反演中也用到了本文的2種機器學習模型,發現RF模型效果最好,SVM模型次之。馬國林等[12]基于多光譜與DEM數據反演艾比湖濕地土壤鹽分時,同樣發現機器學習RF模型的反演效果更優,上述研究均與本文的模型研究結果一致。對表4進一步分析發現,在RF與SVM 2種模型中,基于方案11、12、16構建的模型均具有較高精度且模型精度相差較小,其中方案11引入了波段反射率、植被指數、鹽分指數、極化雷達參數、地表溫度和地形因子,方案12引入了地表溫度和地形因子,方案16引入了極化雷達參數、地表溫度和地形因子,但方案11包含6類變量,信息豐富,且較單一變量最佳模型R2提升0.117,RMSE降低2.556個百分點,綜合考慮認為基于方案11構建的RF模型為研究區土壤鹽分預測的最優模型。
利用選取的最優變量集和RF模型反演居延澤地區土壤鹽分,結果如圖4所示。由圖4可見,土壤鹽分由研究區東北向西南增高,西南方向的鹽分又以東、西居延澤古湖盆區為最。研究區地勢由東北向西南傾斜,東北地區地勢高,為山區與戈壁,鹽分含量較低。由于古湖的退縮,致使西南地區成為水鹽匯聚中心,在高溫干燥環境下,湖水迅速蒸發并干涸,地下水毛細上升使得地表鹽分大量析出,形成低洼地帶鹽分高的分布特征。天鵝湖附近灘涂眾多,旱生植被集中生長,涵養了一定水分,導致含鹽量較低。在西居延澤中部出現斑塊狀的低鹽分區,可能是因為古湖在干涸前涵養了周邊大量的植被。圖中出現的條帶狀低鹽區,為東、西居延澤之間表面覆蓋小礫石的古湖岸線,零星斑塊狀低鹽區為雅丹地貌和灌叢沙堆分布區。東居延澤中存在的大塊低鹽區,是由于該區為地勢較高的沙土區域,含鹽量較低。從分布格局而言,模型預測結果符合野外調查結果,與實地情況較為接近。

圖4 土壤鹽分預測結果
(1)眾多鹽分指示變量中,短波紅外波段(B11)、冠層鹽度響應植被指數(CRSI)、擴展比值植被指數(ERVI)、紅邊鹽分指數(S2re3)、單次散射(FOdd)、地表溫度(LST)與匯水面積(CA)等變量對土壤鹽分監測具有較強的普適性。
(2)單一變量模型的鹽分預測精度根據R2和RMSE綜合判斷,從高到低依次為地形因子、極化雷達參數、地表溫度、鹽分指數、植被指數和波段反射率,其中地形因子構建的RF模型預測精度最高,其驗證集R2為0.704,RMSE為20.257%。
(3)多變量聯合可以進一步提升模型的預測精度與穩定性,隨著環境變量的逐步加入,當6類變量均參與模型構建時,RF模型預測精度最高,其驗證集R2為0.821,RMSE為17.701%,與單一變量最佳模型相比,R2提升0.117,RMSE降低2.556個百分點。
(4)RF模型較SVM模型更適于干旱區土壤鹽分反演,優選全變量組構建的RF模型具有最佳預測精度,其驗證集R2為0.821,RMSE為17.701%。基于該模型的反演結果表明,區域東北及天鵝湖附近鹽漬化程度較低,西南部古湖盆區鹽漬化程度較高。