[R] nls problem
joerg van den hoff
j.van_den_hoff at hzdr.de
Tue Jul 3 13:54:11 CEST 2012
hi list,
used versions: 2.12.1 and 2.14.0 under ubuntu and macosx.
I recently stumbled over a problem with `nls', which occurs if the model
is not specified explicitly but via an evaluation of a 'call' object.
simple example:
8<--------------------------------------------------------------------------------------
nlsProblem <- function (len = 5) {
#=======================================================================
# purpose: to demonstrate an apparent problem with `nls' which occurs,
# if the model is specified by passing th lhs as an evaled 'call'
# object. The problem is related to the way `nls' tries to compute
# its internal variable `varIndex' which rests on the assumption that
# the dependent ("y") and, possibly, the independent ("x") variable
# are identified by having a length equal to the `nls' variable
# `respLength'. the problem arises when there are `varNames'
# components accidentally having this length, too.
# in the present example, setting the number of data points to
# len=2 triggers the error since the `call' object `fifu' has this
# length, too and `nls' constructs an erroneous `mf$formula' internally.
#=======================================================================
#generate some data
x <- seq(0, 4, len = len)
y <- exp(-x)
y <- rnorm(y, y, .001*y)
data <- list(x = x, y = y)
#define suitable model
model <- y ~ exp(-k*x)
fifu <- model[[3]]
#this fit is fine:
fit1 <- nls(model, data = data, start = c(k = 1))
print(summary(fit1))
#this fit crashes `nls' if len = 2:
fit2 <- nls(y ~ eval(fifu), data = data, start = c(k = 1))
print(summary(fit2))
}
8<--------------------------------------------------------------------------------------
to see the problem call `nlsProblem(2)'.
as explained in the above comments in the example function, I tracked it
down to the way
`nls' identifies x and y in the model expression. the problem surfaces in
the line
varIndex <- n%%respLength == 0
(line 70 in the function listing from within R) which, in the case of
`fit2' in the above
example always returns a single TRUE index as long as `len != 2' (which
seems fine for the
further processing) but returns a TRUE value for the index of `fifu' as
well if `len == 2'.
question1: I'm rather sure it worked about 6 months ago with (ca.) 2.11.x
under ubuntu. have there been changes in this area?
question2: is something like the `fit2' line in the example expected to
work or not?
qeustion3: if it is not expected to work, should not the manpage include a
corresponding caveat?
question4: is there a a substitute/workaround for the `fit2' line which
still allows to specify the rhs of the model via a variable instead of a
constant (explicit) expression or function call?
the above example is of course construed but in my use case I actually
need this sort of thing. is there any chance that
the way `nls' analyzes its `model' argument can be changed to parse the
`eval(fifu)' construct correctly in all cases?
since I'm currently not subscribed to the list I'd appreciate if responses
could be Cc'ed to me directly.
thanks in advance,
joerg
More information about the R-help
mailing list