# R: Model selection with categorical variables using leaps and glmnet

I have a linear model containing a few continuous variables and four categorical variables, each represented by 12, 3, 4, and 5 dummy variables respectively. When using model selection criteria such as PRESS, Mallows's Cp, and BIC using the leaps package, the best model returned for each criterion contains only some of the dummy variables for each categorical variable. It is my understanding that this is not good practice, and either all or none of the dummies must be included in the model. Is there a way to have leaps treat the dummy variables for each categorical variable as one variable?

Also, could this method be extended to use with the glmnet package? I'm having the same issue with lasso and ridge regression.

EDIT: Is there a way to specify an lm object with a subset of independent variables to be treated as one?

• Whether all the categories for a dummy should or may not be included in the model depends on what you want to use the model for. Do you want to make statistical inferences about effects, or do you just want to be able to make accurate predictions? – Matthew Drury May 13 '15 at 15:34
• For the purposes of my experiment I want to include all dummies for each categorical variable. I know it might have more predictive power if this were not the case but that's not my goal at the moment. – user3623888 May 13 '15 at 15:37
• Ok, cool. Have you checked out the fused lasso? I believe it is meant to address this kind of thing, but I don't know much (any) of the details: stanford.edu/group/SOL/papers/fused-lasso-JRSSB.pdf – Matthew Drury May 13 '15 at 15:38
• @MatthewDrury: To spare me reading all that: did you mean the group lasso? – Scortchi May 13 '15 at 15:56
• @user3623888: With shrinkage methods you're penalizing coefficients based on their magnitude, so it's not a simple fix. – Scortchi May 13 '15 at 16:00

regsubsets (a function in the leaps package that also performs exhaustive model searches) can accept categorical variables that are not split out into dummy variables and, thus, treats them as groups of variables that are either all part of a model or not.

For example, if Year has levels 2013, 2014 and Treatment has levels C,N,O I can run the following statement:

> search_output<-regsubsets(y~Year+Treatment,data=stats_df, method="exhaustive")


Output:

Subset selection object
Call: regsubsets.formula(mu_ln ~ Year + Treatment, data = SS_stats_df,
nbest = 1, method = "exhaustive")
3 Variables  (and intercept)
Forced in Forced out
Year2014       FALSE      FALSE
TreatmentN     FALSE      FALSE
TreatmentO     FALSE      FALSE
1 subsets of each size up to 3
Selection Algorithm: exhaustive

> summary(search_output)$which (Intercept) Year2014 TreatmentN TreatmentO 1 TRUE FALSE TRUE FALSE 2 TRUE FALSE TRUE TRUE 3 TRUE TRUE TRUE TRUE  When faced with this same problem I found this post very helpful (my answer here is essentially an abbreviated version of the pertinent portion): http://rstudio-pubs-static.s3.amazonaws.com/2897_9220b21cfc0c43a396ff9abf122bb351.html And for recoding or converting to factors or renaming factors these posts are helpful: https://stackoverflow.com/questions/5372896/recoding-variables-with-r http://www.cookbook-r.com/Manipulating_data/Recoding_data/ • This does not appear to answer the question by OP. In your example, line 1 of the summary reads > summary(search_output)$which (Intercept) Year2014 TreatmentN TreatmentO 1 TRUE FALSE TRUE FALSE So, the level TreatmentN is being kept, but Treatment0 is not. Is there a way to force all of Treatment to be kept or dropped? Particularly with categorical variables with many levels, this would be useful. – AmeliaMN Dec 12 '16 at 19:59
• Yes, you should be able to use force.in (cran.r-project.org/web/packages/leaps/leaps.pdf) to requires whichever variables you want to be kept in all models, though I haven't tried that myself. If runtime is not an issue, you can also just run without forcing any, then sort or subset the output based on presence of desired variables. That approach is needed, for example, when dealing with interactions in regsubsets because it cannot be set to require the appropriate main effects when interactions are present. – DirtStats Dec 13 '16 at 1:02
• Thanks @DirtStats. That's not what I or the OP are looking for. Instead, we want individual variables to be selected or not selected based on the overall p-value. I'm not trying to force.in=Treatment, but rather I want it to be forced in or out all together (so Treatment0 cannot be forced out without TreatmentN going along). Do you know of any methods that do this? It's almost like regsubsets() for ANOVA. – AmeliaMN Dec 13 '16 at 19:23