1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301
|
.. _timeseries-analysis:
Manipulation and Analysis of Time Series
****************************************
Combining Time Series
=====================
The :func:`~astropy.table.vstack` and :func:`~astropy.table.hstack` functions
from the :mod:`astropy.table` module can be used to stack time series in
different ways.
Examples
--------
.. EXAMPLE START: Stacking Time Series Row-Wise Using table.vstack
Time series can be stacked "vertically" or row-wise using the
:func:`~astropy.table.vstack` function (although note that sampled time
series cannot be combined with binned time series and vice versa)::
>>> from astropy.table import vstack
>>> from astropy import units as u
>>> from astropy.timeseries import TimeSeries
>>> ts_a = TimeSeries(time_start='2016-03-22T12:30:31',
... time_delta=3 * u.s,
... data={'flux': [1, 4, 5, 3, 2] * u.mJy})
>>> ts_b = TimeSeries(time_start='2016-03-22T12:50:31',
... time_delta=3 * u.s,
... data={'flux': [4, 3, 1, 2, 3] * u.mJy})
>>> ts_ab = vstack([ts_a, ts_b])
>>> ts_ab
<TimeSeries length=10>
time flux
mJy
Time float64
----------------------- -------
2016-03-22T12:30:31.000 1.0
2016-03-22T12:30:34.000 4.0
2016-03-22T12:30:37.000 5.0
2016-03-22T12:30:40.000 3.0
2016-03-22T12:30:43.000 2.0
2016-03-22T12:50:31.000 4.0
2016-03-22T12:50:34.000 3.0
2016-03-22T12:50:37.000 1.0
2016-03-22T12:50:40.000 2.0
2016-03-22T12:50:43.000 3.0
Note that :func:`~astropy.table.vstack` does not automatically sort, nor get rid
of duplicates — this is something you would need to do explicitly afterwards.
.. EXAMPLE END
.. EXAMPLE START: Stacking Time Series Column-Wise Using table.vstack
Time series can also be combined "horizontally" or column-wise with other tables
using the :func:`~astropy.table.hstack` function, though these should not be
time series (as having multiple time columns would be confusing)::
>>> from astropy.table import Table, hstack
>>> data = Table(data={'temperature': [40., 41., 40., 39., 30.] * u.K})
>>> ts_a_data = hstack([ts_a, data])
>>> ts_a_data
<TimeSeries length=5>
time flux temperature
mJy K
Time float64 float64
----------------------- ------- -----------
2016-03-22T12:30:31.000 1.0 40.0
2016-03-22T12:30:34.000 4.0 41.0
2016-03-22T12:30:37.000 5.0 40.0
2016-03-22T12:30:40.000 3.0 39.0
2016-03-22T12:30:43.000 2.0 30.0
.. EXAMPLE END
Sorting Time Series
===================
.. EXAMPLE START: Sorting Time Series
Sorting time series in place can be done using the
:meth:`~astropy.table.Table.sort` method, as for |Table|::
>>> ts = TimeSeries(time_start='2016-03-22T12:30:31',
... time_delta=3 * u.s,
... data={'flux': [1., 4., 5., 3., 2.]})
>>> ts
<TimeSeries length=5>
time flux
Time float64
----------------------- -------
2016-03-22T12:30:31.000 1.0
2016-03-22T12:30:34.000 4.0
2016-03-22T12:30:37.000 5.0
2016-03-22T12:30:40.000 3.0
2016-03-22T12:30:43.000 2.0
>>> ts.sort('flux')
>>> ts
<TimeSeries length=5>
time flux
Time float64
----------------------- -------
2016-03-22T12:30:31.000 1.0
2016-03-22T12:30:43.000 2.0
2016-03-22T12:30:40.000 3.0
2016-03-22T12:30:34.000 4.0
2016-03-22T12:30:37.000 5.0
.. EXAMPLE END
Resampling
==========
We provide a :func:`~astropy.timeseries.aggregate_downsample` function
that can be used to bin values from a time series into equal-size or uneven bins,
and contiguous and non-contiguous bins, using a custom function (mean, median, etc.).
This operation returns a |BinnedTimeSeries|. Note that this is a basic function in
the sense that it does not, for example, know how to treat columns with uncertainties
differently from other values, and it will blindly apply the custom function
specified to all columns.
Example
-------
.. EXAMPLE START: Creating a BinnedTimeSeries with even contiguous bins
The following example shows how to use
:func:`~astropy.timeseries.aggregate_downsample` to bin a light curve from the
Kepler mission into 20 minute contiguous bins using a median function. First,
we read in the data using:
.. plot::
:include-source:
:context: reset
:nofigs:
from astropy.timeseries import TimeSeries
from astropy.utils.data import get_pkg_data_filename
example_data = get_pkg_data_filename('timeseries/kplr010666592-2009131110544_slc.fits')
kepler = TimeSeries.read(example_data, format='kepler.fits')
(See :ref:`timeseries-io` for more details about reading in data). We can then
downsample using:
.. plot::
:context:
:nofigs:
import warnings
warnings.filterwarnings('ignore', message='All-NaN slice encountered')
.. plot::
:include-source:
:context:
:nofigs:
import numpy as np
from astropy import units as u
from astropy.timeseries import aggregate_downsample
kepler_binned = aggregate_downsample(kepler, time_bin_size=20 * u.min, aggregate_func=np.nanmedian)
We can take a look at the results:
.. plot::
:include-source:
:context:
import matplotlib.pyplot as plt
plt.plot(kepler.time.jd, kepler['sap_flux'], 'k.', markersize=1)
plt.plot(kepler_binned.time_bin_start.jd, kepler_binned['sap_flux'], 'r-', drawstyle='steps-pre')
plt.xlabel('Julian Date')
plt.ylabel('SAP Flux (e-/s)')
.. EXAMPLE END
.. EXAMPLE START: Creating a BinnedTimeSeries with uneven contiguous bins
The :func:`~astropy.timeseries.aggregate_downsample` can also be used
to bin the light curve into custom bins. The following example shows
the case of uneven-size contiguous bins:
.. plot::
:context: reset
:nofigs:
import numpy as np
from astropy import units as u
import matplotlib.pyplot as plt
from astropy.timeseries import TimeSeries
from astropy.timeseries import aggregate_downsample
from astropy.utils.data import get_pkg_data_filename
example_data = get_pkg_data_filename('timeseries/kplr010666592-2009131110544_slc.fits')
kepler = TimeSeries.read(example_data, format='kepler.fits')
import warnings
warnings.filterwarnings('ignore', message='All-NaN slice encountered')
.. plot::
:include-source:
:context:
kepler_binned = aggregate_downsample(kepler, time_bin_size=[1000, 125, 80, 25, 150, 210, 273] * u.min,
aggregate_func=np.nanmedian)
plt.plot(kepler.time.jd, kepler['sap_flux'], 'k.', markersize=1)
plt.plot(kepler_binned.time_bin_start.jd, kepler_binned['sap_flux'], 'r-', drawstyle='steps-pre')
plt.xlabel('Julian Date')
plt.ylabel('SAP Flux (e-/s)')
To learn more about the custom binning functionality in
:func:`~astropy.timeseries.aggregate_downsample`, see
:ref:`timeseries-binned-initializing`.
Folding
=======
.. EXAMPLE START: Phase Folding a Time Series
The |TimeSeries| class has a
:meth:`~astropy.timeseries.TimeSeries.fold` method that can be used to
return a new time series with a relative and folded time axis. This method
takes the period as a :class:`~astropy.units.Quantity`, and optionally takes
an epoch as a :class:`~astropy.time.Time`, which defines a zero time offset:
.. plot::
:context: reset
:nofigs:
import numpy as np
from astropy import units as u
import matplotlib.pyplot as plt
from astropy.timeseries import TimeSeries
from astropy.utils.data import get_pkg_data_filename
example_data = get_pkg_data_filename('timeseries/kplr010666592-2009131110544_slc.fits')
kepler = TimeSeries.read(example_data, format='kepler.fits')
.. plot::
:include-source:
:context:
kepler_folded = kepler.fold(period=2.2 * u.day, epoch_time='2009-05-02T20:53:40')
plt.plot(kepler_folded.time.jd, kepler_folded['sap_flux'], 'k.', markersize=1)
plt.xlabel('Time from midpoint epoch (days)')
plt.ylabel('SAP Flux (e-/s)')
Note that in this example we happened to know the period and midpoint from a
previous periodogram analysis. See the example in :doc:`index` for how you
might do this.
.. EXAMPLE END
Arithmetic
==========
.. EXAMPLE START: Arithmetic with Time Series
Since |TimeSeries| objects are subclasses of |Table|, they naturally support
arithmetic on any of the data columns. As an example, we can take the folded
Kepler time series we have seen in previous examples, and normalize it to the
sigma-clipped median value.
.. plot::
:context: reset
:nofigs:
import numpy as np
from astropy import units as u
import matplotlib.pyplot as plt
from astropy.timeseries import TimeSeries
from astropy.utils.data import get_pkg_data_filename
example_data = get_pkg_data_filename('timeseries/kplr010666592-2009131110544_slc.fits')
kepler = TimeSeries.read(example_data, format='kepler.fits')
kepler_folded = kepler.fold(period=2.2 * u.day, epoch_time='2009-05-02T20:53:40')
.. plot::
:context:
:nofigs:
import warnings
warnings.filterwarnings('ignore', message='Input data contains invalid values')
.. plot::
:include-source:
:context:
from astropy.stats import sigma_clipped_stats
mean, median, stddev = sigma_clipped_stats(kepler_folded['sap_flux'])
kepler_folded['sap_flux_norm'] = kepler_folded['sap_flux'] / median
plt.plot(kepler_folded.time.jd, kepler_folded['sap_flux_norm'], 'k.', markersize=1)
plt.xlabel('Time from midpoint epoch (days)')
plt.ylabel('Normalized flux')
.. EXAMPLE END
|