Re: [R] nls() newbie convergence problem

From: Katharine Mullen <kate_at_few.vu.nl>
Date: Thu, 05 Jun 2008 02:59:21 +0200 (CEST)

dydata

   x1 x2 y
1 9 27 248
2 9 22 213
3 4 23 190
4 11 16 183
5 -6 25 144
6 11 14 169
7 -4 13 72
8 2 8 73
9 10 13 156
10 8 30 263
11 12 10 147
12 -7 5 -2
13 0 10 75
14 12 0 77
15 9 8 115
16 12 24 245
17 34 23 370
18 12 1 84
19 10 37 324
20 26 30 371

weight <- function(y,x1,x2,b0,b1,b2)
{
  pred <- b0+b1*x1 + b2*x2
  parms <- abs(b1*b2)^(1/3)
  (y-pred)/parms
}

gmfit <- nls(~weight(y,x1,x2,b0,b1,b2),data=dydata,trace=TRUE,

              start=list(b0=1,b1=1,b2=1), control=list(warnOnly=TRUE))

library(minpack.lm)

weight2 <- function(p, y,x1,x2)
{
  pred <- p[1]+p[2]*x1 + p[3]*x2
  parms <- abs(p[2]*p[3])^(1/3)
  (y-pred)/parms
}

gmfit2 <- nls.lm(par = c(1,1,1),

                  fn = weight2, y=dydata$y,x=dydata$x1,x2=dydata$x2,
                  control=list(nprint=1))


dydata is from Norman R. Draper, Yonghong (Fred) Yang, Generalization of the geometric mean functional relationship, Computational Statistics & Data AnalysisVolume 23, Issue 3, 9 January 1997, Pages 355-372; I don't think the attachment got through.

nls uses an orthogonality convergence criterium. Bates talks about this in the post at
https://stat.ethz.ch/pipermail/r-devel/2000-August/021059.html and it's also described in the book by Bates and Watts (Nonlinear regression analysis its applications). Apparently here the residuals are so small they are effectively 0, and this criteria does not work; there is a warning about using zero-residual data in the help page for nls.

nls.lm from the package minpack.lm stops if any of a few different criteria are met; in this case it stops at the time you think is right.

On Wed, 4 Jun 2008, Bernard Leemon wrote:

> I'm sure this must be a nls() newbie question, but I'm stumped.
> I'm trying to do the example from Draper
> and Yang (1997). They give this snippet of S-Plus code:
>
> Specify the weight function:
> weight < - function(y,x1,x2,b0,b1,b2)
> {
> pred <- b0+b1*x1 + b2*x2
> parms <- abs(b1*b2)^(1/3)
> (y-pred)/parms
> }
> Fit the model
> gmfit < -nls(~weight(y,x1,x2,b0,b1,b2), observe,list("starting value"))
>
> in converting this to R, I left the weight function alone and replaced the
> nls() with
>
> gmfit <-
> nls(~weight(y,x1,x2,b0,b1,b2),data=dydata,trace=TRUE,start=list(b0=1,b1=1,b2=1))
>
> where dydata is the appropriate data.frame.
>
> nls() quickly (6 iterations) finds the exact values from Draper & Yang for
> b0, b1, and b2 but
> despite reporting a discrepancy of only 3.760596e-29 by the 7th iteration,
> it merrily goes on
> to 50 iterations and thinks it never converged. how do I tell nls() that
> I'm actually quite
> happy with 3.760596e-29 and it need not work further?
>
> thanks for any suggestions.
>
> gary mcclelland (aka bernie)
> univ of colorado
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help_at_r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Thu 05 Jun 2008 - 01:30:22 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 05 Jun 2008 - 03:30:37 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive