From 62229b2fe72a76325a19e6ff0561024eeb8d65c1 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 18:04:06 -0400 Subject: [PATCH 01/12] Refactor estimate_incidence_rate --- R/incidence_rate.R | 146 +++++++++++++++++-------- R/utils_default_stats_formats_labels.R | 2 +- man/incidence_rate.Rd | 34 +++++- 3 files changed, 129 insertions(+), 53 deletions(-) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 7c86c448e9..58a9945c80 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -6,6 +6,7 @@ #' as incidence rate. Primary analysis variable is the person-years at risk. #' #' @inheritParams argument_convention +#' @param id_var (`string`)\cr name of variable used as patient identifier if `"n_unique"` is included in `.stats`. #' @param control (`list`)\cr parameters for estimation details, specified by using #' the helper function [control_incidence_rate()]. Possible parameter options are: #' * `conf_level` (`proportion`)\cr confidence level for the estimated incidence rate. @@ -14,7 +15,7 @@ #' * `input_time_unit` (`string`)\cr `day`, `week`, `month`, or `year` (default) #' indicating time unit for data input. #' * `num_pt_year` (`numeric`)\cr time unit for desired output (in person-years). -#' @param n_events (`integer(1)`)\cr number of events observed. +#' @param n_events (`string`)\cr name of integer variable indicating whether an event has been observed (1) or not (0). #' @param .stats (`character`)\cr statistics to select for the table. Run `get_stats("estimate_incidence_rate")` #' to see available statistics for this function. #' @@ -33,30 +34,29 @@ NULL #' - `n_events`: Total number of events observed. #' - `rate`: Estimated incidence rate. #' - `rate_ci`: Confidence interval for the incidence rate. +#' - `n_unique`: Total number of patients with at least one event observed. +#' - `n_rate`: Total number of events observed & estimated incidence rate. #' #' @keywords internal s_incidence_rate <- function(df, .var, n_events, - is_event, + is_event = lifecycle::deprecated(), + id_var = "USUBJID", control = control_incidence_rate()) { - if (!missing(is_event)) { - warning("argument is_event will be deprecated. Please use n_events.") - - if (missing(n_events)) { - assert_df_with_variables(df, list(tte = .var, is_event = is_event)) - checkmate::assert_string(.var) - checkmate::assert_logical(df[[is_event]], any.missing = FALSE) - checkmate::assert_numeric(df[[.var]], any.missing = FALSE) - n_events <- is_event - } - } else { - assert_df_with_variables(df, list(tte = .var, n_events = n_events)) - checkmate::assert_string(.var) - checkmate::assert_numeric(df[[.var]], any.missing = FALSE) - checkmate::assert_integer(df[[n_events]], any.missing = FALSE) + if (lifecycle::is_present(is_event)) { + lifecycle::deprecate_warn( + "0.9.6", "s_incidence_rate(is_event)", "s_incidence_rate(n_events)" + ) + n_events <- as.numeric(is_event) } + assert_df_with_variables(df, list(tte = .var, n_events = n_events)) + checkmate::assert_string(.var) + checkmate::assert_numeric(df[[.var]], any.missing = FALSE) + checkmate::assert_integerish(df[[n_events]], any.missing = FALSE) + + n_unique <- n_available(unique(df[[id_var]][df$n_events == 1])) input_time_unit <- control$input_time_unit num_pt_year <- control$num_pt_year conf_level <- control$conf_level @@ -77,7 +77,12 @@ s_incidence_rate <- function(df, person_years = formatters::with_label(person_years, "Total patient-years at risk"), n_events = formatters::with_label(n_events, "Number of adverse events observed"), rate = formatters::with_label(result$rate, paste("AE rate per", num_pt_year, "patient-years")), - rate_ci = formatters::with_label(result$rate_ci, f_conf_level(conf_level)) + rate_ci = formatters::with_label(result$rate_ci, f_conf_level(conf_level)), + n_unique = formatters::with_label(n_unique, "Total number of patients with at least one adverse event"), + n_rate = formatters::with_label( + c(n_events, result$rate), + paste("Number of adverse events observed (AE rate per", num_pt_year, "patient-years)") + ) ) } @@ -88,15 +93,57 @@ s_incidence_rate <- function(df, #' * `a_incidence_rate()` returns the corresponding list with formatted [rtables::CellValue()]. #' #' @keywords internal -a_incidence_rate <- make_afun( - s_incidence_rate, - .formats = c( - "person_years" = "xx.x", - "n_events" = "xx", - "rate" = "xx.xx", - "rate_ci" = "(xx.xx, xx.xx)" +a_incidence_rate <- function(df, + labelstr = "", + .var, + .df_row, + n_events, + id_var = "USUBJID", + control = control_incidence_rate(), + .stats = NULL, + .formats = c( + "person_years" = "xx.x", + "n_events" = "xx", + "rate" = "xx.xx", + "rate_ci" = "(xx.xx, xx.xx)", + "n_unique" = "xx", + "n_rate" = "xx (xx.x)" + ), + .labels = NULL, + .indent_mods = NULL, + na_str = default_na_str(), + label_fmt = "%s - %.labels") { + x_stats <- s_incidence_rate( + df = df, .var = .var, n_events = n_events, id_var = id_var, control = control + ) + if (is.null(unlist(x_stats))) { + return(NULL) + } + if (is.null(.formats)) .formats <- formals()$.formats %>% eval() + if (is.null(.labels)) .labels <- sapply(x_stats, \(x) attributes(x)$label) + if (nzchar(labelstr) > 0) { + .labels <- sapply(.labels, \(x) gsub("%.labels", x, gsub("%s", labelstr, label_fmt))) + } + + # Fill in with formatting defaults if needed + .stats <- get_stats("estimate_incidence_rate", stats_in = .stats) + .formats <- get_formats_from_stats(.stats, .formats) + .labels <- get_labels_from_stats(.stats, .labels) + .indent_mods <- get_indents_from_stats(.stats, .indent_mods) + + x_stats <- x_stats[.stats] + + # Auto format handling + .formats <- apply_auto_formatting(.formats, x_stats, .df_row, .var) + + in_rows( + .list = x_stats, + .formats = .formats, + .labels = .labels, + .indent_mods = .indent_mods, + .format_na_strs = na_str ) -) +} #' @describeIn incidence_rate Layout-creating function which can take statistics function arguments #' and additional format arguments. This function is a wrapper for [rtables::analyze()]. @@ -115,8 +162,7 @@ a_incidence_rate <- make_afun( #' AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), #' ARM = factor(c("A", "A", "A", "B", "B", "B")) #' ) %>% -#' mutate(is_event = CNSR == 0) %>% -#' mutate(n_events = as.integer(is_event)) +#' mutate(n_events = 1 - CNSR) #' #' basic_table() %>% #' split_cols_by("ARM") %>% @@ -136,9 +182,12 @@ a_incidence_rate <- make_afun( estimate_incidence_rate <- function(lyt, vars, n_events, + id_var = "USUBJID", control = control_incidence_rate(), na_str = default_na_str(), nested = TRUE, + summarize = FALSE, + label_fmt = "%s - %.labels", ..., show_labels = "hidden", table_names = vars, @@ -146,26 +195,31 @@ estimate_incidence_rate <- function(lyt, .formats = NULL, .labels = NULL, .indent_mods = NULL) { - extra_args <- list(n_events = n_events, control = control, ...) - - afun <- make_afun( - a_incidence_rate, - .stats = .stats, - .formats = .formats, - .labels = .labels, - .indent_mods = .indent_mods + extra_args <- c( + list(.stats = .stats, .formats = .formats, .labels = .labels, .indent_mods = .indent_mods, na_str = na_str), + list(n_events = n_events, id_var = id_var, control = control, label_fmt = label_fmt, ...) ) - analyze( - lyt, - vars, - show_labels = show_labels, - table_names = table_names, - afun = afun, - na_str = na_str, - nested = nested, - extra_args = extra_args - ) + if (!summarize) { + analyze( + lyt, + vars, + show_labels = show_labels, + table_names = table_names, + afun = a_incidence_rate, + na_str = na_str, + nested = nested, + extra_args = extra_args + ) + } else { + summarize_row_groups( + lyt, + vars, + cfun = a_incidence_rate, + na_str = na_str, + extra_args = extra_args + ) + } } #' Helper functions for incidence rate diff --git a/R/utils_default_stats_formats_labels.R b/R/utils_default_stats_formats_labels.R index 824f946061..2b744d5a0e 100644 --- a/R/utils_default_stats_formats_labels.R +++ b/R/utils_default_stats_formats_labels.R @@ -382,7 +382,7 @@ tern_default_stats <- list( count_patients_with_flags = c("n", "count", "count_fraction", "count_fraction_fixed_dp", "n_blq"), count_values = c("n", "count", "count_fraction", "count_fraction_fixed_dp", "n_blq"), coxph_pairwise = c("pvalue", "hr", "hr_ci", "n_tot", "n_tot_events"), - estimate_incidence_rate = c("person_years", "n_events", "rate", "rate_ci"), + estimate_incidence_rate = c("person_years", "n_events", "rate", "rate_ci", "n_unique", "n_rate"), estimate_multinomial_response = c("n_prop", "prop_ci"), estimate_odds_ratio = c("or_ci", "n_tot"), estimate_proportion = c("n_prop", "prop_ci"), diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index d0d0c4d337..1ea9ce8bd6 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -11,9 +11,12 @@ estimate_incidence_rate( lyt, vars, n_events, + id_var = "USUBJID", control = control_incidence_rate(), na_str = default_na_str(), nested = TRUE, + summarize = FALSE, + label_fmt = "\%s - \%.labels", ..., show_labels = "hidden", table_names = vars, @@ -27,16 +30,26 @@ s_incidence_rate( df, .var, n_events, - is_event, + is_event = lifecycle::deprecated(), + id_var = "USUBJID", control = control_incidence_rate() ) a_incidence_rate( df, + labelstr = "", .var, + .df_row, n_events, - is_event, - control = control_incidence_rate() + id_var = "USUBJID", + control = control_incidence_rate(), + .stats = NULL, + .formats = c(person_years = "xx.x", n_events = "xx", rate = "xx.xx", rate_ci = + "(xx.xx, xx.xx)", n_unique = "xx", n_rate = "xx (xx.x)"), + .labels = NULL, + .indent_mods = NULL, + na_str = default_na_str(), + label_fmt = "\%s - \%.labels" ) } \arguments{ @@ -44,7 +57,9 @@ a_incidence_rate( \item{vars}{(\code{character})\cr variable names for the primary analysis variable to be iterated over.} -\item{n_events}{(\code{integer(1)})\cr number of events observed.} +\item{n_events}{(\code{string})\cr name of integer variable indicating whether an event has been observed (1) or not (0).} + +\item{id_var}{(\code{string})\cr name of variable used as patient identifier if \code{"n_unique"} is included in \code{.stats}.} \item{control}{(\code{list})\cr parameters for estimation details, specified by using the helper function \code{\link[=control_incidence_rate]{control_incidence_rate()}}. Possible parameter options are: @@ -87,6 +102,12 @@ unmodified default behavior. Can be negative.} by a statistics function.} \item{is_event}{(\code{flag})\cr \code{TRUE} if event, \code{FALSE} if time to event is censored.} + +\item{labelstr}{(\code{string})\cr label of the level of the parent split currently being summarized +(must be present as second argument in Content Row Functions). See \code{\link[rtables:summarize_row_groups]{rtables::summarize_row_groups()}} +for more information.} + +\item{.df_row}{(\code{data.frame})\cr data frame across all of the columns for the given row split.} } \value{ \itemize{ @@ -102,6 +123,8 @@ the statistics from \code{s_incidence_rate()} to the table layout. \item \code{n_events}: Total number of events observed. \item \code{rate}: Estimated incidence rate. \item \code{rate_ci}: Confidence interval for the incidence rate. +\item \code{n_unique}: Total number of patients with at least one event observed. +\item \code{n_rate}: Total number of events observed & estimated incidence rate. } } @@ -136,8 +159,7 @@ df <- data.frame( AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), ARM = factor(c("A", "A", "A", "B", "B", "B")) ) \%>\% - mutate(is_event = CNSR == 0) \%>\% - mutate(n_events = as.integer(is_event)) + mutate(n_events = 1 - CNSR) basic_table() \%>\% split_cols_by("ARM") \%>\% From d636d2f3f7c1955fa02411bbcf15905e922374b7 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 18:19:58 -0400 Subject: [PATCH 02/12] Update docs --- R/incidence_rate.R | 9 +++++++++ man/incidence_rate.Rd | 7 +++++++ 2 files changed, 16 insertions(+) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 58a9945c80..71c6a542d2 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -18,6 +18,11 @@ #' @param n_events (`string`)\cr name of integer variable indicating whether an event has been observed (1) or not (0). #' @param .stats (`character`)\cr statistics to select for the table. Run `get_stats("estimate_incidence_rate")` #' to see available statistics for this function. +#' @param summarize (`flag`)\cr whether the function should act as an analyze function (`summarize = FALSE`), or a +#' summarize function (`summarize = TRUE`). Defaults to `FALSE`. +#' @param label_fmt (`string`)\cr how labels should be formatted after a row split occurs if `summarize = TRUE`. The +#' string should use `"%s"` to represent row split levels, and `"%.labels"` to represent labels supplied to the +#' `.labels` argument. Defaults to `"%s - %.labels"`. #' #' @seealso [control_incidence_rate()] and helper functions [h_incidence_rate]. #' @@ -53,6 +58,8 @@ s_incidence_rate <- function(df, assert_df_with_variables(df, list(tte = .var, n_events = n_events)) checkmate::assert_string(.var) + checkmate::assert_string(n_events) + checkmate::assert_string(id_var) checkmate::assert_numeric(df[[.var]], any.missing = FALSE) checkmate::assert_integerish(df[[n_events]], any.missing = FALSE) @@ -113,6 +120,8 @@ a_incidence_rate <- function(df, .indent_mods = NULL, na_str = default_na_str(), label_fmt = "%s - %.labels") { + checkmate::assert_string(label_fmt) + x_stats <- s_incidence_rate( df = df, .var = .var, n_events = n_events, id_var = id_var, control = control ) diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index 1ea9ce8bd6..3050766ace 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -78,6 +78,13 @@ indicating time unit for data input. possible (\code{TRUE}, the default) or as a new top-level element (\code{FALSE}). Ignored if it would nest a split. underneath analyses, which is not allowed.} +\item{summarize}{(\code{flag})\cr whether the function should act as an analyze function (\code{summarize = FALSE}), or a +summarize function (\code{summarize = TRUE}). Defaults to \code{FALSE}.} + +\item{label_fmt}{(\code{string})\cr how labels should be formatted after a row split occurs if \code{summarize = TRUE}. The +string should use \code{"\%s"} to represent row split levels, and \code{"\%.labels"} to represent labels supplied to the +\code{.labels} argument. Defaults to \code{"\%s - \%.labels"}.} + \item{...}{additional arguments for the lower level functions.} \item{show_labels}{(\code{string})\cr label visibility: one of "default", "visible" and "hidden".} From 12db74109333fba77b8fcbf4f0134b80171a09cf Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 18:22:54 -0400 Subject: [PATCH 03/12] Update docs --- R/incidence_rate.R | 6 +++--- man/incidence_rate.Rd | 6 +++--- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 71c6a542d2..348d8bb0c6 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -6,7 +6,6 @@ #' as incidence rate. Primary analysis variable is the person-years at risk. #' #' @inheritParams argument_convention -#' @param id_var (`string`)\cr name of variable used as patient identifier if `"n_unique"` is included in `.stats`. #' @param control (`list`)\cr parameters for estimation details, specified by using #' the helper function [control_incidence_rate()]. Possible parameter options are: #' * `conf_level` (`proportion`)\cr confidence level for the estimated incidence rate. @@ -16,6 +15,8 @@ #' indicating time unit for data input. #' * `num_pt_year` (`numeric`)\cr time unit for desired output (in person-years). #' @param n_events (`string`)\cr name of integer variable indicating whether an event has been observed (1) or not (0). +#' @param id_var (`string`)\cr name of variable used as patient identifier if `"n_unique"` is included in `.stats`. +#' Defaults to `"USUBJID"`. #' @param .stats (`character`)\cr statistics to select for the table. Run `get_stats("estimate_incidence_rate")` #' to see available statistics for this function. #' @param summarize (`flag`)\cr whether the function should act as an analyze function (`summarize = FALSE`), or a @@ -93,8 +94,7 @@ s_incidence_rate <- function(df, ) } -#' @describeIn incidence_rate Formatted analysis function which is used as `afun` -#' in `estimate_incidence_rate()`. +#' @describeIn incidence_rate Formatted analysis function which is used as `afun` in `estimate_incidence_rate()`. #' #' @return #' * `a_incidence_rate()` returns the corresponding list with formatted [rtables::CellValue()]. diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index 3050766ace..a80c2c2ace 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -59,7 +59,8 @@ a_incidence_rate( \item{n_events}{(\code{string})\cr name of integer variable indicating whether an event has been observed (1) or not (0).} -\item{id_var}{(\code{string})\cr name of variable used as patient identifier if \code{"n_unique"} is included in \code{.stats}.} +\item{id_var}{(\code{string})\cr name of variable used as patient identifier if \code{"n_unique"} is included in \code{.stats}. +Defaults to \code{"USUBJID"}.} \item{control}{(\code{list})\cr parameters for estimation details, specified by using the helper function \code{\link[=control_incidence_rate]{control_incidence_rate()}}. Possible parameter options are: @@ -153,8 +154,7 @@ and additional format arguments. This function is a wrapper for \code{\link[rtab \item \code{s_incidence_rate()}: Statistics function which estimates the incidence rate and the associated confidence interval. -\item \code{a_incidence_rate()}: Formatted analysis function which is used as \code{afun} -in \code{estimate_incidence_rate()}. +\item \code{a_incidence_rate()}: Formatted analysis function which is used as \code{afun} in \code{estimate_incidence_rate()}. }} \examples{ From ea23ac93cc0b706577669c63cfa447b211acf778 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 18:32:07 -0400 Subject: [PATCH 04/12] Move helper functions to different file --- DESCRIPTION | 1 + R/h_incidence_rate.R | 135 +++++++++++++++++ R/incidence_rate.R | 136 ------------------ man/h_incidence_rate.Rd | 2 +- tests/testthat/_snaps/h_incidence_rate.md | 60 ++++++++ tests/testthat/_snaps/incidence_rate.md | 69 +++++++++ tests/testthat/test-h_incidence_rate.R | 38 +++++ ...incidence_rate.R => test-incidence_rate.R} | 39 ----- 8 files changed, 304 insertions(+), 176 deletions(-) create mode 100644 R/h_incidence_rate.R create mode 100644 tests/testthat/_snaps/h_incidence_rate.md create mode 100644 tests/testthat/_snaps/incidence_rate.md create mode 100644 tests/testthat/test-h_incidence_rate.R rename tests/testthat/{test-estimate_incidence_rate.R => test-incidence_rate.R} (65%) diff --git a/DESCRIPTION b/DESCRIPTION index cfe90fd1c6..c72c1c9923 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -130,6 +130,7 @@ Collate: 'h_adsl_adlb_merge_using_worst_flag.R' 'h_biomarkers_subgroups.R' 'h_cox_regression.R' + 'h_incidence_rate.R' 'h_km.R' 'h_logistic_regression.R' 'h_map_for_count_abnormal.R' diff --git a/R/h_incidence_rate.R b/R/h_incidence_rate.R new file mode 100644 index 0000000000..da5a0a5ea8 --- /dev/null +++ b/R/h_incidence_rate.R @@ -0,0 +1,135 @@ +#' Helper functions for incidence rate +#' +#' @description `r lifecycle::badge("stable")` +#' +#' @param control (`list`)\cr parameters for estimation details, specified by using +#' the helper function [control_incidence_rate()]. Possible parameter options are: +#' * `conf_level`: (`proportion`)\cr confidence level for the estimated incidence rate. +#' * `conf_type`: (`string`)\cr `normal` (default), `normal_log`, `exact`, or `byar` +#' for confidence interval type. +#' * `input_time_unit`: (`string`)\cr `day`, `week`, `month`, or `year` (default) +#' indicating time unit for data input. +#' * `num_pt_year`: (`numeric`)\cr time unit for desired output (in person-years). +#' @param person_years (`numeric(1)`)\cr total person-years at risk. +#' @param alpha (`numeric(1)`)\cr two-sided alpha-level for confidence interval. +#' @param n_events (`integer(1)`)\cr number of events observed. +#' +#' @return Estimated incidence rate, `rate`, and associated confidence interval, `rate_ci`. +#' +#' @seealso [incidence_rate] +#' +#' @name h_incidence_rate +NULL + +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated confidence interval based on the normal approximation for the +#' incidence rate. Unit is one person-year. +#' +#' @examples +#' h_incidence_rate_normal(200, 2) +#' +#' @export +h_incidence_rate_normal <- function(person_years, + n_events, + alpha = 0.05) { + checkmate::assert_number(person_years) + checkmate::assert_number(n_events) + assert_proportion_value(alpha) + + est <- n_events / person_years + se <- sqrt(est / person_years) + ci <- est + c(-1, 1) * stats::qnorm(1 - alpha / 2) * se + + list(rate = est, rate_ci = ci) +} + +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated confidence interval based on the normal approximation for the +#' logarithm of the incidence rate. Unit is one person-year. +#' +#' @examples +#' h_incidence_rate_normal_log(200, 2) +#' +#' @export +h_incidence_rate_normal_log <- function(person_years, + n_events, + alpha = 0.05) { + checkmate::assert_number(person_years) + checkmate::assert_number(n_events) + assert_proportion_value(alpha) + + rate_est <- n_events / person_years + rate_se <- sqrt(rate_est / person_years) + lrate_est <- log(rate_est) + lrate_se <- rate_se / rate_est + ci <- exp(lrate_est + c(-1, 1) * stats::qnorm(1 - alpha / 2) * lrate_se) + + list(rate = rate_est, rate_ci = ci) +} + +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated exact confidence interval. Unit is one person-year. +#' +#' @examples +#' h_incidence_rate_exact(200, 2) +#' +#' @export +h_incidence_rate_exact <- function(person_years, + n_events, + alpha = 0.05) { + checkmate::assert_number(person_years) + checkmate::assert_number(n_events) + assert_proportion_value(alpha) + + est <- n_events / person_years + lcl <- stats::qchisq(p = (alpha) / 2, df = 2 * n_events) / (2 * person_years) + ucl <- stats::qchisq(p = 1 - (alpha) / 2, df = 2 * n_events + 2) / (2 * person_years) + + list(rate = est, rate_ci = c(lcl, ucl)) +} + +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated Byar's confidence interval. Unit is one person-year. +#' +#' @examples +#' h_incidence_rate_byar(200, 2) +#' +#' @export +h_incidence_rate_byar <- function(person_years, + n_events, + alpha = 0.05) { + checkmate::assert_number(person_years) + checkmate::assert_number(n_events) + assert_proportion_value(alpha) + + est <- n_events / person_years + seg_1 <- n_events + 0.5 + seg_2 <- 1 - 1 / (9 * (n_events + 0.5)) + seg_3 <- stats::qnorm(1 - alpha / 2) * sqrt(1 / (n_events + 0.5)) / 3 + lcl <- seg_1 * ((seg_2 - seg_3)^3) / person_years + ucl <- seg_1 * ((seg_2 + seg_3) ^ 3) / person_years # styler: off + + list(rate = est, rate_ci = c(lcl, ucl)) +} + +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated confidence interval. +#' +#' @keywords internal +h_incidence_rate <- function(person_years, + n_events, + control = control_incidence_rate()) { + alpha <- 1 - control$conf_level + est <- switch(control$conf_type, + normal = h_incidence_rate_normal(person_years, n_events, alpha), + normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), + exact = h_incidence_rate_exact(person_years, n_events, alpha), + byar = h_incidence_rate_byar(person_years, n_events, alpha) + ) + + num_pt_year <- control$num_pt_year + list( + rate = est$rate * num_pt_year, + rate_ci = est$rate_ci * num_pt_year + ) +} diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 348d8bb0c6..ea234878f4 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -230,139 +230,3 @@ estimate_incidence_rate <- function(lyt, ) } } - -#' Helper functions for incidence rate -#' -#' @description `r lifecycle::badge("stable")` -#' -#' @param control (`list`)\cr parameters for estimation details, specified by using -#' the helper function [control_incidence_rate()]. Possible parameter options are: -#' * `conf_level`: (`proportion`)\cr confidence level for the estimated incidence rate. -#' * `conf_type`: (`string`)\cr `normal` (default), `normal_log`, `exact`, or `byar` -#' for confidence interval type. -#' * `input_time_unit`: (`string`)\cr `day`, `week`, `month`, or `year` (default) -#' indicating time unit for data input. -#' * `num_pt_year`: (`numeric`)\cr time unit for desired output (in person-years). -#' @param person_years (`numeric(1)`)\cr total person-years at risk. -#' @param alpha (`numeric(1)`)\cr two-sided alpha-level for confidence interval. -#' @param n_events (`integer(1)`)\cr number of events observed. -#' -#' @return Estimated incidence rate, `rate`, and associated confidence interval, `rate_ci`. -#' -#' @seealso [incidence_rate] -#' -#' @name h_incidence_rate -NULL - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated confidence interval based on the normal approximation for the -#' incidence rate. Unit is one person-year. -#' -#' @examples -#' h_incidence_rate_normal(200, 2) -#' -#' @export -h_incidence_rate_normal <- function(person_years, - n_events, - alpha = 0.05) { - checkmate::assert_number(person_years) - checkmate::assert_number(n_events) - assert_proportion_value(alpha) - - est <- n_events / person_years - se <- sqrt(est / person_years) - ci <- est + c(-1, 1) * stats::qnorm(1 - alpha / 2) * se - - list(rate = est, rate_ci = ci) -} - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated confidence interval based on the normal approximation for the -#' logarithm of the incidence rate. Unit is one person-year. -#' -#' @examples -#' h_incidence_rate_normal_log(200, 2) -#' -#' @export -h_incidence_rate_normal_log <- function(person_years, - n_events, - alpha = 0.05) { - checkmate::assert_number(person_years) - checkmate::assert_number(n_events) - assert_proportion_value(alpha) - - rate_est <- n_events / person_years - rate_se <- sqrt(rate_est / person_years) - lrate_est <- log(rate_est) - lrate_se <- rate_se / rate_est - ci <- exp(lrate_est + c(-1, 1) * stats::qnorm(1 - alpha / 2) * lrate_se) - - list(rate = rate_est, rate_ci = ci) -} - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated exact confidence interval. Unit is one person-year. -#' -#' @examples -#' h_incidence_rate_exact(200, 2) -#' -#' @export -h_incidence_rate_exact <- function(person_years, - n_events, - alpha = 0.05) { - checkmate::assert_number(person_years) - checkmate::assert_number(n_events) - assert_proportion_value(alpha) - - est <- n_events / person_years - lcl <- stats::qchisq(p = (alpha) / 2, df = 2 * n_events) / (2 * person_years) - ucl <- stats::qchisq(p = 1 - (alpha) / 2, df = 2 * n_events + 2) / (2 * person_years) - - list(rate = est, rate_ci = c(lcl, ucl)) -} - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated Byar's confidence interval. Unit is one person-year. -#' -#' @examples -#' h_incidence_rate_byar(200, 2) -#' -#' @export -h_incidence_rate_byar <- function(person_years, - n_events, - alpha = 0.05) { - checkmate::assert_number(person_years) - checkmate::assert_number(n_events) - assert_proportion_value(alpha) - - est <- n_events / person_years - seg_1 <- n_events + 0.5 - seg_2 <- 1 - 1 / (9 * (n_events + 0.5)) - seg_3 <- stats::qnorm(1 - alpha / 2) * sqrt(1 / (n_events + 0.5)) / 3 - lcl <- seg_1 * ((seg_2 - seg_3)^3) / person_years - ucl <- seg_1 * ((seg_2 + seg_3) ^ 3) / person_years # styler: off - - list(rate = est, rate_ci = c(lcl, ucl)) -} - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated confidence interval. -#' -#' @keywords internal -h_incidence_rate <- function(person_years, - n_events, - control = control_incidence_rate()) { - alpha <- 1 - control$conf_level - est <- switch(control$conf_type, - normal = h_incidence_rate_normal(person_years, n_events, alpha), - normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), - exact = h_incidence_rate_exact(person_years, n_events, alpha), - byar = h_incidence_rate_byar(person_years, n_events, alpha) - ) - - num_pt_year <- control$num_pt_year - list( - rate = est$rate * num_pt_year, - rate_ci = est$rate_ci * num_pt_year - ) -} diff --git a/man/h_incidence_rate.Rd b/man/h_incidence_rate.Rd index 08d974419d..0c896522e9 100644 --- a/man/h_incidence_rate.Rd +++ b/man/h_incidence_rate.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/incidence_rate.R +% Please edit documentation in R/h_incidence_rate.R \name{h_incidence_rate} \alias{h_incidence_rate} \alias{h_incidence_rate_normal} diff --git a/tests/testthat/_snaps/h_incidence_rate.md b/tests/testthat/_snaps/h_incidence_rate.md new file mode 100644 index 0000000000..4ab6826d47 --- /dev/null +++ b/tests/testthat/_snaps/h_incidence_rate.md @@ -0,0 +1,60 @@ +# h_incidence_rate_normal works as expected with healthy input + + Code + res + Output + $rate + [1] 0.01 + + $rate_ci + [1] -0.001630872 0.021630872 + + +# h_incidence_rate_normal_log works as expected with healthy input + + Code + res + Output + $rate + [1] 0.01 + + $rate_ci + [1] 0.003125199 0.031997963 + + +# h_incidence_rate_exact works as expected with healthy input + + Code + res + Output + $rate + [1] 0.01 + + $rate_ci + [1] 0.001776808 0.031478968 + + +# h_incidence_rate_byar works as expected with healthy input + + Code + res + Output + $rate + [1] 0.01 + + $rate_ci + [1] 0.002820411 0.027609866 + + +# h_incidence_rate works as expected with healthy input + + Code + res + Output + $rate + [1] 1 + + $rate_ci + [1] 0.3125199 3.1997963 + + diff --git a/tests/testthat/_snaps/incidence_rate.md b/tests/testthat/_snaps/incidence_rate.md new file mode 100644 index 0000000000..e7597f98fd --- /dev/null +++ b/tests/testthat/_snaps/incidence_rate.md @@ -0,0 +1,69 @@ +# control_incidence_rate works with customized parameters + + Code + res + Output + $conf_level + [1] 0.9 + + $conf_type + [1] "exact" + + $input_time_unit + [1] "month" + + $num_pt_year + [1] 100 + + +# s_incidence_rate works as expected with healthy input + + Code + res + Output + $person_years + [1] 9.058333 + attr(,"label") + [1] "Total patient-years at risk" + + $n_events + [1] 4 + attr(,"label") + [1] "Number of adverse events observed" + + $rate + [1] 44.15823 + attr(,"label") + [1] "AE rate per 100 patient-years" + + $rate_ci + [1] 19.40154 100.50487 + attr(,"label") + [1] "90% CI" + + $n_unique + [1] 4 + attr(,"label") + [1] "Total number of patients with at least one adverse event" + + $n_rate + [1] 4.00000 44.15823 + attr(,"label") + [1] "Number of adverse events observed (AE rate per 100 patient-years)" + + +# estimate_incidence_rate works as expected with healthy input + + Code + res + Output + A B + (N=3) (N=3) + ———————————————————————————————————————————————————————————————————————————————————————————————————— + Total patient-years at risk 3.8 5.2 + Number of adverse events observed 1 3 + AE rate per 100 patient-years 26.20 57.23 + 90% CI (5.06, 135.73) (22.14, 147.94) + Total number of patients with at least one adverse event 1 3 + Number of adverse events observed (AE rate per 100 patient-years) 1 (26.2) 3 (57.2) + diff --git a/tests/testthat/test-h_incidence_rate.R b/tests/testthat/test-h_incidence_rate.R new file mode 100644 index 0000000000..6b4603572c --- /dev/null +++ b/tests/testthat/test-h_incidence_rate.R @@ -0,0 +1,38 @@ +testthat::test_that("h_incidence_rate_normal works as expected with healthy input", { + result <- h_incidence_rate_normal(200, 2, 0.1) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("h_incidence_rate_normal_log works as expected with healthy input", { + result <- h_incidence_rate_normal_log(200, 2, 0.1) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("h_incidence_rate_exact works as expected with healthy input", { + result <- h_incidence_rate_exact(200, 2, 0.1) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("h_incidence_rate_byar works as expected with healthy input", { + result <- h_incidence_rate_byar(200, 2, 0.1) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("h_incidence_rate works as expected with healthy input", { + result <- h_incidence_rate( + 200, + 2, + control_incidence_rate(conf_level = 0.9, conf_type = "normal_log", num_pt_year = 100) + ) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) diff --git a/tests/testthat/test-estimate_incidence_rate.R b/tests/testthat/test-incidence_rate.R similarity index 65% rename from tests/testthat/test-estimate_incidence_rate.R rename to tests/testthat/test-incidence_rate.R index f83173e1a2..e584d2efe8 100644 --- a/tests/testthat/test-estimate_incidence_rate.R +++ b/tests/testthat/test-incidence_rate.R @@ -17,45 +17,6 @@ testthat::test_that("control_incidence_rate fails with wrong input", { testthat::expect_error(control_incidence_rate(num_pt_year = "one")) }) -testthat::test_that("h_incidence_rate_normal works as expected with healthy input", { - result <- h_incidence_rate_normal(200, 2, 0.1) - - res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) -}) - -testthat::test_that("h_incidence_rate_normal_log works as expected with healthy input", { - result <- h_incidence_rate_normal_log(200, 2, 0.1) - - res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) -}) - -testthat::test_that("h_incidence_rate_exact works as expected with healthy input", { - result <- h_incidence_rate_exact(200, 2, 0.1) - - res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) -}) - -testthat::test_that("h_incidence_rate_byar works as expected with healthy input", { - result <- h_incidence_rate_byar(200, 2, 0.1) - - res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) -}) - -testthat::test_that("h_incidence_rate works as expected with healthy input", { - result <- h_incidence_rate( - 200, - 2, - control_incidence_rate(conf_level = 0.9, conf_type = "normal_log", num_pt_year = 100) - ) - - res <- testthat::expect_silent(result) - testthat::expect_snapshot(res) -}) - testthat::test_that("s_incidence_rate works as expected with healthy input", { df <- data.frame( USUBJID = as.character(seq(6)), From b0f265ed0de47b170a586ebe5031c246a75f1573 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 18:56:22 -0400 Subject: [PATCH 05/12] Improve examples --- R/incidence_rate.R | 51 +++++++++++++++++++++++++++++++++---------- man/incidence_rate.Rd | 38 ++++++++++++++++++++++++++------ 2 files changed, 70 insertions(+), 19 deletions(-) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index ea234878f4..8d6bfa7aab 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -27,6 +27,16 @@ #' #' @seealso [control_incidence_rate()] and helper functions [h_incidence_rate]. #' +#' @examples +#' df <- data.frame( +#' USUBJID = as.character(seq(6)), +#' CNSR = c(0, 1, 1, 0, 0, 0), +#' AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), +#' ARM = factor(c("A", "A", "A", "B", "B", "B")), +#' STRATA1 = factor(c("X", "Y", "Y", "X", "X", "Y")) +#' ) +#' df$n_events <- 1 - df$CNSR +#' #' @name incidence_rate #' @order 1 NULL @@ -43,6 +53,13 @@ NULL #' - `n_unique`: Total number of patients with at least one event observed. #' - `n_rate`: Total number of events observed & estimated incidence rate. #' +#' @examples +#' s_incidence_rate( +#' df, +#' .var = "AVAL", +#' n_events = "n_events" +#' ) +#' #' @keywords internal s_incidence_rate <- function(df, .var, @@ -99,6 +116,14 @@ s_incidence_rate <- function(df, #' @return #' * `a_incidence_rate()` returns the corresponding list with formatted [rtables::CellValue()]. #' +#' @examples +#' a_incidence_rate( +#' df, +#' .var = "AVAL", +#' .df_row = df, +#' n_events = "n_events" +#' ) +#' #' @keywords internal a_incidence_rate <- function(df, labelstr = "", @@ -163,19 +188,8 @@ a_incidence_rate <- function(df, #' the statistics from `s_incidence_rate()` to the table layout. #' #' @examples -#' library(dplyr) -#' -#' df <- data.frame( -#' USUBJID = as.character(seq(6)), -#' CNSR = c(0, 1, 1, 0, 0, 0), -#' AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), -#' ARM = factor(c("A", "A", "A", "B", "B", "B")) -#' ) %>% -#' mutate(n_events = 1 - CNSR) -#' -#' basic_table() %>% +#' basic_table(show_colcounts = TRUE) %>% #' split_cols_by("ARM") %>% -#' add_colcounts() %>% #' estimate_incidence_rate( #' vars = "AVAL", #' n_events = "n_events", @@ -186,6 +200,19 @@ a_incidence_rate <- function(df, #' ) %>% #' build_table(df) #' +#' # summarize = TRUE +#' basic_table(show_colcounts = TRUE) %>% +#' split_cols_by("ARM") %>% +#' split_rows_by("STRATA1", child_labels = "visible") %>% +#' estimate_incidence_rate( +#' vars = "AVAL", +#' n_events = "n_events", +#' .stats = c("n_unique", "n_rate"), +#' summarize = TRUE, +#' label_fmt = "%.labels" +#' ) %>% +#' build_table(df) +#' #' @export #' @order 2 estimate_incidence_rate <- function(lyt, diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index a80c2c2ace..fd33fe561d 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -158,19 +158,17 @@ associated confidence interval. }} \examples{ -library(dplyr) - df <- data.frame( USUBJID = as.character(seq(6)), CNSR = c(0, 1, 1, 0, 0, 0), AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), - ARM = factor(c("A", "A", "A", "B", "B", "B")) -) \%>\% - mutate(n_events = 1 - CNSR) + ARM = factor(c("A", "A", "A", "B", "B", "B")), + STRATA1 = factor(c("X", "Y", "Y", "X", "X", "Y")) +) +df$n_events <- 1 - df$CNSR -basic_table() \%>\% +basic_table(show_colcounts = TRUE) \%>\% split_cols_by("ARM") \%>\% - add_colcounts() \%>\% estimate_incidence_rate( vars = "AVAL", n_events = "n_events", @@ -181,6 +179,32 @@ basic_table() \%>\% ) \%>\% build_table(df) +# summarize = TRUE +basic_table(show_colcounts = TRUE) \%>\% + split_cols_by("ARM") \%>\% + split_rows_by("STRATA1", child_labels = "visible") \%>\% + estimate_incidence_rate( + vars = "AVAL", + n_events = "n_events", + .stats = c("n_unique", "n_rate"), + summarize = TRUE, + label_fmt = "\%.labels" + ) \%>\% + build_table(df) + +s_incidence_rate( + df, + .var = "AVAL", + n_events = "n_events" +) + +a_incidence_rate( + df, + .var = "AVAL", + .df_row = df, + n_events = "n_events" +) + } \seealso{ \code{\link[=control_incidence_rate]{control_incidence_rate()}} and helper functions \link{h_incidence_rate}. From 602e5e4bccff8f98258badb7d60c0c4c29a2cb06 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 19:34:23 -0400 Subject: [PATCH 06/12] Update tests --- R/incidence_rate.R | 13 +-- man/incidence_rate.Rd | 2 +- tests/testthat/_snaps/incidence_rate.md | 100 ++++++++++++++++++-- tests/testthat/test-incidence_rate.R | 121 ++++++++++++++++++++---- 4 files changed, 198 insertions(+), 38 deletions(-) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 8d6bfa7aab..38f48848da 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -153,8 +153,12 @@ a_incidence_rate <- function(df, if (is.null(unlist(x_stats))) { return(NULL) } - if (is.null(.formats)) .formats <- formals()$.formats %>% eval() - if (is.null(.labels)) .labels <- sapply(x_stats, \(x) attributes(x)$label) + + # Fill in with defaults + formats_def <- formals()$.formats %>% eval() + .formats <- c(.formats, formats_def)[!duplicated(names(c(.formats, formats_def)))] + labels_def <- sapply(x_stats, \(x) attributes(x)$label) + .labels <- c(.labels, labels_def)[!duplicated(names(c(.labels, labels_def)))] if (nzchar(labelstr) > 0) { .labels <- sapply(.labels, \(x) gsub("%.labels", x, gsub("%s", labelstr, label_fmt))) } @@ -167,9 +171,6 @@ a_incidence_rate <- function(df, x_stats <- x_stats[.stats] - # Auto format handling - .formats <- apply_auto_formatting(.formats, x_stats, .df_row, .var) - in_rows( .list = x_stats, .formats = .formats, @@ -227,7 +228,7 @@ estimate_incidence_rate <- function(lyt, ..., show_labels = "hidden", table_names = vars, - .stats = NULL, + .stats = c("person_years", "n_events", "rate", "rate_ci"), .formats = NULL, .labels = NULL, .indent_mods = NULL) { diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index fd33fe561d..de7572d4fd 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -20,7 +20,7 @@ estimate_incidence_rate( ..., show_labels = "hidden", table_names = vars, - .stats = NULL, + .stats = c("person_years", "n_events", "rate", "rate_ci"), .formats = NULL, .labels = NULL, .indent_mods = NULL diff --git a/tests/testthat/_snaps/incidence_rate.md b/tests/testthat/_snaps/incidence_rate.md index e7597f98fd..44854523b7 100644 --- a/tests/testthat/_snaps/incidence_rate.md +++ b/tests/testthat/_snaps/incidence_rate.md @@ -52,18 +52,98 @@ [1] "Number of adverse events observed (AE rate per 100 patient-years)" -# estimate_incidence_rate works as expected with healthy input +# a_incidence_rate works with default arguments Code res Output - A B - (N=3) (N=3) - ———————————————————————————————————————————————————————————————————————————————————————————————————— - Total patient-years at risk 3.8 5.2 - Number of adverse events observed 1 3 - AE rate per 100 patient-years 26.20 57.23 - 90% CI (5.06, 135.73) (22.14, 147.94) - Total number of patients with at least one adverse event 1 3 - Number of adverse events observed (AE rate per 100 patient-years) 1 (26.2) 3 (57.2) + RowsVerticalSection (in_rows) object print method: + ---------------------------- + row_name formatted_cell indent_mod + 1 person_years 108.7 0 + 2 n_events 4 0 + 3 rate 3.68 0 + 4 rate_ci (0.07, 7.29) 0 + 5 n_unique 4 0 + 6 n_rate 4 (3.7) 0 + row_label + 1 Total patient-years at risk + 2 Number of adverse events observed + 3 AE rate per 100 patient-years + 4 95% CI + 5 Total number of patients with at least one adverse event + 6 Number of adverse events observed (AE rate per 100 patient-years) + +# a_incidence_rate works with customized arguments + + Code + res + Output + RowsVerticalSection (in_rows) object print method: + ---------------------------- + row_name formatted_cell indent_mod + 1 n_rate 4.00 (44.16) 3 + 2 n_unique 4.00 0 + row_label + 1 Total number of applicable adverse events (rate) + 2 Total number of patients with at least one adverse event + +# estimate_incidence_rate works as expected with default input + + Code + res + Output + A B + (N=3) (N=3) + —————————————————————————————————————————————————————————————————— + Total patient-years at risk 45.8 62.9 + Number of adverse events observed 1 3 + AE rate per 100 patient-years 2.18 4.77 + 95% CI (-2.10, 6.46) (-0.63, 10.17) + +# estimate_incidence_rate works as expected with custom input + + Code + res + Output + A B + (N=3) (N=3) + —————————————————————————————————————————————————————————————————————————————————————— + Total number of applicable adverse events (rate) 1.00 (26.20) 3.00 (57.23) + Total number of patients with at least one adverse event 1.00 3.00 + +# estimate_incidence_rate works with default arguments with summarize = TRUE + + Code + res + Output + A B + (N=3) (N=3) + ——————————————————————————————————————————————————————————————————————— + X - Total patient-years at risk 10.1 39.5 + X - Number of adverse events observed 1 2 + X - AE rate per 100 patient-years 9.90 5.06 + X - 95% CI (-9.50, 29.31) (-1.95, 12.08) + Y - Total patient-years at risk 35.7 23.4 + Y - Number of adverse events observed 0 1 + Y - AE rate per 100 patient-years 0.00 4.27 + Y - 95% CI (0.00, 0.00) (-4.10, 12.65) + +# estimate_incidence_rate works with custom arguments with summarize = TRUE + + Code + res + Output + A B + (N=3) (N=3) + ——————————————————————————————————————————————————————————————————————————————————————— + Total patient-years at risk 45.8 62.9 + Number of adverse events observed 1 3 + AE rate per 100 patient-years 2.18 4.77 + X + Total number of patients with at least one adverse event 1 2 + Number of adverse events observed (AE rate per 100 patient-years) 1 (9.9) 2 (5.1) + Y + Total number of patients with at least one adverse event 0 1 + Number of adverse events observed (AE rate per 100 patient-years) 0 (0.0) 1 (4.3) diff --git a/tests/testthat/test-incidence_rate.R b/tests/testthat/test-incidence_rate.R index e584d2efe8..045531971c 100644 --- a/tests/testthat/test-incidence_rate.R +++ b/tests/testthat/test-incidence_rate.R @@ -1,3 +1,12 @@ +df <- data.frame( + USUBJID = as.character(seq(6)), + CNSR = c(0, 1, 1, 0, 0, 0), + AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), + ARM = factor(c("A", "A", "A", "B", "B", "B")), + STRATA1 = factor(c("X", "Y", "Y", "X", "X", "Y")) +) %>% + dplyr::mutate(n_events = 1 - CNSR) + testthat::test_that("control_incidence_rate works with customized parameters", { result <- control_incidence_rate( conf_level = 0.9, @@ -18,14 +27,6 @@ testthat::test_that("control_incidence_rate fails with wrong input", { }) testthat::test_that("s_incidence_rate works as expected with healthy input", { - df <- data.frame( - USUBJID = as.character(seq(6)), - CNSR = c(0, 1, 1, 0, 0, 0), - AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), - ARM = factor(c("A", "A", "A", "B", "B", "B")) - ) %>% - dplyr::mutate(is_event = CNSR == 0) %>% - dplyr::mutate(n_events = as.integer(is_event)) result <- s_incidence_rate( df, .var = "AVAL", @@ -42,19 +43,56 @@ testthat::test_that("s_incidence_rate works as expected with healthy input", { testthat::expect_snapshot(res) }) -testthat::test_that("estimate_incidence_rate works as expected with healthy input", { - df <- data.frame( - USUBJID = as.character(seq(6)), - CNSR = c(0, 1, 1, 0, 0, 0), - AVAL = c(10.1, 20.4, 15.3, 20.8, 18.7, 23.4), - ARM = factor(c("A", "A", "A", "B", "B", "B")) - ) %>% - dplyr::mutate(is_event = CNSR == 0) %>% - dplyr::mutate(n_events = as.integer(is_event)) - - result <- basic_table() %>% +testthat::test_that("a_incidence_rate works with default arguments", { + result <- a_incidence_rate( + df, + .df_row = df, + .var = "AVAL", + n_events = "n_events" + ) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("a_incidence_rate works with customized arguments", { + result <- a_incidence_rate( + df, + .var = "AVAL", + n_events = "n_events", + control = control_incidence_rate( + conf_level = 0.9, + conf_type = "normal_log", + input_time_unit = "month", + num_pt_year = 100 + ), + .df_row = df, + .stats = c("n_rate", "n_unique"), + .formats = c(n_rate = "xx.xx (xx.xx)", n_unique = "xx.xx"), + .labels = c(n_rate = "Total number of applicable adverse events (rate)"), + .indent_mods = c(n_rate = 3L) + ) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("estimate_incidence_rate works as expected with default input", { + result <- basic_table(show_colcounts = TRUE) %>% + split_cols_by("ARM") %>% + estimate_incidence_rate( + vars = "AVAL", + n_events = "n_events" + ) %>% + build_table(df) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("estimate_incidence_rate works as expected with custom input", { + result <- basic_table(show_colcounts = TRUE) %>% split_cols_by("ARM") %>% - add_colcounts() %>% estimate_incidence_rate( vars = "AVAL", n_events = "n_events", @@ -63,7 +101,48 @@ testthat::test_that("estimate_incidence_rate works as expected with healthy inpu conf_type = "normal_log", input_time_unit = "month", num_pt_year = 100 - ) + ), + .stats = c("n_rate", "n_unique"), + .formats = c(n_rate = "xx.xx (xx.xx)", n_unique = "xx.xx"), + .labels = c(n_rate = "Total number of applicable adverse events (rate)"), + .indent_mods = c(n_rate = 3L) + ) %>% + build_table(df) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("estimate_incidence_rate works with default arguments with summarize = TRUE", { + result <- basic_table(show_colcounts = TRUE) %>% + split_cols_by("ARM") %>% + split_rows_by("STRATA1") %>% + estimate_incidence_rate( + vars = "AVAL", + n_events = "n_events", + summarize = TRUE + ) %>% + build_table(df) + + res <- testthat::expect_silent(result) + testthat::expect_snapshot(res) +}) + +testthat::test_that("estimate_incidence_rate works with custom arguments with summarize = TRUE", { + result <- basic_table(show_colcounts = TRUE) %>% + split_cols_by("ARM") %>% + estimate_incidence_rate( + vars = "AVAL", + n_events = "n_events", + .stats = c("person_years", "n_events", "rate") + ) %>% + split_rows_by("STRATA1", child_labels = "visible") %>% + estimate_incidence_rate( + vars = "AVAL", + n_events = "n_events", + .stats = c("n_unique", "n_rate"), + summarize = TRUE, + label_fmt = "%.labels" ) %>% build_table(df) From 0390cd243411d0b0e3afc19290e2ce6864837e61 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 19:41:41 -0400 Subject: [PATCH 07/12] Update NEWS --- NEWS.md | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/NEWS.md b/NEWS.md index e63b36871d..286976e098 100644 --- a/NEWS.md +++ b/NEWS.md @@ -4,6 +4,8 @@ * Reworking of `summarize_glm_count()` documentation and all its associated functions to better describe the results and the functions' purpose. * Added the `.formats` argument to `tabulate_rsp_subgroups` and `tabulate_survival_subgroups` to allow users to specify formats. * Added the `riskdiff` argument to `tabulate_rsp_subgroups` and `tabulate_survival_subgroups` to allow users to add a risk difference table column, and function `control_riskdiff` to specify settings for the risk difference column. +* Added `n_unique` statistic as a non-default option to `estimate_incidence_rate` which returns total number of patients with at least one event observed. +* Refactored `estimate_incidence_rate` to work as both an analyze function and a summarize function, controlled by the added `summarize` parameter. When `summarize = TRUE`, labels can be fine-tuned via the new `label_fmt` argument to the same function. ### Bug Fixes * Fixed a bug in `a_surv_time` that threw an error when split only has `"is_event"`. @@ -15,9 +17,9 @@ * Fixed bug in `decorate_grob` that did not handle well empty strings or `NULL` values for title and footers. * Fixed bug in `g_km` that caused an error when multiple records in the data had estimates at max time. - ### Miscellaneous * Began deprecation of the confusing functions `summary_formats` and `summary_labels`. +* Moved incidence rate helper functions into a separate `h_incidence_rate.R` file. # tern 0.9.5 From 036d63713f59e45814e8b268f6b192bee0242e46 Mon Sep 17 00:00:00 2001 From: github-actions <41898282+github-actions[bot]@users.noreply.github.com> Date: Wed, 4 Sep 2024 23:45:11 +0000 Subject: [PATCH 08/12] [skip style] [skip vbump] Restyle files --- R/h_incidence_rate.R | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/R/h_incidence_rate.R b/R/h_incidence_rate.R index da5a0a5ea8..26558d5fc0 100644 --- a/R/h_incidence_rate.R +++ b/R/h_incidence_rate.R @@ -121,10 +121,10 @@ h_incidence_rate <- function(person_years, control = control_incidence_rate()) { alpha <- 1 - control$conf_level est <- switch(control$conf_type, - normal = h_incidence_rate_normal(person_years, n_events, alpha), - normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), - exact = h_incidence_rate_exact(person_years, n_events, alpha), - byar = h_incidence_rate_byar(person_years, n_events, alpha) + normal = h_incidence_rate_normal(person_years, n_events, alpha), + normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), + exact = h_incidence_rate_exact(person_years, n_events, alpha), + byar = h_incidence_rate_byar(person_years, n_events, alpha) ) num_pt_year <- control$num_pt_year From 231d00c5d9f23d9344278a274201a72f088f4db0 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 19:52:19 -0400 Subject: [PATCH 09/12] Empty commit From a1ffc88f591785317f43f0bede7564acee177271 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 20:03:55 -0400 Subject: [PATCH 10/12] No examples for internal funcs --- R/incidence_rate.R | 7 ------- man/incidence_rate.Rd | 6 ------ 2 files changed, 13 deletions(-) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 38f48848da..7c62993d05 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -53,13 +53,6 @@ NULL #' - `n_unique`: Total number of patients with at least one event observed. #' - `n_rate`: Total number of events observed & estimated incidence rate. #' -#' @examples -#' s_incidence_rate( -#' df, -#' .var = "AVAL", -#' n_events = "n_events" -#' ) -#' #' @keywords internal s_incidence_rate <- function(df, .var, diff --git a/man/incidence_rate.Rd b/man/incidence_rate.Rd index de7572d4fd..3fd8a131c8 100644 --- a/man/incidence_rate.Rd +++ b/man/incidence_rate.Rd @@ -192,12 +192,6 @@ basic_table(show_colcounts = TRUE) \%>\% ) \%>\% build_table(df) -s_incidence_rate( - df, - .var = "AVAL", - n_events = "n_events" -) - a_incidence_rate( df, .var = "AVAL", From 2dc2405056b785bbef9b0922786e1215206a7c87 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Wed, 4 Sep 2024 20:16:51 -0400 Subject: [PATCH 11/12] Export a_incidence_rate --- NAMESPACE | 1 + R/incidence_rate.R | 2 +- 2 files changed, 2 insertions(+), 1 deletion(-) diff --git a/NAMESPACE b/NAMESPACE index f65577812e..e5fbc85f5c 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -28,6 +28,7 @@ export(a_count_patients_with_event) export(a_count_patients_with_flags) export(a_count_values) export(a_coxreg) +export(a_incidence_rate) export(a_length_proportion) export(a_odds_ratio) export(a_proportion) diff --git a/R/incidence_rate.R b/R/incidence_rate.R index 7c62993d05..af5de53112 100644 --- a/R/incidence_rate.R +++ b/R/incidence_rate.R @@ -117,7 +117,7 @@ s_incidence_rate <- function(df, #' n_events = "n_events" #' ) #' -#' @keywords internal +#' @export a_incidence_rate <- function(df, labelstr = "", .var, From 2445da82fa80cd6002a91b695b608bc1e68432b6 Mon Sep 17 00:00:00 2001 From: Emily de la Rua Date: Thu, 5 Sep 2024 16:55:50 -0400 Subject: [PATCH 12/12] Rearrange h_incidence_rate, styler --- R/h_incidence_rate.R | 46 ++++++++++++++++++++--------------------- man/h_incidence_rate.Rd | 22 ++++++++++---------- 2 files changed, 34 insertions(+), 34 deletions(-) diff --git a/R/h_incidence_rate.R b/R/h_incidence_rate.R index 26558d5fc0..2a1e88284a 100644 --- a/R/h_incidence_rate.R +++ b/R/h_incidence_rate.R @@ -21,6 +21,28 @@ #' @name h_incidence_rate NULL +#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and +#' associated confidence interval. +#' +#' @keywords internal +h_incidence_rate <- function(person_years, + n_events, + control = control_incidence_rate()) { + alpha <- 1 - control$conf_level + est <- switch(control$conf_type, + normal = h_incidence_rate_normal(person_years, n_events, alpha), + normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), + exact = h_incidence_rate_exact(person_years, n_events, alpha), + byar = h_incidence_rate_byar(person_years, n_events, alpha) + ) + + num_pt_year <- control$num_pt_year + list( + rate = est$rate * num_pt_year, + rate_ci = est$rate_ci * num_pt_year + ) +} + #' @describeIn h_incidence_rate Helper function to estimate the incidence rate and #' associated confidence interval based on the normal approximation for the #' incidence rate. Unit is one person-year. @@ -107,29 +129,7 @@ h_incidence_rate_byar <- function(person_years, seg_2 <- 1 - 1 / (9 * (n_events + 0.5)) seg_3 <- stats::qnorm(1 - alpha / 2) * sqrt(1 / (n_events + 0.5)) / 3 lcl <- seg_1 * ((seg_2 - seg_3)^3) / person_years - ucl <- seg_1 * ((seg_2 + seg_3) ^ 3) / person_years # styler: off + ucl <- seg_1 * ((seg_2 + seg_3)^3) / person_years list(rate = est, rate_ci = c(lcl, ucl)) } - -#' @describeIn h_incidence_rate Helper function to estimate the incidence rate and -#' associated confidence interval. -#' -#' @keywords internal -h_incidence_rate <- function(person_years, - n_events, - control = control_incidence_rate()) { - alpha <- 1 - control$conf_level - est <- switch(control$conf_type, - normal = h_incidence_rate_normal(person_years, n_events, alpha), - normal_log = h_incidence_rate_normal_log(person_years, n_events, alpha), - exact = h_incidence_rate_exact(person_years, n_events, alpha), - byar = h_incidence_rate_byar(person_years, n_events, alpha) - ) - - num_pt_year <- control$num_pt_year - list( - rate = est$rate * num_pt_year, - rate_ci = est$rate_ci * num_pt_year - ) -} diff --git a/man/h_incidence_rate.Rd b/man/h_incidence_rate.Rd index 1b558e981b..292c42dfdd 100644 --- a/man/h_incidence_rate.Rd +++ b/man/h_incidence_rate.Rd @@ -8,6 +8,8 @@ \alias{h_incidence_rate_byar} \title{Helper functions for incidence rate} \usage{ +h_incidence_rate(person_years, n_events, control = control_incidence_rate()) + h_incidence_rate_normal(person_years, n_events, alpha = 0.05) h_incidence_rate_normal_log(person_years, n_events, alpha = 0.05) @@ -15,26 +17,24 @@ h_incidence_rate_normal_log(person_years, n_events, alpha = 0.05) h_incidence_rate_exact(person_years, n_events, alpha = 0.05) h_incidence_rate_byar(person_years, n_events, alpha = 0.05) - -h_incidence_rate(person_years, n_events, control = control_incidence_rate()) } \arguments{ \item{person_years}{(\code{numeric(1)})\cr total person-years at risk.} \item{n_events}{(\code{integer(1)})\cr number of events observed.} -\item{alpha}{(\code{numeric(1)})\cr two-sided alpha-level for confidence interval.} - \item{control}{(\code{list})\cr parameters for estimation details, specified by using the helper function \code{\link[=control_incidence_rate]{control_incidence_rate()}}. Possible parameter options are: \itemize{ -\item \code{conf_level} (\code{proportion})\cr confidence level for the estimated incidence rate. -\item \code{conf_type} (\code{string})\cr \code{normal} (default), \code{normal_log}, \code{exact}, or \code{byar} +\item \code{conf_level}: (\code{proportion})\cr confidence level for the estimated incidence rate. +\item \code{conf_type}: (\code{string})\cr \code{normal} (default), \code{normal_log}, \code{exact}, or \code{byar} for confidence interval type. -\item \code{input_time_unit} (\code{string})\cr \code{day}, \code{week}, \code{month}, or \code{year} (default) +\item \code{input_time_unit}: (\code{string})\cr \code{day}, \code{week}, \code{month}, or \code{year} (default) indicating time unit for data input. -\item \code{num_pt_year} (\code{numeric})\cr time unit for desired output (in person-years). +\item \code{num_pt_year}: (\code{numeric})\cr time unit for desired output (in person-years). }} + +\item{alpha}{(\code{numeric(1)})\cr two-sided alpha-level for confidence interval.} } \value{ Estimated incidence rate, \code{rate}, and associated confidence interval, \code{rate_ci}. @@ -44,6 +44,9 @@ Estimated incidence rate, \code{rate}, and associated confidence interval, \code } \section{Functions}{ \itemize{ +\item \code{h_incidence_rate()}: Helper function to estimate the incidence rate and +associated confidence interval. + \item \code{h_incidence_rate_normal()}: Helper function to estimate the incidence rate and associated confidence interval based on the normal approximation for the incidence rate. Unit is one person-year. @@ -58,9 +61,6 @@ associated exact confidence interval. Unit is one person-year. \item \code{h_incidence_rate_byar()}: Helper function to estimate the incidence rate and associated Byar's confidence interval. Unit is one person-year. -\item \code{h_incidence_rate()}: Helper function to estimate the incidence rate and -associated confidence interval. - }} \examples{ h_incidence_rate_normal(200, 2)