I received several requests to summarize, so I am posting the replies,
although I have not yet gotten this to work myself (undefined real result
trap while compiling, I may move the computation of the C and M to S-Plus as
suggested in the Geobugs manual). Thanks to Devin Johnson, Erik Sauleau, and
Gentry White for helpful replies.
My original question:
I am trying to predict occurrence of a species (survey success, within a
regular grid overlaid on past survey locations) using a CAR model with
environmental covariates. My initial model was adapted from Gelfand et al.
2005. (Modelling species diversity through species level hierarchical
modelling. Royal Stat Soc C Appl. Statist. 54, Part 1, pp. 1-20) but they
use a car.normal and I would like to rewrite it as a car.proper. I have
looked over the Geobugs examples for car.proper (eg, lip cancer) but am new
at Winbugs and unsure how to adapt them to this problem. Could anyone point
me to a car.proper example similar to this problem?
More generally, what are the pros and cons of using car.normal vs.
car.proper for this type of question?
the car.normal model is:
# likelihood
for (i in 1 : N_nonzeroy) {
y[ind[i]] ~ dbin(p[ind[i]], n[ind[i]])
}
for(i in 1:N_LOC){
logit(p[i]) <- rho[i]+xbeta[i]+mu
xbeta[i]<-beta[1]*tri[i] + beta[2]*con[i]
}
# CAR prior distribution for spatial random effects:
rho[1:N_LOC] ~ car.normal(adj[], weights[], num[], tau)
for(k in 1:sumNumNeigh) {
weights[k] <- 1
}
# other priors
mu ~ dflat()
for ( i in 1:2) { beta[i] ~ dnorm(0, 1) }
vrho ~ dnorm(0, 0.2) I(0,)
tau <- 1/vrho
}
Replies:
Devin Johnson:
Here we used the row standardized parameterization.
This means the rows of the C matrix sum to 1 and the M vector has entries
1/(num neighbors). I always like the proper CAR because it has a joint
distribution for the entire random effects vector, whereas the CAR normal
sets the spatial parameter from the CAR proper to 1 which means there is no
joint distribution for the entire vector. You probably won't see much of a
difference between the two, but I could be wrong.
model
{
for(i in 1:Nplots){
wolverine[i] ~ dbin(p[i],1)
logit(p[i]) <- Int + sig*eps[i]
mu.eps[i] <- 0
}# End site loop
eps[1:Nplots] ~ car.proper(mu.eps[], C[], adj[],num[], M[], 1, gamma)
Int ~ dflat()
sig ~ dunif(0,10)
g.up <- max.bound(C[], adj[], num[], M[])
g.lo <- 0 #min.bound(C[], adj[], num[], M[])
gamma ~ dunif(g.lo, g.up)
}# End model loop
Erik Sauleau:
I would use in the model part:
for(i in 1:N_LOC) {m[i] <- 1/n[ind[i]] }
cumsum[1] <- 0
for(i in 2:(N_LOC+1)) {cumsum[i] <- sum(num[1:(i-1)]) }
for(k in 1:sumNumNeigh) {
for(i in 1:N_LOC) {
pick[k,i]<-step(k-cumsum[i]-epsilon)*step(cumsum[i+1]-k) }
C[k]<-1/ inprod(num[], pick[k,]) }
epsilon <- 0.0001
for (i in 1:N_LOC){
[INSERT HERE YOUR MODEL]
theta[i] <- 0
}
rho[1:N_LOC] ~ car.proper(theta[], C[], adj[], num[], m[], prec, gamma)
prec ~ dgamma(0.5, 0.0005)
gamma.min <- min.bound(C[], adj[], num[], m[])
gamma.max <- max.bound(C[], adj[], num[], m[])
gamma ~ dunif(gamma.min, gamma.max)
Gentry White:
I have adapted the model in He and Sun 2000 for use in WinBUGS using the
car.proper prior like so:
Z[1:N]~car.proper(mu.z[],C[],adj[],num[],M[],tau.z,rho)
rho~dunif(-0.3457,0.1756)
tau.z ~ dgamma(2.15,0.09) # prior on precision
where M and C were vectors of 1's of appropriate lengths.
I have also in the past (before car.proper was available used this model:
for(h in 1:N)
{
for( k in 1:8)
{
zsum[h,k]<-Z[adj[h,k]]
}
muz[h]<-rho*sum(zsum[h,])
Z[h]~dnorm(muz[h], tau.z)
}
Z[115] <- 0
rho~dunif(-0.3457,0.1756)
tau.z ~ dgamma(2.15,0.09) # prior on precision
Which should be the same and again comes from He and Sun 2000.
The use of the car.proper versus car.normal question I will try to answer to
the best of my ability, I think that the car.proper prior is that the
car.normal prior says that the conditional variance depends on the number of
counties and doesn't give any measure of correlation among regions.
Additionally under the car.normal prior the joint density doesn't exist. I
am paraphrasing from the He and Sun paper, though there are a number of
other sources on this.
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
|