Add additional lr_stat_df statistic

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

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(

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
 

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)