Bayesian Hierarchical Model - Estimating Rho (correlation between intercepts and slopes)

Dear List

I'm trying to run a hierarchical model to look at the effect of body condition and environmental variables on breeding probability (1/0). The model I'm using is adapted from Gelman and Hill (2006) pg 376, with the modifications to add a group level predictor, for both intercepts and slopes, listed on pg 379. I'm modelling body condition at the individual level, and NAO at the group (year level). These models are being run in JAGS using Denwood's 'runjags' interface for R.

The issue I'm having is that the 95% credible intervals for rho, the correlation among intercepts and slopes, is -.95 to 0.93, so essentially is spanning the entire available range. I would expect that if the model were 'behaving', then this estimate would take some kind of sensible value i.e. perhaps with  a posterior mode of zero, but much more certainty that it's zero. This may be telling me simply that there is no relationship between intercepts and slopes and that the correlated model is unjustified, but I would appreciate any advice on the model (listed below) which I quite likely have coded incorrectly. If anyone would like simulated data to go with this code I can happily provide it.

As an aside, what are the likely effects of fiddling with the the degrees of freedom (set to 3 here)? Is this akin to the 'nu' parameter in MCMCglmm?

Any help is greatly appreciated

Regards

Xavier

#########################################
#Separate Intercepts and Slopes - CORRELATED Scaled inverse-Wishart Distribution AND Group-Level Predictor
#########################################

cfgroupmod<-"model{

#Individual model

for (i in 1:n.samp){
	logit.p[i]<-a[yearcode[i]] + b[yearcode[i]]*mass[i] 
	p[i]<-exp(logit.p[i])/(1+exp(logit.p[i]))
	breed[i]~dbern(p[i])
	}
	
	
#Group Model 

for (j in 1:5){
	a[j]<-xi.a*B.raw[j,1]
	b[j]<-xi.b*B.raw[j,2]
	B.raw[j,1:2]~dmnorm(B.raw.hat[j,],Tau.B.raw[,])
	B.raw.hat[j,1]<-g.a.0.raw+g.a.1.raw*naoj[j]
	B.raw.hat[j,2]<-g.b.0.raw+g.b.1.raw*naoj[j]
		}
		
	#Recover Parameters on Original Scale
	g.a.0<-xi.a*g.a.0.raw
	g.a.1<-xi.a*g.a.1.raw
	g.b.0<-xi.b*g.b.0.raw
	g.b.1<-xi.b*g.b.1.raw
	
	#Priors for scaled Parameters
	g.a.0.raw~dnorm(0,0.0001)
	g.a.1.raw~dnorm(0,0.0001)
	g.b.0.raw~dnorm(0,0.0001)
	g.b.1.raw~dnorm(0,0.0001)
	
	xi.a~dunif(0,100)	
	xi.b~dunif(0,100)
	
	Tau.B.raw[1:2,1:2]~dwish(W[,],df)
	df<-3
	Sigma.B.raw[1:2,1:2]<-inverse(Tau.B.raw[,])
	sigma.a<-xi.a*sqrt(Sigma.B.raw[1,1])
	sigma.b<-xi.b*sqrt(Sigma.B.raw[2,2])
	rho<-Sigma.B.raw[1,2]/sqrt(Sigma.B.raw[1,1]*Sigma.B.raw[2,2])

}"
	
W<-diag(2)

cfgroupdat<-dump.format(list(n.samp=nrow(fem),n.year=5,breed=fem$futurebreed,mass=fem$seasonbci,yearcode=fem$yearcode,naoj=naoframe$junnaoresids,W=W))

cfgroupparams<-c("a","b","g.a.0","g.a.1","g.b.0","g.b.1","sigma.a","sigma.b","rho")
burncfgroup=20000
sampcfgroup=200000
adaptcfgroup=10000
cfgroupthin=200



---------------------------------------------------------
Dr Xavier Harrison
BBSRC Postdoctoral Research Associate
Centre for Ecology and Conservation
University of Exeter
Tremough Campus
Penryn
Cornwall
TR10 9EZ
Tel: (01326) 371872
xav.harrison at gmail.com










	[[alternative HTML version deleted]]