Automating model selection criteria production

Question

I want to perform model comparison according to several criteria using R. My dataframe's name is df

head(df)

   Y  X1  X2   X3   X4
1 18 307 130 3504 12.0
2 15 350 165 3693 11.5
3 18 318 150 3436 11.0
4 16 304 150 3433 12.0
5 17 302 140 3449 10.5
6 15 429 198 4341 10.0

I want to compare all possible combinations with Y as the dependent variable. I'll give you what I've done so far to get the R.squared for all possible models. I would like your comments for better coding, especially on getting all those model formulas

# First get all possible combinations of the 4 inependent variables
comb_list <- lapply(1:4,function(i) combn(4,i))

# now create the formulas
comb_list_forms <- lapply(unlist(
  sapply(1:length(comb_list), function(i)
  sapply(1:dim(comb_list[[i]])[2], function(x) {nam <- names(df[1+comb_list[[i]][,x]])
  formul <- "Y~"
 sapply(1:length(nam), function(y) formul <<- paste(formul,nam[y],sep="+"))
 formul <- sub("~+", "~", formul,fixed = T)
 }))),as.formula)

# finally get r.squared
attach(df)
models.r.sq <- sapply(comb_list_forms, function(i) summary(lm(i))$r.squared)

Thank you

Answer 1

I think package leaps will help you. It is designed specifically to select a "best" subset of variables for linear regression.

As for the formula creation, you can use foreach package like this:

 #Create the list where each list element contains one combination
  xcomb <- foreach(i=1:4, .combine=c) %do% { 
       combn(paste("X", 1:4, sep=""), i, simplify=FALSE)
 }
 #Convert list elements to formula   
 formlist <- lapply(xcomb, function(l) formula(paste("Y", paste(l, collapse="+"), sep="~")))

Answer 2

Ok, here is another way (which is certainly not as elegant or concise than the other ones, but it has the merit of not requiring Emacs to check parenthesis matches :-)

Let's say we have a vector of predictors of interest like this

> Xs <- paste("X", 1:4, sep="")

Then, we can just use

> allXs <- lapply(seq(along=Xs), 
                  function(x) combn(Xs, x, simplify=FALSE, 
                                    FUN=paste, collapse="+"))

where

> unlist(allXs)

gives you all 15 combinations of the X's.

Another option is to just change the right-hand side of a formula, say

 fm <- as.formula(paste("Y ~ ", paste(Xs, collapse= "+")))

so as to reflect the different combinations that you enumerate in your comb_list object. This can be done using the Formula package:

> fm <- as.formula(paste("Y ~ ", paste(Xs, collapse= "|")))
> fm2 <- Formula(fm)
> foo <- function(x) formula(fm2, rhs=x, collapse=TRUE)
> foo(1:2)
Y ~ X1 + X2
<environment: 0x102593040>
> foo(c(1,3,4))
Y ~ X1 + X3 + X4
<environment: 0x1021d02a0>

Answer 3

This is another way to get the formulas

comb_list_forms <- unlist(lapply(1:length(comb_list), function(i)
  lapply(1:dim(comb_list[[i]])[2], function(x) 
    as.formula(paste("Y ~ ", paste(names(df[1+comb_list[[i]][,x]]), collapse= "+"))))))

Answer 4

You could consider using step on the full model. This uses AIC rather than $R^2$ as the criterion for selection as this is considered to be superior.

In any case, the functions add1 and drop1 that step uses may be of use.