李鑫冉 王歆

(1中國科學(xué)院紫金山天文臺 南京 210008) (2中國科學(xué)院空間目標(biāo)與碎片觀測重點實驗室 南京 210008) (3中國科學(xué)院大學(xué) 北京 100049)

基于遺傳算法的極短弧定軌(續(xù))?

李鑫冉1,2,3?王歆1,2?

(1中國科學(xué)院紫金山天文臺 南京 210008) (2中國科學(xué)院空間目標(biāo)與碎片觀測重點實驗室 南京 210008) (3中國科學(xué)院大學(xué) 北京 100049)

采用遺傳算法解決極短弧定軌問題時,由于遺傳算法不同于經(jīng)典方法的計算過程,野值剔除方法不再適用.在遺傳算法中通過在適值函數(shù)中采用不同損失函數(shù)實現(xiàn)了穩(wěn)健估計,解決了極短弧定軌中的野值處理問題.在遺傳算法中不同損失函數(shù)的引入較經(jīng)典方法大大簡化.通過對多種損失函數(shù)的計算比較,表明采用最小中值二乘(LMS,Least M edian Square)和截尾最小二乘(LTS,Least Trimm ed Square)估計可大幅度提高極短弧定軌的穩(wěn)健性,具有極高的崩潰點.

航天器,天體力學(xué),方法:數(shù)值,統(tǒng)計

1 引言

近年來空間目標(biāo)的巡天觀測采集了大量數(shù)據(jù),但不同于傳統(tǒng)的跟蹤觀測,巡天方式采集到的弧段都十分短,往往只有數(shù)十秒.為了有效利用這些數(shù)據(jù),極短弧(TSA,Too Short A rc)定軌越來越受到關(guān)注[1].

由于短弧定軌問題本質(zhì)上是病態(tài)的,極短弧的病態(tài)性更為突出,傳統(tǒng)方法幾乎無法獲得合理的結(jié)果.文獻[2]將遺傳算法用于解決極短弧定軌問題,采用遺傳算法將定軌問題轉(zhuǎn)換為優(yōu)選問題,避免了求解線性代數(shù)方程,從而克服了經(jīng)典方法的一些不足.遺傳算法也不同于經(jīng)典優(yōu)選算法,對初值沒有特殊要求,不容易陷入局部極值,因此是解決極短弧定軌的一種有效手段.文獻[3]采用完全不同的變量和遺傳算子,也給出了一種遺傳算法求解極短弧定軌的具體方法.

上述研究給出了軌道計算方法,但都沒有涉及野值剔除的問題.在實測資料的初軌計算中,野值問題是不可回避的,觀測數(shù)據(jù)中的少量野值會對計算結(jié)果產(chǎn)生很大影響,甚至是破壞性的[4].最常見的野值處理的方法是采用3σ準(zhǔn)則,在迭代過程中對于大偏差的觀測資料進行剔除,使其不參與軌道計算,逐步收斂后得到可靠的結(jié)果.遺傳算法不同于經(jīng)典計算過程,無法采用這樣的野值剔除過程.同時由于野值的判斷與剔除需要在軌道計算過程中逐步確定,而剔除本身影響軌道計算結(jié)果,因此野值剔除并不容易.

野值處理的另一種手段是采用穩(wěn)健估計方法,文獻[5]在軌道計算中采用了M估計,但其求解方法比較復(fù)雜,特別是需要可靠初值,可靠初值獲取本身就是一個難題,受此局限方法的崩潰點不高.此外方法中引入了額外的參數(shù),參數(shù)確定尚無可靠手段.文獻[5]中根據(jù)歷史數(shù)據(jù)的統(tǒng)計特性給出參數(shù),然而對于極短弧數(shù)據(jù)而言,由于數(shù)據(jù)量少,這種做法并不合適.

文獻[6-7]采用了最小一乘(LAD,Least Absolute Deviation)估計,由于最小一乘估計無需額外參數(shù),避免了參數(shù)確定問題,但其求解采用了線性規(guī)劃方法,該方法和一般計算過程完全不同,僅適用于最小一乘這一特例情況.

在初軌計算的實踐中,野值剔除是不可或缺的環(huán)節(jié).本文在文獻[3]的工作基礎(chǔ)上對遺傳算法中野值的處理進行了討論,根據(jù)遺傳算法的特點采用了穩(wěn)健估計方法來克服野值的影響,并對幾種損失函數(shù)進行了分析,建立了基于遺傳算法的穩(wěn)健初軌計算方法.

2 初軌計算的遺傳算法

遺傳算法是通過模擬自然進化優(yōu)勝劣汰的過程來搜索問題的最優(yōu)解,算法將問題的解向量編碼為染色體,通過遺傳操作實現(xiàn)染色體的進化,從而找到全局最優(yōu)解[8].

為了后續(xù)闡述方便,簡要回顧一下文獻[3]的方法.其基本思路是將定軌問題轉(zhuǎn)換為分層優(yōu)化問題,先優(yōu)選(a,e,M0)3個Kep ler根數(shù),然后保持(a,e,M0)不變,再優(yōu)選(i,?,ω) 3個Kep ler根數(shù).與文獻[2]選擇了觀測首末時刻的斜距(ρ1,ρ2)作為優(yōu)選變量不同,采用軌道根數(shù)作為優(yōu)選變量解除了和觀測資料的耦合,便于資料處理,不受特定資料誤差的影響.采用(ρ1,ρ2)作為優(yōu)選變量時,軌道量無法通過優(yōu)選變量直接得到,需要結(jié)合對應(yīng)時刻的觀測量才能得到,對應(yīng)的觀測資料則尤為重要,其誤差將直接影響計算結(jié)果,當(dāng)其為野值時,解將會完全失敗并且無法克服.

與優(yōu)選法相同,采用遺傳算法初軌計算主要是要構(gòu)造對比量,將對比量采用優(yōu)選變量以及優(yōu)選量結(jié)合觀測量分別構(gòu)造.文獻[3]選擇了兩個觀測時刻的真近點角f的差作為對比量,令僅采用優(yōu)選變量構(gòu)造的f的差記為Δfc,而采用優(yōu)選量結(jié)合觀測量構(gòu)造的記為Δfo,兩者的差為殘差?=Δfo?Δfc.對于n個觀測資料共可形成N=C2n個對比量,其殘差記為?i,i=1,2,···,N.軌道計算問題即求下列適值函數(shù)最小的解:

計算中按照適值函數(shù)對每個染色體進行適應(yīng)度評估,然后通過選擇、交叉、變異的遺傳操作,使染色體不斷進化從而獲得適應(yīng)度更好的新種群,最終獲得最優(yōu)解.

不同于經(jīng)典方法逐步逼近的過程,遺傳算法是基于種群的計算方法,在進化(迭代)過程中遺傳操作對原種群是破壞的,3σ剔除方法將前一次迭代結(jié)果作為后一次迭代的參考,對于遺傳算法顯然是不適用的.從遺傳算法計算過程中可看出,只需要調(diào)整(1)式的適值函數(shù),無需對其他計算環(huán)節(jié)進行任何改變即可實現(xiàn)不同適值函數(shù)的初軌計算.較經(jīng)典計算框架而言十分簡便.在遺傳算法中采用穩(wěn)健估計方法處理野值問題是非常合適的.

3 損失函數(shù)

3.1 最小二乘(OLS,O rd inary Least Square)

在遺傳算法中,(1)式被稱為評價函數(shù)或適值函數(shù),是殘差的函數(shù),因此在估計理論中一般稱為損失函數(shù).不同估計方法主要區(qū)別在于不同的損失函數(shù)形式.(1)式即為最常用的OLS估計的損失函數(shù).OLS方法自Gauss提出以來,已成為數(shù)據(jù)處理中使用最為廣泛的方法,甚至成為了一種事實上的標(biāo)準(zhǔn).

OLS估計計算方法簡便,在誤差正態(tài)分布且沒有野值的情況下,解具有最優(yōu)的性質(zhì).但當(dāng)誤差偏離正態(tài)分布以及有野值時,OLS方法是不穩(wěn)健的,主要原因在于以殘差平方和作為損失函數(shù)時,對于大的殘差在平方作用下將會放大,解會更傾向于大殘差的觀測數(shù)據(jù),從而會使解徹底偏離.由于經(jīng)典初軌計算中都采用OLS估計,因此必須對野值進行處理,不然在實踐中幾乎無用.

3.2 最小一乘(LAD)

為了避免平方對殘差放大這一影響,改用LAD方法成為了自然的想法.文獻[6-7]采用LAD估計實現(xiàn)了穩(wěn)健初軌計算,LAD估計的損失函數(shù)為:

LAD方法的出現(xiàn)早于OLS方法,但由于其求解比較困難,一直沒有得到重視.直到20世紀(jì)中期由于其穩(wěn)健性才重新被人們重視.理想情況下LAD解的統(tǒng)計性質(zhì)原不及OLS解,但即使在少量野值存在時,LAD也不易受到野值影響,此時OLS解的優(yōu)良統(tǒng)計性質(zhì)也不再有任何優(yōu)勢[9].

3.3 最小中值二乘(LM S)

LAD估計的穩(wěn)健性是相對于觀測矩陣的,對于設(shè)計矩陣并不穩(wěn)健,因此并不是完全穩(wěn)健的,此外LAD估計的解不唯一.針對這些缺陷,文獻[10]提出了LMS估計,結(jié)合LAD估計和OLS估計各自的優(yōu)點,其損失函數(shù)為:

其中medi表示對所有i取其中位數(shù).顯然,采用該損失函數(shù)時解是唯一的,同時個別大殘差對損失函數(shù)并無貢獻,LMS估計具有理論最高崩潰點,接近50%.LMS解也具有十分明確的幾何意義,即覆蓋一半觀測資料的最窄帶的中線,LMS估計直觀上十分適合軌道計算.

3.4 截尾最小二乘(LTS)

從穩(wěn)健性角度來看,LMS估計具有極好的特性,但其漸近效率較低,這對于定軌問題而言并不重要.文獻[10]在LMS估計基礎(chǔ)上,提出了LTS估計,等價于損失函數(shù):其中?2r:N表示對于N個殘差?2i按照從小到大排序后第r個值,[N/2]表示不大于N/2的最大整數(shù).可見LTS估計與OLS估計十分類似,區(qū)別在于僅考慮最小的[N/2]+1個殘差. LTS估計保留了OLS估計的漸近效率,同時和LMS估計具有相同的高崩潰點.

除OLS估計外,上述3個估計都是穩(wěn)健估計,但與常用的Huber估計以及Hample估計不同,這3個估計準(zhǔn)則中都不需要額外的參數(shù).

4 數(shù)值驗證

根據(jù)前面的分析,在文獻[3]的基礎(chǔ)上,只需對適值函數(shù)做很小的變化即可實現(xiàn)上述損失函數(shù).采用和文獻[3]相同的中國科學(xué)院空間目標(biāo)與碎片光學(xué)觀測網(wǎng)中的一圈實測資料進行了數(shù)值驗證,資料精度優(yōu)于5′′,采樣頻率1 Hz.選取弧段開始的不同弧長進行了計算.歷元選擇為第1個觀測時刻,算法中參數(shù)和文獻[3]完全相同.

首先為了考察不同損失函數(shù)在無野值情況下的計算結(jié)果,對于60 s的弧段分別采用不同損失函數(shù)進行了計算,結(jié)果見表1.計算中采用了缺省的隨機數(shù)種子,其中POD表示由多天多站精密定軌得到的結(jié)果,作為參考標(biāo)準(zhǔn),OLS為最小二乘法的結(jié)果,和文獻[3]中的結(jié)果一致.LAD、LTS、LMS分別為對應(yīng)損失函數(shù)的結(jié)果,在求解(a,e,M0)和(i,?,ω)中采用了相同的損失函數(shù).

表1 采用不同損失函數(shù)時60 s弧長的軌道根數(shù)Tab le 1 The resu lts of 60-second arc length with d ifferen t loss fun ctions

由表1可見,在無野值的情況下,4種方法得到的結(jié)果相當(dāng),在精度范圍內(nèi)是一致的.理論上OLS估計具有最好的統(tǒng)計性質(zhì),但對于極短弧定軌而言,精度要求不高而且較少的觀測資料在統(tǒng)計性質(zhì)上與正態(tài)分布的偏差相對更大,OLS估計的優(yōu)勢并不能顯現(xiàn).

為了比較不同方法的穩(wěn)健性,對觀測資料引入野值.野值選取來自[120′′,240′′]的均勻分布,對不同野值比例分別采用4種損失函數(shù)進行了計算.考慮到隨機數(shù)種子對計算的影響,對每組采用不同隨機數(shù)種子計算了10次.表2給出了10 s弧長情況下對于初軌計算最為關(guān)心的軌道半長徑a的結(jié)果,其中p為野值比例,ˉa、amin和amax分別為10次計算結(jié)果的平均值、最小值和最大值.可明顯看出:當(dāng)p＞10%時,OLS估計已無法保證結(jié)果的可靠性,出現(xiàn)較大偏差,超出了精度范圍;p=40%時,LAD估計出現(xiàn)一定偏差,而LMS估計和LTS估計仍獲得了可接受的結(jié)果.

圖1給出了20 s弧長的計算結(jié)果,在相同野值比例情況下,LTS估計和LMS估計結(jié)果比較接近,OLS估計結(jié)果偏離比較大,LAD估計結(jié)果在p=40%時出現(xiàn)了較大偏離,和10 s弧長結(jié)果類似.圖2給出了30 s弧長下的計算結(jié)果,橫向比較每種方法隨著野值增多后的表現(xiàn),可明顯看出:對于OLS估計,隨野值增加結(jié)果變化很大,LAD方法的結(jié)果變化有一定范圍,但變化緩慢,LMS估計和LTS估計則隨著p的增加變化不大.

表2 10 s弧長不同野值比例下的軌道半長徑aTab le 2 T he sem i-m a jor ax is a of 10-secon d arc length with d ifferen t p rop or tion s of ou tliers

LMS估計和LTS估計在理論上效果非常接近,僅在漸近效率上有所區(qū)別.由于漸近效率在極短弧定軌中并不是主要問題,從計算結(jié)果看,兩者也確實十分類似,從處理的數(shù)據(jù)結(jié)果來看,LMS略微占優(yōu).結(jié)合極短弧定軌問題本身的特點,在日常計算中推薦使用LMS估計或LTS估計代替OLS估計.

5 結(jié)論與討論

短弧定軌問題觀測數(shù)據(jù)較少,即便可采用野值剔除手段,由于剔除標(biāo)準(zhǔn)是在迭代過程中建立的,在前幾次迭代中不當(dāng)剔除反而會對解產(chǎn)生較大的影響,迭代中可用資料數(shù)量大大減少,對解的穩(wěn)定收斂會有較大影響.由于剔除標(biāo)準(zhǔn)總是以當(dāng)前解作為參考,這種做法顯然也會對收斂至全局最優(yōu)解產(chǎn)生影響.遺傳算法過程中無法使用野值剔除方法,采用穩(wěn)健估計方法可有效克服野值帶來的影響,而且在遺傳算法過程中引入不同的損失函數(shù)幾乎不增加額外負(fù)擔(dān).

圖1 不同野值比例下20 s弧長的結(jié)果Fig.1 T he resu lts of 20-second arc length with d ifferen t p rop ortions of ou tliers

圖2 不同野值比例下30 s弧長的結(jié)果Fig.2 T he resu lts of 30-second arc length with d ifferen t p rop ortions of ou tliers

綜合本文和文獻[3],采用遺傳算法后經(jīng)典計算框架下的各種處理手段都可有效利用,并且引入的方法直觀、簡潔,較原有計算框架大大簡化,能夠滿足實測數(shù)據(jù)處理的各種需求.

本文討論不僅針對遺傳算法,由于絕大多數(shù)的現(xiàn)代優(yōu)化方法都有類似的情況,結(jié)果也可應(yīng)用于這些算法.

[1]M ilani A,Knezevic Z.CeM DA,2005,92:1

[2]Ansa lone L,Cu rti F.A d SpR,2013,52:477

[3]李鑫冉,王歆.天文學(xué)報,2016,57:66

[4]吳連大.人造衛(wèi)星和空間碎片的軌道和探測.北京:中國科學(xué)技術(shù)出版社,2011:181-210

[5]賈沛璋,吳連大.天文學(xué)報,2000,41:123

[6]王歆.天文學(xué)報,2013,54:274

[7]W ang X.ChA&A,2013,37:455

[8]De Jong K A.Evo lu tionary Com pu tation—A Un ified A pp roach.Cam b ridge:M IT P ress,2006:26-27

[9]Portnoy S,K oenker R.StaSc,1997,12:279

[10]Rousseeum P,Leroy A.Robust Regression and Outlier Detection.New York:John W iley&Sons, 1987:112-142

Genetic A lgorithm for Initial O rbit Determ ination with Too Short A rc(Continued)

LIXin-ran1,2,3WANG Xin1,2

(1 Pu rp le M oun tain Observatory,Chinese Academ y of Scien ces,Nan jing 210008) (2 K ey Labo rato ry fo r Space O bjec t and D ebris O bserva tion,Pu rp le M oun tain O bserva to ry,Chinese Academ y of Scien ces,Nan jing 210008) (3 Un iversity of Chinese A cadem y of Scien ces,Beijing 100049)

W hen using the genetic algorithm to solve the problem of too-short-arc (TSA)determ ination,due to the difference of com puting processes between the genetic algorithm and classicalmethod,the methods for outliers editing are no longer app licable.In the genetic algorithm,the robust estimation is acquired by means of using different loss functions in the fitness function,then the outlier problem of TSAs is solved.Compared with the classicalmethod,the app lication of loss functions in the genetic algorithm is greatly simplified.Through the comparison of results of different loss functions,it is clear that themethods of least median square and least trimmed square can greatly im prove the robustness of TSAs,and have a high breakdown point.

space vehicles,celestialmechanics,methods:numerical,statistical

P135;

A

10.15940/j.cnki.0001-5245.2016.02.005

2015-07-16收到原稿,2015-08-21收到修改稿

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

?lixr@pm o.ac.cn

?wangxin@pmo.ac.cn