更新时间:2023-02-14 10:21:10
我不知道这是否可以满足您的需求.请注意,model
代码来自于使用您的代码,然后在光标处键入LINE
.其余的只是标准的bug代码,除了我使用tau = rgamma(1,1)
作为初始值,并且不知道它的标准是什么.我不止一次看到tau = 1
用作初始值.也许那会更好.
I do not know whether this will give you what you want. Note that the model
code came from using your code and then typing LINE
at the cursor. The rest is just standard bugs code, except I used tau = rgamma(1,1)
for an initial value and do not know how standard that is. More than once I have seen tau = 1
used as an initial value. Perhaps that would be better.
实际上,我使用与您使用的相同的model
代码创建了rjags
对象,并添加了jags
语句来运行它.我承认这与将尾声输出转换为bugs
对象不同,但这可能会导致您获得所需的plot
.
In effect, I created an rjags
object using the same model
code you were using and added a jags
statement to run it. I admit that is not the same thing as converting coda output to a bugs
object, but it might result in you getting the desired plot
.
如果您只有一个mcmc.list
而不是一个model
代码,而您只想绘制mcmc.list
,那么我的回答将无济于事.
If all you have is an mcmc.list
and no model
code and you simply want to plot the mcmc.list
, then my answer will not help.
library(R2jags)
x <- c(1, 2, 2, 4, 4, 5, 5, 6, 6, 8)
Y <- c(7, 8, 7, 8, 9, 11, 10, 13, 14, 13)
N <- length(x)
xbar <- mean(x)
summary(lm(Y ~ x))
x2 <- x - xbar
summary(lm(Y ~ x2))
# Specify model in BUGS language
sink("model1.txt")
cat("
model {
for( i in 1 : N ) {
Y[i] ~ dnorm(mu[i],tau)
mu[i] <- alpha + beta * (x[i] - xbar)
}
tau ~ dgamma(0.001,0.001)
sigma <- 1 / sqrt(tau)
alpha ~ dnorm(0.0,1.0E-6)
beta ~ dnorm(0.0,1.0E-6)
}
",fill=TRUE)
sink()
win.data <- list(Y=Y, x=x, N=N, xbar=xbar)
# Initial values
inits <- function(){ list(alpha=rnorm(1), beta=rnorm(1), tau = rgamma(1,1))}
# Parameters monitored
params <- c("alpha", "beta", "sigma")
# MCMC settings
ni <- 25000
nt <- 5
nb <- 5000
nc <- 3
out1 <- jags(win.data, inits, params, "model1.txt", n.chains = nc,
n.thin = nt, n.iter = ni, n.burnin = nb)
print(out1, dig = 2)
plot(out1)
#library(R2WinBUGS)
#plot(out1)
基于评论,也许类似的事情会有所帮助. str(new.data)
行表明有大量数据可用.如果您只是尝试创建默认图的变体,则仅需根据需要提取和设置数据即可.这里plot(new.data$sims.list$P1)
只是一个简单的示例.在不确切知道您想要什么绘图的情况下,我将不会尝试更具体的数据提取.如果您发布的图显示了确切类型的绘图示例,则您可能希望有人可以从此处获取它并发布创建它所需的代码.
Based on the comments perhaps something like this will help. The line str(new.data)
suggests that a large amount of data are available. If you are simply trying to create variations of default plots then doing so may only be a matter of extracting and subsetting the data as desired. Here plot(new.data$sims.list$P1)
is just one straight-forward example. Without knowing exactly what plot you want I will not attempt more specific data extractions. If you post a figure showing an example of the exact kind of plot you want perhaps someone can take it from here and post the code needed to create it.
顺便说一句,我建议将示例数据集的大小减少到大约三个链,并且不超过30次迭代,直到您拥有所需的精确绘图所需的确切代码为止:
By the way, I recommend reducing the size of the example data set to perhaps three chains and perhaps no more than 30 iterations until you have the exact code you want for the exact plot you want:
load("C:/Users/mmiller21/simple R programs/test.mcmc.list.Rdata")
class(test.mcmc.list)
library(R2WinBUGS)
plot(as.bugs.array(sims.array = as.array(test.mcmc.list)))
new.data <- as.bugs.array(sims.array = as.array(test.mcmc.list))
str(new.data)
plot(new.data$sims.list$P1)
还请注意:
class(new.data)
[1] "bugs"
而:
class(test.mcmc.list)
[1] "mcmc.list"
这是您的帖子标题所要求的.
which is what the title of your post requests.