張 彬，梁 楠,2，于徐紅，陳宗浩，劉志杰，游善平*

(1. 貴州師范大學(xué)貴州省信息與計(jì)算科學(xué)重點(diǎn)實(shí)驗(yàn)室，貴州 貴陽 550001；2. 中國科學(xué)院云南天文臺(tái)，云南 昆明 650216)

脈沖星是一種快速旋轉(zhuǎn)、高磁化的中子星，是理解超密度物質(zhì)的關(guān)鍵，可以用于各種基礎(chǔ)物理實(shí)驗(yàn)[1]。目前，脈沖星搜索方法主要分為周期性搜索和單脈沖搜索兩大類。其中，周期性搜索通過應(yīng)用快速傅里葉變換(Fast Fourier Transform, FFT)將時(shí)間序列轉(zhuǎn)化到頻域，以識(shí)別周期性信號(hào)，然后在確定的周期內(nèi)對(duì)原始時(shí)間序列數(shù)據(jù)進(jìn)行折疊，以提高周期性信號(hào)的信噪比[2]。而單脈沖搜索尋找強(qiáng)的、非周期的脈沖，并未使用快速傅里葉變換和折疊。單脈沖搜索非常適合發(fā)現(xiàn)周期性搜索中無法發(fā)現(xiàn)的孤立爆發(fā)。研究人員應(yīng)用單脈沖搜索發(fā)現(xiàn)了旋轉(zhuǎn)射電暫現(xiàn)源(Rotating Radio Transients, RRAT)[3]和快速射電暴的(Fast Radio Bursts, FRB)[4]。

2003年，文[5]首次提出了一個(gè)理論框架，運(yùn)用單脈沖搜索方法探測(cè)射電天文觀測(cè)數(shù)據(jù)中的脈沖星信號(hào)。2009年，文[6]把單脈沖搜索應(yīng)用于阿雷西博L波段饋源陣列脈沖星巡天(Pulsar Arecibo L-band Feed Array, PALFA)七波束觀測(cè)數(shù)據(jù)，發(fā)現(xiàn)了7顆脈沖星。2010年，文[7]使用單脈沖搜索方法對(duì)帕克斯多波束脈沖星巡天(Parkes Multibeam Pulsar Survey, PMPS)數(shù)據(jù)進(jìn)行重新分析，發(fā)現(xiàn)10顆旋轉(zhuǎn)射電暫現(xiàn)源。2015年，文[8]在單脈沖理論的基礎(chǔ)上，設(shè)計(jì)針對(duì)旋轉(zhuǎn)射電暫現(xiàn)源的搜索工具RRATtrap，應(yīng)用在綠岸射電望遠(yuǎn)鏡(Green Bank Telescope, GBT)350 MHz觀測(cè)數(shù)據(jù)，最終探測(cè)到18顆旋轉(zhuǎn)射電暫現(xiàn)源。2016年，文[9]首次將單脈沖搜索與機(jī)器學(xué)習(xí)相結(jié)合，提出遞歸峰值識(shí)別算法(Recursive Algorithm for Peak IDentification, RAPID)，實(shí)現(xiàn)了自動(dòng)化篩選脈沖星候選體，并探測(cè)到6顆脈沖星。

2016年9月,位于貴州喀斯特山區(qū)的500 m口徑球面射電望遠(yuǎn)鏡進(jìn)入科學(xué)運(yùn)行階段，望遠(yuǎn)鏡先后配置了超寬帶接收機(jī)和L波段19波束接收機(jī)[10]。為充分利用FAST優(yōu)異的探測(cè)能力，2018年，文[11]設(shè)計(jì)并實(shí)施了FAST多科學(xué)目標(biāo)同時(shí)掃描巡天，即同時(shí)使用多個(gè)數(shù)字終端采集脈沖星、中子氫、分子譜線、旋轉(zhuǎn)射電暫現(xiàn)源、快速射電暴等多個(gè)科學(xué)目標(biāo)的觀測(cè)數(shù)據(jù)。鑒于單脈沖搜索方法在各大射電望遠(yuǎn)鏡的成功應(yīng)用，本文對(duì)CRAFTS超寬帶數(shù)據(jù)文件的單脈沖搜索結(jié)果進(jìn)行研究，發(fā)現(xiàn)CRAFTS單脈沖搜索篩選得到脈沖星候選體存在數(shù)以萬計(jì)的假陽性樣本。因此，我們根據(jù)脈沖星信號(hào)在色散量信噪比曲線大致是高斯曲線以及信噪比數(shù)據(jù)分布是高斯分布的假設(shè)，構(gòu)造3個(gè)顯著區(qū)分脈沖星信號(hào)和干擾的特征，旨在緩解出現(xiàn)大量假陽性樣本的問題。試驗(yàn)結(jié)果表明，添加特征值判斷的單脈沖篩選工具與傳統(tǒng)方法相比，假陽性樣本數(shù)量減少了20%。

1 CRAFTS超寬帶數(shù)據(jù)單脈沖搜索實(shí)驗(yàn)

在射電天文觀測(cè)數(shù)據(jù)中發(fā)現(xiàn)脈沖星通常分為4個(gè)階段：收集、消色散、周期性搜索或單脈沖搜索和人工檢查[1]。第1階段，原始數(shù)據(jù)以電壓時(shí)間序列的形式由射電望遠(yuǎn)鏡收集；第2階段，脈沖輻射通過星際介質(zhì)(Interstellar Medium, ISM)的色散效應(yīng)，導(dǎo)致較低頻率的脈沖比較高頻率晚到達(dá)[12]，消色散可以去除這些與頻率有關(guān)的延遲效應(yīng)的影響；第3階段，使用周期性搜索或單脈沖搜索找出脈沖星候選體；第4階段，對(duì)判斷為脈沖星候選體的數(shù)據(jù)進(jìn)行人工檢查。

我們先使用PRESTO(Pulsar Exploration and Search Toolkit)脈沖星搜索工具處理2017年8月至2018年5月CRAFTS的超寬帶數(shù)據(jù)，處理流程包括去干擾、消色散和單脈沖搜索等步驟；然后，應(yīng)用單脈沖篩選工具RRATtrap從PRESTO的Single_pulse_search.py輸出文件中篩選出脈沖星候選體。

單脈沖篩選工具RRATtrap根據(jù)脈沖星信號(hào)的特性區(qū)分脈沖星信號(hào)和干擾。(1)脈沖星信號(hào)出現(xiàn)在一定色散量(Dispersion Measure, DM)范圍內(nèi)，在最佳色散量時(shí)檢測(cè)到峰值信噪比，而在該色散量上下信噪比下降，這是由于脈沖以不準(zhǔn)確的色散量消色散時(shí)導(dǎo)致脈沖展寬造成的。(2)由于信號(hào)在最佳色散量檢測(cè)到峰值信噪比，預(yù)計(jì)來自地面的信號(hào)(即射頻干擾)在色散量為0 pc·cm-3時(shí)達(dá)到峰值。同時(shí)射頻干擾(Radio Frequency Interference, RFI)不受色散效應(yīng)的影響，因此，脈沖星信號(hào)出現(xiàn)在一個(gè)非常大的色散量范圍內(nèi)[8]。

圖1是單脈沖篩選工具RRATtrap探測(cè)到已知脈沖星PSR B0540+23的診斷圖。左側(cè)子圖縱坐標(biāo)為試驗(yàn)色散量范圍，橫坐標(biāo)為觀測(cè)時(shí)間，CRAFTS超寬帶數(shù)據(jù)每個(gè)觀測(cè)文件為52 s，右側(cè)子圖描述試驗(yàn)色散量對(duì)應(yīng)的信噪比大小。從圖1可以看到呈紡錘狀的單脈沖事件組(在相鄰的色散量與時(shí)間窗口內(nèi)所有信噪比大于5的事件)在最佳色散量78 pc·cm-3時(shí)取得峰值信噪比68，而在該色散量上下，信噪比下降。其中，強(qiáng)的脈沖星信號(hào)標(biāo)記為Excellent；弱的脈沖星信號(hào)表標(biāo)記為Very good；不規(guī)則的脈沖星信號(hào)標(biāo)記為Good或Ok；射頻干擾標(biāo)記為RFI；宇宙噪聲標(biāo)記為Noise。在圖2中，我們看到在較寬的色散量范圍上極強(qiáng)的射頻干擾的實(shí)例。此外，在t=20～40 s，極強(qiáng)的射頻干擾錯(cuò)誤標(biāo)記為脈沖星信號(hào)。

圖1 RRATtrap探測(cè)到脈沖星PSR B0540+23的診斷圖

圖2 RRATtrap篩選出假陽性實(shí)例的診斷圖

本文試驗(yàn)選取2017年8月至2018年5月共計(jì)約32萬個(gè)超寬帶CRAFTS數(shù)據(jù)文件[1]進(jìn)行單脈沖搜索試驗(yàn)。由于處理的數(shù)據(jù)規(guī)模較大，對(duì)所有數(shù)據(jù)進(jìn)行全面的人工檢查工作量非常大。因此，我們先使用RRATtrap初步篩選脈沖星候選體。然后，我們對(duì)數(shù)據(jù)文件中標(biāo)記為Excellent，Very good，Good和Ok共計(jì)約5萬顆脈沖星候選體的診斷圖進(jìn)行手工檢查，發(fā)現(xiàn)僅有772個(gè)診斷圖真正具有脈沖星信號(hào)，對(duì)應(yīng)101顆已知脈沖星。通過

(1)

計(jì)算的假陽性率(False Positive Rate, FPR)達(dá)到98.5%，其中，F(xiàn)P是沒有脈沖星信號(hào)的數(shù)據(jù)文件標(biāo)記為候選脈沖星的數(shù)量；TP是含有脈沖星信號(hào)的數(shù)據(jù)文件標(biāo)記為候選脈沖星的文件數(shù)量。每100個(gè)脈沖星候選體中，僅有1～2個(gè)包含脈沖星信號(hào)(已知或新發(fā)現(xiàn)的)。由于CRAFTS超寬帶數(shù)據(jù)文件中包含的具體脈沖星數(shù)量未知，本文未對(duì)假陰性情況進(jìn)行分析。

2 區(qū)分脈沖星信號(hào)和干擾的特征

CRAFTS超寬帶數(shù)據(jù)單脈沖搜索結(jié)果表明，脈沖星候選體僅有1.5%真正具有脈沖星信號(hào)，這主要是未能成功區(qū)分脈沖星信號(hào)和干擾造成的，如圖2。本節(jié)首先探討不準(zhǔn)確的色散量對(duì)信噪比的影響。然后，構(gòu)造3個(gè)顯著區(qū)分脈沖星信號(hào)和干擾的特征，并在探測(cè)到的101顆脈沖星樣本中選取同時(shí)具有脈沖星信號(hào)、射頻干擾和宇宙噪聲的79顆脈沖星作為樣本，計(jì)算它們?cè)?個(gè)特征的取值情況。最后，根據(jù)脈沖星信號(hào)、射頻干擾和宇宙噪聲在3個(gè)特征取值分布的差異，提出合理閾值應(yīng)用在單脈沖篩選工具，對(duì)脈沖星候選體進(jìn)行進(jìn)一步篩選。

2.1 色散量與信噪比曲線

文[5]探討了試驗(yàn)與真實(shí)色散量之間的偏差對(duì)信噪比的影響。通常，對(duì)一個(gè)特定色散量的時(shí)間序列進(jìn)行多次下采樣并重新搜索，當(dāng)有效的采樣時(shí)間最接近脈沖寬度時(shí)，得到的信噪比最高。在不同的色散量信道中，隨著試驗(yàn)色散量與真實(shí)色散量的偏差越大，信噪比越小，如圖1。測(cè)量信噪比S(δDM)與真實(shí)信噪比S的比值與色散量的偏差δDM滿足[5]

(2)

這里，

(3)

其中，erf為誤差函數(shù)；δDM為試驗(yàn)色散量與真實(shí)色散量的偏差；Δv為總帶寬；v為中心頻率；W為脈沖寬度，單位ms。

利用(2)式和(3)式計(jì)算天體物理脈沖的信噪比和寬度，我們可以計(jì)算預(yù)期的色散量偏差δDM。圖3展示了在中心頻率546 MHz，時(shí)間分辨率為100 μs和匹配濾波使用30的PSR B2000+40一組單脈沖事件預(yù)期和擬合信噪比下降曲線。在這兩種情況下，信噪比峰值處的色散量是真實(shí)的色散量。我們計(jì)算預(yù)期的信噪比下降，利用峰值信噪比和匹配濾波得到脈沖寬度，和實(shí)際脈沖寬度很接近。首先，我們使用非線性最小二乘法(Non-linear Least Squares, NLS)對(duì)原始數(shù)據(jù)進(jìn)行回歸，得到擬合的峰值信噪比和脈沖寬度。然后，根據(jù)擬合結(jié)果，代入(2)式得到擬合的色散與信噪比曲線。由圖3可以看出，預(yù)期和擬合的色散量與信噪比曲線大致是高斯曲線。使用觀測(cè)到的峰值信噪比和匹配濾波得到脈沖寬度，預(yù)期的信噪比下降與擬合值在峰值右側(cè)非常接近，但在左側(cè)出現(xiàn)偏差。這是匹配濾波得到的脈沖寬度與實(shí)際脈沖寬度的偏差造成的。

2.2 3個(gè)顯著特征

文[13]指出，脈沖星的輪廓可以簡(jiǎn)化為高斯型(對(duì)大多數(shù)脈沖星來說是一個(gè)合理的近似)。按照文[5]提出的單脈沖搜索理論，天體物理脈沖的色散量和信噪比曲線大致是高斯曲線，如圖3。由于高斯曲線通常是對(duì)稱的，文[14]提出了兩種對(duì)稱特征(SIDM和SIS/N)分別表征單脈沖事件組的色散量和信噪比的對(duì)稱性。SIDM的計(jì)算公式為

(4)

其中，DMpeak為信噪比峰值對(duì)應(yīng)的色散量；DMmax為單脈沖事件組最大的色散量；DMmin為單脈沖事件組最小的色散量。SIS/N的計(jì)算公式為

(5)

其中，∑S/Nleft為信噪比峰值左側(cè)所有單脈沖事件的信噪比之和；∑S/Nright為信噪比峰值右側(cè)所有單脈沖事件的信噪比之和。

由(4)式和(5)式定義的SIDM和SIS/N的取值范圍在0～1之間，這兩個(gè)特征值越高，色散量與信噪比曲線越對(duì)稱，相反不遵循(2)式描述規(guī)律的射頻干擾，通常在色散量與信噪比空間是單調(diào)遞減(或遞增)。所以，它們通常具有接近于0的對(duì)稱值，可以與脈沖星信號(hào)進(jìn)行區(qū)分。圖4和圖5分別展示了脈沖星樣本SIDM和SIS/N的對(duì)稱值，并對(duì)脈沖星信號(hào)、射頻干擾和噪聲進(jìn)行了對(duì)比。從圖4和圖5可以看出，脈沖星信號(hào)通常具有較大的對(duì)稱值，而射頻干擾和噪聲的對(duì)稱值相對(duì)較小，這表明脈沖星信號(hào)的色散量和信噪比曲線比干擾更對(duì)稱。此外，圖4和圖5標(biāo)記一些對(duì)稱值比較大的噪聲，這些通常是高斯噪聲。

圖4 脈沖星樣本SIDM對(duì)稱值

圖5 脈沖星樣本SIS/N對(duì)稱值

文[5]預(yù)測(cè)在沒有任何宇宙噪聲和射頻干擾的情況下，信噪比的數(shù)據(jù)分布是高斯分布。峰度是描述總體數(shù)據(jù)分布與高斯分布陡緩程度的統(tǒng)計(jì)量。為了判斷脈沖星樣本信噪比的數(shù)據(jù)分布是否符合高斯分布，我們使用峰度值表征單脈沖事件組內(nèi)信噪比的數(shù)據(jù)分布與高斯分布的符合程度。峰度值的計(jì)算公式為

(6)

其中，S/N為單脈沖事件的信噪比；μ為單脈沖事件組內(nèi)平均信噪比；σ為單脈沖事件組內(nèi)信噪比的方差。

當(dāng)(6)式定義的峰度值是3時(shí)，信噪比的數(shù)據(jù)分布服從高斯分布，隨著峰度值與3的差距越大，其分布形態(tài)的陡緩程度與高斯分布的差異程度越大。圖6用箱線圖展示了脈沖星樣本的峰度值分布情況，并對(duì)脈沖星信號(hào)、射頻干擾和噪聲進(jìn)行了對(duì)比。從圖6可以看出，脈沖星信號(hào)峰度值集中在3附近，而射頻干擾和噪聲的峰度值分散，且與3差距較大，表明脈沖星信號(hào)信噪比分布近似服從高斯分布，而大部分干擾與噪聲不具有此規(guī)律。

圖4～圖6可以明顯看出脈沖星信號(hào)、射頻干擾和宇宙噪聲在3個(gè)特征分布的差異。因此，我們根據(jù)脈沖星信號(hào)和干擾特征值的分布區(qū)間(表1)，選擇合適的特征閾值，剔除不滿足閾值的脈沖星候">選體，從而達(dá)到降低假陽性率的目的(比如設(shè)置脈沖星信號(hào)特征SIDM最小值0.1作為閾值，對(duì)低于閾值的候選體視為干擾，那么，特征值在0.02～0.1區(qū)間的射頻干擾將不會(huì)錯(cuò)誤標(biāo)記為脈沖星候選體)。同時(shí)，我們選取不同閾值進(jìn)行試驗(yàn)，并統(tǒng)計(jì)其對(duì)應(yīng)的假陽性率和遺漏脈沖星信號(hào)的結(jié)果。由圖7各個(gè)特征在不同閾值假陽性及遺漏脈沖星的情況可以看出，隨著閾值變大或者閾值區(qū)間縮小，假陽性率不斷減小，與此同時(shí)，遺漏發(fā)現(xiàn)的脈沖星數(shù)量在不斷增加。脈沖星搜尋的首要前提是保證觀測(cè)數(shù)據(jù)中所有脈沖星信號(hào)不遺漏。因此，我們選取表1中脈沖星信號(hào)特征SIDM最小值0.1和SIS/N最小值0.29作為閾值，特征峰度-1.54～27.11作為閾值區(qū)間，對(duì)低于閾值或不在閾值區(qū)間的脈沖星候選體視為射頻干擾或宇宙噪聲。

圖6 脈沖星樣本峰度值

表1 3個(gè)特征在3種類別的分布區(qū)間

圖7 各個(gè)特征在不同閾值假陽性和遺漏脈沖星的情況，其中紅色豎線代表所選擇的閾值

此外，為了證實(shí)所述特征是否在其他觀測(cè)數(shù)據(jù)有效，我們對(duì)帕克斯多波束數(shù)據(jù)進(jìn)行了試驗(yàn)。根據(jù)文[15]公開的帕克斯單脈沖數(shù)據(jù)庫，我們構(gòu)建了一個(gè)PMPS(Parkes Multibeam Pulsar Survey)數(shù)據(jù)庫，其中包含帕克斯望遠(yuǎn)鏡發(fā)現(xiàn)的部分脈沖星信號(hào)。我們計(jì)算得到PMPS數(shù)據(jù)集中脈沖星信號(hào)在3個(gè)特征值的分布情況，并與CRAFTS進(jìn)行對(duì)比，如表1。由表1中CRAFTS及PMPS的脈沖星信號(hào)在3個(gè)特征值分布情況可以發(fā)現(xiàn)，它們的特征值分布十分接近，盡管它們?cè)谙⒕W(wǎng)格、周圍射頻干擾的環(huán)境等存在顯著差異。因此，本文所述的特征閾值同樣能夠區(qū)分帕克斯數(shù)據(jù)庫中脈沖星信號(hào)和干擾。

我們將上述閾值和閾值區(qū)間應(yīng)用于單脈沖篩選工具，對(duì)脈沖星候選體進(jìn)行重新篩選，試驗(yàn)結(jié)果如表2。由表2可以看出，添加單個(gè)特征閾值判斷的單脈沖篩選工具能在一定程度上降低假陽性率，并且3個(gè)特征閾值組合判斷假陽性率從98.5%降低到78.4%，表明添加特征閾值判斷的單脈沖篩選工具提高了脈沖星搜索效率。

表2 添加特征后假陽性率結(jié)果

3 總結(jié)與展望

本文使用PRESTO的單脈沖搜索方法對(duì)CRAFTS超寬帶數(shù)據(jù)文件進(jìn)行的試驗(yàn)表明，現(xiàn)有的單脈沖篩選工具難以區(qū)分真實(shí)信號(hào)與噪聲或射頻干擾信號(hào)，造成數(shù)以萬計(jì)的假陽性樣本出現(xiàn)，顯著增加人工篩選單脈沖候選體的時(shí)間開銷以及候選數(shù)據(jù)存儲(chǔ)壓力。因此，我們提出3個(gè)顯著區(qū)分脈沖星信號(hào)和干擾的特征，并選取同時(shí)具有脈沖星信號(hào)、射頻干擾和宇宙噪聲的79顆脈沖星樣本，計(jì)算它們?cè)?個(gè)特征取值情況。然后，根據(jù)脈沖星信號(hào)、射頻干擾和宇宙噪聲在3個(gè)特征取值分布的差異，提出合理閾值應(yīng)用于單脈沖篩選工具，用于對(duì)脈沖星候選體進(jìn)行嚴(yán)格的判斷。最后，使用添加特征閾值判斷的篩選工具對(duì)脈沖星候選體進(jìn)行重新試驗(yàn)。結(jié)果表明，在保證所有脈沖星信號(hào)不遺漏的情況下，假陽性率從98.5%降低到78.4%。因此，本文所述3個(gè)特征具有實(shí)用性和有效性，有助于單脈沖搜索在CRAFTS巡天數(shù)據(jù)的應(yīng)用。

致謝：本文在500 m口徑球面射電望遠(yuǎn)鏡(FAST)數(shù)據(jù)基礎(chǔ)上完成。FAST是由中國科學(xué)院國家天文臺(tái)運(yùn)行和管理的國家大科學(xué)裝置。感謝中國科學(xué)院天文大科學(xué)研究中心FAST重大成果培育項(xiàng)目對(duì)本文工作的資助。