| Title: | Bayesian Meta-Analytic Mixture Modelling |
|---|---|
| Description: | Fits Bayesian random-effects meta-analytic mixture models in which the true study effects are drawn from a mixture of normal components with component-specific means and between-study heterogeneity. The mixture model can be combined with step-function selection models (Vevea and Hedges, 1995, <doi:10.1007/BF02294384>) to adjust for publication bias, either with one selection model for all studies or with separate selection models for groups of studies. Models are estimated by Hamiltonian Monte Carlo using 'rstan' and returned as 'stanfit' objects, so posterior summaries, posterior predictive checks and marginal likelihoods for model comparison (via 'bridgesampling') are readily available. Also includes a data simulator and the nudging meta-analysis data of Mertens et al. (2022) <doi:10.1073/pnas.2107346118>. |
| Authors: | Maximilian Maier [aut, cre] (ORCID: <https://orcid.org/0000-0002-9873-6096>), Trustees of Columbia University [cph] |
| Maintainer: | Maximilian Maier <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.1.0 |
| Built: | 2026-10-21 17:52:12 UTC |
| Source: | https://github.com/cran/metamixr |
Effect sizes from the meta-analysis of choice-architecture ("nudge") interventions by Mertens, Herberz, Hahnel and Brosch (2022). The data set was the basis of the debate about publication bias in the nudging literature (Maier et al., 2022) and serves as the real-data example of this package. Each row is one effect size; several effect sizes may stem from the same study and several studies from the same publication.
mertens_nudgemertens_nudge
A data frame with 447 rows and 24 columns. The columns needed for
model fitting are cohens_d and variance_d. The variable descriptions
below are taken from the codebook of Mertens et al. (2022), available at
https://osf.io/49wtf/files/dvg3w.
publication_idPublication ID.
study_idStudy ID. Studies were defined based on independent samples; where an experiment compared multiple intervention conditions to the same control condition, a single study ID was assigned.
es_idEffect size ID.
referenceDescriptive reference of publication (i.e., authors + year of publication).
titleTitle of publication.
yearYear of publication.
locationGeographical location of intervention: 0 = outside United States, 1 = inside United States.
domainBehavioral domain of intervention: "health", "food",
"environment", "finance", "pro-social" or "other" (domains not
covered by the other categories).
intervention_categoryIntervention category based on the
taxonomy by Münscher et al. (2016): "information" (decision
information), "structure" (decision structure) or "assistance"
(decision assistance).
intervention_techniqueIntervention technique based on the
taxonomy by Münscher et al. (2016): "translation" (translate
information), "visibility" (make information visible),
"social_reference" (provide social reference point), "default"
(change choice default), "effort" (change option-related effort),
"composition" (change range or composition of options),
"consequence" (change option consequences), "reminder" (provide
reminders) or "commitment" (facilitate commitment).
type_experimentType of experiment as defined by Harrison and
List (2004): "conventional_lab", "artefactual_field",
"framed_field" or "natural_field".
populationTarget population of intervention: 0 = children and/or adolescents, 1 = adults.
n_studyOverall sample size of study.
n_comparisonSample size of control + intervention condition.
n_controlSample size of control condition.
n_interventionSample size of intervention condition.
binary_outcomeScale of outcome variable: 0 = continuous, 1 = binary.
mean_controlMean of outcome variable in control condition.
sd_controlStandard deviation of outcome variable in control condition.
mean_interventionMean of outcome variable in intervention condition.
sd_interventionStandard deviation of outcome variable in intervention condition.
cohens_dExtracted effect size of intervention (Cohen's d).
variance_dVariance around extracted effect size of intervention.
approximationApproximation involved in effect size extraction: 0 = no, 1 = yes.
Mertens, S., Herberz, M., Hahnel, U. J. J., & Brosch, T. (2022). The effectiveness of nudging: A meta-analysis of choice architecture interventions across behavioral domains. Proceedings of the National Academy of Sciences, 119(1), e2107346118. doi:10.1073/pnas.2107346118
Harrison, G. W., & List, J. A. (2004). Field experiments. Journal of Economic Literature, 42(4), 1009–1055.
Maier, M., Bartoš, F., Stanley, T. D., Shanks, D. R., Harris, A. J. L., & Wagenmakers, E.-J. (2022). No evidence for nudging after adjusting for publication bias. Proceedings of the National Academy of Sciences, 119(31), e2200300119. doi:10.1073/pnas.2200300119
Münscher, R., Vetter, M., & Scheuerle, T. (2016). A review and taxonomy of choice architecture techniques. Journal of Behavioral Decision Making, 29(5), 511–524. doi:10.1002/bdm.1897
data(mertens_nudge) str(mertens_nudge[, c("domain", "cohens_d", "variance_d")]) # z-statistics pile up just above the two-sided significance threshold z <- mertens_nudge$cohens_d / sqrt(mertens_nudge$variance_d) hist(z, breaks = 60, main = "Nudging meta-analysis", xlab = "z") abline(v = qnorm(0.975), lty = 2)data(mertens_nudge) str(mertens_nudge[, c("domain", "cohens_d", "variance_d")]) # z-statistics pile up just above the two-sided significance threshold z <- mertens_nudge$cohens_d / sqrt(mertens_nudge$variance_d) hist(z, breaks = 60, main = "Nudging meta-analysis", xlab = "z") abline(v = qnorm(0.975), lty = 2)
Fits a Bayesian random-effects meta-analysis in which the true study effects
arise from a mixture of M normal components, each with its own mean and
between-study heterogeneity (Maier, 2026). With M = 1 the model is the
ordinary random-effects meta-analysis. No publication bias adjustment is
applied; see sel_mix() and sel_flexpb() for the corresponding selection
models.
re_mix( y, sd, M, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, prior_only = FALSE, ... )re_mix( y, sd, M, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, prior_only = FALSE, ... )
y |
numeric vector of observed effect sizes. |
sd |
numeric vector of the standard errors of |
M |
number of mixture components (a positive integer). |
mu_sd |
standard deviation of the normal prior on the component means
(used by both mean-prior settings: as the prior on every ordered mean when
|
tau_sd |
standard deviation of the half-normal prior on the heterogeneity
SD |
mean_prior |
prior family for the component means. |
tau_prior |
prior family for the heterogeneity SD |
gap_meanlog |
log-scale location of the repulsive lognormal gap prior
(only used when |
gap_sdlog |
log-scale standard deviation of the repulsive lognormal gap
prior (only used when |
gap_min |
hard lower floor on the gaps between adjacent means: each gap is
constrained to be at least |
tau_alpha |
shape of the inverse-gamma prior on |
tau_beta |
scale of the inverse-gamma prior on |
prior_only |
if |
... |
further arguments passed to |
See Maier (2026) for details of the model and the priors.
An object of class stanfit as returned by
rstan::sampling(). It contains posterior draws of
muvector of length M: the ordered component means.
tauvector of length M: the between-study standard deviation of
each component.
thetavector of length M: the mixing weights.
mu1, mu_gap
the underlying parameterisation of the means (the
first mean and the M - 1 gaps between adjacent means).
log_taulog(tau), convenient for funnel diagnostics.
posterior_probsK by M matrix: the posterior probability that
study i belongs to component m.
y_rep, sd_rep
posterior predictive draws of K effect sizes and
the standard errors they were generated with. For each draw a component
is sampled from theta, a study is sampled with probability
proportional to its responsibility for that component, and a new effect
size is generated using that study's standard error.
Use e.g. print(fit, pars = c("mu", "tau", "theta")), summary(fit) or
rstan::extract() to access them.
Maier, M. (2026). Addressing heterogeneity with Bayesian meta-analytic mixture modelling. PsyArXiv. doi:10.31234/osf.io/nkyqm_v2
Maier, M., Bartoš, F., & Wagenmakers, E.-J. (2023). Robust Bayesian meta-analysis: Addressing publication bias with model-averaging. Psychological Methods, 28(1), 107–122. doi:10.1037/met0000405
Stan Development Team (2024). RStan: the R interface to Stan. https://mc-stan.org/
sel_mix() and sel_flexpb() for the same mixture model combined
with step-function selection models; sim_mix() to simulate data from
the model.
Other model fitting functions:
sel_flexpb(),
sel_mix()
# Simulate 40 studies from a two-component mixture without publication bias set.seed(1) dat <- sim_mix(K = 40, M = 2, mu = c(0, 0.6), tau = c(0.1, 0.1), steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. fit <- re_mix(dat$y, dat$sds, M = 2, chains = 1, iter = 500, refresh = 0) print(fit, pars = c("mu", "tau", "theta")) # Posterior probability that each study belongs to the second component resp <- rstan::extract(fit, pars = "posterior_probs")$posterior_probs head(colMeans(resp)[, 2]) # Sample from the prior only (no data enter the model) prior <- re_mix(dat$y, dat$sds, M = 2, prior_only = TRUE, chains = 1, iter = 500, refresh = 0) print(prior, pars = c("mu", "tau"))# Simulate 40 studies from a two-component mixture without publication bias set.seed(1) dat <- sim_mix(K = 40, M = 2, mu = c(0, 0.6), tau = c(0.1, 0.1), steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. fit <- re_mix(dat$y, dat$sds, M = 2, chains = 1, iter = 500, refresh = 0) print(fit, pars = c("mu", "tau", "theta")) # Posterior probability that each study belongs to the second component resp <- rstan::extract(fit, pars = "posterior_probs")$posterior_probs head(colMeans(resp)[, 2]) # Sample from the prior only (no data enter the model) prior <- re_mix(dat$y, dat$sds, M = 2, prior_only = TRUE, chains = 1, iter = 500, refresh = 0) print(prior, pars = c("mu", "tau"))
Fits the meta-analytic mixture model with step-function selection of
sel_mix(), but estimates a separate selection model (vector of selection
weights omega) for each group of studies, e.g. for several meta-analyses
that are analysed jointly, or for published versus preregistered studies.
The mixture components (means, heterogeneity and mixing weights) are shared
across groups. Groups that are known to be free of publication bias can be
declared via unselected, in which case their selection weights are fixed
at 1.
sel_flexpb( y, sd, grp, M, steps = c(0.9, 0.95), one_sided = TRUE, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, unselected = NULL, prior_only = FALSE, ... )sel_flexpb( y, sd, grp, M, steps = c(0.9, 0.95), one_sided = TRUE, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, unselected = NULL, prior_only = FALSE, ... )
y |
numeric vector of observed effect sizes. |
sd |
numeric vector of the standard errors of |
grp |
integer vector of the same length as |
M |
number of mixture components (a positive integer). |
steps |
cutoffs of the p-value intervals of the selection model, given as
cumulative probabilities of the standard normal distribution (i.e.
|
one_sided |
if |
mu_sd |
standard deviation of the normal prior on the component means
(used by both mean-prior settings: as the prior on every ordered mean when
|
tau_sd |
standard deviation of the half-normal prior on the heterogeneity
SD |
mean_prior |
prior family for the component means. |
tau_prior |
prior family for the heterogeneity SD |
gap_meanlog |
log-scale location of the repulsive lognormal gap prior
(only used when |
gap_sdlog |
log-scale standard deviation of the repulsive lognormal gap
prior (only used when |
gap_min |
hard lower floor on the gaps between adjacent means: each gap is
constrained to be at least |
tau_alpha |
shape of the inverse-gamma prior on |
tau_beta |
scale of the inverse-gamma prior on |
unselected |
integer vector of group indices assumed NOT subject to
publication bias: their selection weights |
prior_only |
if |
... |
further arguments passed to |
See Maier (2026) for details of the model and the priors.
An object of class stanfit as returned by
rstan::sampling(). It contains posterior draws of
muvector of length M: the ordered component means.
tauvector of length M: the between-study standard deviation of
each component.
thetavector of length M: the mixing weights among the
published studies.
theta_preselectionvector of length M: the mixing weights
corrected for selection, i.e. the estimated proportions of the
components among all conducted studies rather than among the
published ones.
omegaG by length(steps) + 1 matrix: one row of relative
publication probabilities of the p-value intervals per group,
non-decreasing with the last element fixed at 1. Rows of groups listed
in unselected are constant at 1.
omega_rawthe underlying simplex of each group,
omega[g, ] = cumsum(omega_raw[g, ]).
avg_omegavector of length M: the average selection weight of
the studies attributed to each component, i.e. how strongly each
component is affected by selection.
mu1, mu_gap
the underlying parameterisation of the means (the
first mean and the M - 1 gaps between adjacent means).
log_taulog(tau), convenient for funnel diagnostics.
posterior_probsK by M matrix: the posterior probability that
study i belongs to component m.
y_rep, sd_rep
posterior predictive draws of K effect sizes and
the standard errors they were generated with.
Use e.g. print(fit, pars = c("mu", "tau", "omega")), summary(fit) or
rstan::extract() to access them.
Maier, M. (2026). Addressing heterogeneity with Bayesian meta-analytic mixture modelling. PsyArXiv. doi:10.31234/osf.io/nkyqm_v2
Vevea, J. L., & Hedges, L. V. (1995). A general linear model for estimating effect size in the presence of publication bias. Psychometrika, 60(3), 419–435. doi:10.1007/BF02294384
sel_mix() for a single selection model, re_mix() for no
selection, sim_mix() to simulate data.
Other model fitting functions:
re_mix(),
sel_mix()
# Two sets of studies on the same effect: a published literature subject to # selection (group 1) and preregistered replications that are not (group 2) set.seed(1) lit <- sim_mix(K = 15, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.2, 0.5, 1), one_sided = TRUE) rep <- sim_mix(K = 15, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) y <- c(lit$y, rep$y) sd <- c(lit$sds, rep$sds) grp <- rep(1:2, each = 15) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. Group 2 is known to # be unselected, so a selection model is estimated for group 1 only. fit <- sel_flexpb(y, sd, grp, M = 1, steps = c(0.95, 0.975), unselected = 2, chains = 1, iter = 400, refresh = 0) print(fit, pars = c("mu", "tau", "omega"))# Two sets of studies on the same effect: a published literature subject to # selection (group 1) and preregistered replications that are not (group 2) set.seed(1) lit <- sim_mix(K = 15, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.2, 0.5, 1), one_sided = TRUE) rep <- sim_mix(K = 15, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) y <- c(lit$y, rep$y) sd <- c(lit$sds, rep$sds) grp <- rep(1:2, each = 15) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. Group 2 is known to # be unselected, so a selection model is estimated for group 1 only. fit <- sel_flexpb(y, sd, grp, M = 1, steps = c(0.95, 0.975), unselected = 2, chains = 1, iter = 400, refresh = 0) print(fit, pars = c("mu", "tau", "omega"))
Fits the random-effects meta-analytic mixture model of re_mix() combined
with a step-function (weight-function) selection model
to adjust for publication bias (Maier, 2026).
One common selection model is estimated for all studies.
Use sel_flexpb() to estimate separate selection models for groups of
studies.
sel_mix( y, sd, M, steps = c(0.9, 0.95), one_sided = TRUE, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, prior_only = FALSE, ... )sel_mix( y, sd, M, steps = c(0.9, 0.95), one_sided = TRUE, mu_sd = 1, tau_sd = 0.2, mean_prior = c("normal", "gap"), tau_prior = c("half_normal", "inv_gamma"), gap_meanlog = log(0.2), gap_sdlog = 0.3, gap_min = 0, tau_alpha = 2.5, tau_beta = 0.15, prior_only = FALSE, ... )
y |
numeric vector of observed effect sizes. |
sd |
numeric vector of the standard errors of |
M |
number of mixture components (a positive integer). |
steps |
cutoffs of the p-value intervals of the selection model, given as
cumulative probabilities of the standard normal distribution (i.e.
|
one_sided |
if |
mu_sd |
standard deviation of the normal prior on the component means
(used by both mean-prior settings: as the prior on every ordered mean when
|
tau_sd |
standard deviation of the half-normal prior on the heterogeneity
SD |
mean_prior |
prior family for the component means. |
tau_prior |
prior family for the heterogeneity SD |
gap_meanlog |
log-scale location of the repulsive lognormal gap prior
(only used when |
gap_sdlog |
log-scale standard deviation of the repulsive lognormal gap
prior (only used when |
gap_min |
hard lower floor on the gaps between adjacent means: each gap is
constrained to be at least |
tau_alpha |
shape of the inverse-gamma prior on |
tau_beta |
scale of the inverse-gamma prior on |
prior_only |
if |
... |
further arguments passed to |
See Maier (2026) for details of the model and the priors.
An object of class stanfit as returned by
rstan::sampling(). It contains posterior draws of
muvector of length M: the ordered component means.
tauvector of length M: the between-study standard deviation of
each component.
thetavector of length M: the mixing weights among the
published studies.
theta_preselectionvector of length M: the mixing weights
corrected for selection, i.e. the estimated proportions of the
components among all conducted studies rather than among the
published ones.
omegavector of length length(steps) + 1: the relative
publication probabilities of the p-value intervals, non-decreasing with
the last element fixed at 1.
omega_rawthe underlying simplex, omega = cumsum(omega_raw).
avg_omegavector of length M: the average selection weight of
the studies attributed to each component, i.e. how strongly each
component is affected by selection.
mu1, mu_gap
the underlying parameterisation of the means (the
first mean and the M - 1 gaps between adjacent means).
log_taulog(tau), convenient for funnel diagnostics.
posterior_probsK by M matrix: the posterior probability that
study i belongs to component m.
y_rep, sd_rep
posterior predictive draws of K effect sizes and
the standard errors they were generated with.
Use e.g. print(fit, pars = c("mu", "tau", "omega")), summary(fit) or
rstan::extract() to access them.
Maier, M. (2026). Addressing heterogeneity with Bayesian meta-analytic mixture modelling. PsyArXiv. doi:10.31234/osf.io/nkyqm_v2
Vevea, J. L., & Hedges, L. V. (1995). A general linear model for estimating effect size in the presence of publication bias. Psychometrika, 60(3), 419–435. doi:10.1007/BF02294384
Maier, M., Bartoš, F., & Wagenmakers, E.-J. (2023). Robust Bayesian meta-analysis: Addressing publication bias with model-averaging. Psychological Methods, 28(1), 107–122. doi:10.1037/met0000405
re_mix() for the model without selection, sel_flexpb() for
group-specific selection models, sim_mix() to simulate data.
Other model fitting functions:
re_mix(),
sel_flexpb()
# 30 studies with a true mean of 0.3 under one-sided selection: studies that # are not significant in the expected direction (two-sided p > .10) are # published with probability .2, marginally significant ones with .5 set.seed(1) dat <- sim_mix(K = 30, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.2, 0.5, 1), one_sided = TRUE) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. fit <- sel_mix(dat$y, dat$sds, M = 1, steps = c(0.95, 0.975), chains = 1, iter = 400, refresh = 0) print(fit, pars = c("mu", "tau", "omega"))# 30 studies with a true mean of 0.3 under one-sided selection: studies that # are not significant in the expected direction (two-sided p > .10) are # published with probability .2, marginally significant ones with .5 set.seed(1) dat <- sim_mix(K = 30, M = 1, mu = 0.3, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.2, 0.5, 1), one_sided = TRUE) # A single short chain keeps the example fast; use e.g. chains = 4 and # iter = 2000 (the rstan defaults) for a real analysis. fit <- sel_mix(dat$y, dat$sds, M = 1, steps = c(0.95, 0.975), chains = 1, iter = 400, refresh = 0) print(fit, pars = c("mu", "tau", "omega"))
Simulates K published studies from the random-effects mixture model of
re_mix() combined with the step-function selection model of sel_mix().
Studies are generated until K of them have passed the selection step, so
the returned data set has exactly K rows.
sim_mix( K, M, mu, tau, steps, weights, one_sided, thetas = c(rep(1/M, M)), N_low = 25, N_high = 500, N_shape = 2, N_scale = 58 )sim_mix( K, M, mu, tau, steps, weights, one_sided, thetas = c(rep(1/M, M)), N_low = 25, N_high = 500, N_shape = 2, N_scale = 58 )
K |
number of (published) studies to return. |
M |
number of mixture components. |
mu |
numeric vector of length |
tau |
numeric vector of length |
steps |
cutoffs of the p-value intervals, in the same format as in
|
weights |
publication probabilities of the p-value intervals, a vector
of length |
one_sided |
if |
thetas |
mixing weights of the components; defaults to equal weights. |
N_low |
lower bound on the total sample size of a primary study. |
N_high |
upper bound on the total sample size of a primary study. |
N_shape |
shape ( |
N_scale |
scale of the negative binomial distribution used to generate
sample sizes; the success probability is |
See Maier (2026) for details of the model and the simulation design.
A data frame with K rows and the columns y (observed effect
size) and sds (its standard error), ready to be passed to re_mix(),
sel_mix() or sel_flexpb().
Maier, M. (2026). Addressing heterogeneity with Bayesian meta-analytic mixture modelling. PsyArXiv. doi:10.31234/osf.io/nkyqm_v2
Maier, M., Bartoš, F., & Wagenmakers, E.-J. (2023). Robust Bayesian meta-analysis: Addressing publication bias with model-averaging. Psychological Methods, 28(1), 107–122. doi:10.1037/met0000405
re_mix(), sel_mix(), sel_flexpb()
set.seed(1) # Two components (null effect and medium effect) without publication bias dat <- sim_mix(K = 100, M = 2, mu = c(0, 0.5), tau = c(0.05, 0.1), steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) head(dat) # Strong one-sided selection: studies that are not significant in the # expected direction are published with probability .1 only dat_pb <- sim_mix(K = 100, M = 1, mu = 0, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.1, 0.5, 1), one_sided = TRUE) mean(dat_pb$y / dat_pb$sds > qnorm(0.975))set.seed(1) # Two components (null effect and medium effect) without publication bias dat <- sim_mix(K = 100, M = 2, mu = c(0, 0.5), tau = c(0.05, 0.1), steps = c(0.95, 0.975), weights = c(1, 1, 1), one_sided = TRUE) head(dat) # Strong one-sided selection: studies that are not significant in the # expected direction are published with probability .1 only dat_pb <- sim_mix(K = 100, M = 1, mu = 0, tau = 0.1, steps = c(0.95, 0.975), weights = c(0.1, 0.5, 1), one_sided = TRUE) mean(dat_pb$y / dat_pb$sds > qnorm(0.975))