# Non-linear model fitting problem

**URL:** <https://discourse.mc-stan.org/t/non-linear-model-fitting-problem/9038>\
**Category:** brms\
**Created:** [May 30, 2019, 1:58pm UTC](https://discourse.mc-stan.org/t/non-linear-model-fitting-problem/9038 "2019-05-30T13:58:28Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![ikron](https://avatars.discourse-cdn.com/v4/letter/i/ce73a5/32.png) [@ikron](https://discourse.mc-stan.org/u/ikron)\
**Post date:** [May 30, 2019, 1:58pm UTC](https://discourse.mc-stan.org/t/non-linear-model-fitting-problem/9038/1 "2019-05-30T13:58:28Z")

</div>

Please also provide the following information in addition to your question:

- Operating System: Red Hat Enterprise Linux
- brms Version: version 2.8.0 (package)

Hi,

I’m starting to use brms for non-linear modeling. I’m trying to fit a non-linear model, modeling the growth rate of a fungus in response to temperature. Example data from Ryan et al. 1943. Am. J. Bot. 30: 784-799  
`ryandata <- data.frame(temperature = c(42.30, 41.0, 39.0, 37.0, 35.7, 34.80, 34.35, 33.50, 32.80, 31.30, 29.80, 29.10, 28.50, 27.20, 26.10, 25.30, 24.80, 23.40, 22.70, 21.30, 20.0, 19.70, 18.0, 16.50, 15.0, 14.50, 13.10, 11.60, 11.40, 10.50, 10.20, 9.20, 8.30, 7.50, 6.70, 6.10, 5.80, 5.50, 5.10, 4.50, 3.80), growthrate = c(1.30, 2.50, 4.0, 4.90, 5.20, 5.20, 5.20, 5.0, 5.20, 4.60, 4.40, 4.30, 4.30, 3.80, 3.70, 3.70, 3.50, 3.30, 3.10, 2.90, 2.70, 2.70, 2.30, 2.10, 1.90, 1.70, 1.60, 1.30, 1.20, 1.10, 1.0, 0.86, 0.72, 0.54, 0.37, 0.31, 0.36, 0.30, 0.18, 0.07, 0.06))`

I’m fitting a lactin model to the data, which is motivated by some theory and biological considerations:  
$$ y = \exp{\rho T} - \exp{\rho T\_{max} - \frac{T\_{max} - T}{\delta}} + \lambda $$

I can use standard nls in R to fit this model:

```
start.val <- c(rho = 0.05, Tmax = 43, delta = 8, lambda = -1)
model.lactin <- nls(growthrate ~ exp(rho*temperature) - exp(rho*Tmax - (Tmax-temperature)/delta) + lambda, start = start.val, trace = T, data = ryandata)

summary(model.lactin)

Formula: growthrate ~ exp(rho * temperature) - exp(rho * Tmax - (Tmax - 
    temperature)/delta) + lambda

Parameters:
        Estimate Std. Error t value Pr(>|t|)    
rho 0.073298 0.006612 11.086 2.57e-13 ***
Tmax 43.789369 0.148651 294.578 < 2e-16 ***
delta 6.928934 0.858967 8.067 1.13e-09 ***
lambda -0.980015 0.044927 -21.813 < 2e-16 ***
---
Signif. codes: 0 ‘ ***’ 0.001 ‘** ’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.1657 on 37 degrees of freedom

Number of iterations to convergence: 18 
Achieved convergence tolerance: 7.335e-06

```

Plotting reveals that while the model is not ideal, it is not completely unreasonable:

```
plot(ryandata$temperature, ryandata$growthrate)
lines(ryandata$temperature, predict(model.lactin))

```

I then try to use brms to fit the model:

```
fit_lactin <- brm(
    bf(growthrate ~ exp(rho*temperature) - exp(rho*Tmax - (Tmax-temperature)/delta) + lambda,
       rho ~ 1, Tmax ~ 1 , delta ~ 1 , lambda ~ 1 , nl = TRUE),
    data = ryandata, family = gaussian(),
    prior = c(
        prior(normal(0.1,10), nlpar = "rho"),
        prior(normal(42,10), nlpar = "Tmax", lb = 1, ub = 50),
        prior(normal(7,10), nlpar = "delta"),
        prior(normal(-1,10), nlpar = "lambda")
        ),
    warmup = 5000, iter = 10000, chains = 2, cores = 2, control = list(adapt_delta = 0.95)
    )

```

However, so far I’ve not managed to get the model to convergence. The chains seem to get stuck in different values. I’ve read brms non-linear vignette that states that non-linear model are difficult to fit in general. I’ve tried different priors but no luck so far. So my questions are:

1. Are there some guidelines on how to approach the choice of priors in this case?
2. I guess this is an identifiability issue? So, should I just try to reparameterize the model in a different way?
3. Any advice for reparameterization, on how to go about doing this?

So, if anyone has any suggestions or could point me towards useful resources that would be great.

---

<div class="post-metadata">

**Author:** ![wpetry](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/wpetry/32/4852_2.png) [@wpetry](https://discourse.mc-stan.org/u/wpetry)\
**Post date:** [June 3, 2019, 9:36pm UTC](https://discourse.mc-stan.org/t/non-linear-model-fitting-problem/9038/2 "2019-06-03T21:36:02Z")

</div>

I’m hardly expert on this issue given that I’m stuck with my own non-linear model, but I see two improvements that will get yours to fit.

**1) Strengthen the prior on “lambda”**  
\mathcal{N}(-1,10) is pretty weak, and as I understand it non-linear models need extra attention to the prior specification in order to converge. _How_ you go about this might be a bit more challenging. It’s clear how to proceed if, say, `lambda = 10` is biologically unreasonable in your study system. But if those extreme values _are_ plausible, others would have to weigh in.

**2) Rescale model parameters**  
From [another post](https://discourse.mc-stan.org/t/brms-multilevel-nonlinear-model-model-fails-to-converge-when-converted-from-individual-level-to-multilevel/), @paul.buerkner recommends:

> [@Brms multilevel nonlinear model - model fails to converge when converted from individual level to multilevel](https://discourse.mc-stan.org/t/brms-multilevel-nonlinear-model-model-fails-to-converge-when-converted-from-individual-level-to-multilevel/4691/3):
>
> As another recommendation: Try to make sure that your non-linear parameters are roughly of the same order of magnitute. If one is, say, between 0.001 and 0.01 and another is between 100 and 1000, Stan may have trouble sampling from the model.

Putting the two together, I come to:

```
fit_lactin <- brm(
  bf(growthrate ~ exp(rho*temperature) - exp(rho*TmaxSC*10 - (TmaxSC*10-temperature)/delta) + lambdaSC/10,
     rho ~ 1,
     TmaxSC ~ 1,
     delta ~ 1,
     lambdaSC ~ 1,
     nl = TRUE),
  data = ryandata, family = gaussian(),
  prior = c(
    prior(normal(0.1,10), nlpar = "rho"),
    prior(normal(4.2,1.0), nlpar = "TmaxSC", lb = 0.1, ub = 5.0),
    prior(normal(7,10), nlpar = "delta"),
    prior(normal(-10,2), nlpar = "lambdaSC")
  ),
  warmup = 5000, iter = 10000, chains = 2, cores = 2,
  control = list(adapt_delta = 0.95)
)

```

You can see that I’ve just multiplied `Tmax` and divided `lambda`—both in the model formula _and_ in the priors—to bring them onto a similar order of magnitude as the other parameters. I append “SC” as a suffix to remind myself that the posteriors for these parameters need to be back-transformed. Although it’s a bit harder to see because of this scaling, I’ve also dramatically narrowed the prior on `lambda`. This fits reasonably quickly and gives the following:

```
Population-Level Effects: 
                   Estimate Est.Error l-95% CI
rho_Intercept 0.08 0.01 0.07
TmaxSC_Intercept 4.38 0.01 4.35
delta_Intercept 7.07 0.56 5.94
lambdaSC_Intercept -9.52 0.62 -10.52
                   u-95% CI Eff.Sample Rhat
rho_Intercept 0.09 1902 1.00
TmaxSC_Intercept 4.40 3831 1.00
delta_Intercept 8.03 2010 1.00
lambdaSC_Intercept -8.01 2593 1.00

Family Specific Parameters: 
      Estimate Est.Error l-95% CI u-95% CI
sigma 0.17 0.02 0.14 0.22
      Eff.Sample Rhat
sigma 4260 1.00

```

These correspond (after backtransforming the “SC” parameters) quite nicely to the estimates from `nls`:

```
Parameters:
        Estimate Std. Error t value Pr(>|t|)    
rho 0.073298 0.006612 11.086 2.57e-13 ***
Tmax 43.789369 0.148651 294.579 < 2e-16 ***
delta 6.928933 0.858966 8.067 1.13e-09 ***
lambda -0.980015 0.044927 -21.813 < 2e-16 ***

```

I hope this helps in your case, but so far these “tricks” haven’t gotten me out of my issue. I would also be keen to learn about the strategies for re-parameterization. I also haven’t been able to find specific, practical suggestions on this front. Perhaps these strategies could help you make use of weaker priors that you may have to use when fitting other data to this model.

---

<div class="post-metadata">

**Author:** ![ikron](https://avatars.discourse-cdn.com/v4/letter/i/ce73a5/32.png) [@ikron](https://discourse.mc-stan.org/u/ikron)\
**Post date:** [June 4, 2019, 7:00am UTC](https://discourse.mc-stan.org/t/non-linear-model-fitting-problem/9038/3 "2019-06-04T07:00:30Z")

</div>

Hi, thanks a lot for the advice. Rescaling seems to solve the problem in this case. Good to know.
