So I decided to run nested cross-validation using the specialized `mlr` package in R to see what's actually coming from the modelling approach. #Code (it takes a few minutes to run on an ordinary notebook) library(mlr) daf = read.csv("https://pastebin.com/raw/p1cCCYBR", sep = " ", header = FALSE) tsk = list( tsk1110 = makeRegrTask(id = "tsk1110", data = daf, target = colnames(daf)[1]), tsk500 = makeRegrTask(id = "tsk500", data = daf[, c(1,sample(ncol(daf)-1, 500)+1)], target = colnames(daf)[1]), tsk100 = makeRegrTask(id = "tsk100", data = daf[, c(1,sample(ncol(daf)-1, 100)+1)], target = colnames(daf)[1]), tsk50 = makeRegrTask(id = "tsk50", data = daf[, c(1,sample(ncol(daf)-1, 50)+1)], target = colnames(daf)[1]), tsk10 = makeRegrTask(id = "tsk10", data = daf[, c(1,sample(ncol(daf)-1, 10)+1)], target = colnames(daf)[1]) ) lrn = makeLearner("regr.glmnet", alpha = 0) rdesc = makeResampleDesc("CV", iters = 10) msrs = list(mse, rsq) ctrl = makeTuneControlGrid(resolution = 15L) pset = makeParamSet( makeNumericParam("lambda", lower = -10, upper = 10, trafo = function(x) exp(x))) cvl = makeTuneWrapper(learner = lrn, resampling = rdesc, measures = msrs, control = ctrl, par.set = pset) bm = benchmark(learners = cvl, tasks = tsk, resamplings = rdesc, measures = msrs) #Results > > bm > task.id learner.id mse.test.mean rsq.test.mean > 1 tsk10 regr.glmnet.tuned 0.9029739 -0.34385207 > 2 tsk100 regr.glmnet.tuned 0.9860033 -0.20159950 > 3 tsk1110 regr.glmnet.tuned 0.6836474 0.08630208 > 4 tsk50 regr.glmnet.tuned 0.9159190 -0.15177725 > 5 tsk500 regr.glmnet.tuned 0.8874529 0.04095012 They basically got the same slightly better than random performance, but `glmnet` is good enough to squeeze some information the more features you feed it. So, what about the optimal lambdas? getBMRTuneResults(bm, as.df = TRUE, task.ids = "tsk1110") # task.id learner.id iter lambda mse.test.mean rsq.test.mean #1 tsk1110 regr.glmnet.tuned 1 0.239651 0.7234743 0.06932993 #2 tsk1110 regr.glmnet.tuned 2 0.239651 0.6318552 0.10130075 #3 tsk1110 regr.glmnet.tuned 3 0.239651 0.7202225 0.09277869 #4 tsk1110 regr.glmnet.tuned 4 4.172734 0.7169810 -0.08232901 #5 tsk1110 regr.glmnet.tuned 5 0.239651 0.7163615 -0.05643168 #6 tsk1110 regr.glmnet.tuned 6 0.239651 0.6655465 0.13894617 #7 tsk1110 regr.glmnet.tuned 7 0.239651 0.7125394 0.14249063 #8 tsk1110 regr.glmnet.tuned 8 0.239651 0.7312063 -0.56975408 #9 tsk1110 regr.glmnet.tuned 9 0.239651 0.6895064 -0.03308516 #10 tsk1110 regr.glmnet.tuned 10 0.239651 0.6370360 0.02805295 Notice the lambdas are already transformed. Not a single fold picked the minimal lambda $\exp(-10)\approx 4.54 \mathbf E-5$ (this is the same behavior with `cv.glmnet` in R). I fiddled a bit more with `glmnet` and discovered neither there the minimal lambda is picked. Check: cvfit = cv.glmnet(x = x, y = y, alpha = 0, lambda = exp(seq(-10, 10, length.out = 150))) plot(cvfit) [![enter image description here][1]][1] #Conclusion So, basically, $\lambda>0$ really improves the fit. >How is it possible and what does it say about my dataset? Am I missing something obvious or is it indeed counter-intuitive? We are likely nearer the true distribution of the data setting $\lambda$ to a small value larger than zero. There's nothing counter-intuitive about it though. EDIT: The lambda selection does say something more about your data. As larger lambdas decrease performance, it means there are preferential, *i.e.* larger, coefficients in your model, as large lambdas shrink all coefficients towards zero. Though $\lambda\neq0$ means that the effective degrees of freedom in your model is smaller than the apparent degrees of freedom, $p$. >How can there be any qualitative difference between p=100 and p=1000 given that both are larger than n? $p=1000$ invariably contains at least the same of information or even more than $p=100$. [1]: https://i.sstatic.net/mZxU5.png