Exponential regression with nls in R

Question

I'm trying to solve the following problem:

x <- c(0.11557577149788574,2.1552479877306925,2.5505873377321175,1.0995198836006757,3.710225290286669,2.386870541964232,0.11557577149788574,0.11557577149788574,2.1552479877306925,2.5505873377321175,1.0995198836006757,3.710225290286669,2.386870541964232,0.11557577149788574)
y <- c(16500,11500,11500,13630,7000,11995,13490,16500,11500,11500,13630,7000,11995,13490)

df <- data.frame(x, y)

m <- nls(y ~ I(a*exp(-b*x)+c), data=df, start=list(a=14000, b=1, c=100), control=nls.control(maxiter=10000, minFactor=1e-7, tol=1e-5, printEval=F, warnOnly=F))

But, even if I try to change the start values and the nls control no value is returned. What I'm doing wrong? I need more points to solve that problem?

Thank you!

Answer 1

You need better starting values.

First order the data frame in increasing x so that plotting will work out nicely.

If we set c to zero we can fit the simpler model to log(y) ~ A - b*x which is linear in all coefficients so it can be fit via lm and no starting values are needed.

Use the starting value for b given by that simpler model. Also note that a and c enter the full model linearly so we can use the plinear algorithm of nls which eliminates the need to give starting values to those coefficients:

o <- order(df$x)
df_o <- df[o, ] # order it by increasing x

fm0 <- lm(log(y) ~ I(-x), df_o) # simpler model to get better starting values

st <- list(b = coef(fm0)[[2]])
fm <- nls(y ~ cbind(a = exp(-b*x), c = 1), df_o, start = st, alg = "plinear")

plot(df_o, col = "red", pch = 20)
lines(fitted(fm) ~ x, df_o)

The result is the following where .lin.a is a and .lin.c is c:

> fm
Nonlinear regression model
  model: y ~ cbind(a = exp(-b * x), c = 1)
   data: df_o
         b     .lin.a     .lin.c 
   -0.4903 -1529.0253 16509.4421 
 residual sum-of-squares: 10555038

Number of iterations to convergence: 5 
Achieved convergence tolerance: 7.9e-07

Here we show input data as red points and draw lines through the fitted result:

Answer 2

There is not enough information in your dataset to estimate the coefficient of this non-linear model reliably. There is practically no curvature in your data. If you use Levenberg-Marquardt it gets stuck with parameter values that are more or less equal to the linear fit:

plot(y~x)

library(minpack.lm)
m <- nlsLM(y ~ a*exp(-b*x)+c, data=df, start=list(a=14000, b=1, c=100),  
           control = nls.lm.control(maxiter = 1000))
summary(m)
#Formula: y ~ a * exp(-b * x) + c
#
#Parameters:
#    Estimate Std. Error t value Pr(>|t|)
#a  2.308e+06  6.822e+08   0.003    0.997
#b  8.506e-04  2.518e-01   0.003    0.997
#c -2.292e+06  6.822e+08  -0.003    0.997
#
#Residual standard error: 1250 on 11 degrees of freedom

#Number of iterations till stop: 94 
#Achieved convergence tolerance: 1.49e-08
#Reason stopped: Number of calls to `fcn' has reached or exceeded `maxfev' == 400.

lines(x, predict(m), col="green", lwd=2)
abline(lm(y~x), col="red", lty=2, lwd=2)

If the model is based on science, the range of your measured x values is probably too small.

Answer 3

The problem seems to be in your data.

If you try plot( y ~ x ) then it does not look exponential: indeed it seems to be slightly faster than linear, especially if you see that there are in fact two points in your data with the x y values 3.710225 7000 so you are trying to fit a convex curve to data which suggests a concave curve. Two suggestions:

1. Remove those two points:

   df2 <- df[df$y!=7000,]

   m2 <- nls(y ~ I(a*exp(-b*x)+c), data=df2, start=list(a=14000, b=1, c=100), control=nls.control(maxiter=10000, minFactor=1e-7, tol=1e-5, printEval=F, warnOnly=F))

to give

Nonlinear regression model
  model: y ~ I(a * exp(-b * x) + c)
   data: df2
        a         b         c 
1.418e+04 1.202e-01 1.030e+03 
 residual sum-of-squares: 9781328

Number of iterations to convergence: 14 
Achieved convergence tolerance: 2.573e-06

2. Have your start point suggest a concave curve by reversing the signs of a and b

   m <- nls(y ~ I(a*exp(-b*x)+c), data=df, start=list(a=-14000, b=-1, c=100), control=nls.control(maxiter=10000, minFactor=1e-7, tol=1e-5, printEval=F, warnOnly=F))

to give

Nonlinear regression model
  model: y ~ I(a * exp(-b * x) + c)
   data: df
         a          b          c 
-1529.0204    -0.4903 16509.4360 
 residual sum-of-squares: 10555038

Number of iterations to convergence: 7 
Achieved convergence tolerance: 7.471e-07