姜 濤

(廣東省核工業(yè)地質(zhì)局二九二大隊，廣東 河源 517001)

0 引 言

國家連續(xù)運行參考站(continuously operating referencestations，CORS)目前已在國內(nèi)大部分省市建立，且積累了大量的連續(xù)觀測數(shù)據(jù)[1-2]。這些數(shù)據(jù)為各省市的測量、導(dǎo)航、地震研究等提供了重要的數(shù)據(jù)信息，其中CORS站高程坐標(biāo)數(shù)據(jù)的信息最為重要[3-5]。為此，研究CORS站高程坐標(biāo)方向的運動規(guī)律具有一定的實際價值和意義。

從國內(nèi)不同的CORS站高程坐標(biāo)時序圖可以看出，CORS站高程坐標(biāo)運動規(guī)律存在周期性變化，且不同的CORS站其高程坐標(biāo)運動規(guī)律是不同的[6]。造成不同CORS站高程坐標(biāo)運動周期性變化規(guī)律的因素有很多種，其中包括地球潮汐、地殼運動、大氣層變化、GPS測量技術(shù)等[7]。先前對CORS站高程數(shù)據(jù)建立數(shù)學(xué)模型時，常采用線性最小二乘的方法，該方法事先將固定的整年周期值賦給線性最小二乘模型周期項，導(dǎo)致建立的模型周期值與實際周期值存在一定的偏差[8]。這種建立CORS站高程數(shù)據(jù)數(shù)學(xué)模型的方法明顯地認(rèn)為CORS站高程坐標(biāo)運動只存在線性運動，而忽視了其非線性運動特征[9]。

目前GPS測量技術(shù)已經(jīng)應(yīng)用于交通、建筑、導(dǎo)航、定位等領(lǐng)域，起到了越來越重要的作用[10]。為了對CORS站高程坐標(biāo)建立符合其實際運動規(guī)律的數(shù)學(xué)模型，本文選取了國內(nèi)上海CORS站10余年高程數(shù)據(jù)作為研究對象，對其建立非線性速度場模型[11]。建立非線性速度場模型的前提是先給CORS站高程數(shù)據(jù)建立線性最小二乘模型，其次求出線性最小二乘模型的解，最后將其作為非線性最小二乘模型的迭代初值進(jìn)行未知參數(shù)的求解，得到反映CORS站高程坐標(biāo)運動規(guī)律的非線性速度場模型。

1 CORS站高程數(shù)據(jù)建模方法

1.1 線性最小二乘建模原理

通過對CORS站高程坐標(biāo)時間序列的時序圖分析發(fā)現(xiàn)，CORS站高程方向運動規(guī)律存在明顯的年周期運動。由于不同因素的影響，CORS站高程坐標(biāo)運動周期性變化中可能存在的月、半年、兩年變化在時序圖中不能被準(zhǔn)確發(fā)現(xiàn)，并且不同的CORS站其周期性變化規(guī)律也不盡相同[12]。為了得到符合各個CORS站高程運動規(guī)律的數(shù)學(xué)模型，需進(jìn)行多次的數(shù)據(jù)擬合。包含年、半年周期項的線性最小二乘擬合模型如式(1)所示。

Yi=A+B×ti+C1×sin(2π×ti)+C2×sin(4π×ti)

(1)

式(1)中，Yi為模型擬合值，A為常數(shù)項，B為線性速度值，ti代表時間，C1表示一年周期項對應(yīng)的振幅值，C2表示半年周期項對應(yīng)的振幅值，周期值分別為T1=1,T2=0.5，該周期值表示線性模型給定的固定周期值為一年和半年。

A,B,C1,C2為4個未知數(shù)，設(shè)CORS站高程坐標(biāo)序列觀測值為Z。利用間接平差原理得到的誤差方程如式(2)所示。

V=AN-Z

(2)

式(2)中，

(3)

(4)

按照最小二乘原理VTPV=min，P為單位矩陣，求解未知參數(shù)N的值。

1.2 非線性最小二乘建模原理

非線性最小二乘擬合模型是在線性模型(1)的基礎(chǔ)上加入相應(yīng)的初相值，且周期項不是固定的數(shù)值[13]。具體模型表達(dá)式如式(5)所示。

Yi=A+B×ti+C1×sin(2π×f1×ti+φ1)+C2×sin(4π×f2×ti+φ2)

(5)

式(5)中A,B,C1,f1,φ1,C2,f2,φ2為未知參數(shù)，記為S。

將非線性模型(5)在近似值N0處線性化，線性化后的誤差方程如式(6)所示。

V=BS-Z

(6)

式(6)中，S=[ABC1f1φ1C2f2φ2]T，

qi=sin(2×π×ti+φ1)，wi=2×π×ti×cos(2×π×ti+φ1)，ei=cos(2×π×ti+φ1)，

ri=sin(4×π×ti+φ2)，di=4×π×ti×cos(4×π×ti+φ2)，fi=cos(4×π×ti+φ2)，

最后按照最小二乘原理VTPV=min，P為單位矩陣，求解未知參數(shù)S的值。

2 CORS站高程數(shù)據(jù)建模實驗分析

本文從SOPAC網(wǎng)站下載上海CORS站2005～2017年高程坐標(biāo)數(shù)據(jù)，并繪制上海CORS站高程坐標(biāo)時序如圖1所示。其中，縱軸表示的是上海站高程數(shù)據(jù)值，單位為米，該高程數(shù)據(jù)由原始高程數(shù)據(jù)減去平均值得到。橫軸表示的是高程數(shù)據(jù)對應(yīng)的觀測時間為2005～2017年。從圖1上海站高程坐標(biāo)時序圖可以看出：高程坐標(biāo)數(shù)據(jù)隨著時間的變化具有周期性，且不是線性變化，即不存在直線上升或下降的趨勢[14]。這就為高程數(shù)據(jù)建立數(shù)據(jù)模型提供了相應(yīng)的思路。

圖1 上海站高程坐標(biāo)時序圖

2.1 線性最小二乘建模

通過對上海站高程坐標(biāo)時序圖分析發(fā)現(xiàn)，高程坐標(biāo)數(shù)據(jù)隨著時間的變化其周期性可能存在一年、半年變化。為此，本文在建立線性最小二乘模型時，首先假定模型周期值為固定一年、半年數(shù)值，即T1=1,T2=0.5。然后按照線性最小二乘模型(1)進(jìn)行高程數(shù)據(jù)的擬合，得到的上海站高程數(shù)據(jù)線性最小二乘擬合結(jié)果如圖2所示，其中縱軸表示高程坐標(biāo)數(shù)據(jù)(包括觀測高程數(shù)據(jù)和線性最小二乘擬合建模所求解的高程擬合值)黑色的點表示的是原始高程數(shù)據(jù)，紅色曲線表示的是線性最小二乘擬合曲線圖。從圖2可以看出，上海站高程坐標(biāo)運動存在周期性的變化規(guī)律，且周期性變化主要是一年變化為主、半年變化占比較少。表1給出了線性最小二乘建模初值、線性最小二乘解及建模精度評價指標(biāo)均方根誤差的值[15]。

圖2 上海站高程坐標(biāo)線性最小二乘擬合圖

表1 線性最小二乘建模參數(shù)值

由表1線性最小二乘建模參數(shù)值可以看出：線性最小二乘擬合建模所需的周期項為給定的周期值，分別是一年、半年數(shù)值，由計算出的振幅值C1和C2數(shù)值可得，上海CORS站高程數(shù)據(jù)周期性變化一年周期項對應(yīng)的振幅值明顯大于半年周期項對應(yīng)的振幅值，這就說明上海CORS站高程數(shù)據(jù)周期性變化主要為一年變化，半年變化占比非常少。

2.2 非線性最小二乘建模預(yù)測

2.2.1 非線性最小二乘建模

前文中已經(jīng)提到線性最小二乘建模方法中給定固定周期值的方法與實際CORS站高程數(shù)據(jù)周期性變化值存在一定的偏差，為此，本文建立了非線性最小二乘擬合模型(5)。該模型以線性最小二乘模型的解作為迭代初值，其中非線性最小二乘模型迭代初值如表2所示。

表2 非線性最小二乘建模迭代初值

將表2中給出的非線性最小二乘模型迭代初值代入非線性最小二乘模型(5)，進(jìn)行不斷地迭代解算得到上海站高程數(shù)據(jù)非線性最小二乘擬合建模結(jié)果如圖3所示：縱軸代表高程數(shù)據(jù)值，橫軸表示的是高程值對應(yīng)的時間。

圖3 上海CORS站高程坐標(biāo)非線性速度場擬合模型圖

從圖3上海CORS站高程坐標(biāo)非線性速度場擬合模型圖可以看出：擬合曲線能夠較好地描述CORS站高程數(shù)據(jù)的變化情況，即擬合曲線能夠反映出CORS站高程數(shù)據(jù)的周期性變化的規(guī)律。并且從圖3中也可以看出，擬合模型計算出來的值與實際觀測的高程值誤差在2 cm以內(nèi)。為了說明非線性速度場模型可以反映CORS站高程數(shù)據(jù)周期性變化的實際情況，本文計算出非線性最小二乘擬合模型中的未知參數(shù)值，如表3所示。從表3可以看出，非線性最小二乘模型中的未知參數(shù)f1，f2的數(shù)值并不是整數(shù)值，即周期性變化的數(shù)值并不是固定的整年或整半年數(shù)值，這就說明了非線性最小二乘模型將周期項作為未知數(shù)進(jìn)行迭代是符合CORS站高程數(shù)據(jù)周期性變化規(guī)律的。另外，通過對比非線性最小二乘擬合模型殘差均方根誤差與線性最小二乘擬合模型誤差發(fā)現(xiàn)，非線性最小二乘擬合模型相較于線性最小二乘擬合模型建模精度提高了32%。

表3 非線性最小二乘建模解算未知參數(shù)值

2.2.2 非線性最小二乘預(yù)測

為了驗證本文提出的非線性最小二乘速度場建模方法的可行性，作者在本小節(jié)以上海CORS站2005～2014年的高程數(shù)據(jù)利用非線性速度場模型進(jìn)行擬合建模，然后求解出非線性模型的未知參數(shù)，得到非線性速度場模型的數(shù)學(xué)表達(dá)式。將上海CORS站2015～2017年的高程數(shù)據(jù)值作為實際觀測值，并與非線性預(yù)測模型求解出的預(yù)測值進(jìn)行比較，得到的預(yù)測結(jié)果如圖4所示，前一段綠色的點代表的為2005～2014年的CORS站高程數(shù)據(jù)值，黑色的曲線代表的為利用非線性最小二乘方法對該段高程數(shù)據(jù)的擬合。后一段黑色的點代表的是待預(yù)測的2015～2017年的CORS站高程數(shù)據(jù)觀測值，紅色的曲線為非線性速度場模型預(yù)測值。相應(yīng)的非線性模型參數(shù)值和預(yù)測結(jié)果精度值見表4，其中RMSE/m為預(yù)測值于實際值誤差的均方根誤差統(tǒng)計值。從圖4的預(yù)測曲線可以看出，預(yù)測曲線數(shù)據(jù)值的變化與CORS站高程數(shù)據(jù)實際值的變化基本一致，且預(yù)測時間越短，預(yù)測的精度越高，但從表4可看出總體預(yù)測誤差在1.7 cm以內(nèi)，說明本文建立的非線性速度場模型具有一定可行性與正確性。

圖4 上海CORS站高程坐標(biāo)非線性速度場預(yù)測結(jié)果圖

表4 非線性最小二乘預(yù)測結(jié)果參數(shù)值

3 結(jié) 語

本文對上海CORS站2005～2017年高程數(shù)據(jù)進(jìn)行了建模方法的研究，首先分析對比了線性最小二乘建模和非線性最小二乘建模方法的原理，然后采用CORS站高程數(shù)據(jù)進(jìn)行了建模預(yù)測實驗的分析，得出以下結(jié)論。

(1)線性最小二乘建模采用固定周期項的方法與上海CORS站高程數(shù)據(jù)實際周期變化值有一定的誤差，導(dǎo)致了建立的數(shù)學(xué)模型不能恰當(dāng)?shù)胤从掣叱虜?shù)據(jù)實際的運動變化規(guī)律。

(2)上海CORS站高程數(shù)據(jù)線性最小二乘和非線性最小二乘擬合建模結(jié)果圖都反映出高程運動存在的周期性變化以一年變化為主，半年變化占比較少。

(3)以殘差均方根誤差為建模精度評價標(biāo)準(zhǔn)得到，非線性最小二乘擬合建模精度較線性最小二乘擬合建模精度提高了32%，且非線性模型預(yù)測誤差在1.7 cm以內(nèi)，均方根誤差為0.008 7 m，預(yù)測精度高。