रेजेस पैकेज का उपयोग कर बायेसियन मल्टीमोनियल रिग्रेशन

Aug 17 2020

मैं एक फिट करने के लिए कोशिश कर रहा हूँ बहुपद रसद प्रतिगमन मॉडल का उपयोग कर rjagsपरिणाम के लिए एक स्पष्ट (सांकेतिक) चर (है परिणाम ) 3 स्तरों के साथ, और व्याख्यात्मक चर हैं उम्र (निरंतर) और समूह (3 स्तर के साथ स्पष्ट)। ऐसा करने के लिए, मैं पोस्टीरियर साधन और आयु और समूह के लिए 95% मात्रात्मक-आधारित क्षेत्र प्राप्त करना चाहूंगा ।

मैं वास्तव में महान नहीं for loopहूं जिस पर मुझे लगता है यही कारण है कि मॉडल के लिए मेरा लिखित कोड ठीक से काम नहीं कर रहा है।

मेरे बीटा पादरी एक सामान्य वितरण का पालन करते हैं, prij 0 सामान्य (0,100) j 0 {0, 1, 2} के लिए।

रीप्रोड्यूसिएबल आर कोड

library(rjags)

set.seed(1)
data <- data.frame(Age = round(runif(119, min = 1, max = 18)),
                   Group = c(rep("pink", 20), rep("blue", 18), rep("yellow", 81)), 
                   Outcome = c(rep("A", 45), rep("B", 19), rep("C", 55)))

X <- as.matrix(data[,c("Age", "Group")]) 
J <- ncol(X)
N <- nrow(X)

## Step 1: Specify model
cat("
model {
for (i in 1:N){

    ##Sampling model
    yvec[i] ~ dmulti(p[i,1:J], 1)
    #yvec[i] ~ dcat(p[i, 1:J])  # alternative
    for (j in 1:J){
      log(q[i,j]) <- beta0 + beta1*X[i,1] + beta2*X[i,2] 
      p[i,j] <- q[i,j]/sum(q[i,1:J])  
    } 
    
    ##Priors
    beta0 ~ dnorm(0, 0.001)
    beta1 ~ dnorm(0, 0.001)
    beta2 ~ dnorm(0, 0.001)
}
}",
file="model.txt")

##Step 2: Specify data list 
dat.list <- list(yvec = data$Outcome, X=X, J=J, N=N) 

## Step 3: Compile and adapt model in JAGS 
jagsModel<-jags.model(file = "model.txt",
                      data = dat.list,
                      n.chains = 3,
                      n.adapt = 3000
)

त्रुटि संदेश :

सूत्रों की मदद के लिए मैं देख रहा हूँ :

http://people.bu.edu/dietze/Bayes2018/Lesson21_GLM.pdf

श्रेणी X के साथ JAGS में डिरिचलेट बहुराष्ट्रीय मॉडल

संदर्भ सेhttp://www.stats.ox.ac.uk/~nicholls/MScMCMC15/jags_user_manual.pdf, पृष्ठ ३१

मैंने अभी सीखना शुरू किया है कि rjagsपैकेज का उपयोग कैसे करें ताकि किसी भी संकेत / स्पष्टीकरण और प्रासंगिक स्रोतों से लिंक की बहुत सराहना की जाए!

जवाब

3 Duck Aug 17 2020 at 23:20

मैं आपके मुद्दे पर एक दृष्टिकोण शामिल करूंगा। मैंने उन्हीं पादरियों को लिया है जिन्हें आपने गुणांक के लिए परिभाषित किया है। मुझे केवल यह उल्लेख करने की आवश्यकता है कि जैसा कि आपके पास एक कारक है, Groupमैं इसके एक स्तर को संदर्भ (इस मामले में pink) के रूप में उपयोग करूंगा, इसलिए इसका प्रभाव मॉडल में निरंतर द्वारा ध्यान में रखा जाएगा। अगला कोड:

library(rjags)
#Data
set.seed(1)
data <- data.frame(Age = round(runif(119, min = 1, max = 18)),
                   Group = c(rep("pink", 20), rep("blue", 18), rep("yellow", 81)), 
                   Outcome = c(rep("A", 45), rep("B", 19), rep("C", 55)))

#Input Values we will avoid pink because it is used as reference level
#so constant absorbs the effect of that level
r1 <- as.numeric(data$Group=='pink') r2 <- as.numeric(data$Group=='blue')
r3 <- as.numeric(data$Group=='yellow') age <- data$Age
#Output 2 and 3
o1 <- as.numeric(data$Outcome=='A') o2 <- as.numeric(data$Outcome=='B')
o3 <- as.numeric(data$Outcome=='C')
#Dim, all have the same length
N <- length(r2)

## Step 1: Specify model

model.string <- "
model{
for (i in 1:N){ 

## outcome levels B, C
o1[i] ~ dbern(pi1[i])
o2[i] ~ dbern(pi2[i]) 
o3[i] ~ dbern(pi3[i]) 

## predictors
logit(pi1[i]) <- b1+b2*age[i]+b3*r2[i]+b4*r3[i]
logit(pi2[i]) <- b1+b2*age[i]+b3*r2[i]+b4*r3[i]
logit(pi3[i]) <- b1+b2*age[i]+b3*r2[i]+b4*r3[i]

} 
## priors
b1 ~ dnorm(0, 0.001)
b2 ~ dnorm(0, 0.001)
b3 ~ dnorm(0, 0.001)
b4 ~ dnorm(0, 0.001)
}
"
#Model
model.spec<-textConnection(model.string)

## fit model w JAGS
jags <- jags.model(model.spec,
                   data = list('r2'=r2,'r3'=r3,
                               'o1'=o1,'o2'=o2,'o3'=o3,
                               'age'=age,'N'=N),
                   n.chains=3,
                   n.adapt=3000)

#Update the model
#Update
update(jags, n.iter=1000,progress.bar = 'none')
#Sampling
results <- coda.samples(jags,variable.names=c("b1","b2","b3","b4"),n.iter=1000,
                        progress.bar = 'none')
#Results
Res <- do.call(rbind.data.frame, results)

सहेजे गए मापदंडों के लिए जंजीरों के परिणामों के साथ Res, आप अगले कोड का उपयोग करके पीछे के मीडिया और विश्वसनीय अंतराल की गणना कर सकते हैं:

#Posterior means
apply(Res,2,mean)

         b1          b2          b3          b4 
-0.79447801  0.00168827  0.07240954  0.08650250

#Lower CI limit
apply(Res,2,quantile,prob=0.05)

         b1          b2          b3          b4 
-1.45918662 -0.03960765 -0.61027923 -0.42674155

#Upper CI limit
apply(Res,2,quantile,prob=0.95)

         b1          b2          b3          b4 
-0.13005617  0.04013478  0.72852243  0.61216838 

bमापदंडों चर पर विचार किया (में से प्रत्येक के हैं ageऔर के स्तर Group)। मिश्रित श्रृंखलाओं के कारण अंतिम मूल्य बदल सकते हैं!