Running 'bugfixes.R' [0s]
  Comparing 'bugfixes.Rout' to 'bugfixes.Rout.save' ... OK
  Running 'regtest_ctabs.R' [165s]
  Running 'regtest_libcoin.R' [164s]
Running the tests in 'tests/regtest_ctabs.R' failed.
Complete output:
  > 
  > library("libcoin")
  > 
  > x <- 1:6
  > y <- rep(1:3, length.out = 6)
  > w <- rep(2L, length(x))
  > b <- gl(2, 3)
  > s <- 1:4
  > 
  > ctabs(x)
  [1] 0 1 1 1 1 1 1
  > ctabs(x, weights = w)
  [1] 0 2 2 2 2 2 2
  > ctabs(x, subset = s)
  [1] 0 1 1 1 1 0 0
  > ctabs(x, weights = w, subset = s)
  [1] 0 2 2 2 2 0 0
  > ctabs(x, block = b)
       [,1] [,2]
  [1,]    0    0
  [2,]    1    0
  [3,]    1    0
  [4,]    1    0
  [5,]    0    1
  [6,]    0    1
  [7,]    0    1
  > ctabs(x, weights = w, block = b)
       [,1] [,2]
  [1,]    0    0
  [2,]    2    0
  [3,]    2    0
  [4,]    2    0
  [5,]    0    2
  [6,]    0    2
  [7,]    0    2
  > ctabs(x, subset = s, block = b)
       [,1] [,2]
  [1,]    0    0
  [2,]    1    0
  [3,]    1    0
  [4,]    1    0
  [5,]    0    1
  [6,]    0    0
  [7,]    0    0
  > ctabs(x, weights = w, subset = s, block = b)
       [,1] [,2]
  [1,]    0    0
  [2,]    2    0
  [3,]    2    0
  [4,]    2    0
  [5,]    0    2
  [6,]    0    0
  [7,]    0    0
  > 
  > ctabs(x, y)
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    1    0    0
  [3,]    0    0    1    0
  [4,]    0    0    0    1
  [5,]    0    1    0    0
  [6,]    0    0    1    0
  [7,]    0    0    0    1
  > ctabs(x, y, weights = w)
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    2    0    0
  [3,]    0    0    2    0
  [4,]    0    0    0    2
  [5,]    0    2    0    0
  [6,]    0    0    2    0
  [7,]    0    0    0    2
  > ctabs(x, y, subset = s)
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    1    0    0
  [3,]    0    0    1    0
  [4,]    0    0    0    1
  [5,]    0    1    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  > ctabs(x, y, weights = w, subset = s)
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    2    0    0
  [3,]    0    0    2    0
  [4,]    0    0    0    2
  [5,]    0    2    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  > ctabs(x, y, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    1    0    0
  [3,]    0    0    1    0
  [4,]    0    0    0    1
  [5,]    0    0    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    0    0    0
  [3,]    0    0    0    0
  [4,]    0    0    0    0
  [5,]    0    1    0    0
  [6,]    0    0    1    0
  [7,]    0    0    0    1
  
  > ctabs(x, y, weights = w, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    2    0    0
  [3,]    0    0    2    0
  [4,]    0    0    0    2
  [5,]    0    0    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    0    0    0
  [3,]    0    0    0    0
  [4,]    0    0    0    0
  [5,]    0    2    0    0
  [6,]    0    0    2    0
  [7,]    0    0    0    2
  
  > ctabs(x, y, subset = s, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    1    0    0
  [3,]    0    0    1    0
  [4,]    0    0    0    1
  [5,]    0    0    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    0    0    0
  [3,]    0    0    0    0
  [4,]    0    0    0    0
  [5,]    0    1    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  > ctabs(x, y, weights = w, subset = s, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    2    0    0
  [3,]    0    0    2    0
  [4,]    0    0    0    2
  [5,]    0    0    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0    0    0    0
  [2,]    0    0    0    0
  [3,]    0    0    0    0
  [4,]    0    0    0    0
  [5,]    0    2    0    0
  [6,]    0    0    0    0
  [7,]    0    0    0    0
  
  > 
  > w <- round(1:length(x) / length(x), 2)
  > ctabs(x, weights = w)
  [1] 0.00 0.17 0.33 0.50 0.67 0.83 1.00
  > ctabs(x, weights = w, subset = s)
  [1] 0.00 0.17 0.33 0.50 0.67 0.00 0.00
  > ctabs(x, weights = w, block = b)
       [,1] [,2]
  [1,] 0.00 0.00
  [2,] 0.17 0.00
  [3,] 0.33 0.00
  [4,] 0.50 0.00
  [5,] 0.00 0.67
  [6,] 0.00 0.83
  [7,] 0.00 1.00
  > ctabs(x, weights = w, subset = s, block = b)
       [,1] [,2]
  [1,] 0.00 0.00
  [2,] 0.17 0.00
  [3,] 0.33 0.00
  [4,] 0.50 0.00
  [5,] 0.00 0.67
  [6,] 0.00 0.00
  [7,] 0.00 0.00
  > 
  > ctabs(x, y, weights = w)
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00 0.00  0.0
  [2,]    0 0.17 0.00  0.0
  [3,]    0 0.00 0.33  0.0
  [4,]    0 0.00 0.00  0.5
  [5,]    0 0.67 0.00  0.0
  [6,]    0 0.00 0.83  0.0
  [7,]    0 0.00 0.00  1.0
  > ctabs(x, y, weights = w, subset = s)
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00 0.00  0.0
  [2,]    0 0.17 0.00  0.0
  [3,]    0 0.00 0.33  0.0
  [4,]    0 0.00 0.00  0.5
  [5,]    0 0.67 0.00  0.0
  [6,]    0 0.00 0.00  0.0
  [7,]    0 0.00 0.00  0.0
  > ctabs(x, y, weights = w, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00 0.00  0.0
  [2,]    0 0.17 0.00  0.0
  [3,]    0 0.00 0.33  0.0
  [4,]    0 0.00 0.00  0.5
  [5,]    0 0.00 0.00  0.0
  [6,]    0 0.00 0.00  0.0
  [7,]    0 0.00 0.00  0.0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00 0.00    0
  [2,]    0 0.00 0.00    0
  [3,]    0 0.00 0.00    0
  [4,]    0 0.00 0.00    0
  [5,]    0 0.67 0.00    0
  [6,]    0 0.00 0.83    0
  [7,]    0 0.00 0.00    1
  
  > ctabs(x, y, weights = w, subset = s, block = b)
  , , 1
  
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00 0.00  0.0
  [2,]    0 0.17 0.00  0.0
  [3,]    0 0.00 0.33  0.0
  [4,]    0 0.00 0.00  0.5
  [5,]    0 0.00 0.00  0.0
  [6,]    0 0.00 0.00  0.0
  [7,]    0 0.00 0.00  0.0
  
  , , 2
  
       [,1] [,2] [,3] [,4]
  [1,]    0 0.00    0    0
  [2,]    0 0.00    0    0
  [3,]    0 0.00    0    0
  [4,]    0 0.00    0    0
  [5,]    0 0.67    0    0
  [6,]    0 0.00    0    0
  [7,]    0 0.00    0    0
  
  > 
  > proc.time()
     user  system elapsed 
     0.14    0.09    0.21 
Running the tests in 'tests/regtest_libcoin.R' failed.
Complete output:
  > 
  > library("libcoin")
  > library("coin")
  Loading required package: survival
  > set.seed(29)
  > 
  > n <- 100
  > p <- 4
  > q <- 2
  > X <- matrix(runif(p * n), nc = p)
  > Y <- matrix(runif(q * n), nc = q)
  > w <- as.integer(floor(runif(n, max = 4)))
  > s <- sort(sample(1:n, floor(n/2), replace = TRUE))
  > b <- sample(gl(2, 2, length = n))
  > 
  > tol <- if (.Platform$OS.type == "unix") {
  +     sqrt(.Machine$double.eps)
  + } else {
  +    .Machine$double.eps^(1/5)
  + }
  > 
  > cmp <- function(t1, t2) {
  +     if (is.null(t1$Covariance)) {
  +         var1 <- t1$Variance
  +         var2 <- diag(covariance(t2))
  +     } else {
  +         var1 <- t1$Covariance
  +         var2 <- covariance(t2)
  +         var2 <- var2[!upper.tri(var2)]
  +     }
  +     stopifnot(all.equal(
  +         list(t1$LinearStatistic,
  +              t1$Expectation,
  +              var1),
  +         list(as.vector(statistic(t2, type = "linear")),
  +              as.vector(expectation(t2)),
  +              var2),
  +         check.attributes = FALSE, tolerance = tol
  +     ))
  + }
  > 
  > cmp2 <- function(t1, t2) {
  +     nm <- c("LinearStatistic", "Expectation",
  +             if(t1$varonly == 1) "Variance" else "Covariance")
  +     stopifnot(all.equal(t1[nm], t2[nm], tolerance = tol))
  + }
  > 
  > cmp3 <- function(t1, t2, pvalue = FALSE) {
  +     stopifnot(all.equal(statistic(t1), t2$TestStatistic, tolerance = tol))
  +     if (pvalue)
  +         stopifnot(all.equal(unclass(pvalue(t1)), t2$p.value, check.attributes = FALSE, tolerance = tol))
  + }
  > 
  > cmp4 <- function(t1, t2)
  +     stopifnot(all.equal(t1$TestStatistic, sqrt(t2$TestStatistic), tolerance = tol))
  > 
  > 
  > t1 <- LinStatExpCov(X, Y)
  > t1v <- LinStatExpCov(X, Y, varonly = TRUE)
  > t2 <- independence_test(Y ~ X)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > cmp3(t2, doTest(t1, teststat = "maximum"), pvalue = TRUE)
  > t3 <- independence_test(Y ~ X, teststat = "maximum", alternative = "less")
  > cmp3(t3, doTest(t1, teststat = "maximum", alternative = "less"), pvalue = TRUE)
  > t3 <- independence_test(Y ~ X, teststat = "maximum", alternative = "greater")
  > cmp3(t3, doTest(t1, teststat = "maximum", alternative = "greater"), pvalue = TRUE)
  > t3 <- independence_test(Y ~ X, teststat = "quadratic")
  > cmp3(t3, doTest(t1, teststat = "quadratic"), pvalue = TRUE)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w)
  > t1v <- LinStatExpCov(X, Y, weights = w, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, subset = s)
  > t1v <- LinStatExpCov(X, Y, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, subset = s)
  > t1v <- LinStatExpCov(X, Y, weights = w, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, block = b)
  > t1v <- LinStatExpCov(X, Y, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X  | b)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, block = b)
  > t1v <- LinStatExpCov(X, Y, weights = w, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, weights = w, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > n <- 100
  > n1 <- 5
  > n2 <- 4
  > p <- 4
  > q <- 2
  > X <- rbind(0, matrix(runif(p * n1), nc = p))
  > Y <- rbind(0, matrix(runif(q * n2), nc = q))
  > ix <- sample(1:n1, n, replace = TRUE)
  > iy <- sample(1:n2, n, replace = TRUE)
  > w <- as.integer(floor(runif(n, max = 4)))
  > s <- sort(sample(1:n, floor(n/2), replace = TRUE))
  > b <- sample(gl(2, 2, length = n))
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy)
  > t1v <- LinStatExpCov(X, Y, ix, iy, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,])
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,], weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, subset = s)
  > t1v <- LinStatExpCov(X, Y, ix, iy, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,], subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,], weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, block = b, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,]  | b)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, block = b, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,] | b, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,] | b, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y[iy + 1,] ~ X[ix + 1,]| b, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > 
  > #### X factor
  > n <- 10000
  > p <- 40
  > q <- 20
  > X <- diag(p)[fx <- unclass(factor(sample(1:p, n, replace = TRUE))),]
  > Y <- matrix(runif(q * n), nc = q)
  > w <- as.integer(floor(runif(n, max = 4)))
  > s <- sort(sample(1:n, floor(n/2), replace = TRUE))
  > b <- sample(gl(2, 2, length = n))
  > 
  > t1 <- LinStatExpCov(X, Y)
  > t1v <- LinStatExpCov(X, Y, varonly = TRUE)
  > t2 <- independence_test(Y ~ X)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y)
  > t1vf <- LinStatExpCov(fx, Y, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w)
  > t1v <- LinStatExpCov(X, Y, weights = w, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, weights = w)
  > t1vf <- LinStatExpCov(fx, Y, weights = w, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, subset = s)
  > t1v <- LinStatExpCov(X, Y, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, subset = s)
  > t1vf <- LinStatExpCov(fx, Y, subset = s, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, subset = s)
  > t1v <- LinStatExpCov(X, Y, weights = w, subset = s, varonly = TRUE)
  > t2 <- independence_test(Y ~ X, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, weights = w, subset = s)
  > t1vf <- LinStatExpCov(fx, Y, weights = w, subset = s, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, block = b)
  > t1v <- LinStatExpCov(X, Y, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X  | b)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, block = b)
  > t1vf <- LinStatExpCov(fx, Y, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, block = b)
  > t1v <- LinStatExpCov(X, Y, weights = w, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, weights = w, block = b)
  > t1vf <- LinStatExpCov(fx, Y, weights = w, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, subset = s, block = b)
  > t1vf <- LinStatExpCov(fx, Y, subset = s, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, weights = w, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, weights = w, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(Y ~ X | b, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(fx, Y, weights = w, subset = s, block = b)
  > t1vf <- LinStatExpCov(fx, Y, weights = w, subset = s, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > n1 <- 5
  > n2 <- 7
  > X <- rbind(0, diag(n1))
  > Y <- rbind(0, matrix(runif(q * n2), nc = q))
  > ix <- sample(1:n1, n, replace = TRUE)
  > iy <- sample(1:n2, n, replace = TRUE)
  > w <- as.integer(floor(runif(n, max = 4)))
  > s <- sort(sample(1:n, floor(n/2), replace = TRUE))
  > b <- sample(gl(2, 2, length = n))
  > 
  > YY <- Y[iy + 1,]
  > XX <- X[ix + 1,]
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy)
  > t1v <- LinStatExpCov(X, Y, ix, iy, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, subset = s)
  > t1v <- LinStatExpCov(X, Y, ix, iy, subset = s, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, subset = s)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, subset = s, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, subset = s)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, subset = s, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, block = b, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX  | b)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, block = b)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, block = b, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX | b, weights = ~ w)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, block = b)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX | b, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, subset = s, block = b)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, subset = s, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > t1 <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, block = b)
  > t1v <- LinStatExpCov(X, Y, ix, iy, weights = w, subset = s, block = b, varonly = TRUE)
  > t2 <- independence_test(YY ~ XX| b, weights = ~w, subset = s)
  > cmp(t1, t2)
  > cmp(t1v, t2)
  > t1f <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, subset = s, block = b)
  > t1vf <- LinStatExpCov(numeric(0), Y, ix, iy, weights = w, subset = s, block = b, varonly = TRUE)
  > cmp2(t1, t1f)
  > cmp2(t1v, t1vf)
  > 
  > ### and now maximally selected statistics
  > n <- 100
  > B <- 1000
  > x <- round(runif(n), 2)
  > y <- rnorm(n, mean = x < .5, sd = 2.6)
  > y2 <- runif(n)
  > blk <- gl(4, n/4)
  > ux <- sort(unique(x))
  > ix <- unclass(cut(x, breaks = c(-Inf, ux[-length(ux)] + diff(ux) / 2, Inf)))
  > 
  > cmp3 <- function(t1, t2)
  +     stopifnot(all.equal(statistic(t1), t2$TestStatistic, tolerance = tol))
  > 
  > cmp4 <- function(t1, t2)
  +     stopifnot(all.equal(t1$TestStatistic, sqrt(t2$TestStatistic), tolerance = tol))
  > 
  > (mt <- maxstat_test(y ~ x , distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y by x
  maxT = 2.1537, p-value = 0.26
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: <= 0.53
  
  > lev <- LinStatExpCov(ix, y, nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.153664
  
  $p.value
  [1] 0.245
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst2 <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 4.63827
  
  $p.value
  [1] 0.245
  
  $index
  [1] 36
  
  > cmp4(tst, tst2)
  > ux[tst2$index]
  [1] 0.53
  > lev <- LinStatExpCov(ix, y, nresample = 1000, varonly = TRUE)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.153664
  
  $p.value
  [1] 0.253
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst2 <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 4.63827
  
  $p.value
  [1] 0.253
  
  $index
  [1] 36
  
  > cmp4(tst, tst2)
  > ux[tst2$index]
  [1] 0.53
  > 
  > (mt <- maxstat_test(y ~ x | blk, distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y by x 
  	 stratified by blk
  maxT = 2.2764, p-value = 0.221
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: <= 0.53
  
  > lev <- LinStatExpCov(ix, y, block = blk, nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.276446
  
  $p.value
  [1] 0.192
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst2 <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 5.182204
  
  $p.value
  [1] 0.192
  
  $index
  [1] 36
  
  > cmp4(tst, tst2)
  > ux[tst$index]
  [1] 0.53
  > lev <- LinStatExpCov(ix, y, block = blk, nresample = 1000, varonly = TRUE)
  > try(tst <- doTest(lev, teststat = "maximum"))
  Error in doTest(lev, teststat = "maximum") : 
    need covariance for maximally statistics with blocks
  > 
  > (mt <- maxstat_test(y + y2 ~ x , distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y, y2 by x
  maxT = 2.1537, p-value = 0.445
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: <= 0.53
  
  > lev <- LinStatExpCov(ix, cbind(y, y2), nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.153664
  
  $p.value
  [1] 0.471
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 5.720153
  
  $p.value
  [1] 0.412
  
  $index
  [1] 40
  
  > ux[tst$index]
  [1] 0.57
  > lev <- LinStatExpCov(ix, cbind(y, y2), nresample = 1000, varonly = TRUE)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.153664
  
  $p.value
  [1] 0.476
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 5.720153
  
  $p.value
  [1] 0.429
  
  $index
  [1] 40
  
  > ux[tst$index]
  [1] 0.57
  > 
  > (mt <- maxstat_test(y + y2 ~ x | blk, distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y, y2 by x 
  	 stratified by blk
  maxT = 2.2764, p-value = 0.38
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: <= 0.53
  
  > lev <- LinStatExpCov(ix, cbind(y, y2), block = blk, nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum"))
  $TestStatistic
  [1] 2.276446
  
  $p.value
  [1] 0.394
  
  $index
  [1] 36
  
  > cmp3(mt, tst)
  > ux[tst$index]
  [1] 0.53
  > (tst <- doTest(lev, teststat = "quadratic"))
  $TestStatistic
  [1] 5.406605
  
  $p.value
  [1] 0.504
  
  $index
  [1] 40
  
  > ux[tst$index]
  [1] 0.57
  > lev <- LinStatExpCov(ix, cbind(y, y2), block = blk, nresample = 1000, varonly = TRUE)
  > try(tst <- doTest(lev, teststat = "maximum"))
  Error in doTest(lev, teststat = "maximum") : 
    need covariance for maximally statistics with blocks
  > 
  > x <- sample(gl(5, n))
  > y <- rnorm(length(x), mean = x %in% levels(x)[c(1, 3, 5)], sd = 4.5)
  > y2 <- runif(length(x))
  > ix <- unclass(x)
  > blk <- gl(5, n)
  > 
  > (mt <- maxstat_test(y ~ x , distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y by x (1, 2, 3, 4, 5)
  maxT = 1.7088, p-value = 0.443
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: {1, 3} vs. {2, 4, 5}
  
  > lev <- LinStatExpCov(ix, y, nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 1.708775
  
  $p.value
  [1] 0.485
  
  $index
  [1] 0 1 0 1 1
  
  > cmp3(mt, tst)
  > (tst2 <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 2.919912
  
  $p.value
  [1] 0.485
  
  $index
  [1] 0 1 0 1 1
  
  > cmp4(tst, tst2)
  > lev <- LinStatExpCov(ix, y, nresample = 1000, varonly = TRUE)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 1.708775
  
  $p.value
  [1] 0.476
  
  $index
  [1] 0 1 0 1 1
  
  > cmp3(mt, tst)
  > (tst2 <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 2.919912
  
  $p.value
  [1] 0.476
  
  $index
  [1] 0 1 0 1 1
  
  > cmp4(tst, tst2)
  > 
  > (mt <- maxstat_test(y ~ x | blk, distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y by x (1, 2, 3, 4, 5) 
  	 stratified by blk
  maxT = 1.6489, p-value = 0.543
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: {1, 3} vs. {2, 4, 5}
  
  > lev <- LinStatExpCov(ix, y, block = blk, nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 1.648947
  
  $p.value
  [1] 0.464
  
  $index
  [1] 0 1 0 1 1
  
  > cmp3(mt, tst)
  > (tst2 <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 2.719028
  
  $p.value
  [1] 0.464
  
  $index
  [1] 0 1 0 1 1
  
  > cmp4(tst, tst2)
  > lev <- LinStatExpCov(ix, y, block = blk, nresample = 1000, varonly = TRUE)
  > try(tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  Error in doTest(lev, teststat = "maximum", ordered = FALSE) : 
    need covariance for maximally statistics with blocks
  > 
  > (mt <- maxstat_test(y + y2 ~ x , distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y, y2 by x (1, 2, 3, 4, 5)
  maxT = 2.358, p-value = 0.285
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: {1, 3, 5} vs. {2, 4}
  
  > lev <- LinStatExpCov(ix, cbind(y, y2), nresample = 1000)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 2.357973
  
  $p.value
  [1] 0.301
  
  $index
  [1] 0 1 0 1 0
  
  > cmp3(mt, tst)
  > (tst <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 8.373114
  
  $p.value
  [1] 0.143
  
  $index
  [1] 0 1 0 1 1
  
  > lev <- LinStatExpCov(ix, cbind(y, y2), nresample = 1000, varonly = TRUE)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 2.357973
  
  $p.value
  [1] 0.275
  
  $index
  [1] 0 1 0 1 0
  
  > cmp3(mt, tst)
  > (tst <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 8.373114
  
  $p.value
  [1] 0.134
  
  $index
  [1] 0 1 0 1 1
  
  > 
  > (mt <- maxstat_test(y + y2 ~ x | blk, distrib = approximate(nresample = 1000)))
  
  	Approximative Generalized Maximally Selected Statistics
  
  data:  y, y2 by x (1, 2, 3, 4, 5) 
  	 stratified by blk
  maxT = 2.4473, p-value = 0.218
  alternative hypothesis: two.sided
  sample estimates:
    "best" cutpoint: {1, 3, 5} vs. {2, 4}
  
  > lev <- LinStatExpCov(ix, cbind(y, y2), block = blk, nresample = 50)
  > (tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  $TestStatistic
  [1] 2.447272
  
  $p.value
  [1] 0.2
  
  $index
  [1] 0 1 0 1 0
  
  > cmp3(mt, tst)
  > (tst2 <- doTest(lev, teststat = "quadratic", ordered = FALSE))
  $TestStatistic
  [1] 8.115438
  
  $p.value
  [1] 0.12
  
  $index
  [1] 0 1 0 1 0
  
  > 
  > xx <- factor(x %in% levels(x)[tst2$index == 1])
  > (it <- independence_test(y + y2 ~ xx | blk, teststat = "quadratic"))
  
  	Asymptotic General Independence Test
  
  data:  y, y2 by xx (FALSE, TRUE) 
  	 stratified by blk
  chi-squared = 8.1154, df = 2, p-value = 0.01729
  
  > cmp3(it, tst2)
  > 
  > lev <- LinStatExpCov(ix, cbind(y, y2), block = blk, nresample = 1000, varonly = TRUE)
  > try(tst <- doTest(lev, teststat = "maximum", ordered = FALSE))
  Error in doTest(lev, teststat = "maximum", ordered = FALSE) : 
    need covariance for maximally statistics with blocks
  > 
  > proc.time()
     user  system elapsed 
     6.09    0.65    6.73