0
$\begingroup$

I've created this toy example to demonstrate something that is occurring with my real data. I do not understand how to interpret the apparent failure of cv.glmnet to find a solution when it's explicitly "doped" into the data.

Let's create a matrix of random predictors:

xtrain = replicate(100,rnorm(10))

Now let's create a response variable that is explicitly a function of some of the x:

ytrain = 1 + 2*x[,1] - 3*x[,2]

Run cv.glmnet:

cvfit=cv.glmnet(x,y)

If we query the cvfit object, we'll find that lambda.min is not equal to lambda.1se, as expected, nor is lambda.min the largest magnitude lambda value. If we make predictions on the original training data

ypred=predict(cv.glmnet,xnew=xtrain,s="lambda.min")

we'll get quite good results in comparing ypred to ytrain.

HOWEVER, if we do the exact same thing with a much larger set of potential predictors.

xtrain=replicate(100000,rnorm(10))
ytrain=1 + 2*x[,1] - 3*x[,3]
cvfit=cv.glmnet(xtrain,ytrain)

we end up with a cvfit object where lambda.min = lambda.1se and this lambda value is the largest magnitude lambda value. The number of non-zero coefficients at that lambda is zero. And, of course, predictions at this lambda are nothing but the mean of ytrain.

I could use some help understanding this limitation in applying the cv.glmnet to such a problem and if there is a work-around of some sort (e.g. taking random subsets of the xtrain matrix).

$\endgroup$
1
$\begingroup$

I tried this experiment, and it seems like the lambda sequences are not the same. Store the first cv object like:

cv1 <- cv.glmnet(xtrain,ytrain)

and then, the second time around do:

cv.glmnet(xtrain,ytrain,lambda=cv1$lambda)

This should help, but I think additional messing with the lambda sequence is required to get results closer to those in the first case.

Now a question for you: if the y is an exact function of x, why cant it solve for it exactly? Doing so would make the error in the folds 0 and in the training as well (assuming identity link)... Thanks in advance!

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Thanks for the input. I'm not sure I understand your suggestion or your question. I created this example to try and understand why glmnet seemingly can't find a good solution in my real data (which has a dimension ~1,000,000). It produces a cvfit where lamda.min=lambda.1se and no non-zero coefficients. Essentially it fails to find any fit. Though if I take subsets of my data, it can find fits. So, I'm confused. $\endgroup$ – PickledZebra Apr 16 '14 at 2:11
  • $\begingroup$ cv.glmnet chooses a sequence of lambdas, and the choice it makes seems to depend on the dimensions of the X matrix. If the lambdas chosen dont cover enough of the "interesting" interval of lambda (the places where you get predictive betas) then you will get boring betas. Hence, I suggested using the lambda sequence from the first call to cvglmnet as an input to the second call. $\endgroup$ – Bruno Apr 16 '14 at 3:23
  • $\begingroup$ Yes, I see. Thank you. I understand how that would improve the situation if my real work were like the example I proposed. It seems, in general, that glmnet cannot guarantee the globally optimal solution. It is quite fast on smaller datasets, but as the size of the dataset grows, the solutions returned get increasingly "boring" (the trivial case of no betas and only an intercept). $\endgroup$ – PickledZebra Apr 16 '14 at 4:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.