閆香蓉 王曉春 梁薇 穆雅璇 宋雨龍 丁楠 張文淵

0 引言

水汽在大氣中含量很少,僅占大氣總量的0～4%,同時水汽分布極不均勻,在短時間內變化速率快,常規(guī)一日兩次探空氣球獲取的稀疏氣象數據難以滿足實時監(jiān)測需求．臺風、暴雨等極端天氣都與水汽時空變化密切相關,準確快速地獲取高時空分辨率的水汽觀測資料對天氣、氣候變化研究和天氣預報業(yè)務具有重要意義[1]．

目前,常見的大氣水汽探測方法有:無線電探測儀、氣象衛(wèi)星觀測、水汽輻射計、紅外雷達探測、太陽光譜分析儀等[2]．相較傳統(tǒng)的大氣水汽監(jiān)測方法,通過遙感(Remote Sensing,RS)紅外數據和全球導航衛(wèi)星系統(tǒng)(Global Navigation Satellite System,GNSS)數據都可以獲取大氣可降水量(Precipitable Water Vapor,PWV)．遙感紅外數據具有大范圍、高時空分辨率的優(yōu)勢,通過不同時空亮溫和輻射傳輸比直接反演大氣水汽含量,一般應用分裂窗方法來反演,但是,遙感數據對低層的水汽難以測量[3]．GNSS水汽探測具有連續(xù)運行、高精度、高時空分辨率等優(yōu)點,利用GNSS技術進行大氣水汽反演往往需要獲得相應的地面氣壓和地面溫度等氣象數據,但由于大部分GNSS觀測站建設時沒有安裝氣象傳感器,致使這些數據無法應用于氣象變化研究[4],且水汽信息受觀測站網密度影響較大．因此,合適的插值方法可以有效改善GNSS數據空間分辨率受站網密度影響的問題．

諸多學者圍繞大氣可降水量的季節(jié)性、區(qū)域性、變化特征與實際降水量的關系以及在暴雨預報中的應用等方面對各類插值方法進行了具體的研究[5-9]．李穎等[10]應用交叉驗證法對反距離加權法、全局多項式法、局部多項式法、徑向基函數法、普通克里金法5種常用空間插值方法對PWV數據進行插值,對比分析各方法的精度,得出了隨著大氣可降水量平均值的提高,插值的誤差趨于增大而相對誤差趨于減小的結論．其中空間反距離權重法(Inverse Distance Weighting,IDW)的權重冪數和搜索半徑會影響插值精度,對相關參數的估計過分依賴于經驗,缺乏理論支持;對未知點數據估計的好壞依賴于樣本點的布局[11],當采樣點分布不均勻時,IDW方法會帶來較大的插值誤差．

針對上述問題,本文提出了一種顧及GNSS水汽分布特性和站間距離的優(yōu)化IDW插值方法,該方法通過分析GNSS大氣水汽特性與站點間距離對插值結果的影響,來確定最優(yōu)化的插值參數,進而獲得高質量的大氣水汽插值結果．實驗利用徐州地區(qū)的大氣水汽數據進行實例驗證,結果表明新方法反演高分辨率大氣水汽可適用于GNSS站網布設較為稀疏的地區(qū),有效提高了大氣水汽重構精度．

1 常用的插值方法原理

本文使用3種性能良好的傳統(tǒng)插值方法進行交叉驗證,對比分析各方法的插值質量,并對IDW方法優(yōu)化,提高計算得到的PWV數據質量及穩(wěn)定性．

1.1 線性插值三角網法

線性插值三角網法(Triangulation with Linear Interpolation),是指將二維空間中具備參數值的多個空間點通過直線連接,構建若干插值三角形進行計算,建立起來的所有三角形面積覆蓋整個研究區(qū)域,其中任一三角形的邊不能與其他三角形相交,組成由若干個三角形拼接形成的研究網．線性插值三角網法能夠很好體現原始數據的特征,它將整個研究區(qū)域內的數據進行均勻分配,分配完成后可在數據分布較少的區(qū)域得到一個擁有獨特性質的三角面,以此計算得到一個最優(yōu)插值結果．其計算公式為

(1)

式中:已知測站點的坐標為(x0,y0)和(x1,y1),yi為測區(qū)內一點(xi,yi)處的線性插值計算結果．雖然使用線性插值三角網法能較快計算得到插值結果,但其具有一定程度的局限性[12],受大氣可降水量平均值影響較大,同時使用該方法計算的數據量少時難以保證其插值計算結果的平滑性[13]．

1.2 克里金插值法

克里金插值法,即Kriging插值法,也稱空間自協(xié)方差最佳插值方法．Kriging方法考慮觀測點與估計點之間的相對位置信息,利用觀測點之間的空間位置信息對待求點進行無偏和最優(yōu)估計[14]．克里金插值法的公式為

(2)

式中:ki為待定(克里金)權重系數,其值之和等于1,插值計算中觀測站的個數為n,觀測站的位置為xi,Z(xi)為在xi處的觀測值．

1.3 空間反距離加權法

空間反距離權重法,也稱為反距離加權法．該方法主要依賴于反距離的冪值,它的冪參數可以根據已知點和待插值點的距離來控制已知點對插值點的精度影響．該方法用周邊采樣點的值,估計未知點的值,以待插點與實際觀測樣本點之間的距離為權重,離插值點越近的樣本點賦予的權重越大,其權重貢獻與距離成反比．反距離加權法[15]的公式為

(3)

式中S處的預測值為Z(S0),計算過程使用的樣本數量為n,di為插值點與已知點間的距離,Zi為對應樣本點的值,p為冪參數,一般情況下,p取2．

2 顧及GNSS水汽特性和站間距離的優(yōu)化IDW插值方法

2.1 GNSS數據介紹

本文以徐州市2017年的GNSS數據為基礎,考慮到徐州地區(qū)夏季多雨冬季少雨的降水空間分布特性,選取5—7月降水相對集中且充足的時段開展實驗,選取6個站點作為研究對象,測站名分別為:CHG0、DSH0、HJG0、SHG0、TSG0、TUG0,IGS輔助站分別為:YAKT、GUAM、POL2,探空氣象數據為2017-05-11—2017-07-31有完整數據的共75 d的數據．實驗數據GNSS站點分布如圖 1所示．

圖1 徐州地區(qū)GNSS測站點分布Fig.1 Distribution of GNSS stations in Xuzhou

實驗數據利用GAMIT/GLOBK軟件進行解算,以同步單天時段為單位進行,采用GPS和精密星歷組合方式進行基線解算,最終得到全天候各測站的PWV,各參數配置如表 1所示．

表1 GAMIT參數配置

在進行插值運算前,首先對數據進行預處理,對徐州地區(qū)6個GNSS測站觀測數據進行空間自相關分析,探究各GNSS測站的PWV與λ1(該測站到探空站的距離)、λ2(該測站的水汽密度)和λ3(該測站的水汽壓)3個參數的相關性．通過單點定位的方式獲取各測站坐標,根據歐氏距離計算公式解算各測站到探空站的距離．大氣中的水汽絕大部分位于對流層,但它在垂直方向上隨著高度的升高迅速減小[15],在研究地區(qū)GNSS地表測站數據中,統(tǒng)一選擇地面最低起始高度處的水汽壓和水汽密度數據,用相關系數度量2個變量X(各GNSS測站PWV)和Y(λ1,λ2,λ3)之間的相關性．通過表 2發(fā)現各測站PWV與距離相關性最大且干擾因素較少,而與水汽密度和水汽壓相關性接近,且受其他水汽參數影響較大,故本文針對上述結果提出2種對PWV插值方法的優(yōu)化策略．

表2 各采樣站點PWV與λ1,λ2,λ3相關性

2.2 顧及水汽特性及其變化程度的優(yōu)化策略

根據2.1實驗,PWV與水汽壓、水汽密度相關性接近,且受水汽分布影響較大,因此提出顧及水汽特性及其變化程度的優(yōu)化策略．PWV可由天頂濕延遲(Zenith Wet Delay,ZWD)轉換得到[16],ZWD主要受信號傳播路徑上水汽含量的影響,其中包括水汽密度和水汽壓．

PWV=Π×ZWD,

(4)

(5)

式中,Π為無量綱水汽轉換系數,ρv為液態(tài)水的密度,其值為1 g/cm3,Rw為水汽的氣體常數,k′2,k3為轉換系數,ZWD通過GAMIT高精度GNSS數據處理軟件獲?。?/p>

(6)

Tm=a+bTs,

(7)

式中,Tm為大氣加權平均溫度,e為水汽壓,Ts為絕對溫度．

通過上述公式可知,計算PWV需要已知Tm數值,本實驗采用王明華[17]建立的中國區(qū)域統(tǒng)一模型與各單站的Tm-Ts線性模型,但Tm模型計算值與參考值的偏差隨著季節(jié)不同有明顯變化[18]．

本實驗基于徐州地區(qū)2017年5月1日前3天的探空站數據,將水汽密度和水汽壓分別強制約束為某一值,建立先驗水汽場,通過該水汽場給定Tm模型3個初始值進行約束,增加模型解算結果的穩(wěn)定性,從而提高水汽場的反演精度,得到高精度的PWV,為GNSS站間距離優(yōu)化做準備．

2.3 顧及GNSS站間距離的優(yōu)化策略

由相關性分析實驗可知,λ1(測站到探空站的距離)是3種參數中對PWV插值精度影響最大的因子,只受空間分布的影響,因此提出顧及GNSS站間距離的優(yōu)化策略．根據IDW的數學特性,采樣站點距離與插值精度有關[19],且距離權重函數與采樣點到觀測點的距離次冪成反比[20]．本次實驗采樣站點距離探空站由近及遠依次是:TSG0、HJG0、DSH0、TUG0、CHG0、SHG0．為確定哪幾個測站的聯合插值精度最高,依據采樣點距離相關性依次減少測站數來設計實驗,各個實驗測站數量和位置如圖2所示．

圖2 各方案測站數量和站點分布Fig.2 Number of stations and station distribution map for 4 schemes

本次實驗比較不同插值方法計算得出的最大值、最小值、平均值和標準偏差4項指標．根據表 3,表明方案B:TSG0、HJG0、DSH0、TUG0、CHG0這5個測站插值效果最佳．

表3 不同個數采樣站點插值效果統(tǒng)計

本實驗利用徐州地區(qū)采樣站點分布不均的空間特性,對IDW插值參數優(yōu)化,并基于本地區(qū)的先驗信息加以約束,提出了顧及GNSS水汽特性和站間距離的優(yōu)化IDW插值方法．具體過程如下:

1)基于實驗數據前3天的探空站水汽密度和水汽壓數據建立先驗水汽場,對Tm模型附加初值進行

約束,對IDW參數優(yōu)化;

2)對方案B的5個測站進行反距離權重插值,在以格點為中心的360°方位上,改變搜索橢圓的長短半徑,使其全部覆蓋實驗區(qū),實驗區(qū)不規(guī)則時可將搜索橢圓旋轉一定角度;

3)平均劃分成4個象限,在每個象限區(qū)選擇測站點,每個象限都有站點數據參與插值;

4)優(yōu)化冪參數,將p從2到10分別實驗,通過交叉檢驗比較插值誤差大小,基于基礎數據,取p為5.875 474 8時,插值誤差較?。?/p>

3 結果分析

本文插值檢驗主要采用交叉驗證的方法來對比分析插值結果,即假定一個已知氣象站點的數據未知,使用其他氣象站點的實測值來估算該點的數據,并進行誤差分析．研究過程中采用標準差(Standard Deviation,SD)、平均相對誤差(Mean Relative Error,MRE)、平均絕對誤差(Mean Absolute Error,MAE)及均方根誤差(Root Mean Square Error,RMSE)作為評估插值方法精度的標準．

3.1 插值方法精度分析

本文采用線性插值三角網法(方法1)、IDW插值法(方法2)、Kriging插值法(方法3)和顧及GNSS水汽特性和站間距離的優(yōu)化IDW插值方法(方法4)4種插值方法對數據進行計算．

以往的觀測實驗表明,PWV與局地降水存在密切的關系,每次降水過程都對應著PWV的迅速增加[21]．由圖 3可知,徐州地區(qū)大氣可降水量總體呈現西北至東南逐漸增加的空間格局,西北少、東南多,梯度變化較為明顯,PWV值相差較大．不同插值方法所得結果在色彩平滑程度以及局部地區(qū)的空間分布上存在一定差異:方法3和方法4色彩過渡較為強烈,而方法2的插值結果較為平滑,但基于所選數據為突降暴雨前,說明采用方法4對徐州地區(qū)PWV進行插值,能夠減弱采樣站點分布不均及大氣總水汽量激增的誤差．

圖3 4種不同插值方法的色彩Fig.3 Color images of four interpolation methods a.linear interpolation triangulation;b.Kriging interpolation;c.IDW;d.the improved IDW interpolation

3.2 與探空站數據對比

本次實驗采用SD、MAE、MRE、RMSE作為評估插值方法精度的標準,將徐州地區(qū)探空站每天UTC12:00的數據作為參考標準,對4種插值方法進行精度評定．由表 4可知:方法4的SD數值較小,數據集分布較均勻;相較傳統(tǒng)的IDW插值方法(方法2),方法4的MAE和MRE分別提高了19.48%、19.16%．RMSE可以反映誤差分布大致情況以及誤差在基準線的集中程度,表4表明方法4誤差較集中,其RMSE較方法1、方法2、方法3分別提高了14.88%、15.70%、4.12%．

表4 4種不同插值方法計算誤差交叉驗證結果

進一步對比不同插值方案反演PWV誤差的大小,以探空站數據為參考值對各插值結果求殘差值,并繪制箱型圖．由圖 4可以看出,4種方法的誤差對應的箱形圖具有自身的特點,線性插值三角網法和IDW方法具有的異常值數值較大,離散程度也較大,而Kriging和優(yōu)化IDW插值方法的異常值較小,也較集中．箱形圖各個層上優(yōu)化IDW插值方案的誤差分布較其他方案更為集中,殘差數據較小．

圖4 4種不同插值方法的箱形圖Fig.4 Box-plot of four interpolation methods

4 結語

本文提出顧及GNSS水汽特性和站間距離的優(yōu)化IDW插值方法,基于徐州地區(qū)的先驗信息,比較分析GNSS大氣水汽特性與站點間距離對插值結果的影響,對插值參數進行優(yōu)化,主要結論如下:

1)在3種傳統(tǒng)插值方法中,插值精度的RMSE值均超過3.5 g/m3,其中線性插值三角網法與IDW插值法的插值精度較低．與傳統(tǒng)算法相比,本文所提出的顧及GNSS水汽特性和站間距離的優(yōu)化IDW插值方法能夠有效改善插值結果的質量,其RMSE值為3.49 g/m3,而穩(wěn)定性方面,新方法的誤差分布相較3種傳統(tǒng)方法更加集中,且異常值也較少．

2)在大氣總水汽量激增的情況下,該方法反演的三維水汽色彩圖更能反映真實情況,結合表4可知擬合后探空站位置處的插值精度較高．大氣總水汽量激增常伴隨著降水過程的發(fā)生,此時水汽含量不同區(qū)域PWV分界明顯,在圖上表現為色彩過渡強烈,表明該方法能夠減小采樣站點分布不均及大氣總水汽量激增帶來的誤差．

3)新方法反演大氣水汽可適用于GNSS站網布設較為稀疏的地區(qū),采用該方法能夠在計算結果上減弱大氣可降水量平均值的影響,有效提高大氣水汽時空分辨率精度,估值準確率高．該方法基于一個地區(qū)的先驗信息,通過對其地區(qū)GNSS大氣水汽特性和站網間距離統(tǒng)計分析,對IDW參數優(yōu)化,有效提高了插值精度．