 # Difference between add_fitted_draws() and add_predicted_draws() when applied to brms binomial (Bernoulli) models

Dear Experts,

I am a Bayesian newbie. I have a question about the difference between add_fitted_draws() and add_predicted_draws() when applied to brms binomial (Bernoulli) models. Specifically, using a cut-off of 50%, I compared predictions from add_fitted_draws() with those from add_predicted_draws. I expected the 2 sets of values to be similar - but they differed. Below, I provide a minimal example.

Thank you in advance for your kind guidance and advice.

``````rm(list=ls())
library(rms)
library(tidyverse)
library(brms)
library(tidybayes)
library(janitor)

pr_flrm = c(
prior(normal(0, 0.857), class = "b", coef = "vs1"),
prior(normal(0, 0.00857), class = "b", coef = "hp"))

flrm <- brm(am ~ hp + vs,
data=mtcars,
family='bernoulli',
prior=pr_flrm,
sample_prior = TRUE,  ## so that we can visualize the prior distribution
save_pars = save_pars(all = TRUE),
seed=1202)

## predictions from add_fitted_draws
mtcars %>%
dplyr::select(hp, vs) %>%
add_fitted_draws(flrm, seed = 1234, n=1000) %>%
clean_names %>%
mutate(prediction = (value>0.50)*1) %>%   ## using a cut-off of 0.50
group_by(row)  %>% ## each subject
summarise(across(c(prediction), ~ list(smean.sd(.x)))) %>%
unnest_wider(prediction)

# A tibble: 32 x 3
row  Mean    SD
<int> <dbl> <dbl>
1     1 0.342 0.475
2     2 0.342 0.475
3     3 0.479 0.500
4     4 0.382 0.486
5     5 0.073 0.260
6     6 0.414 0.493
7     7 0.042 0.201
8     8 0.631 0.483
9     9 0.468 0.499
10    10 0.328 0.470
# ... with 22 more rows

## predictions from add_predicted_draws
mtcars %>%
dplyr::select(hp, vs) %>%
add_predicted_draws(flrm, seed = 1234, n=1000) %>%
clean_names %>%
group_by(row)  %>% ## each subject
summarise(across(c(prediction), ~ list(smean.sd(.x)))) %>%   ## smean.sd() from Hmsic
unnest_wider(prediction)

# A tibble: 32 x 3
row  Mean    SD
<int> <dbl> <dbl>
1     1 0.461 0.499
2     2 0.45  0.498
3     3 0.491 0.500
4     4 0.467 0.499
5     5 0.322 0.467
6     6 0.495 0.500
7     7 0.295 0.456
8     8 0.528 0.499
9     9 0.481 0.500
10    10 0.416 0.493
# ... with 22 more rows

``````

.

1 Like

I guess I have answered my own question: If I were to use rbinom() of size 1 to simulate a Bernoulli outcome, the 2 sets of predictions are similar.

``````## use rbinom of size 1 to simulate a Bernoulli outcome
mtcars %>%
dplyr::select(hp, vs) %>%
add_fitted_draws(flrm, seed = 1234, n=1000) %>%
clean_names %>%
group_by(row) %>%
mutate(prediction = map(value, ~ rbinom(1, size = 1, prob = .x))) %>%
unnest(prediction) %>%
summarise(across(c(prediction), ~ list(smean.sd(.x)))) %>%
unnest_wider(prediction)

# A tibble: 32 x 3
row  Mean    SD
<int> <dbl> <dbl>
1     1 0.429 0.495
2     2 0.438 0.496
3     3 0.502 0.500
4     4 0.461 0.499
5     5 0.359 0.480
6     6 0.5   0.500
7     7 0.256 0.437
8     8 0.56  0.497
9     9 0.493 0.500
10    10 0.455 0.498
# ... with 22 more rows

``````
1 Like

Great that you were able to move forward. Yes, the main distinction is that `add_predicted_draws` starts with the “fitted” draws but then adds the measurement noise (as you did via the `rbinom`)

Best of luck in further adventures with Stan!

Thank you for taking time to read and respond to my question. Kudos to the Stan team for the software and support.