neurodsp.spectral.compute_spectrum_welch

neurodsp.spectral.compute_spectrum_welch(sig, fs, avg_type='mean', window='hann', nperseg=None, noverlap=None, f_range=None, outlier_percent=None)[source]

Compute the power spectral density using Welch’s method.

Parameters
sig1d or 2d array

Time series.

fsfloat

Sampling rate, in Hz.

avg_type{‘mean’, ‘median’}, optional

Method to average across the windows:

  • ‘mean’ is the same as Welch’s method, taking the mean across FFT windows.

  • ‘median’ uses median across FFT windows instead of the mean, to minimize outlier effects.

windowstr or tuple or array_like, optional, default: ‘hann’

Desired window to use. See scipy.signal.get_window for a list of available windows. If array_like, the array will be used as the window and its length must be nperseg.

npersegint, optional

Length of each segment, in number of samples. If None, and window is str or tuple, is set to 1 second of data. If None, and window is array_like, is set to the length of the window.

noverlapint, optional

Number of points to overlap between segments. If None, noverlap = nperseg // 8.

f_rangelist of [float, float], optional

Frequency range to sub-select from the power spectrum.

outlier_percentfloat, optional

The percentage of outlier values to be removed. Must be between 0 and 100.

Returns
freqs1d array

Frequencies at which the measure was calculated.

spectrum1d or 2d array

Power spectral density.

Notes

  • Welch’s method ([1]) computes a power spectra by averaging over windowed FFTs.

References

1

Welch, P. (1967). The use of fast Fourier transform for the estimation of power spectra: A method based on time averaging over short, modified periodograms. IEEE Transactions on Audio and Electroacoustics, 15(2), 70–73. DOI: https://doi.org/10.1109/TAU.1967.1161901

Examples

Compute the power spectrum of a simulated time series using Welch’s method:

>>> from neurodsp.sim import sim_combined
>>> sig = sim_combined(n_seconds=10, fs=500,
...                    components={'sim_powerlaw': {}, 'sim_oscillation': {'freq': 10}})
>>> freqs, spec = compute_spectrum_welch(sig, fs=500)