Predict polynomial problem
On Tue, 2010-01-19 at 08:27 +0000, Barry Rowlingson wrote:
On Tue, Jan 19, 2010 at 1:36 AM, Charles C. Berry <cberry at 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...
Perhaps this Programmer's Niche article by Bill Venables might also be useful as it discuss how to manipulate formulas to automate model fitting...? Bill Venables. Programmer's Niche. R News, 2(2):24-26, June 2002. http://cran.r-project.org/doc/Rnews/Rnews_2002-2.pdf HTH G
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)) #3 works, gets 'i' from within: predict(m3,newdata=data.frame(x=9:11)) 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
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~% Dr. Gavin Simpson [t] +44 (0)20 7679 0522 ECRC, UCL Geography, [f] +44 (0)20 7679 0565 Pearson Building, [e] gavin.simpsonATNOSPAMucl.ac.uk Gower Street, London [w] http://www.ucl.ac.uk/~ucfagls/ UK. WC1E 6BT. [w] http://www.freshwaters.org.uk %~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%