9
$\begingroup$

The R function cv.glm (library: boot) calculates the estimated K-fold cross-validation prediction error for generalized linear models and returns delta. Does it make sense to use this function for a lasso regression (library: glmnet) and if so, how can it be carried out? The glmnet library uses a cross-validation to get the best turning parameter, but I did not find any example that cross-validates the final glmnet equation.

$\endgroup$
  • 2
    $\begingroup$ It certainly makes sense, & although LASSO only optimizes over one (hyper-)parameter, if you want to get the best estimate you can of out-of-sample performance you need an outer cross-validation loop. See Can you overfit by training machine learning algorithms using CV/Bootstrap? & the cited papers. $\endgroup$ – Scortchi Jun 10 '15 at 17:50
  • $\begingroup$ It's recommended to use glmnet rather than glm esp. if you're doing CV. If I recall, glm package becomes painful to use. Also, cv.glmnet exposes more parameters than cv.glm. $\endgroup$ – smci Feb 1 '17 at 3:03
9
$\begingroup$

An example on how to do vanilla plain cross-validation for lasso in glmnet on mtcars data set.

  1. Load data set.

  2. Prepare features (independent variables). They should be of matrix class. The easiest way to convert df containing categorical variables into matrix is via model.matrix. Mind you, by default glmnet fits intercept, so you'd better strip intercept from model matrix.

  3. Prepare response (dependent variable). Let's code cars with above average mpg as efficient ('1') and the rest as inefficient ('0'). Convert this variable to factor.

  4. Run cross-validation via cv.glmnet. It will pickup alpha=1 from default glmnet parameters, which is what you asked for: lasso regression.

  5. By examining the output of cross-validation you may be interested in at least 2 pieces of information:

    • lambda, that minimizes cross-validated error. glmnet actually provides 2 lambdas: lambda.min and lambda.1se. It's your judgement call as a practicing statistician which to use.

    • resulting regularized coefficients.

Please see the R code per the above instructions:

# Load data set
data("mtcars")

# Prepare data set 
x   <- model.matrix(~.-1, data= mtcars[,-1])
mpg <- ifelse( mtcars$mpg < mean(mtcars$mpg), 0, 1)
y   <- factor(mpg, labels = c('notEfficient', 'efficient'))

library(glmnet)

# Run cross-validation
mod_cv <- cv.glmnet(x=x, y=y, family='binomial')

mod_cv$lambda.1se
[1] 0.108442

coef(mod_cv, mod_cv$lambda.1se)
                     1
(Intercept)  5.6971598
cyl         -0.9822704
disp         .        
hp           .        
drat         .        
wt           .        
qsec         .        
vs           .        
am           .        
gear         .        
carb         .  

mod_cv$lambda.min
[1] 0.01537137

coef(mod_cv, mod_cv$lambda.min)
                      1
(Intercept)  6.04249733
cyl         -0.95867199
disp         .         
hp          -0.01962924
drat         0.83578090
wt           .         
qsec         .         
vs           .         
am           2.65798203
gear         .         
carb        -0.67974620

Final comments:

  • note, the model's output says nothing about statistical significance of the coefficients, only values.

  • l1 penalizer (lasso), which you asked for, is notorious for instability as evidenced in this blog post and this stackexchange question. A better way could be to cross-validate on alpha too, which would let you decide on proper mix of l1 and l2 penalizers.

  • an alternative way to do cross-validation could be to turn to caret's train( ... method='glmnet')

  • and finally, the best way to learn more about cv.glmnet and it's defaults coming from glmnet is of course ?glmnet in R's console )))

$\endgroup$
  • $\begingroup$ Nice answer. .. it is also worth running the CV multiple times and averging the error curve (see?cv.glmnet) to account for the sampling. $\endgroup$ – user20650 Jun 11 '15 at 22:17
  • $\begingroup$ @SergeyBushmanov VERY USEFUL! $\endgroup$ – theforestecologist Oct 5 '16 at 3:45
  • $\begingroup$ Hi, I know that this is an old post, but I wanted to ask you a question. You mention that the model's output says nothing about statistical significance of the coefficients, so how do you determine that they are significant or not? $\endgroup$ – Jun Jang Jun 15 '18 at 12:07
  • $\begingroup$ @JunJang "There is no statistical significance for coefficients" is the statement from authors of the package, not me. This statement is given, I do not remember exactly, either in one of the book of the package authors or in the package's vignette. In such a case, instead of saying coefficients significant or not, you'd rather say they are "useful" or not in explaining target through cross validation. $\endgroup$ – Sergey Bushmanov Jun 15 '18 at 12:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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