Singular gradient erros, NLS in R

I'm trying to fit nls(Mound~ a*kg.bag.collar^b + c, start = list(a = 83, b = -.5, c=100), data=test) using the dataset here. I've fit it without trouble without the c term with no problem. But adding the c term gives "Error in nls(mToMound.stand ~ a * as.numeric(kg.bag.collar)^b + c, start = list(a = 83, : singular gradient." I know the error can mean that poor starting values are given, but I've tried a wide range, and more to the point, a and be here are very close to the estimates that come from the fit without "c" included. The subject-matter dictates that I need the c term to allow the possibility of a concave up decreasing function that doesn't go through the origin.

What are my options here? Thanks!

• You provide a b variable but it doesn't appear in the formula: typographical error? As far as options go, you can help us out by describing your data: what is the range of values of kg.bag.collar and how many observations do you have? Could you post a scatterplot of the (kg.bag.collar, mound) values?
– whuber
May 19, 2019 at 21:01
• Thanks, whuber. Yes. b missing was a typo. I added a scatterplot, which shows the range of each variable. There are 57 observations. May 20, 2019 at 0:14
• I don't think you can get a good fit with that model. You have pretty strong heteroscedasticity. In particular the variance at your low levels is high. This means the parameter estimate for a is quite uncertain. Your data also doesn't show a clear asymptote for large values yet. This means a parameter estimate for c would be quite uncertain. Finally, you can expect a strong covariance between these two parameters. In total, this makes for a challenging optimization problem. May 20, 2019 at 8:59
• This mean, if you can't assume c to be zero (or at least greater than zero), you need to get more data over a wider range of x values. Or maybe use a Bayesian approach if you have a good idea what to expect as a value for c. May 20, 2019 at 9:02
• All you have to do is plot the function $y= 83x^{-0.5} + 100$ and compare it to your scatterplot: the graph isn't anywhere near the points. For instance, near $x=100$ your points range between $0$ and $15$ in height, but $y=83/\sqrt{100} + 100 = 108.3$ is far different from any of them.
– whuber
May 20, 2019 at 14:13

You can get successful convergence by using the Golub-Pereyra algorithm for partially linear least-squares models (which is a more sophisticated way of doing what @whuber suggested):

fit <- nls(mound~ cbind(1, kg.bag.collar^b), start = list(b = -.5),
data=test, algorithm = "plinear")
plot(mound ~ kg.bag.collar, data = test)
curve(predict(fit, newdata = data.frame(kg.bag.collar = x)), add = TRUE)


That looks somewhat reasonable but summary output shows that none of the parameters is significant (note that .lin1 is the parameter c and .lin2 is a).

summary(fit)
#Formula: mound ~ cbind(1, kg.bag.collar^b)
#
#Parameters:
#      Estimate Std. Error t value Pr(>|t|)
#b       0.7857     1.1344   0.693   0.4915
#.lin1  18.7878    10.5973   1.773   0.0819 .
#.lin2  -0.2737     1.7254  -0.159   0.8745
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
#Residual standard error: 5.402 on 54 degrees of freedom
#
#Number of iterations to convergence: 11
#Achieved convergence tolerance: 9.472e-06


If you use starting values close to this solution, the default Gauss-Newton algorithm also converges.

You have strong collinearity, which means there are more parameters in your model than your data supports:

vcov(fit)
#               b     .lin1      .lin2
#b       1.286925 -11.82171   1.956174
#.lin1 -11.821714 112.30234 -18.076169
#.lin2   1.956174 -18.07617   2.977021


You should try and collect more and better data if you need to fit this model (but that might not be possible).

• Thanks, @Roland. This is helpful, in particular the G-P algorithm for finding starting values to use in the default G-N algorithm. I like the G-N algorithm because it allows me to use confint2 to get 95% CIs. Can you explain or point me to a resource to better understand what you mean by 'collinearity' in this context? I know the term with respect to >1 independent variable in a model being highly correlated (e.g., tinyurl.com/y25ld8ow), but I only have 1 predictor here (kg.bag.collar). I don't see any resources pointing to another meaning, either. Thanks again. May 24, 2019 at 15:38