freqanalysis#

syncopy.freqanalysis(data, method='mtmfft', output='pow', keeptrials=True, foi=None, foilim=None, pad='maxperlen', polyremoval=0, taper='hann', demean_taper=False, taper_opt=None, tapsmofrq=None, nTaper=None, keeptapers=False, toi='all', t_ftimwin=None, wavelet='Morlet', width=6, order=None, order_max=None, order_min=1, c_1=3, adaptive=False, out=None, fooof_opt=None, ft_compat=False, parallel=None, select=None, **kwargs)[source]#

Perform (time-)frequency analysis of Syncopy AnalogData objects

Usage Summary

Options available in all analysis methods:

  • output : one of availableOutputs; return power spectra, complex Fourier spectra or absolute values.

  • foi/foilim : frequencies of interest; either array of frequencies or frequency window (not both)

  • keeptrials : return individual trials or grand average

  • polyremoval : de-trending method to use (0 = mean, 1 = linear or None)

List of available analysis methods and respective distinct options:

“mtmfft”(Multi-)tapered Fourier transform

Perform frequency analysis on time-series trial data using either a single taper window (Hanning) or many tapers based on the discrete prolate spheroidal sequence (DPSS) that maximize energy concentration in the main lobe.

  • taper : one of availableTapers

  • tapsmofrq : spectral smoothing box for slepian tapers (in Hz)

  • nTaper : number of orthogonal tapers for slepian tapers

  • keeptapers : return individual tapers or average

  • pad: either pad to an absolute length or set to ‘nextpow2’

Post-processing of the resulting spectra with FOOOF is available via setting output to one of ‘fooof’, ‘fooof_aperiodic’ or ‘fooof_peaks’, see below for details. The returned spectrum represents the full fooofed spectrum for ‘fooof’, the aperiodic fit for ‘fooof_aperiodic’, and the peaks (Gaussians fit to them) for ‘fooof_peaks’. Returned data is in linear scale. Noisy input data will most likely lead to fitting issues with fooof, always inspect your results!

“mtmconvol”(Multi-)tapered sliding window Fourier transform

Perform time-frequency analysis on time-series trial data based on a sliding window short-time Fourier transform using either a single Hanning taper or multiple DPSS tapers.

  • taper : one of availableTapers

  • tapsmofrq : spectral smoothing box for slepian tapers (in Hz)

  • nTaper : number of orthogonal tapers for slepian tapers

  • keeptapers : return individual tapers or average

  • toi : time-points of interest; can be either an array representing analysis window centroids (in sec), a scalar between 0 and 1 encoding the percentage of overlap between adjacent windows or “all” to center a window on every sample in the data.

  • t_ftimwin : sliding window length (in sec)

“welch”Welch’s method for the estimation of power spectra based on

time-averaging over short, modified periodograms. Here, modified means that a taper is applied. See [Welch1967] for details.

  • toi : time-points of interest; a scalar between 0 and 1 encoding the percentage of overlap between adjacent windows.

  • t_ftimwin : sliding window length (in sec)

  • taper : one of availableTapers

  • tapsmofrq : spectral smoothing box for slepian tapers (in Hz)

  • nTaper : number of orthogonal tapers for slepian tapers

  • keeptapers : must be False with Welch. For multi-tapering, taper averaging happens as part of the modified periodogram computation,

    i.e., before the window averaging performed by Welch.

“wavelet”(Continuous non-orthogonal) wavelet transform

Perform time-frequency analysis on time-series trial data using a non-orthogonal continuous wavelet transform.

  • wavelet : one of availableWavelets

  • toi : time-points of interest; can be either an array representing time points (in sec) or “all”(pre-trimming and subsampling of results)

  • width : Nondimensional frequency constant of Morlet wavelet function (>= 6)

  • order : Order of Paul wavelet function (>= 4) or derivative order of real-valued DOG wavelets (2 = mexican hat)

“superlet”Superlet transform

Perform time-frequency analysis on time-series trial data using the super-resolution superlet transform (SLT) from [Moca2021].

  • order_max : Maximal order of the superlet

  • order_min : Minimal order of the superlet

  • c_1 : Number of cycles of the base Morlet wavelet

  • adaptive : If set to True perform fractional adaptive SLT, otherwise perform multiplicative SLT

Full documentation below

The parameters listed below can be provided as is or a via a cfg configuration ‘structure’, see Notes for details.

Parameters:
  • data (~syncopy.AnalogData) – A non-empty Syncopy AnalogData object

  • method (str) – Spectral estimation method, one of availableMethods (see below).

  • output (str) – Output of spectral estimation. One of availableOutputs (see below); use ‘pow’ for power spectrum (numpy.float32), ‘fourier’ for complex Fourier coefficients (numpy.complex64) or ‘abs’ for absolute values (numpy.float32). Use one of ‘fooof’, ‘fooof_aperiodic’ or ‘fooof_peaks’ to request post-processing of the results with FOOOF, also see the ‘fooof_opt’ parameter description.

  • keeptrials (bool) – If True spectral estimates of individual trials are returned, otherwise results are averaged across trials.

  • foi (array-like or None) – Frequencies of interest (Hz) for output. If desired frequencies cannot be matched exactly, the closest possible frequencies are used. If foi is None or foi = "all", all attainable frequencies (i.e., zero to Nyquist / 2) are selected.

  • foilim (array-like (floats [fmin, fmax]) or None or "all") – Frequency-window [fmin, fmax] (in Hz) of interest. Window specifications must be sorted (e.g., [90, 70] is invalid) and not NaN but may be unbounded (e.g., [-np.inf, 60.5] is valid). Edges fmin and fmax are included in the selection. If foilim is None or foilim = "all", all frequencies are selected.

  • pad ('maxperlen', float or 'nextpow2') – For the default maxperlen, no padding is performed in case of equal length trials, while trials of varying lengths are padded to match the longest trial. If pad is a number all trials are padded so that pad indicates the absolute length of all trials after padding in seconds. For instance pad = 2 pads all trials to an absolute length of 2000 samples, if and only if the longest trial contains at maximum 2000 samples and the samplerate is 1kHz. If pad is ‘nextpow2’ all trials are padded to the nearest power of two (in samples) of the longest trial.

  • polyremoval (int or None) – Order of polynomial used for de-trending data in the time domain prior to spectral analysis. A value of 0 corresponds to subtracting the mean (“de-meaning”), polyremoval = 1 removes linear trends (subtracting the least squares fit of a linear polynomial). If polyremoval is None, no de-trending is performed. Note that for spectral estimation de-meaning is very advisable and hence also the default.

  • tapsmofrq (float or None) – Only valid if method is ‘mtmfft’ or ‘mtmconvol’ Enables multi-tapering and sets the amount of one-sided spectral smoothing with slepian tapers in Hz.

  • nTaper (int or None) – Only valid if method is ‘mtmfft’ or ‘mtmconvol’ and tapsmofrq is set. Number of orthogonal tapers to use for multi-tapering. It is not recommended to set the number of tapers manually! Leave at None for the optimal number to be set automatically.

  • taper (str or None, optional) – Only valid if method is 'mtmfft' or 'mtmconvol'. Windowing function, one of availableTapers For multi-tapering with slepian tapers use tapsmofrq directly.

  • demean_taper (bool) – Set to True to perform de-meaning after tapering. Recommended for later Granger analysis with connectivityanalysis(). Only valid for method='mtmfft'.

  • taper_opt (dict or None) – Dictionary with keys for additional taper parameters. For example kaiser() has the additional parameter ‘beta’. For multi-tapering use tapsmofrq directly.

  • keeptapers (bool) – Only valid if method is ‘mtmfft’ or ‘mtmconvol’ and multi-tapering enabled via setting tapsmofrq. If True, return spectral estimates for each taper. Otherwise power spectrum is averaged across tapers, if and only if output is pow.

  • toi (float or array-like or "all") – Mandatory input for time-frequency analysis methods (method is either “mtmconvol” or “wavelet” or “superlet”). If toi is scalar, it must be a value between 0 and 1 indicating the percentage of overlap between time-windows specified by t_ftimwin (only valid if method is ‘mtmconvol’). If toi is an array it explicitly selects the centroids of analysis windows (in seconds), if toi is “all”, analysis windows are centered on all samples in the data for method=”mtmconvol”. For wavelet based methods (“wavelet” or “superlet”) toi needs to be either an equidistant array of time points or “all”.

  • t_ftimwin (positive float) – Only valid if method is ‘mtmconvol’ or ‘welch’. Sliding window length (in seconds).

  • wavelet (str) – Only valid if method is ‘wavelet’. Wavelet function to use, one of availableWavelets (see below).

  • width (positive float) – Only valid if method is ‘wavelet’ and wavelet is ‘Morlet’. Nondimensional frequency constant of Morlet wavelet function. This number should be >= 6, which corresponds to 6 cycles within the analysis window to ensure sufficient spectral sampling.

  • order (positive int) – Only valid if method is ‘wavelet’ and wavelet is ‘Paul’ or ‘DOG’. Order of the wavelet function. If wavelet is ‘Paul’, order should be chosen >= 4 to ensure that the analysis window contains at least a single oscillation. At an order of 40, the Paul wavelet exhibits about the same number of cycles as the Morlet wavelet with a width of 6. All other supported wavelets functions are real-valued derivatives of Gaussians (DOGs). Hence, if wavelet is ‘DOG’, order represents the derivative order. The special case of a second order DOG yields a function known as “Mexican Hat”, “Marr” or “Ricker” wavelet, which can be selected alternatively by setting wavelet to ‘Mexican_hat’, ‘Marr’ or ‘Ricker’. Note: A real-valued wavelet function encodes only information about peaks and discontinuities in the signal and does not provide any information about amplitude or phase.

  • order_max (int) – Only valid if method is ‘superlet’. Maximal order of the superlet set. Controls the maximum number of cycles within a SL together with the c_1 parameter: c_max = c_1 * order_max

  • order_min (int) – Only valid if method is ‘superlet’. Minimal order of the superlet set. Controls the minimal number of cycles within a SL together with the c_1 parameter: c_min = c_1 * order_min Note that for admissability reasons c_min should be at least 3!

  • c_1 (int) – Only valid if method is ‘superlet’. Number of cycles of the base Morlet wavelet. If set to lower than 3 increase order_min as to never have less than 3 cycles in a wavelet!

  • adaptive (bool) – Only valid if method is ‘superlet’. Wether to perform multiplicative SLT or fractional adaptive SLT. If set to True, the order of the wavelet set will increase linearly with the frequencies of interest from order_min to order_max. If set to False the same SL will be used for all frequencies.

  • fooof_opt (dict or None) – Only valid if method is ‘mtmfft’ and output is one of ‘fooof’, ‘fooof_aperiodic’, or ‘fooof_peaks’. Additional keyword arguments passed to the FOOOF constructor. Available arguments include 'peak_width_limits', 'max_n_peaks', 'min_peak_height', 'peak_threshold', and 'aperiodic_mode'. Please refer to the FOOOF docs for the meanings and the defaults. See the FOOOF reference [Donoghue2020] for details.

  • ft_compat (bool, optional) – Set to True to use Field Trip’s spectral normalization for FFT based methods (method='mtmfft' and method='mtmconvol'). So spectral power is NOT independent of the padding size!

  • parallel (None or bool) – If None (recommended), processing is automatically performed in parallel (i.e., concurrently across trials/channel-groups), provided a dask parallel processing client is running and available. Parallel processing can be manually disabled by setting parallel to False. If parallel is True but no parallel processing client is running, computing will be performed sequentially.

  • select (dict or StructDict or str) – In-place selection of subset of input data for processing. Please refer to syncopy.selectdata() for further usage details.

Returns:

spec – (Time-)frequency spectrum of input data. The spec may contain additional metadata, based on the method used to compute it:

  • For method=’mtmfft’ when output is one of ‘fooof’, ‘fooof_aperiodic’, or ‘fooof_peaks’, the spec.metadata property contains the keys listed and explained in metadata_keys.

Return type:

SpectralData

Notes

This function can be either called providing its input arguments directly or via a cfg configuration ‘structure’. For instance, the following function calls are equivalent

>>> spy.freqanalysis(data, method=...)
>>> cfg = spy.StructDict()
>>> cfg.method = ...
>>> spy.freqanalysis(cfg, data)
>>> cfg.data = data
>>> spy.freqanalysis(cfg)

Please refer to Syncopy for FieldTrip Users for further details.

[Moca2021]

Moca, Vasile V., et al. “Time-frequency super-resolution with superlets.” Nature communications 12.1 (2021): 1-18.

[Donoghue2020]

Donoghue et al. 2020, DOI 10.1038/s41593-020-00744-x.

[Welch1967]

Welch. “The use of fast Fourier transform for the estimation of power spectra: A method based on time averaging over short, modified periodograms.”, 1976, DOI 10.1109/TAU.1967.1161901

Options

syncopy.specest.freqanalysis.availableMethods = ('mtmfft', 'mtmconvol', 'wavelet', 'superlet', 'welch')#

Built-in immutable sequence.

If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.

If the argument is a tuple, the return value is the same object.

syncopy.specest.freqanalysis.availableOutputs = ('pow', 'abs', 'fourier', 'real', 'imag', 'angle', 'absreal', 'absimag', 'complex')#

Built-in immutable sequence.

If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.

If the argument is a tuple, the return value is the same object.

syncopy.shared.const_def.availableTapers = ['boxcar', 'triang', 'parzen', 'bohman', 'blackman', 'nuttall', 'blackmanharris', 'flattop', 'bartlett', 'barthann', 'hamming', 'kaiser', 'kaiser_bessel_derived', 'gaussian', 'general_gaussian', 'general_cosine', 'general_hamming', 'chebwin', 'cosine', 'hann', 'tukey', 'taylor', 'lanczos']#

Built-in mutable sequence.

If no argument is given, the constructor creates a new empty list. The argument must be an iterable if specified.

syncopy.specest.freqanalysis.availableWavelets = ('Morlet', 'Paul', 'DOG', 'Ricker', 'Marr', 'Mexican_hat')#

Built-in immutable sequence.

If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.

If the argument is a tuple, the return value is the same object.

Examples

Generate 10 seconds of white noise, sampled at 1000 Hz:

>>> import syncopy.tests.synth_data as synth_data
>>> wn = synth_data.white_noise(nTrials=2, nChannels=3, nSamples=10000, samplerate=1000)

Configure Welch’s method to estimate the power spectral density with a window length of 0.25 seconds and an overlap of 50 percent between the windows:

>>> cfg = spy.get_defaults(spy.freqanalysis)
>>> cfg.method = "welch"
>>> cfg.t_ftimwin = 0.25  # Window length in seconds.
>>> cfg.toi = 0.5         # Overlap between periodograms (0.5 = 50 percent overlap).

Run Welch:

>>> psd = spy.freqanalysis(cfg, wn)

Visualize the result for the first trial:

>>> _, ax = res.singlepanelplot(trials=0, logscale=False)
>>> ax.set_title("Welch result")

See also

syncopy.specest.mtmfft.mtmfft

(multi-)tapered Fourier transform of multi-channel time series data

syncopy.specest.mtmconvol.mtmconvol

time-frequency analysis of multi-channel time series data with a sliding window FFT

syncopy.specest.wavelet.wavelet

time-frequency analysis of multi-channel time series data using a wavelet transform

syncopy.specest.fooofspy.fooofspy

parameterization of neural power spectra with the ‘fitting oscillations & one over f’ method

numpy.fft.fft

NumPy’s reference FFT implementation

scipy.signal.stft

SciPy’s Short Time Fourier Transform