From d3287697c78bf79567ad0de6ff96ce6a5790b6ef Mon Sep 17 00:00:00 2001 From: Daniel Sabanes Bove Date: Fri, 10 Jul 2026 21:03:03 +0800 Subject: [PATCH 1/4] s_proportion_diff: make the standard error available in the results --- R/prop_diff.R | 15 ++++++++++++--- man/prop_diff.Rd | 2 ++ man/tern-package.Rd | 1 + tests/testthat/_snaps/prop_diff.md | 15 ++++++++++++++- 4 files changed, 29 insertions(+), 4 deletions(-) diff --git a/R/prop_diff.R b/R/prop_diff.R index 631331355b..1423f669d0 100644 --- a/R/prop_diff.R +++ b/R/prop_diff.R @@ -47,6 +47,8 @@ NULL #' #' @return #' * `s_proportion_diff()` returns a named list of elements `diff` and `diff_ci`. +#' Depending on the method used, also the standard error of the difference `se_diff` is +#' returned. #' #' @note When performing an unstratified analysis, methods `"cmh"`, `"cmh_sato"`, `"strat_newcombe"`, #' and `"strat_newcombecc"` are not permitted. For stratified analysis, method @@ -138,6 +140,7 @@ s_proportion_diff <- function(df, weights_method <- "cmh" } + y_stats_from_cmh <- c("diff", "diff_ci", "se_diff") y <- switch(method, "wald" = prop_diff_wald(rsp, grp, conf_level, correct = FALSE), "waldcc" = prop_diff_wald(rsp, grp, conf_level, correct = TRUE), @@ -158,14 +161,17 @@ s_proportion_diff <- function(df, conf_level, correct = TRUE ), - "cmh" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "standard")[c("diff", "diff_ci")], - "cmh_sato" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "sato")[c("diff", "diff_ci")], - "cmh_mn" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "miettinen_nurminen")[c("diff", "diff_ci")], + "cmh" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "standard")[y_stats_from_cmh], + "cmh_sato" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "sato")[y_stats_from_cmh], + "cmh_mn" = prop_diff_cmh(rsp, grp, strata, conf_level, diff_se = "miettinen_nurminen")[y_stats_from_cmh], "uncond_exact_diff" = prop_diff_uncond_exact(rsp, grp, conf_level) ) y$diff <- setNames(y$diff * 100, paste0("diff_", method)) y$diff_ci <- setNames(y$diff_ci * 100, paste0("diff_ci_", method, c("_l", "_u"))) + if (!is.null(y$se_diff)) { + y$se_diff <- setNames(y$se_diff * 100, paste0("se_diff_", method)) + } } attr(y$diff, "label") <- "Difference in Response rate (%)" @@ -173,6 +179,9 @@ s_proportion_diff <- function(df, conf_level, method, long = FALSE ) + if (!is.null(y$se_diff)) { + attr(y$se_diff, "label") <- paste0("Standard Error of Difference in Response rate (%)") + } y } diff --git a/man/prop_diff.Rd b/man/prop_diff.Rd index e43d1d2a83..2a0ce77ec7 100644 --- a/man/prop_diff.Rd +++ b/man/prop_diff.Rd @@ -120,6 +120,8 @@ the statistics from \code{s_proportion_diff()} to the table layout. \itemize{ \item \code{s_proportion_diff()} returns a named list of elements \code{diff} and \code{diff_ci}. +Depending on the method used, also the standard error of the difference \code{se_diff} is +returned. } \itemize{ diff --git a/man/tern-package.Rd b/man/tern-package.Rd index 1bbdae1894..b75b0050ec 100644 --- a/man/tern-package.Rd +++ b/man/tern-package.Rd @@ -38,6 +38,7 @@ Authors: Other contributors: \itemize{ + \item David Munoz Tord \email{david.munoztord@mailbox.org} [contributor] \item F. Hoffmann-La Roche AG [copyright holder, funder] } diff --git a/tests/testthat/_snaps/prop_diff.md b/tests/testthat/_snaps/prop_diff.md index d8c6f80b82..30e8dd0f4e 100644 --- a/tests/testthat/_snaps/prop_diff.md +++ b/tests/testthat/_snaps/prop_diff.md @@ -348,6 +348,12 @@ attr(,"label") [1] "90% CI (CMH, without correction)" + $se_diff + se_diff_cmh + 8.978092 + attr(,"label") + [1] "Standard Error of Difference in Response rate (%)" + # s_proportion_diff works with CMH Sato method @@ -366,10 +372,17 @@ attr(,"label") [1] "90% CI (CMH, Sato variance estimator)" + $se_diff + se_diff_cmh_sato + 10.80533 + attr(,"label") + [1] "Standard Error of Difference in Response rate (%)" + # s_proportion_diff works with CMH Miettinen and Nurminen method list(diff = structure(c(diff_cmh_mn = 13.7686601988347), label = "Difference in Response rate (%)"), diff_ci = structure(c(diff_ci_cmh_mn_l = -3.45069418895496, - diff_ci_cmh_mn_u = 30.2144371774115), label = "90% CI (CMH, Miettinen and Nurminen)")) + diff_ci_cmh_mn_u = 30.2144371774115), label = "90% CI (CMH, Miettinen and Nurminen)"), + se_diff = structure(c(se_diff_cmh_mn = 10.4103330371023), label = "Standard Error of Difference in Response rate (%)")) From e05d6733129ec839291dd16668fcebb8e54136d3 Mon Sep 17 00:00:00 2001 From: Daniel Sabanes Bove Date: Fri, 10 Jul 2026 21:28:08 +0800 Subject: [PATCH 2/4] export the underlying prop_* functions because they are lighter to use in loops etc. and make z statistic available from prop_cmh --- NAMESPACE | 4 ++ R/prop_diff.R | 4 +- R/prop_diff_test.R | 62 +++++++++++++++++-- man/h_prop_diff_test.Rd | 46 +++++++++++++- man/prop_diff.Rd | 4 +- tests/testthat/_snaps/test_proportion_diff.md | 26 +++++++- tests/testthat/test-test_proportion_diff.R | 44 ++++++++++--- 7 files changed, 169 insertions(+), 21 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 4effb5ab9a..0afd85824e 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -252,14 +252,18 @@ export(month2day) export(or_clogit) export(or_glm) export(prop_agresti_coull) +export(prop_chisq) export(prop_clopper_pearson) +export(prop_cmh) export(prop_diff_cmh) export(prop_diff_ha) export(prop_diff_nc) export(prop_diff_strat_nc) export(prop_diff_uncond_exact) export(prop_diff_wald) +export(prop_fisher) export(prop_jeffreys) +export(prop_schouten) export(prop_strat_wilson) export(prop_wald) export(prop_wilson) diff --git a/R/prop_diff.R b/R/prop_diff.R index 1423f669d0..3d002a607a 100644 --- a/R/prop_diff.R +++ b/R/prop_diff.R @@ -323,9 +323,9 @@ estimate_proportion_diff <- function(lyt, na_rm = TRUE, .stats = c("diff", "diff_ci"), .stat_names = NULL, - .formats = c(diff = "xx.x", diff_ci = "(xx.x, xx.x)"), + .formats = c(diff = "xx.x", diff_ci = "(xx.x, xx.x)", se_diff = "xx.x"), .labels = NULL, - .indent_mods = c(diff = 0L, diff_ci = 1L)) { + .indent_mods = c(diff = 0L, diff_ci = 1L, se_diff = 1L)) { # Depending on main functions extra_args <- list( "na_rm" = na_rm, diff --git a/R/prop_diff_test.R b/R/prop_diff_test.R index 12f1b97796..f36647379d 100644 --- a/R/prop_diff_test.R +++ b/R/prop_diff_test.R @@ -308,7 +308,18 @@ NULL #' @describeIn h_prop_diff_test Performs Chi-Squared test. Internally calls [stats::prop.test()]. #' -#' @keywords internal +#' @examples +#' # Chi-Squared test +#' tbl <- matrix( +#' c(13, 7, 8, 12), +#' nrow = 2, +#' byrow = TRUE, +#' dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +#' ) +#' +#' prop_chisq(tbl) +#' +#' @export prop_chisq <- function(tbl, alternative = c("two.sided", "less", "greater")) { checkmate::assert_integer(c(ncol(tbl), nrow(tbl)), lower = 2, upper = 2) tbl <- tbl[, c("TRUE", "FALSE")] @@ -329,7 +340,21 @@ prop_chisq <- function(tbl, alternative = c("two.sided", "less", "greater")) { #' @param transform (`string`)\cr either `none` or `wilson_hilferty`; specifies whether to apply #' the Wilson-Hilferty transformation of the chi-squared statistic. #' -#' @keywords internal +#' @examples +#' # Cochran-Mantel-Haenszel test with two strata +#' ary <- array( +#' c(12, 8, 8, 12, 10, 10, 6, 14), +#' dim = c(2, 2, 2), +#' dimnames = list( +#' group = c("A", "B"), +#' response = c("TRUE", "FALSE"), +#' strata = c("Low", "High") +#' ) +#' ) +#' +#' prop_cmh(ary) +#' +#' @export prop_cmh <- function(ary, alternative = c("two.sided", "less", "greater"), diff_se = c("standard", "sato"), @@ -375,18 +400,34 @@ prop_cmh <- function(ary, z_stat <- num / denom * sign(z_stat) # Preserve the direction of effect. } - if (alternative == "two.sided") { + pval <- if (alternative == "two.sided") { 2 * stats::pnorm(-abs(z_stat)) } else { stats::pnorm(z_stat, lower.tail = (alternative == "greater")) } + + structure( + pval, + z_stat = z_stat + ) } #' @describeIn h_prop_diff_test Performs the Chi-Squared test with Schouten correction. #' #' @seealso Schouten correction is based upon \insertCite{Schouten1980-kd;textual}{tern}. #' -#' @keywords internal +#' @examples +#' # Chi-Squared test with Schouten correction +#' tbl <- matrix( +#' c(13, 7, 8, 12), +#' nrow = 2, +#' byrow = TRUE, +#' dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +#' ) +#' +#' prop_schouten(tbl) +#' +#' @export prop_schouten <- function(tbl, alternative = c("two.sided", "less", "greater")) { checkmate::assert_integer(c(ncol(tbl), nrow(tbl)), lower = 2, upper = 2) alternative <- match.arg(alternative) @@ -422,7 +463,18 @@ prop_schouten <- function(tbl, alternative = c("two.sided", "less", "greater")) #' @describeIn h_prop_diff_test Performs the Fisher's exact test. Internally calls [stats::fisher.test()]. #' -#' @keywords internal +#' @examples +#' # Fisher's exact test +#' tbl <- matrix( +#' c(13, 7, 8, 12), +#' nrow = 2, +#' byrow = TRUE, +#' dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +#' ) +#' +#' prop_fisher(tbl) +#' +#' @export prop_fisher <- function(tbl, alternative = c("two.sided", "less", "greater")) { checkmate::assert_integer(c(ncol(tbl), nrow(tbl)), lower = 2, upper = 2) alternative <- match.arg(alternative) # Is needed here, because stats::fisher.test does not handle defaults. diff --git a/man/h_prop_diff_test.Rd b/man/h_prop_diff_test.Rd index eed793efa6..6ba35f0038 100644 --- a/man/h_prop_diff_test.Rd +++ b/man/h_prop_diff_test.Rd @@ -55,9 +55,53 @@ Note that strata with less than two observations are automatically discarded. \item \code{prop_fisher()}: Performs the Fisher's exact test. Internally calls \code{\link[stats:fisher.test]{stats::fisher.test()}}. }} +\examples{ +# Chi-Squared test +tbl <- matrix( + c(13, 7, 8, 12), + nrow = 2, + byrow = TRUE, + dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +) + +prop_chisq(tbl) + +# Cochran-Mantel-Haenszel test with two strata +ary <- array( + c(12, 8, 8, 12, 10, 10, 6, 14), + dim = c(2, 2, 2), + dimnames = list( + group = c("A", "B"), + response = c("TRUE", "FALSE"), + strata = c("Low", "High") + ) +) + +prop_cmh(ary) + +# Chi-Squared test with Schouten correction +tbl <- matrix( + c(13, 7, 8, 12), + nrow = 2, + byrow = TRUE, + dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +) + +prop_schouten(tbl) + +# Fisher's exact test +tbl <- matrix( + c(13, 7, 8, 12), + nrow = 2, + byrow = TRUE, + dimnames = list(group = c("A", "B"), response = c("TRUE", "FALSE")) +) + +prop_fisher(tbl) + +} \seealso{ \code{\link[=prop_diff_test]{prop_diff_test()}} for implementation of these helper functions. Schouten correction is based upon \insertCite{Schouten1980-kd;textual}{tern}. } -\keyword{internal} diff --git a/man/prop_diff.Rd b/man/prop_diff.Rd index 2a0ce77ec7..fc31e622e3 100644 --- a/man/prop_diff.Rd +++ b/man/prop_diff.Rd @@ -25,9 +25,9 @@ estimate_proportion_diff( na_rm = TRUE, .stats = c("diff", "diff_ci"), .stat_names = NULL, - .formats = c(diff = "xx.x", diff_ci = "(xx.x, xx.x)"), + .formats = c(diff = "xx.x", diff_ci = "(xx.x, xx.x)", se_diff = "xx.x"), .labels = NULL, - .indent_mods = c(diff = 0L, diff_ci = 1L) + .indent_mods = c(diff = 0L, diff_ci = 1L, se_diff = 1L) ) s_proportion_diff( diff --git a/tests/testthat/_snaps/test_proportion_diff.md b/tests/testthat/_snaps/test_proportion_diff.md index 9b3cf7690b..e0002a3519 100644 --- a/tests/testthat/_snaps/test_proportion_diff.md +++ b/tests/testthat/_snaps/test_proportion_diff.md @@ -25,6 +25,8 @@ res Output [1] 0.6477165 + attr(,"z_stat") + [1] 0.4569368 --- @@ -32,6 +34,8 @@ res_less Output [1] 0.3238582 + attr(,"z_stat") + [1] 0.4569368 --- @@ -39,6 +43,8 @@ res_greater Output [1] 0.6761418 + attr(,"z_stat") + [1] 0.4569368 # prop_cmh returns right result for odds ratio > 1 @@ -46,6 +52,8 @@ res Output [1] 0.3833285 + attr(,"z_stat") + [1] -0.8717798 --- @@ -53,6 +61,8 @@ res_less Output [1] 0.8083357 + attr(,"z_stat") + [1] -0.8717798 --- @@ -60,6 +70,8 @@ res_greater Output [1] 0.1916643 + attr(,"z_stat") + [1] -0.8717798 # prop_cmh also works when there are strata with just one observation @@ -67,6 +79,8 @@ res Output [1] 0.3325724 + attr(,"z_stat") + [1] -0.9689455 # prop_fisher returns right result @@ -137,6 +151,8 @@ res_less Output [1] 0.6522653 + attr(,"z_stat") + [1] -0.3914435 --- @@ -144,10 +160,12 @@ res_greater Output [1] 0.3477347 + attr(,"z_stat") + [1] -0.3914435 # prop_cmh with Sato variance estimator and Wilson-Hilferty transformation works - 0.0660308573944943 + structure(0.0660308573944943, z_stat = -1.83821413735577) # s_test_proportion_diff and d_test_proportion_diff return right result @@ -160,6 +178,8 @@ $s $s$pval [1] 0.6477165 + attr(,"z_stat") + [1] 0.4569368 attr(,"label") [1] "p-value (Cochran-Mantel-Haenszel Test)" @@ -176,6 +196,8 @@ $s $s$pval [1] 0.6761418 + attr(,"z_stat") + [1] 0.4569368 attr(,"label") [1] "p-value (Cochran-Mantel-Haenszel Test, 1-sided, direction greater)" @@ -192,6 +214,8 @@ $s $s$pval [1] 0.3238582 + attr(,"z_stat") + [1] 0.4569368 attr(,"label") [1] "p-value (Cochran-Mantel-Haenszel Test, 1-sided, direction less)" diff --git a/tests/testthat/test-test_proportion_diff.R b/tests/testthat/test-test_proportion_diff.R index 8baca3d479..ac9de82c53 100644 --- a/tests/testthat/test-test_proportion_diff.R +++ b/tests/testthat/test-test_proportion_diff.R @@ -31,19 +31,27 @@ testthat::test_that("prop_cmh returns right result for odds ratio < 1", { result <- prop_cmh(tbl) res <- testthat::expect_silent(result) testthat::expect_snapshot(res) - expected <- stats::mantelhaen.test(tbl, correct = FALSE)$p.value + mh_result <- stats::mantelhaen.test(tbl, correct = FALSE) + expected_z_stat <- sqrt(unname(mh_result$statistic)) + expected <- structure(mh_result$p.value, z_stat = expected_z_stat) testthat::expect_equal(result, expected, tolerance = 1e-3) result_less <- prop_cmh(tbl, alternative = "less") res_less <- testthat::expect_silent(result_less) testthat::expect_snapshot(res_less) - expected_less <- stats::mantelhaen.test(tbl, correct = FALSE, alternative = "less")$p.value + expected_less <- structure( + stats::mantelhaen.test(tbl, correct = FALSE, alternative = "less")$p.value, + z_stat = expected_z_stat + ) testthat::expect_equal(result_less, expected_less, tolerance = 1e-3) result_greater <- prop_cmh(tbl, alternative = "greater") res_greater <- testthat::expect_silent(result_greater) testthat::expect_snapshot(res_greater) - expected_greater <- stats::mantelhaen.test(tbl, correct = FALSE, alternative = "greater")$p.value + expected_greater <- structure( + stats::mantelhaen.test(tbl, correct = FALSE, alternative = "greater")$p.value, + z_stat = expected_z_stat + ) testthat::expect_equal(result_greater, expected_greater, tolerance = 1e-3) }) @@ -58,19 +66,27 @@ testthat::test_that("prop_cmh returns right result for odds ratio > 1", { result <- prop_cmh(tbl) res <- testthat::expect_silent(result) testthat::expect_snapshot(res) - expected <- stats::mantelhaen.test(tbl, correct = FALSE)$p.value + mh_result <- stats::mantelhaen.test(tbl, correct = FALSE) + expected_z_stat <- -sqrt(unname(mh_result$statistic)) + expected <- structure(mh_result$p.value, z_stat = expected_z_stat) testthat::expect_equal(result, expected, tolerance = 1e-3) result_less <- prop_cmh(tbl, alternative = "less") res_less <- testthat::expect_silent(result_less) testthat::expect_snapshot(res_less) - expected_less <- stats::mantelhaen.test(tbl, correct = FALSE, alternative = "less")$p.value + expected_less <- structure( + stats::mantelhaen.test(tbl, correct = FALSE, alternative = "less")$p.value, + z_stat = expected_z_stat + ) testthat::expect_equal(result_less, expected_less, tolerance = 1e-3) result_greater <- prop_cmh(tbl, alternative = "greater") res_greater <- testthat::expect_silent(result_greater) testthat::expect_snapshot(res_greater) - expected_greater <- stats::mantelhaen.test(tbl, correct = FALSE, alternative = "greater")$p.value + expected_greater <- structure( + stats::mantelhaen.test(tbl, correct = FALSE, alternative = "greater")$p.value, + z_stat = expected_z_stat + ) testthat::expect_equal(result_greater, expected_greater, tolerance = 1e-3) }) @@ -191,8 +207,14 @@ testthat::test_that("prop_cmh with Wilson-Hilferty transformation works", { result_cmh_less <- prop_cmh(tbl, alternative = "less") result_cmh_greater <- prop_cmh(tbl, alternative = "greater") - expect_equal(result_less, result_cmh_less, tolerance = 1e-1) - expect_equal(result_greater, result_cmh_greater, tolerance = 1e-1) + expect_equal(as.numeric(result_less), as.numeric(result_cmh_less), tolerance = 1e-1) + expect_equal(as.numeric(result_greater), as.numeric(result_cmh_greater), tolerance = 1e-1) + expected_z_stat <- { + z_stat <- attr(result_cmh_less, "z_stat") + ((1 - 2 / 9) - (z_stat^2)^(1 / 3)) / sqrt(2 / 9) * sign(z_stat) + } + expect_equal(attr(result_less, "z_stat"), expected_z_stat) + expect_equal(attr(result_greater, "z_stat"), expected_z_stat) }) testthat::test_that("prop_cmh with Sato variance estimator works", { @@ -212,7 +234,8 @@ testthat::test_that("prop_cmh with Sato variance estimator works", { class = "table" ) result1 <- prop_cmh(tbl1, diff_se = "sato") - testthat::expect_equal(result1, 0.035, tolerance = 1e-3) + testthat::expect_equal(as.numeric(result1), 0.035, tolerance = 1e-3) + testthat::expect_equal(attr(result1, "z_stat"), 2.109, tolerance = 1e-3) tbl2 <- structure( c( @@ -230,7 +253,8 @@ testthat::test_that("prop_cmh with Sato variance estimator works", { class = "table" ) result2 <- prop_cmh(tbl2, diff_se = "sato") - testthat::expect_equal(result2, 0.004, tolerance = 1e-2) + testthat::expect_equal(as.numeric(result2), 0.004, tolerance = 1e-2) + testthat::expect_equal(attr(result2, "z_stat"), 2.864, tolerance = 1e-3) }) testthat::test_that("prop_cmh with Sato variance estimator and Wilson-Hilferty transformation works", { From 7d4017451c208e274fdb27f997253b0574aeb205 Mon Sep 17 00:00:00 2001 From: Daniel Sabanes Bove Date: Fri, 10 Jul 2026 21:34:03 +0800 Subject: [PATCH 3/4] add NEWS entries --- NEWS.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/NEWS.md b/NEWS.md index 00023b0fbc..4c9b83d86c 100644 --- a/NEWS.md +++ b/NEWS.md @@ -8,7 +8,9 @@ * Added `range_with_cens_info` statistic to `s_surv_time()`. * Added `lsmean_se`, `lsmean_ci`, and `lsmean_diffci` statistics to `s_ancova()`. * Added `s_ancova()` to exported functions. +* In `estimate_proportion_diff` now also the standard error is available for selected methods in the `se_diff` statistic. * Exported `h_incidence_rate()` to allow to reuse the incidence rate estimation logic. +* Exported `prop_chisq()`, `prop_cmh()`, `prop_fisher()`, and `prop_schouten()` for standalone proportion difference testing. ### Bug Fixes * Fixed bug in `prop_diff_cmh()` which previously failed when strata combinations had 0 observations. From 454398f14ea25b9a44f6120ca2f23c5b429e9d8b Mon Sep 17 00:00:00 2001 From: Daniel Sabanes Bove Date: Fri, 10 Jul 2026 21:37:43 +0800 Subject: [PATCH 4/4] update wordlist --- inst/WORDLIST | 1 + 1 file changed, 1 insertion(+) diff --git a/inst/WORDLIST b/inst/WORDLIST index b581d3c8c0..54dd768174 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -37,6 +37,7 @@ Satterthwaite Schouten TLG TLGs +Tord binom biomarker biomarkers