分層迴歸UScrime R/JAGS

我組建了一個模型中JAGS爲UScrime數據的層次迴歸（從庫（MASS）包運行代碼。犯罪率作爲響應，用15個預測因子。分層迴歸UScrime R/JAGS

我在我的代碼中的一些錯誤，我無法彌補

JAGS型號：

model{ 
    for (i in 1:m){ 
     for (j in 1:47){ 
     y[j]~dnorm(beta[i,1]+beta[i,2]*x[j]+beta[i,3]*x[j]+beta[i,4]*x[j] 
     +beta[i,5]*x[j]+beta[i,6]*x[j]+beta[i,7]*x[j]+beta[i,8]*x[j] 
     +beta[i,9]*x[j]+beta[i,10]*x[j]+beta[i,11]*x[j]+beta[i,12]*x[j] 
     +beta[i,13]*x[j]+beta[i,14]*x[j]+beta[i,15]*x[j],invsig2)  
     } 
     beta[i,1:15]~dmnorm(theta[],invSig[,]) 
     } 
    theta[1:15]~dmnorm(mu0[],Lam.inv[,]) 
    invSig[1:15,1:15]~dwish(Lam[,],30) 
    invsig2~dgamma(.5,.5*sig20) 
    sig2<-1/invsig2 
    Sig<-inverse(invSig) 
}

而在R上的代碼：

crime=read.table("crime.dat",header=T) 
crime=as.matrix(crime) 

y=crime[,1];M=crime[,2];So=crime[,3];Ed=crime[,4];Po1=crime[,5];Po2=crime[,6] 
LF=crime[,7];M.F=crime[,8];Pop=crime[,9];NW=crime[,10];U1=crime[,11] 
U2=crime[,12];GDP=crime[,13];Ineq=crime[,14];Prob=crime[,15] 
Time=crime[,16] 
Crime=crime[,1] 
m <- length(unique(Crime)) # m = length of crimes 
n <- rep(NA,m) # setting up placeholder vectors 

for (j in 1:m){ 
n[j] <- sum(Crime==j) 
} 

beta_lse=matrix(nrow=m,ncol=15) 
sig2_lse=rep(NA,m) 

for (i in 1:m){ 
crime=y[1:47] 
x1=M[1:47];x2=So[1:47];x3=Ed[1:47];x4=Po1[1:47];x5=Po2[1:47];x6=LF[1:47] 
x7=M.F[1:47];x8=Pop[1:47];x9=NW[1:47];x10=U1[1:47];x11=U2[1:47] 
x12=GDP[1:47];x13=Ineq[1:47];x14=Prob[1:47];x15=Time[1:47] 
lse=lm(crime ~ 0+x1+x2+x3+x4+x5+x6+x7+x8+x9+x10+x11+x12+x13+x14+x15) 
temp=anova(lse) 
beta_lse[i,]=c(lse[[1]]) 
sig2_lse[i]=temp[[3]][2]} 

Lam=cov(beta_lse) 
Lam.inv=solve(Lam) 
mu0=colMeans(beta_lse) 
sig20=mean(sig2_lse) 

library(rjags) 

forJags<-list(y=y,m=m,x=x1+x2+x3+x4+x5+x6+x7+x8+x9+x10+x11+x12+x13+x14+x15, 
    Lam=Lam, Lam.inv=Lam.inv, mu0=mu0, sig20=sig20) 

inits<-list(beta=matrix(0,nrow=m,ncol=15),theta=c(0,0), invSig=.1*diag(2), invsig2=1) 

mod<-jags.model(file="CrimeHierarchicalRegression.bug",data=forJags,inits=inits) 
out<-coda.samples(model=mod,variable.names=c("theta","Sig","sig2","beta"),n.iter=10000,thin=5) 
summary(out)

有兩個加粗的錯誤，Error1和Error2，它們各自的行導致錯誤。第一個錯誤告訴我「要替換的項目數量不是替換長度的倍數」。第二個錯誤告訴我，「有一個意外的符號」在該行的某處。

請讓我知道，如果你看到這些錯誤是如何造成的。提前致謝。

來源

2016-04-26 PapaSwankey

我想補充一點，這不是一項任務，而是在準備期末考試時給出的練習題。 – PapaSwankey

對於您的第一個錯誤，您試圖將長度爲15的行的矢量放入16列的矩陣中。對於第二個錯誤，在x12=x12和x13=x13之間沒有逗號。

添加該逗號並將beta_lse=matrix(nrow=m,ncol=16)更改爲beta_lse=matrix(nrow=m,ncol=15)並且您應該很好。這當然假設您不想爲模型擬合截距，這是您目前在模型的線性預測器中指定的截距。這應該解決這兩個具體的錯誤。

來源

2016-04-27 14:21:18

我不再收到這些錯誤。但是，在我使用Lam.inv = solve（Lam）的線上。我現在得到的錯誤是：「Lapack routine dgesv：system is exactly singular：U [1,1] = 0」。我查看了代碼，這基本上意味着我的Lam矩陣是不可逆的。任何想法如何解決？ – PapaSwankey

這意味着您的矩陣是不可逆的，並且已經在堆棧溢出的其他部分進行了討論。 http://stackoverflow.com/questions/6572119/r-solvesystem-is-exactly-singular –

分層迴歸UScrime R/JAGS

回答

相關問題