prcompBoot: Calculate number of non-mixed EOFs (prcomp version)

In marchtaylor/sinkr: Collection of functions with emphasis in multivariate data analysis

Calculate number of non-mixed EOFs (prcomp version)

Description

The prcompBoot function uses a bootstrap randomization approach to calculate distributions of Empirical Orthogonal Function analysis (EOF) singular values with the prcomp function. EOF mode significance is assessed against the distributions of neighboring EOF singular values ("Lambda") calculated by the permutated models. A bootstrap routine follows the procedure of Babamoradi et al. (2013) whereby permutations sample rows (samples) more than once, which is a non-parametric approach does not make assumptions about the distribution of data.

Usage

prcompBoot(
  x,
  retx = TRUE,
  center = TRUE,
  scale. = FALSE,
  tol = NULL,
  nperm = 99
)

Arguments

x, retx, center, scale., tol	See prcomp for argument definitions.
nperm	Numeric. The number of null model permutations to calculate.

References

Babamoradi, H., van den Berg, F., Rinnan, A, 2013. Bootstrap based confidence limits in principal component analysis - A case study, Chemometrics and Intelligent Laboratory Systems, Volume 120, pp. 97-105. doi:10.1016/j.chemolab.2012.10.007.

Examples

# Generate data
m=50
n=100
frac.gaps <- 0.5 # the fraction of data with NaNs
N.S.ratio <- 0.1 # the Noise to Signal ratio for adding noise to data
x <- (seq(m)*2*pi)/m
t <- (seq(n)*2*pi)/n

# True field
Xt <- 
  outer(sin(x), sin(t)) + 
  outer(sin(2.1*x), sin(2.1*t)) + 
  outer(sin(3.1*x), sin(3.1*t)) +
  outer(tanh(x), cos(t)) + 
  outer(tanh(2*x), cos(2.1*t)) + 
  outer(tanh(4*x), cos(0.1*t)) + 
  outer(tanh(2.4*x), cos(1.1*t)) + 
  tanh(outer(x, t, FUN="+")) + 
  tanh(outer(x, 2*t, FUN="+"))

Xt <- t(Xt)

# Noise field
set.seed(1)
RAND <- matrix(runif(length(Xt), min=-1, max=1), nrow=nrow(Xt), ncol=ncol(Xt))
R <- RAND * N.S.ratio * Xt

# True field + Noise field
Xp <- Xt + R

res <- prcompBoot(Xp, center=FALSE, scale=FALSE, nperm=499)
ylim <- range(res$Lambda.orig, res$Lambda)
boxplot(res$Lambda, log="y", col=8, border=2, outpch="", ylim=ylim)
points(res$Lambda.orig)
abline(v=res$n.sig+0.5, lty=2, col=4)
mtext(paste("Non-mixed PCs =", res$n.sig), side=3, line=0.5, col=4)

marchtaylor/sinkr documentation built on July 4, 2022, 5:48 p.m.