georamps | R Documentation |
General function for fitting Bayesian geostatistical models using the reparameterized and marginalized posterior sampling (RAMPS) algorithm of Yan et al. (2007).
georamps(fixed, random, correlation, data, subset, weights,
variance = list(fixed = ~ 1, random = ~ 1, spatial = ~ 1),
aggregate = list(grid = NULL, blockid = ""), kmat = NULL,
control = ramps.control(...), contrasts = NULL, ...)
fixed |
two-sided linear |
random |
optional one-sided formula of the form |
correlation |
|
data |
optional data frame containing the variables named in |
subset |
optional expression indicating the subset of rows in |
weights |
optional numerical vector of measurement error variance (inverse) weights to be used in the fitting process. Defaults to a value of 1 for point-source measurements and the number of grid points for areal measurements (see the |
variance |
optional list of one-sided formulas, each of the form |
aggregate |
optional list of elements: |
kmat |
optional |
control |
list of parameters for controlling the fitting process. See the |
contrasts |
optional list. See the |
... |
further arguments passed to or from other methods. |
An object of class 'ramps'
containing the following elements:
params |
|
z |
|
loglik |
vector of data log-likelihood values at each MCMC iteration. |
evals |
vector of slice sampler evaluations at each MCMC iteration. |
call |
the matched function call to |
y |
response vector. |
xmat |
design matrix for the main effects. |
terms |
the |
xlevels |
list of the factor levels for |
etype |
grouping factor for the measurement error variances. |
weights |
weights used in the fitting process. |
kmat |
matrix for mapping the spatial parameters to the observed data. |
correlation |
specified |
coords |
matrix of unique coordinates for the measurement and grid sites. |
ztype |
grouping factor for the spatial variances. |
wmat |
matrix for mapping the random effects to the observed data. |
retype |
grouping factor for the random effects variances. |
control |
a list of control parameters used in the fitting process. |
Brian Smith brian-j-smith@uiowa.edu, Jun Yan jun.yan@uconn.edu, and Kate Cowles kate-cowles@uiowa.edu
Yan, J., Cowles, M.K., Wang, S., and Armstrong, M. (2007) “Parallelizing MCMC for Bayesian Spatiotemporal Geostatistical Models”, Statistics and Computing, 17(4), 323-335.
Smith, B. J., Yan, J., and Cowles, M. K. (2008) “Unified Geostatistical Modeling for Data Fusion and Spatial Heteroskedasticity with R Package ramps”, Journal of Statistical Software, 25(10), 1-21.
corRClasses
,
ramps.control
,
mcmc
,
DIC.ramps
,
plot.ramps
,
predict.ramps
,
summary.ramps
,
window.ramps
## Not run:
## Load the included uranium datasets for use in this example
data(NURE)
## Geostatistical analysis of areal measurements
NURE.ctrl1 <- ramps.control(
iter = 25,
beta = param(0, "flat"),
sigma2.e = param(1, "invgamma", shape = 2.0, scale = 0.1, tuning = 0.75),
phi = param(10, "uniform", min = 0, max = 35, tuning = 0.50),
sigma2.z = param(1, "invgamma", shape = 2.0, scale = 0.1)
)
NURE.fit1 <- georamps(log(ppm) ~ 1,
correlation = corRExp(form = ~ lon + lat, metric = "haversine"),
weights = area,
data = NURE,
subset = (measurement == 1),
aggregate = list(grid = NURE.grid, blockid = "id"),
control = NURE.ctrl1
)
print(NURE.fit1)
summary(NURE.fit1)
## Analysis of point-source measurements
NURE.ctrl2 <- ramps.control(
iter = 25,
beta = param(0, "flat"),
sigma2.e = param(1, "invgamma", shape = 2.0, scale = 0.1, tuning = 0.75),
phi = param(10, "uniform", min = 0, max = 35, tuning = 0.5),
sigma2.z = param(1, "invgamma", shape = 2.0, scale = 0.1)
)
NURE.fit2 <- georamps(log(ppm) ~ 1,
correlation = corRExp(form = ~ lon + lat, metric = "haversine"),
data = NURE,
subset = (measurement == 2),
control = NURE.ctrl2
)
print(NURE.fit2)
summary(NURE.fit2)
## Joint analysis of areal and point-source measurements with
## prediction only at grid sites
NURE.ctrl <- ramps.control(
iter = 25,
beta = param(rep(0, 2), "flat"),
sigma2.e = param(rep(1, 2), "invgamma", shape = 2.0, scale = 0.1, tuning = 0.75),
phi = param(10, "uniform", min = 0, max = 35, tuning = 0.5),
sigma2.z = param(1, "invgamma", shape = 2.0, scale = 0.1),
z.monitor = NURE.grid
)
NURE.fit <- georamps(log(ppm) ~ factor(measurement) - 1,
correlation = corRExp(form = ~ lon + lat, metric = "haversine"),
variance = list(fixed = ~ measurement),
weights = area * (measurement == 1) + (measurement == 2),
data = NURE,
aggregate = list(grid = NURE.grid, blockid = "id"),
control = NURE.ctrl
)
print(NURE.fit)
summary(NURE.fit)
## Discard initial 5 MCMC samples as a burn-in sequence
fit <- window(NURE.fit, iter = 6:25)
print(fit)
summary(fit)
## Deviance Information Criterion
DIC(fit)
## Prediction at unmeasured sites
ct <- map("state", "connecticut", plot = FALSE)
lon <- seq(min(ct$x, na.rm = TRUE), max(ct$x, na.rm = TRUE), length = 20)
lat <- seq(min(ct$y, na.rm = TRUE), max(ct$y, na.rm = TRUE), length = 15)
grid <- expand.grid(lon, lat)
newsites <- data.frame(lon = grid[,1], lat = grid[,2],
measurement = 1)
pred <- predict(fit, newsites)
plot(pred, func = function(x) exp(mean(x)),
database = "state", regions = "connecticut",
resolution = c(200, 150), bw = 5,
main = "Posterior Mean",
legend.args = list(text = "ppm", side = 3, line = 1))
plot(pred, func = function(x) exp(sd(x)),
database = "state", regions = "connecticut",
resolution = c(200, 150), bw = 5,
main = "Posterior Standard Deviation",
legend.args = list(text = "ppm", side = 3, line = 1))
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.