Learn R Programming

EFDR (version 1.3)

test.efdr.condsim: 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
)

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.

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

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

Run this code
## 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
if (FALSE) out2 <- test.efdr.condsim(Zvec = z_tilde, H = H, n1 = n, n2 = n, 
                                   parallel = parallel)

Run the code above in your browser using DataLab