LombScargleMultiband#

class astropy.timeseries.LombScargleMultiband(t, y, bands, dy=None, normalization='standard', nterms_base=1, nterms_band=1, reg_base=None, reg_band=1e-06, regularize_by_trace=True, center_data=True, fit_mean=True)[source]#

Bases: LombScargle

Compute the Lomb-Scargle Periodogram.

This implementation is based on code presented in [1] and [2]; if you use this functionality in an academic application, citation of those works would be appreciated.

Parameters:
tarray_like or Quantity [:ref: ‘time’]

sequence of observation times

yarray_like or Quantity

sequence of observations associated with times t

bandsarray_like

sequence of passband labels associated with times t, each unique label defines a single band of data.

dyfloat, array_like, or Quantity, optional

error or sequence of observational errors associated with times t

normalization{‘standard’, ‘model’, ‘log’, ‘psd’}, optional

Normalization to use for the periodogram.

nterms_baseint, optional

number of frequency terms to use for the base model common to all bands. In the case of the fast algorithm, this parameter is passed along to the single band LombScargle method as the nterms parameter.

nterms_bandint, optional

number of frequency terms to use for the residuals between the base model and each individual band

reg_basefloat or None (default = None)

amount of regularization to use on the base model parameters

reg_bandfloat or None (default = 1E-6)

amount of regularization to use on the band model parameters

regularize_by_tracebool (default = True)

if True, then regularization is expressed in units of the trace of the normal matrix

center_databool, optional

if True, pre-center the data by subtracting the weighted mean of the input data.

fit_meanbool, optional

if True, include a constant offset as part of the model at each frequency. This can lead to more accurate results, especially in the case of incomplete phase coverage. Only applicable to the “fast” method

References

[1]

Vanderplas, J., Connolly, A. Ivezic, Z. & Gray, A. Introduction to astroML: Machine learning for astrophysics. Proceedings of the Conference on Intelligent Data Understanding (2012)

[2]

VanderPlas, J. & Ivezic, Z. Periodograms for Multiband Astronomical Time Series. ApJ 812.1:18 (2015)

Attributes Summary

available_methods

Methods Summary

autofrequency([samples_per_peak, ...])

Determine a suitable frequency grid for data.

autopower([method, sb_method, ...])

Compute Lomb-Scargle power at automatically-determined frequencies.

design_matrix(frequency[, t_fit, bands_fit])

Compute the design matrix for a given frequency

distribution()

Not Implemented

false_alarm_level()

Not Implemented

false_alarm_probability()

Not Implemented

from_timeseries(timeseries[, signal_column, ...])

Initialize a multiband periodogram from a time series object.

model(t, frequency[, bands_fit])

Compute the Lomb-Scargle model at the given frequency.

model_parameters(frequency[, units])

Compute the best-fit model parameters at the given frequency.

offset([t_fit, bands_fit])

Return the offset array of the model

power(frequency[, method, sb_method, ...])

Compute the Lomb-Scargle power at the given frequencies.

Attributes Documentation

available_methods = ['fast', 'flexible']#

Methods Documentation

autofrequency(samples_per_peak=5, nyquist_factor=5, minimum_frequency=None, maximum_frequency=None, return_freq_limits=False)[source]#

Determine a suitable frequency grid for data. Note that this assumes the peak width is driven by the observational baseline, which is generally a good assumption when the baseline is much larger than the oscillation period. If you are searching for periods longer than the baseline of your observations, this may not perform well. Even with a large baseline, be aware that the maximum frequency returned is based on the concept of “average Nyquist frequency”, which may not be useful for irregularly-sampled data. The maximum frequency can be adjusted via the nyquist_factor argument, or through the maximum_frequency argument.

Parameters:
samples_per_peakfloat, optional

The approximate number of desired samples across the typical peak

nyquist_factorfloat, optional

The multiple of the average nyquist frequency used to choose the maximum frequency if maximum_frequency is not provided.

minimum_frequencyfloat, optional

If specified, then use this minimum frequency rather than one chosen based on the size of the baseline.

maximum_frequencyfloat, optional

If specified, then use this maximum frequency rather than one chosen based on the average nyquist frequency.

return_freq_limitsbool, optional

if True, return only the frequency limits rather than the full frequency grid.

Returns:
frequencyndarray or Quantity [:ref: ‘frequency’]

The heuristically-determined optimal frequency bin

autopower(method='flexible', sb_method='auto', normalization='standard', samples_per_peak=5, nyquist_factor=5, minimum_frequency=None, maximum_frequency=None)[source]#

Compute Lomb-Scargle power at automatically-determined frequencies.

Parameters:
methodstr, optional

specify the multi-band lomb scargle implementation to use. Options are:

  • ‘flexible’: Constructs a common model, and an offset model per individual

    band. Applies regularization to the resulting model to constrain complexity.

  • ‘fast’: Passes each band individually through LombScargle (single-band),

    combines periodograms at the end by weight. Speed depends on single-band method chosen in ‘sb_method’.

sb_methodstr, optional

specify the single-band lomb scargle implementation to use, only in the case of using the ‘fast’ multiband method. Options are:

  • ‘auto’: choose the best method based on the input

  • ‘fast’: use the O[N log N] fast method. Note that this requires evenly-spaced frequencies: by default this will be checked unless assume_regular_frequency is set to True.

  • ‘slow’: use the O[N^2] pure-python implementation

  • ‘cython’: use the O[N^2] cython implementation. This is slightly faster than method=’slow’, but much more memory efficient.

  • ‘chi2’: use the O[N^2] chi2/linear-fitting implementation

  • ‘fastchi2’: use the O[N log N] chi2 implementation. Note that this requires evenly-spaced frequencies: by default this will be checked unless assume_regular_frequency is set to True.

  • ‘scipy’: use scipy.signal.lombscargle, which is an O[N^2] implementation written in C. Note that this does not support heteroskedastic errors.

normalization{‘standard’, ‘model’, ‘log’, ‘psd’}, optional

If specified, override the normalization specified at instantiation.

samples_per_peakfloat, optional

The approximate number of desired samples across the typical peak

nyquist_factorfloat, optional

The multiple of the average nyquist frequency used to choose the maximum frequency if maximum_frequency is not provided.

minimum_frequencyfloat or Quantity [:ref: ‘frequency’], optional

If specified, then use this minimum frequency rather than one chosen based on the size of the baseline. Should be Quantity if inputs to LombScargle are Quantity.

maximum_frequencyfloat or Quantity [:ref: ‘frequency’], optional

If specified, then use this maximum frequency rather than one chosen based on the average nyquist frequency. Should be Quantity if inputs to LombScargle are Quantity.

Returns:
frequency, powerndarray

The frequency and Lomb-Scargle power

design_matrix(frequency, t_fit=None, bands_fit=None)[source]#

Compute the design matrix for a given frequency

Parameters:
frequencyfloat

the frequency for the model

t_fitarray_like, Quantity, or Time (optional)

Times (length n_samples) at which to compute the model. If not specified, then the times and uncertainties of the input data are used.

bands_fitarray_like, or str

Bands to use in fitting, must be subset of bands in input data.

Returns:
ndarray

The design matrix for the model at the given frequency. This should have a shape of (len(t), n_parameters).

distribution()[source]#

Not Implemented

false_alarm_level()[source]#

Not Implemented

false_alarm_probability()[source]#

Not Implemented

classmethod from_timeseries(timeseries, signal_column=None, uncertainty_column=None, band_labels=None, **kwargs)[source]#

Initialize a multiband periodogram from a time series object.

If a binned time series is passed, the time at the center of the bins is used. Also note that this method automatically gets rid of NaN/undefined values when initializing the periodogram.

Parameters:
signal_columnlist

The names of columns containing the signal values to use.

uncertainty_columnlist, optional

The names of columns containing the errors on the signal.

band_labelslist, optional

The labels for each band, matched by index. If none, uses the labels of signal_column as band names.

**kwargs

Additional keyword arguments are passed to the initializer for this periodogram class.

model(t, frequency, bands_fit=None)[source]#

Compute the Lomb-Scargle model at the given frequency.

The model at a particular frequency is a linear model: model = offset + dot(design_matrix, model_parameters)

Parameters:
tarray_like or Quantity [:ref: ‘time’]

Times (length n_samples) at which to compute the model.

frequencyfloat

the frequency for the model

bands_fitlist or array_like

the unique bands to fit in the model

Returns:
ynp.ndarray

The model fit corresponding to the input times. Will have shape (n_bands,``n_samples``).

model_parameters(frequency, units=True)[source]#

Compute the best-fit model parameters at the given frequency.

The model described by these parameters is:

\[y(t; f, \vec{\theta}) = \theta_0 + \sum_{n=1}^{\tt nterms_base} [\theta_{2n-1}\sin(2\pi n f t) + \theta_{2n}\cos(2\pi n f t)] + \theta_0^{(k)} + \sum_{n=1}^{\tt nterms_band} [\theta_{2n-1}^{(k)}\sin(2\pi n f t) +\]

where \(\vec{\theta}\) is the array of parameters returned by this function.

Parameters:
frequencyfloat

the frequency for the model

unitsbool

If True (default), return design matrix with data units.

Returns:
thetanp.ndarray (n_parameters,)

The best-fit model parameters at the given frequency.

offset(t_fit=None, bands_fit=None)[source]#

Return the offset array of the model

The offset array of the model is the (weighted) mean of the y values in each band. Note that if self.center_data is False, the offset is 0 by definition.

Parameters:
t_fitarray_like, Quantity, or Time (optional)

Times (length n_samples) at which to compute the model. If not specified, then the times and uncertainties of the input data are used.

bands_fitarray_like, or str

Bands to use in fitting, must be subset of bands in input data.

Returns:
offsetarray
power(frequency, method='flexible', sb_method='auto', normalization='standard')[source]#

Compute the Lomb-Scargle power at the given frequencies.

Parameters:
frequencyarray_like or Quantity [:ref: ‘frequency’]

frequencies (not angular frequencies) at which to evaluate the periodogram. Note that in order to use method=’fast’, frequencies must be regularly-spaced.

methodstr, optional

specify the multi-band lomb scargle implementation to use. Options are:

  • ‘flexible’: Constructs a common model, and an offset model per individual

    band. Applies regularization to the resulting model to constrain complexity.

  • ‘fast’: Passes each band individually through LombScargle (single-band),

    combines periodograms at the end by weight. Speed depends on single-band method chosen in ‘sb_method’.

sb_methodstr, optional

specify the single-band lomb scargle implementation to use, only in the case of using the ‘fast’ multiband method. Options can be found in LombScargle.

normalization{‘standard’, ‘model’, ‘log’, ‘psd’}, optional

If specified, override the normalization specified at instantiation.

Returns:
powerndarray

The Lomb-Scargle power at the specified frequency