3

私はよく、既知のパラメーターを使用してシミュレートされたデータで JAGS モデルを実行します。オブジェクトのデフォルトのプロット方法が気に入っていますmcmc。abline(v=TRUE_VALUE)ただし、モデル化された for each パラメータを追加したいと思います。これにより、事後分布が妥当かどうかを簡単に確認できます。

もちろん、これを手動で行うことも、おそらく車輪を再発明して独自の関数を作成することもできます。しかし、既存の方法に基づいたエレガントな方法があるかどうか疑問に思っていましたplot。

これが実際の例です：

require(rjags)
require(coda)

# simulatee data
set.seed(4444)
N <- 100 
Mu <- 100
Sigma <- 15
y <- rnorm(n=N, mean=Mu, sd=Sigma)
jagsdata <- list(y=y)

jags.script <- "
model {
    for (i in 1:length(y)) {
        y[i]  ~ dnorm(mu, tau)
    }
    mu  ~ dnorm(0, 0.001)    
    sigma ~ dunif(0, 1000)
    tau <- 1/sigma^2
}"


mod1 <- jags.model(textConnection(jags.script), data=jagsdata, n.chains=4, 
                   n.adapt=1000)
update(mod1, 200) # burn in
mod1.samples <- coda.samples(model=mod1,
                             variable.names=c('mu', 'sigma'),
                             n.iter=1000)                  
plot(mod1.samples)

ここに画像の説明を入力

abline(v=100)mu やabline(v=15)sigmaのようなものを実行したいだけです。もちろん、他の多くの例では、5、10、20、またはそれ以上の関心のあるパラメーターがあります。したがって、名前付きパラメーターに真の値のベクトルを提供できることに興味があります。

私は見てきたgetAnywhere(plot.mcmc)。それを変更することは良い方法でしょうか？

4

1 に答える 1

7

わかった。そこで、次plot.mcmcのように変更しました。

my.plot.mcmc <- function (x, trace = TRUE, density = TRUE, smooth = FALSE, bwf, 
    auto.layout = TRUE, ask = FALSE, parameters, ...) 
{
    oldpar <- NULL
    on.exit(par(oldpar))
    if (auto.layout) {
        mfrow <- coda:::set.mfrow(Nchains = nchain(x), Nparms = nvar(x), 
            nplots = trace + density)
        oldpar <- par(mfrow = mfrow)
    }
    for (i in 1:nvar(x)) {
        y <- mcmc(as.matrix(x)[, i, drop = FALSE], start(x), 
            end(x), thin(x))
        if (trace) 
            traceplot(y, smooth = smooth, ...)
        if (density) {
            if (missing(bwf)) {
                densplot(y, ...); abline(v=parameters[i])
            } else densplot(y, bwf = bwf, ...)
        }
        if (i == 1) 
            oldpar <- c(oldpar, par(ask = ask))
    }
}

次に、コマンドを実行します

my.plot.mcmc(mod1.samples, parameters=c(Mu, Sigma))

これを生成します

ここに画像の説明を入力

parametersは、JAGS が変数をソートするのと同じソート順の値のベクトルでなければならないことに注意してください。

教訓

  • plot.mcmcおそらく名前空間が原因で、単に new を書くだけではデフォルトでは機能しませんでした。だから私はちょうど新しい関数を作成しました
  • おそらく名前空間のためにも変更set.mfrowする必要がありました。coda:::set.mfrow
  • ask=FALSERStudio では Figure の閲覧が許可されているため、ask を に変更しました。

既存の S3 メソッドをオーバーライドまたは適応させるためのより良い方法についての提案をお待ちしております。

于 2012-06-07T06:32:48.050 に答える