tests/testthat/test-tabforest.R

In R4VN: Health Data Analysis and Publication-Ready Reporting

test_that("tabforest fits crude ORs and respects focal rows", {
  set.seed(20260813)
  n <- 220
  d <- data.frame(
    age = rnorm(n, 50, 12),
    sex = factor(sample(c("Female", "Male"), n, TRUE), levels = c("Female", "Male")),
    smoking = factor(sample(c("No", "Yes"), n, TRUE), levels = c("No", "Yes"))
  )
  lp <- -3 + .04 * d$age + .6 * (d$smoking == "Yes") + .3 * (d$sex == "Male")
  d$y <- factor(ifelse(rbinom(n, 1, plogis(lp)) == 1, "Yes", "No"), levels = c("No", "Yes"))
  z <- tabforest(y, predictors = vars(c.age, sex, smoking), data = d,
                 event = "Yes", show = FALSE)
  expect_s3_class(z, "r4vn_tabforest")
  expect_equal(z$effect, "OR")
  expect_true(all(c("estimate", "lower", "upper", "p", "model") %in% names(z$table)))
  expect_true(all(z$table$estimate[!z$table$reference] > 0, na.rm = TRUE))
})

test_that("adjusted and multi have distinct R4VN semantics", {
  set.seed(20260813)
  n <- 260
  d <- data.frame(
    a = rnorm(n), b = rnorm(n), x = rnorm(n),
    c = factor(sample(c("No", "Yes"), n, TRUE), levels = c("No", "Yes"))
  )
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-.4 + .5*d$a + .4*d$b + .3*d$x + .4*(d$c=="Yes"))) == 1,
                       "Yes", "No"), levels = c("No", "Yes"))
  z <- tabforest(y, predictors = vars(c.a, c.b, c), data = d, event = "Yes",
                 adjusted = vars(c.x), multi = vars(c.a, c.b, c, c.x), show = FALSE)
  expect_true(all(c("Crude", "Adjusted", "Multivariable") %in% unique(z$table$model)))
  expect_true(is.list(z$models$Adjusted))
  expect_s3_class(z$models$Multivariable, "glm")
})

test_that("modified Poisson produces PR and RR with positive limits", {
  set.seed(3)
  n <- 300
  d <- data.frame(age = rnorm(n, 45, 10), smoke = factor(sample(c("No", "Yes"), n, TRUE)))
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-1.3 + .03*d$age + .5*(d$smoke=="Yes"))) == 1, "Yes", "No"),
                levels = c("No", "Yes"))
  pr <- tabforest(y, predictors = vars(c.age, smoke), data = d, event = "Yes", pr = TRUE, show = FALSE)
  rr <- tabforest(y, predictors = vars(c.age, smoke), data = d, event = "Yes", rr = TRUE, show = FALSE)
  expect_equal(pr$effect, "PR")
  expect_equal(rr$effect, "RR")
  expect_true(all(pr$table$lower[!pr$table$reference] > 0, na.rm = TRUE))
})

test_that("continuous outcome produces beta with null zero", {
  set.seed(4)
  n <- 180
  d <- data.frame(age = rnorm(n, 50, 10), sex = factor(sample(c("F", "M"), n, TRUE)))
  d$sbp <- 100 + .8*d$age + 5*(d$sex=="M") + rnorm(n, 0, 8)
  z <- tabforest(sbp, predictors = vars(c.age, sex), data = d, multi = TRUE, show = FALSE)
  expect_equal(z$effect, "Beta")
  expect_equal(z$null, 0)
  expect_equal(z$scale, "linear")
})

test_that("count outcome produces IRR", {
  set.seed(5)
  n <- 180
  d <- data.frame(age = rnorm(n, 45, 9), trt = factor(sample(c("C", "T"), n, TRUE)))
  d$count <- rpois(n, exp(.2 + .01*(d$age-45) - .25*(d$trt=="T")))
  z <- tabforest(count, predictors = vars(c.age, trt), data = d, irr = TRUE, multi = TRUE, show = FALSE)
  expect_equal(z$effect, "IRR")
  expect_equal(z$null, 1)
  expect_equal(z$scale, "log")
})

test_that("per scaling changes continuous ratio effect correctly", {
  set.seed(6)
  n <- 400
  d <- data.frame(age = rnorm(n, 50, 10))
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-2 + .04*d$age)) == 1, "Yes", "No"), levels = c("No", "Yes"))
  a <- tabforest(y, predictors = vars(c.age), data = d, event = "Yes", per = c(age = 1), show = FALSE)
  b <- tabforest(y, predictors = vars(c.age), data = d, event = "Yes", per = c(age = 10), show = FALSE)
  e1 <- a$table$estimate[a$table$variable == "age"][1]
  e10 <- b$table$estimate[b$table$variable == "age"][1]
  expect_equal(e10, e1^10, tolerance = 1e-6)
})

test_that("plot limits do not alter exact numeric table", {
  set.seed(7)
  n <- 220
  d <- data.frame(x = factor(sample(c("No", "Yes"), n, TRUE), levels = c("No", "Yes")))
  d$y <- factor(ifelse(rbinom(n, 1, ifelse(d$x=="Yes", .55, .08)) == 1, "Yes", "No"), levels = c("No", "Yes"))
  z1 <- tabforest(y, predictors = vars(x), data = d, event = "Yes", show = FALSE)
  z2 <- tabforest(y, predictors = vars(x), data = d, event = "Yes", xmax = 2, show = FALSE)
  expect_equal(z1$table$estimate, z2$table$estimate)
  expect_equal(z1$table$upper, z2$table$upper)
  expect_equal(z2$settings$xmax, 2)
})

test_that("Vietnamese display text and label overrides are stored", {
  set.seed(8)
  n <- 120
  d <- data.frame(sex = factor(sample(c("Female", "Male"), n, TRUE)))
  d$y <- factor(sample(c("No", "Yes"), n, TRUE), levels = c("No", "Yes"))
  z <- tabforest(y, predictors = vars(sex), data = d, event = "Yes", lang = "vi",
                 labels = c(sex = "Gi\u1edbi t\u00ednh"),
                 level_labels = list(sex = c(Female = "N\u1eef", Male = "Nam")),
                 text = list(reference = "Nh\u00f3m tham chi\u1ebfu"), show = FALSE)
  expect_equal(z$settings$text_resolved$reference, "Nh\u00f3m tham chi\u1ebfu")
  expect_true("Gi\u1edbi t\u00ednh" %in% z$rows$label)
  expect_true(any(z$rows$label %in% c("N\u1eef", "Nam")))
})

test_that("common sample aligns N across model groups", {
  set.seed(9)
  n <- 200
  d <- data.frame(a = rnorm(n), b = rnorm(n), x = rnorm(n))
  d$b[1:20] <- NA
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-.3 + .5*d$a + .3*ifelse(is.na(d$b),0,d$b))) == 1, "Yes", "No"),
                levels = c("No", "Yes"))
  z <- tabforest(y, predictors = vars(c.a, c.b), data = d, event = "Yes",
                 multi = TRUE, sample = "common", show = FALSE)
  ns <- unique(z$table$n[is.finite(z$table$n)])
  expect_length(ns, 1L)
})

test_that("Cox mode returns HR when survival is installed", {
  skip_if_not_installed("survival")
  set.seed(10)
  n <- 180
  d <- data.frame(age = rnorm(n, 55, 10), trt = factor(sample(c("C", "T"), n, TRUE)))
  h <- .12 * exp(.03*(d$age-55) - .4*(d$trt=="T"))
  te <- rexp(n, h); tc <- runif(n, 1, 7)
  d$time <- pmin(te, tc); d$death <- as.integer(te <= tc)
  z <- tabforest(death, time = time, predictors = vars(c.age, trt), data = d,
                 failure = 1, multi = TRUE, show = FALSE)
  expect_equal(z$effect, "HR")
  expect_equal(z$null, 1)
  expect_true(all(z$table$estimate[!z$table$reference] > 0, na.rm = TRUE))
})

test_that("bN reference syntax is retained in the forest", {
  set.seed(11)
  n <- 240
  d <- data.frame(
    sex = factor(sample(c("Female", "Male"), n, TRUE), levels = c("Female", "Male")),
    age = rnorm(n, 48, 11)
  )
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-2 + .04*d$age + .5*(d$sex=="Male"))) == 1,
                       "Yes", "No"), levels = c("No", "Yes"))
  z <- tabforest(y, predictors = vars(b2.sex, c.age), data = d,
                 event = "Yes", show = FALSE)
  ref <- z$table[z$table$variable == "sex" & z$table$reference, , drop = FALSE]
  expect_equal(ref$level[1], "Male")
  expect_equal(ref$estimate[1], 1)
})

test_that("standard glm objects can be converted without refitting", {
  set.seed(12)
  n <- 220
  d <- data.frame(age = rnorm(n, 50, 10), smoke = factor(sample(c("No", "Yes"), n, TRUE)))
  d$y <- rbinom(n, 1, plogis(-2 + .04*d$age + .6*(d$smoke=="Yes")))
  m <- glm(y ~ age + smoke, family = binomial(), data = d)
  z <- tabforest(m, show = FALSE)
  expect_equal(z$effect, "OR")
  expect_s3_class(z$models$Model, "glm")
})

test_that("current R4VN r4vn_stat regression objects are accepted", {
  skip_if_not(exists("logistic", mode = "function"))
  set.seed(13)
  n <- 220
  d <- data.frame(age = rnorm(n, 50, 10), smoke = factor(sample(c("No", "Yes"), n, TRUE)))
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-2 + .04*d$age + .6*(d$smoke=="Yes"))) == 1,
                       "Yes", "No"), levels = c("No", "Yes"))
  m <- logistic(y ~ age + smoke, data = d, event = "Yes", or = TRUE, show = FALSE)
  expect_true(inherits(m, c("r4vn_stat", "r4vn_result")))
  z <- tabforest(m, show = FALSE)
  expect_equal(z$effect, "OR")
  expect_true(all(z$table$estimate[!z$table$reference] > 0, na.rm = TRUE))
})

test_that("publication-ready data has one row per forest display row", {
  set.seed(14)
  n <- 160
  d <- data.frame(
    age = rnorm(n, 50, 10),
    sex = factor(sample(c("Female", "Male"), n, TRUE),
                 levels = c("Female", "Male"))
  )
  d$y <- factor(
    ifelse(rbinom(n, 1, plogis(-2 + .04 * d$age + .5 * (d$sex == "Male"))) == 1,
           "Yes", "No"),
    levels = c("No", "Yes")
  )
  z <- tabforest(y, predictors = vars(c.age, sex), data = d,
                 event = "Yes", show = FALSE)
  expect_equal(nrow(z$data), nrow(z$rows))
  expect_gt(nrow(z$data), 0L)
  expect_true(any(grepl("OR", names(z$data), fixed = TRUE)))
})

test_that("global omnibus p is hidden by default and can be requested", {
  set.seed(21)
  n <- 260
  d <- data.frame(
    age = rnorm(n, 50, 10),
    sex = factor(sample(c("F", "M"), n, TRUE), levels = c("F", "M"))
  )
  d$y <- factor(ifelse(rbinom(n, 1, plogis(-2 + .04*d$age + .5*(d$sex=="M"))) == 1,
                       "Yes", "No"), levels = c("No", "Yes"))
  a <- tabforest(y, predictors = vars(c.age, sex), data = d, event = "Yes", show = FALSE)
  b <- tabforest(y, predictors = vars(c.age, sex), data = d, event = "Yes", global_p = TRUE, show = FALSE)
  expect_false(a$settings$global_p)
  expect_true(b$settings$global_p)
})

test_that("modelrows layout can be redrawn without refitting", {
  set.seed(22)
  n <- 260
  d <- data.frame(
    age = rnorm(n, 50, 10),
    sex = factor(sample(c("F", "M"), n, TRUE), levels = c("F", "M"))
  )
  d$y <- factor(
    ifelse(
      rbinom(n, 1, plogis(-2 + .04*d$age + .5*(d$sex == "M"))) == 1,
      "Yes", "No"
    ),
    levels = c("No", "Yes")
  )

  z <- tabforest(
    y, predictors = vars(c.age, sex), data = d, event = "Yes",
    crude = TRUE, multi = TRUE, show = FALSE
  )

  tf <- tempfile(fileext = ".pdf")
  on.exit(unlink(tf, force = TRUE), add = TRUE)

  # The contract tested here is redraw + successful file output.
  # Graphics-device/font warnings can differ across operating systems.
  expect_no_error(
    suppressWarnings(
      plot(
        z,
        row_layout = "modelrows",
        zebra = TRUE,
        colors = c(Crude = "gray50", Multivariable = "black"),
        pch = c(Crude = 1, Multivariable = 16),
        file = tf
      )
    )
  )

  expect_true(file.exists(tf))
  expect_gt(unname(file.info(tf)$size), 0)
})

test_that("multi-outcome forest combines panels with a shared row structure", {
  set.seed(23)
  n <- 360
  d <- data.frame(
    age = rnorm(n, 50, 10),
    sex = factor(sample(c("F", "M"), n, TRUE), levels = c("F", "M"))
  )
  d$y1 <- factor(ifelse(rbinom(n, 1, plogis(-2 + .04*d$age + .4*(d$sex=="M"))) == 1,
                        "Yes", "No"), levels = c("No", "Yes"))
  d$y2 <- factor(ifelse(rbinom(n, 1, plogis(-1.5 + .02*d$age - .3*(d$sex=="M"))) == 1,
                        "Yes", "No"), levels = c("No", "Yes"))
  z <- tabforest(outcomes = c(Outcome1 = "y1", Outcome2 = "y2"),
                 predictors = vars(c.age, sex), data = d, event = "Yes",
                 crude = FALSE, multi = TRUE, show = FALSE)
  expect_s3_class(z, "r4vn_tabforest_multi")
  expect_length(z$panels, 2L)
  expect_true(all(c("Outcome1", "Outcome2") %in% unique(z$table$outcome_panel)))
  expect_equal(nrow(z$data), nrow(z$rows))
})

test_that("subgroup forest returns level effects and interaction p fields", {
  set.seed(24)
  n <- 700
  d <- data.frame(
    trt = factor(sample(c("Control", "Intervention"), n, TRUE),
                 levels = c("Control", "Intervention")),
    sex = factor(sample(c("Female", "Male"), n, TRUE), levels = c("Female", "Male")),
    smoke = factor(sample(c("No", "Yes"), n, TRUE), levels = c("No", "Yes")),
    age = rnorm(n, 52, 11)
  )
  lp <- -2 + .35*(d$trt=="Intervention") + .03*(d$age-50) + .3*(d$smoke=="Yes")
  d$y <- factor(ifelse(rbinom(n, 1, plogis(lp)) == 1, "Yes", "No"),
                levels = c("No", "Yes"))
  z <- tabforest(y, predictor = vars(trt), subgroup = vars(sex, smoke),
                 data = d, event = "Yes", type = "subgroup",
                 adjusted = vars(c.age), show = FALSE)
  expect_s3_class(z, "r4vn_tabforest_subgroup")
  expect_true(all(c("subgroup_variable", "subgroup_level", "interaction_p") %in% names(z$table)))
  expect_true(all(c("sex", "smoke") %in% unique(z$table$subgroup_variable)))
})

test_that("subgroup Cox mode returns HR", {
  skip_if_not_installed("survival")
  set.seed(25)
  n <- 500
  d <- data.frame(
    trt = factor(sample(c("Control", "Intervention"), n, TRUE),
                 levels = c("Control", "Intervention")),
    sex = factor(sample(c("Female", "Male"), n, TRUE), levels = c("Female", "Male")),
    age = rnorm(n, 55, 10)
  )
  rate <- .12 * exp(-.35*(d$trt=="Intervention") + .02*(d$age-55))
  te <- rexp(n, rate); tc <- runif(n, 1, 7)
  d$time <- pmin(te, tc); d$death <- as.integer(te <= tc)
  z <- tabforest(death, time = time, predictor = vars(trt), subgroup = vars(sex),
                 data = d, failure = 1, type = "subgroup",
                 adjusted = vars(c.age), estimate = "hr", show = FALSE)
  expect_equal(z$effect, "HR")
  expect_true(all(z$table$estimate[is.finite(z$table$estimate)] > 0))
})