53
$\begingroup$

I am fitting a binomial family glm in R, and I have a whole troupe of explanatory variables, and I need to find the best (R-squared as a measure is fine). Short of writing a script to loop through random different combinations of the explanatory variables and then recording which performs the best, I really don't know what to do. And the leaps function from package leaps does not seem to do logistic regression.

Any help or suggestions would be greatly appreciated.

Improve this question
$\endgroup$
5
  • $\begingroup$ There exist functions that perform automatic search. You should have a look at the step function. Section 5.4 illustrates that point: data.princeton.edu/R/glms.html $\endgroup$
    – ocram
    Mar 15, 2011 at 8:43
  • $\begingroup$ I'm sorry but my post has been edited so that it no longer asks my question. I have 35 (26 significant) explanatory variables in my logistic regression model. I need the best possible combination of 8, not the best subset, and at no point was I interested in a stepwise or all subsets style approach. There is no wiggle room in this 8. I just thought someone may have know how I could fit all combinations of 8 explantory variables and it could tell me which maximises the likelihood (sorry about the R-squared brain fart but AIC isnt relevant either since I have a fixed number of parameters, 8). $\endgroup$
    – Leendert
    Mar 16, 2011 at 10:25
  • $\begingroup$ You can revert to the previous version of your post, or combine both edits. I'm sure @mpiktas was of good intention when trying to improve its appearance and just didn't notice the No. parameters. $\endgroup$
    – chl
    Mar 16, 2011 at 18:09
  • $\begingroup$ @ Everyone: Thank you very much. In the end I used many different thing in the hope they would all give similar answers. And they did. I used the BMA, bestglm and glmnet packages as well as the step function. Fitted models with all of them, and there was not discrepancy in what BMA with maxcol = 9 and step deemed the best model. All the experts in the field around me seemed very content with the variables, and felt that it was quite progressive. So thanks for all the input. I really used it all. $\endgroup$
    – Leendert
    Mar 23, 2011 at 8:22
  • $\begingroup$ glmulti is also a good package for best subset selection, and that one allows you to specify the maximum nr of variables in your model, and also allows one to consider all possible 1st order interaction effects $\endgroup$
    – Tom Wenseleers
    Apr 4, 2017 at 18:51

4 Answers 4

Reset to default
31
$\begingroup$

Stepwise and "all subsets" methods are generally bad. See Stopping Stepwise: Why Stepwise Methods are Bad and what you Should Use by David Cassell and myself (we used SAS, but the lesson applies) or Frank Harrell Regression Modeling Strategies. If you need an automatic method, I recommend LASSO or LAR. A LASSO package for logistic regression is available here, another interesting article is on the iterated LASSO for logistic

Improve this answer
$\endgroup$
6
  • 6
    $\begingroup$ (+1) About R packages, there're also glmnet (implementation with coordinate descent algo., Friedman and coll.) and penalized (allows to keep some var. unpenalized). Of note, F. Harrell provides penalized ML estimation for GLMs (see lrm, or his RMS textbook for further info). $\endgroup$
    – chl
    Mar 15, 2011 at 11:31
  • $\begingroup$ (+1) Nice article, it seems I have to start going far beyond the author states in the question (not the first time I didn't). @chl (+1) perfect alternative suggestions too. $\endgroup$
    – Dmitrij Celov
    Mar 15, 2011 at 15:04
  • $\begingroup$ @chl: +1 for glmnet, that's a great package. $\endgroup$
    – Zach
    Mar 16, 2011 at 16:07
  • 1
    $\begingroup$ @chl Thanks! One of the problems with R is keeping track of packages (there are so many!) and which are best. The task views do help $\endgroup$
    – Peter Flom
    Mar 18, 2011 at 10:42
  • 2
    $\begingroup$ If your variables are collinear it's best to use elastic net using glmnet, say with alpha=0.5, as LASSO tends to randomly kick out highly collinear variables out of the model $\endgroup$
    – Tom Wenseleers
    Apr 4, 2017 at 18:48
17
$\begingroup$

First of all $R^2$ is not an appropriate goodness-of-fit measure for logistic regression, take an information criterion $AIC$ or $BIC$, for example, as a good alternative.

Logistic regression is estimated by maximum likelihood method, so leaps is not used directly here. An extension of leaps to glm() functions is the bestglm package (as usually recommendation follows, consult vignettes there).

You may be also interested in the article by David W. Hosmer, Borko Jovanovic and Stanley Lemeshow Best Subsets Logistic Regression // Biometrics Vol. 45, No. 4 (Dec., 1989), pp. 1265-1270 (usually accessible through the university networks).

Improve this answer
$\endgroup$
5
  • 2
    $\begingroup$ While your comment about $R^2$ being worse than $BIC,AIC$ is useful in general, it actually makes no difference unless you are comparing models of different sizes. The OP clearly states that they are only interest in $8$ variable models, so $BIC$ and $AIC$ will revert back to choosing the model with the highest likelihood. This is equivalent to $R^2$ fitting. $\endgroup$
    – probabilityislogic
    Dec 21, 2011 at 13:32
  • $\begingroup$ Thanks for the remark, but comments by chl below explains why fixed number of explanatory variables is dangerous. By the way the answer appeared earlier than the comment regarding $8$ (up to?) variables rstriction ;) $\endgroup$
    – Dmitrij Celov
    Dec 22, 2011 at 10:18
  • $\begingroup$ any academic reference for R2R2 is not an appropriate goodness-of-fit measure for logistic regression take an information criterion AICAIC or BICBIC $\endgroup$
    – SIslam
    Mar 23, 2017 at 11:46
  • $\begingroup$ a note on bestglm, it uses leaps at the backend for calculation! So it will fail if there is NA in the dataset and will come up with message like Error in leaps.setup(x, y, wt = weights, nbest = nbest, nvmax = nvmax, : NA/NaN/Inf in foreign function call (arg 3) But what is interesting, my dataset has no NA rather some zeros nonetheless this function complains and gives that exact above message!! $\endgroup$
    – SIslam
    Mar 23, 2017 at 12:38
  • $\begingroup$ glmnet is also a good one, and can also do models that consider all possible 1st order interaction effects $\endgroup$
    – Tom Wenseleers
    Apr 4, 2017 at 18:50
6
$\begingroup$

One idea would be to use a random forest and then use the variable importance measures it outputs to choose your best 8 variables. Another idea would be to use the "boruta" package to repeat this process a few hundred times to find the 8 variables that are consistently most important to the model.

Improve this answer
$\endgroup$
9
  • $\begingroup$ @Zach Do you suggest to rely on RFs to perform feature selection, and then apply a GLM -- in this case, there's risk of overfitting or over-optimism --, or to use RFs (with standard measures of var. importance, or all-relevant selection) as a standalone tool? $\endgroup$
    – chl
    Mar 16, 2011 at 18:05
  • $\begingroup$ @chl: I was suggesting using the RFs to perform feature selection, and then apply the GLM. I agree there's risk of over fitting, but the OP said he needed exactly 8 variables. $\endgroup$
    – Zach
    Mar 16, 2011 at 20:47
  • 1
    $\begingroup$ @Zach "exactly 8 variables"... Then, you are somewhat soft-thresholding the potential variables of interest based on a measure of variable importance (which thanks to being based on permutation and double resampling is supposed to be bias-free) and then reinject them into an GLM. IMHO, you break the control on overfitting exerted through bagging. This is also described in ESLII from Hastie et al.: feature selection, if any, must be included in the cross-validation procedure (where cross-validation include the assessment of model performance). $\endgroup$
    – chl
    Mar 16, 2011 at 20:57
  • $\begingroup$ @chl: Doesn't the r package "boruta" cross validate by running the random forest several time? Are you basically saying we need "meta-cross-validation," where you do variable selection and and the fit your model on random subsets of the data? $\endgroup$
    – Zach
    Mar 17, 2011 at 2:49
  • 3
    $\begingroup$ @Zach My point was that bagging, or in particular RFs, prevent from overfitting (to a certain extent) provided you remain in the same toolchain. If you use the results from RFs and then go see how another model would perform with the same data, then you break the CV loop. But why not use directly RFs for classification? Another solution would be to use RFs on a training sample and then apply the GLM on an held-out sample (that could potentially be cross-validated too). $\endgroup$
    – chl
    Mar 17, 2011 at 7:56
2
$\begingroup$

stats::step function or the more general MASS::stepAIC function supports lm, glm (including logistic regression) and aov family models.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.