Fit 'nls': singular gradient matrix at initial parameter estimates

Question

I'm new using 'nls' and I'm encountering problems finding the starting parameters. I've read several posts and tried various parameters and formula constructions but I keep getting errors.

This is a small example of what I'm doing and I'd very much appreciate if anyone could give me some tips!

# Data to which I want to fit a non-linear function
x <- c(0,  4, 13, 30, 63, 92)
y <- c(0.00000000, 0.00508822, 0.01103990, 0.02115466, 0.04036655, 0.05865331)
z <- 0.98

# STEPS:
# 1 pool, z fixed. This works.
fit <- nls(y ~ z * ((1 - exp(-k1*x))),
           start=list(k1=0))


# 2 pool model, z fixed
fit2 <- nls(y ~ z * (1 - exp(-k1*x)) + (1 - exp(-k2*x)),
            start=list(k1=0, k2=0)) # Error: singular gradient matrix at initial parameter estimates


# My goal: 2 pool model, z free
fit3 <- nls(y ~ z * (1 - exp(-k1*x)) + (1 - exp(-k2*x)),
            start=list(z=0.5, k1=0, k2=0)) 

Answer 1

It has been a while since you asked the question but maybe you are still interested in some comments:

At least your fit2 works fine when one varies the starting parameters (see code and plots below). I guess that fit3 is then just a "too complicated" model given these data which follow basically just a linear trend. That implies that two parameters are usually sufficient to describe the data reasonable well (see second plot).

So as a general hint: When you obtain

singular gradient matrix at initial parameter estimates

you can

1) vary the starting values/your initial parameter estimates

and/or

2) try to simplify your model by looking for redundant parameters which usually cause troubles.

I also highly recommend to always plot the data first together with your initial guesses (check also this question).

Here is a plot showing the outcome for your fit, fit2 and a third function defined by me which is given in the code below:

As you can see, there is almost no difference between your fit2 and the function which has a variable z and one additional exponential. Two parameters seem pretty much enough to describe the system reasonable well (also one is already quite good represented by the black line in the plot above). If you then want to fit a line through a certain data point, you can also check out this answer.

So how does it now look like when one uses a linear function with two free parameters and a function with variable z, one exponential term and a variable offset? That is shown in the following plot; again there is not much of a difference:

How do the residuals compare?

> fit
Nonlinear regression model
  model: y ~ zfix * ((1 - exp(-k1 * x)))
   data: parent.frame()
       k1 
0.0006775 
 residual sum-of-squares: 1.464e-05

> fit2
Nonlinear regression model
  model: y ~ zfix * (1 - exp(-k1 * x)) + (1 - exp(-k2 * x))
   data: parent.frame()
        k1         k2 
-0.0006767  0.0014014 
 residual sum-of-squares: 9.881e-06

> fit3
Nonlinear regression model
  model: y ~ Z * (1 - exp(-k1 * x))
   data: parent.frame()
       Z       k1 
0.196195 0.003806 
 residual sum-of-squares: 9.59e-06

> fit4
Nonlinear regression model
  model: y ~ a * x + b
   data: parent.frame()
        a         b 
0.0006176 0.0019234 
 residual sum-of-squares: 6.084e-06

> fit5
Nonlinear regression model
  model: y ~ z * (1 - exp(-k1 * x)) + k2
   data: parent.frame()
       z       k1       k2 
0.395106 0.001685 0.001519 
 residual sum-of-squares: 5.143e-06

As one could guess, the fit with only one free parameter gives the worst while the one with three free parameters gives the best result; however, there is not much of a difference (in my opinion).

Here is the code I used:

x <- c(0,  4, 13, 30, 63, 92)
y <- c(0.00000000, 0.00508822, 0.01103990, 0.02115466, 0.04036655, 0.05865331)
zfix <- 0.98

plot(x,y)

# STEPS:
# 1 pool, z fixed. This works.
fit <- nls(y ~ zfix * ((1 - exp(-k1*x))), start=list(k1=0))
xr = data.frame(x = seq(min(x),max(x),len=200))
lines(xr$x,predict(fit,newdata=xr))

# 2 pool model, z fixed
fit2 <- nls(y ~ zfix * (1 - exp(-k1*x)) + (1 - exp(-k2*x)), start=list(k1=0, k2=0.5))
lines(xr$x,predict(fit2,newdata=xr), col='red')

# 3 z variable
fit3 <- nls(y ~ Z * (1 - exp(-k1*x)), start=list(Z=zfix, k1=0.2))
lines(xr$x,predict(fit3,newdata=xr), col='blue')

legend('topleft',c('fixed z, single exp', 'fixed z, two exp', 'variable z, single exp'), 
       lty=c(1,1,1),
       lwd=c(2.5,2.5,2.5),
       col=c('black', 'red','blue'))

#dev.new()
plot(x,y)

# 4 fit linear function a*x + b
fit4 <- nls(y ~ a *x + b, start=list(a=1, b=0.))
lines(xr$x,predict(fit4,newdata=xr), col='blue')

fit5 <- nls(y ~ z * (1 - exp(-k1*x)) + k2,  start=list(z=zfix, k1=0.1, k2=0.5))
lines(xr$x,predict(fit5,newdata=xr), col='red')

legend('topleft',c('linear approach', 'variable z, single exp, offset'), 
       lty=c(1,1),
       lwd=c(2.5,2.5),
       col=c('blue', 'red'))