5
$\begingroup$

I want to tune lambda of an elastic net. It is possible to optimize lambda based on a cross validation with the glmnet package. However, I recognized that the estimated lambda varies a lot, if I run the cross validation several times. See code below, for a reproducible illustration of this situation:

library("glmnet")

set.seed(1234)

# Some example data
N <- 1000
y <- rnorm(N, 5, 10)
x1 <- y + rnorm(N, 2, 10)
x2 <- y + rnorm(N, - 5, 20)
x3 <- y + rnorm(N, 10, 200)
x4 <- rnorm(N, 20, 50)
x5 <- rnorm(N, - 7, 200)
x6 <- rbinom(N, 1, exp(x1) / (exp(x1) + 1))
x7 <- rbinom(N, 1, exp(x2) / (exp(x2) + 1))
x8 <- rbinom(N, 1, exp(x3) / (exp(x3) + 1))
x9 <- rbinom(N, 1, exp(x4) / (exp(x4) + 1))
x10 <- rbinom(N, 1, exp(x5) / (exp(x5) + 1))

data <- data.frame(y, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10)

# Estimate lambda.1se several times
est_lambda.1se <- numeric()

for(i in 1:1000) {
  # Tuning in glmnet
  mod_cv_test <- cv.glmnet(x = as.matrix(data[ , colnames(data) %in% "y" == FALSE]), 
                       y = y, alpha = 0.5, family = "gaussian")
  # Store lambda.1se
  est_lambda.1se[i] <- mod_cv_test$lambda.1se
}

# Distribution of the estimated lambdas
hist(est_lambda.1se, breaks = 50) # Could the median be taken?

enter image description here

In my example, lambda ranges approximately between 1 and 2.3. I am wondering if there is room for improvements. My idea is to take the median of all runs as best lambda. However, I am afraid that this procedure results in any problems which I can not see right now.

Question: How could the tuning of lambda be improved?

$\endgroup$
3
  • 2
    $\begingroup$ Use leave one out cross validation $\endgroup$ Mar 24, 2017 at 15:49
  • 4
    $\begingroup$ Im not sure i understand the comment above. It would seem to me that while leave one out would give a definitive lambda, it only sweeps the issue under the rug by providing a point estimate of lambda with no variance estimate. The variability of lambda is true variance, the op has done themselves a good service by quantifying it. $\endgroup$ Mar 24, 2017 at 16:00
  • $\begingroup$ Thanks a lot for the comments! It seems to me that k-fold is the standard procedure here, which would be another reason for me to stick with it in addition to the points of Matthew Drury. What do you think about my suggestion to use the median of several runs? I can not see any disadvantages of that, but I am relatively new to cross validation, so maybe I am missing something. $\endgroup$ Mar 24, 2017 at 16:07

2 Answers 2

3
$\begingroup$

I think what you're really asking here is about how sensitive the results are to your choice of $\lambda$, so I suggest you investigate exactly that. When I'm tuning an elastic net, if computation time isn't an issue I like to make a heatmap of the CV loss as a function of $\lambda_1$ and $\lambda_2$ (or $\lambda$ and $\alpha$, depending on your software's parameterization). This exactly tells me the stability of my solution and I can see how much better my optimized point is than a 'random' point. You can also 'zoom in' until the noise of cross validation overpowers the actual increase in model performance (this is like doing a heuristic simplex search). Statements like 'this $\lambda$ is large' aren't really meaningful -- it's all data dependent, so for your data the difference between $\lambda=1$ and $\lambda=2$ might be completely negligible, while for other data sets that's a massive difference.

Here's an example of the kind of heatmap that I'm talking about, taken from this paper (which I suggest you read) by our very own @DikranMarsupial:

heatmap

$\endgroup$
4
  • $\begingroup$ Thank you for your comment @ Chaconne! I just had a look at your suggested paper and it helped a lot to get a better insight! Would you just check the stability of the estimated lambda with your proposed method or would you also try to adjust lambda? What do you think about my suggestion to do several cross validation runs and then take the median of all lambda.1se? $\endgroup$ Mar 28, 2017 at 10:40
  • $\begingroup$ @JoachimSchork if you repeat the CV with a particular $\lambda$ multiple times and average the CV losses you're just getting a lower variance estimate of performance for that particular $\lambda$. If you're really concerned about noise in your assessments of $\lambda$ this can help with that. But when you make a heatmap like the above you don't need to repeat the CV for any particular $\lambda$ because you'll have nearby $\lambda$'s and the loss should be continuously varying at large scales. Does that help? $\endgroup$
    – jld
    Mar 28, 2017 at 17:14
  • $\begingroup$ Yes, that actually helped a lot @ Chaconne. Thank you very much! $\endgroup$ Mar 29, 2017 at 7:03
  • $\begingroup$ > "When I'm tuning an elastic net, if computation time isn't an issue I like to make a heatmap of the CV loss as a function of 𝜆1 and 𝜆2." Just to be clear @jld is it nice convex shapes of the contours that tell you about stability? Or do you also plot the contour and the loss function ponts(like the one with yellow +) to see how they cluster in the valleys?. $\endgroup$
    – Sid
    Mar 11, 2021 at 6:20
2
$\begingroup$

Good question!

You can try the train function in the caret package, it can easily tuning the alpha and lambda parameters at the same time.

What you caught at this unstable situation, it telling you maybe the falling in the local minimum.

I have run this example multiple times, this is what i got:

The local minimum of lambda is 0.201, then I use this new lambda seq(0.01,0.4,0.01) which I got the 'global' minimum is 0.07!

Run the code bellow, it's very impressive! the first plot

eGrid <- expand.grid(alpha = (1:10) * 0.1, 
                 lambda = seq(0.01,0.4,0.01))

Control <- trainControl(method = "repeatedcv",repeats = 3,verboseIter =TRUE)

netFit <- train(x = as.matrix(data[ , colnames(data) %in% "y" == FALSE]), 
            y = y,
            method = "glmnet",
            tuneGrid = eGrid,
            trControl = Control)
plot(netFit)
$\endgroup$
7
  • $\begingroup$ Thank you for your answer @ feng! I ran your code several times and it seems like your proposed plot changes a lot, if I run the code several times. In some runs, it chooses a minima of approximately 0.07, but much more often it chooses a minima of approximately 0.2. Try to run your code with the seed set.seed(123). In this case it chooses a minima of 0.2. If you run the code again with a seed of set.seed(555) it chooses a minima of 0.07. How could I decide which of these minima is better? My idea is to take the median of the minima of several runs. What do you think about that idea? $\endgroup$ Mar 28, 2017 at 7:53
  • $\begingroup$ One further problem with caret is that it calculates just lambda.min. I would like to use lambda.1se, which is the largest value of lambda such that the error is within 1 standard error of the minimum. Do you know a way how I could calculate lambda.1se with caret? $\endgroup$ Mar 28, 2017 at 7:56
  • $\begingroup$ @JoachimSchork that's weird. The code above is not set the seed, I have not reproduce it . I run it with different seed, the lambda always the same (0.18), so 0.18 is corresponding the most minima rmse. $\endgroup$
    – wolfe
    Mar 28, 2017 at 9:34
  • $\begingroup$ There is randomness introduced in the train function. If you run it several times there is a change in the output. Most of the time the best lambda is 0.18, but sometimes it changes to other values. I just ran it 20 times with seeds 1:20 and received 8 different lambdas. Try it for example with a seed of 5, the resultant best lambda would be 0.06. $\endgroup$ Mar 28, 2017 at 10:07
  • 1
    $\begingroup$ Another way to removing the randomness, is try to using two lambada sequences which generate around the two local minimum. Then you will get the unchanged and best results $\endgroup$
    – wolfe
    Mar 28, 2017 at 11:04

Your Answer

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

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