spStack (version 1.0.1)

stackedSampler: Sample from the stacked posterior distribution

Description

A helper function to sample from the stacked posterior distribution to obtain final posterior samples that can be used for subsequent analysis. This function applies on outputs of functions spLMstack() and spGLMstack().

Usage

stackedSampler(mod_out, n.samples)

Value

An object of class stacked_posterior, which is a list that includes the following tags -

beta

samples of the fixed effect from the stacked joint posterior.

z

samples of the spatial random effects from the stacked joint posterior.

In case of model output of class spLMstack, the list additionally contains sigmaSq which are the samples of the variance parameter from the stacked joint posterior of the spatial linear model. For model output of class spGLMstack, the list also contains xi which are the samples of the fine-scale variation term from the stacked joint posterior of the spatial generalized linear model.

Arguments

mod_out

an object of class spLMstack or spGLMstack.

n.samples

(optional) If missing, inherits the number of posterior samples from the original output. Otherwise, it specifies number of posterior samples to draw from the stacked posterior. If it exceeds the number of posterior draws used in the original function, then a warning is thrown and the samples are obtained by resampling. It is recommended, to run the original function with enough samples.

Author

Soumyakanti Pan ,
Sudipto Banerjee

Details

After obtaining the optimal stacking weights \(\hat{w}_1, \ldots, \hat{w}_G\), posterior inference of quantities of interest subsequently proceed from the stacked posterior, $$ \tilde{p}(\cdot \mid y) = \sum_{g = 1}^G \hat{w}_g p(\cdot \mid y, M_g), $$ where \(\mathcal{M} = \{M_1, \ldots, M_g\}\) is the collection of candidate models.

See Also

spLMstack(), spGLMstack()

Examples

data(simGaussian)
dat <- simGaussian[1:100, ]

mod1 <- spLMstack(y ~ x1, data = dat,
                  coords = as.matrix(dat[, c("s1", "s2")]),
                  cor.fn = "matern",
                  params.list = list(phi = c(1.5, 3),
                                     nu = c(0.5, 1),
                                     noise_sp_ratio = c(1)),
                  n.samples = 1000, loopd.method = "exact",
                  parallel = FALSE, solver = "ECOS", verbose = TRUE)
print(mod1$solver.status)
print(mod1$run.time)

post_samps <- stackedSampler(mod1)
post_beta <- post_samps$beta
print(t(apply(post_beta, 1, function(x) quantile(x, c(0.025, 0.5, 0.975)))))