要崩溃了,我按照Pdf书上的代码打下来,一运行就错!我已经反复核查无数遍了,和书一摸一样了,可还是出错[s:15]!
求大神帮助,最好能给个代码关于SSVS的,如果没有那介绍哪个好书有代码也可以阿!(附上代码求鉴定)
model {
for(i in 1:n){
X[i,1]<-1.0
for(j in 1:p){ X[i,j+1]<-x[i,j]
}}
B[1]<-beta0
g[1]<-gamma0
for(j in 1:p){
B[j+1]<-beta[j]
g[j+1]<-gamma[j]
}
for (j in 1:(p+1)){gb[j]<-g[j]*B[j] }
for(i in 1:n){
y~dnorm(mu,tau)
mu<-inprod(X[i,],gb[])
}
gamma0~dbern(0.5)
for ( j in 1:p ){ gamma[j]~dbern(0.5) }
beta0~dnorm( mb0,taub0)
mb0<- (1-gamma0)*prop.mean.beta0
taub0<- gamma0*0.01+(1-gamma0)/pow(prop.sd.beta0,2)
for (j in 1:p){
beta[j]~dnorm( mb[j],taub[j])
mb[j]<-prop.mean.beta[j]
taub[j]<-(gamma[j]/n+(1-gamma[j]))/pow(prop.sd.beta[j],2)
}
tau~dgamma(0.01,0.01)
for (j in 1:(p+1)){ g[j]~dbern(0.5)}
B[1:(p+1)]~dmnorm(mean.beta[1:(p+1)], T[1:(p+1),1:(p+1)])
tau~dgamma(0.01,0.01)
for(j in 1:(p+1) ){ mean.beta[j]<- (1-g[j])*prop.mean.beta[j]}
for(j in 1:(p+1) ){ for (k in 1:(p+1) ){
T[j,k]<-g[j]*g[k]*tau*XTX[j,k]/n
+(1-g[j]*g[k])*equals(j,k)*pow(prop.sd.beta[k],-2)
}}
gamma0<-g[1]
beta0<-B[1]
for (j in 1:p){
beta[j]<-B[j+1]
gamma[j]<-g[j+1]
}
for (j in 1:p){ mindex[j]<-pow(2,j)}
model<-gamma0+inprod(gamma[],mindex[])
pmodel[1]<-equals( model,pow(2,4)+pow(2,5))
pmodel[2]<-equals( model,pow(2,4)+pow(2,5)+pow(2,12))}