Plotting conditional effects for sub-models in a non-linear model

Modeling
1

Is there a way to use the conditional_effects() function to plot the effects of factors at all levels of a non-linear model? For example, in the following model

y ~ var_a + var_b + var_a:var_b
var_a ~ var_c + var_d
var_b ~ 1

I’d like to be able to plot the effect of var_c on var_a using the conditional_effects function in brms. While I can plot the effect of any factor (var_a, var_b, var_c or var_d) on y, I haven’t found a way to plot factor effects for any of the “sub” responses. Is there a way to do this using the conditional_effects() function? If not, is there another function that makes this possible, or does this require a more manual approach?

1 Like
2

Hi @jasonm ,
As far as I know, I don’t think there’s a way to plot non-linear parameters in conditional_effects. But it’s not so bad to generate what you want yourself. The fitted function in brms accepts a nlpar argument which you can you use to pull those estimates on the original dataset, a specific data grid, etc.
Here’s a simulation (mostly) following the example model you mentioned:

Summary 
## packages
library(tidyverse)
library(brms)

## simulate some data according to:
  # y ~ var_a + var_b
    # (dropping the var_a:var_b bit for now since we seem to be assuming just an intercept for b)
  # var_a ~ var_c + var_d
  # var_b ~ 1
set.seed(20250111)
dat <- expand_grid(
  var_c = seq(-2, 2, length.out = 30),
  var_d = c("d1", "d2", "d3")
) %>% 
  mutate(
    # adding a dash of measurement error for a and b
    var_a = 1.5 * var_c + as.integer(str_remove(var_d, "d")) + 
      rnorm(n = n(), sd = 0.25),
    var_b = 5 + rnorm(n = n(), sd = 0.15),
    # simulate a set of y for each combo
    y = map2(var_a, var_b, 
             ~rnorm(mean = .x + .y, sd = 1, n = 30)) 
  ) %>% 
  unnest(cols = y)

# at y level:
dat %>% 
  ggplot(aes(var_a, y))+
  geom_point(aes(color = var_b))
dat %>% 
  ggplot(aes(var_b, y))+
  geom_point(aes(color = var_a))
# at var_a level (also note fewer 'observations'):
dat %>% 
  ggplot(aes(var_c, var_a))+
  geom_point(aes(color = var_d))
# note: you wouldn't observe the parameters var_a and var_b in reality, 
# but just variables c and d

bf_1 <- bf(
  # note: have to strip down parameter names at y level lest you get:
    # Error: Parameter names should not contain dots or underscores.
  y ~ a + b, #  + a*b
  a ~ var_c + var_d,
  b ~ 1,
  family = gaussian(),
  nl = T
)
# glimpse default prior:
get_prior(bf_1, data = dat)
# for nonlinear (multilevel) models, you pretty much always have 
# to set priors in my limited experience
priors_1 <- c(
  # a
  prior(normal(0, 2), coef = 'Intercept', nlpar = 'a'),
  prior(normal(1, 2), coef = 'var_c', nlpar = 'a'), # somewhat informative
  prior(normal(0, 2), coef = 'var_dd2', nlpar = 'a'),
  prior(normal(0, 2), coef = 'var_dd3', nlpar = 'a'),
  # b
  prior(normal(1.5, 3), coef = 'Intercept', nlpar = 'b') # say we knew it should likely be positive
  # sigma -- keep default for now
)

# let 'er rip:
fit_1 <- brm(
  bf_1, prior = priors_1, 
  data = dat,
  chains = 4, cores = 4, iter = 2000, backend = 'cmdstanr',
  control = list(adapt_delta = 0.95, max_treedepth = 12)
)
summary(fit_1)
# vanilla conditional_effects will just show y as fxn of c and d
plot(conditional_effects(fit_1), ask = F)
# to look at nonlinear parameters instead:
?fitted.brmsfit # use nlpar argument

dat %>% 
  # get a grid of c and d to get predictions of a on
  distinct(var_c, var_d) %>% # this works here
    # alternative using data_grid:
  # modelr::data_grid(
  #   var_c = seq_range(var_c, n = 100),
  #   var_d = unique(var_d)
  # ) %>% 
  bind_cols(., # bind the `a` estimates onto the grid
            fitted(fit_1, nlpar = 'a', newdata = .)) %>% 
  ggplot(aes(var_c, Estimate, color = var_d))+
  geom_ribbon(aes(ymin = Q2.5, ymax = Q97.5, fill = var_d), alpha = 0.1)+
  geom_line()
1 Like
3

@zacho That’s exactly what I was looking for, and I really appreciate the code. Thanks so much for your help.

1 Like