Inconsistency Between GridSearchCV/Lasso and LassoCV MSE Curve & Optimal Alpha

### Describe the bug

GridSearchCV with Lasso provides a different MSE curve than LassoCV for small values of alpha (less than about 0.01 in the example below). For most random seeds, they provide the same optimal alpha, but for some seeds (like below), the inconsistency leads to different values.

### Steps/Code to Reproduce

```
import matplotlib.pyplot as plt
import numpy as np
from sklearn.linear_model import Lasso, LassoCV
from sklearn.model_selection import GridSearchCV
from sklearn.model_selection import ShuffleSplit

# Generate data
np.random.seed(13)
x = np.random.choice(2, 250, p=[0.8, 0.2]).reshape((-1,1))
A = np.random.normal(0,1,(130,250))
n = np.random.normal(0,5,(130,1))
y = A.dot(x) + n

# Create parameters
cv = ShuffleSplit(n_splits=10, random_state=0)
lasso_alphas=np.geomspace(start=0.000001, stop=10, num=20)

# Create LassoCV model
model_lasso = LassoCV(alphas=lasso_alphas, max_iter=100000, cv=cv, n_jobs=-1, random_state=0)
model_lasso.fit(A, y.ravel())
print(model_lasso.alpha_)

# Create GridSearchCV/Lasso model
param_gs_lasso = {'alpha':lasso_alphas}
model_gs_lasso = GridSearchCV(Lasso(max_iter=100000, random_state=0), param_gs_lasso, n_jobs=-1, scoring='neg_mean_squared_error', cv=cv)
model_gs_lasso.fit(A, y.ravel())
print(model_gs_lasso.best_estimator_.alpha)

# Compare MSE curves
plt.plot(model_lasso.alphas_, np.average(model_lasso.mse_path_, axis=1), label='LassoCV')
plt.plot(lasso_alphas, -model_gs_lasso.cv_results_['mean_test_score'], label='Lasso/GridSearchCV')
plt.semilogx()
plt.legend()
plt.show()
```

### Expected Results

I would expect the MSE curves and optimal alphas to be the same.

### Actual Results

![image](https://user-images.githubusercontent.com/6155964/200897972-8ceda8e9-faf9-4ba4-842b-9e7328aab087.png)


### Versions

```shell
System:
    python: 3.10.4 (tags/v3.10.4:9d38120, Mar 23 2022, 23:13:41) [MSC v.1929 64 bit (AMD64)]
executable: c:\Program Files\Python310\python.exe
   machine: Windows-10-10.0.19045-SP0

Python dependencies:
      sklearn: 1.1.2
          pip: 22.3
   setuptools: 65.5.0
        numpy: 1.23.4
        scipy: 1.9.3
       Cython: None
       pandas: 1.5.1
   matplotlib: 3.6.0
       joblib: 1.2.0
threadpoolctl: 3.1.0

Built with OpenMP: True

threadpoolctl info:
       user_api: blas
   internal_api: openblas
         prefix: libopenblas
       filepath: C:\Program Files\Python310\Lib\site-packages\numpy\.libs\libopenblas.FB5AE2TYXYH2IJRDKGDGQ3XBKLKTF43H.gfortran-win_amd64.dll
        version: 0.3.20
threading_layer: pthreads
   architecture: SkylakeX
    num_threads: 16

       user_api: blas
   internal_api: openblas
         prefix: libopenblas
       filepath: C:\Program Files\Python310\Lib\site-packages\scipy.libs\libopenblas-57db09cfe174768fb409a6bb5a530d4c.dll
        version: 0.3.18
threading_layer: pthreads
   architecture: SkylakeX
    num_threads: 16

       user_api: openmp
   internal_api: openmp
         prefix: vcomp
       filepath: C:\Program Files\Python310\Lib\site-packages\sklearn\.libs\vcomp140.dll
        version: None
    num_threads: 16
```

COMMENT:
Looking at the code, I don't think that we can expect the same results.

In `LassoCV`, we will compute the model for all `alphas` with the outer loop on the cv iterations. For `Lasso` + `GridSearchCV`, we will compute for one `alpha` all the cv iterations. Therefore, we don't consume the `random_state` in the same order so the results will be different.

Your problem is ill-defined because `n_samples < n_features` (this is why you are probably interested in `Lasso`). `alpha=1e-6` imposes not enough regularization that could be translated to a high variance in MSE and a model that overfit.

What is not clear to me is why the `LassoCV` is stable. It looks like it does not really search for those small `alphas` stopping to a higher `alpha` in the path that leads to a more stable model.

I assume that it is probably on purpose but I don't know enough about the code of the coordinate descent to know what is implementation trick there (if my analysis is not wrong).

@agramfort could you use your wisdom and bring some light regarding the internal of the `LassoCV` that makes it works without overfitting for low alpha values?

Labels: Documentation