n番目の「無効な親値」エラー

Question

私はRjagsで今実行している非常に簡単なモデルで私を助けてくれることを願っています。 次のように私が持っているデータは、次のとおりです。

> print(data) 
$R 
225738  184094  66275  24861  11266 
228662  199379  70308  27511  12229 
246808  224814  78255  30447  13425 
254823  236063  83099  33148  13961 
263772  250706  89182  35450  14750 
272844  262707  96918  37116  15715 
280101  271612  102604  38692  16682 
291493  283018  111125  40996  18064 
310474  299315  119354  44552  19707 
340975  322054  126901  47757  21510 
347597  332946  127708  49103  21354 
354252  355994  130561  51925  22421 
366818  393534  140628  56562  23711 
346430  400629  146037  59594  25313 
316438  399545  150733  62414  26720 
303294  405876  161793  67060  29545 

$N 
9597000  8843000  9154000  9956000 11329000 
9854932  9349814  9532373 10195193 11357751 
9908897  9676950  9303113 10263930 11141510 
9981879  9916245  9248586 10270193 10903446 
10086567 10093723  9307104 10193818 10660101 
10242793 10190641  9479080 10041145 10453320 
10434789 10222806  9712544  9835154 10411620 
10597293 10238784 10014422  9611918 10489448 
10731326 10270163 10229259  9559334 10502839 
10805148 10339566 10393532  9625879 10437809 
10804571 10459413 10466871  9800559 10292169 
10696317 10611599 10477448 10030407 10085603 
10540942 10860363 10539271 10245334  9850488 
10411836 11053751 10569913 10435763  9797028 
10336667 11152428 10652017 10613341  9850533 
10283624 11172747 10826549 10719741  9981814 

$n 
[1] 16 

$na 
[1] 5 

$pbeta 
[1] 0.70 0.95 

をし、次のようにモデルがある：

cat('model{ 
## likelihoods ## 
for(k in 1:na){ for(w in 1:n){ R[w,k] ~ dbin(theta[w,k], N[w,k]) }} 

for(k in 1:na){ for(w in 1:n){ theta[w,k] <- 0.5*beta[w,k]*0.5 }} 

for(k in 1:na){ 
    beta[1,k] ~ dunif(pbeta[1], pbeta[2])   
    beta.plus[1,k] <- beta[1,k]        
    for (w in 2:n){        
    beta.plus[w,k] ~ dunif(beta[(w-1),k], 0.95)    
    beta[w,k] <- beta.plus[w,k]} } }', 
file='model1.bug') 

######## initial random values for beta 
bbb=bb.plus=matrix(rep(NA, 16*5), byrow=T, ncol=5); 
for(k in 1:5){ 
    bbb[1,k]=runif(1, 0.7,0.95); 
    for (w in 2:16){ 
     bb.plus[w,k] = runif(1, bbb[w-1,k], 0.95); 
     bbb[w,k]=bb.plus[w,k]} } 

## data & initial values 
inits1 <- list('beta'= bbb) 

jags_mod <- jags.model('model1.bug', data=data, inits=inits1, n.chains=1, n.adapt=1000) 
update(jags_mod, n.iter=1000) 
posts=coda.samples(model=jags_mod,variable.names=c('beta','theta'), n.iter=niter, thin=1000) 

超簡単。これは実際に私がここに来るのとまったく同じエラーメッセージを与えるより複雑なモデルからの縮小モデルです。 このモデルを実行するたびに、何の問題もありません。 パラメータβのプライオリティは、0.7から0.95に増加するように書かれています。

ここで、モデルの最初の行をコメントアウトすることで、Rの可能性を「シャットオフ」したいと思います。どのようにパラメータのテータがどのように推定されるかを見るために、私はそうしたいと思います（基本的に私はtheta = beta/4を見つけなければなりませんが、それは私にとってうまくなります）。

私は、私は、一般にマトリックスの最下行（行15または16）にパラメータベータの「無効な親」エラーが表示されます。 実際これはこれよりも洗練されています。時にはエラーが発生することもありますが、時々私は（ほとんど、私は）しません。

私は理解していませんなぜこれが起こるのですか：ベータの値は、尤度の有無から独立して生成されるべきではありませんか？

ご迷惑をおかけして申し訳ございませんが、ご迷惑をおかけ致します。 ありがとう、最高！ Emanuele

Answer 1

モデルをもう少し遊んだ後で、問題の原因がわかりました。一様分布の1つの必要な態様（すなわち、unif(a,b)）は、a < bである。あなたのモデル内で均一分布を小さくすると、は、bに近づくようになります。時々、それには達しませんが、他の時にはaはbと等しくなり、無効な親値のエラーが発生します。たとえば、モデルに次のものが含まれているとします。

example ~ dunif(0.4,0.4) 

「ノードの例でエラー、無効な親値」が表示されます。

これを解決するために、プリオーサーをランダムに割り当てるのではなく、指定する方法を調整する方が簡単だと思います。これをベータ版で行うことができます。最初のステップでは、beta(23.48, 4.98)が0.7から0.95までの範囲のほとんどをカバーしますが、その範囲内にあることを確認するために切り捨てることができます。その後、nが増加すると、4.98を下げることができ、前者は0.95に縮小されます。以下のモデルがこれを行います。前任者を検査した後、thetaはbeta/4と等しくなることが分かります。

私は本当にあなたがベータ分布は、かどうかはわからないためのモデルを使用しようとしているかわからないので、私たちは

 beta[1,1] theta[1,1] 
[1,] 0.9448125 0.2362031 
[2,] 0.7788794 0.1947198 
[3,] 0.9498806 0.2374702 

0.9448125/4 
[1] 0.2362031 

を取得し、このモデルからの出力のいくつかを見てみると

data.list <- list(n = 16, na = 5, 
        B = rev(seq(0.1, 4.98, length.out = 16))) 

cat('model{ 
## likelihoods ## 
    #for(k in 1:na){ for(w in 1:n){ R[w,k] ~ dbin(theta[w,k], N[w,k]) }} 

    for(k in 1:na){ for(w in 1:n){ theta[w,k] <- 0.5*beta[w,k]*0.5 }} 
    for(k in 1:na){ 
    for(w in 1:n){ 
    beta[w,k] ~ dbeta(23.48, B[w]) T(0.7,0.95)   

    } } }', 
    file='model1.bug') 

jags_mod <- jags.model('model1.bug', data=data.list, 
inits=inits1, n.chains=1, n.adapt=1000) 

update(jags_mod, n.iter=1000) 

posts=coda.samples(model=jags_mod, 
variable.names=c('beta','theta'), n.iter=10000, thin=10) 

あなたのニーズに合っていますが、上記の方法はあなたがしようとしていることを模倣します。