Model{ for(i in 1:N){ y[i]~dnorm(mu[i],tau[1]) ##ppp calcs y.sim[i]~dnorm(mu[i],tau[1]) res[i]<-y[i]-mu[i] res.sim[i]<-y.sim[i]-mu[i] chi[i]<-res[i]*res[i]/(mu[i]*mu[i]) chi.sim[i]<-res.sim[i]*res.sim[i]/(mu[i]*mu[i]) ## mu[i]<- alpha + mu.cov[i] + frag.e[fragid[i]] #linear predictor mu.cov[i]<-mu.jump[i]-mu.jump[N+1] #remove jump intercept #priors for missing covariate values for(v in 1:Q){ Xc[i,v]<-cut(X[i,v]) X[i,v]~dnorm(0,1) } } alpha~dnorm(0,0.0001) # prior for intercept (=mean BC) mu.jump[1:(N+1)]<-jump.lin.pred(Xc[1:(N+1),1:Q],k,tau.beta) # rj mcmc model selection k~dbin(0.5,kmax) tau.beta<-1/pow(sd.beta,2) sd.beta~dunif(0,maxsd) #prior for random effects... #fragment for(f in 1:Nfrag){ frag.e[f]~dnorm(0,tau[2]) } #prior for residual precision (=1/variance) tau[1]~dgamma(0.001,0.001) #priors for random effects SDs for(t in 2:3){ tau[t]<-1/pow(sd[t],2) sd[t]~dunif(0,maxsd2) } maxsd2<-maxsd*2 sd[1]<-1/sqrt(tau[1]) #resid SD ##extract coefficients and probabilities of inclusion from rj mcmc pred[1:(Q+1)] <- jump.lin.pred.pred(mu.jump[1:(N+1)], X.pred[1:(Q+1), 1:Q]) for (i in 1:Q) { X.pred[i, i] <- 1 for (j in 1:(i - 1)) {X.pred[i, j] <- 0} for (j in (i + 1):Q) {X.pred[i, j] <- 0} X.pred[(Q+1), i] <- 0 beta[i]<-pred[i]-pred[Q+1] inccol[i]<-1-equals(beta[i],0) } for(j in 1:Nvars){ inc[j]<-1-equals(sum(inccol[starts[j]:ends[j]]),0) } for(v in 1:Q){ Xc[N+1,v]<-X[N+1,v] } ### ppp calcs chi.sum[1]<-sum(chi[]) chi.sum[2]<-sum(chi.sim[]) ppp<-step(chi.sum[1]-chi.sum[2]) }