# Loo with custom family and censored data

**URL:** <https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792>\
**Category:** brms\
**Created:** [June 9, 2020, 11:36am UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792 "2020-06-09T11:36:30Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![gkreil](https://avatars.discourse-cdn.com/v4/letter/g/76d3ee/32.png) [@gkreil](https://discourse.mc-stan.org/u/gkreil)\
**Post date:** [June 9, 2020, 11:36am UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792/1 "2020-06-09T11:36:30Z")

</div>

Hi all  
I’d like to ask if one needs a further adjustment of the code provided in the custom family vignette  
([https://cran.r-project.org/web/packages/brms/vignettes/brms\_customfamilies.html](https://cran.r-project.org/web/packages/brms/vignettes/brms_customfamilies.html))  
if one is to evaluate model fit with loo on a dataset with censored observations.

Specifically, do I need to change the _log\_like\_custom_ function as below to get a proper elpd with brms:loo()?

```stan
set.seed(343)
dt.sim <- rnbinom(10125, size = 1, prob = rbeta(n = 10125, shape1 = 3, 25))
cens <- sample(1:10125, size = 1700, replace = F)
dt.sim[cens] <- dt.sim[cens] - 5
dt.sim <- ifelse(dt.sim < 0, 0, dt.sim)
censoring <- rep(0, 10125)
censoring[cens] <- 1
test.data <- data.frame(d = dt.sim, censored = censoring)

discrete_weibull <- custom_family(
  "discrete_weibull", dpars = c("mu", "delta"),
  lb = c(0, 0), ub = c(NA, NA), links = c('log', 'log'), 
  type = "int")

stan_funs <- "
  real discrete_weibull_lpmf(int y, real mu, real delta) {
    return log(((1/exp(mu))^y^delta) - (1/exp(mu))^(y+1)^delta);
  }
  real discrete_weibull_lccdf(int y, real mu, real delta){
    return log((1/exp(mu))^(y+1)^delta);
  }
    real discrete_weibull_lcdf(int y, real mu, real delta){
    return log(1 - ((1/exp(mu))^(y+1)^delta));
  }
"

stanvars <- stanvar(scode = stan_funs, block = "functions")

fit.dwei.brms <- brm(form=bf(d|cens(censored)~1, mu ~ 1, delta ~ 1,
                             family=discrete_weibull),
                     stanvars = stanvars,
                     data= test.data, 
                     chains=2,
                     iter = 2000,
                     prior = c(prior('normal(0.1, 0.5)', class = Intercept),
                               prior('normal(1, 0.5)', class = Intercept, dpar = delta)),
                     cores=2,
                     control = list(adapt_delta =0.9))

log_lik_discrete_weibull <- function(i, prep) {
  mu <- prep$dpars$mu[,i]
  delta <- prep$dpars$delta[,i]
  y <- prep$data$Y[i]
  cens <- prep$data$cens[,i]
  if(cens == 0){
    log(((1/exp(mu))^y^delta) - (1/exp(mu))^(y+1)^delta)
  } else{
    log((1/exp(mu))^(y+1)^delta)
  }
}

loo(fit.dwei.brms)

```

The log\_lik\_discrete\_weibull function does not work because there is no prep$data$cens, but any comment if this is in principle the way to go is very much appreciated

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 9, 2020, 12:23pm UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792/2 "2020-06-09T12:23:33Z")

</div>

I suggest you take a look at [https://github.com/paul-buerkner/brms/blob/master/R/log\_lik.R](https://github.com/paul-buerkner/brms/blob/master/R/log_lik.R) to see how censoring is handled for built-in families.

---

<div class="post-metadata">

**Author:** ![gkreil](https://avatars.discourse-cdn.com/v4/letter/g/76d3ee/32.png) [@gkreil](https://discourse.mc-stan.org/u/gkreil)\
**Post date:** [June 9, 2020, 12:30pm UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792/3 "2020-06-09T12:30:50Z")

</div>

So before I go through 1000 lines of code, there is some additional work to do right?

---

<div class="post-metadata">

**Author:** ![paul.buerkner](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/paul.buerkner/32/3303_2.png) [@paul.buerkner](https://discourse.mc-stan.org/u/paul.buerkner)\
**Post date:** [June 9, 2020, 12:41pm UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792/4 "2020-06-09T12:41:36Z")

</div>

Checking out log\_lik\_censor should mostly suffice. It may well be that your implementation is correct already but I didn’t check your lc©df functions

---

<div class="post-metadata">

**Author:** ![gkreil](https://avatars.discourse-cdn.com/v4/letter/g/76d3ee/32.png) [@gkreil](https://discourse.mc-stan.org/u/gkreil)\
**Post date:** [June 9, 2020, 12:59pm UTC](https://discourse.mc-stan.org/t/loo-with-custom-family-and-censored-data/15792/5 "2020-06-09T12:59:45Z")

</div>

Okay, as far as I can follow the code it is enough to pass the log\_likelihood and _log\_lik_ retrieves the lccdf of the custom family via prep$family$env or something similar.