孫萬鵬,張富文
(桂林理工大學理學院,廣西 桂林 541004)
自2007 年快速射電暴被首次發現以來,該領域已經迅速發展。快速射電暴被認為是宇宙中最神秘的天文現象之一,其強大的脈沖信號具有極高的峰值亮度和短暫的持續時間。雖然人們已經對快速射電暴的基本特征有了一定的了解,但仍有許多關鍵問題需要解決,其中一個便是快速射電暴的分類問題。得益于CHIME團隊提供的第一個快速射電暴大樣本,目前已經有部分研究者通過快速射電暴的關鍵觀測參量之間的差異對其進行分類[1-2],還有其他研究者根據其是否重復爆發的特性將其分為重復暴和非重復暴進行討論[3-4]。然而,當前樣本中重復暴和非重復暴的比例與理論預測結果嚴重不符,這意味著可能有大量重復暴未被觀測到而被錯誤地分類為非重復暴群體[5]。因此,尋找快速射電暴中潛在的重復暴,實現對快速射電暴清晰地分類尤為重要。
現有的大多數研究往往只是簡單地根據某一個或兩個觀測參數的差異來判斷重復暴和非重復暴是否屬于不同類別,很少有研究能夠綜合考慮多個觀測性質的差異,并提出一般的區分方法。而無監督機器學習擅長從輸入數據關系中揭示隱藏信息,因此利用這一方法可以有效地幫助我們結合快速射電暴的多維度特征對其進行分類,并發現潛在的重復暴。
在本文中,我們采用了一種名為t-SNE 的無監督機器學習算法,用于將CHIME/FRB 目錄中的重復暴和非重復暴進行分類。通過這種方法,我們找到了在當前CHIME/FRB 樣本中潛在的重復暴候選體,并提供了一個“未受污染”的樣本,可以清晰地區分重復暴和非重復暴。基于分類的樣本,我們進行了更進一步的統計分析,全面討論了重復暴與非重復暴的參數分布以及一些參數之間的關系。基于重新分類后的統計分析有助于我們揭示重復暴與非重復暴背后不同的物理機制。同時,我們的研究結果為進一步搜尋重復暴和科學地分類快速射電暴提供了指導。
本文的第一章介+紹了樣本選擇和數據處理,并簡要介紹了用于無監督算法的特征參數。第二章給出了t-SNE 算法的介紹和主要參數配置。第三章給出了CHIME/FRB 樣本無監督分類的結果,包括識別出的重復暴候選體,同時對比分析了重復暴和非重復暴的各觀測參數分布以及快速射電暴在觀測者坐標系和靜止坐標系中的統計特性。在第四章中總結了我們的工作。
我們選取了CHIME/FRB 發布的536 個快速射電暴作為我們的樣本,其中包括474 個非重復暴和62 個重復暴,這些快速射電暴樣本均在2018 年7 月25 日至2019 年7 月1 日之間探測到,探測頻率范圍為400 MHz至800 MHz。我們首先排除了6 個未觀測到流量和通量的快速射電暴(均為非重復暴),并在樣本中考慮到一些快速射電暴包含的多個子暴。因此,我們得到了一個包含500 個非重復暴和94 個重復暴的數據集,稱之為樣本1。下面是樣本1 中快速射電暴主要觀測參量的簡要概述,但在此節中我們僅選擇其中前7 個參量用于對快速射電暴進行無監督降維分類:
Width of Sub-Burst (s):這些值表示每個子暴的爆發寬度,具體通過擬合算法f itburst 獲得。紅移展寬效應沒有從這個參數中去除,這些值在快速射電暴的每個子暴之間是不同的。
Scattering Time (s):這個值表示在600 MHz 時由于散射導致的脈沖展寬效應時間。紅移展寬效應沒有從這個參數中去除。這些值對于快速射電暴的每一個子暴都是相同的。
Flux (Jy):表示平均帶寬中的峰值流量。該值對于快速射電暴的每一個子暴都是相同的。
Fluence (Jy·ms):這個值表示所有子脈沖的累計亮度之和。該值對于快速射電暴的每一個子暴都是相同的。
Spectral Index:這些值表征了每個子暴的光譜形狀。在快速射電暴的每個子暴之間光譜指數是不同的。
Spectral Running:這些值代表了頻譜形狀與頻率間的依賴性。這些值在快速射電暴的每個子暴之間是不同的。
Lowest Frequency (MHz):這些值代表的是在脈沖全寬十分之一最大峰值處子脈沖探測的最低頻帶。來自快速射電暴的每個子暴之間的數值是不同的。
Boxcar Width (s):箱式濾波器寬度是結合所有子脈沖的持續時間,包括儀器、散射和紅移增寬效應。這些值是綜合每個子脈沖的平均取值,因此對于快速射電暴的每一個子暴都是相同的。
Highest Frequency (MHz):這些值代表的是在脈沖全寬十分之一最大峰值處子脈沖探測的最高頻帶。來自快速射電暴的每個子暴之間的數值是不同的。
Peak Frequency (MHz):這些值表示每個子暴的峰值頻率。來自快速射電暴的每個子暴之間的數值是不同的。
Dispersion Measure (DMs) (pc cm-3):探測器觀測所得到的色散量,對于每個子暴來說這些值都是相同的。
t-SNE(t-distributed Stochastic Neighbor Embedding)是一種用于降維的無監督機器學習算法,通常用于將高維數據降到三維或者二維進行可視化[6]。從高維映射到低維空間可視化時會綜合考慮保留數據的局部和全局結構。它依賴于隨機初始化將數據點的相似程度轉換為概率,在同一數據集上運行t-SNE 可以產生具有相似拓撲的各種映射。因此對于原本在高維空間中高相似度的數據點在映射到低維空間數據間的距離會更近一些;而對于低相似度的點,映射到低維空間距離會更遠。所以通過t-SNE 得出的數據降維分布圖的坐標軸和標簽沒有特別的意義,只能夠告訴我們聚集在一起的點相似度更大,而相互遠離的點則表明相似度更低。
t-SNE 初始化過程中比較重要的幾個參數有perplexity、early exaggerate、learning rate、number of steps。其中困惑度(perplexity)簡單的可以理解為t-SNE 在生成條件概率時考慮的最近鄰數,通常數據集越大需要使用的困惑度也越高同時所考慮的最近鄰數也越多。也就是說高的困惑度通常會降低對小結構的敏感度,相反較低的困惑度則會考慮更少的最近鄰數,因此會更加注重局部而忽略全局結構。在我們的工作當中困惑度的大小為63.0。早期放大因子(early exaggerate)的默認值是12,此參數的調整會影響到原始空間中的高維數據在嵌入低維空間的緊密程度,有時一些隱約可見的差距可以通過這個參數的調整來觀察映射到低維空間數據的真實差別。但如果數據間本就沒有顯著差別,那么無論如何調整此參數也不會對結果有變化。在我們的工作當中早期放大因子的數值為22.0。另一個關鍵參數是學習率(learning rate),默認值為200,通常在[10.0,1000.0]區間內進行調整。如果學習率太高,數據可視化圖像會變得過于均勻,任何一點與其他相鄰點的距離基本相等;但如果太低則會造成可視化數據被密集的壓縮,幾乎無法判斷真實結構。因此以上這些參數的調整十分重要,需要仔細嘗試,在我們的工作當中學習率的大小設置為265.0。有關超參數的具體介紹可以在sklearn.manifold.t-SNE in python 中查看。
圖圖1 顯示了t-SNE 算法通過結合快速射電暴的子暴寬度、散射時間、光子流量、時間積分通量、譜指數、譜斜率變化率和最低頻率這幾個關鍵特征,對樣本1 進行降維并可視化的結果。我們可以看到快速射電暴明顯分為了兩部分,左側群體與右側群體之間存在清晰界限。而有趣的是,幾乎所有觀察到的重復暴(實心三角)都聚集在右側,而左側幾乎僅聚集了非重復暴(空心圓點),只有1 例重復暴混入其中。在1.2 節中,我們已經解釋了t-SNE 可視化圖像的坐標不代表任何物理意義,只能表示彼此接近的快速射電暴具有更高的相似性。因此,無監督學習算法的分類結果表明,重復暴和非重復暴之間確實存在顯著的特征差異。

圖1 經t-SNE 算法實現的CHIME/FRB 目錄快速射電暴的二維投影,536 個暴被清晰地分為左右兩個簇,在左邊,幾乎所有的快速射電暴都是非重復暴(空心圓點),只有一個重復暴混在其中,觀測到的重復暴(實心三角)幾乎全部分布在右邊的簇中,并且有大量的非重復暴也混在其中
根據無監督降維的結果,可以觀察到圖1 右側區域的一些非重復暴與幾乎所有的重復暴混合在一起,表現出明顯的相似性。因此,我們推測這些與重復暴混合在一起的非重復暴是潛在的重復暴候選體。在圖2 中,這些重復暴候選體被用十字號標出。另外,對于混合在圖1 左側群體中的那例重復暴,我們沒有在其頻譜圖中發現任何特殊特征,因此我們推測可能是測量誤差造成了這種情況。因此,對于聚集在圖1 左側的快速射電暴我們將其全部稱為非重復暴,在圖2 中用空心圓點表示。在接下來的章節中,我們將以t-SNE 算法分類識別的重復暴、重復暴候選體和非重復暴作為新的研究樣本,稱為樣本2,并將重復候選體和重復暴合并為一類進行討論。

圖2 去除污染CHIME/FRB 目錄的分類結果,隱藏在重復暴組中的非重復暴被識別為重復暴候選體,以十字號表示
此前,在CHIME/FRB 目錄中觀測到的重復暴及其子暴有94 個,約占CHIME/FRB 目錄的94/(94+500)=15.8%。如果考慮t-SNE 算法識別的165 個重復暴候選體,當前重復暴及子暴的比例約為(94+165)/(94+500)=43.6%。這表明重復暴的真實比例應該比當前觀測報告的要大得多。這意味著CHIME/FRB 目錄中重復暴與非重復暴間嚴重的“污染”狀況阻礙了對它們相應物理機制的研究。在接下來的小節中,我們將重新研究樣本2 中實際重復暴和非重復暴之間各物理量的分布統計和相關性。
2.1.1 參數分布
為了更直觀地解釋重復暴和非重復暴中參數的分布,我們繪制了樣本1 和樣本2 的參數分布直方圖,如圖3 和圖4 所示。研究的物理量是1.1 節介紹的與快速射電暴特性相關的11 個物理量。除譜指數和譜斜率變化率外,其他參數均采用對數處理。從圖4 可以看出,重復暴和非重復暴的譜指數、譜斜率變化率和頻率帶寬分布存在顯著差異。而其他觀測到的物理量的分布大致相似,無顯著差異。

圖3 CHIME/FRB 目錄中一些重要觀測參數的分布,頻率帶寬指一次爆發中最高頻率和最低頻率之間的差異,重復FRB 和非重復FRB 的分布分別用黑色陰影和黑色階梯線表示

圖4 經t-SNE 重新分類的快速射電暴參數分布;重復暴的分布用黑色陰影表示,非重復暴的分布用黑色階梯線表示
此前,有多位研究者表示重復和非重復暴之間的頻率帶寬存在顯著差異[7],我們的研究再次證實了這一點。從圖4 可以看出,重復暴的帶寬明顯更窄,主要分布在左側。相比之下,右側有大量帶寬為400 MHz 的非重復暴,它們的帶寬跨越了整個CHIME 頻段,甚至無法準確確定它們的真實帶寬。未來需要使用多個超帶寬探測器對非重復暴進行檢測才能解開這個謎團。
如圖4 所示,重復暴與非重復暴的譜指數和譜斜率變化率兩參量的數值也明顯不同。此前,也有研究人員發現重復暴與非重復暴之間的光譜指數和譜斜率變化率存在一定差異[2]。根據CHIME/FRB Collaboration的工作[7],光譜形狀由光譜指數(γ),譜斜率變化率(r)和頻率表示:
Fi表示光束強度的頻率相關值,表示光束強度的頻率相關值,通常代表不同的光譜形狀。其中,f和f0分別表示爆發頻率和CHIME 觀測的下限頻率(400.2 MHz)。因此,光譜指數和譜斜率變化率值的差異反映了快速射電暴光譜形狀的差異。我們的研究證實了重復暴和非重復暴之間光譜形態存在顯著差異,這也意味著它們可能具有不同的物理起源或發射機制。對于圖4 中的其他觀測參數分布,重復暴和非重復暴之間沒有明顯差異。
之前的多項工作都表明,無論是在觀測脈沖寬度還是固有寬度方面,重復暴通常比非重復暴具有更寬的脈沖[8]。但值得注意的是,在我們的“去污染”樣本2 中,我們發現這在兩個群體之間沒有顯著差異(見圖4)。換句話說,我們發現脈沖寬度并不是區分重復暴和非重復暴的重要指標。也有其他研究者注意到了這一點,Connor 等人通過光束發射的選擇效應來解釋這種差異[9]。我們認為,之前研究中表現出爆發寬度存在差異的原因是嚴重污染,因為可能存在被錯誤分類的潛在重復爆發,從而導致先前研究中未考慮的統計差異。這一結論可能需要通過未來發布更大的觀測樣本來逐步證實。
此外,在圖4 中我們還發現,樣本2 中重復暴與非重復暴的色散量和散射時間分布基本一致,這與其他幾位研究者的結論相同。這些反映出對于重復暴和非重復暴來說,爆發環境和宿主星系的性質也應該是相似的。對于這個推論,Bhandari 等人和Ravi 等人通過對已識別的快速射電暴源星系的比較研究提供了更直接的證據[10-11]。但要注意也有不同的討論結果[12]。
為了進一步研究靜止坐標系中快速射電暴的固有特性,我們根據CHIME/FRB 目錄中觀測到的色散值估計了紅移z的上限。由于我們的樣本中幾乎沒有測量到快速射電暴的紅移值,也沒有對宿主星系進行相應的色散測量,因此在我們的計算中統一忽略了宿主星系色散測量的貢獻,這樣會使導出的紅移值偏高。然而,這并不影響我們研究重復暴和非重復暴之間的統計差異。紅移的推導由FRUITBAT 包(https://fruitbat.readthedocs.io.)實現,其中色散測度—紅移關系采用Inoue2004,宇宙學模型為Planck2018。
有了z的上限,就可以推導出快速射電暴的各向同性能量Eiso、各向同性峰值光度Liso和亮溫TB度的上限。根據Zhang 提到的計算方法[13]:
其中Fobs是時間積分通量(單位erg cm-2Hz-1為或Jy·ms),Sv,p是峰值流量(單位erg s-1cm-2Hz-1為或Jy),DL表示光度距離。
亮溫度TB上限可近似表示為:
其中κB代表玻爾茲曼常數,W代表爆發的持續時間。值得注意的是,在Petroff等人和Zhang 的研究中分別采用望遠鏡的帶寬B和中心頻率vc推導出Eiso、Liso和TB等物理量的值[13-14]。但是,本工作中有大量爆發已完全跨越或超過了CHIME 的接收帶寬。因此,對于這一部分爆發在計算中不適合繼續使用設備帶寬,而是使用中心頻率v=600 MHz。對于沒有超過接收帶寬的暴,我們仍采用v=400 MHz。我們分析了Liso、Eiso、TB和紅移z的分布(見圖5)。

圖5 圖中展示了在靜止參考系中,Liso、Eiso、TB 和z 的分布;其中,實線階梯線和虛線階梯線分別表示重復暴和非重復暴,曲線則代表高斯擬合曲線
重復暴樣本和非重復暴樣本分別用實線和虛線表示。顯然,這些參數的分布在重復暴和非重復暴群組中基本一致,均服從對數正態分布,但總的來說,重復暴的四個參數分布范圍更廣。我們發現重復暴的Liso、Eiso和TB值平均小于非重復暴的值,并且數值下限更低。為了給出定量結果,我們進行了高斯擬合,結果列于表1。對于重復暴(非重復暴)樣本,Liso、Eiso、TB和紅移z的中值分別為8.47×1042(2.14×1043erg/s),σ~1.08(0.71);2.22×1040(5.23×1040erg),σ~0.97(0.73);1.08×1038(9.90×1037K),σ~1.54(1.05);和0.60(0.58),σ~0.32(0.31)。

表1 分布特性
2.1.2 相關性分析
為了探究識別出的重復暴、重復暴候選體和非重復暴之間參數分布的差異以及參數之間的關系,我們繪制了密度矩陣散點圖。如圖6 所示,研究的物理量是1.1 節介紹的反應快速射電暴主要特性的11 個參數。除Spectral Index 和Spectral Running 外,其余參數均取對數處理。同時,我們還討論了靜止坐標系中物理參數的關系(見圖7)。這些相同的參數也進行對數處理。

圖6 譜指數、光子流量、時間積分通量、箱式濾波器寬度、子暴寬度、散射時間、譜斜率變化率、最低頻率、最高頻率、峰值頻率和色散量各參數的散點密度矩陣

圖7 靜止坐標系下各向同性光度Liso、各向同性能量Eiso 和亮度溫度TB 的散點密度矩陣
從圖6 中可以清楚地看出,光子流量(Sv)與時間積分通量(Fobs)之間存在一定的相關性,散射時間(τ)與脈沖寬度(W)之間存在明顯的相關性。流量和時間積分通量還與脈沖寬度和散射時間有弱相關性,這將在下一小節中詳細討論。圖7 顯示了靜止系中Liso、Eiso和TB之間的明顯相關性。直觀上可以發現,Liso和Eiso、Liso和TB之間的相關性非常緊密,Liso和TB之間的相關性很強,但是彌散比較大。此外,從圖7 的分布可以看出,這些相關性對于重復暴和非重復暴都滿足,并且它們之間沒有顯著差異。因此,在2.2 節中,我們將重復暴和非重復暴作為一個整體來討論這種關系。
2.2.1 觀測者坐標系
為了探究識別出的重復暴在圖6 中我們發現快速射電暴部分物理量之間存在明顯依賴關系,包括log Svlog Fobs,log Sv-log W,log τ-log Sv,logτ-log Fobs,log τ-log W和log τ-log B W。為了定量描述它們之間的關系,我們對其中具有顯著相關性的兩參數進行了回歸分析,如圖8 所示。我們的模型是y=a+bx,其中y和x代表對數參數。對于沒有明顯線性趨勢的相關性,我們只給出兩個量之間的相關系數,不進行回歸分析。回歸分析結果如表2 所示。

表2 在觀測者坐標系下相關性回歸分析的結果,其中r 表示皮爾森相關系數,p 表示偶然概率,σ 表示彌散度

圖8 快速射電暴主要觀測參量間的關系圖,黑色實線表示對CHIME/FRB 目錄中的所有536 個快速射電暴進行的線性擬合,實心三角,空心圓點和十字號分別代表重復暴、非重復暴和重復暴候選體,從淺到深的灰色陰影區域分別表示1σ、2σ 和3σ 的置信區間
Sv和Fobs、W之間的分布關系在圖8 的頂部給出。我們發現重復暴和非重復暴具有相同的性質,因此我們將它們作為一個整體來進行定量分析并得到,皮爾森相關系數r=0.66,偶然概率p<10-4和彌散σ=0.33;Sv∝W-0.38(r=-0.36,p<10-4和σ=0.38)。此外,我們分析了散射時間τ與其他量之間的關系,并將結果顯示在表2 中。τ和Fobs之間存在弱相關性,r=0.21,但τ和Sv是反相關的,r=-0.41。有趣的是,τ和W、BW之間存在緊密相關性,它們分別遵循W∝τ0.47(r=0.63,p<10-4和σ=0.32)和BW∝τ0.59(r=0.70,p<10-4和σ=0.33)。
此外,我們也分別研究了重復暴和非重復暴之間的關系,但它們基本上遵循相似的關系,因此我們只給出它們的皮爾森相關系數、偶然概率和彌散值(見表2)。
2.2.2 靜止坐標系
靜止坐標系下的統計分析有助于研究快速射電暴的內在特性。Eiso-Liso和TB-Liso的相關性分析顯示在圖9 的頂部。此外,重復暴和非重復暴之間的相關性基本相同,因此我們也對其整體做回歸分析。我們發現Eiso和Liso之間存在非常緊密的相關性。確切地說,我們得到,其中皮爾森相關系數r=0.93,偶然概率p<10-4,彌散σ=0.33。TB和Liso之間也有很強的相關性,但是彌散比較大,具體為(r=0.77,p<10-4和σ=0.85)。然而,TB和Eiso之間的相關性較弱,且彌散較大,(r=0.58,p <10-4和σ=1.08)。

圖9 靜止坐標下Liso、Eiso 和TB 之間的雙參數關系圖,其他符號與圖8 相同
同樣,表3也給出了重復暴和非重復暴的具體回歸分析結果。需要注意的是,這些依賴關系在重復暴中表現得比非重復暴更加顯著,尤其是TB-Eiso關系,對于重復暴來說(r=0.68,p <10-4和σ=1.14),但在非重復暴中(r=0.35,p <10-4和σ=0.99)。重復Eiso暴的Eiso-Liso和TB-Liso的相關系數高于非重復暴,但沒有顯著差異。具體值見表3。

表3 在靜止坐標系下相關性回歸分析的結果,其中r 表示皮爾森相關系數,p 表示偶然概率,σ 表示彌散度
我們使用t-SNE 無監督降維算法對CHIME/FRB 樣本中的500 個非重復暴和94 個重復暴(包含子暴)進行分類,并通過重復暴和非重復暴觀測參數之間的隱藏特性來找到更多可能的重復暴候選者。無監督降維算法結果明顯顯示快速射電暴可分為兩個類別,分別對應于重復暴和非重復暴。重復暴候選體是那些通過機器學習算法分類后與重復暴緊密混合在一起的快速射電暴。我們通過這種方法在CHIME/FRB 目錄中識別出了165 個重復暴候選體,從而將重復暴的比例從15.8%提高到43.6%。這一結果表明,當前樣本中可能存在大量沒有被觀測發現或遺漏的重復暴。這些候選體目前屬于非重復暴類別,但它們表現出明顯的重復爆發特征,與重復暴具有高度相似性。
我們對分類后的快速射電暴進行了統計分析,發現CHIME/FRB 樣本中重復暴和非重復暴的光譜指數和譜斜率變化率值存在顯著差異。從光譜形狀參數化的等式來看,這意味著重復暴和非重復暴之間的光譜形狀存在很大差異。兩個群體之間的另一個顯著差異是重復暴具有更窄的爆發頻率帶寬,而非重復暴往往具有更大的帶寬。表明重復暴和非重復暴可能有不同的起源。
通過研究觀測者坐標系和靜止坐標系中的參數分布,我們發現重復暴和非重復暴的觀測參數分布范圍基本相同。特別是在爆發持續時間方面沒有發現明顯差異,這表明快速射電暴的持續時間可能不是區分兩者的代表性參量。在靜止坐標系中,重復暴的Liso、Eiso和TB的值普遍低于非重復暴的對應值。特別地,Liso的差異較為顯著,這表明重復暴與非重復暴可能具有不同的內在屬性或產生機制。此外,我們還發現重復暴和非重復暴基本遵循相同的兩參數依賴關系。
本文使用的CHIME/FRB 目錄是目前同一探測器在同一固定觀測頻帶中觀測到的第一個快速射電暴大樣本。研究結果為尋找更多的重復暴候選體提供了新的思路和方法,可以彌補尋找隱藏的重復暴的低效觀測方法,對于揭示重復暴和非重復暴之間的物理本質具有重要意義。這些結論還需要進一步的觀測來證實。