Help with error in posterior_epred for AR1 model

Interfaces brms
, ,
1

Hi folks,

I’m fitting an AR1 model with a single binary condition, and I want to get an individual-level estimate of the effect of condition, shrunk hierarchically as is normal. The model structure is something like:

DV ~ Condition + ar(time = TrialNumber, gr = Subject, p = 1, cov = T) +(1+Condition|Study/Subject))

The model fits great, but then when I feed it a dataframe with one row per subject/study combination with a column for Condition, it complains about not having TrialNumber; when I give each subject a dummy trial number (1 for one level of Condition, 2 for the other) it says trial number must be unique within grouping (which it is, in the new dataframe). When I throw up my hands and just feed it the original dataframe it was trained on to try to get something to work, it gives the error “Error in eval(substitute(expr), data, enclos = parent.frame()) : object ‘begin_tg’ not found” (note that begin_tg isn’t part of my analysis, and I don’t know where it comes from).
Would anyone be willing to give me advice on how to get hierarchically-shrunk estimates of the effect of Condition on each subject out of this model? Currently I’ve resorted to manually grabbing the samples, combining the by-subject-by-study random slope of Condition and the main effect of Condition and re-constructing a linear predictor that way, but that yields very wide estimates which don’t seem terribly helpful (I have reason to suspect the estimates should be at least a little shrunk, since there’s a paper in the literature that uses this technique on the same sort of data to pretty nice effect). Any suggestions would be greatly appreciated! (also I apologize for not uploading my dataset / an anonymized dataset here; it’s not sharable at all).

Best,
Canaan

Edit: tagging @Solomon and @martinmodrak here since it’s basically a new glitch arising from a previous non-AR version of this same model that you guys had very helpful comments about; this one: Errors with BLUPs with complex random effects in brms using posterior_epred)

2

This is an interesting one. Could potentially also be a bug in brms. Can you provide a minimial reproducible example, please?

1 Like
3

Hi Paul,

Sure; the data provided here is blinded, but should work.

Data is here:
anon_data.csv (72.9 KB)
; the “newdat” dataframe I’m trying to use is here:
anon_newdat.csv (64.4 KB)

Even though the “newdat” has experiment and trial order info in it, I’d really just like one BLUP for each subject’s effect of Condition, so I’d ideally not use it.

brms code is here:


bf_autoreg <- bf(DependentVariable ~ Condition+ar(time = TrialNumber, gr = Subject, p = 1, cov = T) +(1+Condition|Experiment/Subject))


b_autoreg <- brm(
  formula = bf_autoreg,
  data = anon_tidy_dat,
  family = lognormal(),
  prior = c(prior(normal(2, 1), class = Intercept),
            #prior(normal(0, 1), class = sd)
            prior(normal(0, 0.5), class = b)
            ),
  iter = 10000, 
  warmup = 1000, 
  cores = 4,
  chains = 4,
  seed = 14,
  sample_prior = T,
  save_all_pars = T,
  save_model = "./Models/b_autoreg_stan",
  file = "./Models/b_autoreg",
  control = list(adapt_delta=0.8, max_treedepth = 10)
  )


preds <- posterior_epred(b_autoreg, newdata = newdata, scale = "linear", re_formula = ~ (1+Condition|Experiment/Subject))

Is this enough information for you to work with?

Thanks a lot for your help,
Canaan

4

@paul.buerkner Sorry to bug you about this, just wasn’t sure if you needed any more data from me (if you’re already working on this / aware of this, I apologize for the nag!)

5

Yes, I am aware, but I was on vacation and didn’t have time to look at this.

6

Your example works for me. I suggest you try with the latest github version of brms.

1 Like
7

Hi Paul,

Thanks a lot for your help - the update fixed the strange error I’d been getting. However, there’s still the problem of the “points must be identical within groups” and also not being able to make predictions for the effect of condition on a by-subject level without the entire timecourse of data in there. Is there a way to do this in the current set-up that I’m not aware of, or is my best hope just to manually extract the samples from the fitted model? (Or am I fundamentally misunderstanding the situation…)

Thanks very much,
Canaan

8

I don’t know to be honest. for AR terms you likely need the full time series. if you know what you want you can of course also do it manually with the posterior samples as you suggested but I am not sure I can be of much help with this problem right now unfortunately.

1 Like
9

No worries at all - I appreciate it a lot. I’ll just go for the manual approach then.

10

Hi,
got the same message this morning under brms 2.15. Atfter updating brms to brms 2.16 I got:

Error: Neither ‘ar’ nor ‘ma’ were supplied. Please report a bug.

Some month ago it works. So after installing an older version of brms (brms 2.14.4) it works for me.

11

Can you please provide a minimial reproducible example?

12

I have sent an email to you.

13

thanks! can you check if the problem is already solved when using the current dev version 2.16.1 from github?

2 Likes
14

Yes, problem is solved when using dev version!

1 Like