[R] error using nls with logistic derivative

Francisco Mora Ardila fmora at oikos.unam.mx
Tue Apr 17 18:42:18 CEST 2012

```It seems to work!

Thanks and sorry for the personal message

Francisco

On Tue, 17 Apr 2012 09:09:24 +0200, peter dalgaard wrote
> On Apr 17, 2012, at 06:23 , Francisco Mora Ardila wrote:
>
> > Hi
> >
> > I´m trying to fit a nonlinear model to a derivative of the logistic function
> >
> > y = a/(1+exp((b-x)/c)) (this is the parametrization for the SSlogis function with
nls)
> >
> > The derivative calculated with D function is:
> >
> >> logis<- expression(a/(1+exp((b-x)/c)))
> >> D(logis, "x")
> > a * (exp((b - x)/c) * (1/c))/(1 + exp((b - x)/c))^2
> >
> > So I enter this expression in the nls function:
> >
> > ratelogis <- nls(Y ~ a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,
> > start=list(a = 21.16322, b = 8.83669, c = 2.957765),
> > )
> >
> > The data is:
> >
> >> Y
> > [1]  5.5199668  1.5234525  3.3557000  6.7211704  7.4237955  1.9703127
> > [7]  4.3939336 -1.4380091  3.2650180  3.5760906  0.2947972  1.0569417
> >> X
> > [1]   1   0   0   4   3   5  12  10  12 100 100 100
> >
> > The problem is that I got the next error:
> >
> > Error en nls(Y ~ a * (exp((b - X)/c) * (1/c))/(1 + exp((b - X)/c))^2,  :
> >  step factor 0.000488281 reduced below 'minFactor' of 0.000976563
> >
> > I trien to change the minFactor using the control argument inside nls
> >
> > control=nls.control(maxiter=50, tol=1e-5, minFactor = 1/2048
> >
> > but got a new error message:
> >
> >
> > Error en nls(Y ~ a * (exp((b - X)/c) * (1/c))/(1 + exp((b - X)/c))^2,  :
> >  step factor 0.000244141 reduced below 'minFactor' of 0.000488281
> >
> > So it seems that as I modify minFactor, the step factor reduces also and I can never
> > reach a solution.
> >
> > Does anybody Know what am I doing wrong? Is there a problem with the formula? How
can I
> > solve it? I tried some suggestions on R-help related topics but did not work.
>
> (Please don't send private messages. I don't do free consulting outside of the
> mailing lists.)
>
> With nls(), there's really no substitute for good data and good starting
> values. Did you actually try plotting those data? At best, they are extremely
> noisy. How did you obtain the starting values? They seem remarkably accurate
> for something that fits data so poorly.
>
> What is happening is that the algorithm shoots off into a region of the
> parameter space where it can't make any progress:
>
> > ratelogis <- nls(Y ~ a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,+ start=list(a =
21.16322, b = 8.83669, c = 2.957765),+ trace=TRUE)
> 151.098 :  21.163220  8.836690  2.957765
> 127.1149 :  103.49326  11.43274  20.29663
> 111.344 :  833.02386 -45.86244 140.32985
> 111.3375 :  978.97214 -76.20571 159.90833
> 111.3374 :  1097.1376 -101.6771  174.2037
> 111.3228 :  1179.8451 -119.7416  183.3794
>
> Notice that the "b" parameter which is supposed to be the maximum point starts
> off well to the right of the position of the largest Y values, then shoots
> into large negative values.
>
> With a little eyeballing, I can do better:
>
> > ratelogis <- nls(Y ~ a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,
> + start=list(a = 40, b = 3.5, c = 10),trace=T)
> 139.5173 :  40.0  3.5 10.0
> 97.28614 :  26.513424  4.052639  2.267709
> 81.53127 :  32.384910  3.473411  2.307504
> 65.11387 :  39.53542  3.01097  2.07678
> 50.66328 :  36.223529  2.590102  1.133965
> 50.50921 :  35.984466  2.565698  1.067731
> 50.50162 :  36.149993  2.573420  1.081673
> 50.50145 :  36.129962  2.572195  1.079504
> 50.50144 :  36.133656  2.572402  1.079862
> 50.50144 :  36.133061  2.572368  1.079803
>
> However, the fit isn't exactly stellar. Try this:
>
> X <- c(1, 0, 0, 4, 3, 5, 12, 10, 12, 100, 100, 100)
> Y <- c(5.5199668, 1.5234525, 3.3557, 6.7211704, 7.4237955, 1.9703127,
>    4.3939336, -1.4380091, 3.265018, 3.5760906, 0.2947972, 1.0569417)
> plot(X,Y)
> ratelogis <- nls(Y ~ a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,
>  start=list(a = 40, b = 3.5, c = 10),trace=T)
> x0 <-  seq(0,100,,1000)
> lines(x0,predict(ratelogis,newdata=data.frame(X=x0)))
> lines(x0,evalq(a*(exp((b-X)/c)*(1/c))/(1 + exp((b-X)/c))^2,envir=list(a =
> 21.16322, b = 8.83669, c = 2.957765, X=x0)), lty="dashed")
>
> --
> Peter Dalgaard, Professor,
> Center for Statistics, Copenhagen Business School
> Solbjerg Plads 3, 2000 Frederiksberg, Denmark
> Phone: (+45)38153501
> Email: pd.mes at cbs.dk  Priv: PDalgd at gmail.com

----------------------
Francisco Mora Ardila