test.efdr.condsim: Test for anomalies in wavelet space via conditional...

In EFDR: Wavelet-Based Enhanced FDR for Detecting Signals from Complete or Incomplete Spatially Aggregated Data

Test for anomalies in wavelet space via conditional simulation

Description

Test for anomalies using EFDR and conditional simulation. The noisy image can be partially observed, or/and aggregated at different resolutions

Usage

test.efdr.condsim(
  Zvec,
  H,
  n1,
  n2,
  rho_est_method = c("CPL", "MOM"),
  iter.cs = 100,
  wf = "la8",
  J = 2,
  alpha = 0.05,
  n.hyp = 100,
  b = 11,
  iteration = 200,
  parallel = 1L
)

Arguments

Zvec	vector of observations such that Ztilde = H.Z
H	matrix mapping the fine-resolution image Z in vector form to Ztllde. Must have as many rows as Ztilde and n1 x n2 columns
n1	number of rows in fine-resolution image
n2	number of columns in fine-resolution image
rho_est_method	method with which to estimate the level of exchangeability rho; can be either "CPL" (copula model) or "MOM" (method of moments)
iter.cs	number of conditional simulations to carry out
wf	type of wavelet to employ. Defaults to ‘la8’, the Daubechies orthonormal compactly supported wavelet of length L = 8 (Daubechies, 1992), least asymmetric family. Other options include ‘haar’ (Haar wavelet), ‘fk8’ (Fejer-Korovkin wavelet with L=8) and ‘mb8’ (minimum-bandwidth wavelet with L=8). Please type 'waveslim::wave.filter' in the console for a full list of wavelet names
J	number of resolutions to employ in wavelet decomposition
alpha	significance level at which tests are carried out
n.hyp	number of hypotheses tests to carry out with EFDR. If a vector is supplied, the optimal one from the set of proposed number of tests is chosen
b	the number of neighbours to consider in EFDR
iteration	number of Monte Carlo iterations to employ when determining which of the proposed number of tests in n.hyp is the optimal number of tests
parallel	number of cores to use with parallel backend; needs to be an integer less than or equal to the number of available cores

Value

List with three fields:

filtered

the discrete wavelet transform containing the anomalous wavelet coefficients in the signal

Z

the image containing the anomalous wavelets in the signal

reject_coeff

indices of wavelets under which the null hypothesis of no anomaly was rejected

pvalue_ordered

ordered p-values under the null hypothesis. The column names indicate the wavelet to which the p-value belongs

nhat

the number of tests carried out.

References

Daubechies, I. (1992) Ten Lectures on Wavelets, CBMS-NSF Regional Conference Series in Applied Mathematics, SIAM: Philadelphia.

Shen, X., Huang, H.-C., and Cressie, N. 'Nonparametric hypothesis testing for a spatial signal.' Journal of the American Statistical Association 97.460 (2002): 1122-1140.

Examples

## Set up experiment
n <- 32       # 32 x 32 images
r <- 10       # signal of size 10 x 10
h <- 5        # intensity of 5
grid <- 8     # aggregated to 8 x 8 image
parallel <- 4 # use 4 cores

## Simulate the pixel-level data
raw_grid <- expand.grid(x = seq(1, n), y = seq(1, n))
df <- data.frame(raw_grid)                        # spatial grid
dd <- as.matrix(dist(raw_grid, diag = TRUE))      # distance matrix
Sigma <- exp(-dd/5)                               # cov. fn.
diag(Sigma) <- 1                                  # fix diagonal
L <- t(chol(Sigma))                               # lower Cholesky factor
mu <- matrix(0, n, n)                             # zero mean
mu[(n/2-r/2):(n/2+r/2), (n/2-r/2):(n/2+r/2)] <- h # add signal
Z <- mu + matrix(L %*% rnorm(n^2), n, n)          # simulate data

## Construct H (aggregation) matrix
H <- matrix(0, grid^2, n^2)
for(i in 1:grid^2) {
  ind <- rep(rep(c(0L,1L,0L),
             c((n/grid)*((i-1)%%grid),n/grid,(n-n/grid-n/grid*((i-1)%%grid)))),
          n/grid)
  H[i,which(c(rep(0L,(ceiling(i/grid)-1)*n^2/grid),ind) == TRUE)] <- 1/(n/grid)^2
}

## Aggregate the signal
z_tilde <- c(H %*% c(Z))

## Run EFDR using conditional simulation
## Not run: out2 <- test.efdr.condsim(Zvec = z_tilde, H = H, n1 = n, n2 = n, 
                                   parallel = parallel)
## End(Not run)

EFDR documentation built on Aug. 23, 2023, 1:07 a.m.