welchPSD: Power spectral density estimation using Welch's method.

In bspec: Bayesian Spectral Inference

Power spectral density estimation using Welch's method.

Description

Estimates a time series' power spectral density using Welch's method, i.e., by subdividing the data into segments, computing spectra for each, and averaging over the results.

Usage

welchPSD(x, seglength, two.sided = FALSE, windowfun = tukeywindow,
         method = c("mean", "median"), windowingPsdCorrection = TRUE, ...)

Arguments

x	a time series (ts object).
seglength	the length of the subsegments to be used (in units of time relative to x).
two.sided	a logical flag indicating whether the result is supposed to be one-sided (default) or two-sided.
windowfun	The windowing function to be used.
method	The "averaging" method to be used – either "mean" or "median".
windowingPsdCorrection	a logical flag indicating whether an overall correction for the windowing is supposed to be applied to the result – this essentially specifies whether the result is supposed to reflect the PSD of the windowed or of the "un-windowed" time series.
...	other parameters passed to the windowing function.

Details

The time series will be divided into overlapping sub-segments, each segment is windowed and its "empirical" spectrum is computed. The average of these spectra is the resulting PSD estimate. For robustness, the median may also be used instead of the mean.

Value

A list containing the following elements:

frequency	the Fourier frequencies.
power	the estimated spectral power.
kappa	the number of (by definition) non-zero imaginary components of the Fourier series.
two.sided	a logical flag indicating one- or two-sidedness.
segments	the number of (overlapping) segments used.

Author(s)

Christian Roever, [e-mail]

References

Welch, P. D. 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, AU-15(2):70–73, 1967. doi: 10.1109/TAU.1967.1161901.

Press, W. H., Teukolsky, S. A., Vetterling, W. T., Flannery, B. P. Numerical recipes in C. Cambridge University Press, 1992.

See Also

empiricalSpectrum, tukeywindow, spectrum

Examples

# load example data:
data(sunspots)
# compute and plot the "plain" spectrum:
spec1 <- empiricalSpectrum(sunspots)
plot(spec1$frequency, spec1$power, log="y", type="l")

# plot Welch spectrum using segments of length 10 years:
spec2 <- welchPSD(sunspots, seglength=10)
lines(spec2$frequency, spec2$power, col="red")

# use 20-year segments and a flatter Tukey window:
spec3 <- welchPSD(sunspots, seglength=20, r=0.2)
lines(spec3$frequency, spec3$power, col="blue")

bspec documentation built on April 20, 2022, 5:08 p.m.