Plotting predictions based on the linear combination of predictors

Interfaces brms
, ,
1

I am modelling a forced-choice task with two alternatives where each alternative is characterized by two parameters that are supposed to reflect some underlying signals. Since these parameters have the same effect on each alternative, I want to plot the model predictions using the average of the coefficients for each alternative. How to do it?

Below I show how I compute the hypothesis but I struggle to plot the predictions for the same hypothesis.

library(data.table)
library(brms)

# example data
object_set <- data.table(object_id = 1:800, mean_val = rnorm(800), unc = rgamma(800, 1, 2))
pairs <- object_set[,  CJ(id1 = object_id, id2 = object_id)][id1 < id2]
pairs <- pairs[object_set,on = .(id1=object_id)][object_set, on = .(id2 = object_id)][!is.na(mean_val)]
setnames(pairs, 3:6, c('mean1','unc1','mean2','unc2'))
sample_data <- pairs[sample(1000)]

sample_data[,evidence:=rnorm(.N, mean1, unc1)-rnorm(.N, mean2, unc2), by  = .(id1, id2)]
sample_data[,choice:=as.numeric(evidence>0)]

brm_simple <- brm(choice~mean1+mean2+mean1*unc1+mean2*unc2, sample_data, family = bernoulli())

hypothesis(brm_simple, '(mean1-mean2)/2>0')
hypothesis(brm_simple, '(mean1:unc1-mean2:unc2)/2<0')

So I want to plot, for example, the choice probability as a function of the mean for or the interaction of the mean and uncertainty regardless of alternative number. This would be similar to an average of the plots given by conditional effects here:

conditional_effects(brm_simple, 'mean1:unc1')
conditional_effects(brm_simple, 'mean2:unc2')

with the x-axis reversed for the second plot.

2

For more flexibility in plotting see Extracting and visualizing tidy draws from brms models • tidybayes

1 Like
3

Thanks for the tip! It’s convenient to manipulate the draws with tidybayes but I still struggle to do what I want. I thought it could do something like this:

preds <- sample_data %>%
  data_grid(mean1 = seq_range(c(-3,3), n = 50),
            mean2 = seq_range(c(-3,3), n = 50),
            unc1 = c(0.5),
            unc2 = c(0.5)) %>%
  add_epred_draws(brm_simple, ndraws = 1000) %>%  
  mutate(mean_avg = (mean1-mean2)/2)

preds %>% 
  pivot_longer(starts_with('mean')) %>%
  ggplot(aes(x = value, y = .epred )) + facet_grid(~name) + 
  stat_lineribbon()

But then the results I get for the average effect of the mean does not look right:

image
image1158×526 59.5 KB

I guess that’s because there are different predictions for symmetric values of mean1 and mean2 (e.g., mean1 = -1 and mean2 = 1 would be different from mean1 = 1 and mean2 = -1 though the average is the same). How can I get the smooth prediction curve?

I know I can draw from the predictive posterior and then use glm to get predictions but I doesn’t sound right either.