Why does brms remove AR-effects in hierarchical structure?

Interfaces brms
,
1

It’s great that we can add autoregressive terms in brms with ar(), and it generally works well, even with a grouping. However, I’m confused about how this generates stancode when we have a hierarchical structure.

If we want the autoregressive parameters to have a hierarchical structure based on some grouping variable, it seems like the autoregressive part is removed. Here is an example:

  # Minimal reproducible example
  library(tidyverse)
  library(brms)

  # Make simple data.
  ed <- tsibble::filter_index(tsibble::as_tsibble(EuStockMarkets), ~ "1991-07")

  # Run both models.
  gr_ar1 <- make_stancode(value ~ ar(time = index, gr = key, p = 1), data = ed)
  hi_ar1 <- make_stancode(value ~ (1 + ar(time = index, p = 1) | key), data = ed)

  # Show different stancodes.
  gr_ar1
  # ...
  # // initialize linear predictor term
  # vector[N] mu = Intercept + rep_vector(0.0, N);
  # // include ARMA terms
  # for (n in 1:N) {
  #   err[n] = Y[n] - mu[n];
  #   for (i in 1:J_lag[n]) {
  #     Err[n + 1, i] = err[n + 1 - i];
  #   }
  #   mu[n] += Err[n, 1:Kar] * ar;
  # }
  # ...
  hi_ar1
  # ...
  # // initialize linear predictor term
  # vector[N] mu = Intercept + rep_vector(0.0, N);
  # for (n in 1:N) {
  #   // add more terms to the linear predictor
  #   mu[n] += r_1_1[J_1[n]] * Z_1_1[n];
  # }
  # ...

It seems clear that the AR effect (or ARMA, where I’ve seen the same) disappears if we give a hierarchical structure. However, conceptually I don’t see why we could not give the autoregressive parameters a second-level Normal distribution, so that they shrink towards each other, etc.

My question is, does this difference happen because:

  • brms does not know how to handle both at the same time, and prefers to keep the hierarchical part of the model,
  • The hierarchical grouping part of the model already has taken into account the AR-part of the model (which I can not find in the generated model for the life of me),
  • I am stretching what fits within the brms syntax, therefore should just write my own Stan code? (Should not be a problem, but I like the brms helper functions a lot.)

–

  • Operating System: Windows 10
  • brms Version: 2.14.4
2 
(1 + ar(time = index, p = 1) | key)

is invalid brms syntax and will not work. In fact it should return an error but for some reason it currently does not and just ignores the term.

3

Perhaps what you want is ar(time = index, gr = key, p = 1)?

4

OK, thanks Paul! I agree that it makes sense to be invalid syntax; I was indeed surprised that I could just run the model.

Yes, using the grouping function with ar(time = index, gr = key, p = 1) was one of the suggestions I give in the example, but I am not sure if it is what I want.

This is because I thought it would not put a hierarchical structure on the AR coefficients, as I can’t see it in the model code. It doesn’t right?

5

Indeed, it does not put hierarchical structure on the AR coefficients. Just groups the data so that they are not considered a single long time series.

6

OK, clear, thanks a lot! I think it would be better for me to write it directly in Stan, then.