diff --git a/NEWS.md b/NEWS.md index c632365bbb..22fbdb0d55 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,7 @@ # tern 0.9.10.9006 +* Added `alternative` argument to `s_coxph_pairwise()` to allow one-sided hypothesis testing and added `lr_stat_df` to the parameters return list + # tern 0.9.10 ### Enhancements diff --git a/R/survival_coxph_pairwise.R b/R/survival_coxph_pairwise.R index 742700ec9b..d02511249b 100644 --- a/R/survival_coxph_pairwise.R +++ b/R/survival_coxph_pairwise.R @@ -47,6 +47,7 @@ s_coxph_pairwise <- function(df, strata = NULL, strat = lifecycle::deprecated(), control = control_coxph(), + alternative = c("two.sided", "less", "greater"), ...) { if (lifecycle::is_present(strat)) { lifecycle::deprecate_warn("0.9.4", "s_coxph_pairwise(strat)", "s_coxph_pairwise(strata)") @@ -57,6 +58,8 @@ s_coxph_pairwise <- function(df, checkmate::assert_numeric(df[[.var]]) checkmate::assert_logical(df[[is_event]]) assert_df_with_variables(df, list(tte = .var, is_event = is_event)) + alternative <- match.arg(alternative) + pval_method <- control$pval_method ties <- control$ties conf_level <- control$conf_level @@ -65,6 +68,7 @@ s_coxph_pairwise <- function(df, return( list( pvalue = formatters::with_label(numeric(), paste0("p-value (", pval_method, ")")), + lr_stat_df = formatters::with_label(numeric(), "Log-rank Degrees of freedom"), hr = formatters::with_label(numeric(), "Hazard Ratio"), hr_ci = formatters::with_label(numeric(), f_conf_level(conf_level)), hr_ci_3d = formatters::with_label(numeric(), paste0("Hazard Ratio (", f_conf_level(conf_level), ")")), @@ -99,19 +103,36 @@ s_coxph_pairwise <- function(df, ties = ties ) sum_cox <- summary(cox_fit, conf.int = conf_level, extend = TRUE) - orginal_survdiff <- survival::survdiff( - formula_cox, - data = df_cox - ) - log_rank_pvalue <- 1 - pchisq(orginal_survdiff$chisq, length(orginal_survdiff$n) - 1) + original_survdiff <- survival::survdiff(formula_cox, data = df_cox) + log_rank_stat <- original_survdiff$chisq + log_rank_df <- length(original_survdiff$n) - 1 + log_rank_pvalue <- original_survdiff$pvalue pval <- switch(pval_method, "wald" = sum_cox$waldtest["pvalue"], "log-rank" = log_rank_pvalue, # pvalue from original log-rank test survival::survdiff() "likelihood" = sum_cox$logtest["pvalue"] ) + + # Handle one-sided alternatives. + if (alternative != "two.sided") { + # Need to calculate the signed log-rank statistic, which is not included + # in the original survdiff output. + otmp <- rowSums(original_survdiff$obs) + etmp <- rowSums(original_survdiff$exp) + signed_lr_stat <- (otmp[2] - etmp[2]) / sqrt(original_survdiff$var[2, 2]) + + pval <- if (alternative == "less") { + stats::pnorm(signed_lr_stat) + } else { + stats::pnorm(signed_lr_stat, lower.tail = FALSE) + } + } + + list( pvalue = formatters::with_label(unname(pval), paste0("p-value (", pval_method, ")")), + lr_stat_df = formatters::with_label(unname(c(log_rank_stat, log_rank_df)), "Log-rank Degrees of freedom"), hr = formatters::with_label(sum_cox$conf.int[1, 1], "Hazard Ratio"), hr_ci = formatters::with_label(unname(sum_cox$conf.int[1, 3:4]), f_conf_level(conf_level)), hr_ci_3d = formatters::with_label( diff --git a/man/survival_coxph_pairwise.Rd b/man/survival_coxph_pairwise.Rd index e03fc6c7b5..02555aeb32 100644 --- a/man/survival_coxph_pairwise.Rd +++ b/man/survival_coxph_pairwise.Rd @@ -34,6 +34,7 @@ s_coxph_pairwise( strata = NULL, strat = lifecycle::deprecated(), control = control_coxph(), + alternative = c("two.sided", "less", "greater"), ... ) @@ -107,6 +108,9 @@ by a statistics function.} \item{is_event}{(\code{flag})\cr \code{TRUE} if event, \code{FALSE} if time to event is censored.} \item{strat}{\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#deprecated}{\figure{lifecycle-deprecated.svg}{options: alt='[Deprecated]'}}}{\strong{[Deprecated]}} Please use the \code{strata} argument instead.} + +\item{alternative}{(\code{string})\cr whether \code{two.sided}, or one-sided \code{less} or \code{greater} p-value +should be displayed.} } \value{ \itemize{ diff --git a/tests/testthat/_snaps/survival_coxph_pairwise.md b/tests/testthat/_snaps/survival_coxph_pairwise.md index fc27a5f849..26f0589398 100644 --- a/tests/testthat/_snaps/survival_coxph_pairwise.md +++ b/tests/testthat/_snaps/survival_coxph_pairwise.md @@ -8,6 +8,11 @@ attr(,"label") [1] "p-value (log-rank)" + $lr_stat_df + [1] 2.865544 1.000000 + attr(,"label") + [1] "Log-rank Degrees of freedom" + $hr [1] 0.7108557 attr(,"label") @@ -44,6 +49,11 @@ attr(,"label") [1] "p-value (wald)" + $lr_stat_df + [1] 2.865544 1.000000 + attr(,"label") + [1] "Log-rank Degrees of freedom" + $hr [1] 0.7108557 attr(,"label") @@ -80,6 +90,11 @@ attr(,"label") [1] "p-value (log-rank)" + $lr_stat_df + [1] 4.390702 1.000000 + attr(,"label") + [1] "Log-rank Degrees of freedom" + $hr [1] 0.6251817 attr(,"label") @@ -116,6 +131,11 @@ attr(,"label") [1] "p-value (wald)" + $lr_stat_df + [1] 4.390702 1.000000 + attr(,"label") + [1] "Log-rank Degrees of freedom" + $hr [1] 0.6251817 attr(,"label") @@ -183,21 +203,22 @@ Code res Output - ARM A ARM B ARM C - ————————————————————————————————————————————————————————————— - Stratified Analysis - Hazard Ratio 1.600 2.049 - 99% CI (0.894, 2.863) (1.092, 3.844) + ARM A ARM B ARM C + ———————————————————————————————————————————————————————————————— + Stratified Analysis + Hazard Ratio 1.600 2.049 + 99% CI (0.894, 2.863) (1.092, 3.844) + p-value (likelihood) 0.0181 0.0014 # coxph_pairwise works with NA values Code result Output - ARM A ARM B ARM C - ——————————————————————————————————————————————— - Unstratified Analysis - p-value (log-rank) 1.0000 1.0000 - Hazard Ratio empty empty - 95% CI empty empty + ARM A ARM B ARM C + ————————————————————————————————————————————— + Unstratified Analysis + p-value (log-rank) empty empty + Hazard Ratio empty empty + 95% CI empty empty diff --git a/tests/testthat/test-survival_coxph_pairwise.R b/tests/testthat/test-survival_coxph_pairwise.R index 087c0eb64d..214b6ea03b 100644 --- a/tests/testthat/test-survival_coxph_pairwise.R +++ b/tests/testthat/test-survival_coxph_pairwise.R @@ -80,6 +80,44 @@ testthat::test_that("s_coxph_pairwise works with customized arguments and strati testthat::expect_snapshot(res) }) +testthat::test_that("s_coxph_pairwise works with stratification factors for Log-Rank test", { + adtte_f <- tern_ex_adtte %>% + dplyr::filter(PARAMCD == "OS") %>% + dplyr::mutate(is_event = CNSR == 0) + df <- adtte_f %>% dplyr::filter(ARMCD == "ARM A") + df_ref <- adtte_f %>% dplyr::filter(ARMCD == "ARM B") + + # default control uses pval_method = "log-rank" + result <- s_coxph_pairwise( + df = df, + .ref_group = df_ref, + .in_ref_col = FALSE, + .var = "AVAL", + is_event = "is_event", + strata = c("SEX", "RACE") + ) + + testthat::expect_silent(result) + testthat::expect_true("lr_stat_df" %in% names(result)) + testthat::expect_type(result$lr_stat_df, "double") + testthat::expect_length(result$lr_stat_df, 2) + testthat::expect_identical(attr(result$lr_stat_df, "label"), "Log-rank Degrees of freedom") + + # Check the consistency of the d.f. with the p-value returned by survival::survdiff. + log_rank_pvalue <- stats::pchisq( + result$lr_stat_df[1], + result$lr_stat_df[2], + lower.tail = FALSE + ) + original_survdiff <- survival::survdiff( + survival::Surv(AVAL, is_event) ~ ARMCD + strata(SEX, RACE), + data = adtte_f %>% + dplyr::filter(ARMCD %in% c("ARM A", "ARM B")) %>% + droplevels() + ) + testthat::expect_equal(log_rank_pvalue, original_survdiff$pvalue) +}) + testthat::test_that("coxph_pairwise works with default arguments and no stratification factors", { adtte_f <- tern_ex_adtte %>% dplyr::filter(PARAMCD == "OS") %>% @@ -163,7 +201,8 @@ testthat::test_that("coxph_pairwise works with customized arguments and stratifi var_labels = c("Stratified Analysis"), control = control_coxph(pval_method = "likelihood", conf_level = 0.99), strata = c("SEX", "RACE"), - .stats = c("hr", "hr_ci"), + alternative = "greater", + .stats = c("hr", "hr_ci", "pvalue"), .formats = c(hr = "xx.xxx", hr_ci = "(xx.xxx, xx.xxx)") ) %>% build_table(df = adtte_f)