I recently organized a workshop on non-linear regression with R, covering the fundamentals and a few extras. When preparing the handout (which is available here as a PDF), I realized it is actually quite cumbersome to plot a fitted model with `nls`

, together with the data. To make this process easier, I wrote a simple package - `nlshelper`

- which does just that and a few other things. Read on for a simple demonstration, or visit the Github repository here, or read the full chapter in the link above. The package is also on CRAN.

The main motivation for writing this package was to quickly be able to plot a fitted curve on top of the data, but has since grown to include a few extra tools. The `nlstools`

package can already be used for plotting, with `plotfit`

, but it is a little inflexible, and it also does not accept models fit with `nlsList`

(more on that below).

The following example fits the Chapman-Richards equation to the `Loblolly`

dataset, and plot the fitted curve with the data. With the `extrapolate=TRUE`

argument, `plot_nls`

plots the fitted curve to the X-axis limits (as specified with `xlim`

).

```
# Define function to fit
chapm <- function(x,Asym,b,c)Asym*(1-exp(-b*x))^c
# Fit model
# Note the specification of the starting values.
nls_lob <- nls(height ~ chapm(age, Asym, b,c),
data=Loblolly,
start=list(Asym=100, b=0.1, c=2.5))
# Plot fitted curve
library(nlshelper)
plot_nls(nls_lob, pch=16, points.col="cornflowerblue",
xlim=c(0,40), ylim=c(0,80), extrapolate=TRUE)
```

It is straightforward to fit a non-linear regression by each group in a dataframe with `nlsList`

from `nlme`

, but it is much harder than it should be to plot the data with the fitted curves plotted on top. The `plot_nls`

function also accepts objects returned by `nlsList`

. The following example fits the Gompertz growth model to the `ChickWeight`

data, for each `Diet`

, and plots the curves with the data.

```
library(nlshelper)
library(nlme)
fit1 <- nlsList(weight ~ SSgompertz(Time, Asym, b2, b3)|Diet,
data=ChickWeight)
palette(c("#EAC435","#345995","#E40066","#03CEA4"))
plot_nls(fit1, pch=16, cex=0.5, lwd=2)
```

After having implemented the simple `plot_nls`

function, it turned out that it can be used (or at least was easily modified) for some other objects, including those fit with `loess`

and non-linear quantile regression (`nlrq`

from `quantreg`

).

The following example fits a non-linear quantile regression model to the `ChickWeight`

data, for three quantiles, and plots the fitted curves with the data. To plot multiple quantiles, it is convenient to wrap the fit in a function, accepting only `.tau`

as an argument.

```
library(quantreg)
fit_quan <- function(.tau){
nlrq(weight ~ SSgompertz(Time, Asym, b2, b3),
tau=.tau,
data=ChickWeight)
}
f <- fit_quan(0.1)
plot_nls(fit_quan(0.1), pch=16, cex=0.5, points.col="dimgrey", lines.col="firebrick", lty=3)
plot_nls(fit_quan(0.9), add=TRUE, lines.col="firebrick", lty=3)
plot_nls(fit_quan(0.5), add=TRUE, lines.col="black", lwd=2)
```

A question that frequently arises in non-linear regression is whether the fitted curve differs ‘somehow’ between groups. There are some rather cumbersome solutions to this problem, using indexing, but we’d like a simple function that gives an overall p-value for the contribution of the group to the fit. The `nlshelper`

package provides `anova_nlslist`

for this purpose. It simply performs an F-test of a model fit with `nlsList`

, which includes the grouping variable, and the equivalent model fit with `nls`

. It calculates the appropriate degrees of freedom, and prints a familiar `anova`

-like table.

The following example uses the built-in `Puromycin`

data. The data include reaction velocity () versus subtrate concentration in an enzymatic reaction for cells treated with the antibiotic Puromycin, or an untreated control. For enzymatic reactions that depend on the concentration of the substrate, the Michaelis-Menten model is often used, and follows from simple assumptions on the reaction rate versus the concentration of the substrate and enzyme.

The following example performs an F-test of a model including a grouping variable, versus one without.

```
# Fit the vanilla model without grouping
pur0 <- nls(rate ~ SSmicmen(conc, Vm, K), data=Puromycin)
# Fit a model that includes a grouping variable with nlsList
pur1 <- nlsList(rate ~ SSmicmen(conc, Vm, K)|state, data=Puromycin)
# F-test : does the fitted model differ significantly by 'state'?
# The first argument is the full model, the second one the reduced model.
anova_nlslist(pur1, pur0)
## Analysis of Variance Table
##
## Model 1: rate ~ SSmicmen(conc, Vm, K)
## Model 2: rate ~ SSmicmen(conc, Vm, K) | state
## Res.Df Res.Sum Sq Df Sum Sq F value Pr(>F)
## 1 21 7276.5
## 2 19 2055.1 2 5221.5 24.138 3.608e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

Clearly from this we learn that `state`

affects the relationship very (very) much. Of course, in this example, we don’t need p-values since it is obvious from a plot of the fitted model. Also useful in this context is a quick table of the fitted coefficients, separately for each group, with confidence intervals. The `nlshelper`

package adds a method to the `tidy`

function from the `broom`

package, as the following example shows.

```
# Only works if the nlshelper package is loaded
tidy(pur1, conf.int=TRUE)
## group term estimate std.error statistic p.value
## 1 treated Vm 212.68370711 6.947153189 30.614512 3.241160e-11
## 2 treated K 0.06412123 0.008280943 7.743227 1.565136e-05
## 3 untreated Vm 160.28011589 6.480250310 24.733630 1.384618e-09
## 4 untreated K 0.04770829 0.007781889 6.130683 1.727047e-04
## conf.low conf.high
## 1 197.30212762 229.29006416
## 2 0.04692517 0.08615995
## 3 145.63695367 176.54440095
## 4 0.03137460 0.07006131
# A plot of the fitted model
palette(c("#3BB8C4","#FF3D2F"))
plot_nls(pur1, pch=19, xlim=c(0,1.2))
# A future version of the package may include an automatic legend...
legend("bottomright", names(pur1), pch=19, col=palette())
```

I have a few developments planned for the `nlshelper`

package including:

- Methods for non-linear mixed-effects models
- Computation and plotting of confidence intervals for the fitted line, based on the bootstrap (much like the
`fitplc`

package does for very specific non-linear regressions, see our paper here)