0
$\begingroup$

I am trying to write jags code for the following scenario

  1. Toss a coin with unknown probability of heads (p)
  2. If heads, then draw a random integer from 1 to 12
  3. If tails, then draw a random integer from 1 to 6

Based on data of the integers drawn, compute posterior distribution of p, assuming uniform prior on (0,1).

Here is my code

model = "model
{
  p ~ dunif(0, 1) # prior for coin
  for (i in 1 : n) {
    coin[i] ~ dbern(p) 
    temp[i] = (coin[i]) * (-1/12) + (1/6) # if coin=0(tails),then temp=1/6, if coin = 1(heads), then temp = 1/12
    pi = rep(temp[i], 1/temp[i]) # vector of probabilities, i.e. rep(1/12, 12) or rep(1/6, 6)
    number[i] ~ dcat(pi)
}

}"

When I run this model JAGS throws the following error: RUNTIME ERROR: Expected parameters with fixed values in function rep

If passing a stochastic node to rep is not allowed, then how can I simulate the above scenario?

$\endgroup$
0
$\begingroup$

One way to do this is to think about your vector of probabilities as being one of a pair of fixed possibilities that are known in advance, and to pass them in as data: specifically a matrix with 2 columns where the first column is for the tails situation and the second column is for the heads situation. Also note that both columns can have 12 probabilities - just set the last 6 values to zero for the first column. Then all you need to do is use nested indexing to select the appropriate column within dcat. For example:

m <- 'model{

    for(i in 1:N){      
        coin[i] ~ dbern(p)
        Number[i] ~ dcat(Weights[,coin[i]+1])               
    }

    p ~ dbeta(1,1)

    #data# N, Number, Weights
    #monitor# p

}'


Weights <- matrix(ncol=2, nrow=12)
Weights[1:6,1] <- 1/6
Weights[7:12,1] <- 0
Weights[1:12,2] <- 1/12

N <- 100
p <- 0.75
coin <- rbinom(N, 1, p)
Number <- sapply(c(6,12)[coin+1], function(x) sample(1:x, 1))

library('runjags')
run.jags(m)

This is effectively a type of mixture model if it helps to think about it in those terms (specifically a mixture of categorical distributions).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.