Error while comparing custom family brms models

I have two custom family (reparametrized) lognormal brms models that I would like to compare. However, while trying to add “loo” criterion I get an error as follows:

Error in checkForRemoteErrors(val) : 8 nodes produced errors; first error: object ‘log_lik_lognormal7’ not found

REPRODUCIBLE EXAMPLE

Data

set.seed(0)
pi <- 0
mu_log <- 2
sigma_log <- 0.99
N = 500
y = (1 - rbinom(N, 1, prob = pi)) * rlnorm(N, mu_log, sigma_log)
x1 = seq(1,10, 1)
x2 = seq(5,15, 1)
df = data.frame(y=y, x1 = x1, x2 = x2)

Custom family
I define a custom family, giving standard deviations in natural scale. Reference for re-parametrisation can be found under paragraph 3.3: ProbOnto_2.5_13012017.pdf

custom_family(name = "lognormal7", dpars = c("mu", "sigma"), 
              links = c("log", "log"), 
              lb = c(0, 0), type = "real") ->
  lognormal7

stan_funs <- "
  real lognormal7_lpdf(real y, real mu, real sigma) {
    return lognormal_lpdf(y | log(mu/sqrt(1+sigma^2/mu^2)), sqrt(log(1+sigma^2/mu^2)));
  }
  real lognormal7_rng(real mu, real sigma) {
    return lognormal_rng(log(mu/sqrt(1+sigma^2/mu^2)), sqrt(log(1+sigma^2/mu^2)));
  }
"

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

log_lik_lognormal7 <- function(i, prep) {
  mu <- prep$dpars$mu[, i]
  if(NCOL(prep$dpars$sigma)==1){sigma <- prep$dpars$sigma}else
  {sigma <- prep$dpars$sigma[, i]}   ## [, i] if sigma is modelled, without otherwise
  y <- prep$data$Y[i]
  lognormal7_lpdf(y, mu, sigma)
}


posterior_predict_lognormal7 <- function(i, prep, ...) {
  mu <- prep$dpars$mu[, i]
  if(NCOL(prep$dpars$sigma)==1){sigma <- prep$dpars$sigma}else
  {sigma <- prep$dpars$sigma[, i]}   ## [, i] if sigma is modelled, without otherwise
  lognormal7_rng(mu, sigma)
}

posterior_epred_lognormal7 <- function(prep) {
  mu <- prep$dpars$mu
  return(mu)
}

Models

m1 = brm(bf(y ~ x1, sigma ~ x1), data = df, family = lognormal7, stanvars = stanvars, inits = "0")
m2 = brm(bf(y ~ x2, sigma ~ x2), data = df, family = lognormal7, stanvars = stanvars, inits = "0")

loo() comparison

m1 = m1 %>% add_criterion("loo")
m2 = m2 %>% add_criterion("loo")
loo = loo_compare(m1, m2)
loo

Error traceback

Error in checkForRemoteErrors(val) : 8 nodes produced errors; first error: object 'log_lik_lognormal7' not found
40.
stop(count, " nodes produced errors; first error: ", firstmsg, domain = NA)
39.
checkForRemoteErrors(val)
38.
staticClusterApply(cl, fun, length(x), argfun)
37.
clusterApply(cl = cl, x = splitList(X, nchunks), fun = lapply, FUN = fun, ...)
36.
do.call(c, clusterApply(cl = cl, x = splitList(X, nchunks), fun = lapply, FUN = fun, ...), quote = TRUE)
35.
parallel::parLapply(cl = cl, X = X, fun = FUN, ...)
34.
plapply(seq_len(N), log_lik_fun, cores = cores, prep = object)
33.
log_lik.brmsprep(prep, combine = combine, cores = cores)
32.
log_lik(prep, combine = combine, cores = cores)
31.
log_lik.brmsfit(object = .x1, newdata = .x2, resp = .x3, pointwise = .x4, save_psis = .x5)
30.
.fun(object = .x1, newdata = .x2, resp = .x3, pointwise = .x4, save_psis = .x5)
29.
eval(expr, envir, ...)
28.
eval(expr, envir, ...)
27.
eval2(call, envir = args, enclos = parent.frame())
26.
do_call(log_lik, ll_args)
25.
prepare_loo_args(x, newdata = newdata, resp = resp, pointwise = pointwise, save_psis = save_psis, ...)
24.
.loo(x = .x1, newdata = .x2, resp = .x3, model_name = .x4, pointwise = .x5, k_threshold = .x6, save_psis = .x7, moment_match = .x8, reloo = .x9, moment_match_args = .x10, reloo_args = .x11)
23.
eval(expr, envir, ...)
22.
eval(expr, envir, ...)
21.
eval2(call, envir = args, enclos = parent.frame())
20.
do_call(paste0(".", criterion), args)
19.
.fun(criterion = .x1, pointwise = .x2, resp = .x3, k_threshold = .x4, save_psis = .x5, moment_match = .x6, reloo = .x7, moment_match_args = .x8, reloo_args = .x9, x = .x10, model_name = .x11, use_stored = .x12)
18.
eval(expr, envir, ...)
17.
eval(expr, envir, ...)
16.
eval2(call, envir = args, enclos = parent.frame())
15.
do_call(compute_loo, args)
14.
.fun(models = .x1, criterion = .x2, pointwise = .x3, compare = .x4, resp = .x5, k_threshold = .x6, save_psis = .x7, moment_match = .x8, reloo = .x9, moment_match_args = .x10, reloo_args = .x11)
13.
eval(expr, envir, ...)
12.
eval(expr, envir, ...)
11.
eval2(call, envir = args, enclos = parent.frame())
10.
do_call(compute_loolist, args)
9.
loo.brmsfit(.x1, model_names = .x2)
8.
loo(.x1, model_names = .x2)
7.
eval(expr, envir, ...)
6.
eval(expr, envir, ...)
5.
eval2(call, envir = args, enclos = parent.frame())
4.
do_call(fun, args)
3.
add_criterion.brmsfit(., "loo")
2.
add_criterion(., "loo")
1.
m1 %>% add_criterion("loo")
1 Like

The problem seems to be that loo is computed in multiple threads which do not see the custom log_lik function. I am on my phone now, so can’t easily check the docs, but passing cores = 1 somewhere should probably do the trick.

Best of luck with your model!

@paul.buerkner: maybe brms could default to single core for custom families? Or is this once again a Windows-only thing?

1 Like

I think the best way to avoid that is to store the post-processing functions directly in the applied custom family. See ?custom_family for details. That way, they should be available also when using multiple cores and a forking cluster. Alternatively, explicitely set cores = 1 in add criterion to avoid any parallel processing. That should also fix it.

I will see whether there is also a fix directly from the brms side, but I cannot promise anything right now.

2 Likes

Thanks! I tried but got a slightly different error.

m1 = m1 %>% add_criterion(“loo”, cores = 1)

(I also run the model m1 with cores = 1 argument)

How to solve this? Where should I add cores argument in custom_family() function?

Error in lognormal7_lpdf(y, mu, sigma) : could not find function "lognormal7_lpdf"
37.
log_lik_fun(i, prep)
36.
FUN(X[[i]], ...)
35.
lapply(X, FUN, ...)
34.
plapply(seq_len(N), log_lik_fun, cores = cores, prep = object)
33.
log_lik.brmsprep(prep, combine = combine, cores = cores)
32.
log_lik(prep, combine = combine, cores = cores)
31.
log_lik.brmsfit(object = .x1, newdata = .x2, resp = .x3, pointwise = .x4, save_psis = .x5)
30.
.fun(object = .x1, newdata = .x2, resp = .x3, pointwise = .x4, save_psis = .x5)
29.
eval(expr, envir, ...)
28.
eval(expr, envir, ...)
27.
eval2(call, envir = args, enclos = parent.frame())
26.
do_call(log_lik, ll_args)
25.
prepare_loo_args(x, newdata = newdata, resp = resp, pointwise = pointwise, save_psis = save_psis, ...)
24.
.loo(x = .x1, newdata = .x2, resp = .x3, model_name = .x4, pointwise = .x5, k_threshold = .x6, save_psis = .x7, moment_match = .x8, reloo = .x9, moment_match_args = .x10, reloo_args = .x11)
23.
eval(expr, envir, ...)
22.
eval(expr, envir, ...)
21.
eval2(call, envir = args, enclos = parent.frame())
20.
do_call(paste0(".", criterion), args)
19.
.fun(criterion = .x1, pointwise = .x2, resp = .x3, k_threshold = .x4, save_psis = .x5, moment_match = .x6, reloo = .x7, moment_match_args = .x8, reloo_args = .x9, x = .x10, model_name = .x11, use_stored = .x12)
18.
eval(expr, envir, ...)
17.
eval(expr, envir, ...)
16.
eval2(call, envir = args, enclos = parent.frame())
15.
do_call(compute_loo, args)
14.
.fun(models = .x1, criterion = .x2, pointwise = .x3, compare = .x4, resp = .x5, k_threshold = .x6, save_psis = .x7, moment_match = .x8, reloo = .x9, moment_match_args = .x10, reloo_args = .x11)
13.
eval(expr, envir, ...)
12.
eval(expr, envir, ...)
11.
eval2(call, envir = args, enclos = parent.frame())
10.
do_call(compute_loolist, args)
9.
loo.brmsfit(.x1, model_names = .x2)
8.
loo(.x1, model_names = .x2)
7.
eval(expr, envir, ...)
6.
eval(expr, envir, ...)
5.
eval2(call, envir = args, enclos = parent.frame())
4.
do_call(fun, args)
3.
add_criterion.brmsfit(., "loo")
2.
add_criterion(., "loo")
1.
m1 %>% add_criterion("loo")

Did you run

expose_functions(m1, vectorize = TRUE)

beforehand?

1 Like