How to interpret variables that are excluded from or included in the lasso model?

Question

I got from other posts that one cannot attribute 'importance' or 'significance' to predictor variables that enter a lasso model because calculating those variables' p-values or standard deviations is still a work in progress.

Under that reasoning, is it correct to assert that one CANNOT say that variables that were EXCLUDED from the lasso model are 'irrelevant' or 'insignificant'?

If so, what can I actually claim about the variables that are either excluded or included in a lasso model? In my specific case, I selected the tuning parameter lambda by repeating 10-fold cross-validation 100 times in order to reduce randonmess and to average the error curves.

UPDATE1: I followed a suggestion below and re-ran lasso using bootstrap samples. I had it a go with 100 samples (that amount was what my computer power could manage overnight) and some patterns emerged. 2 of my 41 variables entered the model more then 95% of times, 3 variables more than 90% and 5 variables more than 85%. Those 5 variables are among the 9 that entered the model when I had run it with the original sample and were the ones with the highest coefficient values then. If I run lasso with say 1000 bootstrap samples and those patterns are maintained, what would be the best way to present my results?

• Does 1000 bootstrap samples sound enough? (My sample size is 116)

• Should I list all the variables and how frequently they enter the model, and then argue that those that enter more frequently are more likely to be significant?

• Is that as far as I can go with my claims? Because it is a work in progress (see above) I cannot use a cut-off value, right?

UPDATE2: Following a suggestion below, I have calculated the following: on average, 78% of the variables in original model entered the models generated for the 100 bootstrap samples. On the other hand, only 41% for the other way around. This has to do in great part with the fact that the models generated for the bootstrap samples tended to include much more variables (17 on average) than the original model (9).

UPDATE3: If you could help me in interpreting the results I got from bootstrapping and Monte Carlo simulation, please have a look at this other post.

Answer 1

Your conclusion is correct. Think of two aspects:

1. Statistical power to detect an effect. Unless the power is very high, one can miss even large real effects.
2. Reliability: having a high probability of finding the right (true) features.

There are at least 4 major considerations:

1. Is the method reproducible by you using the same dataset?
2. Is the method reproducible by others using the same dataset?
3. Are the results reproducible using other datasets?
4. Is the result reliable?

When one desires to do more than prediction but to actually draw conclusions about which features are important in predicting the outcome, 3. and 4. are crucial.

You have addressed 3. (and for this purpose, 100 bootstraps is sufficient), but in addition to individual feature inclusion fractions we need to know the average absolute 'distance' between a bootstrap feature set and the original selected feature set. For example, what is the average number of features detected from the whole sample that were found in the bootstrap sample? What is the average number of features selected from a bootstrap sample that were found in the original analysis? What is the proportion of times that a bootstrap found an exact match to the original feature set? What is the proportion that a bootstrap was within one feature of agreeing exactly with the original? Two features?

It would not be appropriate to say that any cutoff should be used in making an overall conclusion.

Regarding part 4., none of this addresses the reliability of the process, i.e., how close the feature set is to the 'true' feature set. To address that, you might do a Monte-Carlo re-simulation study where you take the original sample lasso result as the 'truth' and simulate new response vectors several hundred times using some assumed error structure. For each re-simulation you run the lasso on the original whole predictor matrix and the new response vector, and determine how close the selected lasso feature set is to the truth that you simulated from. Re-simulation conditions on the entire set of candidate predictors and uses coefficient estimates from the initially fitted model (and in the lasso case, the set of selected predictors) as a convenient 'truth' to simulate from. By using the original predictors one automatically gets a reasonable set of co-linearities built into the Monte Carlo simulation.

To simulate new realizations of $Y$ given the original $X$ matrix and now true regression coefficients, one can use the residual variance and assume normality with mean zero, or to be even more empirical, save all the residuals from the original fit and take a bootstrap sample from them to add residuals to the known linear predictor $X\beta$ for each simulation. Then the original modeling process is run from scratch (including selection of the optimum penalty) and a new model is developed. For each of 100 or so iterations compare the new model to the true model you are simulating from.

Again, this is a good check on the reliability of the process -- the ability to find the 'true' features and to get good estimates of $\beta$.

When $Y$ is binary, instead of dealing with residuals, re-simulation involves computing the linear predictor $X\beta$ from the original fit (e.g., using the lasso), taking the logistic transformation, and generating for each Monte Carlo simulation a new $Y$ vector to fit afresh. In R one can say for example

lp <- predict(...) # assuming suitable predict method available, or fitted()
probs <- plogis(lp)
y <- ifelse(runif(n) <= probs, 1, 0)