Jurkevich方法的誤差分析

楊海燕，李孝攀，黨興菊

（昭通學(xué)院，物理與信息工程學(xué)院，云南 昭通 657000）

0 引言

為分析部分活動星系核可能存在的周期性光變，天文學(xué)家們提出了一些可行的周期分析算法。這些算法都充分考慮了天文觀測數(shù)據(jù)非均勻、大誤差的特點(diǎn)。例如：Jurkevich于1971年提出的一種專門用于分析非等間距數(shù)據(jù)的期望均方誤差統(tǒng)計(jì)方法(Jurkevich方法)[1]；Scargle于1981年在Lomb算法[2]基礎(chǔ)上提出的著名的Lomb-Scargle(LS)周期圖法[3]。Stellingwerf于1978年提出的狀態(tài)彌散最小化算法(Phase Dispersion Minimization,PDM)算法[4]；，F(xiàn)erraz-Mello于1981年提出的時間補(bǔ)償離散傅里葉變換(Date-compensated Discrete Fourier Transform, DCDFT)算 法[5]；Simonetti于1985年提出的結(jié)構(gòu)函數(shù)（Structure Functions, SF）法[6]；Edelson 和 Krolik于1988年提出的離散相關(guān)函 數(shù) (Discrete Correlation Function, DCF)算 法[7]；Alexander于1997年對DCF算法進(jìn)行了改進(jìn)和發(fā)展，得到的Z變換離散相關(guān)函數(shù)(Z-transformed Discrete Correlation Function, ZDCF)算法[8]。Foster于1996年提出的對非平穩(wěn)信號分析有獨(dú)特的優(yōu)勢的加權(quán)小波Z變換(Weighted Wavelet Z-transform,WWZ)算法[9]。以上方法被廣泛應(yīng)用于活動星系核周期性光變的研究中，對天體輻射機(jī)制的研究、輻射區(qū)的限定、天體質(zhì)量的估算等起到了重要的促進(jìn)作用。

天文觀測的局限性會導(dǎo)致觀測數(shù)據(jù)不連續(xù)性。Jurkevich方法在不等間距數(shù)據(jù)處理方面有著獨(dú)特的優(yōu)勢，被廣泛用于活動星系核長周期光變的研究。1998年，樊等[10]使用Jurkevich方法研究了NGC4151自有數(shù)據(jù)記載以來的光學(xué)資料，發(fā)現(xiàn)NGC4151的光學(xué)波段可能有（14.04±0.8）a的周期性光變；2000年，蘇[11]也應(yīng)用Jurkevich方法對PKS2251+158近100年光變數(shù)據(jù)進(jìn)行研究，得到該目標(biāo)在光學(xué)B波段的光變周期約為12.39 a；2001年，樊等[12]再次使用Jurkevich方法研究了3C273自1887年至1997年以來的光學(xué)數(shù)據(jù)，發(fā)現(xiàn)光變曲線具有明顯的2.0 a光變周期。近年，Li等[13]使用LS方法和Jurkevich方法研究了3C 454.3的多波段光變曲線，發(fā)現(xiàn)該目標(biāo)在紅外、光學(xué)、軟X射線和γ射線波段的光變均存在(210.8±12)d的周期；使用Jurkevich方法，Li等[14]還得出Mrk421射電15 GHz、X射線和γ射線波段的光變曲線可能具有(287.6±4.4) d、(309.5±5.8) d和(283.4±4.7) d的光變周期；Fan等[15]使用ACF、SF和Jurkevich方法研究了AO 0235+164 R波段2004-2012年光變曲線，得到了該目標(biāo)在R波段可能具有(0.55±0.03) a的光變周期；Xiong等[16]使用REDFIT，LS和Jurkevich方法研究了3C 273 V波段1968-2016年的光變曲線，得到了該目標(biāo)在V波段可能具有(3918±1121) d的光變周期。

本文在介紹Jurkevich方法原理的基礎(chǔ)上，推導(dǎo)并討論了該方法的誤差公式，最后使用該方法研究了耀變體PKS 2155-304光學(xué)I波段光變曲線的光變周期。

1 Jurkevich方法

Jurkevich方法[1]是基于期望均方差的一種統(tǒng)計(jì)學(xué)方法，被用于處理非等間距數(shù)據(jù)。數(shù)據(jù)樣本是由N個數(shù)據(jù)構(gòu)成的，其中樣本的平均值、離差平方和、方差分別為：

Jurkevich方法將數(shù)據(jù)根據(jù)試驗(yàn)周期來進(jìn)行折疊分組處理，然后對得到的每一批分組數(shù)據(jù)計(jì)算平均值、離差平方和再對每組的離差平方和求和得到總離差平方和其中第l分組的平均值、離差平方和及方差分別為：

以上三個式子中的Nl為第l組觀測數(shù)據(jù)的個數(shù)，則m組各自離差平方和的總和為：

對給定的數(shù)據(jù)樣本，V2不受試驗(yàn)周期影響，但和對試驗(yàn)周期的變化很敏感。當(dāng)試驗(yàn)周期不是真實(shí)周期時，很?。划?dāng)試驗(yàn)周期與真實(shí)周期很接近時，由于各組平均值的變化差異較大，此時大很多，即減小。因此當(dāng)取得最小值時，滿足試驗(yàn)周期和真實(shí)周期相等。令試驗(yàn)周期為P，可通過作找到圖中的極小值點(diǎn)對應(yīng)的試驗(yàn)周期值來確定數(shù)據(jù)樣本的周期。

2 Jurkevich方法的誤差

2.1 誤差公式的推導(dǎo)

在物理測量中，設(shè)被測量y和各直接測量量x1，x2，…，xn之間的函數(shù)關(guān)系為：

若xi的誤差為△xi，則y的誤差可表示為

上式為誤差傳遞公式。根據(jù)以上公式，在Jurkevich方法中，我們將xi視為獨(dú)立變量，且只考慮xi測量的誤差△xi，則V2的誤差可表示為

在Jurkevich方法中，我們將所有的數(shù)據(jù)按照相位進(jìn)行分組，分組后每一個分組的和未分組的V2函數(shù)形式是一樣的。因此，的誤差為

2.2 誤差公式的驗(yàn)證和討論

為研究Jurkevich方法誤差公式，我們構(gòu)造了一個信號

添加均值為0.3、方差為0.05的高斯白噪聲，如圖1所示，使用Jurkevich方法隨機(jī)取100個不等間距數(shù)據(jù)進(jìn)行分析，取分組數(shù)為9，如圖2所示，由于誤差較小，Jurkevich圖中不能直觀顯示，這里截取圖中第一個周期成分所在位置附近進(jìn)行放大，從放大圖中可以看出Jurkevich方法的誤差，但該誤差對Jurkevich方法的結(jié)果影響并不大。使用(15)式得到的Jurkevich方法誤差的平均值為0.021，標(biāo)準(zhǔn)差為0.001.當(dāng)取分組數(shù)為5時，使用(15)式得到的Jurkevich方法誤差的平均值為0.012，標(biāo)準(zhǔn)差為0.008.

圖1 添加高斯白噪聲的信號圖

圖2 添加均值為0.3、方差為0.05的高斯白噪聲后的信號及分組數(shù)為9的Jurkevich方法結(jié)果

天文觀測數(shù)據(jù)，特別是X射線和γ射線的觀測數(shù)據(jù)常常具有較大的誤差。為進(jìn)一步研究Jurkevich方法的誤差，我們采用(16)式構(gòu)建的信號，添加均值為1.0、方差為0.1的高斯白噪聲，使用Jurkevich方法隨機(jī)取100個不等間距數(shù)據(jù)進(jìn)行分析，取分組數(shù)為9，如圖3所示，由于誤差較小，Jurkevich圖中不能直觀顯示，因此截取圖3中第一個周期成分所在位置附近進(jìn)行放大，從放大圖中可以看出Jurkevich方法的誤差，但該誤差對Jurkevich方法的結(jié)果影響并不大。根據(jù)(15)式得到分組數(shù)為9時的Jurkevich方法誤差的平均值為0.040，標(biāo)準(zhǔn)差為0.013.當(dāng)取分組數(shù)為5時，Jurkevich方法誤差的平均值為0.032，標(biāo)準(zhǔn)差為0.015.

經(jīng)過以上分析可以看出，當(dāng)測量值xi的誤差較大時，Jurkevich方法的誤差也會隨之變大。分組數(shù)的選擇對Jurkevich方法的誤差也有影響，且誤差會隨分組數(shù)增大而增大；盡管如此，Jurkevich方法仍能直觀高效地反映出待測信號的周期特征，該方法的誤差較小，且遠(yuǎn)遠(yuǎn)小于觀測數(shù)據(jù)的誤差，并不會影響值的判斷。

3 PKS 2155-304光學(xué)I波段的光變周期

PKS 2155-304是南天區(qū)一個很好的高頻峰BL Lac研究對象，紅移量z為0.116，是在X射線和遠(yuǎn)紫外波段最亮的BL Lac天體之一 。PKS 2155-304是自從上世紀(jì)70年代發(fā)現(xiàn)以來最受關(guān)注的研究目標(biāo)之一，目前可以使用的觀測數(shù)據(jù)覆蓋了射電波段到TeV能級的γ射線波段[18-19]。2000年，F(xiàn)an & Lin[20]收集了近17年的UBVRI的觀測數(shù)據(jù)，使用Jurkevich方法等得到了光學(xué)V波段具有約4.16 a和7.0 a的準(zhǔn)周期。2012年，畢等人[21]通過對2004到2011年的觀測數(shù)據(jù)應(yīng)用Jurkevich方法進(jìn)行分析，得到PKS 2155-304光學(xué)R波段具有314 d的中短時標(biāo)光變周期。Zhang等人[22]搜集到了約35 a的R波段光變曲線，使用Jurkevich方法等三種方法得到了PKS 2155-304在光學(xué)R波段具有317 d的準(zhǔn)周期。

本文從文獻(xiàn)[23]中搜集了2005年4月11日到2012年6月30日的281個非均勻觀測數(shù)據(jù)，可以看到，在整個觀測區(qū)間，PKS 2155-304表現(xiàn)出劇烈的光變，其中亮度的最大值11.332星等，最小值為13.492星等，且表現(xiàn)出一定周期性。圖4為BL Lac天體PKS 2155-304在光學(xué)I波段的光變曲線。

圖3 添加均值為1.0、方差為0.1的高斯白噪聲后的信號及分組數(shù)為9的Jurkevich方法結(jié)果

圖4 PKS 2155-304光學(xué)I波段光變曲線

使用Jurkevich方法對PKS 2155-304光學(xué)I波段的光變曲線進(jìn)行周期分析的結(jié)果如圖2所示。使用(15)式得到的Jurkevich方法誤差的平均值為4.51×10-4，標(biāo)準(zhǔn)差為1.96×10-4。由圖5可以看出該天體在I波段的光變曲線中有兩個明顯的低谷，對這兩個低谷做高斯擬合，以FWHM作為這兩個周期的誤差。擬合結(jié)果為P1=(312.00±16.17) d，P2=(644.00±9.77)d，f=2.16＞0.5，上述結(jié)果都符合F判據(jù)中具有明顯周期性的條件，因此對應(yīng)有兩個明顯的周期成分。使用Jurkevich方法對周期可靠性判斷應(yīng)當(dāng)滿足以下兩個條件：(1)數(shù)據(jù)時間跨度大于準(zhǔn)周期的6倍；(2)有明顯振幅。結(jié)果中P1=(312.00±16.17)d滿足Jurkevich方法對周期可靠性的判斷，所以P1=(312.00±16.17) d是PKS 2155-304在光學(xué)I波段具有的準(zhǔn)周期，而P2=(644.00±9.77) d則不滿足時間跨度大于6個周期的限制條件，故舍去。

圖5 PKS2155-304光學(xué)I波段Jurkevich圖

使用Jurkevich方法對PKS 2155-304光學(xué)I波段2005年到2012年近7年的光變曲線進(jìn)行周期分析的結(jié)果表明，PKS 2155-304在光學(xué)I波段可能具有(314.17±11.91) d的光變周期。這一結(jié)果與Sandrinelli等[24]發(fā)現(xiàn)的PKS2155-304在光學(xué)B波段約315 d的光變周期一致，也同Zhang等[22]發(fā)現(xiàn)的PKS2155-304在光學(xué)B波段約317 d的光變周期相吻合。表明以上兩個波段的輻射有著良好的相關(guān)性。

4 結(jié)語

本文在介紹Jurkevich方法原理的基礎(chǔ)上，推導(dǎo)并討論了該方法的誤差公式，最后使用該方法研究了耀變體PKS 2155-304光學(xué)I波段光變曲線的光變周期。通過對平均誤差為0.01為耀變體PKS 2155-304光學(xué)I波段的光變曲線進(jìn)行了周期分析，根據(jù)(15)式得到的Jurkevich方法誤差的平均值為4.51×10-4，標(biāo)準(zhǔn)差為1.96×10-4，以及PKS 2155-304在光學(xué)I波段可能具有（314.17±11.91）d的光變周期。進(jìn)一步驗(yàn)證了Jurkevich方法能夠高效地分析不等間距數(shù)據(jù)中的周期成分，并且是一種研究光變周期誤差極小的方法。同時也證明即使在待分析數(shù)據(jù)誤差極大的情況下，該方法的誤差對低谷的選取和F判據(jù)的使用都沒有影響。