zero_one_inflated_beta2.R

library(brms)
zero_one_inflated_beta2 <- custom_family(
  "zero_one_inflated_beta2", dpars = c("mu", "phi","zero","one"),
  links = c("logit", "log","logit","logit"), lb = c(NA, 0, 0, 0),
  ub = c(NA, NA, 1, 1), type = "real"
)
  


stan_funs <- "
     real zero_one_inflated_beta2_lpdf(real y, real mu, real phi,
                                      real zero, real one) {
       row_vector[2] shape = [mu * phi, (1 - mu) * phi]; 
       if (y == 0) { 
         return bernoulli_lpmf(1 | zero); 
       } else if (y == 1) {
         return bernoulli_lpmf(1 | one);
       } else { 
         return bernoulli_lpmf(0 | zero) + bernoulli_lpmf(0 | one) + beta_lpdf(y | shape[1], shape[2]);
       } 
      }
      
      real beta_binomial2_lpmf(int y, real mu, real phi, int T) {
        return beta_binomial_lpmf(y | T, mu * phi, (1 - mu) * phi);
      }
  
      real zero_one_inflated_beta2_rng(real mu, real phi, real zero, real one) {
  
        if (bernoulli_rng(zero) == 1) { 
         return 0; 
        } else if (bernoulli_rng(one) == 1) {
         return 1;
        } else { 
         return beta_rng(mu * phi, (1 - mu) * phi);
        }
  }
"

log_lik_zero_one_inflated_beta2 <- function(i, draws) {
  mu <- draws$dpars$mu[, i]
  phi <- draws$dpars$phi[, i]
  zero <- draws$dpars$zero[, i]
  one <- draws$dpars$one[, i]
  y <- draws$data$Y[i]
  zero_one_inflated_beta2_lpdf(y, mu, phi, zero, one)
}


predict_zero_one_inflated_beta2 <- function(i, draws, ...) {
  mu <- draws$dpars$mu[, i]
  phi <- draws$dpars$phi[, i]
  zero <- draws$dpars$zero[, i]
  one <- draws$dpars$one[, i]
  zero_one_inflated_beta2_rng(mu, phi, zero, one)
}

fitted_zero_one_inflated_beta2 <- function(draws) {
  with(draws$dpars, mu * (1 - zero - one) + one)
}