基于机器学习方法的快速射电暴分类研究＊

孙万鹏，张富文

（桂林理工大学理学院，广西 桂林 541004）

0 引言

自2007 年快速射电暴被首次发现以来，该领域已经迅速发展。快速射电暴被认为是宇宙中最神秘的天文现象之一，其强大的脉冲信号具有极高的峰值亮度和短暂的持续时间。虽然人们已经对快速射电暴的基本特征有了一定的了解，但仍有许多关键问题需要解决，其中一个便是快速射电暴的分类问题。得益于CHIME团队提供的第一个快速射电暴大样本，目前已经有部分研究者通过快速射电暴的关键观测参量之间的差异对其进行分类[1-2]，还有其他研究者根据其是否重复爆发的特性将其分为重复暴和非重复暴进行讨论[3-4]。然而，当前样本中重复暴和非重复暴的比例与理论预测结果严重不符，这意味着可能有大量重复暴未被观测到而被错误地分类为非重复暴群体[5]。因此，寻找快速射电暴中潜在的重复暴，实现对快速射电暴清晰地分类尤为重要。

现有的大多数研究往往只是简单地根据某一个或两个观测参数的差异来判断重复暴和非重复暴是否属于不同类别，很少有研究能够综合考虑多个观测性质的差异，并提出一般的区分方法。而无监督机器学习擅长从输入数据关系中揭示隐藏信息，因此利用这一方法可以有效地帮助我们结合快速射电暴的多维度特征对其进行分类，并发现潜在的重复暴。

在本文中，我们采用了一种名为t-SNE 的无监督机器学习算法，用于将CHIME/FRB 目录中的重复暴和非重复暴进行分类。通过这种方法，我们找到了在当前CHIME/FRB 样本中潜在的重复暴候选体，并提供了一个“未受污染”的样本，可以清晰地区分重复暴和非重复暴。基于分类的样本，我们进行了更进一步的统计分析，全面讨论了重复暴与非重复暴的参数分布以及一些参数之间的关系。基于重新分类后的统计分析有助于我们揭示重复暴与非重复暴背后不同的物理机制。同时，我们的研究结果为进一步搜寻重复暴和科学地分类快速射电暴提供了指导。

本文的第一章介+绍了样本选择和数据处理，并简要介绍了用于无监督算法的特征参数。第二章给出了t-SNE 算法的介绍和主要参数配置。第三章给出了CHIME/FRB 样本无监督分类的结果，包括识别出的重复暴候选体，同时对比分析了重复暴和非重复暴的各观测参数分布以及快速射电暴在观测者坐标系和静止坐标系中的统计特性。在第四章中总结了我们的工作。

1 样本选择与数据处理

1.1 样本选择

我们选取了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)：探测器观测所得到的色散量，对于每个子暴来说这些值都是相同的。

1.2 t-SNE 算法介绍

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 中查看。

2 结果分析

2.1 无监督识别下的快速射电暴类别

图图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 快速射电暴在观测者/静止坐标系下的统计特性

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 表示偶然概率，σ 表示弥散度

3 结论

我们使用t-SNE 无监督降维算法对CHIME/FRB 样本中的500 个非重复暴和94 个重复暴（包含子暴）进行分类，并通过重复暴和非重复暴观测参数之间的隐藏特性来找到更多可能的重复暴候选者。无监督降维算法结果明显显示快速射电暴可分为两个类别，分别对应于重复暴和非重复暴。重复暴候选体是那些通过机器学习算法分类后与重复暴紧密混合在一起的快速射电暴。我们通过这种方法在CHIME/FRB 目录中识别出了165 个重复暴候选体，从而将重复暴的比例从15.8%提高到43.6%。这一结果表明，当前样本中可能存在大量没有被观测发现或遗漏的重复暴。这些候选体目前属于非重复暴类别，但它们表现出明显的重复爆发特征，与重复暴具有高度相似性。

我们对分类后的快速射电暴进行了统计分析，发现CHIME/FRB 样本中重复暴和非重复暴的光谱指数和谱斜率变化率值存在显著差异。从光谱形状参数化的等式来看，这意味着重复暴和非重复暴之间的光谱形状存在很大差异。两个群体之间的另一个显著差异是重复暴具有更窄的爆发频率带宽，而非重复暴往往具有更大的带宽。表明重复暴和非重复暴可能有不同的起源。

通过研究观测者坐标系和静止坐标系中的参数分布，我们发现重复暴和非重复暴的观测参数分布范围基本相同。特别是在爆发持续时间方面没有发现明显差异，这表明快速射电暴的持续时间可能不是区分两者的代表性参量。在静止坐标系中，重复暴的Liso、Eiso和TB的值普遍低于非重复暴的对应值。特别地，Liso的差异较为显著，这表明重复暴与非重复暴可能具有不同的内在属性或产生机制。此外，我们还发现重复暴和非重复暴基本遵循相同的两参数依赖关系。

本文使用的CHIME/FRB 目录是目前同一探测器在同一固定观测频带中观测到的第一个快速射电暴大样本。研究结果为寻找更多的重复暴候选体提供了新的思路和方法，可以弥补寻找隐藏的重复暴的低效观测方法，对于揭示重复暴和非重复暴之间的物理本质具有重要意义。这些结论还需要进一步的观测来证实。