Enable different tests per cutting for sim_gs_n()

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: simtrial
 Type: Package
 Title: Clinical Trial Simulation
-Version: 0.4.1.1
+Version: 0.4.1.2
 Authors@R: c(
     person("Keaven", "Anderson", email = "[email protected]", role = c("aut")),
     person("Yujie", "Zhao", email = "[email protected]", role =  c("ctb","cre")),

R/maxcombo.R

@@ -35,7 +35,7 @@
 #'
 #' @return A list containing the test method (`method`),
 #' parameters of this test method (`parameter`),
-#' point estimation of the treatment effect (`estimation`),
+#' point estimate of the treatment effect (`estimate`),
 #' standardized error of the treatment effect (`se`),
 #' Z-score of each test of the MaxCombo (`z`),
 #' p-values (`p_value`)
@@ -107,10 +107,10 @@ maxcombo <- function(
 
     res$rho[i] <- rg_unique$rho[i]
     res$gamma[i] <- rg_unique$gamma[i]
-    res$estimation[i] <- weighted_o_minus_e_total
+    res$estimate[i] <- weighted_o_minus_e_total
     res$se[i] <- sqrt(weighted_var_total)
     res$var[i] <- weighted_var_total
-    res$z[i] <- res$estimation[i] / res$se[i]
+    res$z[i] <- res$estimate[i] / res$se[i]
   }
 
   # Merge back to full set of pairs ----
@@ -128,7 +128,7 @@ maxcombo <- function(
   temp <- data.frame(rho = rho, gamma = gamma)
   temp$x <- paste0("FH(", temp$rho, ", ", temp$gamma, ")")
   ans$parameter <- paste(temp$x, collapse = " + ")
-  ans$estimation <- NULL
+  ans$estimate <- NULL
   ans$se <- NULL
 
   # Get z statistics for input rho, gamma combinations

R/milestone.R

@@ -32,7 +32,7 @@
 #' @return A list frame containing:
 #'   - `method` - The method, always `"milestone"`.
 #'   - `parameter` - Milestone time point.
-#'   - `estimation` - Survival difference between the experimental and control arm.
+#'   - `estimate` - Survival difference between the experimental and control arm.
 #'   - `se` - Standard error of the control and experimental arm.
 #'   - `z` - Test statistics.
 #'
@@ -53,7 +53,6 @@
 #' cut_data |>
 #'   milestone(10, test_type = "naive")
 milestone <- function(data, ms_time, test_type = c("log-log", "naive")) {
-
   test_type <- match.arg(test_type)
 
   # Fit into KM curves
@@ -83,7 +82,7 @@ milestone <- function(data, ms_time, test_type = c("log-log", "naive")) {
   if (na_ctrl + na_exp == 2) {
     z <- -Inf
   } else {
-    if (test_type == "naive"){
+    if (test_type == "naive") {
       z_numerator <- surv_diff
       z_denominator <- surv_exp * sqrt(sigma2_exp) + surv_ctrl * sqrt(sigma2_ctrl)
     } else if (test_type == "log-log") {
@@ -92,7 +91,7 @@ milestone <- function(data, ms_time, test_type = c("log-log", "naive")) {
     }
   }
 
-  ans$estimation <- z_numerator
+  ans$estimate <- z_numerator
   ans$se <- z_denominator
   ans$z <- z_numerator / z_denominator
 

R/rmst.R

@@ -104,7 +104,7 @@ rmst <- function(
   ans <- list()
   ans$method <- "RMST"
   ans$parameter <- tau
-  ans$estimation <- res$rmst_diff$rmst_diff
+  ans$estimate <- res$rmst_diff$rmst_diff
   ans$se <- res$rmst_diff$std
   ans$z <- res$rmst_diff$rmst_diff / res$rmst_diff$std
 

R/sim_fixed_n.R

@@ -377,7 +377,7 @@ sim_fixed_n <- function(
   return(results)
 }
 
-# Build a function to calculate test related statistics (e.g., z, estimation, se, etc.) and log-hr
+# Build a function to calculate test related statistics (e.g., z, estimate, se, etc.) and log-hr
 doAnalysis <- function(d, rho_gamma, n_stratum) {
   if (nrow(rho_gamma) == 1) {
     res <- d |>
@@ -386,7 +386,7 @@ doAnalysis <- function(d, rho_gamma, n_stratum) {
     ans <- data.frame(
       method = res$method,
       parameter = res$parameter,
-      estimation = res$estimation,
+      estimate = res$estimate,
       se = res$se,
       z = res$z
     )
@@ -397,7 +397,7 @@ doAnalysis <- function(d, rho_gamma, n_stratum) {
     ans <- data.frame(
       method = rep(res$method, nrow(rho_gamma)),
       parameter = rep(res$parameter, nrow(rho_gamma)),
-      estimation = rep("-", nrow(rho_gamma)),
+      estimate = rep("-", nrow(rho_gamma)),
       se = rep("-", nrow(rho_gamma)),
       z = res$z,
       p_value = rep(res$p_value, nrow(rho_gamma))

R/sim_gs_n.R

@@ -22,15 +22,15 @@
 #' arguments will change as we add additional features.
 #'
 #' @inheritParams sim_fixed_n
-#' @param test One or more test functions such as [wlr()], [maxcombo()], or
-#'   [rmst()]. If a single test function is provided, it will be applied at each
-#'   cut. Alternatively a list of functions created by [create_test()]. The list
-#'   form is experimental and currently limited. It only accepts one test per
-#'   cutting (in the future multiple tests may be accepted), and all the tests
-#'   must consistently return the same exact results (again this may be more
-#'   flexible in the future). Importantly, note that the simulated data set is
-#'   always passed as the first positional argument to each test function
-#'   provided.
+#' @param test One or more test functions such as [wlr()], [rmst()], or
+#'   [milestone()] ([maxcombo()] can only be applied by itself). If a single
+#'   test function is provided, it will be applied at each cut. Alternatively a
+#'   list of functions created by [create_test()]. The list form is experimental
+#'   and currently limited. It only accepts one test per cutting (in the future
+#'   multiple tests may be accepted), and all the tests must consistently return
+#'   the same exact results (again this may be more flexible in the future).
+#'   Importantly, note that the simulated data set is always passed as the first
+#'   positional argument to each test function provided.
 #' @param cut A list of cutting functions created by [create_cut()], see
 #'   examples.
 #' @param seed Random seed.
@@ -81,7 +81,7 @@
 #' # - At least 20 months have elapsed after enrolling 200/400 subjects, with a
 #' #   minimum of 20 months follow-up.
 #' # However, if events accumulation is slow, we will wait for a maximum of 24 months.
-#' ia1 <- create_cut(
+#' ia1_cut <- create_cut(
 #'   planned_calendar_time = 20,
 #'   target_event_overall = 100,
 #'   max_extension_for_target_event = 24,
@@ -95,7 +95,7 @@
 #' # - At least 250 events have occurred.
 #' # - At least 10 months after IA1.
 #' # However, if events accumulation is slow, we will wait for a maximum of 34 months.
-#' ia2 <- create_cut(
+#' ia2_cut <- create_cut(
 #'   planned_calendar_time = 32,
 #'   target_event_overall = 200,
 #'   max_extension_for_target_event = 34,
@@ -106,133 +106,125 @@
 #' # The final analysis will occur at the later of the following 2 conditions:
 #' # - At least 45 months have passed since the start of the study.
 #' # - At least 300 events have occurred.
-#' fa <- create_cut(
+#' fa_cut <- create_cut(
 #'   planned_calendar_time = 45,
 #'   target_event_overall = 350
 #' )
 #'
-#' # Test 1: regular logrank test
+#' # Example 1: regular logrank test at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = wlr,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   weight = fh(rho = 0, gamma = 0)
 #' )
 #'
-#' # Test 2: weighted logrank test by FH(0, 0.5)
+#' # Example 2: weighted logrank test by FH(0, 0.5) at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = wlr,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   weight = fh(rho = 0, gamma = 0.5)
 #' )
 #'
-#' # Test 3: weighted logrank test by MB(3)
+#' # Example 3: weighted logrank test by MB(3) at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = wlr,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   weight = mb(delay = 3)
 #' )
 #'
-#' # Test 4: weighted logrank test by early zero (6)
+#' # Example 4: weighted logrank test by early zero (6) at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = wlr,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   weight = early_zero(6)
 #' )
 #'
-#' # Test 5: RMST
+#' # Example 5: RMST at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = rmst,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   tau = 20
 #' )
 #'
-#' # Test 6: Milestone
+#' # Example 6: Milestone at all 3 analyses
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
 #'   test = milestone,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
 #'   seed = 2024,
 #'   ms_time = 10
 #' )
 #'
-#' # Test 7: MaxCombo (WLR-FH(0,0) + WLR-FH(0, 0.5))
-#' # for all analyses
+#' # Example 7: WLR with fh(0, 0.5) test at IA1,
+#' # WLR with mb(6, Inf) at IA2, and milestone test at FA
+#' ia1_test <- create_test(wlr, weight = fh(rho = 0, gamma = 0.5))
+#' ia2_test <- create_test(wlr, weight = mb(delay = 6, w_max = Inf))
+#' fa_test <- create_test(milestone, ms_time = 10)
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
-#'   test = maxcombo,
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
-#'   seed = 2024,
-#'   rho = c(0, 0),
-#'   gamma = c(0, 0.5)
+#'   test = list(ia1 = ia1_test, ia2 = ia2_test, fa = fa_test),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
+#'   seed = 2024
 #' )
 #'
-#' # Test 8: MaxCombo (WLR-FH(0,0.5) + milestone(10))
-#' # for all analyses
-#' \dontrun{
-#' sim_gs_n(
-#'   n_sim = 3,
-#'   sample_size = 400,
-#'   enroll_rate = enroll_rate,
-#'   fail_rate = fail_rate,
-#'   test = maxcombo(test1 = wlr, test2 = milestone),
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
-#'   seed = 2024,
-#'   test1_par = list(weight = fh(rho = 0, gamma = 0.5)),
-#'   test2_par = list(ms_time = 10)
+#' # WARNING: Multiple tests per cut will be enabled in a future version.
+#' #          Currently does not work.
+#' # Example 8: At IA1, we conduct 3 tests, LR, WLR with fh(0, 0.5), and RMST test.
+#' # At IA2, we conduct 2 tests, LR and WLR with early zero (6).
+#' # At FA, we conduct 2 tests, LR and milestone test.
+#' ia1_test <- list(
+#'   test1 = create_test(wlr, weight = fh(rho = 0, gamma = 0)),
+#'   test2 = create_test(wlr, weight = fh(rho = 0, gamma = 0.5)),
+#'   test3 = create_test(rmst, tau = 20)
+#' )
+#' ia2_test <- list(
+#'   test1 = create_test(wlr, weight = fh(rho = 0, gamma = 0)),
+#'   test2 = create_test(wlr, weight = early_zero(6))
+#' )
+#' fa_test <- list(
+#'   test1 = create_test(wlr, weight = fh(rho = 0, gamma = 0)),
+#'   test3 = create_test(milestone, ms_time = 20)
 #' )
-#' }
-#'
-#' # Test 9: MaxCombo (WLR-FH(0,0) at IAs
-#' # and WLR-FH(0,0) + milestone(10) + WLR-MB(4,2) at FA)
 #' \dontrun{
 #' sim_gs_n(
 #'   n_sim = 3,
 #'   sample_size = 400,
 #'   enroll_rate = enroll_rate,
 #'   fail_rate = fail_rate,
-#'   test = list(ia1 = wlr, ia2 = wlr, fa = maxcombo),
-#'   cut = list(ia1 = ia1, ia2 = ia2, fa = fa),
-#'   seed = 2024,
-#'   test_par = list(
-#'     ia1 = list(weight = fh(rho = 0, gamma = 0)),
-#'     ia2 = list(weight = fh(rho = 0, gamma = 0)),
-#'     ia3 = list(
-#'       test1_par = list(weight = fh(rho = 0, gamma = 0)),
-#'       test2_par = list(ms_time = 10),
-#'       test3_par = list(delay = 4, w_max = 2)
-#'     )
-#'   )
+#'   test = list(ia1 = ia1_test, ia2 = ia2_test, fa = fa_test),
+#'   cut = list(ia1 = ia1_cut, ia2 = ia2_cut, fa = fa_cut),
+#'   seed = 2024
 #' )
 #' }
 sim_gs_n <- function(
@@ -294,11 +286,18 @@ sim_gs_n <- function(
 
       # Test
       ans_1sim_new <- test[[i_analysis]](simu_data_cut, ...)
-      ans_1sim_new$analysis <- i_analysis
-      ans_1sim_new$cut_date <- cut_date[i_analysis]
-      ans_1sim_new$sim_id <- sim_id
-      ans_1sim_new$n <- nrow(simu_data_cut)
-      ans_1sim_new$event <- sum(simu_data_cut$event)
+      ans_1sim_new <- c(sim_id = sim_id, ans_1sim_new)
+      ans_1sim_new <- append(
+        x = ans_1sim_new,
+        values = c(
+          analysis = i_analysis,
+          cut_date = cut_date[i_analysis],
+          n = nrow(simu_data_cut),
+          event = sum(simu_data_cut$event)
+        ),
+        after = 3
+      )
+      ans_1sim_new <- convert_list_to_df_w_list_cols(ans_1sim_new)
 
       # rbind simulation results for all IA(s) and FA in 1 simulation
       ans_1sim <- rbind(ans_1sim, ans_1sim_new)
@@ -422,3 +421,25 @@ multitest <- function(data, ...) {
   }
   return(output)
 }
+
+# Convert a list to a one row data frame using list columns
+convert_list_to_df_w_list_cols <- function(x) {
+  stopifnot(is.list(x), !is.data.frame(x))
+
+  new_list <- vector(mode = "list", length = length(x))
+  names(new_list) <- names(x)
+
+  for (i in seq_along(x)) {
+    if (length(x[[i]]) > 1) {
+      new_list[[i]] <- I(list(x[[i]]))
+    } else {
+      new_list[i] <- x[i]
+    }
+  }
+
+  # Convert the list to a data frame with one row
+  df_w_list_cols <- do.call(data.frame, new_list)
+  stopifnot(nrow(df_w_list_cols) == 1)
+
+  return(df_w_list_cols)
+}

R/wlr.R

@@ -27,7 +27,7 @@
 #'
 #' @return A list containing the test method (`method`),
 #' parameters of this test method (`parameter`),
-#' point estimation of the treatment effect (`estimation`),
+#' point estimate of the treatment effect (`estimate`),
 #' standardized error of the treatment effect (`se`),
 #' Z-score (`z`), p-values (`p_value`).
 #'
@@ -95,9 +95,9 @@ wlr <- function(data, weight, return_variance = FALSE) {
     x <- x |> fh_weight(rho = weight$rho, gamma = weight$gamma)
 
     ans$parameter <- paste0("FH(rho=", weight$rho, ", gamma=", weight$gamma, ")")
-    ans$estimation <- sum(x$weight * x$o_minus_e)
+    ans$estimate <- sum(x$weight * x$o_minus_e)
     ans$se <- sqrt(sum(x$weight^2 * x$var_o_minus_e))
-    ans$z <- ans$estimation / ans$se
+    ans$z <- ans$estimate / ans$se
   } else if (inherits(weight, "mb")) {
     x <- x |> mb_weight(delay = weight$delay, w_max = weight$w_max)