2011-08-03 14 views
3

私はJAGSを初めてお使いになり、打ち切りデータ用のJAGSでdinterval()がどのように動作するかを理解しようとしています。私は、各データポイント(真の値ではない)の上限と下限を持つ粗いデータをモデリングしています。ここで私はそれが動作するはずだと思うかの簡単な例です:間隔検閲データのdinterval()?

各点のいくつかの上限と下限は:

> head(lim) 
     L  U 
[1,] 14.98266 15.68029 
[2,] 21.21827 21.91590 
[3,] 18.34953 19.04716 
[4,] 19.00186 19.69949 
[5,] 15.39891 16.09654 
[6,] 17.81705 18.51468 

モデルを記述するための機能(データを想定して共通で通常から来ます平均と分散):ジャグため

playmodel <- function(){ 
      for (i in 1:50){ 
       is.censored[i] ~ dinterval(t[i], lim[i,]) 
       t[i] ~ dnorm(mu,tau) 
       } 
      mu ~ dnorm(0,.001) 
      tau ~ dgamma(.01,.01) 
      } 
      filename <- "toymod.bug" 
      write.model(toymod,filename) 

いくつかの機能と割り当ては呼び出し:

data <- list("lim"=lim) 
inits <- list(mu=rnorm(1),tau=rgamma(1,.01,.01),t=as.vector(apply(lim,1,mean))) 
#last part is to ensure the starting value is between the upper and lower limit 
#each chain will start at the same place for t but this is just for this case 
params <- c("mu","tau") 

モデルを実行してください:

playmodel.jags <- jags(data,inits, params, model.file="toymod.bug", n.chains=3, 
        n.iter=50000,n.burnin=30000, n.thin=1, DIC=TRUE, 
        working.directory=NULL,refresh = 50000/50, progress.bar = "text") 

これを実行するとどうなりますか?

1)ムーの私の推定値は、それが実行されません15

2)であるべき時に右周りに0に置くとDIC = TRUE:

error: "Error in jags.samples(model, variable.names, n.iter, thin, type = "trace", : Failed to set trace monitor for node deviance

私は愚かな何かをやっている確信していると誰かが私を軌道に乗せるのを助けることができれば感謝しています。

+1

また、I(Lower、Upper)関数を使用してOpenBUGSでモデル化した場合、正常に動作するようです。 – scottyaz

+0

私は統計物理学をどのくらい扱っているかわかりません。これは、stats.stackexchange.comでのクロス投稿が良いアイデアになることはほとんどありません。あるいは、パッケージの著者に連絡するか、パッケージのメーリングリストがあるかどうかを確認することができます。 –

答えて

2

は、以下のマーティン・プラマーからの応答である:

書かれたとおり、あなたのモデルは、任意の観察結果を持っていません。あなたはおそらく、それが本当に速く走っていることに気づいたでしょう。これは、前回のサンプリングから順方向サンプリングであるためです。これはmuの事後平均が以前の平均0と同じ理由です。変数名 "is.censored"は、生存分析で見られるように、左検診または右検診のデータに適していますが、あなたの問題ではありません。だから私はそれを "y"と改名するつもりです。あなたは

y[j] ~ dinterval(t[j], lim[j,]) 

とLIM [j]を持っている場合は2つの列があり、その後、Y [j]の区間打ち切りデータをモデル化するための3つの値

y[j] = 0 if t[j] <= lim[j,1] 
y[j] = 1 if lim[j,1] < t[j] <= lim[j,2] 
y[j] = 2 if lim[j,2] < t[j] 

を取ることができ、あなたはY [j]を供給する必要がありますあなたのモデルのデータとして。あなたの場合、t [j]は常にlim [j、1]とlim [j、2]の間にあるので、あなたのデータはそうでなければなりません。

data <- list("lim"=lim, "y"=rep(1,nrow(lim))) 

DICの問題はかなり深いです。モデルに成果データがないため、偏差は定義されていません。しかし、結果データを提供しても、希望する偏差統計(pDを含む)は得られません。逸脱はゼロになり、 "ジャグ"機能はpDのGelmanヒューリスティックに戻ってしまいます(私はこれを書いていないので、説明しないようにしています)。これもゼロになります。あなたが本当にしたい可能性が

p(lim[j,1] < t[j] <= lim[j,2] | mu, tau) 

ですがぎざぎざはあなたに常に1 DICの「フォーカス」で間違っている

p(y[j] | t[j]) 

を与えています。 WinBUGSがこのような状況で何をするのか分かりません。おそらく、それは検閲された変数のための特別な規則を持っています。