[R] Predict polynomial problem

Barry Rowlingson b.rowlingson at lancaster.ac.uk
Tue Jan 19 00:31:07 CET 2010


 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?

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



More information about the R-help mailing list