OpenBUGS и JAGS не могут читать данные из столбца фрейма данных? Вам обязательно использовать вектор? - PullRequest
Новая
       64

OpenBUGS и JAGS не могут читать данные из столбца фрейма данных? Вам обязательно использовать вектор?

0 голосов
/ 08 мая 2020

Мне было любопытно, почему мой код openBUGS не работает. Я думаю, это из-за того, как я читаю данные. У меня есть фрейм данных, и я пишу data = list (death = com.s $ death, hday = com.s $ hday, age = com.s $ age). Разве нельзя помещать в столбец фрейма данных?

JAGS выдает ошибку «Ошибка в jags.model (file =" cctry.txt ", data = data, inits = Inits, n.chain = 1 ): Ошибка смерти узла [2] Узел несовместим с родительскими ". Хотя openBUGS говорит, что «переменная hday не определена в модели или в наборе данных. Модель должна быть инициализирована перед использованием мониторов, модель должна быть инициализирована до того, как DI C можно будет отслеживать, модель должна быть инициализирована перед обновлением».

Вот мой код. 

library(coda) 
library(R2OpenBUGS)
library(dummies);library(robustHD);require(Hmisc); library(MASS)

#open data 
getHdata(support2);View(support2)


my.s=data.frame(support2$hospdead,support2$death,support2$age,support2$hrt,
support2$hday,support2$temp, support2$pafi, support2$crea,support2$ca, 
support2$wblc, support2$resp,support2$sod,support2$scoma, support2$meanbp, 
support2$alb,support2$bili, support2$dzclass,support2$slos,support2$sex)
names(my.s)=c("hospdead","death","age","hrt","hday","temp","pafi","crea","ca",
"wblc","resp","sod","scoma","meanbp","alb","bili","dzclass","slos","sex")
my.s$hospdead=factor(my.s$hospdead);my.s$death=factor(my.s$death);my.s$sex=factor(my.s$sex);
my.s$dzclass=factor(my.s$dzclass);my.s$ca=factor(my.s$ca)


which(is.na(my.s$meanbp));which(is.na(my.s$scoma));which(is.na(my.s$sod))
which(is.na(my.s$resp));which(is.na(my.s$temp));which(is.na(my.s$hrt))
my.s=my.s[-c(5441,5394),]

#Make dummy variables
my.s=cbind(my.s,dummy(my.s$ca),dummy(my.s$dzclass),dummy(my.s$sex))
colnames(my.s)[20:28]=c("No","Yes","Meta","ARF","COPD","Coma","Cancer","F","M")


#Complete cases
com.s=my.s[complete.cases(my.s), ]

#Logistic regression 
#Standardizing coefficients
for (i in c(3,4,5,6,7,8,10,11,12,13,14,15,16)){
  com.s[,i]=standardize(as.numeric(com.s[,i]))}

#openBUGS
data=list(N=100,death=com.s$death,age=com.s$age,hday=com.s$hday)
Inits1=list(alpha = 0, b.age=0,b.hday=0)
Inits=list(Inits1)

out <- bugs(data = data, inits = Inits, debug=TRUE,
            parameters.to.save = c("alpha","b.age","b.hday"), 
            n.burnin=10,model.file ="cctry.txt", n.chains = 1, n.iter = 100)

И мой cctry.txt - 

model {
for (i in 1:N) {
logit(theta[i]) <- alpha + b.age*age[i] 
+ b.hday*hday[i] 

# Likelihood function for each data point
death[i] ~ dbern(theta[i])
}

#Priors
alpha ~ dnorm(0,.001) # Prior for intercept
b.age~dnorm(0,.00001)
b.hday~dnorm(0,.00001)
}

Я был бы очень признателен за любую помощь!

...