tests/testthat/helper-workflows.R

In MDgof: Various Methods for the Goodness-of-Fit Problem in D>1 Dimensions

expect_gof_contract <- function(x) {
  expect_type(x, "list")
  expect_true(all(c("statistics", "p.values") %in% names(x)))
  expect_true(is.numeric(x$statistics))
  expect_true(is.numeric(x$p.values))
  expect_true(all(is.finite(x$statistics)))
  expect_true(all(is.finite(x$p.values)))
  expect_true(all(x$p.values >= 0 & x$p.values <= 1))
}

expect_probability_output <- function(x) {
  expect_true(is.numeric(x) || is.matrix(x) || is.data.frame(x))
  y <- as.numeric(x)
  expect_true(length(y) > 0L)
  expect_true(all(is.finite(y)))
  expect_true(all(y >= 0 & y <= 1))
}

# Continuous two-dimensional uniform model, matching tests_continuous.Rmd.
pnull_u2 <- function(x) {
  if (!is.matrix(x)) x <- rbind(x)
  apply(x, 1L, function(z) prod(stats::punif(z)))
}

dnull_u2 <- function(x) {
  if (!is.matrix(x)) return(1)
  rep(1, nrow(x))
}

rnull_u2 <- function() {
  matrix(stats::runif(2L * 250L), ncol = 2L)
}

ralt_u2 <- function(a) {
  matrix(stats::rbeta(2L * 250L, a, a), ncol = 2L)
}

# Continuous normal model with estimated mean.
pnull_ne <- function(x, mu) {
  if (!is.matrix(x)) x <- rbind(x)
  apply(x, 1L, function(z) {
    as.numeric(mvtnorm::pmvnorm(rep(-Inf, length(z)), z,
      mean = mu, sigma = diag(length(z))))
  })
}

dnull_ne <- function(x, mu) {
  if (!is.matrix(x)) x <- rbind(x)
  apply(x, 1L, mvtnorm::dmvnorm, mean = mu)
}

rnull_ne <- function(mu) {
  mvtnorm::rmvnorm(250L, mean = mu)
}

ralt_ne <- function(mu) {
  mvtnorm::rmvnorm(250L, mean = c(mu, mu))
}

phat_mean <- function(x) colMeans(x)

# Discrete binomial model, matching tests_disc.Rmd.
make_binomial_table <- function(p1 = 0.5, p2 = 0.6, n = 1000L) {
  nb <- c(5L, 5L)
  x <- stats::rbinom(n, nb[1L], p1)
  y <- stats::rbinom(n, nb[2L], p2)
  z <- matrix(0, prod(nb + 1L), 3L)
  k <- 0L
  for (i in seq_len(nb[1L] + 1L)) {
    for (j in seq_len(nb[2L] + 1L)) {
      k <- k + 1L
      z[k, 1:2] <- c(i, j) - 1L
      z[k, 3L] <- sum(x == i - 1L & y == j - 1L)
    }
  }
  z
}

pnull_binom <- function(x) {
  if (!is.matrix(x)) x <- rbind(x)
  stats::pbinom(x[, 1L], 5L, 0.5) * stats::pbinom(x[, 2L], 5L, 0.6)
}

rnull_binom <- function() make_binomial_table()
ralt_binom <- function(p) make_binomial_table(p2 = p)

pnull_binom_est <- function(x, p) {
  if (!is.matrix(x)) x <- rbind(x)
  stats::pbinom(x[, 1L], 5L, 0.5) * stats::pbinom(x[, 2L], 5L, p)
}

rnull_binom_est <- function(p) make_binomial_table(p2 = p)
ralt_binom_est <- function(p) make_binomial_table(p1 = p, p2 = 0.6)

phat_binom <- function(x) {
  nb <- c(5, 5)
  z <- stats::aggregate(x[, 3L], by = list(x[, 2L]), FUN = sum)
  sum(z[, 1L] * z[, 2L]) / nb[2L] / sum(z[, 2L])
}