[R] BRugs

N S Ha neungsoo.ha at gmail.com
Mon Apr 19 18:04:54 CEST 2010


Thanks for the reply Bob, but it still does not work, you see. I ran this
model, just with the main effects and it ran fine.

n=length(bi.bmi)

Lgen=2
Lrace=5
Lagegp=13
Lstra=15
Lpsu=2

bi.bmi.model=function(){
# likelihood
for(i in 1:n){
	bi.bmi[i]~ dbern(p[i])
	logit(p[i])<- a0 + a1[agegp[i]]+a2[gen[i]]+a3[race[i]]
			    + g[stra[i]]+ u[psu[i],stra[i]]
		}
# constraints for a1, a2, a3
a1[1]<-0.0
a2[1]<-0.0
a3[1]<-0.0
# priors
a0~ dnorm(0.0, 1.0E-4)
for(j in 2:Lagegp){a1[j]~ dnorm(0.0, 1.0E-4)}
for(j in 2:Lgen){ a2[j]~ dnorm(0.0, 1.0E-4)}
for(k in 2:Lrace){ a3[k]~ dnorm(0.0, 1.0E-4)}

for(l in 1:Lstra){
	g[l]~dunif(0, 100)
	}
for( m in 1:Lpsu){
	for(l in 1:Lstra){
		u[m,l]~ dnorm(0.0, tau.u)
		}}
tau.u<-pow(sigma.u, -2)
sigma.u~ dunif(0.0,100)
}
 
library(BRugs)
writeModel(bi.bmi.model, con='bi.bmi.model.txt')
model.data=list( 'n','Lagegp', 'Lgen', 'Lrace', 'Lstra', 'Lpsu',				 
				'bi.bmi','agegp', 'gen', 'race','stra', 'psu')
model.init=function(){
	list( sigma.u=runif(1),
		a0=rnorm(1), a1=c(NA, rep(0,12)),
		a2=c(NA, rep(0, 1)),
		a3=c(NA, rep(0, 4)), 
 		g=rep(0,Lstra), u=matrix(rep(0, 30), nrow=2)	)
	}
model.parameters=c( 'a0', 'a1', 'a2', 'a3')
model.bugs=BRugsFit(modelFile='bi.bmi.model.txt',
				   data=model.data,
				   inits=model.init,
				   numChains=1, 
				   para=model.parameters,
				   nBurnin=50, nIter=100)

This is just with the main effects, and this does not give me any problems,
and I also ran the following model with interaction term between gen and
race, and it also ran fine.
for (i in 1:n){
	bi.bmi[i]~ dbern(p[i])
	logit(p[i])<- a0 + a1[agegp[i]]+a2[gen[i]]+a3[race[i]]
			   + a23[gen[i], race[i]]
			  + gam[stra[i]]+ u[psu[i],stra[i]]
		}
# constraints for a2, a3, a12 and a13
a1[1]<-0.0
a2[1]<-0.0
a3[1]<-0.0
a23[1,1]<-0.0

#gen x race
for(j in 2:Lrace){ a23[1,j]<-0.0}
for(k in 2:Lgen){ a23[k,1]<-0.0}
# priors
a0~ dnorm(0.0, 1.0E-4)
for(i in 2:Lagegp){a1[i]~dnorm(0.0, 1.0E-4)}
for(i in 2:Lgen){ a2[i]~ dnorm(0.0, 1.0E-4)}
for(i in 2:Lrace){ a3[i]~ dnorm(0.0, 1.0E-4)}
for(i in 2:Lgen){
	for(j in 2:Lrace){
		a23[i,j]~ dnorm(0.0, 1.0E-4)
	}}
for(i in 1:Lstra){
	gam[i]~dunif(0, 1000)
	}
for( i in 1:Lpsu){
	for(j in 1:Lstra){
		u[i,j]~ dnorm(0.0, tau.u)
		}}
tau.u<-pow(sigma.u, -2)
sigma.u~ dunif(0.0,100)
}

So, the error happens only when I try to plug in interaction with the agegp.
I still don't know how to correct it.
Thanks

-- 
View this message in context: http://n4.nabble.com/BRugs-tp2015395p2016164.html
Sent from the R help mailing list archive at Nabble.com.



More information about the R-help mailing list