R2OpenBUGSをためす [統計]
Windows版OpenBUGSとDarwine 1.0で、R2OpenBUGSをためしてみた。マニュアルの例をやってみる。
library(R2OpenBUGS) # use Darwine 1.0 WINE <- "/Applications/Darwine/Wine.bundle/Contents/bin/wine" WINEPATH <- "/Applications/Darwine/Wine.bundle/Contents/bin/winepath" OpenBUGS <- paste(Sys.getenv("HOME"), ".wine/drive_c/Program\ Files/OpenBUGS", "OpenBUGS321/OpenBUGS.exe", sep = "/") model.file <- system.file(package="R2OpenBUGS", "model", "schools.txt") data(schools) schools J <- nrow(schools) y <- schools$estimate sigma.y <- schools$sd data <- list ("J", "y", "sigma.y") inits <- function(){ list(theta = rnorm(J, 0, 100), mu.theta = rnorm(1, 0, 100), sigma.theta = runif(1, 0, 100)) } parameters <- c("theta", "mu.theta", "sigma.theta") schools.sim <- bugs(data, inits, parameters, model.file, n.chains = 3, n.iter = 5000, working.directory = getwd(), useWINE = TRUE, WINE = WINE, WINEPATH = WINEPATH, OpenBUGS.pgm = OpenBUGS) print(schools.sim) png("Rplot.png") plot(schools.sim) dev.off()
Inference for Bugs model at "/Library/Frameworks/R.framework/Versions/2.13/Resources/library/R2OpenBUGS/model/schools.txt", Current: 3 chains, each with 5000 iterations (first 2500 discarded) Cumulative: n.sims = 7500 iterations saved mean sd 2.5% 25% 50% 75% 97.5% Rhat n.eff theta[1] 11.2 8.9 -2.8 5.5 9.8 15.7 32.9 1 410 theta[2] 7.5 6.5 -4.9 3.4 7.4 11.5 20.8 1 420 theta[3] 5.8 8.0 -12.3 1.3 6.4 10.5 21.0 1 1200 theta[4] 7.1 6.7 -6.3 3.0 7.1 11.2 20.9 1 610 theta[5] 4.9 6.3 -8.8 1.0 5.5 9.1 16.4 1 680 theta[6] 5.8 6.8 -9.2 1.7 6.1 10.1 18.3 1 490 theta[7] 10.4 7.2 -2.1 5.6 9.7 14.6 26.3 1 300 theta[8] 8.1 8.0 -7.1 3.4 7.8 12.4 25.7 1 420 mu.theta 7.6 5.4 -2.7 4.1 7.6 11.0 18.6 1 310 sigma.theta 6.7 5.9 0.2 2.3 5.3 9.4 21.8 1 540 deviance 60.5 2.3 57.0 59.2 60.1 61.6 66.2 1 1600 For each parameter, n.eff is a crude measure of effective sample size, and Rhat is the potential scale reduction factor (at convergence, Rhat=1). DIC info (using the rule, pD = Dbar-Dhat) pD = 63.4 and DIC = 2.9 DIC is an estimate of expected predictive error (lower deviance is better).
コメント 0