Hello! I need help on how to get the results I’m interested in from a `brms`

model I fit.

Once again, my background in statistics is lackluster, so pardon me if I make some mistake while writing this post.

The data I’m working with is structured as follows:

```
> head(dat_teste)
sample phase TS gamma tau
1 A02010 2 65 69.940 3033
2 A02010 2 65 38.470 2861
3 A02010 2 65 20.480 2726
4 A02010 2 65 10.870 2611
```

Where `phase`

and `TS`

are numerical discrete grouping factors, `gamma`

and `tau`

are my (x,y) data.

Here’s the model I fit using `brms`

.

```
modelo_HB <- brms::bf(tau ~ a + b * (gamma^c),
a ~ 1 + (1|TS) + (1|TS:phase),
b ~ 1 + (1|TS) + (1|TS:phase),
c ~ 1 + (1|TS) + (1|TS:phase),
sigma ~ gamma + (1|TS),
nl = TRUE)
prior <- c(prior(uniform(10, 8000), class = 'b', nlpar = 'a', lb = 10, ub = 8000),
prior(uniform(1, 1000), class = 'b', nlpar = 'b', lb = 1, ub = 1000),
prior(uniform(10^-3, 1), class = 'b', nlpar = 'c', lb = 10^-3, ub = 1),
prior(student_t(3, 0, 2.5), class = 'Intercept', dpar = 'sigma'))
```

As you can see, I’m trying to fit the following nonlinear equation to my data:

\tau = a + b\gamma^c

With a, b and c as the model parameters. The fit went decently well.

What I’m specially interested in is defining the credible intervals for the fit, **for a given value of TS, and unknown value of phase**, and defining the boundaries of the credible interval region in the (\gamma, \tau) plane using curves of the \tau = a + b\gamma^c format.

Simply put, for a given value of `TS`

(a.k.a. for a given group) and a non specified value of `phase`

, I want to have two (a,b,c) parameter vectors (one for the bottom limit, one for the upper limit) which define the boundaries of the credible interval region using the equation \tau = a + b\gamma^c.

How could I do this? Any help is greatly appreciated.