Help with predictions with brms and including errors

General
1

Hi,

I am new to brms and trying to use it to fit linear and non-linear models to my data. I have been able to make some progress following the documentation and other sources, but currently stuck with a couple of queries :

I am using the brms version 2.20.4 on Mac os 14.2.1

  1. I have fit my data with the following model -
fit <- brm(y | se(y_err, sigma=TRUE) ~ 1 + me(x-11.5, x_err), 
           data = poly_data, 
           prior = prior,
           chains = 4, cores = 4, iter = 5000, warmup = 2000, thin = 1)
prior <- prior(normal(8,12), class = Intercept) +
            prior(normal(0, 2), class = b) +
            prior(cauchy(0, 1), class = sigma)

and I have obtained the fit parameters. But now i am unable to use the posterior_predict() to get predictions for new values of x. The posterior_predict(fit) works fine, but whenever I input the newdata parameter - for e.g. posterior_predict(fit, newdata=poly_data) - I get the following error :

Error in terms.formula(formula, ...) : 
  invalid model formula in ExtractVars

Please let me know what am i doing wrong.

To skip this issue, I was thinking of estimating the distribution of slope and intercept from the best fit values+uncertainties and generating the y values in the equation y=bx+c. But I am not sure how to incorporate sigma.

rnorm_multi(n=1000, mu=fixef(fit)[,1], sd=fixef(fit)[,2], r = as.vector(cov2cor(vcov(fit))) )
'rnorm_multi' is from the faux package
  1. Similar to the poly_data, I have another dataset to which I wish to fit a double power law kind of model using the non-linear method,
fit <- brm(bf(y | mi(abs(y_err)) ~ c - log10(10^(-a*(x-d) )+10^(-b*(x-d) )),
              a+b+c+d~1, nl = TRUE), 
              data = poly_data,
              prior = prior,
              chains = 4, cores = 4, iter = 5000, warmup = 2000, thin = 1)

prior = prior(normal(0,1), nlpar="a") +
             prior(normal(1,2), nlpar="b") +
             prior(normal(9,11), nlpar="c") +
             prior(normal(10,13), nlpar="d")

In this model, I could incorporate the errors in y using mi(). But I could not figure how to include the errors in x (me() gives error). Let me know if it is possible to incorporate x errors.

Any help regarding the above two queries is really appreciated.
Thanks.

poly_data.txt (94.2 KB)

EDIT (Aki): added ticks for code blocks

2

For this problem specifically, if I fit your model with me(x, x_err) instead of me(x-11.5, x_err) (i.e. without the -11.5) then posterior_predict(fit, newdata = poly_data) works. So a workaround could be to create a new column in your data frame with that manual offset of -11.5 applied to it and use that in the model formula instead.

(This could be a bug in brms::posterior_predict() so it might be worth filing an issue on the github).

2 Likes