使用R2WinBUGS拟合后怎么绘制后验概率密度

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

使用R2WinBUGS进行线性拟合后，
如何绘制验概率密度（Posterior probability density）
请问大神用哪个语句？
非常感谢，急用
代码如下

1. library(R2WinBUGS)
   
2. library(lattice)
   
3. library(coda)
   
4. dat<-read.csv("c:/data.csv")
   
5. attach(dat)
   
6. setwd("c:/")
   
7. sink("lme.model2.txt")
   
8. cat("
   
9. model {
   
10. 
    
11. # Priors
    
12. for (i in 1:6){
    
13. alpha[i] ~ dnorm(mu.int, tau.int) # Random intercepts
    
14. beta[i] ~ dnorm(mu.slope, tau.slope)# Random slopes
    
15. }
    
16. 
    
17. mu.int ~ dnorm(0, 0.001) # Mean hyperparameter for random intercepts
    
18. tau.int <- 1 / (sigma.int * sigma.int)
    
19. sigma.int ~ dunif(0, 100) # SD hyperparameter for random intercepts
    
20. 
    
21. mu.slope ~ dnorm(0, 0.001) # Mean hyperparameter for random slopes
    
22. tau.slope <- 1 / (sigma.slope * sigma.slope)
    
23. sigma.slope ~ dunif(0, 100) # SD hyperparameter for slopes
    
24. 
    
25. tau <- 1 / ( sigma * sigma) # Residual precision
    
26. sigma ~ dunif(0, 100) # Residual standard deviation
    
27. 
    
28. # Likelihood
    
29. for (i in 1:n) {
    
30. mass[i] ~ dnorm(mu[i], tau)
    
31. mu[i] <- alpha[pop[i]] + beta[pop[i]]* length[i]
    
32. }
    
33. }
    
34. ",fill=TRUE)
    
35. sink()
    
36. 
    
37. # Bundle data
    
38. win.data <- list(mass = as.numeric(mass), pop = pop, length = length, n =length(mass))
    
39. 
    
40. # Inits function
    
41. inits <- function(){ list(alpha = rnorm(6, 0, 2), beta = rnorm(6, 10, 2),
    
42. mu.int = rnorm(1, 0, 1), sigma.int = rlnorm(1), mu.slope = rnorm(1, 0, 1),
    
43. sigma.slope = rlnorm(1), sigma = rlnorm(1))}
    
44. 
    
45. # Parameters to estimate
    
46. parameters <- c("alpha", "beta", "mu.int", "sigma.int", "mu.slope", "sigma.slope", "sigma")
    
47. 
    
48. ni <-10000
    
49. nb <- 1000
    
50. nt <- 2
    
51. nc <- 3
    
52. 
    
53. # Start Gibbs sampling
    
54. #out <- bugs(win.data, inits, parameters, "lme.model1.txt", n.thin=nt, n.chains=nc, n.burnin=nb, n.iter=ni, debug = TRUE，bugs.directory="c:/Program Files/WinBUGS14/")
    
55. 
    
56. out1 <- bugs(win.data, inits, parameters, "lme.model2.txt", n.thin=nt,
    
57. n.chains=nc, n.burnin=nb, n.iter=ni, debug = TRUE)

复制代码

数据在附件中

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：R2WinBUGS winbugs WINBUG 后验概率 概率密度 library Random 如何

http://www.r-tutor.com/bayesian-statistics/openbugs

hope this can help you