extract_posterior_summaries: Extract posterior summaries from a hierarchical model fit

In statdivlab/paramedic: Predicting Absolute and Relative Abundance by Modeling Efficiency to Derive Intervals and Concentrations

Extract posterior summaries from a hierarchical model fit

Description

Return point estimates and credible intervals for the true concentration, and point estimates and prediction intervals for estimated qPCR obtained through a Stan fit.

Usage

extract_posterior_summaries(
  stan_mod,
  stan_samps,
  taxa_of_interest,
  mult_num = 1,
  level = 0.95,
  interval_type = "wald"
)

Arguments

stan_mod	the model summary object from Stan.
stan_samps	the list of MCMC samples from Stan.
taxa_of_interest	the indices of the taxa for which point estimates and posterior summaries are desired.
mult_num	the number to multiply the resulting estimates and standard deviations by (defaults to 1).
level	the alpha level for prediction intervals (defaults to 0.95, for a nominal 95% prediction interval).
interval_type	the type of prediction interval desired (defaults to "wald", but "quantile" is also acceptable).

Details

A paramedic object is a list containing the following elements:

• estimates - the point estimates of qPCR (a matrix with dimension sample size by number of taxa).

• pred_intervals - predction intervals for qPCR (an array with dimension sample size by 2 by number of taxa).

• est_efficiency - point estimates for estimated varying efficiency, if varying efficiency was modeled (a vector of length number of taxa); otherwise, NA.

• efficiency_intervals - posterior level level\times100% confidence intervals for the true efficiency, if efficiency was modeled (a matrix of dimension number of taxa by 2); otherwise, NA.

Value

An object of class paramedic. See Details for more information

Examples

# load the package, read in example data
library("paramedic")
data(example_16S_data)
data(example_qPCR_data)

# run paramedic (with an extremely small number of iterations, for illustration only)
# on only the first 10 taxa
mod <- run_paramedic(W = example_16S_data[, 1:10], V = example_qPCR_data,
n_iter = 30, n_burnin = 25, 
n_chains = 1, stan_seed = 4747)

# get summary, samples
mod_summ <- rstan::summary(mod, probs = c(0.025, 0.975))$summary
mod_samps <- rstan::extract(mod$stan_fit)

# extract relevant summaries
summs <- extract_posterior_summaries(stan_mod = mod_summ, stan_samps = mod_samps, 
taxa_of_interest = 1:3,
mult_num = 1, level = 0.95, interval_type = "wald")

statdivlab/paramedic documentation built on May 3, 2024, 7:08 p.m.