lomb_scargle (pyleoclim.utils.spectral.lomb_scargle)
- pyleoclim.utils.spectral.lomb_scargle(ys, ts, freq=None, freq_method='lomb_scargle', freq_kwargs=None, n50=3, window='hann', detrend=None, sg_kwargs=None, gaussianize=False, standardize=False, average='mean')[source]
Return the computed periodogram using lomb-scargle algorithm
Uses the lombscargle implementation from scipy.signal: https://scipy.github.io/devdocs/generated/scipy.signal.lombscargle.html#scipy.signal.lombscargle
- Parameters
ys (array) – a time series
ts (array) – time axis of the time series
freq (str or array) – vector of frequency. If string, uses the following method:
freq_method (str) –
- Method to generate the frequency vector if not set directly. The following options are avialable:
log
lomb_scargle (default)
welch
scale
nfft
See utils.wavelet.make_freq_vector for details
freq_kwargs (dict) – Arguments for the method chosen in freq_method. See specific functions in utils.wavelet for details By default, uses dt=median(ts), ofac=4 and hifac=1 for Lomb-Scargle
n50 (int) – The number of 50% overlapping segment to apply
window (str or tuple) –
- Desired window to use. Possible values:
boxcar
triang
blackman
hamming
hann (default)
bartlett
flattop
parzen
bohman
blackmanharris
nuttail
barthann
kaiser (needs beta)
gaussian (needs standard deviation)
general_gaussian (needs power, width)
slepian (needs width)
dpss (needs normalized half-bandwidth)
chebwin (needs attenuation)
exponential (needs decay scale)
tukey (needs taper fraction)
If the window requires no parameters, then window can be a string. If the window requires parameters, then window must be a tuple with the first argument the string name of the window, and the next arguments the needed parameters. If window is a floating point number, it is interpreted as the beta parameter of the kaiser window.
- detrendstr
- If None, no detrending is applied. Available detrending methods:
None - no detrending will be applied (default);
linear - a linear least-squares fit to ys is subtracted;
constant - the mean of ys is subtracted
savitzy-golay - ys is filtered using the Savitzky-Golay filters and the resulting filtered series is subtracted from y.
emd - Empirical mode decomposition
- sg_kwargsdict
The parameters for the Savitzky-Golay filters. see pyleoclim.utils.filter.savitzy_golay for details.
- gaussianizebool
If True, gaussianizes the timeseries
- standardizebool
If True, standardizes the timeseriesprep_args : dict
- average{‘mean’,’median’}
Method to use when averaging periodograms. Defaults to ‘mean’.
- Returns
res_dict – the result dictionary, including - freq (array): the frequency vector - psd (array): the spectral density vector
- Return type
dict
See also
pyleoclim.utils.spectral.periodogram
Estimate power spectral density using a periodogram
pyleoclim.utils.spectral.mtm
Retuns spectral density using a multi-taper method
pyleoclim.utils.spectral.welch
Returns power spectral density using the Welch method
pyleoclim.utils.spectral.wwz_psd
Return the psd of a timeseries using wwz method.
pyleoclim.utils.filter.savitzy_golay
Filtering using Savitzy-Golay
pyleoclim.utils.tsutils.detrend
Detrending method
References
Lomb, N. R. (1976). Least-squares frequency analysis of unequally spaced data. Astrophysics and Space Science 39, 447-462.
Scargle, J. D. (1982). Studies in astronomical time series analysis. II. Statistical aspects of spectral analyis of unvenly spaced data. The Astrophysical Journal, 263(2), 835-853.
Scargle, J. D. (1982). Studies in astronomical time series analysis. II. Statistical aspects of spectral analyis of unvenly spaced data. The Astrophysical Journal, 263(2), 835-853.