游高沖,郭 杭,羅孝文,尹海博,王朝陽
(1.南昌大學(xué)信息工程學(xué)院,江西 南昌 330031; 2.自然資源部第二海洋研究所,浙江 杭州 310012)
監(jiān)測海面高度變化對研究全球氣候變化、構(gòu)建海洋高程基準(zhǔn)、預(yù)防沿海臺風(fēng)等自然災(zāi)害具有重要意義。全球?qū)Ш蕉ㄎ幌到y(tǒng)多徑反射(GNSS Multipath Reflectometry,GNSS-MR)技術(shù)是基于全球定位系統(tǒng)多路經(jīng)效應(yīng)發(fā)展起來的一種新型地基遙感技術(shù),利用地表反射信號與導(dǎo)航衛(wèi)星直射信號的相干信號獲取地表環(huán)境信息[1]。相比于全球定位系統(tǒng)干涉(GNSS Reflectometry,GNSS-R)遙感技術(shù)[2],GNSS-MR不需要另外架設(shè)左旋圓極化天線[3],具有設(shè)備要求低、信號源穩(wěn)定、全天候、覆蓋范圍廣等優(yōu)勢,被廣泛應(yīng)用于地表土壤水分監(jiān)測[4]、積雪監(jiān)測[5]和海洋潮汐監(jiān)測[6-7]等領(lǐng)域中。
LARSONKM et al[8-11]提出基于載波SNR(Signal-to-Noise Ratio)值的GNSS-MR技術(shù),利用該技術(shù)對土壤水分和雪深進(jìn)行了反演研究,并推廣到潮位反演。由于GNSS-MR僅使用低高度角信噪比數(shù)據(jù),在進(jìn)行潮位監(jiān)測時(shí)常存在時(shí)間分辨率和反演精度低的問題,為此,國內(nèi)外學(xué)者開展了一系列的改進(jìn)研究。針對GNSS-MR中低階多項(xiàng)式提取反射信號存在信號混雜的問題,王瑞芳[12]、ZHANG et al[13]等采用經(jīng)驗(yàn)?zāi)B(tài)分解提取反射信號的方法,提純了海面反射信號,提高了潮位反演精度。為提高GNSS-MR潮位監(jiān)測的時(shí)間分辨率,王杰 等[14]采用小波分析的方法,提取反射信號瞬時(shí)頻率進(jìn)行潮位反演,有效地增加了潮位反演點(diǎn)數(shù),但該方法提高潮位監(jiān)測時(shí)間分辨率有限,仍難以滿足高精度、高時(shí)間分辨率潮位監(jiān)測的要求。陳發(fā)德 等[15]利用GPS、BDS、Galileo多系統(tǒng)數(shù)據(jù)進(jìn)行融合潮位反演,有效提高了潮位反演的連續(xù)性和時(shí)間分辨率,但該方法將不同系統(tǒng)的反演結(jié)果無差別地組合在一起作為反演結(jié)果,聯(lián)合反演潮位的均方根誤差(Root Mean Squared Error,RMSE)為0.323 m,反演精度偏低。何秀鳳 等[16]在進(jìn)行風(fēng)暴潮反演時(shí),采用滑動窗口最小二乘法,對多系統(tǒng)反演數(shù)據(jù)進(jìn)行融合,在提高了風(fēng)暴潮監(jiān)測的時(shí)間分辨率的同時(shí),也提高了監(jiān)測的精度,但反演精度仍為分米級。
為進(jìn)一步提高GNSS-MR潮位反演的時(shí)間分辨率和精度,本文提出基于LS-SVM的多系統(tǒng)數(shù)據(jù)融合潮位反演方法,利用多個(gè)系統(tǒng)的潮位反演數(shù)據(jù),對LS-SVM進(jìn)行訓(xùn)練,得到多系統(tǒng)數(shù)據(jù)融合的LS-SVM潮位反演模型,進(jìn)而估計(jì)出融合后的潮位值,實(shí)現(xiàn)高時(shí)間分辨率、高精度潮位監(jiān)測。
圖1為利用GNSS-MR技術(shù)反演潮位變化的示意圖。只考慮1次反射的情況下,載波多路徑相位延遲φ可以表示成:
(1)
式中:λ表示載波波長,h表示天線相位中心到海面的高度,θ表示衛(wèi)星高度角。受多路徑效應(yīng)的影響,海面高度信息以直、反射相干涉的方式保存在混合信號的信噪比中。海面多路徑效應(yīng)下SNR滿足如下等式:
(2)
式中:Ad是直射信號幅值,Ar是反射信號幅值,φ是多路徑效應(yīng)相位延遲量。SNR可視為直射信號作用的趨勢項(xiàng)與反射信號作用的周期項(xiàng)的疊加,圖2給出了GPS 01號衛(wèi)星的SNR序列和去趨勢后的殘差序列SNRm。由圖2a可知,SNR整體呈拋物線形式,通過二次擬合去趨勢,即可得到反射信號作用的周期項(xiàng)(圖2b),周期項(xiàng)SNRm滿足如下等式:
SNRm=Acos(4πhλ-1sinθ)
(3)
式中:A是信號振幅。記f=2hλ-1,t=sinθ,則上式可進(jìn)一步整理為如下等式:
SNRm=Acos(2πft)
(4)
周期信號SNRm的頻率包含著與潮位反演相關(guān)的h值,SNRm是以高度角正弦值為變量的不等間距信號,無法通過傅里葉變換獲取。因此,本文采用Lomb-Scargle(L-S)頻譜分析[17],提取SNRm的頻率f,再由h=f·λ/2、hs=H-h,求出海面高度值hs。圖3給出了2020年214年積日部分GPS衛(wèi)星反射信號SNRm及相應(yīng)的L-S頻譜圖,L-S頻譜圖(圖3b)中曲線峰值所對應(yīng)的橫坐標(biāo)為有效反射高h(yuǎn)。
圖1 GNSS-MR監(jiān)測潮位變化示意圖
圖2 GPS 01號衛(wèi)星SNR序列(a)和二次擬合去趨勢后的殘差SNRm序列(b)Fig.2 SNR sequence diagram (a) and residual SNRm sequence chart after two fitting detrended trend (b) of GPS 01 satellite
圖3 HKQT站GPS部分衛(wèi)星海面反射信號SNRm(a)和相應(yīng)的Lomb-Scargle頻譜(b)Fig.3 Partial satellite sea surface reflection signal SNRm (a) and corresponding Lomb-Scargle spectrum (b) of GPS at HKQT station
根據(jù)SVM思想,在高維空間的回歸方程為
f(X)=〈ω,φ(x)〉+b
(5)
式中:〈·,·〉表示點(diǎn)積;ω∈Rnk,是原權(quán)重空間中的權(quán)重向量;φ(·):Rn→Rnk,是將輸入數(shù)據(jù)映射到高維特征空間的非線性函數(shù);b是偏置量,b∈R。
對于LS-SVM可將優(yōu)化問題轉(zhuǎn)化為[18]
(6)
式中:ek為樣本誤差;γ表示正則化參數(shù),在噪聲較小的情況下,較小的γ可以避免過度擬合。此時(shí)的約束條件為
yk=〈ω,φ(xk)〉+b+ek,k=1,…,N
(7)
拉格朗日函數(shù)為
(8)
式中:αk為拉格朗日乘子,αk∈R。根據(jù)KKT條件消去ω和ek,得到如下線性方程式:
(9)
式中:y=[y1,y2…,yN]T;1v=[1,1,…,1]T;α=[α1,α2,…,αN]T;Ωk,l=〈φ(xk),φ(xl)〉=K(xk,xl),k,l=1,…,N。
利用最小二乘法對公式進(jìn)行求解,即可得到LS-SVM 非線性回歸模型:
(10)
基于LS-SVM多系統(tǒng)數(shù)據(jù)融合方法的潮位反演性能與核函數(shù)K、核參數(shù)σ和正則化參數(shù)γ的選取密切相關(guān)[18],本文選用能很好反映模型復(fù)雜程度和普適性較好的徑向基函數(shù)為核函數(shù),通過網(wǎng)格搜索法優(yōu)選σ和γ參數(shù)?;贚S-SVM的多系統(tǒng)融合模型潮位反演流程如圖4所示。
圖4 基于LS-SVM的多系統(tǒng)融合潮位反演流程圖Fig.4 Flow chart of multi system fusion tidal level inversion based on LS-SVM
為驗(yàn)證基于LS-SVM的多系統(tǒng)融合潮位監(jiān)測方法,本文選取2020年香港北岸HKQT站(22.29103351°N,114.21322098°E)214至244年積日的衛(wèi)星觀測數(shù)據(jù)進(jìn)行潮位反演實(shí)驗(yàn)。該測站配備有TrimbleNETR9型號接收機(jī),可以接收GPS、BDS、Galileo、GLONASS四個(gè)系統(tǒng)的衛(wèi)星信號。此外,該站同一位置建設(shè)有驗(yàn)潮站,提供1次/6 min的潮位數(shù)據(jù),可以用于對LS-SVM多系統(tǒng)融合潮位監(jiān)測方法的精度進(jìn)行評估。
潮位數(shù)據(jù)由國際海平面實(shí)施網(wǎng)( http://www.ioc-sealevelmonitoring.org/)提供,該潮位數(shù)據(jù)以平均海平面(Mean Sea Level, MSL)為基準(zhǔn)。為使反演潮位值與驗(yàn)潮站的潮位數(shù)據(jù)具有統(tǒng)一的基準(zhǔn),本文依據(jù)該網(wǎng)站提供的驗(yàn)潮站修訂本地參考圖(Revised Local Reference, RLR)進(jìn)行相應(yīng)的高程轉(zhuǎn)換,可得到HKQT站天線相位中心到MSL的高。圖5中RLR的橢球高為-8.274±0.004 m,經(jīng)精密單點(diǎn)定位解算出天線相位中心的橢球高為5.683 m,將天線橢球高減去RLR橢球高,再減去7.01 m即可得到天線相位中心到MSL的高度值,經(jīng)解算高度值為6.407 m。
潮位反演所需的載波波長、高度角、方向角、天線高等基本參數(shù)設(shè)定如表1所示。
圖5 驗(yàn)潮站修訂本地參考圖Fig.5 Revising local reference diagram of tide gauge station
表1 GNSS-MR潮位反演參數(shù)表Tab.1 Inversion parameters for GNSS-MR tide level
利用香港HKQT站提供的2020年214年積日至2020年244年積日連續(xù)30 d的GPS、BDS、Galileo、GLONASS衛(wèi)星數(shù)據(jù)對HKQT站附近潮位進(jìn)行反演,結(jié)果如圖6所示。圖6a、6c、6e和6g分別為GPS、BDS、Galileo和GLONASS反演結(jié)果及驗(yàn)潮站實(shí)測的潮位變化曲線,從中可以看出,HKQT站附近潮位在 ±1 m 區(qū)間范圍內(nèi)周期性波動,各系統(tǒng)反演點(diǎn)均落在曲線上或曲線附近,反演點(diǎn)的變化趨勢同驗(yàn)潮站實(shí)測潮位變化趨勢相同,反演潮位值和實(shí)測值具有明顯的一致性,表明4個(gè)系統(tǒng)均能有效實(shí)現(xiàn)潮位監(jiān)測。將反演潮位值減去同時(shí)刻驗(yàn)潮站實(shí)測潮位值,即可得到反演潮位的偏差值,圖6b、6d、6f和6h分別GPS、BDS、Galileo和GLONASS系統(tǒng)反演潮位的偏差序列,各系統(tǒng)反演偏差均在±0.5 m范圍內(nèi)變化。本文采用均方根誤差、相關(guān)系數(shù)和時(shí)間分辨率作為潮位反演的性能評定指標(biāo),評定結(jié)果如表2所示。
表2 單系統(tǒng)潮位反演性能對比Tab.2 Comparison of inversion performance of single systems tidal level
由表2可知,4個(gè)系統(tǒng)均能用于HKQT站附近的潮位監(jiān)測,且各系統(tǒng)反演潮位的精度相當(dāng),反演精度在0.120~0.139 m之間,相關(guān)性系數(shù)在0.953 6~0.963 2 之間,反演潮位值同驗(yàn)潮站實(shí)測值之間具有極強(qiáng)的相關(guān)性。各系統(tǒng)反演潮位的時(shí)間分辨率差異較大,其中GLONASS系統(tǒng)潮位監(jiān)測時(shí)間分辨率最高,平均1.23 h一個(gè)反演值;BDS系統(tǒng)潮位監(jiān)測時(shí)間分辨率最低,平均5.56 h一個(gè)反演值。總體而言,單系統(tǒng)GNSS-MR潮位反演精度和時(shí)間分辨率均較差,難以滿足高精度、高時(shí)間分辨率潮位反演性能需求。
為提高潮位監(jiān)測的時(shí)間分辨率和精度,考慮到不同系統(tǒng)之間的互補(bǔ)性,本文采用LS-SVM對GPS、BDS、Galileo、GLONASS反演潮位數(shù)據(jù)進(jìn)行融合。首先,將4個(gè)系統(tǒng)的反演潮位值組合在一起,作為 LS-SVM 的訓(xùn)練集;然后對訓(xùn)練集進(jìn)行LS-SVM訓(xùn)練,并采用十倍交叉驗(yàn)證方法分析確定正則優(yōu)化參數(shù)γ和核參數(shù)σ,得到潮位LS-SVM模型;最后,由LS-SVM潮位反演模型估計(jì)出潮位值。
圖6 GPS、BDS、Galileo、GLONASS反演潮位結(jié)果、反演偏差及驗(yàn)潮站實(shí)測潮位曲線Fig.6 GPS, BDS, Galileo and GLONASS respectively inversion tidal level results, inversion deviations and measuring tide level curves at tide gauge station
圖7a是GPS、BDS、Galileo、GLONASS反演結(jié)果組合得到的訓(xùn)練集序列,4個(gè)系統(tǒng)組合后反演點(diǎn)數(shù)目明顯增多,反演點(diǎn)數(shù)目達(dá)到1 487個(gè),平均每0.5 h一個(gè)反演值,相比于GLONASS系統(tǒng),時(shí)間分辨率提高了59.3%以上,極大地增強(qiáng)了GNSS-MR潮位監(jiān)測的時(shí)間分辨率。圖7b是與訓(xùn)練集相對應(yīng)的反演偏差序列,反演偏差在±0.5 m區(qū)間內(nèi)變化,訓(xùn)練集的RMSE值為0.140 m,由此可見,直接將不同系統(tǒng)反演結(jié)果疊加在一起,并不能提高潮位反演精度。圖7c 是訓(xùn)練集經(jīng)LS-SVM訓(xùn)練后得到的LS-SVM潮位反演模型估計(jì)的潮位時(shí)間序列。由圖7d可知,LS-SVM潮位模型估計(jì)值的偏差幅值多集中在-0.2~0.2 m區(qū)間內(nèi),較訓(xùn)練集偏差明顯減小。經(jīng)分析,LS-SVM潮位反演模型估計(jì)潮位值的RMSE值為 0.053 m,較訓(xùn)練集RMSE值減小了62.1%,較單系統(tǒng)潮位反演RMSE最小的BDS系統(tǒng),RMSE值減小了55.8%,較單系統(tǒng)潮位反演時(shí)間分辨率最高的GLONASS系統(tǒng),時(shí)間分辨率提高了59.3%。由此可見,LS-SVM融合模型在提高潮位反演時(shí)間分辨率的同時(shí),也提高了潮位反演精度。
圖7 4種系統(tǒng)組合的訓(xùn)練集和經(jīng)LS-SVM融合后的潮位反演結(jié)果及相應(yīng)的反演偏差序列Fig.7 The training set of the four system combination and the inversion result of tidal level after LS-SVM fusion and the corresponding inversion deviation sequence
圖8是LS-SVM融合模型反演潮位值與實(shí)測潮位值的相關(guān)性散點(diǎn)圖,經(jīng)LS-SVM融合后,反演潮位值與實(shí)測潮位值更緊密地分布在直線附近,反演潮位值與實(shí)測值的相關(guān)系數(shù)為0.992 5,較訓(xùn)練集提高了4.2%;較單系統(tǒng)潮位反演相關(guān)性系數(shù)最大的GLONASS系統(tǒng),相關(guān)性系數(shù)提高了3.0%,有效提高了潮位監(jiān)測的精度。
圖8 LS-SVM多系統(tǒng)融合反演潮位值與驗(yàn)潮站實(shí)測潮位值相關(guān)性圖Fig.8 Correlation diagram of LS-SVM multi system fusion inversion of tidal level and measured tidal level at tide gauge station
為了分析LS-SVM潮位融合算法的性能,本節(jié)將對比LS-SVM、SVR、滑動窗口最小二乘法3種潮位融合算法?;赟VR的多系統(tǒng)融合方法同LS-SVM方法類似,對多系統(tǒng)組合訓(xùn)練集進(jìn)行SVR訓(xùn)練,得到潮位反演模型,進(jìn)而估計(jì)出潮位值。本文在利用SVR進(jìn)行多系統(tǒng)數(shù)據(jù)融合時(shí),選用高斯徑向核函數(shù)RBF和網(wǎng)格尋優(yōu)法進(jìn)行參數(shù)優(yōu)選,不敏感損失系數(shù)ε設(shè)定為0.1,SVR算法詳見文獻(xiàn)[19]。
圖9a和9b分別為基于SVR的數(shù)據(jù)融合潮位反演結(jié)果圖和反演偏差圖。經(jīng)SVR融合后,反演潮位的RMSE值為0.111 m,較訓(xùn)練集RMSE值減小了20.7%;較BDS潮位反演RMSE值減小了7.5%,反演精度有所提高。圖9c和9d分別為基于滑動窗口最小二乘法的多系統(tǒng)融合潮位反演結(jié)果和偏差圖,同訓(xùn)練集相比,偏差變化幅值明顯減小。經(jīng)滑動窗口最小二乘法融合后,RMSE值為0.090 m,較訓(xùn)練集RMSE減小了35.7%,較BDS系統(tǒng)RMSE值減小了25%。圖10為基于SVR算法和滑動窗口最小二乘法的多系統(tǒng)融合潮位反演結(jié)果與驗(yàn)潮站實(shí)測值的相關(guān)性圖,經(jīng)分析,SVR反演結(jié)果與驗(yàn)潮站實(shí)測值之間的相關(guān)系數(shù)為0.971 1,較訓(xùn)練集相關(guān)性系數(shù)提高了2.0%;滑動窗口最小二乘法反演結(jié)果與驗(yàn)潮站實(shí)測值的相關(guān)系數(shù)為0.980 4,較訓(xùn)練集提高了2.9%。由此可知,這兩種方法均能提高潮位反演精度。此外,基于滑動窗口最小二乘法的多系統(tǒng)融合潮位反演方法的反演精度為厘米級,優(yōu)于基于SVR的多系統(tǒng)融合潮位反演方法。
圖9 基于SVR算法和滑動窗口最小二乘法的多系統(tǒng)數(shù)據(jù)融合潮位反演結(jié)果及偏差Fig.9 The inversion result and deviation of multi system data fusion tide level based on SVR algorithm and sliding window least square method
圖10 基于SVR算法和滑動窗口最小二乘法反演潮位與實(shí)測潮位的相關(guān)性圖Fig.10 Correlation diagram between the inversion tide level based on the SVR algorithm and the sliding window least square method and the measured tide level
由表3可知,3種多系統(tǒng)融合模型均能提高GNSS-MR反演潮位的精度和時(shí)間分辨率,且時(shí)間分辨率相同,均為0.5 h。LS-SVM模型較SVR模型,RMSE值減小了52.3%,相關(guān)系數(shù)提高2.2%;較滑動窗口最小二乘法,RMSE值減小了41.1%,相關(guān)系數(shù)提高了1.2%。LS-SVM融合模型優(yōu)于SVR融合模型和滑動窗口最小二乘法融合模型。
表3 3種融合模型潮位反演性能對比Tab.3 Comparison of inversion performance of three kinds of fusion models
本文針對單系統(tǒng)GNSS-MR潮位監(jiān)測中時(shí)間分辨率和反演精度低的問題,提出一種基于LS-SVM的多系統(tǒng)融合潮位反演方法,利用香港HKQT站提供的衛(wèi)星觀測數(shù)據(jù),對比了滑動窗口最小二乘法、SVR算法和LS-SVM算法在多系統(tǒng)融合潮位反演中的表現(xiàn),結(jié)果表明:
(1)單系統(tǒng)GNSS-MR潮位反演的RMSE最小值為0.120 m、最大值為0.139 m;相關(guān)系數(shù)最大值為 0.963 2,最小值為0.953 6;GLONASS潮位反演的時(shí)間分辨率最高,BDS潮位反演的時(shí)間分辨率最低。單系統(tǒng)潮位監(jiān)測的時(shí)間分辨率和反演精度偏差無法滿足高精度、高時(shí)間分辨率潮位監(jiān)測的要求。
(2)在多系統(tǒng)融合GNSS-MR潮位反演方面,基于LS-SVM的多系統(tǒng)融合潮位反演比滑動窗口最小二乘法、SVR算法具有更優(yōu)的潮位反演性能。相比于SVR模型,LS-SVM模型的RMSE值減小了52.3%,相關(guān)系數(shù)提高了2.2%;相比于滑動窗口最小二乘法,LS-SVM模型的RMSE值減小了41.1%,相關(guān)系數(shù)提高了1.2%。
本文所用方法較單系統(tǒng)而言,時(shí)間分辨率提升到0.5 h,為進(jìn)一步提高時(shí)間分辨率,在未來可以采用小波變換、NTFT等視頻變換方法,提取瞬時(shí)反射信號的頻率值,達(dá)到提升時(shí)間分辨率的效果。