This function sets up a spike-and-slab model for variable selection and model
choice in generalized additive models and samples its posterior. It uses a
blockwise Metropolis-within-Gibbs sampler and the redundant multiplicative
parameter expansion described in the reference. This routine is not meant to
be called directly by the user -- spikeSlabGAM
provides a
formula-based interface for specifying models and takes care of (most of) the
housekeeping. Sampling of the chains is done in parallel using package
parallel
. A "SOCK" cluster is set up under Windows to do so (and
closed after computations are done, I try to clean up after myself), see
makeCluster
etc. Use options(mc.cores =<foo>)
to set the (maximal) number of processes forked by the parallelization. If
options()$mc.cores
is unspecified, it is set to 2.
spikeAndSlab(
y,
X,
family = c("gaussian", "binomial", "poisson"),
hyperparameters = list(),
model = list(),
mcmc = list(),
start = list()
)
response
design matrix
(character) the family of the response, defaults to normal/Gaussian response
a list of hyperparameters controlling the priors (see details)
a list with information about the grouping structure of the model (see details)
(optional) list setting arguments for the sampler (see details)
(optional) list containing the starting values for
a list with components:
formula
see arguments
data
see arguments
family
see arguments
y
see arguments
X
see arguments
hyperparameters
see arguments
model
see arguments
mcmc
see arguments
start
see arguments
posteriorPred
a list with entries mu
and
y
containing samples of the expected values and realizations of the
response from the posterior predictive
postMeans
a list containing the posterior means of the parameters:
beta
the regression coefficients
alpha
ksi
tau
hypervariances of the penalized model terms
gamma
inclusion indicator variables of the model terms
pV1
w
hyperparameter
for gamma
sigma2
error variance (for Gaussian data)
logLik
log likelihood
logPost
log of (unnormalized) posterior
samples
a list containing the posterior samples of the parameters, see above for explanation of the entries
DIC
a vector with
fitted
a matrix with the posterior mean of the linear predictor in the first column and the posterior mean of the expected response in the second.
runTime
of the sampler, in seconds
Details for model specification:
hyperparameters
w
hyperparameters for the c(alphaW = 1, betaW = 1)
, i.e. a uniform distribution.
tau2
hyperparameters for the c(a1 = 5, a2 = 25)
gamma
sets c(v0 = 0.00025)
sigma2
hyperparameters for c(b1 = 1e-4, b2 = 1e-4)
. Only relevant for Gaussian
response.
varKsi
variance for prior of
ksiDF
defaults to 0 for a gaussian prior for
model
groupIndicators
a factor that maps the columns of X to the different model terms
H
a matrix containing the hierarchy of the penalized model terms
n
number of observations
q
length of
scale
scale/weights of
the response, defaults to rep(1, n)
, use this to specify number of
trials for binomial data
offset
defaults to rep(0,
n)
mcmc
nChains
how many parallel chains to run: defaults to 3
chainLength
how many samples should be generated per chain, defaults to 500
burnin
how many initial iterations should be discarded, defaults to 100
thin
save only every thin
-th
iteration, defaults to 5
verbose
verbose output and report progress? defaults to TRUE
returnSamples
defaults to TRUE
sampleY
generate samples of y and its conditional expectation from posterior predictive? defaults to FALSE
useRandomStart
use random draw or ridge estimate for beta as starting value? defaults to TRUE, i.e. random starting values.
blocksize
approx. blocksizes of the updates for
scalemode
how to do term-wise rescaling of
subvectors of
modeSwitching
probability to do P-IWLS
with the mode of the proposal set to the current value, which is useful if
the chain gets stuck. Defaults to
reduceRet
don't return data and samples for
start
gamma
starting values for mcmc$useRandomStart
is FALSE
, otherwise drawn from the
prior.
tau2
starting values for mcmc$useRandomStart
is FALSE
, otherwise
drawn from the prior.
sigma2
starting values for
mcmc$useRandomStart
is FALSE
, otherwise drawn from the prior.
w
starting value for mcmc$useRandomStart
is FALSE
, otherwise drawn from the
prior.
seed
Sets RNG seed for reproducible results. Parallel chains are seeded with this seed incremented by the number of the chain.
Scheipl, F. (2010) Normal-Mixture-of-Inverse-Gamma Priors for Bayesian Regularization and Model Selection in Structured Additive Regression Models. LMU Munich, Department of Statistics: Technical Reports, No.84 (https://epub.ub.uni-muenchen.de/11785/)