I have a set of True/False observations. Each of those observations is the result of a process that generates a set of items, and yields True if *any* of those individual items is true. Which of the items that are in the set changes for each observation.

I have modelled this with something that looks like:

```
data {
// Which sub_items were present for observation i:
int<lower=0, upper=1> item_occurences[n_observations, n_items]
...
}
parameters {
// Individual sub_item probabilities for each sub-item:
real<lower=0, upper=1> item_probabilities[n_items]
}
model {
real p;
for (i in 1:n_observations) {
// Probability that at least 1 item is true
// i.e. probability that all items are not false
// sub_items holds the individual success
// probabilities for each observation
p = 1 - prod(1 - (item_occurences[i] .* item_probabilities));
observation[i] ~ bernoulli(p);
}
}
```

However, this is struggling (failing) to recover parameters in a simulation test.

Is there a better way to do this?