# Singular gradient error in nls with correct starting values

I'm trying to fit a line+exponential curve to some data. As a start, I tried to do this on some artificial data. The function is: $$y=a+b\cdot r^{(x-m)}+c\cdot x$$ It is effectively an exponential curve with a linear section, as well as an additional horizontal shift parameter (m). However, when I use R's nls() function I get the dreaded "singular gradient matrix at initial parameter estimates" error, even if I use the same parameters that I used to generate the data in the first place.
I've tried the different algorithms, different starting values and tried to use optim to minimise the residual sum of squares, all to no avail. I've read that a possible reason for this could be an over-parametrisation of the formula, but I don't think it is (is it?)
Does anyone have a suggestion for this problem? Or is this just an awkward model?

A short example:

#parameters used to generate the data
reala=-3
realb=5
realc=0.5
realr=0.7
realm=1
x=1:11 #x values - I have 11 timepoint data
#linear+exponential function
y=reala + realb*realr^(x-realm) + realc*x
#add a bit of noise to avoid zero-residual data
jitter_y = jitter(y,amount=0.2)
testdat=data.frame(x,jitter_y)

#try the regression with similar starting values to the the real parameters
linexp=nls(jitter_y~a+b*r^(x-m)+c*x, data=testdat, start=list(a=-3, b=5, c=0.5, r=0.7, m=1), trace=T)


Thanks!

• Hint: look at the coefficient of $r^x$ (for a fixed $r$) and note that $b r^{-m} = \text{constant}$ has a one-dimensional family of solutions $(b,m)$ with $b = r^m \cdot \text{constant}$. – whuber Jul 14 '11 at 17:35
• This isn't an identified model, unless $b$ or $r$ are somehow constrained. I think requiring the $r \in (0,1)$ would do the job. – Macro Jul 14 '11 at 19:02

I've got bitten by this recently. My intentions were the same, make some artificial model and test it. The main reason is the one given by @whuber and @marco. Such model is not identified. To see that, remember that NLS minimizes the function:

$$\sum_{i=1}^n(y_i-a-br^{x_i-m}-cx_i)^2$$

Say it is minimized by the set of parameters $$(a,b,m,r,c)$$. It is not hard to see that the the set of parameters $$(a,br^{-m},0,r,c)$$ will give the same value of the function to be minimized. Hence the model is not identified, i.e. there is no unique solution.

It is also not hard to see why the gradient is singular. Denote

$$f(a,b,r,m,c,x)=a+br^{x-m}+cx$$

Then

$$\frac{\partial f}{\partial b}=r^{x-m}$$

$$\frac{\partial f}{\partial m}=-b\ln rr^{x-m}$$

and we get that for all $$x$$

$$b\ln r\frac{\partial f}{\partial b}+\frac{\partial f}{\partial m}=0.$$

Hence the matrix

\begin{align} \begin{pmatrix} \nabla f(x_1)\\\\ \vdots\\\\ \nabla f(x_n) \end{pmatrix} \end{align}

will not be of full rank and this is why nls will give the the singular gradient message.

I've spent over a week looking for bugs in my code elsewhere till I noticed that the main bug was in the model :)

• This is ages old I know but just wondering, does this mean nls cannot be used on models that are not identifiable? For example, a neural network? – Count Zero Jul 20 '15 at 19:47
• fat chance, I know, but could you break this down for less calc-remembering folks? :). also, what is the solution to the OP's problem, then? Give up and go home? – theforestecologist May 4 '18 at 19:37
• The solution to OP problem is to use one parameter instead of two, i.e. instead of $b\cdot r^{x-m}$ use $\beta \cdot r^x$. The $m$ parameter is superflous, since it merges into $\beta$, i.e. $\beta = b\cdot r^{-m}$. – mpiktas May 8 '18 at 6:44
• @CountZero, basically yes, usual optimisation methods would fail if the parameters are unindentified. Neural networks sidestep this problem however, by adding additional contraints and using other interesting tricks. – mpiktas May 8 '18 at 6:47
• I would think $\frac {\partial f}{\partial m} = -b \ln{r}\ r^{x-m}$ ? there is a lnr missing in the answer ? – wiswit Feb 19 at 5:54

The answers above are, of course, correct. For what its worth, in addition to the explanations given, if you are trying this on an artificial data set, according to the nls help page found at: http://stat.ethz.ch/R-manual/R-patched/library/stats/html/nls.html

R's nls wont be able to handle it. The help page specifically states:

Warning

Do not use nls on artificial "zero-residual" data.

The nls function uses a relative-offset convergence criterion that compares the numerical imprecision at the current parameter estimates to the residual sum-of-squares. This performs well on data of the form

y = f(x, θ) + eps

(with var(eps) > 0). It fails to indicate convergence on data of the form

y = f(x, θ)

because the criterion amounts to comparing two components of the round-off error. If you wish to test nls on artificial data please add a noise component, as shown in the example below.

So, no noise==no good for R's nls.

• Welcome to the site, @B_D_Dubbya. I took the liberty of formatting your answer, I hope you don't mind. You can find more info about editing your answers on CV here. – gung - Reinstate Monica Sep 13 '12 at 1:54
• I am aware of this issue - hence using the "jitter" function to add some noise – steiny Aug 30 '17 at 13:41