From 99493361b7ecc7895d7239fe08fc64ff1d254635 Mon Sep 17 00:00:00 2001 From: Nan Xiao Date: Sun, 26 Nov 2023 23:55:18 -0500 Subject: [PATCH 1/2] Use base pipe operator --- DESCRIPTION | 7 ++- NAMESPACE | 2 - R/MBdelayed.R | 20 ++++----- R/counting_process.R | 4 +- R/cut_data_by_date.R | 2 +- R/cut_data_by_event.R | 4 +- R/early_zero_weight.R | 22 +++++----- R/fh_weight.R | 24 +++++----- R/get_cut_date_by_event.R | 2 +- R/mb_weight.R | 12 ++--- R/pvalue_maxcombo.R | 6 +-- R/randomize_by_fixed_block.R | 6 +-- R/sim_fixed_n.R | 20 ++++----- R/simfix2simpwsurv.R | 10 ++--- R/utils-pipe.R | 14 ------ man/MBdelayed.Rd | 20 ++++----- man/counting_process.Rd | 4 +- man/cut_data_by_date.Rd | 2 +- man/cut_data_by_event.Rd | 2 +- man/early_zero_weight.Rd | 22 +++++----- man/fh_weight.Rd | 22 +++++----- man/get_cut_date_by_event.Rd | 2 +- man/mb_weight.Rd | 12 ++--- man/pipe.Rd | 20 --------- man/pvalue_maxcombo.Rd | 6 +-- man/randomize_by_fixed_block.Rd | 6 +-- man/sim_fixed_n.Rd | 20 ++++----- man/simfix2simpwsurv.Rd | 10 ++--- man/simtrial-package.Rd | 2 +- .../test-double_programming_sim_pw_surv.R | 4 +- tests/testthat/test-double_programming_wMB.R | 12 ++--- .../test-independent_test_counting_process.R | 24 +++++----- ...t-independent_test_get_cut_date_by_event.R | 8 ++-- .../test-independent_test_pvalue_maxcombo.R | 4 +- tests/testthat/test-independent_test_wlr.R | 8 ++-- vignettes/arbitraryhazard.Rmd | 12 ++--- vignettes/modestWLRTVignette.Rmd | 44 +++++++++---------- vignettes/parallel.Rmd | 8 ++-- vignettes/pvalue_maxcomboVignette.Rmd | 38 ++++++++-------- vignettes/simtrialroutines.Rmd | 16 +++---- 40 files changed, 219 insertions(+), 264 deletions(-) delete mode 100644 R/utils-pipe.R delete mode 100644 man/pipe.Rd diff --git a/DESCRIPTION b/DESCRIPTION index 95de5327..a9841880 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: simtrial Type: Package Title: Clinical Trial Simulation -Version: 0.3.0.6 +Version: 0.3.0.7 Authors@R: c( person("Keaven", "Anderson", email = "keaven_anderson@merck.com", role = c("aut")), person("Yilong", "Zhang", email = "elong0527@gmail.com", role = c("aut")), @@ -20,7 +20,7 @@ Authors@R: c( person("John", "Blischak", role = c("ctb")), person("Merck & Co., Inc., Rahway, NJ, USA and its affiliates", role = "cph") ) -Description: simtrial provides some basic routines for simulating a +Description: Provides some basic routines for simulating a clinical trial. The primary intent is to provide some tools to generate trial simulations for trials with time to event outcomes. Piecewise exponential failure rates and piecewise constant @@ -34,14 +34,13 @@ BugReports: https://github.com/Merck/simtrial/issues Encoding: UTF-8 LazyData: true VignetteBuilder: knitr -Depends: R (>= 3.5.0) +Depends: R (>= 4.1.0) Imports: Rcpp, data.table, doFuture, foreach, future, - magrittr, methods, mvtnorm, stats, diff --git a/NAMESPACE b/NAMESPACE index ff79e821..ed287d0c 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,6 +1,5 @@ # Generated by roxygen2: do not edit by hand -export("%>%") export(counting_process) export(cut_data_by_date) export(cut_data_by_event) @@ -34,7 +33,6 @@ importFrom(data.table,setorderv) importFrom(data.table,uniqueN) importFrom(doFuture,"%dofuture%") importFrom(future,plan) -importFrom(magrittr,"%>%") importFrom(methods,is) importFrom(mvtnorm,GenzBretz) importFrom(mvtnorm,pmvnorm) diff --git a/R/MBdelayed.R b/R/MBdelayed.R index 0a44a48d..46469e0d 100644 --- a/R/MBdelayed.R +++ b/R/MBdelayed.R @@ -48,36 +48,36 @@ #' #' # Set up time, event, number of event dataset for testing #' # with arbitrary weights -#' ten <- MBdelayed %>% counting_process(arm = "experimental") +#' ten <- MBdelayed |> counting_process(arm = "experimental") #' head(ten) #' #' # MaxCombo with logrank, FH(0,1), FH(1,1) -#' ten %>% -#' fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) %>% +#' ten |> +#' fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) |> #' pvalue_maxcombo() #' #' # Magirr-Burman modestly down-weighted rank test with 6 month delay #' # First, add weights -#' ten <- ten %>% mb_weight(6) +#' ten <- ten |> mb_weight(6) #' head(ten) #' #' # Now compute test based on these weights -#' ten %>% +#' ten |> #' summarize( #' S = sum(o_minus_e * mb_weight), #' V = sum(var_o_minus_e * mb_weight^2), #' Z = S / sqrt(V) -#' ) %>% +#' ) |> #' mutate(p = pnorm(Z)) #' #' # Create 0 weights for first 6 months -#' ten <- ten %>% mutate(w6 = 1 * (tte >= 6)) -#' ten %>% +#' ten <- ten |> mutate(w6 = 1 * (tte >= 6)) +#' ten |> #' summarize( #' S = sum(o_minus_e * w6), #' V = sum(var_o_minus_e * w6^2), #' Z = S / sqrt(V) -#' ) %>% +#' ) |> #' mutate(p = pnorm(Z)) #' #' # Generate another dataset @@ -100,5 +100,5 @@ #' ) #' ) #' # Cut data at 24 months after final enrollment -#' MBdelayed2 <- ds %>% cut_data_by_date(max(ds$enroll_time) + 24) +#' MBdelayed2 <- ds |> cut_data_by_date(max(ds$enroll_time) + 24) "MBdelayed" diff --git a/R/counting_process.R b/R/counting_process.R index 4ff023cc..ccfeb492 100644 --- a/R/counting_process.R +++ b/R/counting_process.R @@ -62,8 +62,6 @@ #' @export #' #' @examples -#' library(magrittr) -#' #' # Example 1 #' x <- data.frame( #' stratum = c(rep(1, 10), rep(2, 6)), @@ -75,7 +73,7 @@ #' #' # Example 2 #' x <- sim_pw_surv(n = 400) -#' y <- cut_data_by_event(x, 150) %>% counting_process(arm = "experimental") +#' y <- cut_data_by_event(x, 150) |> counting_process(arm = "experimental") #' # Weighted logrank test (Z-value and 1-sided p-value) #' z <- sum(y$o_minus_e) / sqrt(sum(y$var_o_minus_e)) #' c(z, pnorm(z)) diff --git a/R/cut_data_by_date.R b/R/cut_data_by_date.R index 052512b7..d8224c76 100644 --- a/R/cut_data_by_date.R +++ b/R/cut_data_by_date.R @@ -31,7 +31,7 @@ #' @examples #' # Use default enrollment and event rates and #' # cut at calendar time 5 after start of randomization -#' sim_pw_surv(n = 20) %>% cut_data_by_date(5) +#' sim_pw_surv(n = 20) |> cut_data_by_date(5) cut_data_by_date <- function(x, cut_date) { ans <- as.data.table(x) ans <- ans[enroll_time <= cut_date, ] diff --git a/R/cut_data_by_event.R b/R/cut_data_by_event.R index e4ae8f6c..c17c6d15 100644 --- a/R/cut_data_by_event.R +++ b/R/cut_data_by_event.R @@ -31,10 +31,10 @@ #' #' @examples #' # Use default enrollment and event rates at cut at 100 events -#' x <- sim_pw_surv(n = 200) %>% cut_data_by_event(100) +#' x <- sim_pw_surv(n = 200) |> cut_data_by_event(100) #' table(x$event, x$treatment) cut_data_by_event <- function(x, event) { cut_date <- get_cut_date_by_event(x, event) - ans <- x %>% cut_data_by_date(cut_date = cut_date) + ans <- x |> cut_data_by_date(cut_date = cut_date) return(ans) } diff --git a/R/early_zero_weight.R b/R/early_zero_weight.R index 8302cbb4..5117271f 100644 --- a/R/early_zero_weight.R +++ b/R/early_zero_weight.R @@ -40,23 +40,23 @@ #' library(gsDesign2) #' #' # Example 1: Unstratified -#' sim_pw_surv(n = 200) %>% -#' cut_data_by_event(125) %>% -#' counting_process(arm = "experimental") %>% -#' early_zero_weight(early_period = 2) %>% +#' sim_pw_surv(n = 200) |> +#' cut_data_by_event(125) |> +#' counting_process(arm = "experimental") |> +#' early_zero_weight(early_period = 2) |> #' filter(row_number() %in% seq(5, 200, 40)) #' #' # Example 2: Stratified #' n <- 500 #' # Two strata #' stratum <- c("Biomarker-positive", "Biomarker-negative") -#' prevelance_ratio <- c(0.6, 0.4) +#' prevalence_ratio <- c(0.6, 0.4) #' #' # Enrollment rate #' enroll_rate <- define_enroll_rate( #' stratum = rep(stratum, each = 2), #' duration = c(2, 10, 2, 10), -#' rate = c(c(1, 4) * prevelance_ratio[1], c(1, 4) * prevelance_ratio[2]) +#' rate = c(c(1, 4) * prevalence_ratio[1], c(1, 4) * prevalence_ratio[2]) #' ) #' enroll_rate$rate <- enroll_rate$rate * n / sum(enroll_rate$duration * enroll_rate$rate) #' @@ -80,16 +80,16 @@ #' sim_pw_surv( #' n = n, # Sample size #' # Stratified design with prevalence ratio of 6:4 -#' stratum = tibble(stratum = stratum, p = prevelance_ratio), +#' stratum = tibble(stratum = stratum, p = prevalence_ratio), #' # Randomization ratio #' block = c("control", "control", "experimental", "experimental"), #' enroll_rate = enroll_rate, # Enrollment rate #' fail_rate = temp$fail_rate, # Failure rate #' dropout_rate = temp$dropout_rate # Dropout rate -#' ) %>% -#' cut_data_by_event(125) %>% -#' counting_process(arm = "experimental") %>% -#' early_zero_weight(early_period = 2, fail_rate = fail_rate) %>% +#' ) |> +#' cut_data_by_event(125) |> +#' counting_process(arm = "experimental") |> +#' early_zero_weight(early_period = 2, fail_rate = fail_rate) |> #' filter(row_number() %in% seq(5, 200, 40)) early_zero_weight <- function(x, early_period = 4, fail_rate = NULL) { ans <- as.data.table(x) diff --git a/R/fh_weight.R b/R/fh_weight.R index d46f29b7..8cf26fe3 100644 --- a/R/fh_weight.R +++ b/R/fh_weight.R @@ -88,8 +88,8 @@ #' #' # Example 1 #' # Use default enrollment and event rates at cut at 100 events -#' x <- sim_pw_surv(n = 200) %>% -#' cut_data_by_event(100) %>% +#' x <- sim_pw_surv(n = 200) |> +#' cut_data_by_event(100) |> #' counting_process(arm = "experimental") #' #' # Compute logrank FH(0, 1) @@ -104,9 +104,9 @@ #' # Example 2 #' # Use default enrollment and event rates at cut of 100 events #' set.seed(123) -#' x <- sim_pw_surv(n = 200) %>% -#' cut_data_by_event(100) %>% -#' counting_process(arm = "experimental") %>% +#' x <- sim_pw_surv(n = 200) |> +#' cut_data_by_event(100) |> +#' counting_process(arm = "experimental") |> #' fh_weight(rho_gamma = data.frame(rho = c(0, 0), gamma = c(0, 1)), return_corr = TRUE) #' #' # Compute p-value for MaxCombo @@ -118,11 +118,11 @@ #' )[1] #' #' # Check that covariance is as expected -#' x <- sim_pw_surv(n = 200) %>% -#' cut_data_by_event(100) %>% +#' x <- sim_pw_surv(n = 200) |> +#' cut_data_by_event(100) |> #' counting_process(arm = "experimental") #' -#' x %>% fh_weight( +#' x |> fh_weight( #' rho_gamma = data.frame( #' rho = c(0, 0), #' gamma = c(0, 1) @@ -131,7 +131,7 @@ #' ) #' #' # Off-diagonal element should be variance in following -#' x %>% fh_weight( +#' x |> fh_weight( #' rho_gamma = data.frame( #' rho = 0, #' gamma = .5 @@ -140,10 +140,10 @@ #' ) #' #' # Compare off diagonal result with fh_weight() -#' x %>% fh_weight(rho_gamma = data.frame(rho = 0, gamma = .5)) +#' x |> fh_weight(rho_gamma = data.frame(rho = 0, gamma = .5)) fh_weight <- function( - x = sim_pw_surv(n = 200) %>% - cut_data_by_event(150) %>% + x = sim_pw_surv(n = 200) |> + cut_data_by_event(150) |> counting_process(arm = "experimental"), rho_gamma = data.frame( rho = c(0, 0, 1, 1), diff --git a/R/get_cut_date_by_event.R b/R/get_cut_date_by_event.R index de9af1c2..1c6c311d 100644 --- a/R/get_cut_date_by_event.R +++ b/R/get_cut_date_by_event.R @@ -56,7 +56,7 @@ #' ) #' ) #' -#' d <- get_cut_date_by_event(x %>% filter(stratum == "Positive"), event = 50) +#' d <- get_cut_date_by_event(x |> filter(stratum == "Positive"), event = 50) #' #' y <- cut_data_by_date(x, cut_date = d) #' table(y$stratum, y$event) diff --git a/R/mb_weight.R b/R/mb_weight.R index 61f4d2fe..387e53f7 100644 --- a/R/mb_weight.R +++ b/R/mb_weight.R @@ -75,14 +75,14 @@ #' #' # Use default enrollment and event rates at cut at 100 events #' # For transparency, may be good to set either `delay` or `w_max` to `Inf` -#' x <- sim_pw_surv(n = 200) %>% -#' cut_data_by_event(125) %>% +#' x <- sim_pw_surv(n = 200) |> +#' cut_data_by_event(125) |> #' counting_process(arm = "experimental") #' #' # Example 1 #' # Compute Magirr-Burman weights with `delay = 6` -#' ZMB <- x %>% -#' mb_weight(delay = 6, w_max = Inf) %>% +#' ZMB <- x |> +#' mb_weight(delay = 6, w_max = Inf) |> #' summarize( #' S = sum(o_minus_e * mb_weight), #' V = sum(var_o_minus_e * mb_weight^2), @@ -94,8 +94,8 @@ #' #' # Example 2 #' # Now compute with maximum weight of 2 as recommended in Magirr, 2021 -#' ZMB2 <- x %>% -#' mb_weight(delay = Inf, w_max = 2) %>% +#' ZMB2 <- x |> +#' mb_weight(delay = Inf, w_max = 2) |> #' summarize( #' S = sum(o_minus_e * mb_weight), #' V = sum(var_o_minus_e * mb_weight^2), diff --git a/R/pvalue_maxcombo.R b/R/pvalue_maxcombo.R index 19285bca..71308428 100644 --- a/R/pvalue_maxcombo.R +++ b/R/pvalue_maxcombo.R @@ -62,9 +62,9 @@ #' head(xx) #' #' # MaxCombo power estimate for cutoff at max of targeted events, minimum follow-up -#' p <- xx %>% -#' group_by(sim) %>% -#' group_map(~ pvalue_maxcombo(.x)) %>% +#' p <- xx |> +#' group_by(sim) |> +#' group_map(~ pvalue_maxcombo(.x)) |> #' unlist() #' mean(p < .025) pvalue_maxcombo <- function( diff --git a/R/randomize_by_fixed_block.R b/R/randomize_by_fixed_block.R index 4885bba0..5d0b0b67 100644 --- a/R/randomize_by_fixed_block.R +++ b/R/randomize_by_fixed_block.R @@ -37,12 +37,12 @@ #' #' # Example 1 #' # 2:1 randomization with block size 3, treatments "A" and "B" -#' data.frame(x = 1:10) %>% mutate(Treatment = randomize_by_fixed_block(block = c("A", "B", "B"))) +#' data.frame(x = 1:10) |> mutate(Treatment = randomize_by_fixed_block(block = c("A", "B", "B"))) #' #' # Example 2 #' # Stratified randomization -#' data.frame(stratum = c(rep("A", 10), rep("B", 10))) %>% -#' group_by(stratum) %>% +#' data.frame(stratum = c(rep("A", 10), rep("B", 10))) |> +#' group_by(stratum) |> #' mutate(Treatment = randomize_by_fixed_block()) randomize_by_fixed_block <- function(n = 10, block = c(0, 0, 1, 1)) { length_block <- length(block) diff --git a/R/sim_fixed_n.R b/R/sim_fixed_n.R index 7b27bd0a..867e705d 100644 --- a/R/sim_fixed_n.R +++ b/R/sim_fixed_n.R @@ -94,24 +94,24 @@ #' rho_gamma = data.frame(rho = 0, gamma = c(0, 1)) #' ) #' # Get power approximation for FH, data cutoff combination -#' xx %>% -#' group_by(cut, rho, gamma) %>% +#' xx |> +#' group_by(cut, rho, gamma) |> #' summarize(mean(z <= qnorm(.025))) #' #' # MaxCombo power estimate for cutoff at max of targeted events, minimum follow-up -#' p <- xx %>% -#' filter(cut != "Targeted events") %>% -#' group_by(sim) %>% -#' group_map(~ pvalue_maxcombo(.x)) %>% +#' p <- xx |> +#' filter(cut != "Targeted events") |> +#' group_by(sim) |> +#' group_map(~ pvalue_maxcombo(.x)) |> #' unlist() #' #' mean(p < .025) #' #' # MaxCombo estimate for targeted events cutoff -#' p <- xx %>% -#' filter(cut == "Targeted events") %>% -#' group_by(sim) %>% -#' group_map(~ pvalue_maxcombo(.x)) %>% +#' p <- xx |> +#' filter(cut == "Targeted events") |> +#' group_by(sim) |> +#' group_map(~ pvalue_maxcombo(.x)) |> #' unlist() #' #' mean(p < .025) diff --git a/R/simfix2simpwsurv.R b/R/simfix2simpwsurv.R index bef075b8..d8088747 100644 --- a/R/simfix2simpwsurv.R +++ b/R/simfix2simpwsurv.R @@ -41,8 +41,6 @@ #' @export #' #' @examples -#' library(magrittr) -#' #' # Example 1 #' # Convert standard input #' simfix2simpwsurv() @@ -75,10 +73,10 @@ #' ) #' #' # Cut after 200 events and do a stratified logrank test -#' dat <- sim %>% -#' cut_data_by_event(200) %>% # cut data -#' counting_process(arm = "experimental") %>% # convert format for tenFH -#' fh_weight(rho_gamma = data.frame(rho = 0, gamma = 0)) # stratified logrank +#' dat <- sim |> +#' cut_data_by_event(200) |> # Cut data +#' counting_process(arm = "experimental") |> # Convert format for fh_weight() +#' fh_weight(rho_gamma = data.frame(rho = 0, gamma = 0)) # Stratified logrank simfix2simpwsurv <- function( # Failure rates as in sim_fixed_n() fail_rate = data.frame( diff --git a/R/utils-pipe.R b/R/utils-pipe.R deleted file mode 100644 index fd0b1d13..00000000 --- a/R/utils-pipe.R +++ /dev/null @@ -1,14 +0,0 @@ -#' Pipe operator -#' -#' See \code{magrittr::\link[magrittr:pipe]{\%>\%}} for details. -#' -#' @name %>% -#' @rdname pipe -#' @keywords internal -#' @export -#' @importFrom magrittr %>% -#' @usage lhs \%>\% rhs -#' @param lhs A value or the magrittr placeholder. -#' @param rhs A function call using the magrittr semantics. -#' @return The result of calling `rhs(lhs)`. -NULL diff --git a/man/MBdelayed.Rd b/man/MBdelayed.Rd index 4e667ea6..7dd631bd 100644 --- a/man/MBdelayed.Rd +++ b/man/MBdelayed.Rd @@ -35,36 +35,36 @@ legend("topright", legend = c("control", "experimental"), lty = 1:2) # Set up time, event, number of event dataset for testing # with arbitrary weights -ten <- MBdelayed \%>\% counting_process(arm = "experimental") +ten <- MBdelayed |> counting_process(arm = "experimental") head(ten) # MaxCombo with logrank, FH(0,1), FH(1,1) -ten \%>\% - fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) \%>\% +ten |> + fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) |> pvalue_maxcombo() # Magirr-Burman modestly down-weighted rank test with 6 month delay # First, add weights -ten <- ten \%>\% mb_weight(6) +ten <- ten |> mb_weight(6) head(ten) # Now compute test based on these weights -ten \%>\% +ten |> summarize( S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), Z = S / sqrt(V) - ) \%>\% + ) |> mutate(p = pnorm(Z)) # Create 0 weights for first 6 months -ten <- ten \%>\% mutate(w6 = 1 * (tte >= 6)) -ten \%>\% +ten <- ten |> mutate(w6 = 1 * (tte >= 6)) +ten |> summarize( S = sum(o_minus_e * w6), V = sum(var_o_minus_e * w6^2), Z = S / sqrt(V) - ) \%>\% + ) |> mutate(p = pnorm(Z)) # Generate another dataset @@ -87,7 +87,7 @@ ds <- sim_pw_surv( ) ) # Cut data at 24 months after final enrollment -MBdelayed2 <- ds \%>\% cut_data_by_date(max(ds$enroll_time) + 24) +MBdelayed2 <- ds |> cut_data_by_date(max(ds$enroll_time) + 24) } \references{ Magirr, Dominic, and Carlā€Fredrik Burman. 2019. diff --git a/man/counting_process.Rd b/man/counting_process.Rd index 7c7ac7ef..6056bdec 100644 --- a/man/counting_process.Rd +++ b/man/counting_process.Rd @@ -53,8 +53,6 @@ The function only considered two group situation. The tie is handled by the Breslow's Method. } \examples{ -library(magrittr) - # Example 1 x <- data.frame( stratum = c(rep(1, 10), rep(2, 6)), @@ -66,7 +64,7 @@ counting_process(x, arm = 1) # Example 2 x <- sim_pw_surv(n = 400) -y <- cut_data_by_event(x, 150) \%>\% counting_process(arm = "experimental") +y <- cut_data_by_event(x, 150) |> counting_process(arm = "experimental") # Weighted logrank test (Z-value and 1-sided p-value) z <- sum(y$o_minus_e) / sqrt(sum(y$var_o_minus_e)) c(z, pnorm(z)) diff --git a/man/cut_data_by_date.Rd b/man/cut_data_by_date.Rd index c86c3795..f5c15d3a 100644 --- a/man/cut_data_by_date.Rd +++ b/man/cut_data_by_date.Rd @@ -21,5 +21,5 @@ Cut a dataset for analysis at a specified date \examples{ # Use default enrollment and event rates and # cut at calendar time 5 after start of randomization -sim_pw_surv(n = 20) \%>\% cut_data_by_date(5) +sim_pw_surv(n = 20) |> cut_data_by_date(5) } diff --git a/man/cut_data_by_event.Rd b/man/cut_data_by_event.Rd index 23b4a25e..733760ba 100644 --- a/man/cut_data_by_event.Rd +++ b/man/cut_data_by_event.Rd @@ -21,6 +21,6 @@ event count is reached. } \examples{ # Use default enrollment and event rates at cut at 100 events -x <- sim_pw_surv(n = 200) \%>\% cut_data_by_event(100) +x <- sim_pw_surv(n = 200) |> cut_data_by_event(100) table(x$event, x$treatment) } diff --git a/man/early_zero_weight.Rd b/man/early_zero_weight.Rd index afbb77d2..db22d012 100644 --- a/man/early_zero_weight.Rd +++ b/man/early_zero_weight.Rd @@ -26,23 +26,23 @@ library(dplyr) library(gsDesign2) # Example 1: Unstratified -sim_pw_surv(n = 200) \%>\% - cut_data_by_event(125) \%>\% - counting_process(arm = "experimental") \%>\% - early_zero_weight(early_period = 2) \%>\% +sim_pw_surv(n = 200) |> + cut_data_by_event(125) |> + counting_process(arm = "experimental") |> + early_zero_weight(early_period = 2) |> filter(row_number() \%in\% seq(5, 200, 40)) # Example 2: Stratified n <- 500 # Two strata stratum <- c("Biomarker-positive", "Biomarker-negative") -prevelance_ratio <- c(0.6, 0.4) +prevalence_ratio <- c(0.6, 0.4) # Enrollment rate enroll_rate <- define_enroll_rate( stratum = rep(stratum, each = 2), duration = c(2, 10, 2, 10), - rate = c(c(1, 4) * prevelance_ratio[1], c(1, 4) * prevelance_ratio[2]) + rate = c(c(1, 4) * prevalence_ratio[1], c(1, 4) * prevalence_ratio[2]) ) enroll_rate$rate <- enroll_rate$rate * n / sum(enroll_rate$duration * enroll_rate$rate) @@ -66,16 +66,16 @@ set.seed(2023) sim_pw_surv( n = n, # Sample size # Stratified design with prevalence ratio of 6:4 - stratum = tibble(stratum = stratum, p = prevelance_ratio), + stratum = tibble(stratum = stratum, p = prevalence_ratio), # Randomization ratio block = c("control", "control", "experimental", "experimental"), enroll_rate = enroll_rate, # Enrollment rate fail_rate = temp$fail_rate, # Failure rate dropout_rate = temp$dropout_rate # Dropout rate -) \%>\% - cut_data_by_event(125) \%>\% - counting_process(arm = "experimental") \%>\% - early_zero_weight(early_period = 2, fail_rate = fail_rate) \%>\% +) |> + cut_data_by_event(125) |> + counting_process(arm = "experimental") |> + early_zero_weight(early_period = 2, fail_rate = fail_rate) |> filter(row_number() \%in\% seq(5, 200, 40)) } \references{ diff --git a/man/fh_weight.Rd b/man/fh_weight.Rd index 7e3cb3d1..06d7ac8c 100644 --- a/man/fh_weight.Rd +++ b/man/fh_weight.Rd @@ -5,7 +5,7 @@ \title{Fleming-Harrington weighted logrank tests} \usage{ fh_weight( - x = sim_pw_surv(n = 200) \%>\% cut_data_by_event(150) \%>\% counting_process(arm = + x = counting_process(cut_data_by_event(sim_pw_surv(n = 200), 150), arm = "experimental"), rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_variance = FALSE, @@ -86,8 +86,8 @@ library(dplyr) # Example 1 # Use default enrollment and event rates at cut at 100 events -x <- sim_pw_surv(n = 200) \%>\% - cut_data_by_event(100) \%>\% +x <- sim_pw_surv(n = 200) |> + cut_data_by_event(100) |> counting_process(arm = "experimental") # Compute logrank FH(0, 1) @@ -102,9 +102,9 @@ fh_weight(x, rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 0)), retur # Example 2 # Use default enrollment and event rates at cut of 100 events set.seed(123) -x <- sim_pw_surv(n = 200) \%>\% - cut_data_by_event(100) \%>\% - counting_process(arm = "experimental") \%>\% +x <- sim_pw_surv(n = 200) |> + cut_data_by_event(100) |> + counting_process(arm = "experimental") |> fh_weight(rho_gamma = data.frame(rho = c(0, 0), gamma = c(0, 1)), return_corr = TRUE) # Compute p-value for MaxCombo @@ -116,11 +116,11 @@ library(mvtnorm) )[1] # Check that covariance is as expected -x <- sim_pw_surv(n = 200) \%>\% - cut_data_by_event(100) \%>\% +x <- sim_pw_surv(n = 200) |> + cut_data_by_event(100) |> counting_process(arm = "experimental") -x \%>\% fh_weight( +x |> fh_weight( rho_gamma = data.frame( rho = c(0, 0), gamma = c(0, 1) @@ -129,7 +129,7 @@ x \%>\% fh_weight( ) # Off-diagonal element should be variance in following -x \%>\% fh_weight( +x |> fh_weight( rho_gamma = data.frame( rho = 0, gamma = .5 @@ -138,5 +138,5 @@ x \%>\% fh_weight( ) # Compare off diagonal result with fh_weight() -x \%>\% fh_weight(rho_gamma = data.frame(rho = 0, gamma = .5)) +x |> fh_weight(rho_gamma = data.frame(rho = 0, gamma = .5)) } diff --git a/man/get_cut_date_by_event.Rd b/man/get_cut_date_by_event.Rd index 6586f921..2f6f65c3 100644 --- a/man/get_cut_date_by_event.Rd +++ b/man/get_cut_date_by_event.Rd @@ -46,7 +46,7 @@ x <- sim_pw_surv( ) ) -d <- get_cut_date_by_event(x \%>\% filter(stratum == "Positive"), event = 50) +d <- get_cut_date_by_event(x |> filter(stratum == "Positive"), event = 50) y <- cut_data_by_date(x, cut_date = d) table(y$stratum, y$event) diff --git a/man/mb_weight.Rd b/man/mb_weight.Rd index b67d27ea..eda2c380 100644 --- a/man/mb_weight.Rd +++ b/man/mb_weight.Rd @@ -55,14 +55,14 @@ library(dplyr) # Use default enrollment and event rates at cut at 100 events # For transparency, may be good to set either `delay` or `w_max` to `Inf` -x <- sim_pw_surv(n = 200) \%>\% - cut_data_by_event(125) \%>\% +x <- sim_pw_surv(n = 200) |> + cut_data_by_event(125) |> counting_process(arm = "experimental") # Example 1 # Compute Magirr-Burman weights with `delay = 6` -ZMB <- x \%>\% - mb_weight(delay = 6, w_max = Inf) \%>\% +ZMB <- x |> + mb_weight(delay = 6, w_max = Inf) |> summarize( S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), @@ -74,8 +74,8 @@ pnorm(ZMB$z) # Example 2 # Now compute with maximum weight of 2 as recommended in Magirr, 2021 -ZMB2 <- x \%>\% - mb_weight(delay = Inf, w_max = 2) \%>\% +ZMB2 <- x |> + mb_weight(delay = Inf, w_max = 2) |> summarize( S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), diff --git a/man/pipe.Rd b/man/pipe.Rd deleted file mode 100644 index a648c296..00000000 --- a/man/pipe.Rd +++ /dev/null @@ -1,20 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/utils-pipe.R -\name{\%>\%} -\alias{\%>\%} -\title{Pipe operator} -\usage{ -lhs \%>\% rhs -} -\arguments{ -\item{lhs}{A value or the magrittr placeholder.} - -\item{rhs}{A function call using the magrittr semantics.} -} -\value{ -The result of calling \code{rhs(lhs)}. -} -\description{ -See \code{magrittr::\link[magrittr:pipe]{\%>\%}} for details. -} -\keyword{internal} diff --git a/man/pvalue_maxcombo.Rd b/man/pvalue_maxcombo.Rd index e083f9cc..d216be6a 100644 --- a/man/pvalue_maxcombo.Rd +++ b/man/pvalue_maxcombo.Rd @@ -53,9 +53,9 @@ xx <- sim_fixed_n( head(xx) # MaxCombo power estimate for cutoff at max of targeted events, minimum follow-up -p <- xx \%>\% - group_by(sim) \%>\% - group_map(~ pvalue_maxcombo(.x)) \%>\% +p <- xx |> + group_by(sim) |> + group_map(~ pvalue_maxcombo(.x)) |> unlist() mean(p < .025) } diff --git a/man/randomize_by_fixed_block.Rd b/man/randomize_by_fixed_block.Rd index b603780e..95317d32 100644 --- a/man/randomize_by_fixed_block.Rd +++ b/man/randomize_by_fixed_block.Rd @@ -27,11 +27,11 @@ library(dplyr) # Example 1 # 2:1 randomization with block size 3, treatments "A" and "B" -data.frame(x = 1:10) \%>\% mutate(Treatment = randomize_by_fixed_block(block = c("A", "B", "B"))) +data.frame(x = 1:10) |> mutate(Treatment = randomize_by_fixed_block(block = c("A", "B", "B"))) # Example 2 # Stratified randomization -data.frame(stratum = c(rep("A", 10), rep("B", 10))) \%>\% - group_by(stratum) \%>\% +data.frame(stratum = c(rep("A", 10), rep("B", 10))) |> + group_by(stratum) |> mutate(Treatment = randomize_by_fixed_block()) } diff --git a/man/sim_fixed_n.Rd b/man/sim_fixed_n.Rd index d9c85c51..cf97481e 100644 --- a/man/sim_fixed_n.Rd +++ b/man/sim_fixed_n.Rd @@ -102,24 +102,24 @@ xx <- sim_fixed_n( rho_gamma = data.frame(rho = 0, gamma = c(0, 1)) ) # Get power approximation for FH, data cutoff combination -xx \%>\% - group_by(cut, rho, gamma) \%>\% +xx |> + group_by(cut, rho, gamma) |> summarize(mean(z <= qnorm(.025))) # MaxCombo power estimate for cutoff at max of targeted events, minimum follow-up -p <- xx \%>\% - filter(cut != "Targeted events") \%>\% - group_by(sim) \%>\% - group_map(~ pvalue_maxcombo(.x)) \%>\% +p <- xx |> + filter(cut != "Targeted events") |> + group_by(sim) |> + group_map(~ pvalue_maxcombo(.x)) |> unlist() mean(p < .025) # MaxCombo estimate for targeted events cutoff -p <- xx \%>\% - filter(cut == "Targeted events") \%>\% - group_by(sim) \%>\% - group_map(~ pvalue_maxcombo(.x)) \%>\% +p <- xx |> + filter(cut == "Targeted events") |> + group_by(sim) |> + group_map(~ pvalue_maxcombo(.x)) |> unlist() mean(p < .025) diff --git a/man/simfix2simpwsurv.Rd b/man/simfix2simpwsurv.Rd index 1698b2d8..22fd7314 100644 --- a/man/simfix2simpwsurv.Rd +++ b/man/simfix2simpwsurv.Rd @@ -31,8 +31,6 @@ to analyze or otherwise evaluate individual simulations in ways that \code{\link[=sim_fixed_n]{sim_fixed_n()}} does not. } \examples{ -library(magrittr) - # Example 1 # Convert standard input simfix2simpwsurv() @@ -65,8 +63,8 @@ sim <- sim_pw_surv( ) # Cut after 200 events and do a stratified logrank test -dat <- sim \%>\% - cut_data_by_event(200) \%>\% # cut data - counting_process(arm = "experimental") \%>\% # convert format for tenFH - fh_weight(rho_gamma = data.frame(rho = 0, gamma = 0)) # stratified logrank +dat <- sim |> + cut_data_by_event(200) |> # Cut data + counting_process(arm = "experimental") |> # Convert format for fh_weight() + fh_weight(rho_gamma = data.frame(rho = 0, gamma = 0)) # Stratified logrank } diff --git a/man/simtrial-package.Rd b/man/simtrial-package.Rd index c8bd6e1a..5e01517a 100644 --- a/man/simtrial-package.Rd +++ b/man/simtrial-package.Rd @@ -6,7 +6,7 @@ \alias{simtrial-package} \title{simtrial: Clinical Trial Simulation} \description{ -simtrial provides some basic routines for simulating a clinical trial. The primary intent is to provide some tools to generate trial simulations for trials with time to event outcomes. Piecewise exponential failure rates and piecewise constant enrollment rates are the underlying mechanism used to simulate a broad range of scenarios. However, the basic generation of data is done using pipes to allow maximum flexibility for users to meet different needs. +Provides some basic routines for simulating a clinical trial. The primary intent is to provide some tools to generate trial simulations for trials with time to event outcomes. Piecewise exponential failure rates and piecewise constant enrollment rates are the underlying mechanism used to simulate a broad range of scenarios. However, the basic generation of data is done using pipes to allow maximum flexibility for users to meet different needs. } \seealso{ Useful links: diff --git a/tests/testthat/test-double_programming_sim_pw_surv.R b/tests/testthat/test-double_programming_sim_pw_surv.R index 9d0dbb59..277066a3 100644 --- a/tests/testthat/test-double_programming_sim_pw_surv.R +++ b/tests/testthat/test-double_programming_sim_pw_surv.R @@ -31,8 +31,8 @@ x <- sim_pw_surv( ) # prepare to test block -block1 <- x %>% dplyr::filter(stratum == "Low") -block2 <- x %>% dplyr::filter(stratum == "High") +block1 <- x |> dplyr::filter(stratum == "Low") +block2 <- x |> dplyr::filter(stratum == "High") bktest1 <- c() j <- 1 for (i in seq(1, floor(nrow(block1) / 4))) { diff --git a/tests/testthat/test-double_programming_wMB.R b/tests/testthat/test-double_programming_wMB.R index ca6ae55e..413aaa2a 100644 --- a/tests/testthat/test-double_programming_wMB.R +++ b/tests/testthat/test-double_programming_wMB.R @@ -16,8 +16,8 @@ test_mb_weight <- function(x, delay = 4) { # Test 1: for the situation of single stratum #### test_that("Validation passed for the situation of single stratum", { - x <- sim_pw_surv(n = 200) %>% - cut_data_by_event(125) %>% + x <- sim_pw_surv(n = 200) |> + cut_data_by_event(125) |> counting_process(arm = "experimental") out1 <- test_mb_weight(x, delay = 3) @@ -49,8 +49,8 @@ test_that("Validation passed for the situation of multiple stratum", { duration = rep(1, 4), rate = rep(.001, 4) ) - ) %>% - cut_data_by_event(125) %>% + ) |> + cut_data_by_event(125) |> counting_process(arm = "experimental") out1 <- test_mb_weight(x, delay = 3) @@ -82,8 +82,8 @@ test_that("Validation passed for the situation of a stratum with no records", { duration = rep(1, 4), rate = rep(.001, 4) ) - ) %>% - cut_data_by_event(125) %>% + ) |> + cut_data_by_event(125) |> counting_process(arm = "experimental") out <- mb_weight(x, delay = 1) diff --git a/tests/testthat/test-independent_test_counting_process.R b/tests/testthat/test-independent_test_counting_process.R index 4ba6cd08..4bd51f63 100644 --- a/tests/testthat/test-independent_test_counting_process.R +++ b/tests/testthat/test-independent_test_counting_process.R @@ -12,8 +12,8 @@ surv_to_count <- function(time, status, trt, strats) { surv = c(1, .survfit$surv[-n]) ) # ensure left continuous } - km <- db %>% - dplyr::group_by(strats) %>% + km <- db |> + dplyr::group_by(strats) |> dplyr::do(tidy_survfit(Surv(time, status) ~ 1, data = .)) # KM estimator by stratum and treatment Group Predicted at Specified Time @@ -25,24 +25,24 @@ surv_to_count <- function(time, status, trt, strats) { .x1 <- data.frame(time = pred_time, n.risk) # Number of Event - .x2 <- data.frame(time = .survfit$time, n.event = .survfit$n.event) %>% subset(n.event > 0) + .x2 <- data.frame(time = .survfit$time, n.event = .survfit$n.event) |> subset(n.event > 0) - merge(.x1, .x2, all = TRUE) %>% dplyr::mutate(n.event = dplyr::if_else(is.na(n.event), 0, n.event)) + merge(.x1, .x2, all = TRUE) |> dplyr::mutate(n.event = dplyr::if_else(is.na(n.event), 0, n.event)) } - km_by_trt <- db %>% - dplyr::group_by(strats, trt) %>% + km_by_trt <- db |> + dplyr::group_by(strats, trt) |> dplyr::do(pred_time = pred_survfit(km[km$strats == .$strats[1], ]$time, Surv(time, status) ~ 1, data = . - )) %>% - tidyr::unnest(cols = pred_time) %>% + )) |> + tidyr::unnest(cols = pred_time) |> dplyr::rename(tn.risk = n.risk, tn.event = n.event) # Log Rank Expectation Difference and Variance - res <- merge(km, km_by_trt, all = TRUE) %>% - dplyr::arrange(trt, strats, time) %>% + res <- merge(km, km_by_trt, all = TRUE) |> + dplyr::arrange(trt, strats, time) |> dplyr::mutate( OminusE = tn.event - tn.risk / n.risk * n.event, Var = (n.risk - tn.risk) * tn.risk * n.event * (n.risk - n.event) / n.risk^2 / (n.risk - 1) @@ -61,7 +61,7 @@ testthat::test_that("Counting Process Format without ties", { res_counting_process <- simtrial::counting_process(x, arm) res_test <- surv_to_count(time = x$tte, status = x$event, trt = x$treatment, strats = x$stratum) - res_test <- data.frame(subset(res_test, trt == 1)) %>% + res_test <- data.frame(subset(res_test, trt == 1)) |> subset(n.event > 0 & n.risk - tn.risk > 0 & tn.risk > 0) testthat::expect_equal(res_counting_process$o_minus_e, res_test$OminusE) @@ -80,7 +80,7 @@ testthat::test_that("Counting Process Format with ties", { res_counting_process <- counting_process(x, arm) res_test <- surv_to_count(time = x$tte, status = x$event, trt = x$treatment, strats = x$stratum) - res_test <- data.frame(subset(res_test, trt == 1)) %>% + res_test <- data.frame(subset(res_test, trt == 1)) |> subset(n.event > 0 & n.risk - tn.risk > 0 & tn.risk > 0) testthat::expect_equal(res_counting_process$o_minus_e, res_test$OminusE) diff --git a/tests/testthat/test-independent_test_get_cut_date_by_event.R b/tests/testthat/test-independent_test_get_cut_date_by_event.R index 5d0b4d6e..a958146e 100644 --- a/tests/testthat/test-independent_test_get_cut_date_by_event.R +++ b/tests/testthat/test-independent_test_get_cut_date_by_event.R @@ -3,10 +3,10 @@ test_that("get_cut_date_by_event returns the correct cut date", { y <- sim_pw_surv(n = 200) ycutdate <- get_cut_date_by_event(y, event) - x <- y %>% - dplyr::ungroup() %>% - dplyr::filter(fail == 1) %>% - dplyr::arrange(cte) %>% + x <- y |> + dplyr::ungroup() |> + dplyr::filter(fail == 1) |> + dplyr::arrange(cte) |> dplyr::slice(event) expect_equal(x$cte, ycutdate) }) diff --git a/tests/testthat/test-independent_test_pvalue_maxcombo.R b/tests/testthat/test-independent_test_pvalue_maxcombo.R index 06895bf0..563a31a6 100644 --- a/tests/testthat/test-independent_test_pvalue_maxcombo.R +++ b/tests/testthat/test-independent_test_pvalue_maxcombo.R @@ -1,7 +1,7 @@ testthat::test_that("the p-values correspond to pvalue_maxcombo", { set.seed(2022) # this part is a double programming - y <- sim_pw_surv(n = 300) %>% cut_data_by_event(30) + y <- sim_pw_surv(n = 300) |> cut_data_by_event(30) adjust.methods <- "asymp" wt <- list(a1 = c(0, 0), a2 = c(0, 1), a3 = c(1, 0), a4 = c(1, 1)) ties.method <- "efron" @@ -76,7 +76,7 @@ testthat::test_that("the p-values correspond to pvalue_maxcombo", { } p1 <- pval - a2 <- y %>% counting_process(arm = "experimental") + a2 <- y |> counting_process(arm = "experimental") aa <- fh_weight(a2, rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_corr = TRUE) p2 <- pvalue_maxcombo(z = aa) diff --git a/tests/testthat/test-independent_test_wlr.R b/tests/testthat/test-independent_test_wlr.R index fc0c00d9..97f4e86f 100644 --- a/tests/testthat/test-independent_test_wlr.R +++ b/tests/testthat/test-independent_test_wlr.R @@ -1,6 +1,6 @@ testthat::test_that("the z values match with the correspondings in fh_weight", { set.seed(1234) - y <- sim_pw_surv(n = 300) %>% cut_data_by_event(30) + y <- sim_pw_surv(n = 300) |> cut_data_by_event(30) adjust.methods <- "asymp" wt <- list(a1 = c(0, 0), a2 = c(0, 1), a3 = c(1, 0), a4 = c(1, 1)) ties.method <- "efron" @@ -19,7 +19,7 @@ testthat::test_that("the z values match with the correspondings in fh_weight", { })) tst.rslt <- attr(fit, "lrt") z1 <- tst.rslt$Z - a2 <- y %>% counting_process(arm = "experimental") + a2 <- y |> counting_process(arm = "experimental") aa <- fh_weight(a2, rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1))) z2 <- aa$z expect_equal(c(z1[1], z1[7:9]), z2, tolerance = 0.00001) @@ -27,7 +27,7 @@ testthat::test_that("the z values match with the correspondings in fh_weight", { testthat::test_that("fh_weight calculated correct correlation value when input a sequence of rho and gamma", { set.seed(123) - y <- sim_pw_surv(n = 300) %>% cut_data_by_event(30) + y <- sim_pw_surv(n = 300) |> cut_data_by_event(30) adjust.methods <- "asymp" wt <- list(a1 = c(0, 0), a2 = c(0, 1), a3 = c(1, 0), a4 = c(1, 1)) ties.method <- "efron" @@ -101,7 +101,7 @@ testthat::test_that("fh_weight calculated correct correlation value when input a pval <- pval2 / 2 } corr1 <- cor.tst[2:5, 2:5] - a2 <- y %>% counting_process(arm = "experimental") + a2 <- y |> counting_process(arm = "experimental") corr2 <- fh_weight(a2, rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_corr = TRUE) corr2 <- rbind(corr2$v1, corr2$v2, corr2$v3, corr2$v4) expect_equal(corr1, corr2, tolerance = 0.00001) diff --git a/vignettes/arbitraryhazard.Rmd b/vignettes/arbitraryhazard.Rmd index 12816009..59a4e5c4 100644 --- a/vignettes/arbitraryhazard.Rmd +++ b/vignettes/arbitraryhazard.Rmd @@ -35,14 +35,14 @@ dloglogis <- function(x, alpha = 1, beta = 4) { 1 / (1 + (x / alpha)^beta) } times <- (1:150) / 50 -xx <- data.frame(Times = times, Survival = dloglogis(times, alpha = .5, beta = 4)) %>% +xx <- data.frame(Times = times, Survival = dloglogis(times, alpha = .5, beta = 4)) |> mutate( duration = Times - lag(Times, default = 0), H = -log(Survival), rate = (H - lag(H, default = 0)) / duration / 3 - ) %>% + ) |> select(duration, rate) -ggplot(data = xx %>% mutate(Time = lag(cumsum(duration), default = 0)), aes(x = Time, y = rate)) + +ggplot(data = xx |> mutate(Time = lag(cumsum(duration), default = 0)), aes(x = Time, y = rate)) + geom_line() ``` @@ -58,7 +58,7 @@ x <- sim_pw_surv( n = 250, # sample size block = block, enroll_rate = enroll_rate, - fail_rate = xx %>% mutate(stratum = "All", treatment = tx, period = 1:n(), stratum = "All"), + fail_rate = xx |> mutate(stratum = "All", treatment = tx, period = 1:n(), stratum = "All"), dropout_rate = dropout_rate ) ``` @@ -66,7 +66,7 @@ x <- sim_pw_surv( We assume the entire study lasts 3 years ```{r,fig.height=4,fig.width=7.5} -y <- x %>% cut_data_by_date(3) +y <- x |> cut_data_by_date(3) head(y) ``` @@ -83,5 +83,5 @@ We overlay the actual rates in red. ```{r,fig.height=4,fig.width=7.5} fit <- bshazard(Surv(tte, event) ~ 1, data = y, nk = 120) plot(fit, conf.int = TRUE, xlab = "Time", xlim = c(0, 3), ylim = c(0, 2.5), lwd = 2) -lines(x = times, y = (xx %>% mutate(Time = lag(cumsum(duration), default = 0)))$rate, col = 2) +lines(x = times, y = (xx |> mutate(Time = lag(cumsum(duration), default = 0)))$rate, col = 2) ``` diff --git a/vignettes/modestWLRTVignette.Rmd b/vignettes/modestWLRTVignette.Rmd index 878e232e..f9a9365d 100644 --- a/vignettes/modestWLRTVignette.Rmd +++ b/vignettes/modestWLRTVignette.Rmd @@ -69,7 +69,7 @@ MBdelay <- sim_pw_surv( enroll_rate = enroll_rate, fail_rate = xpar$fail_rate, dropout_rate = xpar$dropout_rate -) %>% +) |> cut_data_by_date(studyDuration) fit <- survfit(Surv(tte, event) ~ treatment, data = MBdelay) plot(fit, col = 1:2, mark = "|", xaxt = "n") @@ -86,9 +86,9 @@ This requires generating weights and then computing the test. We begin with the default of `w_max=Inf` which corresponds to the original @MagirrBurman test and set the time until maximum weight $\tau$ with `delay = 6`. ```{r} -ZMB <- MBdelay %>% - counting_process(arm = "experimental") %>% - mb_weight(delay = 6) %>% +ZMB <- MBdelay |> + counting_process(arm = "experimental") |> + mb_weight(delay = 6) |> summarize(S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), z = S / sqrt(V)) # Compute p-value of modestly weighted logrank of Magirr-Burman pnorm(ZMB$z) @@ -97,9 +97,9 @@ pnorm(ZMB$z) Now we set the maximum weight to be 2 as in @Magirr2021 and set the `delay=Inf` so that the maximum weight begins at the observed median of the observed combined treatment Kaplan-Meier curve. ```{r} -ZMB <- MBdelay %>% - counting_process(arm = "experimental") %>% - mb_weight(delay = Inf, w_max = 2) %>% +ZMB <- MBdelay |> + counting_process(arm = "experimental") |> + mb_weight(delay = Inf, w_max = 2) |> summarize(S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), z = S / sqrt(V)) # Compute p-value of modestly weighted logrank of Magirr-Burman pnorm(ZMB$z) @@ -113,9 +113,9 @@ and let $\gamma=0, \rho = -1/2.$ ```{r} w_max <- 2 -Z_modified_FH <- MBdelay %>% - counting_process(arm = "experimental") %>% - mutate(w = pmin(w_max, 1 / s)) %>% +Z_modified_FH <- MBdelay |> + counting_process(arm = "experimental") |> + mutate(w = pmin(w_max, 1 / s)) |> summarize(S = sum(o_minus_e * w), V = sum(var_o_minus_e * w^2), z = S / sqrt(V)) # Compute p-value of modestly weighted logrank of Magirr-Burman pnorm(Z_modified_FH$z) @@ -156,7 +156,7 @@ FHwn <- sim_pw_surv( enroll_rate = enroll_rate, fail_rate = xpar$fail_rate, dropout_rate = xpar$dropout_rate -) %>% +) |> cut_data_by_date(studyDuration) fit <- survfit(Surv(tte, event) ~ treatment, data = FHwn) plot(fit, col = 1:2, mark = "|", xaxt = "n") @@ -166,9 +166,9 @@ axis(1, xaxp = c(0, 36, 6)) We perform a logrank and weighted logrank tests as suggested for more limited downweighting by follows: ```{r} -xx <- FHwn %>% - counting_process(arm = "experimental") %>% - fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) %>% +xx <- FHwn |> + counting_process(arm = "experimental") |> + fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1)), return_corr = TRUE) |> mutate(p = pnorm(z)) xx ``` @@ -176,16 +176,16 @@ xx Now for a MaxCombo test with the above component tests, we have p-value of ```{r} -xx %>% pvalue_maxcombo() +xx |> pvalue_maxcombo() ``` Next, we consider the @MagirrBurman modestly weighted logrank test with down-weighting specified for the first 6 months but a maximum weight of 2. This requires generating weights and then computing the test. ```{r} -ZMB <- FHwn %>% - counting_process(arm = "experimental") %>% - mb_weight(delay = 6, w_max = 2) %>% +ZMB <- FHwn |> + counting_process(arm = "experimental") |> + mb_weight(delay = 6, w_max = 2) |> summarize(S = sum(o_minus_e * mb_weight), V = sum(var_o_minus_e * mb_weight^2), z = S / sqrt(V)) # Compute p-value of modestly weighted logrank of Magirr-Burman @@ -196,9 +196,9 @@ Finally, we consider weighted logrank tests with less down-weighting. Results are quite similar to the results with greater down-weighting. ```{r} -xx <- FHwn %>% - counting_process(arm = "experimental") %>% - fh_weight(rho_gamma = data.frame(rho = c(0, 0, .5), gamma = c(0, .5, .5)), return_corr = TRUE) %>% +xx <- FHwn |> + counting_process(arm = "experimental") |> + fh_weight(rho_gamma = data.frame(rho = c(0, 0, .5), gamma = c(0, .5, .5)), return_corr = TRUE) |> mutate(p = pnorm(z)) xx ``` @@ -206,7 +206,7 @@ xx Now for a MaxCombo test with the above component tests, we have p-value of ```{r} -xx %>% pvalue_maxcombo() +xx |> pvalue_maxcombo() ``` Thus, with less down-weighting the MaxCombo test appears less problematic. diff --git a/vignettes/parallel.Rmd b/vignettes/parallel.Rmd index 838de33a..1a33cc01 100644 --- a/vignettes/parallel.Rmd +++ b/vignettes/parallel.Rmd @@ -134,11 +134,11 @@ over that of enrollment 1, which is We also see that there is a distinction between the duration of the study under the proposed enrollment strategies. ```{r sequential-display-results, eval=FALSE, echo=FALSE} -seq_result1 %>% - head(5) %>% +seq_result1 |> + head(5) |> kable(digits = 2) -seq_result2 %>% - head(5) %>% +seq_result2 |> + head(5) |> kable(digits = 2) ``` diff --git a/vignettes/pvalue_maxcomboVignette.Rmd b/vignettes/pvalue_maxcomboVignette.Rmd index 2cfd3298..a3d007ca 100644 --- a/vignettes/pvalue_maxcomboVignette.Rmd +++ b/vignettes/pvalue_maxcomboVignette.Rmd @@ -50,7 +50,7 @@ library(dplyr) ```{r} x <- sim_fixed_n(n_sim = 1, timing_type = 5, rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1))) -x %>% kable(digits = 2) +x |> kable(digits = 2) ``` Once you have this format, the MaxCombo p-value per @Karrison2016, @NPHWGDesign can be computed as follows (note that you will need to have the package **mvtnorm** installed): @@ -66,25 +66,25 @@ Again, we use defaults for that routine. ```{r,message=FALSE,warning=FALSE,cache=FALSE} s <- sim_pw_surv(n = 100) -head(s) %>% kable(digits = 2) +head(s) |> kable(digits = 2) ``` Once generated, we need to cut the data for analysis. Here we cut after 75 events. ```{r,warning=FALSE,message=FALSE} -x <- s %>% cut_data_by_event(75) -head(x) %>% kable(digits = 2) +x <- s |> cut_data_by_event(75) +head(x) |> kable(digits = 2) ``` Now we can analyze this data. We begin with `s` to show how this can be done in a single line. In this case, we use the 4 test combination suggested in @NPHWGSimulation, @NPHWGDesign. ```{r,warning=FALSE,message=FALSE} -z <- s %>% - cut_data_by_event(75) %>% - counting_process(arm = "experimental") %>% +z <- s |> + cut_data_by_event(75) |> + counting_process(arm = "experimental") |> fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_corr = TRUE) -z %>% kable(digits = 2) +z |> kable(digits = 2) ``` Now we compute our p-value as before: @@ -97,7 +97,7 @@ Suppose we want the p-value just based on the logrank and FH(0,1) and FH(1,0) as We remove the rows and columns associated with FH(0,0) and FH(1,1) and then apply `pvalue_maxcombo()`. ```{r,warning=FALSE,message=FALSE} -pvalue_maxcombo(z %>% select(-c(v1, v4)) %>% filter((rho == 0 & gamma == 1) | (rho == 1 & gamma == 0))) +pvalue_maxcombo(z |> select(-c(v1, v4)) |> filter((rho == 0 & gamma == 1) | (rho == 1 & gamma == 0))) ``` ### Using survival data in another format @@ -108,21 +108,21 @@ In this case we use the small `aml` dataset from the **survival** package. ```{r,warning=FALSE,message=FALSE} library(survival) -head(aml) %>% kable() +head(aml) |> kable() ``` We rename variables and create a stratum variable as follows: ```{r,warning=FALSE,message=FALSE} -x <- aml %>% transmute(tte = time, event = status, stratum = "All", treatment = as.character(x)) -head(x) %>% kable() +x <- aml |> transmute(tte = time, event = status, stratum = "All", treatment = as.character(x)) +head(x) |> kable() ``` Now we analyze the data with a MaxCombo with the logrank and FH(0,1) and compute a p-value. ```{r,warning=FALSE,message=FALSE} -x %>% - counting_process(arm = "Maintained") %>% - fh_weight(rho_gamma = data.frame(rho = 0, gamma = c(0, 1)), return_corr = TRUE) %>% +x |> + counting_process(arm = "Maintained") |> + fh_weight(rho_gamma = data.frame(rho = 0, gamma = c(0, 1)), return_corr = TRUE) |> pvalue_maxcombo() ``` @@ -134,7 +134,7 @@ We now consider the example simulation from the `pvalue_maxcombo()` help file to # Only use cut events + min follow-up xx <- sim_fixed_n(n_sim = 100, timing_type = 5, rho_gamma = data.frame(rho = c(0, 0, 1), gamma = c(0, 1, 1))) # MaxCombo power estimate for cutoff at max of targeted events, minimum follow-up -p <- unlist(xx %>% group_by(sim) %>% group_map(~ pvalue_maxcombo(.x))) +p <- unlist(xx |> group_by(sim) |> group_map(~ pvalue_maxcombo(.x))) mean(p < .001) ``` @@ -145,14 +145,14 @@ The latter can be done without having to re-run all simulations as follows, demo ```{r,cache=FALSE,warning=FALSE,message=FALSE} # Only use cuts for events and events + min follow-up xx <- sim_fixed_n(n_sim = 100, timing_type = c(2, 5), rho_gamma = data.frame(rho = 0, gamma = c(0, 1))) -head(xx) %>% kable(digits = 2) +head(xx) |> kable(digits = 2) ``` Now we compute a p-value separately for each cut type, first for targeted event count. ```{r,warning=FALSE,message=FALSE} # Subset to targeted events cutoff tests -p <- unlist(xx %>% filter(cut == "Targeted events") %>% group_by(sim) %>% group_map(~ pvalue_maxcombo(.x))) +p <- unlist(xx |> filter(cut == "Targeted events") |> group_by(sim) |> group_map(~ pvalue_maxcombo(.x))) mean(p < .025) ``` @@ -160,7 +160,7 @@ Now we use the later of targeted events and minimum follow-up cutoffs. ```{r,warning=FALSE,message=FALSE} # Subset to targeted events cutoff tests -p <- unlist(xx %>% filter(cut != "Targeted events") %>% group_by(sim) %>% group_map(~ pvalue_maxcombo(.x))) +p <- unlist(xx |> filter(cut != "Targeted events") |> group_by(sim) |> group_map(~ pvalue_maxcombo(.x))) mean(p < .025) ``` diff --git a/vignettes/simtrialroutines.Rmd b/vignettes/simtrialroutines.Rmd index ee80ca67..bdc1ef51 100644 --- a/vignettes/simtrialroutines.Rmd +++ b/vignettes/simtrialroutines.Rmd @@ -131,7 +131,7 @@ x <- sim_pw_surv( fail_rate = fail_rate, dropout_rate = dropout_rate ) -head(x) %>% kable(digits = 2) +head(x) |> kable(digits = 2) ``` ## Cutting data for analysis @@ -143,7 +143,7 @@ Observations enrolled after the input `cut_date` are deleted and events and cens ```{r} y <- cut_data_by_date(x, cut_date = 5) -head(y) %>% kable(digits = 2) +head(y) |> kable(digits = 2) ``` For instance, if we wish to cut the entire dataset when 50 events are observed in the Positive stratum we can use the `get_cut_date_by_event` function as follows: @@ -170,7 +170,7 @@ The counting process format is further discussed in the next section where we co ```{r} ten150 <- counting_process(y150, arm = "experimental") -head(ten150) %>% kable(digits = 2) +head(ten150) |> kable(digits = 2) ``` ## Logrank and weighted logrank testing @@ -196,7 +196,7 @@ c(z, pnorm(z)) For Fleming-Harrington tests, a routine has been built to do these tests for you: ```{r} -fh_weight(x = ten150, rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1))) %>% kable(digits = 2) +fh_weight(x = ten150, rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1))) |> kable(digits = 2) ``` If we wanted to take the minimum of these for a MaxCombo test, we would first use `fh_weight()` to compute a correlation matrix for the above z-statistics as follows. @@ -204,8 +204,8 @@ Note that the ordering of `rho_gamma` and `g` in the argument list is opposite o The correlation matrix for the `z`-values is now in `V1`-`V4`. ```{r,message=FALSE} -x <- ten150 %>% fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_corr = TRUE) -x %>% kable(digits = 2) +x <- ten150 |> fh_weight(rho_gamma = data.frame(rho = c(0, 0, 1, 1), gamma = c(0, 1, 0, 1)), return_corr = TRUE) +x |> kable(digits = 2) ``` We can compute a p-value for the MaxCombo as follows using the `pmvnorm` function from the `mvtnorm` package. Note the arguments for `GenzBretz` which are more stringent than the defaults; we have also used these more stringent parameters in the example in the help file. @@ -258,14 +258,14 @@ sim_fixed_n( block = block, # Block for treatment timing_type = 1:5, # Use all possible data cutoff methods rho_gamma = rho_gamma # FH test(s) to use; in this case, logrank -) %>% kable(digits = 2) +) |> kable(digits = 2) ``` If you look carefully, you should be asking why the cutoff with the planned number of events is so different than the other data cutoff methods. To explain, we note that generally you will want `sample_size` above to match the enrollment specified in `enroll_rate`: ```{r} -enroll_rate %>% summarize("Targeted enrollment based on input enrollment rates" = sum(duration * rate)) +enroll_rate |> summarize("Targeted enrollment based on input enrollment rates" = sum(duration * rate)) ``` The targeted enrollment takes, on average, 30 months longer than the sum of the enrollment durations in `enroll_rate` (14 months) at the input enrollment rates. To achieve the input `sample_size` of 500, the final enrollment rate is assumed to be steady state and extends in each simulation until the targeted enrollment is achieved. The planned duration of the trial is taken as 30 months as specified in `totalDuration`. The targeted minimum follow-up is From c6e5e194363012c3291df71c2928409a09dfb801 Mon Sep 17 00:00:00 2001 From: Nan Xiao Date: Mon, 27 Nov 2023 00:02:45 -0500 Subject: [PATCH 2/2] Fix typo on prevalence --- tests/testthat/test-unvalidated-data.table.R | 6 +++--- tests/testthat/test-unvalidated-early_zero_weight.R | 12 ++++++------ 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/tests/testthat/test-unvalidated-data.table.R b/tests/testthat/test-unvalidated-data.table.R index 1ee1d0ad..2f3c7155 100644 --- a/tests/testthat/test-unvalidated-data.table.R +++ b/tests/testthat/test-unvalidated-data.table.R @@ -88,12 +88,12 @@ test_that("functions that use data.table do not modify input data table", { n <- 500 # Two strata stratum <- c("Biomarker-positive", "Biomarker-negative") - prevelance_ratio <- c(0.6, 0.4) + prevalence_ratio <- c(0.6, 0.4) # Enrollment rate enroll_rate <- gsDesign2::define_enroll_rate( stratum = rep(stratum, each = 2), duration = c(2, 10, 2, 10), - rate = c(c(1, 4) * prevelance_ratio[1], c(1, 4) * prevelance_ratio[2]) + rate = c(c(1, 4) * prevalence_ratio[1], c(1, 4) * prevalence_ratio[2]) ) enroll_rate$rate <- enroll_rate$rate * n / sum(enroll_rate$duration * enroll_rate$rate) # Failure rate @@ -114,7 +114,7 @@ test_that("functions that use data.table do not modify input data table", { x <- sim_pw_surv( n = n, # Sample size # Stratified design with prevalence ratio of 6:4 - stratum = data.frame(stratum = stratum, p = prevelance_ratio), + stratum = data.frame(stratum = stratum, p = prevalence_ratio), # Randomization ratio block = c("control", "control", "experimental", "experimental"), enroll_rate = enroll_rate, # Enrollment rate diff --git a/tests/testthat/test-unvalidated-early_zero_weight.R b/tests/testthat/test-unvalidated-early_zero_weight.R index 18a5d73f..9ccf04fd 100644 --- a/tests/testthat/test-unvalidated-early_zero_weight.R +++ b/tests/testthat/test-unvalidated-early_zero_weight.R @@ -16,12 +16,12 @@ test_that("early_zero_weight() with stratified data", { n <- 500 # Two strata stratum <- c("Biomarker-positive", "Biomarker-negative") - prevelance_ratio <- c(0.6, 0.4) + prevalence_ratio <- c(0.6, 0.4) # Enrollment rate enroll_rate <- gsDesign2::define_enroll_rate( stratum = rep(stratum, each = 2), duration = c(2, 10, 2, 10), - rate = c(c(1, 4) * prevelance_ratio[1], c(1, 4) * prevelance_ratio[2]) + rate = c(c(1, 4) * prevalence_ratio[1], c(1, 4) * prevalence_ratio[2]) ) enroll_rate$rate <- enroll_rate$rate * n / sum(enroll_rate$duration * enroll_rate$rate) # Failure rate @@ -42,7 +42,7 @@ test_that("early_zero_weight() with stratified data", { input <- sim_pw_surv( n = n, # Sample size # Stratified design with prevalence ratio of 6:4 - stratum = data.frame(stratum = stratum, p = prevelance_ratio), + stratum = data.frame(stratum = stratum, p = prevalence_ratio), # Randomization ratio block = c("control", "control", "experimental", "experimental"), enroll_rate = enroll_rate, # Enrollment rate @@ -63,12 +63,12 @@ test_that("early_zero_weight() fails with bad input", { n <- 500 # Two strata stratum <- c("Biomarker-positive", "Biomarker-negative") - prevelance_ratio <- c(0.6, 0.4) + prevalence_ratio <- c(0.6, 0.4) # Enrollment rate enroll_rate <- gsDesign2::define_enroll_rate( stratum = rep(stratum, each = 2), duration = c(2, 10, 2, 10), - rate = c(c(1, 4) * prevelance_ratio[1], c(1, 4) * prevelance_ratio[2]) + rate = c(c(1, 4) * prevalence_ratio[1], c(1, 4) * prevalence_ratio[2]) ) enroll_rate$rate <- enroll_rate$rate * n / sum(enroll_rate$duration * enroll_rate$rate) # Failure rate @@ -89,7 +89,7 @@ test_that("early_zero_weight() fails with bad input", { input <- sim_pw_surv( n = n, # Sample size # Stratified design with prevalence ratio of 6:4 - stratum = data.frame(stratum = stratum, p = prevelance_ratio), + stratum = data.frame(stratum = stratum, p = prevalence_ratio), # Randomization ratio block = c("control", "control", "experimental", "experimental"), enroll_rate = enroll_rate, # Enrollment rate