[R] Dropping terms from regression w/ poly()

Joshua Stults joshua.stults at gmail.com
Thu Jun 4 16:07:00 CEST 2009


Hello r-help,

I'm fitting a model with lm() and using the orthogonal polynomials
from poly() as my basis:

dat <- read.csv("ConsolidatedData.csv", header=TRUE)
attach(dat)
nrows <- 1925
Rad <- poly(Radius, 2)
ntheta <- 14
Theta <- poly(T.Angle..deg., ntheta)
nbeta <- 4
Beta <- poly(B.Beta..deg., nbeta)

model.1 <- lm( Measurement ~ Block + Rad + Theta + Beta + Rad:Theta +
              Rad:Beta + Theta:Beta)

Which works splendidly, and for my data set shows that the odd orders
in Theta and Beta are not significant (expected because it is a
symmetric measurement that I'm fitting). Now I want to simplify my
model by dropping all the odd-order terms.  I couldn't figure out a
way to do that without defining new design matrices like this:

T2 <- matrix(0,nrows,ntheta)
for(i in 1:ntheta){
  T2[,i] <- Theta[,i]
}
B2 <- matrix(0,nrows,nbeta)
for(i in 1:nbeta){
  B2[,i] <- Beta[,i]
}
R2 <-matrix(0,nrows,2)
R2[,1] <- Rad[,1]
R2[,2] <- Rad[,2]

And then fitting a model using the individual columns of T2, R2 and
B2.  That lets me drop all the odd terms, but it is very cumbersome:

model.2 <- lm( Measurement ~ Block + (R2[,1] + R2[,2]) + (T2[,1] +
              T2[,2] + T2[,3] + T2[,4] + T2[,5] + T2[,6] + T2[,7] +
              T2[,8] + T2[,9] + T2[,10] + T2[,11]+ T2[,12] + T2[,13] +
              T2[,14]) + (B2[,1] + B2[,2] + B2[,3] + B2[,4]) + (R2[,1]
              + R2[,2]) : (T2[,1] + T2[,2] + T2[,3] + T2[,4] + T2[,5]
              + T2[,6] + T2[,7] + T2[,8] + T2[,9] + T2[,10] + T2[,11]+
              T2[,12] + T2[,13] + T2[,14]) + (R2[,1] + R2[,2]) :
              (B2[,1] + B2[,2] + B2[,3] + B2[,4]) + (T2[,1] + T2[,2] +
              T2[,3] + T2[,4] + T2[,5] + T2[,6] + T2[,7] + T2[,8] +
              T2[,9] + T2[,10] + T2[,11]+ T2[,12] + T2[,13] + T2[,14])
              : (B2[,1] + B2[,2] + B2[,3] + B2[,4]))

because I have to call out each column explicitly.  I'm new to R so I
have a suspicion that there is a much better way to accomplish this
(dropping odd terms from a model based on poly), but I haven't figured
it out.  Any hints / suggestions?  Thanks.


-- 
Joshua Stults
Website: http://j-stults.blogspot.com




More information about the R-help mailing list