趙建輝,張晨陽,閔 林,李 寧,王穎琳
(1.河南大學計算機與信息工程學院,開封 475004;2.河南省大數據分析與處理重點實驗室,開封 475004;3.河南省智能技術與應用工程技術研究中心,開封 475004;4.河南大學信息化管理辦公室,開封 475004)
土壤水分對于植被生長、農業生產和生態系統循環具有重要影響。土壤水分監測在氣象、水文、農業等多學科中發揮著重要的作用[1-2]。在現代農業中,土壤水分更是不可或缺的重要指標,全面監測土壤水分可以對農作物產量、旱情墑情和農作物長勢起到指導和決定性作用。因此,研究農作物覆蓋下的地表土壤水分分布情況具有重要意義[3-5]。
在現有土壤水分監測方法中,遙感技術已成為一種重要的監測手段,光學遙感和合成孔徑雷達(Synthetic Aperture Radar,SAR)微波遙感在土壤水分監測中得到越來越多的應用[6-8]。韓玲等[9]使用ASAR雙極化雷達數據,利用高級積分方程模型(Advanced Integral Equation Model,AIEM)模擬地表后向散射,提出一種組合粗糙度,構建土壤水分反演模型,并通過粒子群算法求解得到較為準確的土壤水分值。楊貴軍等[10]使用 Radarsat-2全極化數據,將水云模型中植被參數改為雷達植被指數,構建改進的水云模型,得到精度較高的冬小麥覆蓋下的地表土壤水分反演結果。王樹果等[11]采用Sentinel-1 SAR數據和Landsat8多光譜數據,結合Oh模型、水云模型以及復型洗牌全局優化算法,發展了一種土壤水分、植被含水量和地表粗糙度協同反演方案,可以得到與地面實測情況較一致的反演結果。地表參數與植被指數、雷達后向散射系數之間存在不易描述的非線性關系,為了提高反演精度,往往導致土壤水分反演模型參數眾多、結構復雜。神經網絡具有極強的非線性擬合能力,并且可以自主學習,在解決土壤水分反演過程中的非線性問題時得到了越來越多的應用。韓穎娟等[12]從風云衛星中獲取溫度和波段信息,構建卷積神經網絡反演地表土壤濕度,試驗結果精度較高。郭交等[13]使用 Sentinel-1和Sentinel-2多源遙感數據,利用支持向量回歸(Support Vector Regression,SVR)和廣義神經網絡(Generalized Regression Neural Network,GRNN)反演土壤水分,收斂速度快,精度高。余凡等[14]使用Envisat ASAR數據和TM數據,將雙極化數據與光學遙感數據相結合,采用遺傳算法(Genetic Algorithm,GA)優化的BP神經網絡反演作物覆蓋區土壤水分,反演結果與實際情況基本吻合。
在利用神經網絡進行土壤水分反演時,選擇合適規模的訓練數據非常重要,過多的網絡輸入數據會導致收斂速度變慢,對反演精度產生負面影響。對特征集合進行降維可以有效提高預測準確性,以便于構造效率更高、消耗更低的預測模型。實現特征降維有特征選擇和特征抽取 2種方式[15]。差分進化特征選擇(Differential Evolution Feature Selection,DEFS)算法是一種行之有效的特征選擇方法,可以從數據集中篩選出最優特征子集。Khushaba等[16]將DEFS技術應用于人類腦電圖的分類,并與其他降維技術進行了比較,結果表明該算法在求解最優性、內存需求和計算代價等方面具有重要意義。主成分分析(Principal Component Analysis,PCA)是一種常見的特征抽取方法,可以將數據集映射到低維空間,對數據進一步壓縮,剔除冗余數據。王雅婷[17]從Radarsat-2數據中提取12個極化特征參數,利用PCA方法進行降維處理,然后再進行支持向量回歸建模,取得了較好的反演精度。
本文結合Sentinel-1微波遙感數據和Sentinel-2光學遙感數據,經過預處理后從中提取 21個特征參數,使用DEFS算法對特征數據進行篩選,得到一定數量的最優特征子集,之后使用PCA算法再次進行降維;然后構建經過GA算法優化的BP神經網絡(簡稱GA-BP神經網絡),使用經過特征降維后的數據集和部分實測數據對網絡進行訓練,并使用訓練好的網絡對研究區的土壤水分進行反演;最后利用實測數據對反演結果精度進行對比驗證。
研究區位于河南省開封市祥符區,面積約900 km2,大致范圍為 34°36'~34°51'N,114°30'~114°45'E,如圖1所示。研究區為黃河沖積平原的一部分,地勢平坦,屬溫帶大陸性季風氣候,冬季寒冷干燥,夏季高溫多雨,平均海拔65 m,年平均氣溫14 ℃,年降水量628 mm,無霜期214 d。主要種植冬小麥、玉米、棉花等作物,其中冬小麥的生長周期為8個月,一般10月播種,次年6月收獲。冬小麥在不同物候期的生物量和生物特征不同,地表土壤水分在不同時間、不同地點也有較大差異。試驗處于冬小麥的出苗期、分蘗期和越冬期,這 3個物候期內冬小麥植株較小,植被高度及覆蓋度均較低,地面植被覆蓋情況變化不大,且農田已不再發生犁地、播種等影響地表粗糙度的田間活動。因此,本文針對這 3個相近的物候期進行統一建模和分析,使用這 3個物候期的遙感數據和地面實測數據進行神經網絡訓練、土壤水分反演與結果精度評價。
1.2.1 地面實測數據
土壤水分遙感反演研究中使用的地面實測土壤水分數據主要有兩種來源,一種是來源于研究區內的地面觀測站點或自動觀測網絡[11-12,17-18],此類實測數據獲取方便且通常采集頻率較高、采集數量較大,以此為基礎開展的土壤水分反演研究中反演算法或反演模型的可選種類較多、選擇空間較大;另一種是來源于傳統的人工測量方法[5-6,13-14],在衛星過境日期依靠人工進行地面采樣測量,此類實測數據往往獲取困難且通常采集次數有限、采集數量較少,以此為基礎開展的土壤水分反演研究中反演算法或反演模型的可選種類相對較少、選擇空間相對較小。本文研究區內沒有地面觀測站點和自動觀測網絡,因此采用人工測量方法來獲取地面實測數據,并開展基于小樣本量實測數據的土壤水分反演研究。
在Sentinel-1A衛星過境的時間,同步進行3次野外實地考察和采樣,現場采集土壤水分值和經緯度坐標。研究區地面共設置20個采樣點,采樣點分布如圖1b和圖1c所示。
3次采樣共采集60組實測數據,去除2個異常點后,共有58組有效實測數據用于后續試驗。野外采樣過程中使用TDR350土壤水分儀測量農田表層土壤體積含水量,探針長度為3.8 cm,在每個采樣點以十字測量法測5個點的土壤水分值,以 5個測量點的土壤水分平均值作為該采樣點的最終土壤水分實測值。使用室外手持式集思寶UG905定位儀(定位精度1~3 m)定位采樣點,選用WGS84坐標系記錄采樣點坐標。
1.2.2 遙感數據
本文使用歐空局(European Space Agency,ESA)提供的遙感數據,如表1所示。使用歐空局研發的哨兵應用平臺(SentiNel Application Platform,SNAP)軟件對所獲取的SAR圖像進行輻射定標、多視、Refined Lee濾波和地形校正等預處理操作。根據Sentinel-1A SAR圖像獲取日期和是否出現云霧雨等影響土壤水分大幅度波動的天氣因素,選擇相近日期的3景準同步光學圖像作為試驗數據。在SNAP軟件中使用Sen2Cor插件對所獲取的多光譜成像儀(Multi Spectral Image,MSI)圖像進行大氣校正等預處理操作。

表1 遙感數據Table 1 Remote sensing data
SAR通過向地物發射微波波束和接收回波信號來探測地物特性,波長、入射角和極化方式等雷達系統參數和目標地物的介電常數、物理結構等特征參數對雷達信息具有直接的影響,從雷達數據中提取可以表征地物特性的特征信息是雷達遙感反演的基礎。
2.1.1 后向散射系數
土壤水分反演時主動微波遙感主要通過后向散射系數所反映的信息進行反演。依據采樣點的經緯度,從預處理之后的 SAR數據中提取相應位置的入射角(θ)、VV極化后向散射系數()、VH極化后向散射系數()作為后續試驗的特征參數。由于cos(θ)和sin(θ)與土壤濕度也存在一定的關系[18],并且在雷達入射角一定的情況下其后向散射系數僅與地表粗糙度有關[19],所以將cos(θ)、sin(θ)和也作為特征參數,并將加入其中,從 SAR數據中提取共計 9個與雷達后向散射系數相關的特征參數。
2.1.2 極化特征參數
極化分解可以將地物較為復雜的散射過程分解成若干簡單的散射機理。通過極化分解的方式,可以從 SAR遙感數據中提取更多的特征參數[20]。對雙極化Sentine1-1A數據采用H/A/α分解,對目標地物的相干矩陣或者協方差矩陣進行特征值分解,可以從中提取出表征目標散射極化程度的極化熵H,表征目標散射機理的平均散射角α,極化熵的補充參數—反熵A,以及可以表示該散射機制強度的特征值λ1和λ2[21-22],從SAR數據中提取共計5個極化特征參數。
2.1.3 植被指數及地表粗糙度
植被和地表粗糙度對地表散射信息具有直接的影響,土壤水分反演需要有效抑制植被覆蓋和表面粗糙度的影響。光學遙感數據中可以提取的特征參數主要為植被指數。植被指數(Vegetation Index, VI)是2個或多個波長范圍內的地物反射率組合運算,以增強植被某一特性或者細節。目前遙感領域提出的植被指數有100多種[23],受限于傳感器類型和所用波段組合,不同的植被指數有不同的波段適用范圍和應用領域。基于Sentinel-2A攜帶的多光譜成像儀所提供的多波段數據,結合本研究區的實際植被覆蓋情況,本試驗最終選取了土壤水分反演研究中常用的 6個植被指數[24-26],包括歸一化差異植被指數(Normalized Difference Vegetation Index, NDVI)、歸一化差異水分指數(Normalized Difference Water Index,NDWI)、比值植被指數(Ratio Vegetation Index, RVI)、水分脅迫指數(Moisture Stress Index, MSI)、水波段指數(Water Band Index, WBI)和融合植被指數(Fusion Vegetation Index, FVI)等,其計算公式如公式(1)至公式(6)所示。
式中ρ842、ρ665、ρ1610、ρ865、ρ945分別表示 Sentinel-2 數據中對應中心波長為842、665、1 610、865、945 nm的波段值。
使用Sentinel-1A數據計算地表組合粗糙度[27],如公式(7)所示,其中使用非線性最小二乘法和線性回歸法擬合獲得C波段的Av和Bv,如公式(8)和公式(9)所示。
式中Zs為組合粗糙度,Av和Bv是只與入射角有關的系數。
在使用神經網絡進行預測時,過多的輸入數據可能會造成數據災難,影響神經網絡的學習,通常會采用特征選擇和特征抽取這兩種方法將數據降到更低維度,達到去除冗余數據的目的。特征選擇是單純地從提取到的所有特征中有依據地選擇較為重要的部分特征作為訓練集特征,被選擇的特征可以大概表征全體數據的信息,特征在選擇前后并不改變本身值的大小。而特征抽取的本質上是從一個維度空間映射到另一個維度空間,映射過程相當于一個黑盒,沒有確定的篩選依據,而是借助數學工具來進行降維,特征抽取后的特征失去了本身的物理意義,相應特征值也會改變。
本文使用 DEFS算法[16]對所提取的特征參數進行初步篩選。算法的第一步是從初始種群生成新的種群向量;第二步是變異操作,隨機選取群體中的兩個不同個體,通過差分進化實現個體突變,將其向量差進行縮放后,對要變異的個體執行向量合成;第三步是交叉操作,變異向量與原始矩陣中占據該位置的原向量交叉,這個操作的結果稱為試驗向量。新群體中的相應位置將包含試驗向量或原始目標向量,這取決于其中哪一個達到了更高的適應度;第四步是去除冗余,使用輪盤賭算法[28],通過與每個特征相關的分布因子計算單個特征的概率,來去除多余向量,獲得新種群。特征分布因子fi由公式(10)給出:
式中α1、α2是常數,ε是為了避免出現分母為零的情況而設置的一個極小的數,PDi是從高于平均精度的子集中計算到的正分布因子,NDi則恰恰相反,是從低于平均精度的子集中計算得到的負分布因子。
本文利用PCA方法對特征參數集合進行降維。PCA是數據分析中的一種重要方式,可以將數據中主要的特征變量抽取出來,常用于機器學習中對高維數據的降維。PCA將高維數據集映射到低維數據集,同時盡可能地保留原有信息。高維數據集中,各向量具有相關性,低維數據集中則線性無關,這樣便可以消除掉高維數據集中相互重疊的信息[29-30]。DEFS和PCA算法流程如圖2所示。
BP神經網絡已在許多領域得到廣泛應用,但存在易陷入局部極小、依賴于設計結構等缺陷,有時無法找到全局最優值。遺傳算法雖然不具備自學習能力,但具有全局尋優能力。因此,利用遺傳算法對神經網絡進行優化可以改善神經網絡的缺點,既發揮了神經網絡非線性映射能力和遺傳算法的全局尋優能力,又加快了神經網絡的學習速度,綜合提高了整個預測模型的精確度和擬合能力[31]。GA-BP神經網絡的構建流程大致分為2步:
1)構建BP神經網絡模型
建立包括輸入層、隱藏層和輸出層的BP網絡。確定每層神經元的個數、傳遞函數、BP網絡的具體參數和訓練次數。隱含層節點數s參照Kolmogorov定理[25]由公式(11)確定。
式中m為輸入層個數,n為輸出層個數。
2)使用GA算法優化BP神經網絡
遺傳算法優化 BP神經網絡的目的是通過遺傳算法得到更合適的網絡初始權值和閾值。首先對需要選擇的特征進行編號,將每個特征視為一個基因個體,使用初始化的BP神經網絡的預測誤差作為該個體的適應度值,通過選擇、交叉、變異操作尋找最優個體,即最優的BP神經網絡初始權值和閾值[32]。設置遺傳算法的初始參數時,交叉概率和變異概率一般在0到1之間取值。
本文使用遺傳算法優化BP神經網絡具體流程如圖3所示。
為了消除冗余特征對土壤水分反演結果精度造成的影響,本文提出了一種基于特征選擇和GA-BP神經網絡的土壤水分反演方法,技術路線如圖4所示,主要步驟如下:
1)對Sentinel-1和Sentinel-2數據進行預處理,從中提取前文所述的21個特征參數,并對所提取的特征參數進行編號,如表2所示。

表2 特征參數Table 2 Characteristic parameters
土壤濕度越大,后向散射系數與sin(θ)之間的相關性越高;土壤濕度越小,后向散射系數與cos(θ)之間的相關性越高[18]。在入射角一定的情況下,后向散射系數隨體積含水量的增大而增大,為不同極化后向散射系數的4種組合方式。除反熵A與土壤水分呈負相關外,極化熵H、平均散射角α、特征值λ1和λ2與土壤水分都呈正相關[17]。6個植被指數均與植被覆蓋程度呈正相關。
2)使用DEFS算法從21個特征參數中選出10個特征參數作為最優特征子集,包括NDVI、NDWI、WBI、FVI、α和H等參數。
在進行特征選擇時,所選特征參數所含信息量要大,并且特征參數之間的相關性要小。在使用 DEFS算法選擇最優特征子集的過程中,相關性較強的特征參數會因為相互之間存在重復的信息量而被去除掉,信息含量多且相關性較小的特征參數會予以保留。
3)使用 PCA方法對最優特征子集進一步降維。由于無法確定最合適的最優特征子集的個數,所以最優特征子集中仍可能冗余。使用PCA方法對特征值的累計貢獻率進行計算,發現前 8個主成分可以包含原有數據集信息的99.99%左右,所以選擇前8個主成分組成的特征矩陣作為后續的神經網絡輸入。
4)構建并使用 GA-BP神經網絡對土壤水分進行反演。本文的樣本數據量較小,而GA-BP神經網絡可基于小樣本數據進行訓練和反演,經過調參后可以保證模型的訓練和反演精度。將58個采樣點實測土壤水分值隨機分為兩組,其中50個用于神經網絡訓練,8個用于有效性驗證和精度評價。設置神經網絡的各個參數:輸入層個數為8;輸出層為1;根據公式(11)計算出隱含層的個數為5;學習率為0.1;迭代次數為5 000。設置遺傳算法的初始參數:迭代次數為 100;種群規模為 60;交叉概率為0.4;變異概率為0.1。
為了驗證本文所提方法的有效性,以河南省開封市祥符區冬小麥農田為研究區域,開展土壤水分反演試驗,并以GA-BP土壤水分反演模型[14]作為參照,對方法性能進行了對比分析。
試驗設置了3種方案:方案一采用GA-BP土壤水分反演模型,直接使用21個特征參數作為神經網絡輸入參數進行土壤水分反演;方案二在方案一的基礎上增加了DEFS算法,先對21個特征參數進行特征選擇,再使用GA-BP土壤水分反演模型進行土壤水分反演;方案三為本文所提方法,在方案二的基礎上再增加PCA方法,對DEFS特征選擇后的最優特征子集進行進一步降維優化,再使用GA-BP土壤水分反演模型進行土壤水分反演。采用偏差(Bias)、均方根誤差(Root Mean Square Error,RMSE)、無偏均方根誤差(unbiased Root Mean Square Error,ubRMSE)、決定系數R2共4個評價指標對反演精度進行評價。為了降低試驗結果偶發性,試驗結果均為多次試驗后求取的平均數。
表3為不同試驗方案反演結果精度對比,圖5為不同試驗方案土壤水分反演結果與實測值對比。由表3和圖5試驗結果可以看出,本文所提方法的反演值和實測值更為接近,Bias、RMSE和 ubRMSE比方案一、方案二更小,同時R2更高。并且在方案二中,僅使用 DEFS算法進行特征選擇后也比方案一的反演結果精度更高。本文方法反演結果的決定系數為 0.789 3,均方根誤差為0.028 7 cm3/cm3,相比單純使用GA-BP神經網絡,加入DEFS和PCA之后決定系數提高了0.215 7,同時均方根誤差降低了0.029 5 cm3/cm3。試驗結果表明,本文所提方法可以有效去除多余特征參數,提高土壤水分反演精度。

表3 反演結果精度對比Table 3 Accuracy comparison of inversion results
使用本文所提方法獲得的研究區農田土壤水分反演結果如圖6所示,其中為了去除非農田區域對土壤水分反演的影響、更好地顯示土壤水分分布情況,試驗中濾除了建筑、道路、河流等非農田區域,如圖6中白色區域所示。結果顯示2019年10月18日土壤水分反演值整體較高,均值為0.155 cm3/cm3,主要是由于10月上旬研究區多次降雨,土壤比較濕潤。10月30日研究區土壤水分反演值均值為0.136 cm3/cm3,比10月18日略低,與10月下旬天氣晴朗關系較大。2019年12月29日反演結果整體較干旱,均值為0.070 cm3/cm3,主要是因為入冬后溫度有時會降到 0 ℃以下,低溫會導致土壤含水量降低,寒風也會助長土壤水分的蒸發。經分析可知,這 3個日期的土壤水分反演結果與實際天氣情況比較吻合。此外,3個日期的采樣點土壤水分實測數據均值分別為0.162、0.136和0.065 cm3/cm3,反演結果與采樣點實測土壤水分值頻率分布較為一致,進一步驗證了本文所提方法的有效性。
本文基于Sentinel-1和Sentinel-2多源遙感數據,提取了21個與土壤含水量相關的特征參數,經過差分進化特征選擇(DEFS)和主成分分析(PCA)算法對特征參數進行篩選和降維后,結合地面實測數據,使用 GA-BP神經網絡反演土壤水分,并探討了DEFS和PCA算法對土壤水分反演精度的影響,主要結論如下:
1)不同特征參數所含的信息存在不同程度的重復和冗余,DEFS算法可以去除掉相關性較大、重復性較高的特征參數,保留信息含量多且相關性較小的特征參數。
2)在使用GA-BP神經網絡反演土壤水分的過程中,多余的特征參數會影響土壤水分反演結果的精度,組合使用DEFS和PCA算法可以剔除冗余特征參數,有效提高反演精度。本文方法反演結果的決定系數為 0.789 3,均方根誤差為0.028 7 cm3/cm3,相比單純使用GA-BP神經網絡,加入DEFS和PCA之后決定系數提高了0.215 7,同時均方根誤差降低了0.029 5 cm3/cm3。
雖然本文試驗考慮了地表粗糙度,但由于缺乏地面實測粗糙度數據,所以試驗中參考已有研究成果使用了從SAR數據中提取的地表粗糙度參數,這可能會影響反演精度。此外,相比于本文所用的雙極化SAR數據,全極化SAR數據中包含更多與土壤濕度相關的信息,可以提取更多的特征參數,使用全極化SAR數據有可能會進一步提高反演精度。在今后的研究中,可以考慮針對以上因素對試驗進行改進。