## Separation through the origin: the check itself, then what the fitters do
## with it.

## Classes split by the sign of the first covariate, which is centred on zero,
## so the covariates alone separate them with no help from the intercept.
sep_binary <- function(n = 60, seed = 1) {
  set.seed(seed)
  x1 <- c(-seq_len(n / 2), seq_len(n / 2)) / 10
  d <- data.frame(y = as.integer(x1 > 0), x1 = x1,
                  x2 = stats::rnorm(n), x3 = stats::rnorm(n))
  list(data = d, formula = y ~ x1 + x2 + x3)
}

sep_warnings <- function(expr) {
  n <- 0L
  value <- withCallingHandlers(expr, bbqr_separation_warning = function(w) {
    n <<- n + 1L
    invokeRestart("muffleWarning")
  })
  list(n = n, value = value)
}

test_that("the check decides separable and non-separable designs", {
  X <- cbind(c(-2, -1, 1, 2), c(0.3, -0.1, 0.2, -0.4))
  expect_true(bbqr:::.bbqr_separated_origin(X, c(0, 0, 1, 1)))
  # the same points with interleaved labels: no direction can split them
  expect_false(bbqr:::.bbqr_separated_origin(X[, 1, drop = FALSE], c(0, 1, 0, 1)))
  # a threshold away from zero needs the intercept, so not through the origin
  expect_false(bbqr:::.bbqr_separated_origin(X[, 1, drop = FALSE] + 10,
                                             c(0, 0, 1, 1)))
  # an observation at the origin is on every such hyperplane
  expect_false(bbqr:::.bbqr_separated_origin(rbind(X, 0), c(0, 0, 1, 1, 1)))
})

test_that("the check agrees with linear programming", {
  skip_if_not_installed("lpSolve")
  by_lp <- function(X, y) {
    A <- (2 * y - 1) * X
    lpSolve::lp("min", rep(0, 2 * ncol(A)), cbind(A, -A),
                rep(">=", nrow(A)), rep(1, nrow(A)))$status == 0
  }
  set.seed(101)
  ours <- theirs <- logical(150)
  for (i in seq_along(ours)) {
    # twenty points in ten dimensions: about half of such samples separate
    X <- matrix(stats::rnorm(20 * 10), 20, 10)
    y <- stats::rbinom(20, 1, 0.5)
    ours[i]   <- bbqr:::.bbqr_separated_origin(X, y)
    theirs[i] <- by_lp(X, y)
  }
  expect_false(anyNA(ours))
  expect_identical(ours, theirs)
  expect_gt(sum(theirs), 30)
  expect_gt(sum(!theirs), 30)
})

test_that("a separated sample is flagged and warned about once", {
  s <- sep_binary()
  for (pen in c("alasso", "lasso")) {
    expect_warning(
      fit <- bbqr(s$formula, s$data, penalty = pen, ndraw = 200, burn = 50),
      class = "bbqr_separation_warning")
    expect_true(fit$separated_origin)
    expect_output(print(fit), "separated through the origin")
    expect_output(print(summary(fit)), "separated through the origin")
    expect_true(summary(fit)$separated_origin)
  }

  # several quantile levels share one design: one warning, every fit flagged
  w <- sep_warnings(bbqr(s$formula, s$data, quantile = c(0.25, 0.5, 0.75),
                         ndraw = 200, burn = 50))
  expect_identical(w$n, 1L)
  expect_true(all(vapply(w$value, function(f) f$separated_origin, logical(1))))
  expect_output(print(w$value), "separated through the origin")

  # the exported fitters check too when called directly
  expect_warning(bbqr_alasso(s$formula, s$data, ndraw = 200, burn = 50),
                 class = "bbqr_separation_warning")
  expect_warning(bbqr_lasso(s$formula, s$data, ndraw = 200, burn = 50),
                 class = "bbqr_separation_warning")
})

test_that("nothing is said outside the configuration the result covers", {
  s <- sep_binary()
  # normal slope prior: the fact is recorded, nothing is claimed
  expect_warning(
    none <- bbqr(s$formula, s$data, penalty = "none", ndraw = 200, burn = 50),
    NA)
  expect_true(none$separated_origin)
  expect_false(any(grepl("separated", capture.output(print(none)))))
  # a sampled sigma is outside the stated case
  expect_warning(
    free <- bbqr(s$formula, s$data, anchor = "free", ndraw = 200, burn = 50),
    NA)
  expect_true(free$separated_origin)
  # so is the published hierarchy
  expect_warning(
    v3 <- bbqr(s$formula, s$data, prior = prior("alasso", model = "v3"),
               ndraw = 200, burn = 50),
    NA)
  expect_true(v3$separated_origin)

  # overlapping classes
  o <- sim_binary(n = 200, seed = 3)
  expect_warning(fit <- bbqr(o$formula, o$data, ndraw = 200, burn = 50), NA)
  expect_false(fit$separated_origin)
  expect_false(any(grepl("separated", capture.output(print(summary(fit))))))
})

test_that("the check runs on the design the sampler sees", {
  set.seed(5)
  x0 <- c(-(1:20), 1:20) / 10
  d <- data.frame(y = as.integer(x0 > 0), x = x0 + 10)
  # as supplied, x is positive throughout: the split needs the intercept
  expect_warning(raw <- bbqr(y ~ x, d, ndraw = 200, burn = 50), NA)
  expect_false(raw$separated_origin)
  # centred, the same split passes through the origin
  expect_warning(std <- bbqr(y ~ x, d, standardize = TRUE,
                             ndraw = 200, burn = 50),
                 class = "bbqr_separation_warning")
  expect_true(std$separated_origin)
})

test_that("the check can be switched off and never touches the draws", {
  s <- sep_binary()
  set.seed(9)
  on <- suppressWarnings(bbqr(s$formula, s$data, ndraw = 300, burn = 50))
  set.seed(9)
  expect_warning(
    off <- bbqr(s$formula, s$data, ndraw = 300, burn = 50,
                check_separation = FALSE),
    NA)
  expect_identical(off$separated_origin, NA)
  expect_identical(on$beta, off$beta)
  expect_identical(on$beta0, off$beta0)
  expect_false(any(grepl("separated", capture.output(print(off)))))

  old <- options(bbqr.check_separation = FALSE)
  on.exit(options(old), add = TRUE)
  expect_warning(opt <- bbqr(s$formula, s$data, ndraw = 200, burn = 50), NA)
  expect_identical(opt$separated_origin, NA)
  expect_warning(
    direct <- bbqr_lasso(s$formula, s$data, ndraw = 200, burn = 50),
    NA)
  expect_identical(direct$separated_origin, NA)

  expect_error(bbqr(s$formula, s$data, ndraw = 200, check_separation = NA),
               "check_separation")
})

test_that("a fit saved before the check existed prints as it always did", {
  s <- sim_binary()
  fit <- bbqr(s$formula, s$data, ndraw = 200, burn = 50)
  fit$separated_origin <- NULL
  expect_false(any(grepl("separated", capture.output(print(fit)))))
  expect_false(any(grepl("separated", capture.output(print(summary(fit))))))
})