# Loo with brms custom family

**URL:** <https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791>\
**Category:** General\
**Tags:** loo\
**Created:** [June 9, 2020, 10:30am UTC](https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791 "2020-06-09T10:30:51Z")\
**Posts on this page:** 4\
**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, 10:30am UTC](https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791/1 "2020-06-09T10:30:51Z")

</div>

Hi all  
There’s something wrong with my R code, I’m trying to fit a custom family in brms to  
a right censored dataset and do loo(). When I try to compute the log likelihood, I run out of memory. It says

_Error: vector memory exhausted (limit reached?)_

Below the reproducible example:

```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
  y <- prep$data$Y[i]
  log(((1/exp(mu))^y^delta) - (1/exp(mu))^(y+1)^delta)
}

loo(fit.dwei.brms)

```

My session info:  
R version 4.0.0 (2020-04-24)  
Platform: x86\_64-apple-darwin17.0 (64-bit)  
Running under: macOS Catalina 10.15.4

Matrix products: default  
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib  
LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib

locale:  
[1] en\_US.UTF-8/en\_US.UTF-8/en\_US.UTF-8/C/en\_US.UTF-8/en\_US.UTF-8

attached base packages:  
[1] stats graphics grDevices utils datasets methods base

other attached packages:  
[1] brms\_2.13.0 Rcpp\_1.0.4.6

---

<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, 10:41am UTC](https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791/2 "2020-06-09T10:41:07Z")

</div>

As always, things became clearer when posting the issue. So because there is censored data,  
I might need some lccdf function in the log\_lik… However, also without using censoring, the error persists.

---

<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, 11:11am UTC](https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791/3 "2020-06-09T11:11:57Z")

</div>

The line

```
delta <- prep$dpars$delta

```

is the problem as you predict delta as well. Please replace by

```
delta <- prep$dpars$delta[, i]

```

---

<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:12am UTC](https://discourse.mc-stan.org/t/loo-with-brms-custom-family/15791/4 "2020-06-09T11:12:26Z")

</div>

Thanks, I saw it as well just now.