Non-linear Regression

Dr. Peng Zhao (✉ peng.zhao@xjtlu.edu.cn)

Department of Health and Environmental Sciences
Xi’an Jiaotong-Liverpool University

1 Learning objectives

1. Understand what is non-linear regression
2. Understand how non-linear regression is performed
3. Evaluate non-linear regression models
4. Choose possible non-linear models

2 Principle

2.1 Definition

Non-linear regression:

A model for the relationship between the dependent variable(s) and the independent variable(s) is not linear but a curve

Two categories:

1. One that can be transformed into a linear model

2. One that cannot be transformed into a linear model

Examples:

1. $y=a+b{e}^x$: Set $x^{\prime }=e^x$, then $y=a+b{x}^{\prime }$.
2. $y=a+b_{1}x+b_2{x}^2+...+b_m{x}^m$: Set $x_i=x^i,i=1,2,...,m$, then $y=a+b_1{x}_1+b_2{x}_2+...+b_m{x}_m$
3. $y=a{e}^{bx}$: Set $log(y)=log(a)+bx$ and $y^{\prime }=log(y)$, then $y^{\prime }=log(a)+bx$.
4. $z=a+b{e}^{xy}$: cannot be transformed into linear.

Name	Equation
Asymptotic functions
Michaelis-Menten	$y=\frac{ax}{b+cx}$ or $y=\frac{ax}{b+x}$ or…
2-parameter asymptotic exponential	$y=a(1-{\mathrm{e}}^{-bx})$
3-parameter asymptotic exponential	$y=a-b{\mathrm{e}}^{-cx}$
S-shaped functions
2-parameter logistic	$y=\frac{{\mathrm{e}}^{a+bx}}{1+{\mathrm{e}}^{a+bx}}$
3-parameter logistic	$y=\frac{a}{1+b{\mathrm{e}}^{-cx}}$
4-parameter logistic	$y=a+\frac{b-a}{1+{\mathrm{e}}^{(c-x)/d}}$
Weibull	$y=a-b{\mathrm{e}}^{-\left(c{x}^d\right)}$
Gompertz	$y=a{\mathrm{e}}^{-b{\mathrm{e}}^{-cx}}$
Humped curves
Ricker curve	$y=ax{\mathrm{e}}^{-bx}$
First-order compartment	$y=k\exp (-\exp (a)x)-\exp (-\exp (b)x)$
Bell-shaped	$y=a\exp (-|bx{|}^2)$
Biexponential	$y=a{\mathrm{e}}^{bx}-c{\mathrm{e}}^{-dx}$

par(mfrow = c(2,2))
curve(x/(1+x), 0, 10)
curve(1 - exp(-x), -100, -90)
curve(exp(x)/(1+exp(x)), -5, 5)
curve(x * exp(x), 95, 100)

2.2 Model

$$y_i=f(x_i,\theta )+{\epsilon }_i,i=1,2,...,n$$

• $y_i$: The response variable.
• $x_i=(x_{i1},x_{i2},...,x_{ik}{)}^{\prime }$: The explanatory variable.
• $\theta =({\theta }_0,{\theta }_1,...{\theta }_p{)}^{\prime }$: Unknown parameter vector.
• ${\epsilon }_i$: Random error term.

Assumptions:

• $E({\epsilon }_i)=0,i=1,2,...,n$: The mean value of random error term is 0
• $cor({\epsilon }_i,{\epsilon }_j)=0,i,j=1,2,...,n$, where $i\ne j$: No correlation.
• $cov({\epsilon }_i,{\epsilon }_j)={\sigma }^2,i,j=1,2,...,n$, where $i=j$: Equal variance
• The explanatory variable is non-random variable
• $f(x_i,\theta )$ is assumed to be a continuous differentiable function

2.3 Methods

Non-linear least squares (NLS):

$$Q(\theta )=\sum _{i=1}^n(y_i-f(x_i,\theta ){)}^2$$ For the minimal $Q(\theta )$:

$$\frac{\partial Q}{\partial {\theta }_j}=-2\sum _{i=1}^n(y_i-f(x_i,\theta ))\frac{\partial f}{\partial {\theta }_j}=0$$

Partial coefficient of determination: Measure the fitting of non-linear regression model.

$$R^2=1-\frac{SSE}{SST}$$

3 Workflow

3.1 Data

Reaction rate ~ concentration

dtf <- read.table("data/mm.txt", header = TRUE)
library(ggplot2)
ggplot(dtf) + geom_point(aes(conc, rate))

3.2 Fit the Model

Michaelis-Menten model $y=\frac{ax}{b+x}$:

• Biochemistry: the relationship between the reaction rate and the substrate concentration
• Ecology: (Holling’s disc equation) the relationship between the feeding rate of predator and the prey density.

Features:

1. The curve passes through the origin.
2. Rising rate of the curve diminishes with the increasing of x.
3. Function has an asymptote $y=a$.
4. $y=a/2$ when $x=b$.

m1 <- nls(rate ~ a * conc / (b + conc), data = dtf, start = list(a = 200, b = 0.03))
m1
summary(m1)

$$r=\frac{212.7c}{0.06412+c}$$

3-parameter asymptotic exponential model $y=a-b{\mathrm{e}}^{-cx}$:

m2 <- nls(rate ~ a - b * exp(-c * conc), data = dtf, start = list(a = 200, b = 150, c = 5))
m2
summary(m2)

$$r=201-153e^{-6.38c}$$

xv <- seq(0, 1.2, .01)
y1 <- predict(m1, list(conc = xv))
y2 <- predict(m2, list(conc = xv))
dtf_predicted <- data.frame(xv, y1, y2)
ggplot(dtf) + 
  geom_point(aes(conc, rate)) +
  geom_line(aes(xv, y1), data = dtf_predicted, color = 'blue') + 
  geom_line(aes(xv, y2), data = dtf_predicted, color = 'red')

Partial coefficient of determination:

$$R^2=1-\frac{SSE}{SST}$$

calc_R2 <- function(model) {
  ms <- summary(model)
  sse <- as.vector((ms[[3]])^2 * ms$df[2])
  null <- lm(dtf$rate ~ 1)
  sst <- as.vector(unlist(summary.aov(null)[[1]][2]))
  1 - sse/sst
}

calc_R2(m1)

[1] 0.9612608

calc_R2(m2)

[1] 0.9672987

3.3 Results

The Michaelis-Menten model can explain 96.1% of the total variation in the reaction rate, while the 3-parameter asymptotic exponential model can explain 96.4%.

3.4 Automatic starting values

m3 <- nls(rate ~ SSmicmen(conc, a, b), data = dtf)
summary(m3)

m4 <- nls(rate ~ SSasympOff(conc, a, b, c), data = dtf)
summary(m4)

Function	Model
SSasymp()	asymptotic regression model
SSasympOff()	asymptotic regression model with an offset
SSasympOrig()	asymptotic regression model through the origin
SSbiexp()	biexponential model
SSfol()	first-order compartment model
SSfpl()	four-parameter logistic model
SSgompertz()	Gompertz growth model
SSlogis()	logistic model
SSmicmen()	Michaelis–Menten model
SSweibull()	Weibull growth curve model

4 Further readings

• Nonlinear Regression with R. Christian Ritz & Jens Carl Streibig.
• First steps with Non-Linear Regression in R. R bloggers.