蔣伊琳， 張芳園， 鄭　輝

(哈爾濱工程大學(xué) 信息與通信工程學(xué)院， 哈爾濱 150001)

Bootstrap方法對晶振穩(wěn)定度估計*

蔣伊琳， 張芳園， 鄭輝

(哈爾濱工程大學(xué) 信息與通信工程學(xué)院， 哈爾濱 150001)

針對小樣本條件下晶振能夠快速穩(wěn)定實現(xiàn)時鐘同步的問題，提出了在蒙特卡洛模擬實驗下結(jié)合傳統(tǒng)和改進(jìn)后Bootstrap方法與卡爾曼濾波方法消除相關(guān)噪聲，從而維持晶振穩(wěn)定性的一種方法.由于傳統(tǒng)Bootstrap方法限制再生樣本的生成范圍，使計算結(jié)果遠(yuǎn)遠(yuǎn)偏離真實情況，而改進(jìn)Bootstrap方法可以拓展再生樣本的生成范圍，避免了只抽取原始數(shù)據(jù)的情況，使得頻偏和時偏的估計分布盡可能地接近真實分布.仿真結(jié)果表明，改進(jìn)Bootstrap對頻偏和時偏的校正優(yōu)于傳統(tǒng)Bootstrap方法，利用改進(jìn)Bootstrap與卡爾曼相結(jié)合的方法能夠?qū)崿F(xiàn)對晶振的快速穩(wěn)定.

Bootstrap方法； 卡爾曼濾波； 時鐘同步； 蒙特卡洛； 小樣本； 再生樣本； 原始數(shù)據(jù)； 估計分布

高精度的時間同步技術(shù)在通信、國防和深空探索等領(lǐng)域是一項關(guān)鍵性技術(shù).運用該技術(shù)可以大大提高系統(tǒng)中時間的一致性和準(zhǔn)確性，從而提高系統(tǒng)的穩(wěn)定性和可靠性，滿足越來越高的工程應(yīng)用和科學(xué)研究的要求.

時鐘在運行過程中，在不同的運行環(huán)境、溫度、濕度和運行時間等條件以及晶振頻率本身偏差影響，晶振的穩(wěn)定度會逐步變差，時鐘時間也將各不相同.為了快速實現(xiàn)時鐘同步，對時鐘的頻率和偏差進(jìn)行同時校正，并對晶振進(jìn)行快速校正，縮短穩(wěn)定時間，通常采用卡爾曼濾波的方法對晶振的時偏和頻偏進(jìn)行估計，但卡爾曼濾波在運行過程中需要足夠多的原始數(shù)據(jù)維持系統(tǒng)的穩(wěn)定性.在小樣本條件下，卡爾曼濾波過程不能夠?qū)崿F(xiàn)對晶振的精確校正，因此，可通過Bootstrap方法獲取大量的數(shù)據(jù)，再結(jié)合卡爾曼濾波，改善晶振的穩(wěn)定度.

Bootstrap方法也被稱為自助法，是一種通過等概率有放回模擬再抽樣來分析不確定性的方法，適用于任何分布.在實際工程中，Bootstrap克服了數(shù)據(jù)有限這一困難，但在確定分布函數(shù)時只依賴于原始數(shù)據(jù)，從而自主分布產(chǎn)生偏離，不能接近真實情形.在這種情況下，通過擴(kuò)展生成原始樣本范圍改進(jìn)Bootstrap方法，避免了再抽樣時再生數(shù)據(jù)只為原始數(shù)據(jù)這一情況，使估計分布盡可能地接近真實分布，獲得高精度估計.

1　卡爾曼濾波

1.1時偏和頻偏估計

卡爾曼濾波[1]是在信號和噪聲的統(tǒng)計特性已知的條件下，采用高效率遞歸方式獲得最小方差的最優(yōu)估計濾波器.在IEEE1588協(xié)議[2]條件下，采用卡爾曼濾波估計[3]對晶振進(jìn)行校正，提高同步精度.

IEEE1588協(xié)議的高同步精度建立在硬件時間戳的基礎(chǔ)上，通過主從時鐘之間時間戳信息的交換而達(dá)到同步的目的，其中，假設(shè)時間數(shù)據(jù)包的往返時延對稱.圖1為主從時鐘的交互過程.協(xié)議把時間同步過程分為兩個階段：1)偏移測量階段；2)延遲測量階段.

在本文中假設(shè)往返路徑對稱，即數(shù)據(jù)包在傳播過程中的時延相等，則有

(1)

(2)

式中：offset為主從時鐘之間的時偏；TTrans1和TTrans2分別為從時鐘接收到主時鐘發(fā)送Sync數(shù)據(jù)包的傳播時延和主時鐘接收從時鐘發(fā)送Delay-Req數(shù)據(jù)包的傳播時延；Tm1為主時鐘發(fā)送Sync數(shù)據(jù)包的時間；Ts2為從時鐘發(fā)送Delay-Req數(shù)據(jù)

圖1　主從時鐘之間的交互過程

包的時間；Ts1為從時鐘接收到主時鐘發(fā)送的Sync數(shù)據(jù)包的時間；Tm2為主時鐘接收到從時鐘發(fā)送的Delay-Req數(shù)據(jù)包的時間.

將式(1)、(2)進(jìn)行離散化，且認(rèn)為時間偏差的測量誤差大小為0，則從時鐘的時間偏差和頻率偏差的關(guān)系式為

θ(k+1)=θ(k)-uθ(k)+(s(k)-

us(k))τ+wθ(k)

(3)

s(k+1)=s(k)-us(k)+ws(k)

(4)

式中：θ(k)和s(k)為k時刻時鐘的時偏和頻偏；uθ(k)和us(k)為k時刻主從時鐘實現(xiàn)時鐘同步時，時偏和頻偏的修正值；wθ(k)和ws(k)為k時刻時鐘的時偏和頻偏誤差；τ為k時刻采樣樣本的采樣時間，因采樣頻率一定，則采樣周期恒定.

設(shè)Z(k)=[θ(k)，s(k)]T為k時刻系統(tǒng)的當(dāng)前狀態(tài)向量，Y(k)=[uθ(k)，us(k)]T為k時刻時鐘時間和頻率的調(diào)整值.在同步間隔較短時認(rèn)為時鐘的頻偏變化平穩(wěn)，假定在同步間隔內(nèi)頻偏的變化率保持不變，即用上次同步的頻偏作為同步間隔內(nèi)的頻偏計算下一狀態(tài)的時間偏差.卡爾曼濾波的狀態(tài)方程和測量方程關(guān)系式為

Z(k)=AZ(k-1)+BY(k-1)+W(k-1)

(5)

Y(k)=HZ(k)+V(k)

(6)

式中：A和B為系統(tǒng)參數(shù)；H為觀測系統(tǒng)的參數(shù)；W(k-1)為k-1時刻的系統(tǒng)噪聲；V(k)為k時刻的測量噪聲.

根據(jù)卡爾曼濾波的遞推關(guān)系式，則有如下的求解過程.

狀態(tài)預(yù)測：

Zk/k-1=Ak/k-1Zk-1/k-1

(7)

狀態(tài)估計：

Zk/k=Zk/k-1+Kk(Yk-HZk/k-1)

(8)

濾波增益：

Kk=Pk/k-1HT(HPk/k-1HT+R)-1

(9)

均方誤差預(yù)測：

Pk/k-1=Ak/k-1Pk-1/k-1Ak/k-1+Q

(10)

均方誤差估計：

Pk/k=(I-KkH)Pk/k-1

(11)

式中：R為觀測噪聲方差；Q為狀態(tài)轉(zhuǎn)移協(xié)方差矩陣；I為單位矩陣；K為系統(tǒng)增益；P為狀態(tài)協(xié)方差矩陣.為了盡可能消除頻偏和時偏的相關(guān)噪聲，通過上述迭代計算過程可提高時鐘的同步精度.

1.2仿真結(jié)果

通過MATLAB仿真驗證上述過程，觀測時間為200 s，采樣間隔為0.1 s，觀測點為2 000個.本文中頻偏為晶振的實際頻率和標(biāo)準(zhǔn)頻率之間每秒產(chǎn)生的偏差值，通過上述求解過程對系統(tǒng)的頻偏和時偏進(jìn)行校正，提高本地時鐘的穩(wěn)定性和準(zhǔn)確性，仿真圖如圖2所示.

圖2　時鐘的時偏和頻偏的狀態(tài)估計值

由圖2可知，在觀測點為2 000個的條件下，卡爾曼濾波在約50 s時，頻偏和時偏狀態(tài)估計值趨于穩(wěn)定，時偏估計精度為10-5，頻偏估計精度為10-7.

2　Bootstrap方法

2.1傳統(tǒng)Bootstrap方法

Bootstrap方法[4-8]的核心思想是用已知的經(jīng)驗分布函數(shù)代替未知分布函數(shù)，即利用從原始樣本中抽取的再生樣本計算相關(guān)統(tǒng)計量，從而來估計總體分布特性.Bootstrap再抽樣有3種方法，分別為參數(shù)法、半?yún)?shù)法和非參數(shù)法.

已知一組隨機(jī)觀測樣本X=(x1，x2，…，xN)，N為總觀測點的個數(shù)，F(xiàn)為未知分布函數(shù)，R(X，F(xiàn))是X、F的函數(shù)，根據(jù)原始樣本X估計隨機(jī)變量R(X，F(xiàn))的分布特征.

Bootstrap方法基本步驟如下：

1) 根據(jù)已知的觀測樣本X構(gòu)造其經(jīng)驗分布函數(shù)Fn，每個觀測樣本被抽取的概率為1/N；

3) 根據(jù)得到的Bootstrap子樣本計算相應(yīng)的統(tǒng)計量Rb=Rb(Xb，F(xiàn)n)；

4) 重復(fù)步驟2)、3)過程共n次(n≥1 000)，獲得相應(yīng)的n個統(tǒng)計量Rb；

5) 利用Rb分布去逼近R(X，F(xiàn))分布，Rb的分布被稱為Bootstrap分布.

2.2改進(jìn)Bootstrap方法

Bootstrap方法作為一種統(tǒng)計方法被應(yīng)用于許多場合，但是其也有許多的不足之處.由于再生樣本是由原始樣本直接生成，當(dāng)原始樣本數(shù)據(jù)較小時，再生樣本會非常相似于原始樣本，將會導(dǎo)致計算結(jié)果遠(yuǎn)遠(yuǎn)偏離真實情況，并且當(dāng)真實分布為連續(xù)分布時，無法獲取未在樣本觀測點處的分布特性，從而不能估計出更加接近于真實分布的分布函數(shù).

針對上述不足，對Bootstrap方法做出了如下改進(jìn)[9-10]：

1) 將原始觀測樣本X=(x1，x2，…，xN)按從小到大的順序排列，結(jié)果仍記為X=(x1，x2，…，xN).

2) 對變化后的樣本做出相應(yīng)的鄰域，即A1=[x1-(x2-x1)/m，x1+(x2-x1)/m]；Ai=[xi-(xi-xi-1)/m，xi+(xi+1-xi)/m]；AN=[xN-(xN-xN-1)/m，xN+(xN-xN-1)/m]，i=2，3，…，N-1.m≥2，通常取m=2，當(dāng)樣本容量過小時，即N<10時，為避免再生樣本過大偏離真實分布，m的取值應(yīng)大一點.

3) 確定L為鄰域指標(biāo)，其分布特性為P(L=i)=1/N，i=1，2，…，N.根據(jù)該分布確定鄰域指標(biāo)樣本(l1，l2，…，lN).

改進(jìn)后的Bootstrap方法擴(kuò)展了觀測數(shù)據(jù)，避免只抽取原始數(shù)據(jù)這一情況，使估計分布盡可能地接近真實分布，提高了估計的準(zhǔn)確性.

2.3基于Bootstrap的卡爾曼對時偏和頻偏的估計

蒙特卡洛方法[11-12]也被稱為統(tǒng)計模擬方法，被應(yīng)用在許多計算領(lǐng)域且有廣泛的影響.它利用隨機(jī)數(shù)構(gòu)建模型對各種因素進(jìn)行評估，通過設(shè)定的隨機(jī)過程，反復(fù)生成隨機(jī)序列，計算參數(shù)估計量和統(tǒng)計量以及其分布特征.理論上，通過計算機(jī)的高速率和大容量計算，蒙特卡洛模擬方法隨著增加的模擬次數(shù)可以獲得任何精度水平的計算，并且計算誤差不會隨著計算維度的增加而增加.

在小樣本(N=10)條件下，分別采用傳統(tǒng)Bootstrap方法和改進(jìn)Bootstrap方法對原始樣本進(jìn)行重采樣，獲得Bootstrap樣本，文中假定重采樣的次數(shù)為2 000次.將Bootstrap樣本作為新的觀測數(shù)據(jù)，采樣間隔為0.01 s，卡爾曼濾波器的參數(shù)設(shè)置與上文相同.將Bootstrap樣本送入卡爾曼濾波器，消除噪聲影響，對晶振的時偏和頻偏進(jìn)行校正.在蒙特卡洛模擬實驗條件下重復(fù)上述過程100次，得到傳統(tǒng)Bootstrap和改進(jìn)Bootstrap方法對晶振進(jìn)行校正后時偏和頻偏估計精度的集中范圍.采用MATLAB仿真，得到傳統(tǒng)和改進(jìn)后的Bootstrap仿真圖，如圖3、4所示.

圖3　傳統(tǒng)Bootstrap采樣后時鐘時偏和頻偏估計

圖4　改進(jìn)Bootstrap采樣后時鐘時偏和頻偏估計

由圖3、4可知，經(jīng)過Bootstrap 2 000次重采樣，蒙特卡洛100次模擬實驗后，傳統(tǒng)Bootstrap方法和改進(jìn)Bootstrap方法的頻偏和時偏狀態(tài)估計值都在約10 s時趨于穩(wěn)定，而傳統(tǒng)Bootstrap方法的時偏估計精度和頻偏估計精度分別為10-4和10-6，改進(jìn)后的Bootstrap方法的時偏和頻偏估計精度少數(shù)為10-4和10-6，時偏和頻偏估計精度絕大多數(shù)集中于10-5和10-7.

3　結(jié)　論

通過驗證，改進(jìn)后的Bootstrap方法雖然擴(kuò)展了樣本范圍，使采樣樣本不再局限于原始數(shù)據(jù)，盡可能地使估計分布接近真實分布，提高了估計的準(zhǔn)確性，但是其改進(jìn)程度有限，不能使其分布與真實分布完全一致.經(jīng)過蒙特卡洛模擬實驗后，頻偏和時偏的估計精度絕大多數(shù)集中于10-5和10-7，少數(shù)的估計精度為10-4和10-6.在小樣本條件下，改進(jìn)后的Bootstrap與卡爾曼相結(jié)合的方法在校正精度上與大樣本時的卡爾曼濾波效果接近一致，且在相等的運行時間，改進(jìn)Bootstrap方法能夠更加快速達(dá)到穩(wěn)定，縮短穩(wěn)定時間，同時小樣本所需要的采樣時間更短.通過仿真分析可知，在觀測數(shù)據(jù)為小樣本時，改進(jìn)Bootstrap方法與卡爾曼相結(jié)合的方法能夠快速實現(xiàn)晶振的頻偏和時偏校正，在實際應(yīng)用中具有重要意義.

[1]王慶輝，魯婷婷，劉孟龍，等.基于TOA減小非視距誤差的方案設(shè)計 [J].沈陽工業(yè)大學(xué)學(xué)報，2014，36(2)：204-209.

(WANG Qing-hui，LU Ting-ting，LIU Meng-long，et al.Scheme design of TOA-based system with mitigating NLOS error [J].Journal of Shenyang University of Technology，2014，36(2)：204-209.)

[2]李曉珍.基于IEEE1588的網(wǎng)絡(luò)時間同步系統(tǒng)研究 [D].西安：中國科學(xué)院研究生院，2011.

(LI Xiao-zhen.Research on the network time synchronization system based on IEEE1588 [D].Xi’an：Graduate School of Chinese Academy of Sciences，2011.)

[3]周思捷.基于IEEE1588無線網(wǎng)絡(luò)時間同步技術(shù)研究 [D].上海：上海交通大學(xué)，2013.

(ZHOU Si-jie.The technology research on time synchronization in wireless network based on IEEE1588 [D].Shanghai：Shanghai Jiaotong University，2013.)

[4]吳慶平，林素仙，黃飛.非參數(shù)Bootstrap方法及其MATLAB實現(xiàn) [J].麗水學(xué)院學(xué)報，2012，34(2)：14-18.

(WU Qing-ping，LIN Su-xian，HUANG Fei.Non-parametric Bootstrap method and its realiztion in MATLAB [J].Journal of Lishui University，2012，34(2)：14-18.)

[5]袁修開，呂震宙，岳珠峰.小樣本下分位數(shù)函數(shù)的Bootstrap置信區(qū)間估計 [J].航空學(xué)報，2012，33(10)：1842-1849.

(YUAN Xiu-kai，Lü Zhen-zhou，YUE Zhu-feng.Bootstrap confidence interval of quantile function estimation for small samples [J].Acta Aeronautica et Astronutica Sinica，2012，33(10)：1842-1849.)

[6]張萍.基于Bootstrap方法的統(tǒng)計分析 [J].宜賓學(xué)院學(xué)報，2011，11(12)：31-33.

(ZHANG Ping.Statistical analysis based on Bootstrap method [J].Journal of Yibin University，2011，11(12)：31-33.)

[7]孫國童，唐碩.仿真試驗中Bootstrap方法的應(yīng)用 [J].系統(tǒng)仿真學(xué)報，2010，22(6)：1347-1349.

(SUN Guo-tong，TANG Shuo.Application of Bootstrap estimators in simulation experiment [J].Journal of System Simulation，2010，22(6)：1347-1349.)

[8]韓冬，馬進(jìn)，賀仁睦.基于Bootstrap 的負(fù)荷模型的小樣本不確定性分析 [J].電力系統(tǒng)保護(hù)與控制，2012，40(18)：95-100.

(HAN Dong，MA Jin，HE Ren-mu.Small sample uncertainty analysis of load model based on Bootstrap [J].Power System Protection and Control，2012，40(18)：95-100.)

[9]劉建，吳翊，譚璐.對 Bootstrap 方法的自助抽樣的改進(jìn) [J].數(shù)學(xué)理論與應(yīng)用，2006，26(1)：69-72.

(LIU Jian，WU Yi，TAN Lu.An improvement to the resampling procedure of Bootstrap method [J].Ma-thematical Theory and Applications，2006，26(1)：69-72.)

[10]鄒艷，羅文強(qiáng).改進(jìn)的 Bootstrap 方法對比及應(yīng)用研究 [J].應(yīng)用數(shù)學(xué)，2008，21(增刊1)：62-66.

(ZOU Yan，LUO Wen-qiang.Comparative study anel application of improved Bootstrap method [J].Ma-thematica Applicata，2008，21(Sup1)：62-66.)

[11]邵偉.蒙特卡洛方法及在一些統(tǒng)計模型中的應(yīng)用 [D].濟(jì)南：山東大學(xué)，2012.

(SHAO Wei.Monte Carlo methods and their applications in some statistical model [D].Jinan：Shandong University，2012.)

[12]黃殿勛，張文，郭萍，等.發(fā)輸電系統(tǒng)可靠性評估的蒙特卡洛改進(jìn)算法 [J].電力系統(tǒng)保護(hù)與控制，2010，38(21)：179-183.

(HUANG Dian-xun，ZHANG Wen，GUO Ping，et al.The Monte Carlo improved method for reliability eva-luation of generation and transmission systems [J].Power System Protection and Control，2010，38(21)：179-183.)

(責(zé)任編輯：鐘媛英文審校：尹淑英)

Estimation of crystal oscillator stability with Bootstrap method

JIANG Yi-lin, ZHANG Fang-yuan, ZHENG Hui

(College of Information and Communication Engineering, Harbin Engineering University, Harbin 150001, China)

In order to solve the problem that the crystal oscillator can realize the fast stability and clock synchronization under the small sample condition, a method which could eliminate the correlated noise to maintain the crystal oscillator stability through combining with the traditional Bootstrap method, improved Bootstrap method and Kalman filtering method under the condition of Monte Carlo simulation experiment, was proposed. Because the traditional Bootstrap method limits the formation range of the regeneration sample, the calculated results largely deviate from the true situation. However, the improved Bootstrap method can expand the formation range of the regeneration sample, and such situation that only the original data are extracted can be avoided, which makes the estimation distribution of frequency offset and time offset be close to the real distribution. The simulated results show that the improved Bootstrap method is superior to the traditional Bootstrap method in the calibration of frequency offset and time offset. The method combined the improved Bootstrap method with the Kalman filtering method can realize the fast stability of crystal oscillator.

Bootstrap method; Kalman filtering; clock synchronization; Monte Carlo; small sample; regeneration sample; original data; estimation distribution

2015-08-31.

國家自然科學(xué)基金資助項目(61202410).

蔣伊琳(1980-)，男，浙江鎮(zhèn)江人，講師，博士，主要從事寬帶信號檢測與估值、無線時間同步等方面的研究.

10.7688/j.issn.1000-1646.2016.03.12

TM 343

A

1000-1646(2016)03-0304-05

*本文已于2015-12-07 16∶18在中國知網(wǎng)優(yōu)先數(shù)字出版. 網(wǎng)絡(luò)出版地址： http：∥www.cnki.net/kcms/detail/21.1189.T.20151207.1618.040.html