Non-linear specification (variable is not a valid distributional or non-linear parameter)

I am trying to use non-linear syntax to estimate several parameters in a model with brms where I am treating the outcome as a Beta proportion and trying to infer several other parameters.

The general specification is:

\begin{align*} S_{r,t} & \sim Beta(\mu, \kappa) \\ \mu & = I \\ I_{p,d} & = 1 + B_{p,d} \\ B_{p,d} & = a * R^b + c \\ \end{align*}

S_{r,t}, C, and R are observed, so I want to estimate I, B, a, b, and c

I first tried

bform <- bf(
Srt ~ I,
I ~ 1 + (C * B )+ (1|r) + (1|t),
B ~ 1 + (a * R^b + c),
a ~ 1,
b ~ 1,
c ~ 1,
nl = TRUE
)


But, that gives me a power error:
Error in terms.formula(as.formula(x)) : invalid power in formula

I assumed this was do to how R is interpreting the ^

So I eliminated b and assumed a 2

bform <- bf(
seroprev ~ I,
I ~ 1 + (C * B + (1|r) + (1|t),
B ~ 1 + (a * R^2 + c),
a ~ 1,
c ~ 1,
nl = TRUE
)


This gives a different error:
Error: The parameter 'B' is not a valid distributional or non-linear parameter. Did you forget to set 'nl = TRUE'?

I have set nl = TRUE.

Is it possible to estimate multiple parameters this way? Or should I be “injecting” these parameters into the transformed parameters block somehow?

When you eliminate the b and assume a 2, the model is no longer non-linear (R is observed, so R^2 is a constant).

You may well be right that the original problem has to do with how R interprets ^ in a model formula (i.e. as a degree of interaction). However,^ works in the example A real-world non-linear model in the relevant brms vignette here:
https://cran.r-project.org/web/packages/brms/vignettes/brms_nonlinear.html

If you could post a reproducible example of the original failure that yields Error in terms.formula(as.formula(x)) : invalid power in formula (running the call to bf from your post doesn’t reproduce the error for me) then I might be able to help. But I’ll also go ahead and pass you up the food chain by paging @paul.buerkner.

Thank you very much!

Here’s is an (unprincipled) reprex:

# load packaged
library(dplyr)
library(brms)

# Some fake data
Data <- expand.grid(r = factor(c(1:5)),
t = factor(c(1:8))
) %>%
mutate(C = round(rnorm(nrow(.),500, 150)),
D = round(rnorm(nrow(.),10000, 1000)),
R = C/D,
Srt = rnorm(nrow(.),0.5, 0.15))

# Non-linear formula
bform <- bf(
Srt ~ I,
I ~ 1 + (C * B) + (1|r) + (1|t),
B ~ 1 + (a * R(...^b) + c),
a ~ 1,
b ~ 1,
c ~ 1,
nl = TRUE
)

# Get available priors
get_prior(bform,
data = Data,
family = Beta())


This gives the following error:
Error: The parameter 'B' is not a valid distributional or non-linear parameter. Did you forget to set 'nl = TRUE'?

I see 2 problems:

(a) you can’t mix linear and non-linear formula syntax, that is I ~ B + (1 | r) for example is invalid if B is supposed to be a non-linear parameter.

(b) non-linear formula except for the main one need to be wrapped in nlf() to work. that is for example nlf(B ~ b^2)

Thank you, @paul.buerkner !

Where is the correct place to enter grouping variables for the non-linear parameters? In this case I want to group I, C, and B by r and t and leave a, b, and c independent of r and t.

I tried:

bform <- bf(
Srt ~ I + (1|r) + (1|t),
nlf(I ~ C * B),
nlf(B ~ a * R^b + c),
nlf(a ~ 1),
nlf(b ~ 1),
nlf(c ~ 1),
nl = TRUE
)

get_prior(bform,
data = Data,
family = Beta())


But this gives me:
Error in vector("list", nrow(re)) : invalid 'length' argument

group level terms need to be entered in linear formulas only.

and you need to distinguish between a non-linear parameter (e.g. c) and a non-linear formula. for example you don’t need nlf() around c ~ 1 because the formula itself is linear.

Thank you very much, @paul.buerkner !

This seems to be approaching the correct syntax (at least no immediate errors):

bform <- bf(
Srt ~ I + (1|r) + (1|t),
nlf(I ~ C * B),
nlf(B ~ a * R^b + c),
a ~ 1,
b ~ 1,
c ~ 1,
nl = TRUE
)

get_prior(bform,
data = Data,
family = Beta())



Is it possible to set priors on group-level effects in non-linear models?

Make some fake data:

Data <- expand.grid(r = factor(c(1:5)),
t = factor(c(1:8))
) %>%
mutate(C = round(rnorm(nrow(.),500, 150)),
D = round(rnorm(nrow(.),10000, 1000)),
R = C/D,
Srt = rnorm(nrow(.),0.5, 0.15))



Model formula:

bform <- bf(
Srt ~ I,
I ~ (1|r) + (1|t),
nlf(I ~ C * B),
nlf(B ~ a * R^b + c),
a ~ 1,
b ~ 1,
c ~ 1,
nl = TRUE
)


Get available priors:

prior_list <- get_prior(bform,
data = Data,
family = Beta())


Set priors on nlpars:

nl_priors <- c(
set_prior("normal(0, 10)",
class = "b",
nlpar = "a"),
set_prior("normal(0, 10)",
class = "b",
nlpar = "b"),
set_prior("normal(0, 10)",
class = "b",
nlpar = "c")
)


Compile and run model:

NL_Mod <- brm(bform,
Data,
family = Beta(),
prior = nl_priors,
inits = "random",
iter = 2000,
warmup = 1000,
chains = 4,
cores = ncores,
backend = "cmdstan",
normalize = FALSE,
save_pars = save_pars(all = TRUE),

However, I only get estimates for population level-parameters (a, b, and c) and I is a primary estimate of interest. How might I specify the model to return the estimates of I?
You can use posterior_epred(..., nlpar = "I")