posterior_g_comp {beaver} | R Documentation |
Compute Posterior G-Computation Estimate
Description
Calculate the estimated effect for each observation at each dose and average over all observations. This function calculates the posterior marginal treatment effect at each dose.
Usage
posterior_g_comp(
x,
doses = attr(x, "doses"),
reference_dose = NULL,
prob = c(0.025, 0.975),
return_stats = TRUE,
return_samples = FALSE,
new_data = NULL,
reference_type = c("difference", "ratio")
)
Arguments
x |
an object output from |
doses |
doses at which to obtain the posterior. |
reference_dose |
dose to which to compare as either a difference or ratio. |
prob |
the percentiles of the posterior to calculate for each dose. |
return_stats |
logical indicating if the posterior mean and quantiles should be returned. |
return_samples |
logical indicating if posterior mean samples should be returned. |
new_data |
a dataframe containing all the variables used in the
covariate adjustments to the model used to obtain |
reference_type |
whether to provide the posterior of the difference or the ratio between each dose and the reference dose. |
Value
A list with the elements stats
and samples
. When using this
function with default settings, samples
is NULL and stats
is a
dataframe summarizing the posterior samples. stats
contains, at a
minimum, the columns "dose", "value", and variables corresponding to the
values passed in prob
("2.50%" and "97.50%" by default). When
return_stats
is set to FALSE
, stats
is NULL. When return_samples
is set to TRUE
, samples
is a dataframe with the posterior samples for
each iteration of the MCMC.
- When x is of class 'beaver_mcmc_bma':
-
The dataframe will have, at a minimum, the columns "iter" and "model", indicating the MCMC iteration and the model that was used in the calculations, as well as the columns "dose" and "value". The functions used for each model are defined within the
model_negbin_XYZ()
functions and used in thebeaver_mcmc()
function. - When x is of class 'beaver_mcmc':
-
The dataframe will have, at a minimum, the column "iter", indicating the MCMC iteration, as well as the columns "dose" and "value". The functions used for each model are defined within the
model_negbin_XYZ()
functions and used in therun_mcmc()
function.
See Also
Other posterior calculations:
beaver_mcmc()
,
posterior.beaver_mcmc_bma()
,
posterior.beaver_mcmc()
,
pr_eoi_g_comp()
,
pr_eoi()
Examples
# The {beaver} package, by definition, performs MCMC for multiple models.
# Even with a small number of chains/burn-ins/samples, a minimally illustrative
# example requires >5s to run.
library(dplyr)
# No covariates----
set.seed(100)
df <- data_negbin_emax(
n_per_arm = 10,
doses = 0:3,
b1 = 0,
b2 = 2.5,
b3 = 0.5,
ps = 0.75
)
df %>%
group_by(dose) %>%
summarize(
mean = mean(response),
se = sd(response) / sqrt(n()),
.groups = "drop"
)
mcmc <- beaver_mcmc(
emax = model_negbin_emax(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 1.5,
sigma_b3 = 3,
w_prior = 1 / 4
),
linear = model_negbin_linear(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
w_prior = 1 / 4
),
quad = model_negbin_quad(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 1.5,
sigma_b3 = 3,
w_prior = 1 / 4
),
exp = model_negbin_exp(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 0,
sigma_b3 = 3,
w_prior = 1 / 4
),
formula = ~ 1,
data = df,
n_iter = 1e2,
n_chains = 1,
quiet = TRUE
)
mcmc$w_post
draws <- try(draws(mcmc)) #draws() is intended for single model fits only
draws_emax <- draws(mcmc$models$emax$mcmc)
draws_linear <- draws(mcmc$models$linear$mcmc)
draws_quad <- draws(mcmc$models$quad$mcmc)
draws_exp <- draws(mcmc$models$exp$mcmc)
post <- posterior(
mcmc,
contrast = matrix(1, 1, 1),
doses = 0:3,
reference_dose = 0,
reference_type = "difference"
)
pr_eoi(
mcmc,
eoi = c(5, 8),
contrast = matrix(1, 1, 1),
reference_dose = 0,
reference_type = "difference"
)
post_g_comp <- posterior_g_comp(
mcmc,
new_data = df,
reference_dose = 0,
reference_type = "difference"
)
pr_eoi_g_comp(
mcmc,
eoi = c(5, 8),
new_data = df,
reference_dose = 0,
reference_type = "difference"
)
plot(mcmc, contrast = matrix(1, 1, 1))
# With covariates----
set.seed(1000)
x <-
data.frame(
gender = factor(sample(c("F", "M"), 40, replace = TRUE))
) %>%
model.matrix(~ gender, data = .)
df_cov <-
data_negbin_emax(
n_per_arm = 10,
doses = 0:3,
b1 = c(0, 0.5),
b2 = 2.5,
b3 = 0.5,
ps = 0.75,
x = x
) %>%
mutate(
gender = case_when(
genderM == 1 ~ "M",
TRUE ~ "F"
),
gender = factor(gender)
) %>%
select(subject, dose, gender, response)
df_cov %>%
group_by(dose, gender) %>%
summarize(
mean = mean(response),
se = sd(response) / sqrt(n()),
.groups = "drop"
)
mcmc_cov <- beaver_mcmc(
emax = model_negbin_emax(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 1.5,
sigma_b3 = 3,
w_prior = 1 / 4
),
linear = model_negbin_linear(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
w_prior = 1 / 4
),
quad = model_negbin_quad(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 1.5,
sigma_b3 = 3,
w_prior = 1 / 4
),
exp = model_negbin_exp(
mu_b1 = 0,
sigma_b1 = 10,
mu_b2 = 0,
sigma_b2 = 10,
mu_b3 = 0,
sigma_b3 = 3,
w_prior = 1 / 4
),
formula = ~ gender,
data = df_cov,
n_iter = 1e2,
n_chains = 1,
quiet = TRUE
)
mcmc_cov$w_post
draws_cov <- try(draws(mcmc_cov)) #draws() is intended for single model fits only
draws_cov_emax <- draws(mcmc_cov$models$emax$mcmc)
draws_cov_linear <- draws(mcmc_cov$models$linear$mcmc)
draws_cov_quad <- draws(mcmc_cov$models$quad$mcmc)
draws_cov_exp <- draws(mcmc_cov$models$exp$mcmc)
post_cov <- posterior(
mcmc_cov,
contrast = matrix(c(1, 1, 0, 1), 2, 2),
doses = 0:3,
reference_dose = 0,
reference_type = "difference"
)
pr_eoi(
mcmc_cov,
eoi = c(5, 8),
contrast = matrix(c(1, 1, 0, 1), 2, 2),
reference_dose = 0,
reference_type = "difference"
)
post_g_comp_cov <- posterior_g_comp(
mcmc_cov,
new_data = df_cov,
reference_dose = 0,
reference_type = "difference"
)
pr_eoi_g_comp(
mcmc_cov,
eoi = c(5, 8),
new_data = df_cov,
reference_dose = 0,
reference_type = "difference"
)
plot(mcmc_cov, new_data = df_cov, type = "g-comp")