Calculate posterior mean (and quantiles for specific doses for each MCMC iteration of the model.

posterior(
  x,
  doses,
  times,
  probs,
  reference_dose,
  predictive,
  return_samples,
  iter,
  return_stats
)

# S3 method for dreamer_mcmc
posterior(
  x,
  doses = attr(x, "doses"),
  times = attr(x, "times"),
  probs = c(0.025, 0.975),
  reference_dose = NULL,
  predictive = 0,
  return_samples = FALSE,
  iter = NULL,
  return_stats = TRUE
)

# S3 method for dreamer_mcmc_independent
posterior(
  x,
  doses = attr(x, "doses"),
  times = attr(x, "times"),
  probs = c(0.025, 0.975),
  reference_dose = NULL,
  predictive = 0,
  return_samples = FALSE,
  iter = NULL,
  return_stats = TRUE
)

# S3 method for dreamer_bma
posterior(
  x,
  doses = x$doses,
  times = x$times,
  probs = c(0.025, 0.975),
  reference_dose = NULL,
  predictive = 0,
  return_samples = FALSE,
  iter = NULL,
  return_stats = TRUE
)

Arguments

x

output from a call to dreamer_mcmc.

doses

doses at which to estimate posterior quantities.

times

a vector of times at which to calculate the posterior response (for longitudinal models only).

probs

quantiles of the posterior to be calculated.

reference_dose

the dose at which to adjust the posterior plot. Specifying a dose returns the plot of pr(trt_dose - trt_reference_dose | data).

predictive

An integer. If greater than 0, the return values will be from the predictive distribution of the mean of predictive observations. If 0 (default), the posterior on the dose response mean is returned.

return_samples

logical indicating if the weighted raw MCMC samples from the Bayesian model averaging used to calculate the mean and quantiles should be returned.

iter

an index on which iterations of the MCMC should be used in the calculations. By default, all MCMC iterations are used.

return_stats

logical indicating whether or not the posterior statistics should be calculated.

Value

A named list with the following elements:

  • stats: a tibble the dose, posterior mean, and posterior quantiles.

  • samps: the weighted posterior samples. Only returned if return_samples = TRUE.

Methods (by class)

  • posterior(dreamer_mcmc): posterior summary for linear model.

  • posterior(dreamer_mcmc_independent): posterior summary for independent model.

  • posterior(dreamer_bma): posterior summary for Bayesian model averaging fit.

Examples

set.seed(888)
data <- dreamer_data_linear(
  n_cohorts = c(20, 20, 20),
  dose = c(0, 3, 10),
  b1 = 1,
  b2 = 3,
  sigma = 5
)

# Bayesian model averaging
output <- dreamer_mcmc(
 data = data,
 n_adapt = 1e3,
 n_burn = 1e3,
 n_iter = 1e4,
 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-04-01 17:24:44.889
#> 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-04-01 17:24:44.957
#> total : 0.07 secs
#> mod_quad
#> start : 2024-04-01 17:24:44.957
#> 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-04-01 17:24:45.030
#> total : 0.07 secs

posterior(output)
#> $stats
#> # A tibble: 3 × 4
#>    dose   mean `2.50%` `97.50%`
#>   <dbl>  <dbl>   <dbl>    <dbl>
#> 1     0  0.572  -0.803     1.95
#> 2     3  9.56    7.92     11.0 
#> 3    10 30.8    28.6      33.1 
#> 

# return posterior samples of the mean
post <- posterior(output, return_samples = TRUE)
head(post$samps)
#>   iter dose mean_response    model
#> 1    1    0     0.6844290 mod_quad
#> 2    6    0    -0.4293629 mod_quad
#> 3    8    0     1.5144070 mod_quad
#> 4   12    0     0.6374851 mod_quad
#> 5   13    0     1.5117433 mod_quad
#> 6   19    0     0.8239574 mod_quad

# from a single model
posterior(output$mod_quad)
#> $stats
#> # A tibble: 3 × 4
#>    dose   mean `2.50%` `97.50%`
#>   <dbl>  <dbl>   <dbl>    <dbl>
#> 1     0  0.612  -0.808     2.06
#> 2     3  9.37    7.22     11.5 
#> 3    10 30.9    28.6      33.3 
#> 

# posterior of difference of doses
posterior(output, reference_dose = 0)
#> $stats
#> # A tibble: 3 × 5
#>    dose reference_dose  mean `2.50%` `97.50%`
#>   <dbl>          <dbl> <dbl>   <dbl>    <dbl>
#> 1     0              0  0       0         0  
#> 2     3              0  8.98    7.16     10.4
#> 3    10              0 30.3    27.5      33.0
#>