时间序列 (astropy.timeseries

介绍

从在固定时间对连续变量进行采样,到将事件合并到时间窗口中,天体物理学的许多不同领域都需要处理一维时间序列数据。为了满足这一需求 astropy.timeseries 子包提供表示和操作时间序列的类。

下面介绍的时间序列类是 QTable 具有特殊列的子类使用 Time 班级。因此,中描述的大部分功能 数据表 (astropy.table ) 适用于此。但是新类的主要目的是提供时间序列特定的功能 QTable .

入门

在本节中,我们将快速了解如何读取时间序列、访问数据并执行一些基本分析。有关创建和使用时间序列的更多详细信息,请参阅中的完整文档 使用 timeseries .

最基本的时间序列类是 TimeSeries -它将时间序列表示为特定时间点的值集合。如果您对将时间序列表示为离散时间容器中的度量值感兴趣,您可能会对 BinnedTimeSeries 我们在中显示的子类 使用 timeseries

首先,我们检索一个包含源的开普勒光曲线的FITS文件:

>>> from astropy.utils.data import get_pkg_data_filename
>>> filename = get_pkg_data_filename('timeseries/kplr010666592-2009131110544_slc.fits')  

注解

此处提供的光曲线是为示例目的而手动选择的。有关开普勒拟合格式的详细信息,请参见 Kepler Data Validation Document 还有开普勒科学中心 Light Curve Files 文档。要使用Python获取其他用于科学目的的开普勒光曲线,请参见 astroquery 附属套餐。

然后我们可以使用 TimeSeries 要在该文件中读取的类::

>>> from astropy.timeseries import TimeSeries
>>> ts = TimeSeries.read(filename, format='kepler.fits')  

时间序列是一种特殊的类型 Table 物体::

>>> ts  
<TimeSeries length=14280>
          time             timecorr   ...   pos_corr1      pos_corr2
                              d       ...      pix            pix
         object            float32    ...    float32        float32
----------------------- ------------- ... -------------- --------------
2009-05-02T00:41:40.338  6.630610e-04 ...  1.5822421e-03 -1.4463664e-03
2009-05-02T00:42:39.188  6.630857e-04 ...  1.5743829e-03 -1.4540013e-03
2009-05-02T00:43:38.045  6.631103e-04 ...  1.5665225e-03 -1.4616371e-03
2009-05-02T00:44:36.894  6.631350e-04 ...  1.5586632e-03 -1.4692718e-03
2009-05-02T00:45:35.752  6.631597e-04 ...  1.5508028e-03 -1.4769078e-03
2009-05-02T00:46:34.601  6.631844e-04 ...  1.5429436e-03 -1.4845425e-03
2009-05-02T00:47:33.451  6.632091e-04 ...  1.5350844e-03 -1.4921773e-03
2009-05-02T00:48:32.291  6.632337e-04 ...  1.5272264e-03 -1.4998110e-03
2009-05-02T00:49:31.149  6.632584e-04 ...  1.5193661e-03 -1.5074468e-03
                    ...           ... ...            ...            ...
2009-05-11T17:58:22.526  1.014493e-03 ...  3.6121816e-03  3.1950327e-03
2009-05-11T17:59:21.376  1.014518e-03 ...  3.6102540e-03  3.1872767e-03
2009-05-11T18:00:20.225  1.014542e-03 ...  3.6083264e-03  3.1795206e-03
2009-05-11T18:01:19.065  1.014567e-03 ...  3.6063993e-03  3.1717657e-03
2009-05-11T18:02:17.923  1.014591e-03 ...  3.6044715e-03  3.1640085e-03
2009-05-11T18:03:16.772  1.014615e-03 ...  3.6025438e-03  3.1562524e-03
2009-05-11T18:04:15.630  1.014640e-03 ...  3.6006160e-03  3.1484952e-03
2009-05-11T18:05:14.479  1.014664e-03 ...  3.5986886e-03  3.1407392e-03
2009-05-11T18:06:13.328  1.014689e-03 ...  3.5967610e-03  3.1329831e-03
2009-05-11T18:07:12.186  1.014713e-03 ...  3.5948332e-03  3.1252259e-03

以同样的方式 Table 对象的各个列和行 TimeSeries 可以使用索引符号访问和切片对象:

>>> ts['sap_flux']  
<Quantity [1027045.06, 1027184.44, 1027076.25, ..., 1025451.56, 1025468.5 ,
           1025930.9 ] electron / s>

>>> ts['time', 'sap_flux']  
<TimeSeries length=14280>
          time             sap_flux
                         electron / s
         object            float32
----------------------- --------------
2009-05-02T00:41:40.338  1.0270451e+06
2009-05-02T00:42:39.188  1.0271844e+06
2009-05-02T00:43:38.045  1.0270762e+06
2009-05-02T00:44:36.894  1.0271414e+06
2009-05-02T00:45:35.752  1.0271569e+06
2009-05-02T00:46:34.601  1.0272296e+06
2009-05-02T00:47:33.451  1.0273199e+06
2009-05-02T00:48:32.291  1.0271497e+06
2009-05-02T00:49:31.149  1.0271755e+06
                    ...            ...
2009-05-11T17:58:22.526  1.0234769e+06
2009-05-11T17:59:21.376  1.0234574e+06
2009-05-11T18:00:20.225  1.0238128e+06
2009-05-11T18:01:19.065  1.0243234e+06
2009-05-11T18:02:17.923  1.0244257e+06
2009-05-11T18:03:16.772  1.0248654e+06
2009-05-11T18:04:15.630  1.0250156e+06
2009-05-11T18:05:14.479  1.0254516e+06
2009-05-11T18:06:13.328  1.0254685e+06
2009-05-11T18:07:12.186  1.0259309e+06

>>> ts[0:4]  
<TimeSeries length=4>
          time             timecorr   ...   pos_corr1      pos_corr2
                              d       ...      pix            pix
         object            float32    ...    float32        float32
----------------------- ------------- ... -------------- --------------
2009-05-02T00:41:40.338  6.630610e-04 ...  1.5822421e-03 -1.4463664e-03
2009-05-02T00:42:39.188  6.630857e-04 ...  1.5743829e-03 -1.4540013e-03
2009-05-02T00:43:38.045  6.631103e-04 ...  1.5665225e-03 -1.4616371e-03
2009-05-02T00:44:36.894  6.631350e-04 ...  1.5586632e-03 -1.4692718e-03

如前面的例子所示, TimeSeries 对象有一个 time 列,它始终是第一列。也可以使用 .time 属性:

>>> ts.time  
<Time object: scale='tdb' format='isot' value=['2009-05-02T00:41:40.338' '2009-05-02T00:42:39.188'
  '2009-05-02T00:43:38.045' ... '2009-05-11T18:05:14.479'
  '2009-05-11T18:06:13.328' '2009-05-11T18:07:12.186']>

第一列总是 Time 对象(参见 Times and Dates ),因此支持转换为不同时间尺度和格式的能力:

>>> ts.time.mjd  
array([54953.0289391 , 54953.02962023, 54953.03030145, ...,
       54962.7536398 , 54962.75432093, 54962.75500215])

>>> ts.time.unix  
array([1.24122483e+09, 1.24122489e+09, 1.24122495e+09, ...,
       1.24206505e+09, 1.24206511e+09, 1.24206517e+09])

我们还可以检查时间定义的时间刻度:

>>> ts.time.scale  
'tdb'

这是重心动力学时间尺度(参见 时间和日期 (astropy.time ) 更多细节)。我们可以用我们目前所看到的来制作一个情节:

import matplotlib.pyplot as plt
plt.plot(ts.time.jd, ts['sap_flux'], 'k.', markersize=1)
plt.xlabel('Julian Date')
plt.ylabel('SAP Flux (e-/s)')

(png _, svgpdf

../_images/index-22.png

看起来好像有几个过境点!我们可以使用 BoxLeastSquares 类来使用“盒最小二乘法”(BLS)算法估计周期:

>>> import numpy as np
>>> from astropy import units as u
>>> from astropy.timeseries import BoxLeastSquares
>>> periodogram = BoxLeastSquares.from_timeseries(ts, 'sap_flux')  

要运行周期图分析,我们使用一个持续时间为0.2天的框:

>>> results = periodogram.autopower(0.2 * u.day)  
>>> best = np.argmax(results.power)  
>>> period = results.period[best]  
>>> period  
<Quantity 2.20551724 d>
>>> transit_time = results.transit_time[best]  
>>> transit_time  
<Time object: scale='tdb' format='isot' value=2009-05-02T20:51:16.338>

有关可用周期图算法的详细信息,请参阅 周期图算法 .

我们现在可以使用上面找到的时间段来折叠时间序列 fold() 方法:

>>> ts_folded = ts.fold(period=period, epoch_time=transit_time)  

现在我们可以看看折叠的时间序列:

plt.plot(ts_folded.time.jd, ts_folded['sap_flux'], 'k.', markersize=1)
plt.xlabel('Time (days)')
plt.ylabel('SAP Flux (e-/s)')

(png _, svgpdf

../_images/index-6.png

使用 天体统计工具 (astropy.stats ) 模块,我们可以通过sigma剪裁数据来规范化通量,以确定基线通量:

>>> from astropy.stats import sigma_clipped_stats
>>> mean, median, stddev = sigma_clipped_stats(ts_folded['sap_flux'])  
>>> ts_folded['sap_flux_norm'] = ts_folded['sap_flux'] / median  

我们可以通过将这些点划分为相等时间的单元来降低时间序列的采样率-这将返回一个 BinnedTimeSeries ::

>>> from astropy.timeseries import aggregate_downsample
>>> ts_binned = aggregate_downsample(ts_folded, time_bin_size=0.03 * u.day)  
>>> ts_binned  
<BinnedTimeSeries length=74>
   time_bin_start     time_bin_size    ...   sap_flux_norm
                            s          ...
       object            float64       ...       float64
------------------- ------------------ ... ------------------
-1.1022116370482966             2592.0 ... 0.9998741745948792
-1.0722116370482966             2592.0 ... 0.9999074339866638
-1.0422116370482966             2592.0 ...  0.999972939491272
-1.0122116370482965             2592.0 ... 1.0000077486038208
-0.9822116370482965             2592.0 ... 0.9999921917915344
-0.9522116370482965             2592.0 ... 1.0000101327896118
-0.9222116370482966             2592.0 ... 1.0000121593475342
-0.8922116370482965             2592.0 ... 0.9999905228614807
-0.8622116370482965 2592.0000000000023 ... 1.0000263452529907
                ...                ... ...                ...
 0.8177883629517035 2591.9999999999977 ... 1.0000624656677246
 0.8477883629517035 2592.0000000000014 ... 1.0000633001327515
 0.8777883629517035  2592.000000000019 ... 1.0000433921813965
 0.9077883629517037 2591.9999999999814 ...  1.000024676322937
 0.9377883629517034   2592.00000000002 ... 1.0000224113464355
 0.9677883629517037  2591.999999999981 ... 1.0000698566436768
 0.9977883629517035             2592.0 ... 0.9999606013298035
 1.0277883629517035             2592.0 ... 0.9999635815620422
 1.0577883629517035             2592.0 ... 0.9999105930328369
 1.0877883629517036 2592.0000000000095 ... 0.9998687505722046

现在我们可以看看最终结果:

plt.plot(ts_folded.time.jd, ts_folded['sap_flux_norm'], 'k.', markersize=1)
plt.plot(ts_binned.time_bin_start.jd, ts_binned['sap_flux_norm'], 'r-', drawstyle='steps-post')
plt.xlabel('Time (days)')
plt.ylabel('Normalized flux')

(png _, svgpdf

../_images/index-10.png

进一步了解 astropy.timeseries 模块,您可以在下一节中找到指向完整文档的链接。

参考/API

astropy.timeseries公司包裹

此子包包含用于处理时间序列的类和函数。

功能

aggregate_downsample \(时间系列), * [, ...] )

通过将值分到固定大小的存储单元中,使用单个函数组合存储单元中的值,对时间序列进行降采样。

autocheck_required_columns \(cls)

这是一个decorator,它确保表包含由u required_columns属性指示的特定方法。

Classes

BasePeriodogram \(t,y[, dy] )

BaseTimeSeries \ [data, masked, names, dtype, ...] )

BinnedTimeSeries \ [data, time_bin_start, ...] )

以表格形式表示二进制时间序列数据的类。

BoxLeastSquares \(t,y[, dy] )

计算盒最小二乘周期图

BoxLeastSquaresResults * ARGs)

无框正方形搜索的结果

LombScargle \(t,y[, dy, fit_mean, ...] )

计算Lomb擦伤周期图。

TimeSeries \ [data, time, time_start, ...] )

以表格形式表示时间序列数据的类。

类继承图

Inheritance diagram of astropy.timeseries.periodograms.base.BasePeriodogram, astropy.timeseries.core.BaseTimeSeries, astropy.timeseries.binned.BinnedTimeSeries, astropy.timeseries.periodograms.bls.core.BoxLeastSquares, astropy.timeseries.periodograms.bls.core.BoxLeastSquaresResults, astropy.timeseries.periodograms.lombscargle.core.LombScargle, astropy.timeseries.sampled.TimeSeries

astropy.timeseries.io包裹

功能

kepler_fits_reader \(文件名)

这可以作为astropy时间序列中开普勒或TESS文件的FITS阅读器。