Description Usage Arguments Examples
this function runs Markov Chain Monte Carlo to estimate the Bayesian multi-resolution spatial regression model.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 |
y |
is a n vector of Gaussian data. |
X |
is a n x p matrix of fixed effects (like latitude, elevation, etc) |
locs |
is a n x 2 matrix of observation locations. |
params |
is a list of parameter settings. The list
|
priors |
is the list of prior settings. |
M |
The number of resolutions. |
n_neighbors |
The expected number of neighbors for each interior basis function. This determines the basis radius parameter. |
n_coarse_grid |
The number of basis functions in one direction (e.g. |
n_padding |
The number of additional boundary points to add on each boundary. For example, n_padding = 5 will add 5 boundary knots to the both the left and right side of the grid). |
n_cores |
is the number of cores for parallel computation using openMP. |
inits |
is the list of initial values if the user wishes to specify initial values. If these values are not specified, then the initial values will be randomly sampled from the prior. |
config |
is the list of configuration values if the user wishes to specify initial values. If these values are not specified, then default a configuration will be used. |
verbose |
Should verbose output be printed? Typically this is only useful for troubleshooting. |
use_spam |
is a boolean flag to determine whether the output is a list of spam matrix objects ( |
n_chain |
is the MCMC chain id. The default is 1. |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 | set.seed(111)
## genereate the locations
locs <- matrix(runif(20), 10, 2)
## generate some covariates and regression coefficients
X <- cbind(1, matrix(rnorm(30), 10, 3))
beta <- rnorm(ncol(X))
## simulate the MRA process
M <- 2
MRA <- mra_wendland_2d(locs, M = M, n_coarse_grid = 4)
W <- do.call(cbind, MRA$W)
n_dims <- rep(NA, length(MRA$W))
dims_idx <- NULL
for (i in 1:M) {
n_dims[i] <- ncol(MRA$W[[i]])
dims_idx <- c(dims_idx, rep(i, n_dims[i]))
}
## set up the process precision matrices
Q_alpha <- make_Q_alpha_2d(sqrt(n_dims), c(0.9, 0.8))
Q_alpha_tau2 <- make_Q_alpha_tau2(Q_alpha, tau2 = c(2, 4))
## add in constraints so each resolution has random effects that sum to 0
A_constraint <- sapply(1:M, function(i){
tmp = rep(0, sum(n_dims))
tmp[dims_idx == i] <- 1
return(tmp)
})
a_constraint <- rep(0, M)
alpha <- as.vector(spam::rmvnorm.prec.const(
n = 1,
mu = rep(0, nrow(W)),
Q = Q_alpha_tau2,
A = t(A_constraint),
a = a_constraint))
## define the data
y <- as.vector(X %*% beta + W %*% alpha + rnorm(10))
## define the params for MCMC fitting
params <- list(
n_mcmc = 5,
n_adapt = 5,
n_thin = 1,
n_message = 5)
## define the model priors
priors <- list(
alpha_tau2 = 1,
beta_tau2 = 1,
alpha_sigma2 = 1,
beta_sigma2 = 1,
mu_beta = rep(0, ncol(X)),
Sigma_beta = 5 * diag(ncol(X)))
## Fit the MRA model using MCMC
out <- mcmc_mra(
y = y,
X = X,
locs = locs,
params = params,
priors = priors,
M = 2,
n_coarse_grid = 4,
n_cores = 1L,
verbose = FALSE
)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.