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

Joerg van den Hoff j.van_den_hoff at fzd.de
Mon Apr 28 14:56:52 CEST 2008


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



More information about the R-help mailing list