盒最小二乘(BLS)周期图#
“盒最小二乘法”(BLS)周期图 [1] 是一种统计工具,用于在时间序列光度数据中检测凌日系外行星和日食双星。此实现的主接口是 BoxLeastSquares
班级。
数学背景#
BLS方法通过将一个过境点建模为一个具有四个参数的周期性倒置的礼帽来寻找候选过境点:周期、持续时间、深度和参考时间。在这个实现中,参考时间被选为观测基线中第一次过境的中间过境时间。这些参数如下图所示:

假设测量磁通量的不确定度是已知的、独立的和高斯的,则传输通量的最大似然可以计算为
在哪里? \(y_n\) 是亮度测量值, \(\sigma_n\) 是相关的不确定性,两个总和都是在运输数据点上计算的。
同样地,最大可能的过境流量是
这些金额超过了过境观测。利用这些结果,我们得到了在给定时间段内运输模型(最大化深度)的对数似然 \(P\) ,持续时间 \(\tau\) ,和参考时间 \(t_0\) 是
这个等式可能很熟悉,因为它与“卡平方”成正比 \(\chi^2\) 对于这个模型,这是高斯不确定性假设的直接结果。
这个 \(\chi^2\) 被称为“信号残留物” [1], 因此,在持续时间和参考时间内最大化对数似然相当于从 [1].
在实践中,这是通过找到网格上持续时间和参考时间的最大似然模型来实现的 durations
和 oversample
的参数 power
方法。
基本用法#
经纬周期图以时间序列观测值作为输入,其中时间戳 t
以及观察结果 y
(通常为亮度)存储为 numpy
数组或 Quantity
物体。如果已知,误差线 dy
也可以选择提供。
例子#
要评估模拟数据集的周期图:
>>> import numpy as np
>>> import astropy.units as u
>>> from astropy.timeseries import BoxLeastSquares
>>> rng = np.random.default_rng(42)
>>> t = rng.uniform(0, 20, 2000)
>>> y = np.ones_like(t) - 0.1*((t%3)<0.2) + 0.01*rng.standard_normal(len(t))
>>> model = BoxLeastSquares(t * u.day, y, dy=0.01)
>>> periodogram = model.autopower(0.2)
的输出 astropy.timeseries.BoxLeastSquares.autopower
方法是 BoxLeastSquaresResults
通常情况下,最有用的属性是 period
和 power
属性。
可以使用matplotlib绘制此结果:
>>> import matplotlib.pyplot as plt
>>> plt.plot(periodogram.period, periodogram.power)

在这个图中,您可以看到峰值出现在正确的三天内。
目标#
默认情况下, power
方法计算模型拟合的对数似然,并在参考时间和持续时间内最大化。也可以使用测量渡越深度的信噪比作为目标函数。
例子#
要计算模型拟合的对数可能性,请调用 power
或 autopower
具有 objective='snr'
如下:
>>> model = BoxLeastSquares(t * u.day, y, dy=0.01)
>>> periodogram = model.autopower(0.2, objective="snr")

这个目标通常会产生一个周期图,在性质上类似于对数似然谱,但它已被用来提高在相关噪声存在下过境搜索的可靠性。
周期网格#
过境周期图总是在周期网格上计算的,其结果对采样很敏感。如中所述 [1], 经纬周期图法的性能对周期网格的敏感性要比周期图法高 LombScargle
周期图。
公交周期图的这个实现包括一个保守的启发式算法,用于估计由 autoperiod
和 autopower
方法和此方法的详细信息在 autoperiod
.
例子#
可以提供具体的时段网格,如下所示:
>>> model = BoxLeastSquares(t * u.day, y, dy=0.01)
>>> periods = np.linspace(2.5, 3.5, 1000) * u.day
>>> periodogram = model.power(periods, 0.2)

但是,如果时段网格过于粗糙,则可能会错过正确的时段。
>>> model = BoxLeastSquares(t * u.day, y, dy=0.01)
>>> periods = np.linspace(0.5, 10.5, 15) * u.day
>>> periodogram = model.power(periods, 0.2)

峰值统计#
为了帮助传输审核过程并调试候选峰值的问题,可以使用 compute_stats
该方法可用于计算一个候选公交的若干统计信息。
其中许多统计数据都基于中描述的VARTOOLS包 [2]. 这通常用于计算周期图中最大点的统计信息:
>>> model = BoxLeastSquares(t * u.day, y, dy=0.01)
>>> periodogram = model.autopower(0.2)
>>> max_power = np.argmax(periodogram.power)
>>> stats = model.compute_stats(periodogram.period[max_power],
... periodogram.duration[max_power],
... periodogram.transit_time[max_power])
这将计算包含有关该候选人的统计信息的词典。本词典中的每个条目都在 compute_stats
.