On Mon, 18 Jan 2010, Barry Rowlingson wrote:
I have a function that fits polynomial models for the orders in n:
lmn <- function(d,n){
models=list()
for(i in n){
models[[i]]=lm(y~poly(x,i),data=d)
}
return(models)
}
My data is:
> d=data.frame(x=1:10,y=runif(10))
So first just do it for a cubic:
> mmn = lmn(d,3)
> predict(mmn[[3]])
1 2 3 4 5 6 7 8
0.6228353 0.5752811 0.5319524 0.4957381 0.4695269 0.4562077 0.4586691 0.4798001
9 10
0.5224893 0.5896255
and lets extrapolate a bit:
> predict(mmn[[3]],newdata=data.frame(x=c(9,10,11)))
1 2 3
0.5224893 0.5896255 0.6840976
now let's to it for cubic to quintic:
> mmn = lmn(d,3:5)
check the cubic:
> predict(mmn[[3]])
1 2 3 4 5 6 7 8
0.6228353 0.5752811 0.5319524 0.4957381 0.4695269 0.4562077 0.4586691 0.4798001
9 10
0.5224893 0.5896255
- thats the same as last time. Extrapolate?
> predict(mmn[[3]],newdata=data.frame(x=c(9,10,11)))
Error: variable 'poly(x, i)' was fitted with type "nmatrix.3" but type
"nmatrix.5" was supplied
In addition: Warning message:
In Z/rep(sqrt(norm2[-1L]), each = length(x)) :
longer object length is not a multiple of shorter object length
it falls over. I can't see the difference between the objects,
summary() looks the same. Is something wrapped up in an environment
somewhere, or some lazy evaluation thing, or have I just done
something stupid?
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.
HTH,
Chuck
Here's a complete example you can paste in - R --vanilla < this.R
gives the error above - R 2.10.1 on Ubuntu, and also on R 2.8.1 I had
lying around on a Windows box:
d = data.frame(x=1:10,y=runif(10))
lmn <- function(d,n){
models=list()
for(i in n){
models[[i]]=lm(y~poly(x,i),data=d)
}
return(models)
}
mmn = lmn(d,3)
predict(mmn[[3]])
predict(mmn[[3]],newdata=data.frame(x=c(9,10,11)))
mmn2 = lmn(d,3:5)
predict(mmn2[[3]])
predict(mmn2[[3]],newdata=data.frame(x=c(9,10,11)))
Barry
______________________________________________
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.
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.
