Model stability when dealing with large $p$, small $n$ problem

Question

Intro:

I have a dataset with a classical "large p, small n problem". The number available samples n=150 while the number of possible predictors p=400. The outcome is a continuous variable.

I want to find the most "important" descriptors, i.e., those that are best candidates for explaining the outcome and helping to build a theory.

After research on this topic I found LASSO and Elastic Net are commonly used for the case of large p, small n. Some of my predictors are highly correlated and I want to preserve their groupings in the importance assessment, therefore, I opted for Elastic Net. I suppose that I can use absolute values of regression coefficients as a measure of importance (please correct me if I am wrong; my dataset is standardized).

Problem:

As my number of samples is small, how can I achieve a stable model?

My current approach is to find best tuning parameters (lambda and alpha) in a grid search on 90% of the dataset with 10-fold cross-validation averaging MSE score. Then I train the model with the best tuning parameters on the whole 90% of dataset. I am able to evaluate my model using R squared on the holdout 10% of the dataset (which account to only 15 samples).

Running repeatedly this procedure, I found a large variance in R squared assessments. As well, the number of non-zeroed predictors varies as well as their coefficients.

How can I get a more stable assessment of predictors' importance and more stable assessment of final model performance?

Can I repeatedly run my procedure to create a number of models, and then average regression coefficients? Or should I use the number of occurrences of a predictor in the models as its importance score?

Currently, I get around 40-50 non-zeroed predictors. Should I penalize number of predictors harder for better stability?

Answer 1

"Sparse Algorithms are not Stable: A No-free-lunch Theorem"

I guess the title says a lot, as you pointed out.

[...] a sparse algorithm can have non-unique optimal solutions, and is therefore ill-posed

Check out randomized lasso, and the talk by Peter Buhlmann.

Update:

I found this paper easier to follow than the paper by Meinshausen and Buhlmann called "Stability Selection".

In "Random Lasso", the authors consider the two important drawbacks of the lasso for large $p$, small $n$ problems, that is,

1. In the case where there exist several correlated variables, lasso only picks one or a few, thus leading to the instability that you talk about
2. Lasso cannot select more variables than the sample size $n$ which is a problem for many models

The main idea for random lasso, that is able to deal with both drawbacks of lasso is the following

If several independent data sets were generated from the same distribution, then we would expect lasso to select nonidentical subsets of those highly correlated important variables from different data sets, and our final collection may be most, or perhaps even all, of those highly correlated important variables by taking a union of selected variables from different data sets. Such a process may yield more than $n$ variables, overcoming the other limitation of lasso.

Bootstrap samples are drawn to simulate multiple data sets. The final coefficients are obtained by averaging over the results of each bootstrap sample.

It would be great if somebody could elaborate on and explain this algorithm further in the answers.

Answer 2

My current approach is to find best tuning parameters (lambda and alpha) in a grid search on 90% of the dataset with 10-fold cross-validation averaging MSE score. Then I train the model with the best tuning parameters on the whole 90% of dataset. I am able to evaluate my model using R squared on the holdout 10% of the dataset (which account to only 15 samples).

How stable are the tuning parameters?

Do you see large differences between goodness-of-fit (e.g. MSE of the optimal parameter's cross validation) and the 10% independent test performance?

That would be a symptom of overfitting:

The problem with the grid search (and many other parameter optimization strategies) is that you basically assume a rather smooth behaviour of $MSE = f (grid parameters)$. But for small test sets the variance due to the small test set size (= 135 samples total in 10 c.v.-folds) can be larger than the actual differences of $MSE = f (grid parameters)$. In that case already the parameters are rather unstable.

Can I repeatedly run my procedure to create a number of models, and then average regression coefficients? Or should I use the number of occurrences of a predictor in the models as its importance score?

There are several possibilities to build such aggregated models:

• linear models can be averaged by averaging the coefficients
• more generally, you can predict a sample by each of the $m$ different models, and average the $m$ predictions (you could also derive an idea of the uncertainty looking at the distribution of the predictions).

Search terms would be "aggregated models", "bootstrap aggregating", "bagging".

Side thought: some types of data have expected and interpretable collinearity that can cause variable selectio to "jump" between more or less equal solutions.

Answer 3

There's no way out of it. As some said, models are unstable by nature (otherwise statistics would not be needed).

But instability itself brings information. So instead of trying to get rid of it I tried to analyze it.

I run cross validation simulations many times and then get the coefficients for the best selected parameters in each run and put them together.

In the case of the elastic net I run a cross validation test for each alpha (0..1 by 0.1) with the same k folded data (you should compare alphas on the same data set) and select the $\lambda$/$\alpha$ pair associated with less test error... Than I repeat it for n times with different randomly selected k folded data and pick the best pair for each iteration.

Then I extract the regression coefficients for each parameters pair and this gives me a distribution of values for each parameter. This way I can use the mean/median value to describe the strength of the predictor and its standard deviation/IQR to describe its variability, that is its stability.

A very stabile predictor means you can expect its effect to be similar also with new data; a predictor which is unstable even in your data, would probably be very unstable even with new data.