Mixture model fails on cmdstanr when conditioning theta on participant, but not on rstan

I am playing with the data in https://doi.org/10.1371/journal.pcbi.1005833 , where different kinds of primates play a “matching pennies” game (guess in which hand is the tasty treat).

I define a handful of possible strategies and want to build a mixture model trying to identify which strategy each primate is mostly likely using.

The two simplest ones (as an example here):

Strategy 1: Simple Bias: the primate just picks the left hand with a fixed rate. E.g. the primate is left handed and tends to pick left 80% of the time.
Strategy 2: Biased Nash: the primate keeps track of how many times the food was in the left hand so far, and picks accordingly (more likely to pick left, the higher the proportion of times food was in the left hand)

I build the models on a single primate and all works nicely. I then try to scale up to all the primates in a given species. The single strategy models work nicely. The mixture model fails, plausibly due to how I specify that theta (the probability of picking strategy 1) could differ by primate.

Data: https://www.dropbox.com/s/dfjaxgoha1x5pw4/sampledata.csv?dl=0

# Setting the family as a mixture of two bernoulli
mix <- mixture(bernoulli, bernoulli, order = "none")

# Defining the two bernoullis and conditioning theta on individual
mix_f1 <- bf(
  left ~ 1 ,
  mu1 ~ 1 + (1 | individual),
  mu2 ~ 1 + LeftProb + (1 + LeftProb | individual),
  theta1 ~ 1 + (1|individual)
)

# defining priors
mix_p1 <- c(
  prior(normal(0,0.3), class=Intercept, dpar=mu1),
  prior(normal(0,0.3), class=Intercept, dpar=mu2),
  prior(normal(0,0.3), class=b, dpar=mu2),
  prior(normal(0, .3), class = sd, dpar=mu1),
  prior(normal(0, .3), class = sd, dpar=mu2),
  prior(normal(0, .3), class = sd, dpar=theta1),
  prior(lkj(3), class = cor)
)

# Fittng the model
mix_m1 <- brm(
  mix_f1,
  d,
  family = mix,
  mix_p1,
  sample_prior = T,
  chains = 1,
  cores = 2,
  backend = "cmdstanr",
  threads = threading(2),
  file = paste0("mix_m1"),
  file_refit = "on_change",
  control = list(
    adapt_delta = 0.99,
    max_treedepth = 20
  )
)

Commenting away the backend all works. but with cmstanr (no matter if threading or not), I get

Compiling Stan program…
Start sampling
Running MCMC with 1 chain, with 2 thread(s) per chain…

Chain 1 num_threads = 2
Warning: Chain 1 finished unexpectedly!

Error in rstan::read_stan_csv(out$output_files()) :
csvfiles does not contain any CSV file name
In addition: Warning message:
No chains finished successfully. Unable to retrieve the fit.

I am unable to figure out why this is. Any suggestion?

Hi,
we will probably need a bit more details on what happened. brms doesn’t directly support detailed output from cmdstanr, but you can get it via:

scode <- make_stancode(... parameters to brm ....)
sdata <- make_standata(... parameters to brm ....)
m <- cmdstan_model(write_stan_file(scode))
class(sdata) <- NULL # Necessary for stupid reasons
fit <- m$sample(data = as.list(sdata))
fit$output()

I tried running your model, but the data in the .csv file seem to be in a very different format than what the model expects so I was not able to run the model directly.

Best of luck with the model

2 Likes

thanks. I load the data as

d <- readr::read_csv("SampleData.csv")

and it works, sorry for not having specified that.

When following your code I get a bunch of:

Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:

Chain 1 Exception: Exception: bernoulli_logit_lpmf: Logit transformed probability parameter is nan, but must be not nan! (in ‘/var/folders/3m/f039n0x549vfxhdtj55yykzhfjr0d6/T/Rtmprg7Jbb/model-4b44ac60106.stan’, line 68, column 6 to column 63) (in ‘/var/folders/3m/f039n0x549vfxhdtj55yykzhfjr0d6/T/Rtmprg7Jbb/model-4b44ac60106.stan’, line 68, column 6 to column 63)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.

but then the model runs. When calling things via brms, the model still fails as in the original post. So it seems to be an issue in how brms and cmdstanr interact?

1 Like

So this one is on me - I looked at preview of the data and it looked very different, but in fact it runs, sorry for the confusion.

Anyway, I confirm that I get the same behaviour on my computer. It’s super weird - brms should almost directly just call make_standata and make_stancode. It would be great if you could try to see how much you can simplify the model to still exhibit the problem and report this as an issue at Issues · paul-buerkner/brms · GitHub.

Do you need a workaround for this or can you use rstan backend for now?

1 Like

I can make do with rstan, but yeah, I’d like to understand what’s going on.
If I remember correctly the model without random effects did work in cmdstanr via brms, so this is the simplest model w the problem . But I’ll double check and open the issue. Thanks again!

For the record, the github issue is posted here: Mixture model fails on brms using cmdstanr, but not using cmdstanr directly · Issue #1150 · paul-buerkner/brms · GitHub

1 Like

the bug has now been fixed: Mixture model fails on brms using cmdstanr, but not using cmdstanr directly · Issue #1150 · paul-buerkner/brms · GitHub

1 Like