sppm {ppmSuite} | R Documentation |
Function that fits spatial product partition model with Gaussian likelihood
Description
sppm
is the main function used to fit model with Guassian likelihood and spatial PPM as prior on partitions.
Usage
sppm(y,s,
s.pred=NULL,
cohesion,
M=1,
modelPriors=c(0, 100^2, 10, 10),
cParms=c(1, 1.5, 0, 1, 2, 2),
mh=c(0.5, 0.5),
draws=1100,burn=100,thin=1)
Arguments
y |
numeric vector containing response variable |
s |
Two-column matrix containing spatial locations (i.e., longitude and lattitude). |
s.pred |
Two-column matrix containing spatial locations at which out-of-sample predictions will be collected. |
cohesion |
Scalar that indicates which cohesion to use.
|
M |
Parameter related to Dirichlet process scale or dispersion parameter. |
modelPriors |
Vector containing model prior values (see below for more details) |
cParms |
Vector containing partition model prior values (see below for more details) |
mh |
Tuning standard deviations for metropolis updates for sigma2 and sigma20 |
draws |
Number of MCMC samples to collect |
burn |
Number of the MCMC samples discarded in the burn-in phase of the sampler |
thin |
The amount of thinning desired for the chain |
Details
The vector modelPriors = c(m0, s20, ms, ms0)
where each prior parameter is listed in the model description below.
The cParm vector contains values associated with the cohesion function.
cParm = c(
epsilon value - cohesion 1 only,
distance bound - cohesion 2 only,
mu0 - center of NNIG for cohesion 3 and 4
k0 - scale parm of gaussian part of NNIG for cohesion 3 and 4
v0 - degrees of freedom IG part of NNIG for cohesion 3 and 4
L0 - scale parm (scalar of identity matrix) IG part of NNIG for cohesion 3 and 4).
The model this function fits is Gaussian likelihood model using the sPPM prior on partitions (Page and Quintana, 2016). Specific model details are
y_i | \mu^*, \sigma^{2*}, c_i \sim N(\mu_{c_i}^*, \sigma^{2*}_{c_i}), i=1,\ldots,n
\mu^*_j | \mu_0, \sigma^2_0 \sim N(\mu_0,\sigma^2_0)
\sigma^*_j | A \sim UN(0, ms)
\rho|M,\xi \sim sPPM
To complete the model specification, the folloing hyperpriors are assumed,
\mu_0 | m, s^2 \sim N(m0,s0^2)
\sigma_0 | B \sim UN(0,ms0)
Note that we employ uniform prior distributions on variance components as suggest in Gelman's 2006 Bayesian paper. "sPPM" in the model specificaiton denotes the the spatial product partition model. The computational implementation of the model is based algorithm 8 found in Neal's 2000 JCGS paper.
Value
This function returns in a list all MCMC interates for each model parameter, posterior predictive, and fitted values. In addition the LPML model fit metric is provided.
Examples
data(scallops)
Y<-log(scallops[,5]+1)
s_coords <- scallops[,3:4] #lat and long
m <- dim(s_coords)[1]
# standardize spatial coordinates
smn <- apply(s_coords,2,mean)
ssd <- apply(s_coords,2,sd)
s_std <- t((t(s_coords) - smn)/ssd)
# Create a grid of prediction locations
np <- 10
sp <- expand.grid(seq(min(s_coords[,1]), max(s_coords[,1]),length=np),
seq(min(s_coords[,2]), max(s_coords[,2]), length=np))
sp_std <- t((t(sp) - smn)/ssd) # standardized prediction spatial coordinates
niter <- 20000
nburn <- 10000
nthin <- 10
nout <- (niter - nburn)/nthin
out <- sppm(y=Y,s=s_std,s.pred=sp_std,cohesion=4, M=1, draws=niter, burn=nburn, thin=nthin)
# fitted values
fitted.values <- out$fitted
fv.mn <- apply(fitted.values, 2,mean)
mean((Y - fv.mn)^2) # MSE
out$lpml #lpml value
ppred <- out$ppred
predmn <- apply(ppred,2,mean)
# The first partition iterate is used for plotting
# purposes only. We recommended using the salso
# R-package to estimate the partition based on Si
Si <- out$Si
plot(s_coords[,1], s_coords[,2], col=Si[1,])