孫 為 朱明晨,2

1 銅陵學(xué)院建筑工程學(xué)院，安徽省銅陵市翠湖四路1335號，2440612 東南大學(xué)交通學(xué)院，南京市東南大學(xué)路2號，218889

全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)能夠增加可觀測衛(wèi)星的數(shù)量和多樣性[1]，極大推動了衛(wèi)星導(dǎo)航定位技術(shù)的發(fā)展。GNSS信號穿過大氣時，會受到電離層和對流層折射的影響而產(chǎn)生延遲，是GNSS的一個主要誤差源。電離層屬于彌散性介質(zhì)，其延遲效應(yīng)可以通過頻率間差分消除，而對流層屬于中性大氣，信號傳播速度與頻率無關(guān)，無法采用類似方法消除延遲效應(yīng)[2]。對流層天頂延遲(zenith total delay，ZTD)一般能達(dá)到2.3 m左右[3]，可分為靜力學(xué)延遲(zenith hydrostatic delay，ZHD)和濕延遲(zenith wet delay，ZWD)兩個部分。靜力學(xué)延遲較為穩(wěn)定，但濕延遲極易受天氣效應(yīng)影響，隨機(jī)性極強(qiáng)，難以用模型估計[4]，所以對流層延遲的精確估計是提高GNSS實(shí)時導(dǎo)航定位精度的關(guān)鍵技術(shù)之一。

目前國內(nèi)外常用的對流層模型可分為氣象和無氣象參數(shù)模型兩大類。經(jīng)典的氣象模型有Hopfield模型[5]、Saastamoinen模型[6]和Black模型[7]等，這些模型能夠?qū)Ω裳舆t進(jìn)行高精度估計，但對濕延遲的估計效果并不理想[8]。此類模型的應(yīng)用依賴于地表氣象參數(shù)，這極大限制了模型的應(yīng)用范圍。無氣象參數(shù)模型或是對地表氣象參數(shù)進(jìn)行經(jīng)驗(yàn)建模，再代入第一類模型中計算，或是直接對ZTD進(jìn)行時空分析，這兩種方法都能擺脫對氣象參數(shù)的依賴，應(yīng)用范圍更廣泛。無氣象參數(shù)模型常用的有EGNOS模型[9]、GPT系列模型[10-11]和CTrop模型[12]等，這些模型僅需輸入測站的概略坐標(biāo)和年積日即可獲取經(jīng)驗(yàn)ZTD延遲。由于對流層延遲變化復(fù)雜，難以通過一個固定的數(shù)學(xué)模型進(jìn)行表達(dá)，人工神經(jīng)網(wǎng)絡(luò)技術(shù)可為解決這一問題提供途徑[13]，其強(qiáng)大的自適應(yīng)學(xué)習(xí)能力和非線性逼近能力可以避免參數(shù)化建模的缺點(diǎn)，從而通過權(quán)值的自動調(diào)整達(dá)到最佳的逼近效果。

中國大陸構(gòu)造環(huán)境監(jiān)測網(wǎng)(Crustal Movement Observation Network of China，CMONOC)提供260多個連續(xù)監(jiān)測站，用于監(jiān)測中國區(qū)域地殼運(yùn)動、重力場變化以及電離層與對流層參數(shù)。CMONOC提供逐小時由GAMIT/GLOBK軟件解算的高精度ZTD產(chǎn)品(下文簡稱ZTDC)，可以進(jìn)行對流層延遲的建模[14]。陳俊平等[3]通過分析ZTDC的時空變化規(guī)律，構(gòu)建2.5°×2.0°分辨率的SHAtrop模型。該模型考慮到高程改正，僅需要輸入經(jīng)緯度與時間即可滿足中國區(qū)域GNSS用戶實(shí)時導(dǎo)航ZTD改正的需求，優(yōu)于常見的EGONS模型、UNB3模型和GPT2模型。采用CMONOC提供的對流層延遲產(chǎn)品，基于BP神經(jīng)網(wǎng)絡(luò)技術(shù)，建立多個弱神經(jīng)網(wǎng)絡(luò)非線性映射，集成構(gòu)造BP-Adaboost模型。用戶僅需輸入測站的概略坐標(biāo)、年積日和時間，即可獲得高精度的ZTD估計值。將新模型與GPT3模型、CTrop模型和SHAtrop模型進(jìn)行比較，可證明建模方法的有效性。

1 建模數(shù)據(jù)與驗(yàn)證數(shù)據(jù)

選取2014～2019年具有連續(xù)觀測數(shù)據(jù)的CMONOC測站，剔除其中誤差大于5 mm的ZTD值，共得到測站217個。均勻選取其中155個測站2014～2018年的數(shù)據(jù)進(jìn)行建模，利用剩余62個測站2019年的數(shù)據(jù)進(jìn)行驗(yàn)證，得到訓(xùn)練樣本1 883 076個、檢驗(yàn)樣本136 958個。為進(jìn)一步驗(yàn)證模型的適用性，同時選取由IGRA(integrated global radiosonde archive)提供的中國區(qū)域內(nèi)86個探空站的數(shù)據(jù)，剔除探空數(shù)據(jù)中觀測層數(shù)少于12或頂層水汽壓大于0.015 hPa的數(shù)據(jù)，得到檢驗(yàn)樣本48 910個，計算各個探空站的ZTD精確值(下文簡稱ZTDR)。建模測站和驗(yàn)證測站的分布如圖1所示。

圖1 CMONOC與探空測站分布Fig.1 Distribution of CMONOC and radiosonde stations

由于探空數(shù)據(jù)僅給出各個等壓面的氣溫、水汽壓、位勢高度等參數(shù)，可以利用Saastamoinen公式計算干延遲、利用積分方法計算濕延遲，并相加得到測站最終的ZTD值[8]：

(1)

式中，P為大氣壓(單位hPa)；φ為緯度；h為大地高；Nwi為探空數(shù)據(jù)第i層大氣折射率指數(shù)的濕分量，可由公式(2)計算得到：

(2)

式中，Rd和Rw為干空氣和濕空氣常數(shù)，k1、k2和k3為折射率常數(shù)，e為水汽壓(單位hPa)，T為熱力學(xué)溫度[15-16]。探空數(shù)據(jù)使用的位勢高由位勢除以重力常數(shù)(gn= 9.806 65 m/s2)得到，并不嚴(yán)格等于幾何高度。因此為統(tǒng)一探空站與CMONOC測站的高程系統(tǒng)，聯(lián)合公式(3)和公式(4)將位勢高轉(zhuǎn)換成幾何高度并最終轉(zhuǎn)換為大地高[2]：

(3)

g(φ,h)=gn(1-0.002 637 3cos(2φ)+

5.9×10-6cos2(2φ))·(1-3.14×10-7h)

(4)

式中，hd為位勢高(單位m)，C為位勢，h為正高(單位m)。采用平均偏差(bias)和均方根誤差(RMS)作為評定標(biāo)準(zhǔn)。平均偏差表示準(zhǔn)確度，即模型與真值的偏離程度；均方根誤差表示精度，用于衡量模型的可靠性和穩(wěn)定性[17]，其計算公式為：

(5)

2 BP-Adaboost模型的建立

BP-Adaboost算法的思想是合并多個“弱”神經(jīng)網(wǎng)絡(luò)的輸出，進(jìn)而得到更加精確的估計值[18-19]。其算法步驟如下：

1)神經(jīng)網(wǎng)絡(luò)初始化及訓(xùn)練。弱預(yù)測器選取前饋神經(jīng)網(wǎng)絡(luò)(feedforward neural network，F(xiàn)NN)。FNN網(wǎng)絡(luò)由輸入層、隱含層、輸出層以及模型參數(shù)組成，模型的學(xué)習(xí)誤差ε取0.1 mm。FNN網(wǎng)絡(luò)輸入層節(jié)點(diǎn)數(shù)為5，分別代表測站經(jīng)度、緯度、高程、年積日(day of year，DOY)和小時；輸出層節(jié)點(diǎn)數(shù)為1，代表ZTD(單位mm)。隱含層節(jié)點(diǎn)的選取無規(guī)律可循，一般通過試算確定。節(jié)點(diǎn)數(shù)過多可能導(dǎo)致學(xué)習(xí)時間過長且結(jié)果不穩(wěn)定，節(jié)點(diǎn)數(shù)過少會導(dǎo)致網(wǎng)絡(luò)不收斂。單個FNN網(wǎng)絡(luò)基本結(jié)構(gòu)如圖2所示。

圖2 FNN結(jié)構(gòu)Fig.2 Feed forward neural network structure

激活函數(shù)選取Tan-Sigmoid，其表達(dá)式為：

(6)

FNN訓(xùn)練采用最速下降法，調(diào)用MATLAB的toolbox。

2)累積誤差計算。共選取20個FNN網(wǎng)絡(luò)，將其合并成BP-Adaboost強(qiáng)預(yù)測器，用于精確計算ZTD。將n個訓(xùn)練樣本的分布權(quán)值初始化為D0(i)=1/n,i=1,2,3,…,n，每個網(wǎng)絡(luò)的累積誤差初始化為Error(t)=0,t=1,2,3,…, 20。若第t個FNN(FNNt)得到的第i個樣本的預(yù)測誤差大于50 mm，則認(rèn)為該樣本是需要加強(qiáng)學(xué)習(xí)的樣本，需要按照公式(7)更新該網(wǎng)絡(luò)的累積誤差Error(t)及樣本的分布誤差Dt，否則FNNt對應(yīng)的累積誤差及樣本的分布權(quán)值不變：

(7)

具體的累積誤差計算流程如圖 3所示。

圖3 累積誤差計算流程Fig.3 Cumulative error calculation flow chart

3)各弱預(yù)測器權(quán)重計算。第t個FNN網(wǎng)絡(luò)的權(quán)重Pt可由公式(8)計算并歸一化得到：

(8)

最后利用BP-Adaboost計算ZTD時，取這20個FNN網(wǎng)絡(luò)預(yù)測結(jié)果的加權(quán)平均值作為最終的ZTD輸出。算法流程如圖4所示。

圖4 BP-Adaboost算法流程Fig.4 The steps of BP-Adaboost algorithm

為確定第1步中各FNN網(wǎng)絡(luò)隱含層節(jié)點(diǎn)的個數(shù)，試算隱含層節(jié)點(diǎn)數(shù)為7～36共30種BP-Adaboost模型，由CMONOC驗(yàn)證站2019年ZTDC進(jìn)行檢驗(yàn)，得到其bias和RMS(圖5)。

圖5 不同隱含層節(jié)點(diǎn)數(shù)結(jié)構(gòu)的bias和RMSFig.5 Bias and RMS with different hidden layer node number structure

由圖5可見，當(dāng)隱含層節(jié)點(diǎn)數(shù)少于27時，其bias絕對值和RMS均呈現(xiàn)出下降趨勢；當(dāng)節(jié)點(diǎn)數(shù)超過27時，bias與RMS趨于平緩，精度提升有限。所以隱含層節(jié)點(diǎn)數(shù)確定為27，選取的FNN網(wǎng)絡(luò)結(jié)構(gòu)為5×27×1。

3 精度驗(yàn)證

為討論BP-Adaboost模型的適用性，利用2019年CMONOC的141個參與建模測站(CMONOC建模站)、62個未參與建模測站(CMONOC驗(yàn)證站)和86個探空站的ZTD精確值分別比較CTrop模型、GPT3模型以及SHAtrop模型在中國區(qū)域的精度，結(jié)果見表1(單位mm)。

由表1可見，BP-Adaboost模型的精度最優(yōu)。當(dāng)用CMONOC建模站和驗(yàn)證站進(jìn)行數(shù)據(jù)評定時，BP-Adaboost的bias僅為0.62 mm和-1.16 mm，表明該模型沒有明顯的系統(tǒng)誤差，體現(xiàn)出和GNSS數(shù)據(jù)相當(dāng)程度的自洽性。其RMS比CTrop模型、GPT3模型和SHAtrop模型分別減小9.6%、2.8%和3.3%。當(dāng)用探空站進(jìn)行數(shù)據(jù)驗(yàn)證時，新模型bias為-12.32 mm，其絕對值為4種模型中最小，其RMS與GPT3模型相當(dāng)，但比其他2種模型分別提高3.9%和2.5%。

表1 4種模型的bias與RMS

3.1 模型精度的空間分析

為進(jìn)一步分析BP-Adaboost模型的空間適用性，繪制出各個CMONOC站與探空站的bias與RMS分布(圖6)。

圖6 BP-Adaboost模型各測站的bias與RMS分布Fig.6 Bias and RMS distribution of each station in BP-Adaboost model

由圖6可見，BP-Adaboost模型的RMS在整體上呈現(xiàn)出東南區(qū)域精度低于西北區(qū)域的分布規(guī)律，這可能是因?yàn)橹袊臇|南區(qū)域水汽變化較活躍，難以建模。進(jìn)一步定量討論模型精度的水平分布，將中國區(qū)域按緯度劃分為41°～50° N、31°～40° N和21°～30° N三個區(qū)間段，在3個緯度內(nèi)的bias和RMS如表2(單位mm)所示。

由表2可見，4種模型均在31°～40° N內(nèi)取得最高精度。在高緯度地區(qū)，BP-Adaboost模型均取得最優(yōu)精度，3類驗(yàn)證數(shù)據(jù)表明，在41°～50° N區(qū)間段內(nèi)其RMS分別優(yōu)于GPT3模型3.8%、4.7%和1.7%；而在31°～40° N內(nèi)BP-Adaboost模型精度分別優(yōu)于GPT3模型2.1%、1.5%和3.2%，充分說明BP-Adaboost模型在不同緯度的適用性。新模型在低緯度精度較差的原因可能是低緯度的CMONOC測站較少、建模數(shù)據(jù)有限。

表2 不同緯度區(qū)間4種模型的bias與RMS

為進(jìn)一步分析這4種模型在不同海拔高度的適用性，按照測站高程，將各個測站按照0～2 km、2～4 km以及大于4 km這3個區(qū)間分別統(tǒng)計，如表3(單位mm)所示。

表3 3個高程區(qū)間的bias和RMS

由表3可見，各個模型精度均隨海拔高程升高而增加。在高程0～2 km區(qū)間內(nèi)，兩組CMONOC測站數(shù)據(jù)驗(yàn)證結(jié)果中BP-Adaboost模型都表現(xiàn)出最優(yōu)的精度，在探空數(shù)據(jù)驗(yàn)證結(jié)果中BP-Adaboost與GPT3模型相當(dāng)，均優(yōu)于其他2種模型；在2～4 km區(qū)間內(nèi)，兩組CMONOC測站數(shù)據(jù)驗(yàn)證中BP-Adaboost模型均取得最優(yōu)精度，但探空數(shù)據(jù)的結(jié)果顯示，SHAtrop模型精度最高，BP-Adaboost模型精度略低；在大于4 km的區(qū)間內(nèi)沒有CMONOC驗(yàn)證站，其余2組數(shù)據(jù)的結(jié)果表明，BP-Adaboost模型精度仍然優(yōu)于其他3種模型。以上結(jié)論充分說明了BP-Adaboost模型在不同高度區(qū)間的適用性。

3.2 模型精度的季節(jié)性分析

為分析模型在不同時間段的適用精度，計算這4種模型在每個季度的RMS(圖7)。

圖7 4種模型的季節(jié)性RMSFig.7 Seasonal RMS of four models

由圖7可見，4種模型在冬季精度高而在夏季精度較低，在春秋兩季精度相當(dāng)，原因可能是夏季水汽活動較為劇烈。除夏季和秋季的探空數(shù)據(jù)驗(yàn)證結(jié)果外，BP-Adaboost模型都取得最優(yōu)的精度。 CMONOC建模測站驗(yàn)證的春季和夏季結(jié)果表明，BP-Adaboost模型精度比GPT3模型在2個季節(jié)分別提高5.2%和3.8%，表現(xiàn)出更佳的適用性和與GNSS觀測數(shù)據(jù)的自洽性。

4 結(jié) 語

本文利用CMONOC提供的155個測站2014～2018年的ZTDC產(chǎn)品，采用BP-Adaboost算法將多個弱神經(jīng)網(wǎng)絡(luò)集成為強(qiáng)預(yù)測器，從而建立無氣象參數(shù)的ZTD模型，并對模型進(jìn)行精度評定，得出以下結(jié)論：

1)利用2019年CMONOC參與建模的155個測站以及未參與建模的62個測站數(shù)據(jù)進(jìn)行驗(yàn)證，BP-Adaboost模型在中國區(qū)域內(nèi)的總體bias為0.62 mm和-1.16 mm，優(yōu)于常用的CTrop模型、GPT3模型和SHAtrop模型，表明該模型沒有明顯的系統(tǒng)誤差；其RMS為25.30 mm和26.72 mm，比上述3種模型精度分別提高9.6%、2.8%和3.3%。

2)利用IGRA提供的中國區(qū)域內(nèi)86個無線電探空站的ZTDR進(jìn)行驗(yàn)證，BP-Adaboost模型的bias為-12.32 mm，RMS為46.29 mm；其RMS略高于GPT3模型，但明顯優(yōu)于其他2種常用模型。該數(shù)據(jù)的驗(yàn)證結(jié)果中4種模型均出現(xiàn)一定的負(fù)偏差。

3)BP-Adaboost模型的精度分布呈現(xiàn)出內(nèi)陸高、東南沿海低的特性，在高緯度地區(qū)BP-Adaboost模型的優(yōu)勢更明顯，這是因?yàn)镃MONOC在內(nèi)陸地區(qū)的測站更為密集、建模數(shù)據(jù)更加豐富，證明BP-Adaboost算法建模的有效性，即如果融合更豐富的數(shù)據(jù)，該算法能進(jìn)一步提高建模精度。BP-Adaboost在高海拔地區(qū)也優(yōu)于其他3種模型。

4)4種模型的精度都呈現(xiàn)出夏季低而冬季高的特性。CMONOC兩組數(shù)據(jù)驗(yàn)證結(jié)果表明，BP-Adaboost的精度均優(yōu)于其他3種模型，其精度在春季和夏季比GPT3模型高5%左右。除秋季外，其余3個季節(jié)探空數(shù)據(jù)的驗(yàn)證結(jié)果中4種模型精度相當(dāng)。

本文初步論證集成神經(jīng)網(wǎng)絡(luò)在對流層應(yīng)用中的有效性，新模型相比于已有的模型，精度有較為顯著的提高，可以進(jìn)一步用于GNSS導(dǎo)航及GNSS氣象學(xué)等相關(guān)領(lǐng)域的應(yīng)用研究。