結(jié)合小波降噪與BP神經(jīng)網(wǎng)絡(luò)的GNSS-IR土壤濕度反演

吳昊艦，劉立龍，章傳銀，薛張芳，張志

(1.桂林理工大學(xué) 測繪地理信息學(xué)院，廣西 桂林 541004；2.廣西空間信息與測繪重點實驗室，廣西 桂林 541004；3.中國測繪科學(xué)研究院，北京 100036)

0 引言

土壤濕度即土壤中水分在土壤所占的比重。作為當(dāng)今農(nóng)業(yè)、氣象和水文研究中的一個重要參數(shù)[1]，土壤濕度發(fā)揮著無可替代的作用。傳統(tǒng)的地基測量方法[2](如土壤濕度計法、烘焙法以及時域反射法等)雖然精度較高，但是只適用于小范圍的測量并且容易受到外界環(huán)境因素的影響。通過利用全球衛(wèi)星導(dǎo)航系統(tǒng)(global navigation satellite system，GNSS)接收機來接收衛(wèi)星傳播過程中產(chǎn)生的信噪比觀測值來監(jiān)測地表環(huán)境，能夠有效彌補傳統(tǒng)方法中時空分辨率不足的缺點。基于此建立起來的全球衛(wèi)星導(dǎo)航系統(tǒng)干涉測量(global navigation satellite system interferometry and reflectomentry，GNSS-IR)遙感方法已經(jīng)發(fā)展成了一項全新的技術(shù)手段。

國內(nèi)外學(xué)者對于如何利用GNSS-IR方法獲取土壤濕度開展了大量的研究。Zavorotny等[3]通過建立多路徑電動力模型和裸土模型的反射系數(shù)進行實驗，證明了多徑干涉相位反演結(jié)果相對于幅度而言更加準(zhǔn)確。Jia等[4]提出利用隨機森林與支持向量機對實驗數(shù)據(jù)進行最小似然檢索，構(gòu)造了一個包含不同土壤類型的綜合模擬數(shù)據(jù)集，最終結(jié)果表明隨機森林方法更適合處理反演問題。Martin等[5]基于人工神經(jīng)網(wǎng)絡(luò)算法將輔助數(shù)據(jù)與反射率平均值進行不同組合，探究了地表粗糙度與植被衰減對反射率所造成的影響。萬瑋等[6]依據(jù)GNSS-R接收機的特征和原理將其分為了雙天線模式和單天線模式，并且給出了不同模式下的實驗，證明兩種模式的精度差異。梁月吉等[7]利用美國板塊觀測計劃(plate boundary observation，PBO)數(shù)據(jù)建立了土壤濕度多星線性回歸反演模型，通過提升聯(lián)合反演衛(wèi)星數(shù)量來提高反演精度。孫波等[8]提出了一種自適應(yīng)融合算法，將GPS不同頻段的多顆衛(wèi)星進行加權(quán)融合，相較傳統(tǒng)均值融合方法有效提高了反演精度。

然而上述研究大多采用的是傳統(tǒng)的一元線性回歸模型進行模型構(gòu)建，且實驗條件基本都是接近理想情況的。實驗中常用低階多項式擬合獲取信噪比殘差序列分離效果不夠理想，同時存在一系列問題如未曾考慮植被覆蓋的影響或者是地表粗糙度的散射影響。小波變換[9]作為一種表述信號時間與頻率的分析方法，具備多分辨率分析的優(yōu)點，能夠同時表達時間、頻率的局部特征，因此可以利用這一優(yōu)點來實現(xiàn)衛(wèi)星信號中直射信號的分離從而獲取信號趨勢項。BP(back propagation)神經(jīng)網(wǎng)絡(luò)則實現(xiàn)了一個從輸入到輸出的過程而且具備實現(xiàn)任何復(fù)雜非線性映射的能力，通過反向傳播來調(diào)整網(wǎng)絡(luò)中的閾值和權(quán)值從而使得網(wǎng)絡(luò)實際輸出與期望輸出達到最小值[10]，在不出現(xiàn)過擬合的情況下，對地表植被覆蓋和地表粗糙度引起的噪聲進行抑制?；诖?，本文建立了一種小波變換與BP神經(jīng)網(wǎng)絡(luò)(為便于后續(xù)描述將該模型簡稱為小波神經(jīng)網(wǎng)絡(luò)模型)相結(jié)合的GNSS-IR土壤濕度反演模型，并利用PBO H2O項目組以及積雪遙測系統(tǒng)(snow telemetry,SNOTEL)提供的土壤濕度數(shù)據(jù)作為精度評定的標(biāo)準(zhǔn)。

1 實驗基本原理與方法

1.1 地基GNSS-IR反演經(jīng)典模型

本文下載了PBO監(jiān)測網(wǎng)中P036測站的觀測數(shù)據(jù)以及附近氣候站的土壤濕度實測數(shù)據(jù)進行比較分析。圖1為GNSS接收機同時接收直射信號與反射信號的原理示意圖，其中H為天線高，即接收機天線相位中心到地球表面的垂直距離；AD為衛(wèi)星直射信號；AR為經(jīng)過地表反射后的反射信號；θ為衛(wèi)星高度角。GNSS接收機此時是同時接收衛(wèi)星直射信號AD與經(jīng)由地面反射的反射信號AR。

圖1 GNSS-IR反射原理示意圖

由圖1 可知，信噪比與振幅的關(guān)系如式(1)所示[11]。

(1)

式中：Ac為合成信號振幅；Ad與Ar分別為直射信號與反射信號的振幅；ψ為兩種信號之間相位之差，一般通過多路徑效應(yīng)從所獲得的信噪比(signal to noise ratio，SNR)觀測值中分離出來，并以此來得到SNR中因地表反射引起的多路徑變化信息。結(jié)合式(1)可知，多路徑效應(yīng)對信噪比的影響主要來源于衛(wèi)星高度角，上式中Ad的值通常遠(yuǎn)大于Ar的值，因此采用低階多項式擬合的方法消除趨勢項Ad。圖2為去除趨勢項后的SNR殘差序列。

圖2 去除趨勢項后的信噪比殘差序列

由式(1)可知，合成信號的振幅與直射信號幅度、反射信號幅度和相位延遲三者有關(guān)。其中直射信號與反射信號的相位差ψ如式(2)所示。

(2)

式中：h為反射面高度；λ為載波波長。h與頻率f之間的關(guān)系如式(3)所示。

(3)

由式(3)可知，反射信號頻率與反射面高度h呈線性關(guān)系。又因為sinθ的變化是非線性的，所以通過L-S(lomb-scargle)頻譜變換[12]的方法對頻率進行分析，最終得到反射面高度h，如圖3所示。

圖3 L-S頻譜分析圖

衛(wèi)星干涉信號的信噪比表達如式(4)所示[13]。

SNR=Acos(2πfsinθ+φ)

(4)

式中：A、f、φ分別為干涉信號參量中的振幅、頻率及相位。其中的sinθ可以在導(dǎo)出衛(wèi)星導(dǎo)航文件時一并導(dǎo)出衛(wèi)星高度角來計算具體數(shù)值，在得到去除趨勢項的信噪比之后，再將經(jīng)過L-S頻譜分析后得到的頻率f帶入式(4)，進行最小二乘擬合后獲得干涉信號SNR的振幅和相位。文獻[13]表明上述三個特征參量都可以用于土壤濕度的反演，且與土壤濕度之間呈線性關(guān)系。圖4為小波變換神經(jīng)網(wǎng)絡(luò)模型反演土壤濕度基本流程圖。

圖4 土壤濕度反演流程圖

1.2 小波變換分解理論

小波變換的實質(zhì)是把函數(shù)空間中的任意函數(shù)表示成其在具有不同伸縮因子a和平移因子b在小波基函數(shù)上的投影疊加，小波基函數(shù)是基于Mallat 算法[14]實現(xiàn)時間-頻率的分解與重構(gòu)。通過改變a和b的數(shù)值繼而得到具有不同時間-頻率寬度的小波。

考慮到初始信號的震蕩性與光滑度，本文采用正交對稱小波sym作為基函數(shù)。圖5表示了P036測站2017年DOY 123 PRN 32號衛(wèi)星的信噪比SNR趨勢分離結(jié)果，橫軸為等間隔觀測歷元，縱軸為信噪比指數(shù)變化，單位為dB-Hz，數(shù)據(jù)采樣率為15 s。由圖可見，在衛(wèi)星信號的低頻部分，多項式擬合表現(xiàn)較差，難以準(zhǔn)確地擬合SNR信號變化趨勢，而最終反演所需要的反射信號是通過將原始SNR減去擬合值來獲取，因此準(zhǔn)確擬合原始信噪比至關(guān)重要。小波變換方法能很好地表達該信號的變化趨勢，繼而實現(xiàn)信號高頻與低頻部分的分離，實現(xiàn)高質(zhì)量信噪比殘差序列的提取。

圖5 PRN 07號衛(wèi)星的信號信噪比擬合圖

1.3 BP神經(jīng)網(wǎng)絡(luò)模型

BP神經(jīng)網(wǎng)絡(luò)模型是在1986年由Rumelhart和McCelland 為首的科研團隊提出來的，其核心思想是通過負(fù)梯度下降算法將誤差控制在一個較小的范圍內(nèi)，再將誤差的變化量反向傳播到神經(jīng)網(wǎng)絡(luò)的每一層，進而調(diào)整每層的參數(shù)值，經(jīng)過多次迭代之后誤差就會穩(wěn)定在一個范圍內(nèi)[15]。

2 實例分析

2.1 數(shù)據(jù)獲取

為了驗證傳統(tǒng)相位線性回歸模型、BP神經(jīng)網(wǎng)絡(luò)模型及小波神經(jīng)網(wǎng)絡(luò)模型反演土壤濕度的可靠性與精度，本文下載了PBO中P036測站2017DOY121至2017DOY181該段時間的數(shù)據(jù)。該測站位于美國新墨西哥州天使火機場(36°25′13"N，105°17′37"W)附近，平均海拔2 529.9 m，建于2004年。采用的接收機類型為TRIMBLE NETRS，天線類型為TRM29659.00，可提供時間分辨率為15 s的L2C觀測值，周圍地貌較為平坦沒有大型的阻礙物，因此測站數(shù)據(jù)適合于土壤濕度反演。圖6為P036測站的周圍地表環(huán)境。

圖6 P036周圍地表環(huán)境圖

圖7為P036與SNOTEL中心提供的PALO氣象站點數(shù)據(jù)以及該測段時間范圍內(nèi)的降雨情況。該氣象站提供了地表深處5 cm左右的土壤濕度實測值。由圖7可知，PBO與SNOTEL兩種土壤濕度的數(shù)據(jù)整體變化趨勢基本保持一致，都與降雨事件的發(fā)生保持緊密的聯(lián)系，發(fā)生強降雨事件時土壤濕度數(shù)值都會隨之劇烈上升，但是從具體數(shù)值上來看存在一定的偏差，SNOTEL數(shù)值要明顯小于前者。產(chǎn)生該現(xiàn)象的原因可能是PBO用來反演土壤濕度的模型存在一定的誤差。

圖7 土壤濕度參考值與降雨值

2.2 實驗數(shù)據(jù)處理

當(dāng)衛(wèi)星高度角過小時，信號容易受到周邊建筑或樹木的影響導(dǎo)致干涉現(xiàn)象不明顯，因此本實驗截取了高度角在5°～30°之間的衛(wèi)星信噪比數(shù)據(jù)進行處理。利用TEQC軟件對觀測文件進行處理得到信噪比數(shù)據(jù)、衛(wèi)星高度角和方位角等一系列數(shù)據(jù)。為了去除以直射分量為主的趨勢項，在傳統(tǒng)的低階多項式擬合方法的基礎(chǔ)上，采用了小波變換方法來代替多項式擬合，最終得到實驗所需要的信噪比殘差值即反射信號分量。

將處理后得到的反射信號分量進行L-S頻譜分析變換得到頻率f，再對式(4)進行最小二乘擬合得到振幅A和相位φ，把得到的三個土壤濕度特征向量作為BP神經(jīng)網(wǎng)絡(luò)輸入數(shù)據(jù)，實際的土壤濕度值作為輸出數(shù)據(jù)進行模型構(gòu)建。進行上述處理后一共得到122組數(shù)據(jù)，經(jīng)過多次調(diào)整訓(xùn)練集與測試集之間的比例，最終確定當(dāng)兩者之間的比例為3∶1時反演效果最佳。通過隨機抽取的方法將數(shù)據(jù)劃分為91組訓(xùn)練集與31組測試集。由于傳統(tǒng)的三個特征向量中相位與土壤濕度之間的相關(guān)性最高，因此本文建立了相位與土壤濕度之間的線性相關(guān)模型，如式(5)所示。

y=1.107 1x-0.043 9

(5)

式中：y為土壤濕度值；x為信號相位值。

2.3 實驗結(jié)果分析

為了進一步探究三種模型反演能力的高低，利用P036測站DOY213至DOY274的L2C數(shù)據(jù)進行測試。圖8為三種不同模型反演結(jié)果與實測值分析圖，圖8(a)、圖8(b)、圖8(c)分別代表傳統(tǒng)相位線性回歸、BP神經(jīng)網(wǎng)絡(luò)、小波神經(jīng)網(wǎng)絡(luò)模型分析圖。由圖8可知，在土壤濕度值較大的時候，傳統(tǒng)相位線性回歸模型明顯與實測值有一定的偏移量，反演精度不如另外兩種模型；BP神經(jīng)網(wǎng)絡(luò)模型整體反演效果較為理想，但依然不如小波神經(jīng)網(wǎng)絡(luò)模型，造成這種結(jié)果的原因可能是實驗樣本量較小產(chǎn)生了過擬合現(xiàn)象導(dǎo)致的。

圖8 三種模型結(jié)果分析圖

從圖9能夠明顯地看出,BP神經(jīng)網(wǎng)絡(luò)模型與小波神經(jīng)網(wǎng)絡(luò)模型總體上都能很好地吻合土壤濕度真值的走勢。傳統(tǒng)相位線性回歸模型雖然整體趨勢較好，但是從個別天數(shù)來看，與前者存在較大的偏差，尤其是在DOY214、DOY218、DOY236、DOY246，主要是因為GPS不同衛(wèi)星之間存在軌道差異，導(dǎo)致接收機接收到的信號分別來自不同的地表環(huán)境，進而導(dǎo)致反演結(jié)果發(fā)生跳變現(xiàn)象，而且多半出現(xiàn)在降雨事件發(fā)生后一段時間土壤濕度開始下降過程中。在所選取的實驗天數(shù)內(nèi)發(fā)生了數(shù)次降雨事件，由于受降雨天氣影響，土壤濕度出現(xiàn)急劇上升。而且從圖中能夠明顯看出，在土壤濕度上升的過程中小波神經(jīng)網(wǎng)絡(luò)模型相比另外兩種模型的反演精度在一定程度上有較大的提升。

圖9 三種模型反演結(jié)果與真值對比圖

表1為三種模型反演結(jié)果精度對比。由表1可知，BP神經(jīng)網(wǎng)絡(luò)模型的均方根誤差RMSE與平均絕對偏差MAE分別為0.027和0.024，決定系數(shù)R2為0.865，精度比傳統(tǒng)相位線性回歸模型提高了15.33%；而小波神經(jīng)網(wǎng)絡(luò)模型的均方根誤差與平均絕對偏差分別為0.015和0.013，決定系數(shù)為0.943，精度比傳統(tǒng)相位線性回歸模型提升了25.73%，比單獨的BP神經(jīng)網(wǎng)絡(luò)模型又提升了9.02%。結(jié)果證明，小波神經(jīng)網(wǎng)絡(luò)模型相對前兩者而言，在誤差方面有了明顯的改善，說明改進后結(jié)合了小波變換的神經(jīng)網(wǎng)絡(luò)模型與實測土壤濕度的吻合度得到了進一步的提升。圖10為三種模型與土壤濕度真值的絕對誤差對比圖，也直接證明了小波神經(jīng)網(wǎng)絡(luò)模型的誤差的確要明顯低于前兩者。傳統(tǒng)相位線性回歸模型的誤差控制在±0.08 cm3/cm-3以內(nèi)，BP神經(jīng)網(wǎng)絡(luò)模型的誤差控制在±0.05 cm3/cm-3，而小波神經(jīng)網(wǎng)絡(luò)模型僅僅在±0.03 cm3/cm-3以內(nèi)，相比之下，在精度上得到較大提升。

表1 三種模型反演精度對比

圖10 三種模型誤差對比圖

3 結(jié)束語

針對GNSS-IR反演土壤濕度研究出現(xiàn)的一些問題，如反射信號分離過程不理想、衛(wèi)星信號不穩(wěn)定以及地表粗糙度和植被系數(shù)等帶來的影響，本文構(gòu)建了一種將小波變換分解與BP神經(jīng)網(wǎng)絡(luò)相結(jié)合的新模型，得到以下結(jié)論。

1)土壤濕度與特征向量之間的確存在著一定的線性關(guān)系，其中反演結(jié)果較為理想的相位模型決定系數(shù)達到了0.75，反演誤差控制在±0.08 cm3/cm-3以內(nèi)。

2)BP神經(jīng)網(wǎng)絡(luò)模型總體反演結(jié)果有一定的改善，尤其是在降雨事件發(fā)生后土壤濕度急劇上升的過程中，能夠有效抑制單顆衛(wèi)星反演容易出現(xiàn)的異常跳變現(xiàn)象，而且也有效地抑制了地表粗糙度及植被覆蓋等因素造成的影響。

3)BP神經(jīng)網(wǎng)絡(luò)模型與小波神經(jīng)網(wǎng)絡(luò)模型相對傳統(tǒng)相位線性回歸模型精度分別提高了15.33%、25.73%，誤差同時得到了有效的改正。這說明兩種模型相對于傳統(tǒng)相位線性回歸模型都有很大程度的提升，均能夠有效抑制地表粗糙度和植被覆蓋的影響，而神經(jīng)網(wǎng)絡(luò)模型較后者精度較低的原因是因為在衛(wèi)星反射信號分離的過程中采取的方法存在較大偏差，也從側(cè)面證明了小波變換的頻率局部優(yōu)化分析法對于提取反射信號分量有著較大的優(yōu)勢。

后續(xù)將在本文基礎(chǔ)上引入其他算法與現(xiàn)有模型相結(jié)合，以及聯(lián)合多模多頻GNSS數(shù)據(jù)反演土壤濕度，旨在進一步提高GNSS-IR反演土壤濕度的精度。

致謝：非常感謝UNAVCO提供的GPS衛(wèi)星觀測數(shù)據(jù)以及ISMN項目組提供的土壤濕度參考值。