Thanks as always to Kit Baum, a new package -bayesinits- is now
available from the SSC archive.
To install:
. ssc install bayesinits
. help bayesinits
-bayesinits- (version 1.4.1) generates domain-aware, randomized,
reproducible initial values for one or more Markov chains, for use
with -bayesmh-'s init#() options. Running multiple chains from
dispersed, feasible starting points is the basis of the Gelman-Rubin
convergence diagnostic (-bayesstats grubin-), but constructing such
values by hand is tedious and error-prone when parameters are
constrained: a carelessly drawn variance may be nonpositive, a
correlation may fall outside (-1, 1), and mixture weights may fail to
sum to one.
-bayesinits- samples on domain-appropriate working scales and
transforms back, so every chain starts at a valid point. Six parameter
domains are supported:
unconstrained() normal draws on the raw scale
positive() log-scale draws, exponentiated
corr() Fisher-z draws, mapped back via tanh()
fisherz() draws directly on the z scale
phi() angular parameters in (0, pi)
simplex() normalized log-normal weights (positive, sum 1)
Chain-specific strings are returned in r(init1), ..., r(initK),
together with composite strings r(init_all) and r(init_all_wrap) that
can be passed to -bayesmh- verbatim. Default centers and scales may be
overridden per parameter via center() and scale(), which accept
numeric expressions. The mle option derives centers and scales
automatically from the estimation results in memory: centers at the
point estimates, scales as inflate() (default 2) times the
delta-method standard error on the working scale, with map() linking
Bayesian parameter names to e(b) coefficient names. seed() and
rngstream() provide reproducibility and independent streams for
parallel chains.
A typical workflow:
. bayesinits, nchains(3) unconstrained(b0 b1) positive(sig2)
. bayesmh y = ({b0} + {b1}*x), likelihood(normal({sig2})) ///
prior({b0}, normal(0,100)) prior({b1}, normal(0,100)) ///
prior({sig2}, igamma(0.01,0.01)) ///
nchains(3) `r(init_all)' rseed(12345)
. bayesstats grubin
Or, anchored at a preliminary fit with threefold inflation:
. regress y x
. bayesinits, nchains(4) unconstrained(b0 b1) positive(sig2) ///
mle map(b0=_cons b1=x) inflate(3) center(sig2=e(rmse)^2)
-bayesinits- requires Stata 15 or later. It is particularly useful for
hierarchical models with correlation structures (for example,
bivariate diagnostic test accuracy meta-analysis), mixture models, and
any setting where dispersed multi-chain runs and credible convergence
assessment matter.
Comments, bug reports, and suggestions are welcome.
Ben Adarkwa Dwamena, MD
Division of Nuclear Medicine and Molecular Imaging
Department of Radiology, University of Michigan
[email protected]
available from the SSC archive.
To install:
. ssc install bayesinits
. help bayesinits
-bayesinits- (version 1.4.1) generates domain-aware, randomized,
reproducible initial values for one or more Markov chains, for use
with -bayesmh-'s init#() options. Running multiple chains from
dispersed, feasible starting points is the basis of the Gelman-Rubin
convergence diagnostic (-bayesstats grubin-), but constructing such
values by hand is tedious and error-prone when parameters are
constrained: a carelessly drawn variance may be nonpositive, a
correlation may fall outside (-1, 1), and mixture weights may fail to
sum to one.
-bayesinits- samples on domain-appropriate working scales and
transforms back, so every chain starts at a valid point. Six parameter
domains are supported:
unconstrained() normal draws on the raw scale
positive() log-scale draws, exponentiated
corr() Fisher-z draws, mapped back via tanh()
fisherz() draws directly on the z scale
phi() angular parameters in (0, pi)
simplex() normalized log-normal weights (positive, sum 1)
Chain-specific strings are returned in r(init1), ..., r(initK),
together with composite strings r(init_all) and r(init_all_wrap) that
can be passed to -bayesmh- verbatim. Default centers and scales may be
overridden per parameter via center() and scale(), which accept
numeric expressions. The mle option derives centers and scales
automatically from the estimation results in memory: centers at the
point estimates, scales as inflate() (default 2) times the
delta-method standard error on the working scale, with map() linking
Bayesian parameter names to e(b) coefficient names. seed() and
rngstream() provide reproducibility and independent streams for
parallel chains.
A typical workflow:
. bayesinits, nchains(3) unconstrained(b0 b1) positive(sig2)
. bayesmh y = ({b0} + {b1}*x), likelihood(normal({sig2})) ///
prior({b0}, normal(0,100)) prior({b1}, normal(0,100)) ///
prior({sig2}, igamma(0.01,0.01)) ///
nchains(3) `r(init_all)' rseed(12345)
. bayesstats grubin
Or, anchored at a preliminary fit with threefold inflation:
. regress y x
. bayesinits, nchains(4) unconstrained(b0 b1) positive(sig2) ///
mle map(b0=_cons b1=x) inflate(3) center(sig2=e(rmse)^2)
-bayesinits- requires Stata 15 or later. It is particularly useful for
hierarchical models with correlation structures (for example,
bivariate diagnostic test accuracy meta-analysis), mixture models, and
any setting where dispersed multi-chain runs and credible convergence
assessment matter.
Comments, bug reports, and suggestions are welcome.
Ben Adarkwa Dwamena, MD
Division of Nuclear Medicine and Molecular Imaging
Department of Radiology, University of Michigan
[email protected]
