I am unfortunately very confused about the relationship between parameters and the user supplied formula in a non-linear model under a log normal distribution. I have read the Bürkner, 2017 paper “Bayesian Distributional Non-Linear Multilevel Modeling with the R Package brms” but I fear I may lack the statistical background to understand how the equations therein actually relate to a real world example. I would be extremely grateful if someone can explain in simpleton friendly terms, and/or point me towards some resources to help me get a better grounding. The Package brms is fantastic, and its applications would be extremely useful for my work, so I am trying to get my head around it. Basically, I just want to get a mathematical representation of the equation using the model output.

For example, say I have this model that follows (what I hope is the correct syntax) for a power relationship a*x^b:

```
b <- c(2.12, 0.35)
x <- rnorm(500)
y <- rlnorm(500, mean = b[1] * x^ b[2])
dat1 <- data.frame(x, y) %>% filter(x > 0)
prior1 <- prior(normal(0, 1), nlpar = "b1", lb = 0) +
prior(normal(0, 1), nlpar = "b2", lb = 0)
formula.fit.1.gaus <- brm(bf(y ~ b1 * x ^ b2,
b1 ~ 1,
b2 ~ 1,
nl = TRUE),
data = dat1,
prior = prior1,
iter = 5000,
warmup = 2000,
thin = 5,
cores = 4,
family = gaussian(),
file = "formula.fit.1.gaussian",
sample_prior = "yes",
save_all_pars = T,
seed = 1234)
plot(conditional_effects(formula.fit.1.gaus), points = TRUE)
```

Returning output

```
Family: gaussian
Links: mu = identity; sigma = identity
Formula: y ~ b1 * x^b2
b1 ~ 1
b2 ~ 1
Data: dat1 (Number of observations: 262)
Samples: 4 chains, each with iter = 5000; warmup = 2000; thin = 5;
total post-warmup samples = 2400
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b1_Intercept 6.53 0.82 4.97 8.09 1.00 2352 2220
b2_Intercept 1.47 0.26 0.88 1.94 1.00 2101 2046
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma 12.87 0.62 11.69 14.17 1.00 2377 2333
Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).
```

My understanding is the estimates from the model output cannot simply be put back into the user supplied formula, as each non-linear parameter is a placeholder for a linear equation? If so, how then can I obtain the values that can be used to summaries the model?

However, when I try to put the estimates into the power equation (b1*x^b2) and plot them ( `curve(6.53*x^1.47, from = 0, to = 3.315768)`

), at a glance it appears to be fairly similar to the model described using `plot(conditional_effects(formula.fit.1.gaus), points = TRUE)`

. Can someone explain why this is? Is it just coincidence, or can the estimates sometimes be used, or have I just completely misunderstood?

Following this, if I do the same model, but under a lognormal distribution:

```
prior1 <- prior(normal(0, 1), nlpar = "b1", lb = 0) +
prior(normal(0, 1), nlpar = "b2", lb = 0)
formula.fit.1.log <- brm(bf(y ~ b1 * x ^ b2,
b1 ~ 1,
b2 ~ 1,
nl = TRUE),
data = dat1,
prior = prior1,
iter = 5000,
warmup = 2000,
thin = 5,
cores = 4,
family = lognormal(),
#control = list(adapt_delta = 0.99,
# max_treedepth = 12),
file = "formula.fit.1.log",
sample_prior = "yes",
save_all_pars = T,
seed = 1234)
```

with results

```
Family: lognormal
Links: mu = identity; sigma = identity
Formula: y ~ b1 * x^b2
b1 ~ 1
b2 ~ 1
Data: dat1 (Number of observations: 262)
Samples: 4 chains, each with iter = 5000; warmup = 2000; thin = 5;
total post-warmup samples = 2400
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
b1_Intercept 2.11 0.07 1.97 2.24 1.00 2607 2353
b2_Intercept 0.44 0.05 0.34 0.54 1.00 2138 2292
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma 0.97 0.04 0.89 1.06 1.00 2394 2290
Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).
```

I understand I cannot simply back transform parameters using exp() - but how could I go about getting the parameters from this model to place back into the supplied formula as a mathematical approximation?

Again would appreciate any insight whatsoever!