Non-linear model with integrate_1d in Bmrs

Interfaces brms
1

I am trying to define a non-linear model in Brms using “integrate_1d” but I am having problems dealing with the type of parameter arguments passed to integrate_1d.

Here is a simplified version of the code:

wtpdata <- read.table(".../wtpdata_test.csv", header=TRUE, sep=",")

stan_funs <- "
  real integrand(real x, real xc, real[] theta, real[] x_r, int[] x_i) {
	  real a = theta[1];
      real p = x_r[1];                	
    return exp(-(-log(x*p))^a);
  }
  real integralfun(data real p, real a, data int xi) {
  	  	return integrate_1d(integrand, 0.0, 1.0, {a}, {p}, {xi}, 0.01);
  }
"

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

nlform <- bf(wtp ~ integralfun(prob, a, id), a ~ (1|id), nl = TRUE)

nlprior <- c(prior(lognormal(2, 1), lb=0,  nlpar = "a"))
        
fit2 <- brm(formula = nlform, data = wtpdata, family = zero_one_inflated_beta,
                 prior = nlprior, control = list(adapt_delta = 0.9), stanvars = stanvars)

Here is a link to the test dataset: wtpdata_test

Here is the error message I get:

Error in compileCode(f, code, language = language, verbose = verbose) : ^/Library/Frameworks/R.framework/Versions/4.1/Resources/library/StanHeaders/include/stan/math/prim/arr/functor/integrate_1d.hpp:195:15: note: candidate function template not viable: no known conversion from ‘vector<local_scalar_t__>’ to ‘const vector’ for 4th argument in line double integrate_1d( ^22 warnings and 1 error generated.make: *** [file5c9845fad7bc.o] Error 1
Error in sink(type = “output”) : invalid connection

  • Operating System: MacOS 11.2.3
  • brms Version: 2.15.0
2

Thanks for changing to the forum!

I just tried your example and it does fail for me as well when using rstan as you do. The error you see looked to me as if this should not happen such that I switched to use as backend for brms cmdstan which allows me to use the current 2.27.0. Doing so fixes the problem and your example is sampling on my screen just fine (slow though).

So can you please switch to cmdstan? To do so use as argument backend=“cmdstanr” for the brm call and install cmdstanr. See Getting started with CmdStanR • cmdstanr

As a bonus you will be able to benefit from within-chain parallelisation and use multiple cores per chain (lookup the threading argument of brm and there is also a vignette for this available on CRAN).

2 Likes
3

Thanks a lot for the suggestion and for the tip about using within-chain parallelization! I changed the above code and it works fine as long as I do not use within-chain parallelization. However, if I try threading I get a compilation error.

Specifically, if I run this code the sampling goes fine:

fit2 <- brm(formula = nlform, data = wtpdata, family = zero_one_inflated_beta,
                 prior = nlprior, control = list(adapt_delta = 0.9), stanvars = stanvars, chains = 4, cores = 4, backend = "cmdstanr")

However, if I run this code

it2 <- brm(formula = nlform, data = wtpdata, family = zero_one_inflated_beta,
                 prior = nlprior, control = list(adapt_delta = 0.9), stanvars = stanvars, chains = 2, cores = 2, threads = threading(2), backend = "cmdstanr")

I get the following error:

Ill-typed arguments supplied to function ‘integralfun’. Available signatures:
(data real, real, data int) => real
Instead supplied arguments of incompatible type: real, real, int.make

4

Oh… you are hitting an issue in brms here, so tagging @paul.buerkner .

brms generates this partial sum function:

  // compute partial sums of the log-likelihood
  real partial_log_lik_lpmf(int[] seq, int start, int end, vector Y, matrix X_a, vector b_a, vector C_1, int[] C_2, real phi, real zoi, real coi, int[] J_1, vector Z_1_a_1, vector r_1_a_1) {
    real ptarget = 0;
    int N = end - start + 1;
    // initialize linear predictor term
    vector[N] nlp_a = X_a[start:end] * b_a;
    // initialize non-linear predictor term
    vector[N] mu;
    for (n in 1:N) {
      // add more terms to the linear predictor
      int nn = n + start - 1;
      nlp_a[n] += r_1_a_1[J_1[nn]] * Z_1_a_1[nn];
    }
    for (n in 1:N) {
      int nn = n + start - 1;
      // compute non-linear predictor values
      mu[n] = inv_logit(integralfun(C_1[nn] , nlp_a[n] , C_2[nn]));
    }
    for (n in 1:N) {
      int nn = n + start - 1;
      ptarget += zero_one_inflated_beta_lpdf(Y[nn] | mu[n], phi, zoi, coi);
    }
    return ptarget;
  }

but because you are using the covariate as an argument to the data x_r argument of integrate_1d the newer stanc3 compiler enforce that you have to declare C_1 as argument to the partial sum function also with a “data” qualifier. So this compiles:

  // compute partial sums of the log-likelihood
  real partial_log_lik_lpmf(int[] seq, int start, int end, vector Y, matrix X_a, vector b_a, data vector C_1, int[] C_2, real phi, real zoi, real coi, int[] J_1, vector Z_1_a_1, vector r_1_a_1) {
    real ptarget = 0;
    int N = end - start + 1;
    // initialize linear predictor term
    vector[N] nlp_a = X_a[start:end] * b_a;
    // initialize non-linear predictor term
    vector[N] mu;
    for (n in 1:N) {
      // add more terms to the linear predictor
      int nn = n + start - 1;
      nlp_a[n] += r_1_a_1[J_1[nn]] * Z_1_a_1[nn];
    }
    for (n in 1:N) {
      int nn = n + start - 1;
      // compute non-linear predictor values
      mu[n] = inv_logit(integralfun(C_1[nn] , nlp_a[n] , C_2[nn]));
    }
    for (n in 1:N) {
      int nn = n + start - 1;
      ptarget += zero_one_inflated_beta_lpdf(Y[nn] | mu[n], phi, zoi, coi);
    }
    return ptarget;
  }

note that only “vector C_1” is replaced with “data vector C_1”… the real issue is that this difficult to make consistent with old and new stan (aka stanc3 based vs old compiler).

Can you file an issue for this now? This is a legit problem to sort out.

1 Like
5

Thanks for finding the problem with the declaration of data arguments. I’ll file the issue.