Confidence intervals for multitaper power spectral density estimates
Confidence intervals for multitaper power spectral density estimates
spec_confint(x, ...) ## S3 method for class 'spec' spec_confint(x, ...) ## S3 method for class 'tapers' spec_confint(x, ...) ## Default S3 method: spec_confint(x, ...) .spec_confint(dof, p = 0.95, as.db = FALSE, ...)
x |
object to calculate spectral properties |
... |
additional arguments |
dof |
numeric; the degrees of freedom ν |
p |
numeric; the coverage probability p, bound within [0,1) |
as.db |
logical; should the values be returned as decibels? |
The errors are estimated from the number of degrees of freedom ν by evaluating the χ_{p,ν}^{2}(ν,ν) distribution for an optional coverage probability p (defaulting to p=0.95). Additionally, the p=0.5 values and an approximation from 1/√{ν - 1} are returned.
A more sophisticated (and complicated) approach would be to estimate via jack-knifing (Prieto et al 2007), but this is not yet made available.
Additive uncertainties δ S are returned, such that the spectrum with confidence interval is S \pm δ S.
A data.frame
with the following properties (and names):
lower
: Based on upper tail probabilities (p)
upper
: Based on lower tail probabilities (1-p)
median
: Based on lower tail probabilities (p=0.5)
approx
: Approximation based on 1/√(ν - 1).
A.J. Barbour; some code modified from the spec.ci
function inside stats::plot.spec
spectral_properties
, psd-package
, stats::plot.spec
, dB
## Not run: #REX library(psd) ## ## Confidence intervals from taper numbers ## sp <- spectral_properties(as.tapers(1:50), p=0.95, db.ci=TRUE) # standard errors as a function of tapers par(las=1) plot(stderr.chi.upper ~ taper, sp, type="s", ylim=c(-10,20), yaxs="i", xaxs="i", xlab=expression("number of tapers ("* nu/2 *")"), ylab="dB", main="Spectral uncertainties") mtext("(additive factor)", line=.3) lines(stderr.chi.lower ~ taper, sp, type="s") lines(stderr.chi.median ~ taper, sp, type="s", lwd=2) lines(stderr.chi.approx ~ taper, sp, type="s", col="red",lwd=2) # indicate K needed to reach 3 dB wide confidence interval (p=.95) abline(v=33, lty=3) legend("topright", c(expression("Based on "* chi^2 *"(p,"*nu*") and (1-p,"*nu*")"), expression(""* chi^2 *"(p=0.5,"*nu*")"), "approximation"), lwd=c(1,3,3), col=c("black","black","red"), bg="grey98") ## End(Not run)#REX
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.