[R] Scoping and nls.
Gabor Grothendieck
ggrothendieck at myway.com
Tue Oct 26 05:30:42 CEST 2004
Rolf Turner <rolf <at> math.unb.ca> writes:
:
: A colleague of mine is trying to use nls() to effect an optimization,
: and is encountering a scoping problem. I should know how to solve it
: for him but .... well, I just don't.
:
: I also had a quick scrounge of the archives --- I know I've seen this
: topic addressed before --- but I couldn't track it down.
:
: So here's a toy example that demonstrates the problem:
:
: hhh <- function(y,x) {
: g <- function(a,b,x) {1/(1+a^2 + b^2*exp(x))}
: nls(y~g(a,b,x),data=data.frame(x=x,y=y))
: }
: set.seed(123)
: x <- runif(50,1,10)
: y <- 1/(1+16 + 36*exp(x)) + rnorm(50,0,0.1)
: hhh(y,x)
:
: which results in an error
:
: Error in get(x, envir, mode, inherits) : variable "g" was not found
:
: So can one assign the function ``g'' somewhere where it **can**
: be found? I.e. use something like
:
: assign("g",g,envir=??????)
:
: Or is there some other magic incantation that can be used here
: to get nls to deal with a ``g'' defined inside the function
: which call nsl()?
:
: In the real problem the function g is a complicated gadget, taking
: different forms in different parts of its domain, so it is
: inconvenient-to-impossible to specify it explicitly in the formula in
: the call to nls().
:
You can assign g in the global environment. Alternately, note from
?nls that functions in the formula are looked up in the environment
of the formula so just be sure to set up your formula so that the
its environment is the one inside hhh.
Thus any of these should work:
1. add the statement assign("g", g, .GlobalEnv) after the definition of g.
(This has the drawback of adding variable g to the your Global Environment
possibly clobbering any g that was already there.)
2. replace the nls statement with:
fo <- y ~ g(a,b,x)
nls(fo,data=data.frame(x=x,y=y),start=c(a=0,b=0))
3. replace the nls statement with:
nls(as.formula("y~g(a,b,x)"),data=data.frame(x=x,y=y),start=c(a=0,b=0))
More information about the R-help
mailing list