Calculate the probability a dose being the smallest dose that has at least X% of the maximum efficacy.
pr_medx(
x,
doses = attr(x, "doses"),
ed,
greater = TRUE,
small_bound = 0,
time = NULL
)output from a call to dreamer_mcmc().
the doses for which pr(minimum effective X% dose) is to be calculated.
a number between 0 and 100 indicating the ed% dose that is being sought.
if TRUE, higher responses indicate better efficacy. If
FALSE, lower responses indicate better efficacy.`
the lower (upper) bound of the response variable
when greater = TRUE (FALSE). This is used to calculate the
ed% effect as ed / 100 * (effect_100 - small_bound) + small_bound.
the time (scalar) at which the Pr(MEDX) should be calculated.
A data frame with the following columns:
dose: numeric dose levels.
prob: Prob(EDX | data) for each dose. Note: these probabilities do
not necessarily sum to 1 because the EDX may not exist. In fact,
Pr(EDX does not exist | data) = 1 - sum(prob).
Obtaining the probability of a particular does being the
minimum efficacious dose achieving ed% efficacy is dependent on
the doses specified.
For a given MCMC sample of parameters, the 100% efficacy value is defined
as the highest efficacy of the doses specified. For each posterior draw
of MCMC parameters, the minimum ed% efficacious dose is defined as the
lowest dose what has at least ed% efficacy relative to the 100%
efficacy value.
The ed% effect is calculated as
ed / 100 * (effect_100 - small_bound) + small_bound where effect_100
is the largest mean response among doses for a given MCMC iteration.
set.seed(888)
data <- dreamer_data_linear(
n_cohorts = c(20, 20, 20),
dose = c(0, 3, 10),
b1 = 1,
b2 = .1,
sigma = 5
)
# Bayesian model averaging
output <- dreamer_mcmc(
data = data,
n_adapt = 1e3,
n_burn = 1e3,
n_iter = 1e3,
n_chains = 2,
silent = FALSE,
mod_linear = model_linear(
mu_b1 = 0,
sigma_b1 = 1,
mu_b2 = 0,
sigma_b2 = 1,
shape = 1,
rate = .001,
w_prior = 1 / 2
),
mod_quad = model_quad(
mu_b1 = 0,
sigma_b1 = 1,
mu_b2 = 0,
sigma_b2 = 1,
mu_b3 = 0,
sigma_b3 = 1,
shape = 1,
rate = .001,
w_prior = 1 / 2
)
)
#> mod_linear
#> start : 2024-12-19 14:43:33.027
#> Compiling model graph
#> Resolving undeclared variables
#> Allocating nodes
#> Graph information:
#> Observed stochastic nodes: 6
#> Unobserved stochastic nodes: 3
#> Total graph size: 50
#>
#> Initializing model
#>
#> finish: 2024-12-19 14:43:33.045
#> total : 0.02 secs
#> mod_quad
#> start : 2024-12-19 14:43:33.045
#> Compiling model graph
#> Resolving undeclared variables
#> Allocating nodes
#> Graph information:
#> Observed stochastic nodes: 6
#> Unobserved stochastic nodes: 4
#> Total graph size: 59
#>
#> Initializing model
#>
#> finish: 2024-12-19 14:43:33.069
#> total : 0.02 secs
pr_medx(output, ed = 80)
#> # A tibble: 3 × 3
#> ed dose prob
#> <dbl> <chr> <dbl>
#> 1 80 0 0.0605
#> 2 80 3 0.415
#> 3 80 10 0.522
# single model
pr_medx(output$mod_linear, ed = 80)
#> # A tibble: 3 × 3
#> ed dose prob
#> <dbl> <chr> <dbl>
#> 1 80 0 0.13
#> 2 80 3 0.018
#> 3 80 10 0.846