葉錦強(qiáng)

（1.惠州市國(guó)土資源局 大亞灣經(jīng)濟(jì)技術(shù)開(kāi)發(fā)區(qū)分局，廣東 惠州 516081）

GPS已被廣泛運(yùn)用到地殼形變的研究上，而將GPS資料制成時(shí)間序列后，可以分析GPS的連續(xù)觀測(cè)資料得到測(cè)站點(diǎn)與時(shí)間及空間的關(guān)系，進(jìn)而推估出各種構(gòu)造與地殼形變之間的關(guān)聯(lián)性。

GPS的觀測(cè)數(shù)據(jù)中存在噪聲。數(shù)據(jù)噪聲可以簡(jiǎn)單地分為2種：與時(shí)間無(wú)關(guān)的白噪聲和與時(shí)間相關(guān)的有色噪聲。白噪聲是一種隨機(jī)的信號(hào)，在功率頻譜上呈現(xiàn)水平的功率譜密度。白噪聲在光功率頻譜內(nèi)位于可見(jiàn)光波段，當(dāng)眼睛接受光的三原色分量相等時(shí)為白光。白噪聲之外的噪聲統(tǒng)稱為有色噪聲。

在多數(shù)的GPS時(shí)間序列的研究中，僅將觀測(cè)資料的誤差視為與觀測(cè)時(shí)間無(wú)關(guān)的白噪聲，因?yàn)榘自肼暤臄?shù)值模型和計(jì)算比較容易。但大量研究表明，如果GPS時(shí)間序列分析中不考慮與時(shí)間相關(guān)的有色噪聲，將會(huì)低估地殼形變速度的誤差[1-3]。尤其是在一些地殼活動(dòng)速度比較低的地區(qū)測(cè)得的GPS位移量很小，若低估誤差可能導(dǎo)致地殼變形的大小及方向上的錯(cuò)誤，這時(shí)候考慮有色噪聲的影響就格外重要。本文中將有色噪聲分為閃爍噪聲和隨機(jī)游走噪聲。前者可能來(lái)自GPS系統(tǒng)中未被模型化的大氣效應(yīng)、衛(wèi)星軌道及星歷誤差、多路徑效應(yīng)、衛(wèi)星海拔高度及方位角所造成的效應(yīng)[3-5]；后者來(lái)自于地表潛移或降雨使儀器支柱產(chǎn)生布朗運(yùn)動(dòng)所造成的效應(yīng)[6]。

本文利用時(shí)間序列分析軟件CATS[7]，基于兩種不同的方法來(lái)進(jìn)行資料的噪聲分析，分別是時(shí)間域上的最大似然估計(jì)法和頻率域上的頻譜分析法。

1 GPS資料處理與噪聲分析方法

對(duì)于GPS的時(shí)間序列的分析，通?？紤]長(zhǎng)周期運(yùn)動(dòng)、同震變形以及震后變形。GPS單站、單分量坐標(biāo)序列一般用下列周期性模型來(lái)表示[8]：

式中，ti以年為單位；系數(shù)a、b分別代表地殼位置和線性變換率；系數(shù)c、d、e和f描述了周年和半周年運(yùn)動(dòng)振幅；系數(shù)gj代表地震造成的同震偏移（offset）；H（t）為階躍函數(shù)（heaviside step function）；vi為殘差值，代表觀測(cè)值與預(yù)測(cè)值間的差異，也就是本文需要分析的噪聲信號(hào)。

由于研究目的是與地殼形變無(wú)關(guān)的噪聲信息，所以針對(duì)GPS時(shí)間序列都必須按照式（1）事先扣除已知的信號(hào)，包括變形速率、各種原因造成的偏移量、季節(jié)性的周期變化（在此為年周期和半年周期），研究區(qū)站點(diǎn)分布見(jiàn)圖1。

圖1 研究地區(qū)站點(diǎn)分布圖

圖2 SUIY測(cè)站時(shí)間序列處理事例

圖2a中從上至下分別是南北分量、東西分量、垂直分量的原始時(shí)間序列（紅色線）及其擬合線（綠色線），圖2b中從上至下分別是南北分量、東西分量、垂直分量的殘差時(shí)間序列，也就是本文所要研究的噪聲時(shí)間序列。

1.1 頻譜分析

GPS原始時(shí)間序列扣除擬合值后的噪聲時(shí)間序列，可以用頻譜定性分析噪聲的類型。頻譜分析是一種在頻率域上分析信號(hào)的方法。噪聲時(shí)間序列在頻譜域上的功率譜可以使用Power Law的形式表示[6]：

式中，α為頻譜指數(shù)；P0為常數(shù)。α越大，表示噪聲序列的時(shí)間相關(guān)性越高。

將式（2）取對(duì)數(shù)進(jìn)行曲線擬合，運(yùn)用最小二乘法可以得到P0及α值。α通常是介于-1～3之間的任意實(shí)數(shù)，其中整數(shù)α代表一些特殊的噪聲類型：α=0時(shí)，是標(biāo)準(zhǔn)的白噪聲；α=1時(shí)，是標(biāo)準(zhǔn)的閃爍噪聲；α=2時(shí)，則是標(biāo)準(zhǔn)的隨機(jī)游走噪聲。另外，除了標(biāo)準(zhǔn)的白噪聲以外，其余的部分統(tǒng)稱為有色噪聲。

1.2 最大似然估計(jì)法

最大似然法是一種非線性的最小二乘法，目的在于找到與時(shí)間序列最相近的模型參數(shù)。通過(guò)調(diào)整協(xié)方差矩陣使得似然函數(shù)取得最大值，即可得到與該時(shí)間序列最相近的噪聲模型，這樣就可以通過(guò)協(xié)方差矩陣估算出時(shí)間序列中的噪聲振幅大小。與頻譜分析相比，最大似然估計(jì)法可以定量地計(jì)算噪聲大小。

最大似然估計(jì)法的方程式如下[1]：

或

式中，ln是自然對(duì)數(shù)；lik（x，C）是概略值；n是資料點(diǎn)數(shù)；det是行列式；C是資料的協(xié)方差矩陣。確定協(xié)方差矩陣的方法可以參考文獻(xiàn) [1]、[2]、[9]、[10]。

2 時(shí)間序列分析軟件CATS

CATS軟件包含2個(gè)程序：線性部分使用最小二乘，利用式（1）計(jì)算截距、斜率、偏移以及周期信號(hào)的振幅；非線性部分則利用線性計(jì)算后的殘差求得特定噪聲模型的參數(shù)及振幅大小。

假設(shè)模型中只有一種噪聲，則模型的協(xié)方差矩陣C=σ2J（σ表示噪聲振幅大?。?，則式（4）可以改寫(xiě)為：

此時(shí)，可以由式（5）得到噪聲振幅的關(guān)系式，此關(guān)系式可以令式（5）達(dá)到最大值：

當(dāng)模型中不只一種噪聲模型時(shí)，CATS用下面的方法來(lái)計(jì)算最大似然估計(jì)值。設(shè)2個(gè)噪聲振幅分別σ1和σ2，將振幅以角度φ和一個(gè)比例常數(shù)r表示：

而協(xié)方差矩陣可以表示成：

代入式（4），就可以簡(jiǎn)單地運(yùn)用角度φ計(jì)算最大似然估計(jì)法的數(shù)值。

3 最佳噪聲模型分析

本文使用3個(gè)噪聲模型來(lái)分析噪聲時(shí)間序列：白噪聲（WN）、白噪聲+閃爍噪聲（FL）、白噪聲+隨機(jī)游走噪聲（RW），利用CATS軟件，使用最大似然估計(jì)法在時(shí)間域上求最佳噪聲模型，并使用頻譜分析法從頻率域上找出最佳模型。

表1中針對(duì)3種模型計(jì)算出各測(cè)站噪聲三分量的MLE值。從表中可以統(tǒng)計(jì)出，16個(gè)GPS測(cè)站的48個(gè)分量中有90%以上的FL有較大的MLE值（相比RW和WN的MLE值）。只有極少分量（JLYJ的E分量和SUIY的N、E分量）的RW有較大MLE值，這可能是因?yàn)橛^測(cè)或解算中存在誤差，也不排除與該站的基座有關(guān)。雖然存在上述異常情況，但整體上可以說(shuō)明，白噪聲+閃爍噪聲是描述研究地區(qū)GPS觀測(cè)網(wǎng)噪聲時(shí)間序列的較好模型。

表1 GPS觀測(cè)網(wǎng)中各測(cè)站的MLE值

另外需要注意的是，“陸態(tài)網(wǎng)絡(luò)”的GPS儀器基座為水泥柱，屬于深基及淺基支柱型的較穩(wěn)定的儀器基座，相對(duì)來(lái)說(shuō)不容易產(chǎn)生隨機(jī)游走噪聲，這與計(jì)算的結(jié)果也是一致的。

利用CATS軟件計(jì)算出各測(cè)站三分量的有色噪聲頻譜指數(shù)，將頻譜分析的模型設(shè)定為噪聲加上一個(gè)Power Law的有色噪聲，計(jì)算結(jié)果如表2。從表中可以看出，研究地區(qū)的GPS資料的頻譜指數(shù)基本介于0～-1之間，與前人的研究結(jié)果相同[2,3]，與最大似然估計(jì)法的計(jì)算結(jié)果也相符。

表2 觀測(cè)網(wǎng)中各測(cè)站的頻譜指數(shù)

4 結(jié) 語(yǔ)

本文對(duì)GPS觀測(cè)數(shù)據(jù)進(jìn)行噪聲分析。基于最大似然估計(jì)法計(jì)算結(jié)果表明，白噪聲+閃爍噪聲是較近似的噪聲模型。少數(shù)幾個(gè)噪聲接近白噪聲+隨機(jī)游走噪聲模型的測(cè)站，可能是由于儀器基座不穩(wěn)定或者時(shí)間序列長(zhǎng)度較長(zhǎng)所產(chǎn)生的效應(yīng)。而頻譜分析所得到的結(jié)果與最大似然估計(jì)法一致。不同的噪聲模型將對(duì)地殼形變速度的誤差造成很大的影響，而研究更精確的模型對(duì)于反演地球自身的運(yùn)動(dòng)狀態(tài)尤為重要。

[1]Langbein J,Johnson H.Correlated Errors in Geodetic Timeseries:Implications for Time-dependent Deformation[J].Geophys Res,1997,102(1):591-603

[2]Zhang J,Bock Y,Johnson H,et al.Southern California Permanent GPS Geodetic Array:Error Analysis of Daily Position Estimates and Site Velocities[J].Geophys Res,1997,102(8): 35-55

[3]Mao A,Harrison C G A,Dixon T H.Noise in GPS Coordinate Time Series[J].Geophys Res, 1999,104(1):2 797-2 816

[4]Williams S,Bock Y,Fang P,et al.Error Analysis of Continuous GPS Position Time Series[J].Geophys Res,2004,109(3):412-421

[5]Langbein J.Noise in Two-color Electronic Distance Meter Measurements Revisited[J].Geophys Res,2004,109(4):1-16

[6]Agnew D C.The Time-domain Behavior of Power Law Noises[J].Geophysical Research Letters,1992,19(4):333-336

[7]Williams S D P.The Effect of Coloured Noise on the Uncertainties of Rates Estimated from Geodetic Time Series[J].Geodesy, 2003,76(10):483-494

[8]Nikolaidis R.Observation of Geodetic and Deismic Deformation with the Global Positioning System[D].University of California,2002

[9]黃立人.GPS基準(zhǔn)站坐標(biāo)分量時(shí)間序列的噪聲特性分析[J].大地測(cè)量與地球動(dòng)力學(xué),2006,26(2):31-34

[10]Gardner M.White and Brown Music, Fractal Curves and Oneover-fluctuations[J].Scientific American, 1978,238(4):16-32