Model 15.8

[ASKurz]

In Section 15.3.1, McElreath fit a model using what he called the “weighted average” approach to handling missing values for categorical data (p. 517). Here's his model:

# load
library(rethinking)

# data
set.seed(9) 

N_houses <- 100L 
alpha <- 5
beta <- (-3)
k <- 0.5
r <- 0.2
cat <- rbern( N_houses , k )
notes <- rpois( N_houses , alpha + beta*cat )
R_C <- rbern( N_houses , r )
cat_obs <- cat
cat_obs[R_C==1] <- (-9L)

dat <- list(
  notes = notes,
  cat = cat_obs,
  RC = R_C,
  N = as.integer(N_houses) )

# fit the model
m15.8 <- ulam( alist(
  # singing bird model
  ## cat known present/absent: 
  notes|RC==0 ~ poisson( lambda ), 
  log(lambda) <- a + b*cat,
  ## cat NA:
  notes|RC==1 ~ custom( log_sum_exp(
    log(k) + poisson_lpmf( notes | exp(a + b) ),
    log(1-k) + poisson_lpmf( notes | exp(a) ) 
  ) ),

  # priors
  a ~ normal(0,1), 
  b ~ normal(0,0.5),

  # sneaking cat model 
  cat|RC==0 ~ bernoulli(k),
  k ~ beta(2,2)
), data=dat , chains=4 , cores=4 )

I'm not sure if this is even possible with brms. If you know the answer or, better yet, know how to fit the model, I'd love to see your code.