[R] weighted nonlinear fits: `nls' and `eval'

From: Joerg van den Hoff <j.van_den_hoff_at_fzd.de>
Date: Mon, 28 Apr 2008 14:56:52 +0200


dear list,

my question concerns the use of `eval' in defining the model formula for `nls' (version 2.6.2.).

consider the following simple example, where the same model and data are used to perform unweighted and weighted fits. I intentionally used very uneven weights to guarantee large differences in the results

#================================CUT===========================

ln      <- 100
x       <- 1:ln
y0      <- exp(-.02*x)
y       <- rnorm(y0, y0, .01)

wts     <- rep(1, ln)
y[30]   <- 1.2*y[30]

wts[30] <- 1e3
model <- y ~ exp(-k*x)

xydata <- list(x=x, y=y)

#simple unweighted fit works as expected:
r0 <- nls(model, start = c(k=.01), data = list(x=x, y=y))

#simple weighted fit works as expected:
r1 <- nls(model, start = c(k=.01), data = xydata, weights = wts)

#this actually performs an unweighted fit (issuing a warning):
mdl <- model[[3]]
r2 <- nls(y ~ eval(mdl), start = c(k=.01), data = xydata, weights = wts)

#this, too, actually performs an unweighted fit (issuing a warning):
dv1 <- deriv(model, "k")
r3 <- nls(y ~ eval(dv1), start = c(k=.01), data = xydata, weights = wts)

#weighted fit, works as expected

dv2 <- deriv(model, "k", c("k", "x"))
r4 <- nls(y ~ dv2(k, x), start = c(k=.01), data = xydata, weights = wts)
#================================CUT===========================
as you can see the fits producing `r2' and `r3' do _not_ do what's intended. looking at the source code of `nls' I get some ideas where it is going wrong (in setting up the model frame in `mf'?) but not what exactly is the problem (scoping rules?).

`r2' and `r3' can be convinced to do what's intended by modifying `xydata' to

xydata <- list(x=x, y=y, `(weights)` = wts)

i.e. by adding `(weights)` component here, too. but that probably is not the way to go ...

while it is clear in the case of `r2' that my example is artificial, `r3' is what I have used up to now for unweighted fits without problems. moreover there _are_ circumstances where `eval' constructs for defining the model formula occur quite naturally.

questions:

1.)
is it actually not intended/supported to use `eval' in defining the model formula argument of `nls'? if so, why not? or what am I missing (e.g. necessity to provide explicitely an `envir' argument to `eval': which one?)?

2.)
could/should not `nls' be modified to cleanly support the use of `eval' (remarks along the line of "submit changes" would not help much :-)).

3.)
if the policy of the core group is that the current state of affairs is quite perfect, should not the manpage of `nls' be augmented to warn people very clearly that any use of `eval' in the model formula is calling for trouble. I find the example above especially annoying: `r2' and `r3' run apparently fine, only issuing a warning which sure is cryptic for joe user. so, if this happens in some wrapper function provided to others (which I do) they very well might believe having performed a weighted fit when they did not.

joerg



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 Mon 28 Apr 2008 - 13:02:13 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 Mon 28 Apr 2008 - 13:30:32 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