and the values in those
places are different:
On Tue, 19 Jan 2010, Barry Rowlingson wrote:
On Tue, Jan 19, 2010 at 1:36 AM, Charles C. Berry <cbe...@tajo.ucsd.edu> wrote:
Its the environment thing.
I think you want something like this:
models[[i]]=lm( bquote( y ~ poly(x,.(i)) ), data=d)
Use
terms( mmn[[3]] )
both with and without this change and
ls( env = environment( formula( mmn[[3]] ) ) )
get("i",env=environment(formula(mmn[[3]])))
sapply(mmn,function(x) environment( formula( x ) ) )
to see what gives.
Think I see it now. predict involves evaluating poly, and poly here
needs 'i' for the order. If the right 'i' isn't gotten when predict is
called then I get the error. Your fix sticks the right 'i' into the
environment when predict is called.
I haven't quite got my head round _how_ it does it, and I have no
idea how I could have figured this out for myself. Oh well...
Per ?bquote, "bquote quotes its argument except that terms wrapped in
'.()' are evaluated in the specified 'where' environment.
(which by default is the parent.frame)
Note:
i <- 20
bquote(y ~ poly(x,.(i)))
y ~ poly(x, 20)
So, now 'i' is irrelevant as the expression returned by bquote has '20' as
the 'degree' arg.
The following lines are also illustrative:
d = data.frame(x=1:10,y=runif(10))
i=3
#1 naive model:
m1 = lm(y~poly(x,i),data=d)
#2,3 bquote, without or with i-wrapping:
m2 = lm(bquote(y~poly(x,i)),data=d)
m3 = lm(bquote(y~poly(x,.(i))),data=d)
#1 works, gets 'i' from global i=3 above:
predict(m1,newdata=data.frame(x=9:11))
#2 fails - why?
predict(m2,newdata=data.frame(x=9:11))
Well, the terms() objects are the same:
all.equal(terms(m1),terms(m2))
[1] TRUE
but they will look in different places for 'i':
environment(terms(m2))
<environment: 0x01b7c178>
environment(terms(m1))
<environment: R_GlobalEnv>
and the values in those places are different:
environment(terms(m2))$i
[1] 2
environment(terms(m1))$i
[1] 3
#3 works, gets 'i' from within:
predict(m3,newdata=data.frame(x=9:11))
It doesn't need 'i', because the i was evaluated and substituted by
bquote. That is, it doesn't get("i") as the expression returned by bquote
has no 'i' in it.
HTH,
Chuck
rm(i)
#1 now fails because we removed 'i' from top level:
predict(m1,newdata=data.frame(x=9:11))
#2 still fails:
predict(m2,newdata=data.frame(x=9:11))
#3 still works:
predict(m3,newdata=data.frame(x=9:11))
Thanks
--
blog: http://geospaced.blogspot.com/
web: http://www.maths.lancs.ac.uk/~rowlings
web: http://www.rowlingson.com/
twitter: http://twitter.com/geospacedman
pics: http://www.flickr.com/photos/spacedman
Charles C. Berry (858) 534-2098
Dept of Family/Preventive Medicine
E mailto:cbe...@tajo.ucsd.edu UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
______________________________________________
R-help@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.