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

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?

Any help will be greatly appreciated!

-
Maybe regularized regression is just not enough? Have you tried some machine learning approach? – mbq May 31 '12 at 19:03
which ML approaches suited for continuous or ordinal dependent variable would you suggest? – dimi May 31 '12 at 21:08
Random forest, SVR... For feature selection you can try some all-relevant method (should be more stable) like ACE or Boruta. – mbq May 31 '12 at 23:21

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.

-
Nice links(+1). – jbowman May 31 '12 at 21:13
Thank you for remarks. I have also considered randomized lasso, but would it suit in the case of collinearity? – dimi Jun 1 '12 at 9:28
What happens under approximate collinearity between the explanatory variables? In the conventional forward search algorithm in regression analysis, we are often faced with the situation where two variables x1 and x2 have similar explanatory power. If x1 is in the model, then there is no need to include x2; conversely, if x2 is in the model there is no need to include x1. If I understand your procedure correctly, you will tend to include x1 half the time and x2 half the time, leading to stability probabilities of about 50% each. If so, you might falsely conclude that neither variable is needed. – dimi Jun 1 '12 at 9:36
I imagine, I can similarly use randomized elastic net in this case. – dimi Jun 1 '12 at 9:37
show 1 more comment

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.

-