diff --git a/NEWS.md b/NEWS.md index c4baf42a..c572080e 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,6 @@ # bayesplot (development version) +* Add extrapolation_factor parameter to `ppc_km_overlay()` and `ppc_km_overlay_grouped()` by @Sakuski * Add possibility for left-truncation to `ppc_km_overlay()` and `ppc_km_overlay_grouped()` by @Sakuski * Added `ppc_loo_pit_ecdf()` by @TeemuSailynoja diff --git a/R/ppc-censoring.R b/R/ppc-censoring.R index 7f39a78b..d11312ee 100644 --- a/R/ppc-censoring.R +++ b/R/ppc-censoring.R @@ -24,11 +24,12 @@ #' @section Plot Descriptions: #' \describe{ #' \item{`ppc_km_overlay()`}{ -#' Empirical CCDF estimates of each dataset (row) in `yrep` are overlaid, -#' with the Kaplan-Meier estimate (Kaplan and Meier, 1958) for `y` itself on -#' top (and in a darker shade). This is a PPC suitable for right-censored -#' `y`. Note that the replicated data from `yrep` is assumed to be -#' uncensored. +#' Empirical CCDF estimates of each dataset (row) in `yrep` are overlaid, with +#' the Kaplan-Meier estimate (Kaplan and Meier, 1958) for `y` itself on top +#' (and in a darker shade). This is a PPC suitable for right-censored `y`. +#' Note that the replicated data from `yrep` is assumed to be uncensored. Left +#' truncation (delayed entry) times for `y` can be specified using +#' `left_truncation_y`. #' } #' \item{`ppc_km_overlay_grouped()`}{ #' The same as `ppc_km_overlay()`, but with separate facets by `group`. @@ -40,24 +41,33 @@ #' @template reference-km #' #' @examples +#' \donttest{ #' color_scheme_set("brightblue") -#' y <- example_y_data() +#' #' # For illustrative purposes, (right-)censor values y > 110: +#' y <- example_y_data() #' status_y <- as.numeric(y <= 110) #' y <- pmin(y, 110) +#' #' # In reality, the replicated data (yrep) would be obtained from a #' # model which takes the censoring of y properly into account. Here, #' # for illustrative purposes, we simply use example_yrep_draws(): #' yrep <- example_yrep_draws() #' dim(yrep) -#' \donttest{ +#' +#' # Overlay 25 curves #' ppc_km_overlay(y, yrep[1:25, ], status_y = status_y) -#' } +#' +#' # With extrapolation_factor = 1 (no extrapolation) +#' ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, extrapolation_factor = 1) +#' +#' # With extrapolation_factor = Inf (show all posterior predictive draws) +#' ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, extrapolation_factor = Inf) +#' #' # With separate facets by group: #' group <- example_group_data() -#' \donttest{ #' ppc_km_overlay_grouped(y, yrep[1:25, ], group = group, status_y = status_y) -#' } +#' #' # With left-truncation (delayed entry) times: #' min_vals <- pmin(y, apply(yrep, 2, min)) #' left_truncation_y <- rep(0, length(y)) @@ -66,7 +76,6 @@ #' runif(sum(condition), min = 0.6, max = 0.99) * y[condition], #' min_vals[condition] - 0.001 #' ) -#' \donttest{ #' ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, #' left_truncation_y = left_truncation_y) #' } @@ -78,15 +87,23 @@ NULL #' be a numeric vector of the same length as `y` with values in \{0, 1\} (0 = #' right censored, 1 = event). #' @param left_truncation_y Optional parameter that specifies left-truncation -#' (delayed entry) times for the observations from `y`. This must -#' be a numeric vector of the same length as `y`. If `NULL` (default), -#' no left-truncation is assumed. +#' (delayed entry) times for the observations from `y`. This must be a numeric +#' vector of the same length as `y`. If `NULL` (default), no left-truncation +#' is assumed. +#' @param extrapolation_factor A numeric value (>=1) that controls how far the +#' plot is extended beyond the largest observed value in `y`. The default +#' value is 1.2, which corresponds to 20 % extrapolation. Note that all +#' posterior predictive draws may not be shown by default because of the +#' controlled extrapolation. To display all posterior predictive draws, set +#' `extrapolation_factor = Inf`. +#' ppc_km_overlay <- function( y, yrep, ..., status_y, left_truncation_y = NULL, + extrapolation_factor = 1.2, size = 0.25, alpha = 0.7 ) { @@ -97,15 +114,25 @@ ppc_km_overlay <- function( suggested_package("ggfortify") if (!is.numeric(status_y) || length(status_y) != length(y) || !all(status_y %in% c(0, 1))) { - stop("`status_y` must be a numeric vector of 0s and 1s the same length as `y`.") + stop("`status_y` must be a numeric vector of 0s and 1s the same length as `y`.", call. = FALSE) } if (!is.null(left_truncation_y)) { if (!is.numeric(left_truncation_y) || length(left_truncation_y) != length(y)) { - stop("`left_truncation_y` must be a numeric vector of the same length as `y`.") + stop("`left_truncation_y` must be a numeric vector of the same length as `y`.", call. = FALSE) } } + if (extrapolation_factor < 1) { + stop("`extrapolation_factor` must be greater than or equal to 1.", call. = FALSE) + } + if (extrapolation_factor == 1.2) { + message( + "Note: `extrapolation_factor` now defaults to 1.2 (20%).\n", + "To display all posterior predictive draws, set `extrapolation_factor = Inf`." + ) + } + data <- ppc_data(y, yrep, group = status_y) # Modify the status indicator: @@ -149,6 +176,10 @@ ppc_km_overlay <- function( fsf$is_y_size <- ifelse(fsf$is_y_color == "yrep", size, 1) fsf$is_y_alpha <- ifelse(fsf$is_y_color == "yrep", alpha, 1) + max_time_y <- max(y, na.rm = TRUE) + fsf <- fsf %>% + dplyr::filter(is_y_color != "yrep" | time <= max_time_y * extrapolation_factor) + # Ensure that the observed data gets plotted last by reordering the # levels of the factor "strata" fsf$strata <- factor(fsf$strata, levels = rev(levels(fsf$strata))) @@ -194,6 +225,7 @@ ppc_km_overlay_grouped <- function( ..., status_y, left_truncation_y = NULL, + extrapolation_factor = 1.2, size = 0.25, alpha = 0.7 ) { @@ -207,7 +239,8 @@ ppc_km_overlay_grouped <- function( status_y = status_y, left_truncation_y = left_truncation_y, size = size, - alpha = alpha + alpha = alpha, + extrapolation_factor = extrapolation_factor ) p_overlay + diff --git a/man/PPC-censoring.Rd b/man/PPC-censoring.Rd index 736a15a9..04224863 100644 --- a/man/PPC-censoring.Rd +++ b/man/PPC-censoring.Rd @@ -12,6 +12,7 @@ ppc_km_overlay( ..., status_y, left_truncation_y = NULL, + extrapolation_factor = 1.2, size = 0.25, alpha = 0.7 ) @@ -23,6 +24,7 @@ ppc_km_overlay_grouped( ..., status_y, left_truncation_y = NULL, + extrapolation_factor = 1.2, size = 0.25, alpha = 0.7 ) @@ -49,6 +51,13 @@ right censored, 1 = event).} be a numeric vector of the same length as \code{y}. If \code{NULL} (default), no left-truncation is assumed.} +\item{extrapolation_factor}{A numeric value (>=1) that controls how far the +plot is extended beyond the largest observed value in \code{y}. The default +value is 1.2, which corresponds to 20 \% extrapolation. Note that all +posterior predictive draws may not be shown by default because of +the controlled extrapolation. To display all posterior predictive draws, +set \code{extrapolation_factor = Inf}.} + \item{size, alpha}{Passed to the appropriate geom to control the appearance of the \code{yrep} distributions.} @@ -76,11 +85,12 @@ additional plots at \describe{ \item{\code{ppc_km_overlay()}}{ -Empirical CCDF estimates of each dataset (row) in \code{yrep} are overlaid, -with the Kaplan-Meier estimate (Kaplan and Meier, 1958) for \code{y} itself on -top (and in a darker shade). This is a PPC suitable for right-censored -\code{y}. Note that the replicated data from \code{yrep} is assumed to be -uncensored. +Empirical CCDF estimates of each dataset (row) in \code{yrep} are overlaid, with +the Kaplan-Meier estimate (Kaplan and Meier, 1958) for \code{y} itself on top +(and in a darker shade). This is a PPC suitable for right-censored \code{y}. +Note that the replicated data from \code{yrep} is assumed to be uncensored. Left +truncation (delayed entry) times for \code{y} can be specified using +\code{left_truncation_y}. } \item{\code{ppc_km_overlay_grouped()}}{ The same as \code{ppc_km_overlay()}, but with separate facets by \code{group}. @@ -89,24 +99,33 @@ The same as \code{ppc_km_overlay()}, but with separate facets by \code{group}. } \examples{ +\donttest{ color_scheme_set("brightblue") -y <- example_y_data() + # For illustrative purposes, (right-)censor values y > 110: +y <- example_y_data() status_y <- as.numeric(y <= 110) y <- pmin(y, 110) + # In reality, the replicated data (yrep) would be obtained from a # model which takes the censoring of y properly into account. Here, # for illustrative purposes, we simply use example_yrep_draws(): yrep <- example_yrep_draws() dim(yrep) -\donttest{ + +# Overlay 25 curves ppc_km_overlay(y, yrep[1:25, ], status_y = status_y) -} + +# With extrapolation_factor = 1 (no extrapolation) +ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, extrapolation_factor = 1) + +# With extrapolation_factor = Inf (show all posterior predictive draws) +ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, extrapolation_factor = Inf) + # With separate facets by group: group <- example_group_data() -\donttest{ ppc_km_overlay_grouped(y, yrep[1:25, ], group = group, status_y = status_y) -} + # With left-truncation (delayed entry) times: min_vals <- pmin(y, apply(yrep, 2, min)) left_truncation_y <- rep(0, length(y)) @@ -115,7 +134,6 @@ left_truncation_y[condition] <- pmin( runif(sum(condition), min = 0.6, max = 0.99) * y[condition], min_vals[condition] - 0.001 ) -\donttest{ ppc_km_overlay(y, yrep[1:25, ], status_y = status_y, left_truncation_y = left_truncation_y) } diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-default-2.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-default-2.svg index 792f234d..a04d8745 100644 --- a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-default-2.svg +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-default-2.svg @@ -25,17 +25,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -49,16 +49,18 @@ - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 y diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-default-2.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-default-2.svg index a24e432f..73ef1f39 100644 --- a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-default-2.svg +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-default-2.svg @@ -25,17 +25,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -50,17 +50,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -89,27 +89,31 @@ - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 0.0 0.5 diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-left-truncation-y.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-left-truncation-y.svg index 8b812eda..1d094bde 100644 --- a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-left-truncation-y.svg +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-left-truncation-y.svg @@ -25,17 +25,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -50,17 +50,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -89,27 +89,31 @@ - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 0.0 0.5 diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-max-extrapolation.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-max-extrapolation.svg new file mode 100644 index 00000000..ec169b9e --- /dev/null +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-max-extrapolation.svg @@ -0,0 +1,129 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +1 + + + + + + + + + +2 + + + + + + + + +0 +10 +20 +30 +40 + + + + + + +0 +10 +20 +30 +40 + +0.0 +0.5 +1.0 + + + + + +y +y +r +e +p +ppc_km_overlay_grouped (max extrapolation) + + diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-no-extrapolation.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-no-extrapolation.svg new file mode 100644 index 00000000..c3d814cd --- /dev/null +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-grouped-no-extrapolation.svg @@ -0,0 +1,129 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +1 + + + + + + + + + +2 + + + + + + + + +0 +5 +10 +15 +20 + + + + + + +0 +5 +10 +15 +20 + +0.0 +0.5 +1.0 + + + + + +y +y +r +e +p +ppc_km_overlay_grouped (no extrapolation) + + diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-left-truncation-y.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-left-truncation-y.svg index 9496f17d..7ad60bab 100644 --- a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-left-truncation-y.svg +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-left-truncation-y.svg @@ -25,17 +25,17 @@ - - - - - - - - - - - + + + + + + + + + + + @@ -49,16 +49,18 @@ - - - - - -0 -10 -20 -30 -40 + + + + + + +0 +5 +10 +15 +20 +25 y diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-max-extrapolation.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-max-extrapolation.svg new file mode 100644 index 00000000..84cf51e6 --- /dev/null +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-max-extrapolation.svg @@ -0,0 +1,71 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +0.0 +0.5 +1.0 + + + + + + + + + +0 +10 +20 +30 +40 + + +y +y +r +e +p +ppc_km_overlay (max extrapolation) + + diff --git a/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-no-extrapolation.svg b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-no-extrapolation.svg new file mode 100644 index 00000000..5b220bc5 --- /dev/null +++ b/tests/testthat/_snaps/ppc-censoring/ppc-km-overlay-no-extrapolation.svg @@ -0,0 +1,71 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +0.0 +0.5 +1.0 + + + + + + + + + +0 +5 +10 +15 +20 + + +y +y +r +e +p +ppc_km_overlay (no extrapolation) + + diff --git a/tests/testthat/test-ppc-censoring.R b/tests/testthat/test-ppc-censoring.R index 611d2c3f..d733c0ee 100644 --- a/tests/testthat/test-ppc-censoring.R +++ b/tests/testthat/test-ppc-censoring.R @@ -5,7 +5,8 @@ source(test_path("data-for-ppc-tests.R")) test_that("ppc_km_overlay returns a ggplot object", { skip_if_not_installed("ggfortify") - expect_gg(ppc_km_overlay(y, yrep, status_y = status_y, left_truncation_y = left_truncation_y, size = 0.5, alpha = 0.2)) + expect_gg(ppc_km_overlay(y, yrep, status_y = status_y, left_truncation_y = left_truncation_y, size = 0.5, alpha = 0.2, extrapolation_factor = Inf)) + expect_gg(ppc_km_overlay(y, yrep, status_y = status_y, left_truncation_y = left_truncation_y, size = 0.5, alpha = 0.2, extrapolation_factor = 1)) expect_gg(ppc_km_overlay(y2, yrep2, status_y = status_y2)) }) @@ -60,6 +61,22 @@ test_that("ppc_km_overlay errors if bad left_truncation_y value", { ) }) +test_that("ppc_km_overlay errors if bad extrapolation_factor value", { + skip_if_not_installed("ggfortify") + expect_error( + ppc_km_overlay(y, yrep, status_y = status_y, extrapolation_factor = 0.99), + "`extrapolation_factor` must be greater than or equal to 1." + ) +}) + +test_that("ppc_km_overlay messages if extrapolation_factor left at default value", { + skip_if_not_installed("ggfortify") + expect_message( + ppc_km_overlay(y, yrep, status_y = status_y), + "To display all posterior predictive draws, set `extrapolation_factor = Inf`.", + ) +}) + # Visual tests ----------------------------------------------------------------- test_that("ppc_km_overlay renders correctly", { @@ -82,13 +99,31 @@ test_that("ppc_km_overlay renders correctly", { p_base2 <- ppc_km_overlay(vdiff_y3, vdiff_yrep3, status_y = vdiff_status_y3) vdiffr::expect_doppelganger("ppc_km_overlay (default 2)", p_base2) - p_custom2 <- ppc_km_overlay( + p_custom2_left_truncation <- ppc_km_overlay( vdiff_y3, vdiff_yrep3, status_y = vdiff_status_y3, left_truncation_y = vdiff_left_truncation_y3) vdiffr::expect_doppelganger("ppc_km_overlay (left_truncation_y)", - p_custom2) + p_custom2_left_truncation) + + p_custom2_no_extrapolation <- ppc_km_overlay( + vdiff_y3, + vdiff_yrep3, + status_y = vdiff_status_y3, + extrapolation_factor = 1 + ) + vdiffr::expect_doppelganger("ppc_km_overlay (no extrapolation)", + p_custom2_no_extrapolation) + + p_custom2_max_extrapolation <- ppc_km_overlay( + vdiff_y3, + vdiff_yrep3, + status_y = vdiff_status_y3, + extrapolation_factor = Inf + ) + vdiffr::expect_doppelganger("ppc_km_overlay (max extrapolation)", + p_custom2_max_extrapolation) }) test_that("ppc_km_overlay_grouped renders correctly", { @@ -119,7 +154,7 @@ test_that("ppc_km_overlay_grouped renders correctly", { status_y = vdiff_status_y3) vdiffr::expect_doppelganger("ppc_km_overlay_grouped (default 2)", p_base2) - p_custom2 <- ppc_km_overlay_grouped( + p_custom2_left_truncation <- ppc_km_overlay_grouped( vdiff_y3, vdiff_yrep3, vdiff_group3, @@ -129,6 +164,32 @@ test_that("ppc_km_overlay_grouped renders correctly", { vdiffr::expect_doppelganger( "ppc_km_overlay_grouped (left_truncation_y)", - p_custom2 + p_custom2_left_truncation + ) + + p_custom2_no_extrapolation <- ppc_km_overlay_grouped( + vdiff_y3, + vdiff_yrep3, + vdiff_group3, + status_y = vdiff_status_y3, + extrapolation_factor = 1 + ) + + vdiffr::expect_doppelganger( + "ppc_km_overlay_grouped (no extrapolation)", + p_custom2_no_extrapolation + ) + + p_custom2_max_extrapolation <- ppc_km_overlay_grouped( + vdiff_y3, + vdiff_yrep3, + vdiff_group3, + status_y = vdiff_status_y3, + extrapolation_factor = Inf + ) + + vdiffr::expect_doppelganger( + "ppc_km_overlay_grouped (max extrapolation)", + p_custom2_max_extrapolation ) })