[R] fit nonlinear model with groups

Douglas Bates bates at stat.wisc.edu
Fri Feb 8 15:31:51 CET 2002


Bill Simpson <W.Simpson at gcal.ac.uk> writes:

> I wanted to fit a nonlinear model using nls to two groups. Here is a toy
> example.
> 
> library(nlme)
> x<-rep(1:5,2)
> cond<-c(rep("a",5),rep("b",5))
> cond<-as.factor(cond)
> y[cond=="a"]<-1/x +2
> y[cond=="b"]<-2/x+3
> df<-data.frame(x,y,cond)
> fit<-nls(y~c(a1,b1)[cond]/x+c(a2,b2)[cond],data=df,
> start=list(a1=1,b1=2,a2=2,b2=3))
> 
> I based above on an example posted by Bill Venables.
> 
> I get message:
> maximum number of iterations exceeded
> even if I reset maxiter using nls.control and even though I start the
> params at the right place.
> 
> This example seems to work OK if I do
> fit<-nls(y~c(a1,b1)[cond]/x+a2,data=df,
> start=list(a1=1,b1=2,a2=2))
> So maybe the trick using c()[cond] can be used only once per formula.
> 
> How to do it? Thanks for any help.

In a situation like this the first thing to do is to add the optional
argument trace = TRUE to the call to nls.  This provides some
information on the progress of the iterations.

Doing that will probably show that the iterations are at the right
place but the convergence criterion is not below the threshhold.  The
reason is that the convergence criterion is a relative criterion that
compares the expected change in the sum-of-squares at the next step to
the current residual sum-of-squares.

*When the data are simulated without any random noise* this takes the
form of (round-off error 1)/(round-off error 2).

Just add some random noise to y and you should be ok.

And no I don't regard this as a bug in nls, I regard it as a bug in
the simulation.

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._



More information about the R-help mailing list