Cross-validated prediction error worse for LARS than ordinary linear regression

Question

I am analysing microarray data in order to build a model for predicting cell proliferation (a continuous variable) based on gene expression (also a continuous variable). There are many more genes than samples (p>n), so I can't use ordinary linear regression methods for variable selection.

For this reason I tried to use least angle regression (LARS), using the lars package in R. I performed cross-validation to determine the number of steps at which the cross-validated mean squared error was lowest, and used this number of steps to build the model.

Out of curiosity, I also built a simple linear regression model using the single most correlated gene. To my surprise, the cross-validated mean squared error for this model was lower than that of the LARS model with more predictors.

I then tried making a LARS model with a single predictor, which was the same predictor that was most correlated and used in the linear regression model. This was an attempt to test that I was building the LARS model correctly, thinking that I should get the same model for both methods using the same predictor, as my understanding was that LARS was a way to select variables for a regression model in cases of p>n, but that the resulting coefficients would be the same for the same set of predictors. However, this was not the case, and the LARS model with the single predictor performed worse in cross-validation than the linear regression model. After checking my code, I don't think that I am building the LARS model incorrectly, I think that the simple linear regression model and the LARS model genuinely give different results with the same predictor.

In all cases leave-one-out cross-validation was performed. I chose this method as I only have 7 samples.

My questions are:

• Why did the LARS model perform worse than the simple linear regression model with a single predictor? To me this suggests overfitting in the LARS model, but performing the cross-validation of the LARS model in order to select the parameter for the number of steps should have avoided this - if fewer predictors in this model had been better, surely the cross-validated error would have been lowest at a lower number of steps?
• Why was the LARS model with a single predictor not the same (worse) as the simple linear regression model with a single predictor?

I am fairly inexperienced in statistics, so easy-to-understand explanations would be appreciated!

Answer 1

The key here is to remember that shrinkage methods like LARS, LASSO, and ridge regression are intended to protect you from over-optimism. If you want your results to generalize, not just fit your present samples, you do not want to rely on regression coefficients that might be highly dependent on the particular samples at hand or the randomness particular to these sets of analyses.

With only 7 samples and presumably 15,000-20,000 genes covered in your microarray, this is a substantial issue. With 7 samples, there are only $7!=5040$ different rank orders of sets of results among the samples. So you should not be too surprised if 3 or 4 individual gene-expression values have perfect rank correlations with your outcome variable of cell proliferation, just by chance.

The "single most correlated gene" will have an over-optimistically strong relation to your outcome variable; you chose it as the one most closely related in this data set. Remember that the word "regression" was originally used in the context of "regression to the mean"; the vagaries of sampling mean that the variable most highly correlated to outcome in this data set is not likely to be so in subsequent samples. You could think of your "single most correlated gene" as being overfit in this sample relative to its actual importance.

To reduce this over-optimism, shrinkage methods set regression coefficients to be lower in magnitude than they would be in an ordinary regression. (You can check that for yourself, from your results.) Thus you should not be surprised that cross-validation errors on your data set are larger with your LARS model than with the single-variable ordinary regression. That's the trade-off between fitting your particular set of samples and producing a model that has some hope of generalizing to other samples.

My guess is that the presence of only 1 predictor variable posed problems for the termination step of the LARS algorithm, so that it never reached the ultimate ordinary least squares coefficient when you tried that. I'm not very familiar with the algorithm, but as you are using R the source code is readily available for inspection.