孫波， 梁勇,*， 漢牟田， 楊磊， 荊麗麗， 俞永慶

(1. 山東農(nóng)業(yè)大學(xué)信息科學(xué)與工程學(xué)院, 泰安 271019； 2. 北京航空航天大學(xué)電子信息工程學(xué)院, 北京 100083； 3. 中國石油化工股份有限公司勝利油田分公司海洋采油廠, 東營 257237)

土壤濕度是全球碳水循環(huán)的重要環(huán)節(jié)，是量化陸地及大氣能量交換的關(guān)鍵參數(shù)[1]。及時(shí)、準(zhǔn)確地獲取農(nóng)田土壤濕度數(shù)據(jù)對于在農(nóng)業(yè)生產(chǎn)中進(jìn)行合理灌溉、減少水資源浪費(fèi)、降低生產(chǎn)成本和提高農(nóng)作物產(chǎn)量都是非常重要的。

相比于傳統(tǒng)的烘干稱重法以及時(shí)域反射法(TDR)、頻域反射法(FDR)等接觸方式獲取土壤濕度的方法[2]，利用全球衛(wèi)星導(dǎo)航系統(tǒng)反射測量法(Global Navigation Satellite System Reflectometry, GNSS-R)技術(shù)進(jìn)行土壤濕度探測是一種新興的技術(shù)手段，其具有非接觸、大面積、實(shí)時(shí)性、連續(xù)性的優(yōu)點(diǎn)，近年來受到越來越多的關(guān)注。

1993年，歐洲航天局(ESA)的Martin-Neira[3]首次提出了利用GNSS-R信號可以進(jìn)行海面測高的技術(shù)。隨后，科羅拉多大學(xué)與美國國家航空航天局(NASA)合作研制了GNSS-R接收處理設(shè)備，并開展了一系列機(jī)載實(shí)驗(yàn)，探索該技術(shù)用于土壤濕度探測的可行性。近年來，Rodriguez-Alvarez等[4]在GNSS-R技術(shù)的基礎(chǔ)上提出可以利用衛(wèi)星直射信號和多徑信號的干涉現(xiàn)象進(jìn)行土壤濕度探測的方法，稱為干涉圖技術(shù)(Interference Pattern Technique，IPT)，也稱為全球衛(wèi)星導(dǎo)航系統(tǒng)干涉測量法(Global Navigation Satellite System Interferometry and Reflectometry，GNSS-IR)。Larson等[5]指出GNSS接收機(jī)記錄的信噪比(Signal to Noise Ratio，SNR)數(shù)據(jù)受多徑信號的影響呈現(xiàn)出衰減振蕩特征，該特征是干涉現(xiàn)象的表現(xiàn)，其通過實(shí)驗(yàn)證明信噪比數(shù)據(jù)的振蕩振幅和相位與土壤濕度有關(guān)。Chew等[6]在此基礎(chǔ)上建立了信噪比相位和土壤濕度的經(jīng)驗(yàn)?zāi)Ｐ?，可以用來反演裸土表層的土壤濕度?018年，李飛等[7]對鏡面反射點(diǎn)、空間分辨率和土壤探測深度進(jìn)行了理論分析和公式推導(dǎo)。2018年，嚴(yán)頌華等[8]探討了基于信號強(qiáng)度指示器(SSI)數(shù)據(jù)的土壤水分估算的可行性，并進(jìn)行了SSI相位估算和原位土壤水分?jǐn)?shù)據(jù)記錄的實(shí)驗(yàn)比較。在植被影響研究方面，張思勃等[9]考慮小麥不同生長階段的影響，建立了土壤濕度與信噪比振蕩頻率和相位的經(jīng)驗(yàn)?zāi)Ｐ汀?/p>

然而，上述經(jīng)驗(yàn)反演模型的建立均假設(shè)土壤濕度與信噪比幅度、相位等觀測量呈線性關(guān)系，因此模型建立過程多采用傳統(tǒng)的一元統(tǒng)計(jì)回歸方法，其反演結(jié)果的相關(guān)系數(shù)介于0.76～0.90之間[10-12]，反演精度不高。而機(jī)器學(xué)習(xí)等人工智能算法的興起，使得復(fù)雜回歸模型的建立成為可能。楊磊等[13]探索了機(jī)器學(xué)習(xí)方法在GNSS-R土壤濕度技術(shù)中的應(yīng)用，取得了較好的結(jié)果。機(jī)器學(xué)習(xí)方法具有黑箱特征，不用指定回歸函數(shù)的形式，可以進(jìn)行任意的線性或非線性回歸，并且隱含實(shí)現(xiàn)了特征提取的環(huán)節(jié)，可以對土壤粗糙度等引起的噪聲進(jìn)行抑制，在不出現(xiàn)過擬合的情況下，其效果是傳統(tǒng)的統(tǒng)計(jì)回歸方法不可及的。

本文提出使用支持向量機(jī)(Support Vector Machine，SVM)這一機(jī)器學(xué)習(xí)方法，進(jìn)行GNSS-IR土壤濕度反演模型建立，特別地，采用遺傳算法(Genetic Algorithm，GA)解決SVM的參數(shù)自動優(yōu)選問題，并通過實(shí)測數(shù)據(jù)處理，證明該方法可取得較好的結(jié)果。進(jìn)一步與統(tǒng)計(jì)回歸及其他機(jī)器學(xué)習(xí)方法進(jìn)行比較，證明了GA-SVM方法在同等條件下效果更優(yōu)。

1 GNSS-IR土壤濕度反演原理

相比于傳統(tǒng)的GNSS-R反演模式，GNSS-IR模式只采用一根右旋圓極化(Right Hand Circular Polarization，RHCP)天線同時(shí)接收GNSS的直反射衛(wèi)星信號，GNSS發(fā)射的信號為RHCP信號，根據(jù)電磁波反射理論，當(dāng)?shù)孛鏋椴缓判猿煞值耐寥澜橘|(zhì)時(shí)，在衛(wèi)星高度角θ較低的情況下，反射信號極性仍以RHCP為主[14]，其場景如圖1所示。

圖1 干涉場景Fig.1 Scenario of interference

在傳統(tǒng)的衛(wèi)星信號接收中，通常把多徑反射信號當(dāng)作干擾進(jìn)行抑制，但實(shí)際上反射信號與反射面的物理特征如土壤濕度具有相關(guān)性。由于地基情況下直射跟反射信號頻率近似相同，只是傳播路徑長度不同，因此2個(gè)信號可在接收天線處發(fā)生較為穩(wěn)定的干涉現(xiàn)象，形成干涉信號,如圖2(a)所示?？梢钥闯觯诘透叨冉菚r(shí)，干涉波動現(xiàn)象十分明顯。該信號經(jīng)接收機(jī)接收處理后以信噪比數(shù)據(jù)的形式呈現(xiàn)出來。信噪比的數(shù)學(xué)模型[15]為

(1)

式中：Ad、Am分別為直射、反射信號的幅度；ψ為直射信號和反射信號的相位差。

則有

(2)

式中：φ1為直反射路徑差導(dǎo)致的相位差；φ0為干涉初始相位；H為接收機(jī)天線的等效高度，其隨著電磁波的穿透深度和反射面的介電特性的變化而變化；λ為GNSS衛(wèi)星信號的波長。

通過式(2)，可以得到多徑振蕩的頻率為

(3)

考慮在每天信噪比觀測的幾小時(shí)時(shí)間內(nèi)，天線等效高度的變化率可以忽略不計(jì)，并進(jìn)一步令x=sinθ，式(3)可簡化為

(4)

式(4)說明H與f成線性關(guān)系，天線越高，振蕩頻率越高。結(jié)合式(1)～式(4)，將直射信號剔除，只保留與反射面參數(shù)有關(guān)的多徑信號，如圖2(a)、(b)所示。該多徑信號SNRm可以表示為

SNRm=Amcos(4πHλ-1sinθ+φ)

(5)

式中：φ為多徑信號的初始相位。

通過對SNRm信號進(jìn)行譜分析，可以獲得多徑信號的振蕩頻率f和天線等效高度H，如圖2(c)所示。由于SNRm信號是非均勻采樣信號，常規(guī)的傅里葉分析并不適用，引入Lomb-Scargle方法[16]對SNRm進(jìn)行譜分析， Lomb-Scargle方法能夠提供一種針對不均勻采樣時(shí)間序列的譜幅值計(jì)算方法，可以處理非均勻采樣的信號。

圖2 上升段信噪比分析Fig.2 Ascending SNR analysis

此后通過最小二乘法對SNRm進(jìn)行擬合，可以得到多徑信號的幅度Am和初始相位φ，如圖2(b)所示，并分別建立振蕩頻率、幅度、相位和土壤濕度之間的經(jīng)驗(yàn)?zāi)Ｐ汀?/p>

2 GA-SVM土壤濕度反演模型

文獻(xiàn)[10-12]提出建立一元線性回歸模型進(jìn)行土壤濕度的反演，然而考慮到振蕩頻率、幅度、相位與土壤濕度的相關(guān)性，采用多元回歸模型可以有效利用特征信息，充分利用它們之間的相關(guān)性，提高擬合精度，本文采用非線性SVM模型進(jìn)行土壤濕度的反演研究。

SVM是20世紀(jì)90年代由Vapnik[17]提出來的一種研究小樣本、小概率事件的神經(jīng)網(wǎng)絡(luò)模型,可以根據(jù)有限的樣本信息在模型復(fù)雜性和學(xué)習(xí)能力之間尋找最優(yōu)解, 以獲得最好的預(yù)測精度, 其在解決小樣本、非線性和高維學(xué)習(xí)中表現(xiàn)出許多特有的優(yōu)勢，而土壤濕度反演本質(zhì)上是一個(gè)回歸問題，適合選用SVM算法進(jìn)行回歸分析。

SVM利用核函數(shù)映射的方法，通過非線性映 射將已知空間變換到高維空間，從而利用解線性問題的方法求解非線性問題。本文中采用徑向基(Radial Basis Function,RBF)核函數(shù)來建立非線性 SVM模型。RBF核函數(shù)相比于線性核能夠處理分類標(biāo)注和屬性的非線性關(guān)系，相比于多項(xiàng)式核有更少的參數(shù)，同時(shí)具有簡單實(shí)用、普適性好的優(yōu)點(diǎn)[18]。并且根據(jù)實(shí)踐證明，SVM的性能與核函數(shù)的參數(shù)g以及懲罰因子C有很大關(guān)系，人為規(guī)定和調(diào)整參數(shù)的話，既耗費(fèi)時(shí)間又不能確保預(yù)測結(jié)果的準(zhǔn)確性，因此，本文引入了GA對SVM回歸模型的參數(shù)進(jìn)行自動尋找最優(yōu)解。

遺傳算法是 Holland[19]于1975年提出的一種基于生物遺傳和進(jìn)化機(jī)制的、適用于復(fù)雜系統(tǒng)優(yōu)化的自適應(yīng)概率優(yōu)化技術(shù)。該算法是借鑒進(jìn)化生物學(xué)中遺傳、突變、自然選擇、雜交等現(xiàn)象而探索出的一種搜索最優(yōu)解方法，主要包括編碼、初始種群、計(jì)算適應(yīng)度以及選擇、交叉、變異3個(gè)基本運(yùn)算。與傳統(tǒng)優(yōu)化算法相比，GA尋優(yōu)具有以下特點(diǎn)：①搜索效率高，以群體為單位，對自變量中所有個(gè)體進(jìn)行并行搜索；②可進(jìn)行全局最優(yōu)解搜索，避免陷入局部最優(yōu)解；③不受目標(biāo)函數(shù)影響，適應(yīng)性強(qiáng)。

基于GA-SVM的GNSS-IR土壤濕度反演模型建立過程如圖3所示。

圖3 基于GA-SVM的GNSS-IR土壤濕度反演模型Fig.3 GNSS-IR soil moisture inversion model based on GA-SVM

3 實(shí)驗(yàn)與結(jié)果分析

3.1 實(shí)驗(yàn)數(shù)據(jù)采集

實(shí)驗(yàn)場地位于法國圖盧茲(Toulouse)市拉馬斯克雷(Lamasquere)鎮(zhèn)的一處農(nóng)場(北緯 43°29′14.45″，東經(jīng)1°13′44.11″)，實(shí)驗(yàn)田種植作物為大豆。2014年2月3日—2014年3月21日進(jìn)行了為期47 d的土壤濕度探測實(shí)驗(yàn)。實(shí)驗(yàn)期間該地區(qū)正處于冬季，日平均降雨量僅為1.19 mm，地表植被稀少，可以視為裸土，且場地四周無遮蔽。土壤表面均方根高度約為0.02 m，反射區(qū)域2～5 cm深度的土壤成分為:沙土18%，粉土41%，黏土41%。使用Leica GR25接收機(jī)和AR10基準(zhǔn)站天線，同時(shí)接收GPS的直射、反射信號，天線的架設(shè)高度為1.69～1.70 m，相位中心指向天頂。距離天線相位中心在地表投影約2 m的位置采用填埋方式布置2枚ML3 Theta Probe土壤濕度傳感器采集土壤濕度同比數(shù)據(jù)，2枚傳感器深度分別為2 cm和5 cm，精度為±1%，取樣間隔時(shí)間為2 min。本文選取干涉現(xiàn)象明顯的2°～30°低高度角的L1波段S1C SNR信號進(jìn)行土壤濕度反演。

3.2 反演模型建立

3.2.1 樣本數(shù)據(jù)的選取

選取觀測效果較好的GPS PRN 12的SNR數(shù)據(jù)按照第1節(jié)的原理進(jìn)行處理，得到47組與土壤濕度具有相關(guān)性的多徑信號的初始相位φ、振蕩頻率f、幅度Am，結(jié)合衛(wèi)星高度角θ，與實(shí)測的2 cm與5 cm深度的土壤濕度日平均值形成SVM數(shù)據(jù)集，這47組數(shù)據(jù)對應(yīng)的傳感器土壤濕度實(shí)測平均值范圍為23.24 ～29.45 cm3/cm3。為了消除各影響因素因類型及量綱不同帶來的影響，對輸入變量均做歸一化處理。將這47組數(shù)據(jù)按照約3∶1的比例劃分成訓(xùn)練集與測試集，訓(xùn)練集與測試集互斥，沒有交集，將傳感器測量的日平均土壤濕度作為SVM的輸出。

3.2.2 GA-SVM模型訓(xùn)練

本文SVM模型選用普適度較好的RBF核函數(shù)，并采用GA對其進(jìn)行參數(shù)尋優(yōu)。為了提高反演模型的精度及泛化能力，對核函數(shù)參數(shù)g和懲罰因子C進(jìn)行優(yōu)選是關(guān)鍵。GA的相關(guān)設(shè)置如下：最大進(jìn)化代數(shù)為200，種群數(shù)量為20，懲罰因子C的取值范圍為[0,100]，核函數(shù)參數(shù)g的取值范圍為[0,1 000]，交叉驗(yàn)證次數(shù)V=5。尋優(yōu)過程的適應(yīng)度曲線如圖4所示?？梢钥闯?，經(jīng)過30次進(jìn)化迭代后，種群最佳適應(yīng)度基本保持穩(wěn) 定，可以認(rèn)為近似達(dá)到最優(yōu)解。最終尋優(yōu)結(jié)果為:C=3.23，g=0.08，適應(yīng)度MSE=0.045。

圖4 GA-SVM進(jìn)化代數(shù)曲線Fig.4 Evolution algebra curves of GA-SVM

3.3 結(jié)果分析

將最佳參數(shù)C和g代入SVM模型，選取GPS PRN 12衛(wèi)星的前35 d數(shù)據(jù)和土壤濕度真值作為訓(xùn)練集，對后12 d的土壤濕度進(jìn)行預(yù)測反演，所得結(jié)果如圖5所示。

圖5(a)～(c)顯示訓(xùn)練樣本預(yù)測值與土壤濕度實(shí)測值具有很好的相關(guān)性,誤差控制在±0.5 cm3/cm3之內(nèi)，決定系數(shù)R2達(dá)到了0.882 7。圖5(d)～(f)表明利用訓(xùn)練后的GA-SVM模型可以對測試樣本進(jìn)行較好的反演，反演結(jié)果接近實(shí)測值，反演誤差控制在±0.3 cm3/cm3之內(nèi)，沒有出現(xiàn)偏態(tài)分布，反演值與實(shí)測值的線性回歸方程決定系數(shù)R2達(dá)到了0.956 9。在圖5(d)～(f)中，還進(jìn)行了GA-SVM方法與傳統(tǒng)的利用SNR頻率、幅度、相位的一元線性回歸方法的對比，結(jié)果表明，本文方法的反演結(jié)果更接近實(shí)測值，誤差較小，相比于單變量反演效果最好的相位反演(R2=0.874 7),決定系數(shù)R2提高了約9.40%，說明GA-SVM優(yōu)化模型反演精度高，泛化能力較強(qiáng)。

為了進(jìn)一步驗(yàn)證GA-SVM模型的優(yōu)越性，在使用同一數(shù)據(jù)樣本的條件下，本文還與粒子群優(yōu)化的SVM模型(PSO-SVM)、反向傳播(BP)神經(jīng)網(wǎng)絡(luò)這2種機(jī)器學(xué)習(xí)方法進(jìn)行了對比，表1給出了對比結(jié)果。

圖5 GA-SVM土壤濕度反演模型結(jié)果分析Fig.5 Result analysis of GA-SVM soil moisture inversion model

日期實(shí)測值/(cm3·cm-3)GA-SVMPSO-SVMBP神經(jīng)網(wǎng)絡(luò)反演值/(cm3·cm-3)絕對誤差/(cm3·cm-3)相對誤差/%反演值/(cm3·cm-3)絕對誤差/(cm3·cm-3)相對誤差/%反演值/(cm3·cm-3)絕對誤差/(cm3·cm-3)相對誤差/%2014-03-1025.8325.76-0.070.2726.140.311.2026.520.692.672014-03-1125.5525.41-0.140.5526.180.632.4726.951.405.482014-03-1225.2125.10-0.110.4425.300.090.3627.091.887.462014-03-1324.5524.850.301.2225.721.174.7724.950.401.632014-03-1424.4524.640.190.7825.380.933.8024.990.542.212014-03-1524.2724.480.210.8724.870.602.4724.440.170.702014-03-1624.0724.240.170.7124.240.170.7125.231.164.822014-03-1723.9724.180.210.8823.990.020.0825.471.506.262014-03-1823.8323.990.160.6723.940.110.4624.100.271.132014-03-1923.8823.71-0.170.7124.620.743.1024.660.783.272014-03-2023.6223.58-0.040.1724.060.441.8625.551.938.172014-03-2123.2423.490.251.0823.620.381.6424.361.124.82

表2 土壤濕度反演結(jié)果評價(jià)比較

4 結(jié) 論

本文提出了基于GA-SVM優(yōu)化模型的GNSS-IR土壤濕度反演方法，給出了數(shù)據(jù)處理的相關(guān)流程，并通過實(shí)測數(shù)據(jù)進(jìn)行了理論驗(yàn)證，結(jié)果表明：

1) 在土壤濕度介于23.24～29.45 cm3/cm3范圍內(nèi)，反演值與高精度土壤濕度計(jì)測得的實(shí)測值的相關(guān)度較好，線性回歸方程決定系數(shù)達(dá)到了0.956 9，反演誤差控制在±0.3 cm3/cm3之內(nèi)。

2) 與傳統(tǒng)利用相位等進(jìn)行一元線性回歸的方法相比，決定系數(shù)提高了約9.40%，進(jìn)一步驗(yàn)證了GA-SVM優(yōu)化模型的有效性。

3) 在相同數(shù)據(jù)集下，與PSO-SVM、BP神經(jīng)網(wǎng)絡(luò)的對比結(jié)果表明，本文模型各項(xiàng)誤差較小，說明其更適用于小樣本條件下GNSS單天線的土壤濕度反演，具有較好的普適性。

本文為土壤濕度的遙感監(jiān)測提供了一種新的思路和方法。由于季節(jié)所限，本次實(shí)驗(yàn)及數(shù)據(jù)分析都是基于裸土的土壤濕度分析，地表植被和農(nóng)作物對土壤濕度的影響并未進(jìn)行修正，這也是后續(xù)實(shí)驗(yàn)研究的主要內(nèi)容。

致謝感謝法國CESBIO實(shí)驗(yàn)室的Baup和Boniface博士收集的氣象數(shù)據(jù)，感謝法國圖盧茲第三大學(xué)的Roussel和Frappart博士提供的GNSS觀測數(shù)據(jù),同時(shí)感謝金華職業(yè)技術(shù)學(xué)院阮宏梁老師在實(shí)驗(yàn)中提供的支持。