From 75d4da82a64ad53a7155a437b7936c5f4057b747 Mon Sep 17 00:00:00 2001 From: Sam Date: Tue, 19 Nov 2024 16:12:20 +0000 Subject: [PATCH 01/16] update gitignore to ignore local build paths for testing --- .gitignore | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/.gitignore b/.gitignore index 65d4fd1cc..b9c70bc28 100644 --- a/.gitignore +++ b/.gitignore @@ -25,3 +25,7 @@ data/models/*reparam docs /doc/ /Meta/ +vignettes/**_cache/ +*.pdf +.vscode/ +vignettes/figures/ From 871c82750dd890a15c926bbf560ea41f68e7b87a Mon Sep 17 00:00:00 2001 From: Sam Date: Tue, 19 Nov 2024 16:50:32 +0000 Subject: [PATCH 02/16] breaking rename --- R/latent_lognormal.R | 6 +++--- R/latent_model.R | 5 ++++- man/{log_lik_latent_lognormal.Rd => log_lik_latent.Rd} | 6 +++--- ..._epred_latent_lognormal.Rd => posterior_epred_latent.Rd} | 6 +++--- ...dict_latent_lognormal.Rd => posterior_predict_latent.Rd} | 6 +++--- 5 files changed, 16 insertions(+), 13 deletions(-) rename man/{log_lik_latent_lognormal.Rd => log_lik_latent.Rd} (82%) rename man/{posterior_epred_latent_lognormal.Rd => posterior_epred_latent.Rd} (78%) rename man/{posterior_predict_latent_lognormal.Rd => posterior_predict_latent.Rd} (79%) diff --git a/R/latent_lognormal.R b/R/latent_lognormal.R index a9de72c26..485576496 100644 --- a/R/latent_lognormal.R +++ b/R/latent_lognormal.R @@ -8,7 +8,7 @@ #' @param ... Additional arguments #' @autoglobal #' @keywords internal -posterior_predict_latent_lognormal <- function(i, prep, ...) { # nolint: object_length_linter +posterior_predict_latent <- function(i, prep, ...) { # nolint: object_length_linter mu <- brms::get_dpar(prep, "mu", i = i) sigma <- brms::get_dpar(prep, "sigma", i = i) @@ -42,7 +42,7 @@ posterior_predict_latent_lognormal <- function(i, prep, ...) { # nolint: object_ #' @param prep The result of a call to [`brms::prepare_predictions`] #' @autoglobal #' @keywords internal -posterior_epred_latent_lognormal <- function(prep) { # nolint: object_length_linter +posterior_epred_latent <- function(prep) { # nolint: object_length_linter mu <- brms::get_dpar(prep, "mu") sigma <- brms::get_dpar(prep, "sigma") exp(mu + sigma^2 / 2) @@ -56,7 +56,7 @@ posterior_epred_latent_lognormal <- function(prep) { # nolint: object_length_lin #' @param prep The result of a call to [brms::prepare_predictions()] #' @autoglobal #' @keywords internal -log_lik_latent_lognormal <- function(i, prep) { +log_lik_latent <- function(i, prep) { mu <- brms::get_dpar(prep, "mu", i = i) sigma <- brms::get_dpar(prep, "sigma", i = i) y <- prep$data$Y[i] diff --git a/R/latent_model.R b/R/latent_model.R index 99ccdd5c9..f2939d6aa 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -89,7 +89,10 @@ epidist_family_model.epidist_latent_model <- function( ub = c(NA, as.numeric(lapply(family$other_bounds, "[[", "ub"))), type = family$type, vars = c("pwindow", "swindow", "vreal1"), - loop = FALSE + loop = FALSE, + log_lik = log_lik_latent, + posterior_predict = posterior_predict_latent, + posterior_epred = posterior_epred_latent ) custom_family$reparm <- family$reparm return(custom_family) diff --git a/man/log_lik_latent_lognormal.Rd b/man/log_lik_latent.Rd similarity index 82% rename from man/log_lik_latent_lognormal.Rd rename to man/log_lik_latent.Rd index 4cac6907c..ecadb37a5 100644 --- a/man/log_lik_latent_lognormal.Rd +++ b/man/log_lik_latent.Rd @@ -1,10 +1,10 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_lognormal.R -\name{log_lik_latent_lognormal} -\alias{log_lik_latent_lognormal} +\name{log_lik_latent} +\alias{log_lik_latent} \title{Calculate the pointwise log likelihood of the \code{latent_gamma} family} \usage{ -log_lik_latent_lognormal(i, prep) +log_lik_latent(i, prep) } \arguments{ \item{i}{The index of the observation to calculate the log likelihood of} diff --git a/man/posterior_epred_latent_lognormal.Rd b/man/posterior_epred_latent.Rd similarity index 78% rename from man/posterior_epred_latent_lognormal.Rd rename to man/posterior_epred_latent.Rd index c8fb0d6ac..8bce50dd1 100644 --- a/man/posterior_epred_latent_lognormal.Rd +++ b/man/posterior_epred_latent.Rd @@ -1,11 +1,11 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_lognormal.R -\name{posterior_epred_latent_lognormal} -\alias{posterior_epred_latent_lognormal} +\name{posterior_epred_latent} +\alias{posterior_epred_latent} \title{Draws from the expected value of the posterior predictive distribution of the \code{latent_gamma} family} \usage{ -posterior_epred_latent_lognormal(prep) +posterior_epred_latent(prep) } \arguments{ \item{prep}{The result of a call to \code{\link[brms:prepare_predictions]{brms::prepare_predictions}}} diff --git a/man/posterior_predict_latent_lognormal.Rd b/man/posterior_predict_latent.Rd similarity index 79% rename from man/posterior_predict_latent_lognormal.Rd rename to man/posterior_predict_latent.Rd index bb5612305..997936ae2 100644 --- a/man/posterior_predict_latent_lognormal.Rd +++ b/man/posterior_predict_latent.Rd @@ -1,11 +1,11 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_lognormal.R -\name{posterior_predict_latent_lognormal} -\alias{posterior_predict_latent_lognormal} +\name{posterior_predict_latent} +\alias{posterior_predict_latent} \title{Draws from the posterior predictive distribution of the \code{latent_lognormal} family} \usage{ -posterior_predict_latent_lognormal(i, prep, ...) +posterior_predict_latent(i, prep, ...) } \arguments{ \item{i}{The index of the observation to predict} From d9936d4fc63b46a1bc074ff459f92325d643da93 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 12:13:01 +0000 Subject: [PATCH 03/16] start work on primarycensored port in and generalised predict and epred --- DESCRIPTION | 3 +- NAMESPACE | 2 + R/globals.R | 8 +- R/latent_gamma.R | 4 + R/latent_lognormal.R | 87 ------------- R/latent_model.R | 121 +++++++++++++++++- man/as_epidist_latent_model.Rd | 4 +- ...dist_latent_model.epidist_linelist_data.Rd | 4 +- ...idist_family_model.epidist_latent_model.Rd | 4 +- ...dist_formula_model.epidist_latent_model.Rd | 4 +- man/is_epidist_latent_model.Rd | 4 +- ...ik_latent.Rd => log_lik_epidist_latent.Rd} | 8 +- man/new_epidist_latent_model.Rd | 4 +- man/posterior_epred_epidist_latent.Rd | 39 ++++++ man/posterior_epred_latent.Rd | 16 --- man/posterior_predict_epidist_latent.Rd | 38 ++++++ man/posterior_predict_latent.Rd | 20 --- 17 files changed, 225 insertions(+), 145 deletions(-) rename man/{log_lik_latent.Rd => log_lik_epidist_latent.Rd} (74%) create mode 100644 man/posterior_epred_epidist_latent.Rd delete mode 100644 man/posterior_epred_latent.Rd create mode 100644 man/posterior_predict_epidist_latent.Rd delete mode 100644 man/posterior_predict_latent.Rd diff --git a/DESCRIPTION b/DESCRIPTION index d0b4f6e05..c7162d733 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -34,7 +34,8 @@ Imports: rstan (>= 2.26.0), dplyr, tibble, - lubridate + lubridate, + primarycensored Suggests: bookdown, testthat (>= 3.0.0), diff --git a/NAMESPACE b/NAMESPACE index 8fb43ebb5..6b87480fd 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -48,6 +48,8 @@ export(lognormal) export(new_epidist_latent_model) export(new_epidist_linelist_data) export(new_epidist_naive_model) +export(posterior_epred_epidist_latent) +export(posterior_predict_epidist_latent) export(predict_delay_parameters) export(predict_dpar) export(simulate_exponential_cases) diff --git a/R/globals.R b/R/globals.R index 41593e849..622145bf7 100644 --- a/R/globals.R +++ b/R/globals.R @@ -1,16 +1,10 @@ # Generated by roxyglobals: do not edit by hand utils::globalVariables(c( - ".data", # "samples", # - ".data", # + "family", # "woverlap", # - ".data", # - ".data", # - ".data", # "rlnorm", # - ".data", # - ".data", # <.replace_prior> "prior_new", # <.replace_prior> "source_new", # <.replace_prior> NULL diff --git a/R/latent_gamma.R b/R/latent_gamma.R index 46b7c2d87..698971b7b 100644 --- a/R/latent_gamma.R +++ b/R/latent_gamma.R @@ -27,6 +27,10 @@ posterior_predict_latent_gamma <- function(i, prep, ...) { # nolint: object_leng d_censored <- s_censored - p_censored } return(d_censored) + + .predict <- function(s) { + primarycensored::rpcens() + } } # Within brms this is a helper function called rblapply diff --git a/R/latent_lognormal.R b/R/latent_lognormal.R index 485576496..e69de29bb 100644 --- a/R/latent_lognormal.R +++ b/R/latent_lognormal.R @@ -1,87 +0,0 @@ -#' Draws from the posterior predictive distribution of the `latent_lognormal` -#' family -#' -#' See [brms::posterior_predict()]. -#' -#' @param i The index of the observation to predict -#' @param prep The result of a call to [brms::posterior_predict()] -#' @param ... Additional arguments -#' @autoglobal -#' @keywords internal -posterior_predict_latent <- function(i, prep, ...) { # nolint: object_length_linter - mu <- brms::get_dpar(prep, "mu", i = i) - sigma <- brms::get_dpar(prep, "sigma", i = i) - - relative_obs_time <- prep$data$vreal1[i] - pwindow <- prep$data$vreal2[i] - swindow <- prep$data$vreal3[i] - - .predict <- function(s) { - d_censored <- relative_obs_time + 1 - # while loop to impose the truncation - while (d_censored > relative_obs_time) { - p_latent <- stats::runif(1, 0, 1) * pwindow - d_latent <- stats::rlnorm(1, meanlog = mu[s], sdlog = sigma[s]) - s_latent <- p_latent + d_latent - p_censored <- .floor_mult(p_latent, pwindow) - s_censored <- .floor_mult(s_latent, swindow) - d_censored <- s_censored - p_censored - } - return(d_censored) - } - - # Within brms this is a helper function called rblapply - do.call(rbind, lapply(seq_len(prep$ndraws), .predict)) -} - -#' Draws from the expected value of the posterior predictive distribution of the -#' `latent_gamma` family -#' -#' See [brms::posterior_epred()]. -#' -#' @param prep The result of a call to [`brms::prepare_predictions`] -#' @autoglobal -#' @keywords internal -posterior_epred_latent <- function(prep) { # nolint: object_length_linter - mu <- brms::get_dpar(prep, "mu") - sigma <- brms::get_dpar(prep, "sigma") - exp(mu + sigma^2 / 2) -} - -#' Calculate the pointwise log likelihood of the `latent_gamma` family -#' -#' See [brms::log_lik()]. -#' -#' @param i The index of the observation to calculate the log likelihood of -#' @param prep The result of a call to [brms::prepare_predictions()] -#' @autoglobal -#' @keywords internal -log_lik_latent <- function(i, prep) { - mu <- brms::get_dpar(prep, "mu", i = i) - sigma <- brms::get_dpar(prep, "sigma", i = i) - y <- prep$data$Y[i] - relative_obs_time <- prep$data$vreal1[i] - pwindow <- prep$data$vreal2[i] - swindow <- prep$data$vreal3[i] - - # Generates values of the swindow_raw and pwindow_raw, but really these should - # be extracted from prep or the fitted raws somehow. See: - # https://github.com/epinowcast/epidist/issues/267 - swindow_raw <- stats::runif(prep$ndraws) - pwindow_raw <- stats::runif(prep$ndraws) - - swindow <- swindow_raw * swindow - - # For no overlap calculate as usual, for overlap ensure pwindow < swindow - if (i %in% prep$data$noverlap) { - pwindow <- pwindow_raw * pwindow - } else { - pwindow <- pwindow_raw * swindow - } - - d <- y - pwindow + swindow - obs_time <- relative_obs_time - pwindow - lpdf <- stats::dlnorm(d, meanlog = mu, sdlog = sigma, log = TRUE) - lcdf <- stats::plnorm(obs_time, meanlog = mu, sdlog = sigma, log.p = TRUE) - return(lpdf - lcdf) -} diff --git a/R/latent_model.R b/R/latent_model.R index f2939d6aa..6b92eed70 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -82,7 +82,7 @@ epidist_family_model.epidist_latent_model <- function( data, family, ...) { # Really the name and vars are the "model-specific" parts here custom_family <- brms::custom_family( - paste0("latent_", family$family), + "epidist_latent", dpars = family$dpars, links = c(family$link, family$other_links), lb = c(NA, as.numeric(lapply(family$other_bounds, "[[", "lb"))), @@ -90,14 +90,127 @@ epidist_family_model.epidist_latent_model <- function( type = family$type, vars = c("pwindow", "swindow", "vreal1"), loop = FALSE, - log_lik = log_lik_latent, - posterior_predict = posterior_predict_latent, - posterior_epred = posterior_epred_latent + log_lik = log_lik_epidist_latent(family), + posterior_predict = posterior_predict_epidist_latent(family), + posterior_epred = posterior_epred_epidist_latent(family) ) custom_family$reparm <- family$reparm return(custom_family) } +#' Calculate the pointwise log likelihood of the `latent_gamma` family +#' +#' See [brms::log_lik()]. +#' +#' @param i The index of the observation to calculate the log likelihood of +#' @param prep The result of a call to [brms::prepare_predictions()] +#' @autoglobal +#' @keywords internal +log_lik_epidist_latent <- function(i, prep) { + mu <- brms::get_dpar(prep, "mu", i = i) + sigma <- brms::get_dpar(prep, "sigma", i = i) + y <- prep$data$Y[i] + relative_obs_time <- prep$data$vreal1[i] + pwindow <- prep$data$vreal2[i] + swindow <- prep$data$vreal3[i] + + # Generates values of the swindow_raw and pwindow_raw, but really these should + # be extracted from prep or the fitted raws somehow. See: + # https://github.com/epinowcast/epidist/issues/267 + swindow_raw <- stats::runif(prep$ndraws) + pwindow_raw <- stats::runif(prep$ndraws) + + swindow <- swindow_raw * swindow + + # For no overlap calculate as usual, for overlap ensure pwindow < swindow + if (i %in% prep$data$noverlap) { + pwindow <- pwindow_raw * pwindow + } else { + pwindow <- pwindow_raw * swindow + } + + d <- y - pwindow + swindow + obs_time <- relative_obs_time - pwindow + lpdf <- stats::dlnorm(d, meanlog = mu, sdlog = sigma, log = TRUE) + lcdf <- stats::plnorm(obs_time, meanlog = mu, sdlog = sigma, log.p = TRUE) + return(lpdf - lcdf) +} + +#' Create a function to draw from the posterior predictive distribution for a +#' latent model +#' +#' This function creates a function that draws from the posterior predictive +#' distribution for a latent model using [primarycensored::rpcens()] to handle +#' censoring and truncation. The returned function takes a prep argument from +#' `brms` and returns posterior predictions. This is used internally by +#' [brms::posterior_predict()] to generate predictions for latent models. +#' +#' @inheritParams epidist_family_model +#' +#' @return A function that takes a prep argument from brms and returns a matrix +#' of posterior predictions, with one row per posterior draw and one column +#' per observation. +#' +#' @seealso [brms::posterior_predict()] for details on how this is used within +#' `brms`, [primarycensored::rpcens()] for details on the censoring approach +#' @autoglobal +#' @family latent_model +#' @export +posterior_predict_epidist_latent <- function(i, prep, ...) { + fn_string <- paste0("posterior_predict_", family$name) + dist_fn <- eval(parse(text = fn_string)) + + rdist <- function(n, i, prep, ...) { + purrr::map_dbl( + seq_len(n), + ~ do.call(dist_fn, list(i = i, prep = prep)) + ) + } + + .predict <- function(i, prep) { + relative_obs_time <- prep$data$vreal1[i] + pwindow <- prep$data$vreal2[i] + swindow <- prep$data$vreal3[i] + + primarycensored::rpcens( + n = prep$ndraws, + rdist = rdist, + rprimary = stats::runif, + pwindow = prep$data$vreal2[i], + swindow = prep$data$vreal3[i], + D = prep$data$vreal1[i], + i = i, + prep = prep + ) + } + return(.predict) +} + +#' Create a function to draw from the expected value of the posterior predictive +#' distribution for a latent model +#' +#' This function creates a function that calculates the expected value of the +#' posterior predictive distribution for a latent model. The returned function +#' takes a prep argument (from brms) and returns posterior expected values. +#' This is used internally by [brms::posterior_epred()] to calculate expected +#' values for latent models. +#' +#' @inheritParams epidist_family_model +#' +#' @return A function that takes a prep argument from brms and returns a matrix +#' of posterior expected values, with one row per posterior draw and one column +#' per observation. +#' +#' @seealso [brms::posterior_epred()] for details on how this is used within +#' `brms`. +#' @autoglobal +#' @family latent_model +#' @export +posterior_epred_epidist_latent <- function(family) { + fn_string <- paste0("posterior_epred_", family$name) + eval(parse(text = fn_string)) +} + #' Define the model-specific component of an `epidist` custom formula #' #' @inheritParams epidist_formula_model diff --git a/man/as_epidist_latent_model.Rd b/man/as_epidist_latent_model.Rd index fb2329eda..2c8eb8378 100644 --- a/man/as_epidist_latent_model.Rd +++ b/man/as_epidist_latent_model.Rd @@ -18,6 +18,8 @@ Other latent_model: \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()} +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/as_epidist_latent_model.epidist_linelist_data.Rd b/man/as_epidist_latent_model.epidist_linelist_data.Rd index ca3c6acf0..dc4951d33 100644 --- a/man/as_epidist_latent_model.epidist_linelist_data.Rd +++ b/man/as_epidist_latent_model.epidist_linelist_data.Rd @@ -18,6 +18,8 @@ Other latent_model: \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()} +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/epidist_family_model.epidist_latent_model.Rd b/man/epidist_family_model.epidist_latent_model.Rd index ae3dae26b..dd5defc9c 100644 --- a/man/epidist_family_model.epidist_latent_model.Rd +++ b/man/epidist_family_model.epidist_latent_model.Rd @@ -23,6 +23,8 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()} +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/epidist_formula_model.epidist_latent_model.Rd b/man/epidist_formula_model.epidist_latent_model.Rd index 5811233df..df1298052 100644 --- a/man/epidist_formula_model.epidist_latent_model.Rd +++ b/man/epidist_formula_model.epidist_latent_model.Rd @@ -26,6 +26,8 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()} +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/is_epidist_latent_model.Rd b/man/is_epidist_latent_model.Rd index b9c4115ae..76671a95f 100644 --- a/man/is_epidist_latent_model.Rd +++ b/man/is_epidist_latent_model.Rd @@ -18,6 +18,8 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()} +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/log_lik_latent.Rd b/man/log_lik_epidist_latent.Rd similarity index 74% rename from man/log_lik_latent.Rd rename to man/log_lik_epidist_latent.Rd index ecadb37a5..755f0dfe7 100644 --- a/man/log_lik_latent.Rd +++ b/man/log_lik_epidist_latent.Rd @@ -1,10 +1,10 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_lognormal.R -\name{log_lik_latent} -\alias{log_lik_latent} +% Please edit documentation in R/latent_model.R +\name{log_lik_epidist_latent} +\alias{log_lik_epidist_latent} \title{Calculate the pointwise log likelihood of the \code{latent_gamma} family} \usage{ -log_lik_latent(i, prep) +log_lik_epidist_latent(i, prep) } \arguments{ \item{i}{The index of the observation to calculate the log likelihood of} diff --git a/man/new_epidist_latent_model.Rd b/man/new_epidist_latent_model.Rd index 02841df23..bfb7cfba6 100644 --- a/man/new_epidist_latent_model.Rd +++ b/man/new_epidist_latent_model.Rd @@ -21,6 +21,8 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{is_epidist_latent_model}()} +\code{\link{is_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()}, +\code{\link{posterior_predict_epidist_latent}()} } \concept{latent_model} diff --git a/man/posterior_epred_epidist_latent.Rd b/man/posterior_epred_epidist_latent.Rd new file mode 100644 index 000000000..930f072b6 --- /dev/null +++ b/man/posterior_epred_epidist_latent.Rd @@ -0,0 +1,39 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/latent_model.R +\name{posterior_epred_epidist_latent} +\alias{posterior_epred_epidist_latent} +\title{Create a function to draw from the expected value of the posterior predictive +distribution for a latent model} +\usage{ +posterior_epred_epidist_latent(family) +} +\arguments{ +\item{family}{Output of a call to \code{brms::brmsfamily()} with additional +information as provided by \code{.add_dpar_info()}} +} +\value{ +A function that takes a prep argument from brms and returns a matrix +of posterior expected values, with one row per posterior draw and one column +per observation. +} +\description{ +This function creates a function that calculates the expected value of the +posterior predictive distribution for a latent model. The returned function +takes a prep argument (from brms) and returns posterior expected values. +This is used internally by \code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}} to calculate expected +values for latent models. +} +\seealso{ +\code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}} for details on how this is used within +\code{brms}. + +Other latent_model: +\code{\link{as_epidist_latent_model}()}, +\code{\link{as_epidist_latent_model.epidist_linelist_data}()}, +\code{\link{epidist_family_model.epidist_latent_model}()}, +\code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{is_epidist_latent_model}()}, +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_predict_epidist_latent}()} +} +\concept{latent_model} diff --git a/man/posterior_epred_latent.Rd b/man/posterior_epred_latent.Rd deleted file mode 100644 index 8bce50dd1..000000000 --- a/man/posterior_epred_latent.Rd +++ /dev/null @@ -1,16 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_lognormal.R -\name{posterior_epred_latent} -\alias{posterior_epred_latent} -\title{Draws from the expected value of the posterior predictive distribution of the -\code{latent_gamma} family} -\usage{ -posterior_epred_latent(prep) -} -\arguments{ -\item{prep}{The result of a call to \code{\link[brms:prepare_predictions]{brms::prepare_predictions}}} -} -\description{ -See \code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}}. -} -\keyword{internal} diff --git a/man/posterior_predict_epidist_latent.Rd b/man/posterior_predict_epidist_latent.Rd new file mode 100644 index 000000000..fb4d45190 --- /dev/null +++ b/man/posterior_predict_epidist_latent.Rd @@ -0,0 +1,38 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/latent_model.R +\name{posterior_predict_epidist_latent} +\alias{posterior_predict_epidist_latent} +\title{Create a function to draw from the posterior predictive distribution for a +latent model} +\usage{ +posterior_predict_epidist_latent(i, prep, ...) +} +\arguments{ +\item{...}{Additional arguments passed to \code{fn} method.} +} +\value{ +A function that takes a prep argument from brms and returns a matrix +of posterior predictions, with one row per posterior draw and one column +per observation. +} +\description{ +This function creates a function that draws from the posterior predictive +distribution for a latent model using \code{\link[primarycensored:rprimarycensored]{primarycensored::rpcens()}} to handle +censoring and truncation. The returned function takes a prep argument from +\code{brms} and returns posterior predictions. This is used internally by +\code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}} to generate predictions for latent models. +} +\seealso{ +\code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}} for details on how this is used within +\code{brms}, \code{\link[primarycensored:rprimarycensored]{primarycensored::rpcens()}} for details on the censoring approach + +Other latent_model: +\code{\link{as_epidist_latent_model}()}, +\code{\link{as_epidist_latent_model.epidist_linelist_data}()}, +\code{\link{epidist_family_model.epidist_latent_model}()}, +\code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{is_epidist_latent_model}()}, +\code{\link{new_epidist_latent_model}()}, +\code{\link{posterior_epred_epidist_latent}()} +} +\concept{latent_model} diff --git a/man/posterior_predict_latent.Rd b/man/posterior_predict_latent.Rd deleted file mode 100644 index 997936ae2..000000000 --- a/man/posterior_predict_latent.Rd +++ /dev/null @@ -1,20 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_lognormal.R -\name{posterior_predict_latent} -\alias{posterior_predict_latent} -\title{Draws from the posterior predictive distribution of the \code{latent_lognormal} -family} -\usage{ -posterior_predict_latent(i, prep, ...) -} -\arguments{ -\item{i}{The index of the observation to predict} - -\item{prep}{The result of a call to \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}}} - -\item{...}{Additional arguments} -} -\description{ -See \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}}. -} -\keyword{internal} From ae02c94814e91d9dbb72c5f75a3feafbafeb0e4c Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 12:22:07 +0000 Subject: [PATCH 04/16] check new implementatons against ebola vignette --- NAMESPACE | 4 ++-- R/globals.R | 1 - R/latent_model.R | 20 +++++++++---------- man/as_epidist_latent_model.Rd | 4 ++-- ...dist_latent_model.epidist_linelist_data.Rd | 4 ++-- ...idist_family_model.epidist_latent_model.Rd | 4 ++-- ...dist_formula_model.epidist_latent_model.Rd | 4 ++-- man/is_epidist_latent_model.Rd | 4 ++-- ...ik_epidist_latent.Rd => log_lik_latent.Rd} | 6 +++--- man/new_epidist_latent_model.Rd | 4 ++-- ...st_latent.Rd => posterior_epred_latent.Rd} | 8 ++++---- ..._latent.Rd => posterior_predict_latent.Rd} | 11 +++++----- 12 files changed, 37 insertions(+), 37 deletions(-) rename man/{log_lik_epidist_latent.Rd => log_lik_latent.Rd} (83%) rename man/{posterior_epred_epidist_latent.Rd => posterior_epred_latent.Rd} (89%) rename man/{posterior_predict_epidist_latent.Rd => posterior_predict_latent.Rd} (84%) diff --git a/NAMESPACE b/NAMESPACE index 6b87480fd..2266b4020 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -48,8 +48,8 @@ export(lognormal) export(new_epidist_latent_model) export(new_epidist_linelist_data) export(new_epidist_naive_model) -export(posterior_epred_epidist_latent) -export(posterior_predict_epidist_latent) +export(posterior_epred_latent) +export(posterior_predict_latent) export(predict_delay_parameters) export(predict_dpar) export(simulate_exponential_cases) diff --git a/R/globals.R b/R/globals.R index 622145bf7..bfa8d919d 100644 --- a/R/globals.R +++ b/R/globals.R @@ -2,7 +2,6 @@ utils::globalVariables(c( "samples", # - "family", # "woverlap", # "rlnorm", # "prior_new", # <.replace_prior> diff --git a/R/latent_model.R b/R/latent_model.R index 6b92eed70..883ae0929 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -82,7 +82,7 @@ epidist_family_model.epidist_latent_model <- function( data, family, ...) { # Really the name and vars are the "model-specific" parts here custom_family <- brms::custom_family( - "epidist_latent", + paste0("latent_", family$family), dpars = family$dpars, links = c(family$link, family$other_links), lb = c(NA, as.numeric(lapply(family$other_bounds, "[[", "lb"))), @@ -90,9 +90,9 @@ epidist_family_model.epidist_latent_model <- function( type = family$type, vars = c("pwindow", "swindow", "vreal1"), loop = FALSE, - log_lik = log_lik_epidist_latent(family), - posterior_predict = posterior_predict_epidist_latent(family), - posterior_epred = posterior_epred_epidist_latent(family) + log_lik = log_lik_latent, + posterior_predict = posterior_predict_latent(family), + posterior_epred = posterior_epred_latent(family) ) custom_family$reparm <- family$reparm return(custom_family) @@ -106,7 +106,7 @@ epidist_family_model.epidist_latent_model <- function( #' @param prep The result of a call to [brms::prepare_predictions()] #' @autoglobal #' @keywords internal -log_lik_epidist_latent <- function(i, prep) { +log_lik_latent <- function(i, prep) { mu <- brms::get_dpar(prep, "mu", i = i) sigma <- brms::get_dpar(prep, "sigma", i = i) y <- prep$data$Y[i] @@ -156,8 +156,8 @@ log_lik_epidist_latent <- function(i, prep) { #' @autoglobal #' @family latent_model #' @export -posterior_predict_epidist_latent <- function(i, prep, ...) { - fn_string <- paste0("posterior_predict_", family$name) +posterior_predict_latent <- function(family) { + fn_string <- paste0("brms:::posterior_predict_", family$family) dist_fn <- eval(parse(text = fn_string)) rdist <- function(n, i, prep, ...) { @@ -167,7 +167,7 @@ posterior_predict_epidist_latent <- function(i, prep, ...) { ) } - .predict <- function(i, prep) { + .predict <- function(i, prep, ...) { relative_obs_time <- prep$data$vreal1[i] pwindow <- prep$data$vreal2[i] swindow <- prep$data$vreal3[i] @@ -206,8 +206,8 @@ posterior_predict_epidist_latent <- function(i, prep, ...) { #' @autoglobal #' @family latent_model #' @export -posterior_epred_epidist_latent <- function(family) { - fn_string <- paste0("posterior_epred_", family$name) +posterior_epred_latent <- function(family) { + fn_string <- paste0("brms:::posterior_epred_", family$family) eval(parse(text = fn_string)) } diff --git a/man/as_epidist_latent_model.Rd b/man/as_epidist_latent_model.Rd index 2c8eb8378..0a600d099 100644 --- a/man/as_epidist_latent_model.Rd +++ b/man/as_epidist_latent_model.Rd @@ -19,7 +19,7 @@ Other latent_model: \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/as_epidist_latent_model.epidist_linelist_data.Rd b/man/as_epidist_latent_model.epidist_linelist_data.Rd index dc4951d33..d9620e5fb 100644 --- a/man/as_epidist_latent_model.epidist_linelist_data.Rd +++ b/man/as_epidist_latent_model.epidist_linelist_data.Rd @@ -19,7 +19,7 @@ Other latent_model: \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/epidist_family_model.epidist_latent_model.Rd b/man/epidist_family_model.epidist_latent_model.Rd index dd5defc9c..61a7fad8d 100644 --- a/man/epidist_family_model.epidist_latent_model.Rd +++ b/man/epidist_family_model.epidist_latent_model.Rd @@ -24,7 +24,7 @@ Other latent_model: \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/epidist_formula_model.epidist_latent_model.Rd b/man/epidist_formula_model.epidist_latent_model.Rd index df1298052..dbd10c3b9 100644 --- a/man/epidist_formula_model.epidist_latent_model.Rd +++ b/man/epidist_formula_model.epidist_latent_model.Rd @@ -27,7 +27,7 @@ Other latent_model: \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/is_epidist_latent_model.Rd b/man/is_epidist_latent_model.Rd index 76671a95f..54754bc4d 100644 --- a/man/is_epidist_latent_model.Rd +++ b/man/is_epidist_latent_model.Rd @@ -19,7 +19,7 @@ Other latent_model: \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/log_lik_epidist_latent.Rd b/man/log_lik_latent.Rd similarity index 83% rename from man/log_lik_epidist_latent.Rd rename to man/log_lik_latent.Rd index 755f0dfe7..c58db85f4 100644 --- a/man/log_lik_epidist_latent.Rd +++ b/man/log_lik_latent.Rd @@ -1,10 +1,10 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_model.R -\name{log_lik_epidist_latent} -\alias{log_lik_epidist_latent} +\name{log_lik_latent} +\alias{log_lik_latent} \title{Calculate the pointwise log likelihood of the \code{latent_gamma} family} \usage{ -log_lik_epidist_latent(i, prep) +log_lik_latent(i, prep) } \arguments{ \item{i}{The index of the observation to calculate the log likelihood of} diff --git a/man/new_epidist_latent_model.Rd b/man/new_epidist_latent_model.Rd index bfb7cfba6..9809489cf 100644 --- a/man/new_epidist_latent_model.Rd +++ b/man/new_epidist_latent_model.Rd @@ -22,7 +22,7 @@ Other latent_model: \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_epred_latent}()}, +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/posterior_epred_epidist_latent.Rd b/man/posterior_epred_latent.Rd similarity index 89% rename from man/posterior_epred_epidist_latent.Rd rename to man/posterior_epred_latent.Rd index 930f072b6..b257df5fe 100644 --- a/man/posterior_epred_epidist_latent.Rd +++ b/man/posterior_epred_latent.Rd @@ -1,11 +1,11 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_model.R -\name{posterior_epred_epidist_latent} -\alias{posterior_epred_epidist_latent} +\name{posterior_epred_latent} +\alias{posterior_epred_latent} \title{Create a function to draw from the expected value of the posterior predictive distribution for a latent model} \usage{ -posterior_epred_epidist_latent(family) +posterior_epred_latent(family) } \arguments{ \item{family}{Output of a call to \code{brms::brmsfamily()} with additional @@ -34,6 +34,6 @@ Other latent_model: \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_predict_epidist_latent}()} +\code{\link{posterior_predict_latent}()} } \concept{latent_model} diff --git a/man/posterior_predict_epidist_latent.Rd b/man/posterior_predict_latent.Rd similarity index 84% rename from man/posterior_predict_epidist_latent.Rd rename to man/posterior_predict_latent.Rd index fb4d45190..4d9dba15c 100644 --- a/man/posterior_predict_epidist_latent.Rd +++ b/man/posterior_predict_latent.Rd @@ -1,14 +1,15 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/latent_model.R -\name{posterior_predict_epidist_latent} -\alias{posterior_predict_epidist_latent} +\name{posterior_predict_latent} +\alias{posterior_predict_latent} \title{Create a function to draw from the posterior predictive distribution for a latent model} \usage{ -posterior_predict_epidist_latent(i, prep, ...) +posterior_predict_latent(family) } \arguments{ -\item{...}{Additional arguments passed to \code{fn} method.} +\item{family}{Output of a call to \code{brms::brmsfamily()} with additional +information as provided by \code{.add_dpar_info()}} } \value{ A function that takes a prep argument from brms and returns a matrix @@ -33,6 +34,6 @@ Other latent_model: \code{\link{epidist_formula_model.epidist_latent_model}()}, \code{\link{is_epidist_latent_model}()}, \code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_epidist_latent}()} +\code{\link{posterior_epred_latent}()} } \concept{latent_model} From 0a044e88bf87ac03446e3e955af683c00fae23cb Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 12:49:11 +0000 Subject: [PATCH 05/16] refine posterior prediction and make vignettes faster --- R/latent_gamma.R | 4 ---- R/latent_lognormal.R | 0 R/latent_model.R | 10 ++++------ vignettes/ebola.Rmd | 15 +++++++++------ vignettes/epidist.Rmd | 4 +++- 5 files changed, 16 insertions(+), 17 deletions(-) delete mode 100644 R/latent_lognormal.R diff --git a/R/latent_gamma.R b/R/latent_gamma.R index 698971b7b..46b7c2d87 100644 --- a/R/latent_gamma.R +++ b/R/latent_gamma.R @@ -27,10 +27,6 @@ posterior_predict_latent_gamma <- function(i, prep, ...) { # nolint: object_leng d_censored <- s_censored - p_censored } return(d_censored) - - .predict <- function(s) { - primarycensored::rpcens() - } } # Within brms this is a helper function called rblapply diff --git a/R/latent_lognormal.R b/R/latent_lognormal.R deleted file mode 100644 index e69de29bb..000000000 diff --git a/R/latent_model.R b/R/latent_model.R index 883ae0929..8b7690bbb 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -161,10 +161,8 @@ posterior_predict_latent <- function(family) { dist_fn <- eval(parse(text = fn_string)) rdist <- function(n, i, prep, ...) { - purrr::map_dbl( - seq_len(n), - ~ do.call(dist_fn, list(i = i, prep = prep)) - ) + prep$ndraws <- n + do.call(dist_fn, list(i = i, prep = prep)) } .predict <- function(i, prep, ...) { @@ -172,7 +170,7 @@ posterior_predict_latent <- function(family) { pwindow <- prep$data$vreal2[i] swindow <- prep$data$vreal3[i] - primarycensored::rpcens( + as.matrix(primarycensored::rpcens( n = prep$ndraws, rdist = rdist, rprimary = stats::runif, @@ -181,7 +179,7 @@ posterior_predict_latent <- function(family) { D = prep$data$vreal1[i], i = i, prep = prep - ) + )) } return(.predict) } diff --git a/vignettes/ebola.Rmd b/vignettes/ebola.Rmd index 0e3592597..a41a63fd5 100644 --- a/vignettes/ebola.Rmd +++ b/vignettes/ebola.Rmd @@ -239,8 +239,9 @@ fit <- epidist( formula = mu ~ 1, family = lognormal(), algorithm = "sampling", - refresh = 0, - silent = 2, + chains = 2, + cores = 2, + refresh = as.integer(interactive()), seed = 1, backend = "cmdstanr" ) @@ -264,8 +265,9 @@ fit_sex <- epidist( formula = bf(mu ~ 1 + sex, sigma ~ 1 + sex), family = lognormal(), algorithm = "sampling", - refresh = 0, - silent = 2, + chains = 2, + cores = 2, + refresh = as.integer(interactive()), seed = 1, backend = "cmdstanr" ) @@ -294,8 +296,9 @@ fit_sex_district <- epidist( ), family = lognormal(), algorithm = "sampling", - refresh = 0, - silent = 2, + chains = 2, + cores = 2, + refresh = as.integer(interactive()), seed = 1, backend = "cmdstanr" ) diff --git a/vignettes/epidist.Rmd b/vignettes/epidist.Rmd index 86fa9b254..b9896177a 100644 --- a/vignettes/epidist.Rmd +++ b/vignettes/epidist.Rmd @@ -285,7 +285,9 @@ The parameters of the model are inferred using Bayesian inference. In particular, we use the the No-U-Turn Sampler (NUTS) Markov chain Monte Carlo (MCMC) algorithm via the [`brms`](https://paul-buerkner.github.io/brms/) R package [@brms]. ```{r} -fit <- epidist(data = data, chains = 2, cores = 2, refresh = 0) +fit <- epidist( + data = data, chains = 2, cores = 2, refresh = as.integer(interactive()) +) ``` The `fit` object is a `brmsfit` object containing MCMC samples from each of the parameters in the model, shown in the table below. From 9f68a66a61c97a8eea90928fbc90bc82b360d899 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 17:53:53 +0000 Subject: [PATCH 06/16] rename and move around --- NAMESPACE | 4 +- R/gen.R | 79 ++++++++++ R/latent_gamma.R | 84 ---------- R/latent_model.R | 149 ++++++------------ _pkgdown.yml | 4 + man/as_epidist_latent_model.Rd | 5 +- ...dist_latent_model.epidist_linelist_data.Rd | 5 +- ...idist_family_model.epidist_latent_model.Rd | 5 +- ...dist_formula_model.epidist_latent_model.Rd | 5 +- man/epidist_gen_log_lik_latent.Rd | 41 +++++ ...tent.Rd => epidist_gen_posterior_epred.Rd} | 20 +-- ...nt.Rd => epidist_gen_posterior_predict.Rd} | 27 ++-- man/is_epidist_latent_model.Rd | 5 +- man/log_lik_latent.Rd | 17 -- man/log_lik_latent_gamma.Rd | 17 -- man/new_epidist_latent_model.Rd | 5 +- man/posterior_epred_latent_gamma.Rd | 16 -- man/posterior_predict_latent_gamma.Rd | 19 --- 18 files changed, 209 insertions(+), 298 deletions(-) create mode 100644 R/gen.R delete mode 100644 R/latent_gamma.R create mode 100644 man/epidist_gen_log_lik_latent.Rd rename man/{posterior_epred_latent.Rd => epidist_gen_posterior_epred.Rd} (65%) rename man/{posterior_predict_latent.Rd => epidist_gen_posterior_predict.Rd} (67%) delete mode 100644 man/log_lik_latent.Rd delete mode 100644 man/log_lik_latent_gamma.Rd delete mode 100644 man/posterior_epred_latent_gamma.Rd delete mode 100644 man/posterior_predict_latent_gamma.Rd diff --git a/NAMESPACE b/NAMESPACE index 2266b4020..5fa0b68b8 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -38,6 +38,8 @@ export(epidist_family_prior) export(epidist_family_reparam) export(epidist_formula) export(epidist_formula_model) +export(epidist_gen_posterior_epred) +export(epidist_gen_posterior_predict) export(epidist_model_prior) export(epidist_prior) export(epidist_stancode) @@ -48,8 +50,6 @@ export(lognormal) export(new_epidist_latent_model) export(new_epidist_linelist_data) export(new_epidist_naive_model) -export(posterior_epred_latent) -export(posterior_predict_latent) export(predict_delay_parameters) export(predict_dpar) export(simulate_exponential_cases) diff --git a/R/gen.R b/R/gen.R new file mode 100644 index 000000000..615d96697 --- /dev/null +++ b/R/gen.R @@ -0,0 +1,79 @@ +#' Create a function to draw from the posterior predictive distribution for a +#' latent model +#' +#' This function creates a function that draws from the posterior predictive +#' distribution for a latent model using [primarycensored::rpcens()] to handle +#' censoring and truncation. The returned function takes a prep argument from +#' `brms` and returns posterior predictions. This is used internally by +#' [brms::posterior_predict()] to generate predictions for latent models. +#' +#' @inheritParams epidist_family_model +#' +#' @return A function that takes a prep argument from brms and returns a matrix +#' of posterior predictions, with one row per posterior draw and one column +#' per observation. The prep object must have the following variables: +#' * `vreal1`: relative observation time +#' * `vreal2`: primary event window +#' * `vreal3`: secondary event window +#' +#' @seealso [brms::posterior_predict()] for details on how this is used within +#' `brms`, [primarycensored::rpcens()] for details on the censoring approach +#' @autoglobal +#' @family gen +#' @export +epidist_gen_posterior_predict <- function(family) { + dist_fn <- get( + paste0("posterior_predict_", family$family), + asNamespace("brms") + ) + + rdist <- function(n, i, prep, ...) { + prep$ndraws <- n + do.call(dist_fn, list(i = i, prep = prep)) + } + + .predict <- function(i, prep, ...) { + relative_obs_time <- prep$data$vreal1[i] + pwindow <- prep$data$vreal2[i] + swindow <- prep$data$vreal3[i] + + as.matrix(primarycensored::rpcens( + n = prep$ndraws, + rdist = rdist, + rprimary = stats::runif, + pwindow = prep$data$vreal2[i], + swindow = prep$data$vreal3[i], + D = prep$data$vreal1[i], + i = i, + prep = prep + )) + } + return(.predict) +} + +#' Create a function to draw from the expected value of the posterior predictive +#' distribution for a latent model +#' +#' This function creates a function that calculates the expected value of the +#' posterior predictive distribution for a latent model. The returned function +#' takes a prep argument (from brms) and returns posterior expected values. +#' This is used internally by [brms::posterior_epred()] to calculate expected +#' values for latent models. +#' +#' @inheritParams epidist_family_model +#' +#' @return A function that takes a prep argument from brms and returns a matrix +#' of posterior expected values, with one row per posterior draw and one column +#' per observation. +#' +#' @seealso [brms::posterior_epred()] for details on how this is used within +#' `brms`. +#' @autoglobal +#' @family gen +#' @export +epidist_gen_posterior_epred <- function(family) { + get( + paste0("posterior_epred_", family$family), + asNamespace("brms") + ) +} diff --git a/R/latent_gamma.R b/R/latent_gamma.R deleted file mode 100644 index 46b7c2d87..000000000 --- a/R/latent_gamma.R +++ /dev/null @@ -1,84 +0,0 @@ -#' Draws from the posterior predictive distribution of the `latent_gamma` family -#' -#' See [brms::posterior_predict()]. -#' -#' @param i The index of the observation to predict -#' @param prep The result of a call to [brms::posterior_predict()] -#' @param ... Additional arguments -#' @autoglobal -#' @keywords internal -posterior_predict_latent_gamma <- function(i, prep, ...) { # nolint: object_length_linter - mu <- brms::get_dpar(prep, "mu", i = i) - shape <- brms::get_dpar(prep, "shape", i = i) - - relative_obs_time <- prep$data$vreal1[i] - pwindow <- prep$data$vreal2[i] - swindow <- prep$data$vreal3[i] - - .predict <- function(s) { - d_censored <- relative_obs_time + 1 - # while loop to impose the truncation - while (d_censored > relative_obs_time) { - p_latent <- stats::runif(1, 0, 1) * pwindow - d_latent <- stats::rgamma(1, shape = shape[s], scale = mu[s] / shape[s]) - s_latent <- p_latent + d_latent - p_censored <- .floor_mult(p_latent, pwindow) - s_censored <- .floor_mult(s_latent, swindow) - d_censored <- s_censored - p_censored - } - return(d_censored) - } - - # Within brms this is a helper function called rblapply - do.call(rbind, lapply(seq_len(prep$ndraws), .predict)) -} - -#' Draws from the expected value of the posterior predictive distribution of the -#' `latent_gamma` family -#' -#' See [brms::posterior_epred()]. -#' -#' @param prep The result of a call to [`brms::prepare_predictions`] -#' @autoglobal -#' @keywords internal -posterior_epred_latent_gamma <- function(prep) { # nolint: object_length_linter - brms::get_dpar(prep, "mu") -} - -#' Calculate the pointwise log likelihood of the `latent_gamma` family -#' -#' See [brms::log_lik()]. -#' -#' @param i The index of the observation to calculate the log likelihood of -#' @param prep The result of a call to [brms::prepare_predictions()] -#' @autoglobal -#' @keywords internal -log_lik_latent_gamma <- function(i, prep) { - mu <- brms::get_dpar(prep, "mu", i = i) - shape <- brms::get_dpar(prep, "shape", i = i) - y <- prep$data$Y[i] - relative_obs_time <- prep$data$vreal1[i] - pwindow <- prep$data$vreal2[i] - swindow <- prep$data$vreal3[i] - - swindow_raw <- stats::runif(prep$ndraws) - pwindow_raw <- stats::runif(prep$ndraws) - - swindow <- swindow_raw * swindow - - # For no overlap calculate as usual, for overlap ensure pwindow < swindow - if (i %in% prep$data$noverlap) { - pwindow <- pwindow_raw * pwindow - } else { - pwindow <- pwindow_raw * swindow - } - - d <- y - pwindow + swindow - obs_time <- relative_obs_time - pwindow - lpdf <- stats::dgamma(d, shape = shape, scale = mu / shape, log = TRUE) - lcdf <- stats::pgamma( - obs_time, - shape = shape, scale = mu / shape, log.p = TRUE - ) - return(lpdf - lcdf) -} diff --git a/R/latent_model.R b/R/latent_model.R index 8b7690bbb..cba1eb202 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -90,123 +90,76 @@ epidist_family_model.epidist_latent_model <- function( type = family$type, vars = c("pwindow", "swindow", "vreal1"), loop = FALSE, - log_lik = log_lik_latent, - posterior_predict = posterior_predict_latent(family), - posterior_epred = posterior_epred_latent(family) + log_lik = epidist_gen_log_lik_latent(family), + posterior_predict = epidist_gen_posterior_predict(family), + posterior_epred = epidist_gen_posterior_epred(family) ) custom_family$reparm <- family$reparm return(custom_family) } -#' Calculate the pointwise log likelihood of the `latent_gamma` family +#' Create a function to calculate the pointwise log likelihood of the latent +#' model #' -#' See [brms::log_lik()]. +#' This function creates a log likelihood function that accounts for the latent +#' variables in the model, including primary and secondary event windows and +#' their overlap. The returned function calculates the log likelihood for a +#' single observation by augmenting the data with the latent variables and +#' using the underlying brms log likelihood function. #' -#' @param i The index of the observation to calculate the log likelihood of -#' @param prep The result of a call to [brms::prepare_predictions()] -#' @autoglobal -#' @keywords internal -log_lik_latent <- function(i, prep) { - mu <- brms::get_dpar(prep, "mu", i = i) - sigma <- brms::get_dpar(prep, "sigma", i = i) - y <- prep$data$Y[i] - relative_obs_time <- prep$data$vreal1[i] - pwindow <- prep$data$vreal2[i] - swindow <- prep$data$vreal3[i] - - # Generates values of the swindow_raw and pwindow_raw, but really these should - # be extracted from prep or the fitted raws somehow. See: - # https://github.com/epinowcast/epidist/issues/267 - swindow_raw <- stats::runif(prep$ndraws) - pwindow_raw <- stats::runif(prep$ndraws) - - swindow <- swindow_raw * swindow - - # For no overlap calculate as usual, for overlap ensure pwindow < swindow - if (i %in% prep$data$noverlap) { - pwindow <- pwindow_raw * pwindow - } else { - pwindow <- pwindow_raw * swindow - } - - d <- y - pwindow + swindow - obs_time <- relative_obs_time - pwindow - lpdf <- stats::dlnorm(d, meanlog = mu, sdlog = sigma, log = TRUE) - lcdf <- stats::plnorm(obs_time, meanlog = mu, sdlog = sigma, log.p = TRUE) - return(lpdf - lcdf) -} - -#' Create a function to draw from the posterior predictive distribution for a -#' latent model -#' -#' This function creates a function that draws from the posterior predictive -#' distribution for a latent model using [primarycensored::rpcens()] to handle -#' censoring and truncation. The returned function takes a prep argument from -#' `brms` and returns posterior predictions. This is used internally by -#' [brms::posterior_predict()] to generate predictions for latent models. +#' @seealso [brms::log_lik()] for details on the brms log likelihood interface. #' #' @inheritParams epidist_family_model #' -#' @return A function that takes a prep argument from brms and returns a matrix -#' of posterior predictions, with one row per posterior draw and one column -#' per observation. +#' @return A function that calculates the log likelihood for a single +#' observation. The prep object must have the following variables: +#' * `vreal1`: relative observation time +#' * `vreal2`: primary event window +#' * `vreal3`: secondary event window #' -#' @seealso [brms::posterior_predict()] for details on how this is used within -#' `brms`, [primarycensored::rpcens()] for details on the censoring approach -#' @autoglobal #' @family latent_model -#' @export -posterior_predict_latent <- function(family) { - fn_string <- paste0("brms:::posterior_predict_", family$family) - dist_fn <- eval(parse(text = fn_string)) - - rdist <- function(n, i, prep, ...) { - prep$ndraws <- n - do.call(dist_fn, list(i = i, prep = prep)) - } +#' @autoglobal +epidist_gen_log_lik_latent <- function(family) { + # Get internal brms log_lik function + log_lik_brms <- get( + paste0("log_lik_", family$family), + asNamespace("brms") + ) - .predict <- function(i, prep, ...) { + .log_lik <- function(i, prep) { + y <- prep$data$Y[i] relative_obs_time <- prep$data$vreal1[i] pwindow <- prep$data$vreal2[i] swindow <- prep$data$vreal3[i] - as.matrix(primarycensored::rpcens( - n = prep$ndraws, - rdist = rdist, - rprimary = stats::runif, - pwindow = prep$data$vreal2[i], - swindow = prep$data$vreal3[i], - D = prep$data$vreal1[i], - i = i, - prep = prep - )) + # Generates values of the swindow_raw and pwindow_raw, but really these + # should be extracted from prep or the fitted raws somehow. See: + # https://github.com/epinowcast/epidist/issues/267 + swindow_raw <- stats::runif(prep$ndraws) + pwindow_raw <- stats::runif(prep$ndraws) + + swindow <- swindow_raw * swindow + + # For no overlap calculate as usual, for overlap ensure pwindow < swindow + if (i %in% prep$data$noverlap) { + pwindow <- pwindow_raw * pwindow + } else { + pwindow <- pwindow_raw * swindow + } + + d <- y - pwindow + swindow + obs_time <- relative_obs_time - pwindow + # Create brms truncation upper bound + prep$data$ub <- rep(obs_time, length(prep$data$Y)) + # Update augmented data + prep$data$Y <- rep(y, length(prep$data$Y)) + + # Call internal brms log_lik function with augmented data + lpdf <- log_lik_brms(i, prep) + return(lpdf) } - return(.predict) -} -#' Create a function to draw from the expected value of the posterior predictive -#' distribution for a latent model -#' -#' This function creates a function that calculates the expected value of the -#' posterior predictive distribution for a latent model. The returned function -#' takes a prep argument (from brms) and returns posterior expected values. -#' This is used internally by [brms::posterior_epred()] to calculate expected -#' values for latent models. -#' -#' @inheritParams epidist_family_model -#' -#' @return A function that takes a prep argument from brms and returns a matrix -#' of posterior expected values, with one row per posterior draw and one column -#' per observation. -#' -#' @seealso [brms::posterior_epred()] for details on how this is used within -#' `brms`. -#' @autoglobal -#' @family latent_model -#' @export -posterior_epred_latent <- function(family) { - fn_string <- paste0("brms:::posterior_epred_", family$family) - eval(parse(text = fn_string)) + return(.log_lik) } #' Define the model-specific component of an `epidist` custom formula diff --git a/_pkgdown.yml b/_pkgdown.yml index 0f393849a..51ef0149f 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -72,6 +72,10 @@ reference: desc: Functions for specifying custom Stan code to put into `brms` contents: - has_concept("stan") +- title: Generator model functions + desc: Generator functions for creating family specific model utilities. + contents: + - has_concept("gen") - title: Simulation desc: Tools for simulating datasets contents: diff --git a/man/as_epidist_latent_model.Rd b/man/as_epidist_latent_model.Rd index 0a600d099..2aa6e600b 100644 --- a/man/as_epidist_latent_model.Rd +++ b/man/as_epidist_latent_model.Rd @@ -17,9 +17,8 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{epidist_gen_log_lik_latent}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{new_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/as_epidist_latent_model.epidist_linelist_data.Rd b/man/as_epidist_latent_model.epidist_linelist_data.Rd index d9620e5fb..03e839d66 100644 --- a/man/as_epidist_latent_model.epidist_linelist_data.Rd +++ b/man/as_epidist_latent_model.epidist_linelist_data.Rd @@ -17,9 +17,8 @@ Other latent_model: \code{\link{as_epidist_latent_model}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{epidist_gen_log_lik_latent}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{new_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/epidist_family_model.epidist_latent_model.Rd b/man/epidist_family_model.epidist_latent_model.Rd index 61a7fad8d..6794976a2 100644 --- a/man/epidist_family_model.epidist_latent_model.Rd +++ b/man/epidist_family_model.epidist_latent_model.Rd @@ -22,9 +22,8 @@ Other latent_model: \code{\link{as_epidist_latent_model}()}, \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{epidist_gen_log_lik_latent}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{new_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/epidist_formula_model.epidist_latent_model.Rd b/man/epidist_formula_model.epidist_latent_model.Rd index dbd10c3b9..6a5cdfc46 100644 --- a/man/epidist_formula_model.epidist_latent_model.Rd +++ b/man/epidist_formula_model.epidist_latent_model.Rd @@ -25,9 +25,8 @@ Other latent_model: \code{\link{as_epidist_latent_model}()}, \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, +\code{\link{epidist_gen_log_lik_latent}()}, \code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{new_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/epidist_gen_log_lik_latent.Rd b/man/epidist_gen_log_lik_latent.Rd new file mode 100644 index 000000000..289e96ab8 --- /dev/null +++ b/man/epidist_gen_log_lik_latent.Rd @@ -0,0 +1,41 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/latent_model.R +\name{epidist_gen_log_lik_latent} +\alias{epidist_gen_log_lik_latent} +\title{Create a function to calculate the pointwise log likelihood of the latent +model} +\usage{ +epidist_gen_log_lik_latent(family) +} +\arguments{ +\item{family}{Output of a call to \code{brms::brmsfamily()} with additional +information as provided by \code{.add_dpar_info()}} +} +\value{ +A function that calculates the log likelihood for a single +observation. The prep object must have the following variables: +\itemize{ +\item \code{vreal1}: relative observation time +\item \code{vreal2}: primary event window +\item \code{vreal3}: secondary event window +} +} +\description{ +This function creates a log likelihood function that accounts for the latent +variables in the model, including primary and secondary event windows and +their overlap. The returned function calculates the log likelihood for a +single observation by augmenting the data with the latent variables and +using the underlying brms log likelihood function. +} +\seealso{ +\code{\link[brms:log_lik.brmsfit]{brms::log_lik()}} for details on the brms log likelihood interface. + +Other latent_model: +\code{\link{as_epidist_latent_model}()}, +\code{\link{as_epidist_latent_model.epidist_linelist_data}()}, +\code{\link{epidist_family_model.epidist_latent_model}()}, +\code{\link{epidist_formula_model.epidist_latent_model}()}, +\code{\link{is_epidist_latent_model}()}, +\code{\link{new_epidist_latent_model}()} +} +\concept{latent_model} diff --git a/man/posterior_epred_latent.Rd b/man/epidist_gen_posterior_epred.Rd similarity index 65% rename from man/posterior_epred_latent.Rd rename to man/epidist_gen_posterior_epred.Rd index b257df5fe..d385fe3f6 100644 --- a/man/posterior_epred_latent.Rd +++ b/man/epidist_gen_posterior_epred.Rd @@ -1,11 +1,11 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_model.R -\name{posterior_epred_latent} -\alias{posterior_epred_latent} +% Please edit documentation in R/gen.R +\name{epidist_gen_posterior_epred} +\alias{epidist_gen_posterior_epred} \title{Create a function to draw from the expected value of the posterior predictive distribution for a latent model} \usage{ -posterior_epred_latent(family) +epidist_gen_posterior_epred(family) } \arguments{ \item{family}{Output of a call to \code{brms::brmsfamily()} with additional @@ -27,13 +27,7 @@ values for latent models. \code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}} for details on how this is used within \code{brms}. -Other latent_model: -\code{\link{as_epidist_latent_model}()}, -\code{\link{as_epidist_latent_model.epidist_linelist_data}()}, -\code{\link{epidist_family_model.epidist_latent_model}()}, -\code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_predict_latent}()} +Other gen: +\code{\link{epidist_gen_posterior_predict}()} } -\concept{latent_model} +\concept{gen} diff --git a/man/posterior_predict_latent.Rd b/man/epidist_gen_posterior_predict.Rd similarity index 67% rename from man/posterior_predict_latent.Rd rename to man/epidist_gen_posterior_predict.Rd index 4d9dba15c..fb60d3bd0 100644 --- a/man/posterior_predict_latent.Rd +++ b/man/epidist_gen_posterior_predict.Rd @@ -1,11 +1,11 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_model.R -\name{posterior_predict_latent} -\alias{posterior_predict_latent} +% Please edit documentation in R/gen.R +\name{epidist_gen_posterior_predict} +\alias{epidist_gen_posterior_predict} \title{Create a function to draw from the posterior predictive distribution for a latent model} \usage{ -posterior_predict_latent(family) +epidist_gen_posterior_predict(family) } \arguments{ \item{family}{Output of a call to \code{brms::brmsfamily()} with additional @@ -14,7 +14,12 @@ information as provided by \code{.add_dpar_info()}} \value{ A function that takes a prep argument from brms and returns a matrix of posterior predictions, with one row per posterior draw and one column -per observation. +per observation. The prep object must have the following variables: +\itemize{ +\item \code{vreal1}: relative observation time +\item \code{vreal2}: primary event window +\item \code{vreal3}: secondary event window +} } \description{ This function creates a function that draws from the posterior predictive @@ -27,13 +32,7 @@ censoring and truncation. The returned function takes a prep argument from \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}} for details on how this is used within \code{brms}, \code{\link[primarycensored:rprimarycensored]{primarycensored::rpcens()}} for details on the censoring approach -Other latent_model: -\code{\link{as_epidist_latent_model}()}, -\code{\link{as_epidist_latent_model.epidist_linelist_data}()}, -\code{\link{epidist_family_model.epidist_latent_model}()}, -\code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{is_epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()} +Other gen: +\code{\link{epidist_gen_posterior_epred}()} } -\concept{latent_model} +\concept{gen} diff --git a/man/is_epidist_latent_model.Rd b/man/is_epidist_latent_model.Rd index 54754bc4d..647725836 100644 --- a/man/is_epidist_latent_model.Rd +++ b/man/is_epidist_latent_model.Rd @@ -18,8 +18,7 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{new_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{epidist_gen_log_lik_latent}()}, +\code{\link{new_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/log_lik_latent.Rd b/man/log_lik_latent.Rd deleted file mode 100644 index c58db85f4..000000000 --- a/man/log_lik_latent.Rd +++ /dev/null @@ -1,17 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_model.R -\name{log_lik_latent} -\alias{log_lik_latent} -\title{Calculate the pointwise log likelihood of the \code{latent_gamma} family} -\usage{ -log_lik_latent(i, prep) -} -\arguments{ -\item{i}{The index of the observation to calculate the log likelihood of} - -\item{prep}{The result of a call to \code{\link[brms:prepare_predictions]{brms::prepare_predictions()}}} -} -\description{ -See \code{\link[brms:log_lik.brmsfit]{brms::log_lik()}}. -} -\keyword{internal} diff --git a/man/log_lik_latent_gamma.Rd b/man/log_lik_latent_gamma.Rd deleted file mode 100644 index 5b460412b..000000000 --- a/man/log_lik_latent_gamma.Rd +++ /dev/null @@ -1,17 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_gamma.R -\name{log_lik_latent_gamma} -\alias{log_lik_latent_gamma} -\title{Calculate the pointwise log likelihood of the \code{latent_gamma} family} -\usage{ -log_lik_latent_gamma(i, prep) -} -\arguments{ -\item{i}{The index of the observation to calculate the log likelihood of} - -\item{prep}{The result of a call to \code{\link[brms:prepare_predictions]{brms::prepare_predictions()}}} -} -\description{ -See \code{\link[brms:log_lik.brmsfit]{brms::log_lik()}}. -} -\keyword{internal} diff --git a/man/new_epidist_latent_model.Rd b/man/new_epidist_latent_model.Rd index 9809489cf..fe5e672a7 100644 --- a/man/new_epidist_latent_model.Rd +++ b/man/new_epidist_latent_model.Rd @@ -21,8 +21,7 @@ Other latent_model: \code{\link{as_epidist_latent_model.epidist_linelist_data}()}, \code{\link{epidist_family_model.epidist_latent_model}()}, \code{\link{epidist_formula_model.epidist_latent_model}()}, -\code{\link{is_epidist_latent_model}()}, -\code{\link{posterior_epred_latent}()}, -\code{\link{posterior_predict_latent}()} +\code{\link{epidist_gen_log_lik_latent}()}, +\code{\link{is_epidist_latent_model}()} } \concept{latent_model} diff --git a/man/posterior_epred_latent_gamma.Rd b/man/posterior_epred_latent_gamma.Rd deleted file mode 100644 index 1bc3a4e14..000000000 --- a/man/posterior_epred_latent_gamma.Rd +++ /dev/null @@ -1,16 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_gamma.R -\name{posterior_epred_latent_gamma} -\alias{posterior_epred_latent_gamma} -\title{Draws from the expected value of the posterior predictive distribution of the -\code{latent_gamma} family} -\usage{ -posterior_epred_latent_gamma(prep) -} -\arguments{ -\item{prep}{The result of a call to \code{\link[brms:prepare_predictions]{brms::prepare_predictions}}} -} -\description{ -See \code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}}. -} -\keyword{internal} diff --git a/man/posterior_predict_latent_gamma.Rd b/man/posterior_predict_latent_gamma.Rd deleted file mode 100644 index cc0a1d3dd..000000000 --- a/man/posterior_predict_latent_gamma.Rd +++ /dev/null @@ -1,19 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/latent_gamma.R -\name{posterior_predict_latent_gamma} -\alias{posterior_predict_latent_gamma} -\title{Draws from the posterior predictive distribution of the \code{latent_gamma} family} -\usage{ -posterior_predict_latent_gamma(i, prep, ...) -} -\arguments{ -\item{i}{The index of the observation to predict} - -\item{prep}{The result of a call to \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}}} - -\item{...}{Additional arguments} -} -\description{ -See \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}}. -} -\keyword{internal} From d5bf00f23acea57f907eba1805e6f5b5dcc6f8bf Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 18:07:05 +0000 Subject: [PATCH 07/16] refactor tests --- tests/testthat/test-gen.R | 127 +++++++++++++++++++++++++ tests/testthat/test-latent_gamma.R | 64 ------------- tests/testthat/test-latent_lognormal.R | 64 ------------- tests/testthat/test-latent_model.R | 22 +++++ 4 files changed, 149 insertions(+), 128 deletions(-) create mode 100644 tests/testthat/test-gen.R delete mode 100644 tests/testthat/test-latent_gamma.R delete mode 100644 tests/testthat/test-latent_lognormal.R diff --git a/tests/testthat/test-gen.R b/tests/testthat/test-gen.R new file mode 100644 index 000000000..5a0865381 --- /dev/null +++ b/tests/testthat/test-gen.R @@ -0,0 +1,127 @@ +test_that("epidist_gen_posterior_predict returns a function that outputs positive integers with length equal to draws", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + prep <- brms::prepare_predictions(fit) + i <- 1 + family <- epidist_family(data = prep_obs, family = lognormal()) + predict_fn <- epidist_gen_posterior_predict(family) + pred_i <- predict_fn(i = i, prep) + expect_identical(floor(pred_i), pred_i) + expect_length(pred_i, prep$ndraws) + expect_gte(min(pred_i), 0) + + # Test gamma + prep_gamma <- brms::prepare_predictions(fit_gamma) + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + pred_i_gamma <- predict_fn_gamma(i = i, prep_gamma) + expect_identical(floor(pred_i_gamma), pred_i_gamma) + expect_length(pred_i_gamma, prep_gamma$ndraws) + expect_gte(min(pred_i_gamma), 0) +}) + +test_that("epidist_gen_posterior_predict returns a function that errors for i out of bounds", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + prep <- brms::prepare_predictions(fit) + i_out_of_bounds <- length(prep$data$Y) + 1 + family <- epidist_family(data = prep_obs, family = lognormal()) + predict_fn <- epidist_gen_posterior_predict(family) + expect_error(predict_fn(i = i_out_of_bounds, prep)) + + # Test gamma + prep_gamma <- brms::prepare_predictions(fit_gamma) + i_out_of_bounds_gamma <- length(prep_gamma$data$Y) + 1 + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + expect_error(predict_fn_gamma(i = i_out_of_bounds_gamma, prep_gamma)) +}) + +test_that("epidist_gen_posterior_predict returns a function that can generate predictions with no censoring", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + family <- epidist_family(data = prep_obs, family = lognormal()) + predict_fn <- epidist_gen_posterior_predict(family) + draws <- data.frame(relative_obs_time = 1000, pwindow = 0, swindow = 0) |> + tidybayes::add_predicted_draws(fit, ndraws = 100) + expect_identical(draws$.draw, 1:100) + pred <- draws$.prediction + expect_gte(min(pred), 0) + expect_true(all(abs(pred - round(pred)) > .Machine$double.eps^0.5)) + + # Test gamma + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + draws_gamma <- data.frame( + relative_obs_time = 1000, pwindow = 0, swindow = 0 + ) |> + tidybayes::add_predicted_draws(fit_gamma, ndraws = 100) + expect_identical(draws_gamma$.draw, 1:100) + pred_gamma <- draws_gamma$.prediction + expect_gte(min(pred_gamma), 0) + expect_true( + all(abs(pred_gamma - round(pred_gamma)) > .Machine$double.eps^0.5) + ) +}) + +test_that("epidist_gen_posterior_predict returns a function that predicts delays in the 95% credible interval", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + prep <- brms::prepare_predictions(fit) + prep$ndraws <- 1000 # Down from the 4000 for time saving + family <- epidist_family(data = prep_obs, family = lognormal()) + predict_fn <- epidist_gen_posterior_predict(family) + q <- purrr::map_vec(seq_along(prep$data$Y), function(i) { + y <- predict_fn(i, prep) + ecdf <- ecdf(y) + q <- ecdf(prep$data$Y[i]) + return(q) + }) + expect_lt(quantile(q, 0.1), 0.3) + expect_gt(quantile(q, 0.9), 0.7) + expect_lt(min(q), 0.1) + expect_gt(max(q), 0.9) + expect_lt(mean(q), 0.65) + expect_gt(mean(q), 0.35) + + # Test gamma + prep_gamma <- brms::prepare_predictions(fit_gamma) + prep_gamma$ndraws <- 1000 + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + q_gamma <- purrr::map_vec(seq_along(prep_gamma$data$Y), function(i) { + y <- predict_fn_gamma(i, prep_gamma) + ecdf <- ecdf(y) + q <- ecdf(prep_gamma$data$Y[i]) + return(q) + }) + expect_lt(quantile(q_gamma, 0.1), 0.3) + expect_gt(quantile(q_gamma, 0.9), 0.7) + expect_lt(min(q_gamma), 0.1) + expect_gt(max(q_gamma), 0.9) + expect_lt(mean(q_gamma), 0.65) + expect_gt(mean(q_gamma), 0.35) +}) + +test_that("epidist_gen_posterior_epred returns a function that creates arrays with correct dimensions", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + prep <- brms::prepare_predictions(fit) + family <- epidist_family(data = prep_obs, family = lognormal()) + epred_fn <- epidist_gen_posterior_epred(family) + epred <- epred_fn(prep) + expect_setequal(class(epred), c("matrix", "array")) + expect_identical(nrow(epred), prep$ndraws) + expect_identical(ncol(epred), length(prep$data$Y)) + expect_gte(min(epred), 0) + + # Test gamma + prep_gamma <- brms::prepare_predictions(fit_gamma) + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + epred_fn_gamma <- epidist_gen_posterior_epred(family_gamma) + epred_gamma <- epred_fn_gamma(prep_gamma) + expect_setequal(class(epred_gamma), c("matrix", "array")) + expect_identical(nrow(epred_gamma), prep_gamma$ndraws) + expect_identical(ncol(epred_gamma), length(prep_gamma$data$Y)) + expect_gte(min(epred_gamma), 0) +}) diff --git a/tests/testthat/test-latent_gamma.R b/tests/testthat/test-latent_gamma.R deleted file mode 100644 index e91bbf4a5..000000000 --- a/tests/testthat/test-latent_gamma.R +++ /dev/null @@ -1,64 +0,0 @@ -test_that("posterior_predict_latent_gamma outputs positive integers with length equal to draws", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit_gamma) - i <- 1 - pred_i <- posterior_predict_latent_gamma(i = i, prep) - expect_identical(floor(pred_i), pred_i) - expect_length(pred_i, prep$ndraws) - expect_gte(min(pred_i), 0) -}) - -test_that("posterior_predict_latent_gamma errors for i out of bounds", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit_gamma) - i_out_of_bounds <- length(prep$data$Y) + 1 - expect_error(posterior_predict_latent_gamma(i = i_out_of_bounds, prep)) -}) - -test_that("posterior_predict_latent_gamma can generate predictions with no censoring", { # nolint: line_length_linter. - skip_on_cran() - draws <- data.frame(relative_obs_time = 1000, pwindow = 0, swindow = 0) |> - tidybayes::add_predicted_draws(fit_gamma, ndraws = 100) - expect_identical(draws$.draw, 1:100) - pred <- draws$.prediction - expect_gte(min(pred), 0) - expect_true(all(abs(pred - round(pred)) > .Machine$double.eps^0.5)) -}) - -test_that("posterior_predict_latent_gamma predicts delays for which the data is in the 95% credible interval", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit_gamma) - prep$ndraws <- 1000 # Down from the 4000 for time saving - q <- purrr::map_vec(seq_along(prep$data$Y), function(i) { - y <- posterior_predict_latent_gamma(i, prep) - ecdf <- ecdf(y) - q <- ecdf(prep$data$Y[i]) - return(q) - }) - expect_lt(quantile(q, 0.1), 0.3) - expect_gt(quantile(q, 0.9), 0.7) - expect_lt(min(q), 0.1) - expect_gt(max(q), 0.9) - expect_lt(mean(q), 0.65) - expect_gt(mean(q), 0.35) -}) - -test_that("posterior_epred_latent_gamma creates a array of non-negative numbers with the correct dimensions", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit_gamma) - epred <- posterior_epred_latent_gamma(prep) - expect_setequal(class(epred), c("matrix", "array")) - expect_identical(nrow(epred), prep$ndraws) - expect_identical(ncol(epred), length(prep$data$Y)) - expect_gte(min(epred), 0) -}) - -test_that("log_lik_latent_gamma produces a vector with length ndraws of finite non-NA numbers", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit_gamma) - i <- 1 - log_lik <- log_lik_latent_gamma(i, prep) - expect_length(log_lik, prep$ndraws) - expect_false(anyNA(log_lik)) - expect_true(all(is.finite(log_lik))) -}) diff --git a/tests/testthat/test-latent_lognormal.R b/tests/testthat/test-latent_lognormal.R deleted file mode 100644 index 5ec8d7489..000000000 --- a/tests/testthat/test-latent_lognormal.R +++ /dev/null @@ -1,64 +0,0 @@ -test_that("posterior_predict_latent_lognormal outputs positive integers with length equal to draws", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit) - i <- 1 - pred_i <- posterior_predict_latent_lognormal(i = i, prep) - expect_identical(floor(pred_i), pred_i) - expect_length(pred_i, prep$ndraws) - expect_gte(min(pred_i), 0) -}) - -test_that("posterior_predict_latent_lognormal errors for i out of bounds", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit) - i_out_of_bounds <- length(prep$data$Y) + 1 - expect_error(posterior_predict_latent_lognormal(i = i_out_of_bounds, prep)) -}) - -test_that("posterior_predict_latent_lognormal can generate predictions with no censoring", { # nolint: line_length_linter. - skip_on_cran() - draws <- data.frame(relative_obs_time = 1000, pwindow = 0, swindow = 0) |> - tidybayes::add_predicted_draws(fit, ndraws = 100) - expect_identical(draws$.draw, 1:100) - pred <- draws$.prediction - expect_gte(min(pred), 0) - expect_true(all(abs(pred - round(pred)) > .Machine$double.eps^0.5)) -}) - -test_that("posterior_predict_latent_lognormal predicts delays for which the data is in the 95% credible interval", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit) - prep$ndraws <- 1000 # Down from the 4000 for time saving - q <- purrr::map_vec(seq_along(prep$data$Y), function(i) { - y <- posterior_predict_latent_lognormal(i, prep) - ecdf <- ecdf(y) - q <- ecdf(prep$data$Y[i]) - return(q) - }) - expect_lt(quantile(q, 0.1), 0.3) - expect_gt(quantile(q, 0.9), 0.7) - expect_lt(min(q), 0.1) - expect_gt(max(q), 0.9) - expect_lt(mean(q), 0.65) - expect_gt(mean(q), 0.35) -}) - -test_that("posterior_epred_latent_lognormal creates a array of non-negative numbers with the correct dimensions", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit) - epred <- posterior_epred_latent_lognormal(prep) - expect_setequal(class(epred), c("matrix", "array")) - expect_identical(nrow(epred), prep$ndraws) - expect_identical(ncol(epred), length(prep$data$Y)) - expect_gte(min(epred), 0) -}) - -test_that("log_lik_latent_lognormal produces a vector with length ndraws of finite non-NA numbers", { # nolint: line_length_linter. - skip_on_cran() - prep <- brms::prepare_predictions(fit) - i <- 1 - log_lik <- log_lik_latent_lognormal(i, prep) - expect_length(log_lik, prep$ndraws) - expect_false(anyNA(log_lik)) - expect_true(all(is.finite(log_lik))) -}) diff --git a/tests/testthat/test-latent_model.R b/tests/testthat/test-latent_model.R index 579af69fc..53cd4f40d 100644 --- a/tests/testthat/test-latent_model.R +++ b/tests/testthat/test-latent_model.R @@ -56,3 +56,25 @@ test_that("epidist_stancode.epidist_latent_model produces valid stanvars", { # n ) expect_s3_class(stancode, "stanvars") }) + +test_that("epidist_gen_log_lik_latent returns a function that produces valid log likelihoods", { # nolint: line_length_linter. + skip_on_cran() + # Test lognormal + prep <- brms::prepare_predictions(fit) + i <- 1 + family <- epidist_family(data = prep_obs, family = lognormal()) + log_lik_fn <- epidist_gen_log_lik_latent(family) + log_lik <- log_lik_fn(i = i, prep) + expect_length(log_lik, prep$ndraws) + expect_false(anyNA(log_lik)) + expect_true(all(is.finite(log_lik))) + + # Test gamma + prep_gamma <- brms::prepare_predictions(fit_gamma) + family_gamma <- epidist_family(data = prep_obs, family = gamma()) + log_lik_fn_gamma <- epidist_gen_log_lik_latent(family_gamma) + log_lik_gamma <- log_lik_fn_gamma(i = i, prep_gamma) + expect_length(log_lik_gamma, prep_gamma$ndraws) + expect_false(anyNA(log_lik_gamma)) + expect_true(all(is.finite(log_lik_gamma))) +}) From c2cf2f2f30421e5c2da7620fc572690dff96250c Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 18:13:54 +0000 Subject: [PATCH 08/16] check tests --- tests/testthat/test-gen.R | 10 +++++----- tests/testthat/test-latent_model.R | 2 +- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/tests/testthat/test-gen.R b/tests/testthat/test-gen.R index 5a0865381..39ea95ef2 100644 --- a/tests/testthat/test-gen.R +++ b/tests/testthat/test-gen.R @@ -12,7 +12,7 @@ test_that("epidist_gen_posterior_predict returns a function that outputs positiv # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) pred_i_gamma <- predict_fn_gamma(i = i, prep_gamma) expect_identical(floor(pred_i_gamma), pred_i_gamma) @@ -32,7 +32,7 @@ test_that("epidist_gen_posterior_predict returns a function that errors for i ou # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) i_out_of_bounds_gamma <- length(prep_gamma$data$Y) + 1 - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) expect_error(predict_fn_gamma(i = i_out_of_bounds_gamma, prep_gamma)) }) @@ -50,7 +50,7 @@ test_that("epidist_gen_posterior_predict returns a function that can generate pr expect_true(all(abs(pred - round(pred)) > .Machine$double.eps^0.5)) # Test gamma - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) draws_gamma <- data.frame( relative_obs_time = 1000, pwindow = 0, swindow = 0 @@ -87,7 +87,7 @@ test_that("epidist_gen_posterior_predict returns a function that predicts delays # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) prep_gamma$ndraws <- 1000 - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) q_gamma <- purrr::map_vec(seq_along(prep_gamma$data$Y), function(i) { y <- predict_fn_gamma(i, prep_gamma) @@ -117,7 +117,7 @@ test_that("epidist_gen_posterior_epred returns a function that creates arrays wi # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) epred_fn_gamma <- epidist_gen_posterior_epred(family_gamma) epred_gamma <- epred_fn_gamma(prep_gamma) expect_setequal(class(epred_gamma), c("matrix", "array")) diff --git a/tests/testthat/test-latent_model.R b/tests/testthat/test-latent_model.R index 53cd4f40d..ca8228380 100644 --- a/tests/testthat/test-latent_model.R +++ b/tests/testthat/test-latent_model.R @@ -71,7 +71,7 @@ test_that("epidist_gen_log_lik_latent returns a function that produces valid log # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = gamma()) + family_gamma <- epidist_family(data = prep_obs, family = Gamma()) log_lik_fn_gamma <- epidist_gen_log_lik_latent(family_gamma) log_lik_gamma <- log_lik_fn_gamma(i = i, prep_gamma) expect_length(log_lik_gamma, prep_gamma$ndraws) From 418cb37d9b55b0bcad473aa502235b7f4870ff9b Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 19:59:15 +0000 Subject: [PATCH 09/16] refine tess --- R/gen.R | 10 ++-------- R/latent_model.R | 5 +---- R/utils.R | 19 +++++++++++++++++++ tests/testthat/test-gen.R | 40 ++++++++++++++++++--------------------- 4 files changed, 40 insertions(+), 34 deletions(-) diff --git a/R/gen.R b/R/gen.R index 615d96697..abd86190f 100644 --- a/R/gen.R +++ b/R/gen.R @@ -22,10 +22,7 @@ #' @family gen #' @export epidist_gen_posterior_predict <- function(family) { - dist_fn <- get( - paste0("posterior_predict_", family$family), - asNamespace("brms") - ) + dist_fn <- .get_brms_fn("posterior_predict", family) rdist <- function(n, i, prep, ...) { prep$ndraws <- n @@ -72,8 +69,5 @@ epidist_gen_posterior_predict <- function(family) { #' @family gen #' @export epidist_gen_posterior_epred <- function(family) { - get( - paste0("posterior_epred_", family$family), - asNamespace("brms") - ) + .get_brms_fn("posterior_epred", family) } diff --git a/R/latent_model.R b/R/latent_model.R index cba1eb202..11a94c201 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -121,10 +121,7 @@ epidist_family_model.epidist_latent_model <- function( #' @autoglobal epidist_gen_log_lik_latent <- function(family) { # Get internal brms log_lik function - log_lik_brms <- get( - paste0("log_lik_", family$family), - asNamespace("brms") - ) + log_lik_brms <- .get_brms_fn("log_lik", family) .log_lik <- function(i, prep) { y <- prep$data$Y[i] diff --git a/R/utils.R b/R/utils.R index 00629dab5..ef2d3f73a 100644 --- a/R/utils.R +++ b/R/utils.R @@ -157,3 +157,22 @@ return(df) } + +#' Get a brms function by prefix and family +#' +#' Helper function to get internal brms functions by constructing their name +#' from a prefix and family. Used to get functions like `log_lik_*`, +#' `posterior_predict_*` etc. +#' +#' @param prefix Character string prefix of the brms function to get (e.g. +#' "log_lik") +#' +#' @param family A brms family object +#' @return The requested brms function +#' @keywords internal +.get_brms_fn <- function(prefix, family) { + get( + paste0(prefix, "_", tolower(family$family)), + asNamespace("brms") + ) +} diff --git a/tests/testthat/test-gen.R b/tests/testthat/test-gen.R index 39ea95ef2..3e6a246b5 100644 --- a/tests/testthat/test-gen.R +++ b/tests/testthat/test-gen.R @@ -3,8 +3,7 @@ test_that("epidist_gen_posterior_predict returns a function that outputs positiv # Test lognormal prep <- brms::prepare_predictions(fit) i <- 1 - family <- epidist_family(data = prep_obs, family = lognormal()) - predict_fn <- epidist_gen_posterior_predict(family) + predict_fn <- epidist_gen_posterior_predict(lognormal()) pred_i <- predict_fn(i = i, prep) expect_identical(floor(pred_i), pred_i) expect_length(pred_i, prep$ndraws) @@ -12,8 +11,7 @@ test_that("epidist_gen_posterior_predict returns a function that outputs positiv # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + predict_fn_gamma <- epidist_gen_posterior_predict(Gamma()) pred_i_gamma <- predict_fn_gamma(i = i, prep_gamma) expect_identical(floor(pred_i_gamma), pred_i_gamma) expect_length(pred_i_gamma, prep_gamma$ndraws) @@ -25,23 +23,26 @@ test_that("epidist_gen_posterior_predict returns a function that errors for i ou # Test lognormal prep <- brms::prepare_predictions(fit) i_out_of_bounds <- length(prep$data$Y) + 1 - family <- epidist_family(data = prep_obs, family = lognormal()) - predict_fn <- epidist_gen_posterior_predict(family) - expect_error(predict_fn(i = i_out_of_bounds, prep)) + predict_fn <- epidist_gen_posterior_predict(lognormal()) + expect_warning( + expect_error( + predict_fn(i = i_out_of_bounds, prep) + ) + ) # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) i_out_of_bounds_gamma <- length(prep_gamma$data$Y) + 1 - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) - expect_error(predict_fn_gamma(i = i_out_of_bounds_gamma, prep_gamma)) + predict_fn_gamma <- epidist_gen_posterior_predict(Gamma()) + expect_warning( + expect_error(predict_fn_gamma(i = i_out_of_bounds_gamma, prep_gamma)) + ) }) test_that("epidist_gen_posterior_predict returns a function that can generate predictions with no censoring", { # nolint: line_length_linter. skip_on_cran() # Test lognormal - family <- epidist_family(data = prep_obs, family = lognormal()) - predict_fn <- epidist_gen_posterior_predict(family) + predict_fn <- epidist_gen_posterior_predict(lognormal()) draws <- data.frame(relative_obs_time = 1000, pwindow = 0, swindow = 0) |> tidybayes::add_predicted_draws(fit, ndraws = 100) expect_identical(draws$.draw, 1:100) @@ -50,8 +51,7 @@ test_that("epidist_gen_posterior_predict returns a function that can generate pr expect_true(all(abs(pred - round(pred)) > .Machine$double.eps^0.5)) # Test gamma - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + predict_fn_gamma <- epidist_gen_posterior_predict(Gamma()) draws_gamma <- data.frame( relative_obs_time = 1000, pwindow = 0, swindow = 0 ) |> @@ -69,8 +69,7 @@ test_that("epidist_gen_posterior_predict returns a function that predicts delays # Test lognormal prep <- brms::prepare_predictions(fit) prep$ndraws <- 1000 # Down from the 4000 for time saving - family <- epidist_family(data = prep_obs, family = lognormal()) - predict_fn <- epidist_gen_posterior_predict(family) + predict_fn <- epidist_gen_posterior_predict(lognormal()) q <- purrr::map_vec(seq_along(prep$data$Y), function(i) { y <- predict_fn(i, prep) ecdf <- ecdf(y) @@ -87,8 +86,7 @@ test_that("epidist_gen_posterior_predict returns a function that predicts delays # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) prep_gamma$ndraws <- 1000 - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - predict_fn_gamma <- epidist_gen_posterior_predict(family_gamma) + predict_fn_gamma <- epidist_gen_posterior_predict(Gamma()) q_gamma <- purrr::map_vec(seq_along(prep_gamma$data$Y), function(i) { y <- predict_fn_gamma(i, prep_gamma) ecdf <- ecdf(y) @@ -107,8 +105,7 @@ test_that("epidist_gen_posterior_epred returns a function that creates arrays wi skip_on_cran() # Test lognormal prep <- brms::prepare_predictions(fit) - family <- epidist_family(data = prep_obs, family = lognormal()) - epred_fn <- epidist_gen_posterior_epred(family) + epred_fn <- epidist_gen_posterior_epred(lognormal()) epred <- epred_fn(prep) expect_setequal(class(epred), c("matrix", "array")) expect_identical(nrow(epred), prep$ndraws) @@ -117,8 +114,7 @@ test_that("epidist_gen_posterior_epred returns a function that creates arrays wi # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - epred_fn_gamma <- epidist_gen_posterior_epred(family_gamma) + epred_fn_gamma <- epidist_gen_posterior_epred(Gamma()) epred_gamma <- epred_fn_gamma(prep_gamma) expect_setequal(class(epred_gamma), c("matrix", "array")) expect_identical(nrow(epred_gamma), prep_gamma$ndraws) From 4f21389aabb8a9893fde56562be272de1339e4f5 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 20:50:04 +0000 Subject: [PATCH 10/16] add note about direct usage failing but tidybayes working --- tests/testthat/test-gen.R | 11 ++++------- tests/testthat/test-latent_model.R | 6 ++---- vignettes/ebola.Rmd | 3 --- 3 files changed, 6 insertions(+), 14 deletions(-) diff --git a/tests/testthat/test-gen.R b/tests/testthat/test-gen.R index 3e6a246b5..f51162663 100644 --- a/tests/testthat/test-gen.R +++ b/tests/testthat/test-gen.R @@ -104,13 +104,10 @@ test_that("epidist_gen_posterior_predict returns a function that predicts delays test_that("epidist_gen_posterior_epred returns a function that creates arrays with correct dimensions", { # nolint: line_length_linter. skip_on_cran() # Test lognormal - prep <- brms::prepare_predictions(fit) - epred_fn <- epidist_gen_posterior_epred(lognormal()) - epred <- epred_fn(prep) - expect_setequal(class(epred), c("matrix", "array")) - expect_identical(nrow(epred), prep$ndraws) - expect_identical(ncol(epred), length(prep$data$Y)) - expect_gte(min(epred), 0) + epred <- prep_obs |> + tidybayes::add_epred_draws(fit) + expect_equal(mean(epred$.epred), 5.97, tolerance = 0.1) + expect_gte(min(epred$.epred), 0) # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) diff --git a/tests/testthat/test-latent_model.R b/tests/testthat/test-latent_model.R index ca8228380..c1a64560f 100644 --- a/tests/testthat/test-latent_model.R +++ b/tests/testthat/test-latent_model.R @@ -62,8 +62,7 @@ test_that("epidist_gen_log_lik_latent returns a function that produces valid log # Test lognormal prep <- brms::prepare_predictions(fit) i <- 1 - family <- epidist_family(data = prep_obs, family = lognormal()) - log_lik_fn <- epidist_gen_log_lik_latent(family) + log_lik_fn <- epidist_gen_log_lik_latent(lognormal()) log_lik <- log_lik_fn(i = i, prep) expect_length(log_lik, prep$ndraws) expect_false(anyNA(log_lik)) @@ -71,8 +70,7 @@ test_that("epidist_gen_log_lik_latent returns a function that produces valid log # Test gamma prep_gamma <- brms::prepare_predictions(fit_gamma) - family_gamma <- epidist_family(data = prep_obs, family = Gamma()) - log_lik_fn_gamma <- epidist_gen_log_lik_latent(family_gamma) + log_lik_fn_gamma <- epidist_gen_log_lik_latent(Gamma()) log_lik_gamma <- log_lik_fn_gamma(i = i, prep_gamma) expect_length(log_lik_gamma, prep_gamma$ndraws) expect_false(anyNA(log_lik_gamma)) diff --git a/vignettes/ebola.Rmd b/vignettes/ebola.Rmd index a41a63fd5..1dbf76016 100644 --- a/vignettes/ebola.Rmd +++ b/vignettes/ebola.Rmd @@ -322,7 +322,6 @@ Figure \@ref(fig:epred)B illustrates the higher mean of men as compared with wom ```{r} epred_draws <- obs_prep |> - as.data.frame() |> data_grid(NA) |> mutate(relative_obs_time = NA, pwindow = NA, swindow = NA) |> add_epred_draws(fit, dpar = TRUE) @@ -334,7 +333,6 @@ epred_base_figure <- epred_draws |> theme_minimal() epred_draws_sex <- obs_prep |> - as.data.frame() |> data_grid(sex) |> mutate(relative_obs_time = NA, pwindow = NA, swindow = NA) |> add_epred_draws(fit_sex, dpar = TRUE) @@ -346,7 +344,6 @@ epred_sex_figure <- epred_draws_sex |> theme_minimal() epred_draws_sex_district <- obs_prep |> - as.data.frame() |> data_grid(sex, district) |> mutate(relative_obs_time = NA, pwindow = NA, swindow = NA) |> add_epred_draws(fit_sex_district, dpar = TRUE) From 477d7a6af14b194fda02373ff2d01598b99e1ac7 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 20:59:14 +0000 Subject: [PATCH 11/16] catch lok_lik bug --- R/latent_model.R | 2 +- man/dot-get_brms_fn.Rd | 23 +++++++++++++++++++++++ 2 files changed, 24 insertions(+), 1 deletion(-) create mode 100644 man/dot-get_brms_fn.Rd diff --git a/R/latent_model.R b/R/latent_model.R index 11a94c201..6eba158d5 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -149,7 +149,7 @@ epidist_gen_log_lik_latent <- function(family) { # Create brms truncation upper bound prep$data$ub <- rep(obs_time, length(prep$data$Y)) # Update augmented data - prep$data$Y <- rep(y, length(prep$data$Y)) + prep$data$Y <- rep(d, length(prep$data$Y)) # Call internal brms log_lik function with augmented data lpdf <- log_lik_brms(i, prep) diff --git a/man/dot-get_brms_fn.Rd b/man/dot-get_brms_fn.Rd new file mode 100644 index 000000000..c1ad322e0 --- /dev/null +++ b/man/dot-get_brms_fn.Rd @@ -0,0 +1,23 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/utils.R +\name{.get_brms_fn} +\alias{.get_brms_fn} +\title{Get a brms function by prefix and family} +\usage{ +.get_brms_fn(prefix, family) +} +\arguments{ +\item{prefix}{Character string prefix of the brms function to get (e.g. +"log_lik")} + +\item{family}{A brms family object} +} +\value{ +The requested brms function +} +\description{ +Helper function to get internal brms functions by constructing their name +from a prefix and family. Used to get functions like \verb{log_lik_*}, +\verb{posterior_predict_*} etc. +} +\keyword{internal} From 98960df66617666fc29218c6cd08b70fe8d50924 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 21:08:06 +0000 Subject: [PATCH 12/16] add a summary section to getting started as dev sense check + it feels like gap --- vignettes/epidist.Rmd | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/vignettes/epidist.Rmd b/vignettes/epidist.Rmd index b9896177a..262369036 100644 --- a/vignettes/epidist.Rmd +++ b/vignettes/epidist.Rmd @@ -310,6 +310,14 @@ Any tool that supports `brms` fitted model objects will be compatible with `fit` **Note that here we use the default `rstan` backend but we generally recommend using the `cmdstanr` backend for faster sampling and additional features.** **This can be set using `backend = "cmdstanr"` after following the installing CmdStan instructions in the README.** +For example, we can use the built in `summary()` function + +```{r} +summary(fit) +``` + +to see that the model has converged and recovered the delay distribution distribution paramters reasonably well. + The `epidist` package also provides functions to make common post-processing tasks easy. For example, individual predictions of the lognormal delay parameters can be extracted using: From 510658fbe8430bfff8334e6ee6b872761fcb61b4 Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 21:12:18 +0000 Subject: [PATCH 13/16] add spaghetti plot to getting started --- vignettes/epidist.Rmd | 31 +++++++++++++++++++++---------- 1 file changed, 21 insertions(+), 10 deletions(-) diff --git a/vignettes/epidist.Rmd b/vignettes/epidist.Rmd index 262369036..9293eeab1 100644 --- a/vignettes/epidist.Rmd +++ b/vignettes/epidist.Rmd @@ -325,31 +325,42 @@ For example, individual predictions of the lognormal delay parameters can be ext pred <- predict_delay_parameters(fit) ``` -Figure \@ref(fig:fitted-lognormal) shows the lognormal delay distribution obtained using the average of the `mu` and `sigma` draws. +Figure \@ref(fig:fitted-lognormal) shows the lognormal delay distribution obtained for 100 draws from the posterior distribution of the `mu` and `sigma` parameters. Whereas in Figure \@ref(fig:obs-est) the histogram of censored, truncated, sampled data was substantially different to the underlying delay distribution, using `epidist()` we have obtained a much closer match to the truth. (ref:fitted-lognormal) A fitted delay distribution (in pink) as compared with the true delay distribution (in black). ```{r fitted-lognormal, fig.cap="(ref:fitted-lognormal)"} +# Sample 100 random draws from the posterior +set.seed(123) +sample_draws <- sample(nrow(pred), 100) + ggplot() + + # Plot the true distribution geom_function( data = data.frame(x = c(0, 30)), aes(x = x), fun = dlnorm, + linewidth = 2, args = list( meanlog = secondary_dist[["mu"]], sdlog = secondary_dist[["sigma"]] ) ) + - geom_function( - data = data.frame(x = c(0, 30)), - aes(x = x), fun = dlnorm, - args = list( - meanlog = mean(pred$mu), - sdlog = mean(pred$sigma) - ), - col = "#CC79A7" - ) + + # Plot 100 sampled posterior distributions + lapply(sample_draws, \(i) { + geom_function( + data = data.frame(x = c(0, 30)), + aes(x = x), + fun = dlnorm, + args = list( + meanlog = pred$mu[i], + sdlog = pred$sigma[i] + ), + alpha = 0.1, + color = "#CC79A7" + ) + }) + labs( x = "Delay between primary and secondary event (days)", y = "Probability density" From cbec12ccda19f467375fd5501c84c1225acfc27f Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 21:14:16 +0000 Subject: [PATCH 14/16] update fig caption --- vignettes/epidist.Rmd | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/vignettes/epidist.Rmd b/vignettes/epidist.Rmd index 9293eeab1..40a561d41 100644 --- a/vignettes/epidist.Rmd +++ b/vignettes/epidist.Rmd @@ -328,7 +328,7 @@ pred <- predict_delay_parameters(fit) Figure \@ref(fig:fitted-lognormal) shows the lognormal delay distribution obtained for 100 draws from the posterior distribution of the `mu` and `sigma` parameters. Whereas in Figure \@ref(fig:obs-est) the histogram of censored, truncated, sampled data was substantially different to the underlying delay distribution, using `epidist()` we have obtained a much closer match to the truth. -(ref:fitted-lognormal) A fitted delay distribution (in pink) as compared with the true delay distribution (in black). +(ref:fitted-lognormal) The fitted delay distribution (pink lines, 100 draws from the posterior) compared to the true underlying delay distribution (black line). Each pink line represents one possible delay distribution based on the fitted model parameters. ```{r fitted-lognormal, fig.cap="(ref:fitted-lognormal)"} # Sample 100 random draws from the posterior From 57fb92418039ea171b40ac8d9ceda1c741c0861d Mon Sep 17 00:00:00 2001 From: Sam Date: Wed, 20 Nov 2024 22:04:37 +0000 Subject: [PATCH 15/16] switch gamma test over to tidybayes as well --- tests/testthat/test-gen.R | 11 ++++------- 1 file changed, 4 insertions(+), 7 deletions(-) diff --git a/tests/testthat/test-gen.R b/tests/testthat/test-gen.R index f51162663..5dc4c5686 100644 --- a/tests/testthat/test-gen.R +++ b/tests/testthat/test-gen.R @@ -110,11 +110,8 @@ test_that("epidist_gen_posterior_epred returns a function that creates arrays wi expect_gte(min(epred$.epred), 0) # Test gamma - prep_gamma <- brms::prepare_predictions(fit_gamma) - epred_fn_gamma <- epidist_gen_posterior_epred(Gamma()) - epred_gamma <- epred_fn_gamma(prep_gamma) - expect_setequal(class(epred_gamma), c("matrix", "array")) - expect_identical(nrow(epred_gamma), prep_gamma$ndraws) - expect_identical(ncol(epred_gamma), length(prep_gamma$data$Y)) - expect_gte(min(epred_gamma), 0) + epred_gamma <- prep_obs |> + tidybayes::add_epred_draws(fit_gamma) + expect_equal(mean(epred_gamma$.epred), 6.56, tolerance = 0.1) + expect_gte(min(epred_gamma$.epred), 0) }) From 2e6cd342510932c20ffdb99f47b33db08676dd7f Mon Sep 17 00:00:00 2001 From: Sam Date: Thu, 21 Nov 2024 18:34:40 +0000 Subject: [PATCH 16/16] deal with review comments --- R/gen.R | 14 +++++++------- R/latent_model.R | 2 +- R/utils.R | 2 +- man/dot-get_brms_fn.Rd | 7 ++++++- man/epidist_gen_log_lik_latent.Rd | 8 ++++++-- man/epidist_gen_posterior_epred.Rd | 10 +++++++--- man/epidist_gen_posterior_predict.Rd | 16 ++++++++++------ vignettes/ebola.Rmd | 2 +- vignettes/epidist.Rmd | 8 ++++---- 9 files changed, 43 insertions(+), 26 deletions(-) diff --git a/R/gen.R b/R/gen.R index abd86190f..ebef86b3b 100644 --- a/R/gen.R +++ b/R/gen.R @@ -3,15 +3,15 @@ #' #' This function creates a function that draws from the posterior predictive #' distribution for a latent model using [primarycensored::rpcens()] to handle -#' censoring and truncation. The returned function takes a prep argument from +#' censoring and truncation. The returned function takes a `prep` argument from #' `brms` and returns posterior predictions. This is used internally by #' [brms::posterior_predict()] to generate predictions for latent models. #' -#' @inheritParams epidist_family_model +#' @inheritParams epidist_family #' -#' @return A function that takes a prep argument from brms and returns a matrix -#' of posterior predictions, with one row per posterior draw and one column -#' per observation. The prep object must have the following variables: +#' @return A function that takes a `prep` argument from brms and returns a +#' matrix of posterior predictions, with one row per posterior draw and one +#' column per observation. The `prep` object must have the following variables: #' * `vreal1`: relative observation time #' * `vreal2`: primary event window #' * `vreal3`: secondary event window @@ -53,11 +53,11 @@ epidist_gen_posterior_predict <- function(family) { #' #' This function creates a function that calculates the expected value of the #' posterior predictive distribution for a latent model. The returned function -#' takes a prep argument (from brms) and returns posterior expected values. +#' takes a `prep` argument (from brms) and returns posterior expected values. #' This is used internally by [brms::posterior_epred()] to calculate expected #' values for latent models. #' -#' @inheritParams epidist_family_model +#' @inheritParams epidist_family #' #' @return A function that takes a prep argument from brms and returns a matrix #' of posterior expected values, with one row per posterior draw and one column diff --git a/R/latent_model.R b/R/latent_model.R index 6eba158d5..96215b971 100644 --- a/R/latent_model.R +++ b/R/latent_model.R @@ -109,7 +109,7 @@ epidist_family_model.epidist_latent_model <- function( #' #' @seealso [brms::log_lik()] for details on the brms log likelihood interface. #' -#' @inheritParams epidist_family_model +#' @inheritParams epidist_family #' #' @return A function that calculates the log likelihood for a single #' observation. The prep object must have the following variables: diff --git a/R/utils.R b/R/utils.R index ef2d3f73a..9fd6a419e 100644 --- a/R/utils.R +++ b/R/utils.R @@ -167,7 +167,7 @@ #' @param prefix Character string prefix of the brms function to get (e.g. #' "log_lik") #' -#' @param family A brms family object +#' @inheritParams epidist_family #' @return The requested brms function #' @keywords internal .get_brms_fn <- function(prefix, family) { diff --git a/man/dot-get_brms_fn.Rd b/man/dot-get_brms_fn.Rd index c1ad322e0..414c7a6e0 100644 --- a/man/dot-get_brms_fn.Rd +++ b/man/dot-get_brms_fn.Rd @@ -10,7 +10,12 @@ \item{prefix}{Character string prefix of the brms function to get (e.g. "log_lik")} -\item{family}{A brms family object} +\item{family}{A description of the response distribution and link function to +be used in the model. Every family function has a link argument allowing +users to specify the link function to be applied on the response variable. +If not specified, default links are used. For details of all supported +families see \code{\link[=brmsfamily]{brmsfamily()}}. Commonly used, such as \code{\link[=lognormal]{lognormal()}}, are also +reexported as part of \code{epidist}.} } \value{ The requested brms function diff --git a/man/epidist_gen_log_lik_latent.Rd b/man/epidist_gen_log_lik_latent.Rd index 289e96ab8..5fc1e1dee 100644 --- a/man/epidist_gen_log_lik_latent.Rd +++ b/man/epidist_gen_log_lik_latent.Rd @@ -8,8 +8,12 @@ model} epidist_gen_log_lik_latent(family) } \arguments{ -\item{family}{Output of a call to \code{brms::brmsfamily()} with additional -information as provided by \code{.add_dpar_info()}} +\item{family}{A description of the response distribution and link function to +be used in the model. Every family function has a link argument allowing +users to specify the link function to be applied on the response variable. +If not specified, default links are used. For details of all supported +families see \code{\link[=brmsfamily]{brmsfamily()}}. Commonly used, such as \code{\link[=lognormal]{lognormal()}}, are also +reexported as part of \code{epidist}.} } \value{ A function that calculates the log likelihood for a single diff --git a/man/epidist_gen_posterior_epred.Rd b/man/epidist_gen_posterior_epred.Rd index d385fe3f6..05d072457 100644 --- a/man/epidist_gen_posterior_epred.Rd +++ b/man/epidist_gen_posterior_epred.Rd @@ -8,8 +8,12 @@ distribution for a latent model} epidist_gen_posterior_epred(family) } \arguments{ -\item{family}{Output of a call to \code{brms::brmsfamily()} with additional -information as provided by \code{.add_dpar_info()}} +\item{family}{A description of the response distribution and link function to +be used in the model. Every family function has a link argument allowing +users to specify the link function to be applied on the response variable. +If not specified, default links are used. For details of all supported +families see \code{\link[=brmsfamily]{brmsfamily()}}. Commonly used, such as \code{\link[=lognormal]{lognormal()}}, are also +reexported as part of \code{epidist}.} } \value{ A function that takes a prep argument from brms and returns a matrix @@ -19,7 +23,7 @@ per observation. \description{ This function creates a function that calculates the expected value of the posterior predictive distribution for a latent model. The returned function -takes a prep argument (from brms) and returns posterior expected values. +takes a \code{prep} argument (from brms) and returns posterior expected values. This is used internally by \code{\link[brms:posterior_epred.brmsfit]{brms::posterior_epred()}} to calculate expected values for latent models. } diff --git a/man/epidist_gen_posterior_predict.Rd b/man/epidist_gen_posterior_predict.Rd index fb60d3bd0..b572d4c6e 100644 --- a/man/epidist_gen_posterior_predict.Rd +++ b/man/epidist_gen_posterior_predict.Rd @@ -8,13 +8,17 @@ latent model} epidist_gen_posterior_predict(family) } \arguments{ -\item{family}{Output of a call to \code{brms::brmsfamily()} with additional -information as provided by \code{.add_dpar_info()}} +\item{family}{A description of the response distribution and link function to +be used in the model. Every family function has a link argument allowing +users to specify the link function to be applied on the response variable. +If not specified, default links are used. For details of all supported +families see \code{\link[=brmsfamily]{brmsfamily()}}. Commonly used, such as \code{\link[=lognormal]{lognormal()}}, are also +reexported as part of \code{epidist}.} } \value{ -A function that takes a prep argument from brms and returns a matrix -of posterior predictions, with one row per posterior draw and one column -per observation. The prep object must have the following variables: +A function that takes a \code{prep} argument from brms and returns a +matrix of posterior predictions, with one row per posterior draw and one +column per observation. The \code{prep} object must have the following variables: \itemize{ \item \code{vreal1}: relative observation time \item \code{vreal2}: primary event window @@ -24,7 +28,7 @@ per observation. The prep object must have the following variables: \description{ This function creates a function that draws from the posterior predictive distribution for a latent model using \code{\link[primarycensored:rprimarycensored]{primarycensored::rpcens()}} to handle -censoring and truncation. The returned function takes a prep argument from +censoring and truncation. The returned function takes a \code{prep} argument from \code{brms} and returns posterior predictions. This is used internally by \code{\link[brms:posterior_predict.brmsfit]{brms::posterior_predict()}} to generate predictions for latent models. } diff --git a/vignettes/ebola.Rmd b/vignettes/ebola.Rmd index 1dbf76016..bd1c199fe 100644 --- a/vignettes/ebola.Rmd +++ b/vignettes/ebola.Rmd @@ -315,7 +315,7 @@ ranef(fit_sex_district) To go further than summaries of the fitted model, we recommend using the `tidybayes` package. For example, to obtain the posterior expectation of the delay distribution, under no censoring or truncation, we may use the `modelr::data_grid()` function in combination with the `tidybayes::add_epred_draws()` function. -The `tidybayes::add_epred_draws()` function uses the (internal) `posterior_epred_latent_lognormal()` function implemented in `epidist` for the `latent_lognormal` custom `brms` family. +The `tidybayes::add_epred_draws()` function uses the `epidist_gen_posterior_predict()` function to generate a posterior prediction function for the `lognormal()` distribution. In Figure \@ref(fig:epred) we show the posterior expectation of the delay distribution for each of the three fitted models. Figure \@ref(fig:epred)B illustrates the higher mean of men as compared with women. diff --git a/vignettes/epidist.Rmd b/vignettes/epidist.Rmd index 40a561d41..399652a84 100644 --- a/vignettes/epidist.Rmd +++ b/vignettes/epidist.Rmd @@ -290,6 +290,9 @@ fit <- epidist( ) ``` +**Note that here we use the default `rstan` backend but we generally recommend using the `cmdstanr` backend for faster sampling and additional features.** +**This can be set using `backend = "cmdstanr"` after following the installing CmdStan instructions in the README.** + The `fit` object is a `brmsfit` object containing MCMC samples from each of the parameters in the model, shown in the table below. Many of these parameters (e.g. `swindow` and `pwindow`) are the so called latent variables, and have lengths corresponding to the `sample_size`. @@ -307,9 +310,6 @@ data.frame( Users familiar with Stan and `brms`, can work with `fit` directly. Any tool that supports `brms` fitted model objects will be compatible with `fit`. -**Note that here we use the default `rstan` backend but we generally recommend using the `cmdstanr` backend for faster sampling and additional features.** -**This can be set using `backend = "cmdstanr"` after following the installing CmdStan instructions in the README.** - For example, we can use the built in `summary()` function ```{r} @@ -341,7 +341,7 @@ ggplot() + data = data.frame(x = c(0, 30)), aes(x = x), fun = dlnorm, - linewidth = 2, + linewidth = 1.5, args = list( meanlog = secondary_dist[["mu"]], sdlog = secondary_dist[["sigma"]]