1
$\begingroup$

Package glmnet provides a cross validation function called cv.glmnet that allows us to choose between two suggested models (from the many), labelled "lambda.min" and "lambda.1se". However even "lambda.1se" does not provide enough sparsity for me, as can be seen in the example below:

Here is a 27x27 variance-covariance matrix, from real life, based on foreign exchange market returns (sample size 5 years), from which I have reconstructed some data using the MASS package's mvrnorm function.

rr <- matrix(c(+3.905427e-04,-1.672944e-04,-3.868682e-05,-2.012400e-05,-5.208216e-05,-1.676518e-05,+3.191132e-05,-8.393085e-05,+2.994233e-05,+4.142394e-05,
+2.788218e-05,+2.100455e-05,+8.915135e-06,+3.237941e-05,-1.370915e-05,-1.869792e-05,-2.818437e-05,+4.867242e-06,-1.484061e-06,-9.871850e-07,
-2.301505e-05,+1.239755e-05,+2.707801e-06,-2.730393e-05,-3.439632e-05,-2.850229e-06,-1.212873e-06,-1.672944e-04,+3.148921e-04,+3.981354e-05,
-2.159599e-05,-7.679438e-06,+7.327542e-07,+8.133157e-06,+1.944368e-05,-3.291582e-05,-2.830392e-05,-9.711673e-06,-1.953129e-05,-8.678863e-06,
-1.468040e-05,-2.999372e-05,-8.443994e-06,-3.243783e-06,+8.457205e-06,-3.331677e-05,-2.585784e-05,-9.963355e-06,-1.006820e-05,-1.510451e-05,
-2.297489e-05,-2.822061e-05,-1.058895e-05,-1.523963e-05,-3.868682e-05,+3.981354e-05,+6.136468e-05,+6.385946e-06,+5.402061e-06,+1.340335e-06,
-2.120207e-05,-1.004486e-05,-2.145629e-06,-2.930022e-06,-2.912909e-06,-1.334397e-06,-1.120393e-06,-3.215213e-06,-3.748124e-06,-5.266702e-07,
+1.385250e-06,+2.594281e-06,-2.793294e-06,-3.918135e-06,+3.701831e-06,-1.357614e-06,+1.050664e-06,+1.330178e-06,-1.893086e-06,-8.355664e-07,
-9.429525e-07,-2.012400e-05,-2.159599e-05,+6.385946e-06,+4.514971e-05,+1.487041e-05,+4.051007e-06,-2.579603e-05,-7.896501e-06,+3.735320e-06,
-8.991534e-07,-1.909784e-06,-2.652713e-07,-3.995973e-07,-1.947112e-06,+1.223827e-05,+2.273053e-06,+5.392409e-06,-4.957146e-07,+6.043096e-07,
+1.373694e-06,+7.145170e-06,+4.236828e-08,+2.855058e-06,+1.196767e-05,+6.420234e-06,+3.212209e-06,+5.828089e-06,-5.208216e-05,-7.679438e-06,
+5.402061e-06,+1.487041e-05,+3.968051e-05,+8.180625e-06,-1.679490e-05,+1.438802e-06,-9.333792e-07,-3.856856e-06,-4.196170e-06,-2.545341e-06,
-1.465675e-06,-4.679705e-06,+5.273817e-06,+4.392979e-06,+6.707402e-06,-1.777882e-06,+3.340004e-06,+3.041857e-06,+1.046557e-05,-1.228042e-06,
+3.006442e-06,+1.042182e-05,+6.513678e-06,+2.545036e-06,+3.286945e-06,-1.676518e-05,+7.327542e-07,+1.340335e-06,+4.051007e-06,+8.180625e-06,
+4.515972e-06,-2.384442e-06,+2.346723e-06,-6.908337e-07,-1.428231e-06,-1.410132e-06,-1.144762e-06,-4.898731e-07,-1.704683e-06,+1.394768e-07,
+7.696074e-07,+1.714821e-06,-4.269990e-07,+3.351277e-07,+5.123849e-07,+2.740856e-06,-4.688568e-07,+3.553067e-07,+2.336701e-06,+5.732269e-07,
+3.052545e-07,+2.869335e-07,+3.191132e-05,+8.133157e-06,-2.120207e-05,-2.579603e-05,-1.679490e-05,-2.384442e-06,+1.496400e-04,+2.803930e-05,
-7.596035e-06,-2.577285e-06,-1.057582e-06,-1.404342e-05,-5.620170e-06,-1.767036e-06,-3.168803e-05,-1.800529e-05,-1.292217e-05,-4.453331e-06,
-3.319978e-05,-1.786611e-05,-1.284912e-05,+4.276601e-07,-1.347161e-05,-3.513421e-05,-2.668011e-05,-1.162202e-05,-1.245654e-05,-8.393085e-05,
+1.944368e-05,-1.004486e-05,-7.896501e-06,+1.438802e-06,+2.346723e-06,+2.803930e-05,+1.250962e-04,-8.983686e-06,-1.080656e-05,-5.574779e-06,
-1.244544e-05,-5.308212e-06,-9.088221e-06,-1.134928e-05,-7.996861e-06,-1.762774e-06,-8.772736e-06,-6.480971e-06,-9.015757e-06,-9.206922e-06,
-5.871784e-06,-5.546947e-06,-9.127687e-06,-1.096436e-05,-5.814707e-06,-7.572505e-06,+2.994233e-05,-3.291582e-05,-2.145629e-06,+3.735320e-06,
-9.333792e-07,-6.908337e-07,-7.596035e-06,-8.983686e-06,+1.254965e-05,+7.247157e-06,+1.709803e-06,+2.861585e-06,+1.100077e-06,+2.904772e-06,
+4.547653e-06,+2.510788e-07,-6.168826e-07,-3.769066e-07,+1.891723e-06,+2.750163e-06,-8.707644e-07,+1.931014e-06,+1.641571e-06,+2.066175e-06,
+7.301858e-07,+1.004776e-06,+2.087421e-06,+4.142394e-05,-2.830392e-05,-2.930022e-06,-8.991534e-07,-3.856856e-06,-1.428231e-06,-2.577285e-06,
-1.080656e-05,+7.247157e-06,+1.374557e-05,+2.806961e-06,+3.332886e-06,+1.514102e-06,+3.783409e-06,-1.389411e-06,-1.112583e-06,-2.430049e-06,
+2.748058e-07,+5.139875e-07,+1.300466e-06,-2.462923e-06,+1.697307e-06,+1.043461e-06,-7.175271e-07,-2.794225e-06,-1.512459e-07,+2.815655e-07,
+2.788218e-05,-9.711673e-06,-2.912909e-06,-1.909784e-06,-4.196170e-06,-1.410132e-06,-1.057582e-06,-5.574779e-06,+1.709803e-06,+2.806961e-06,
+5.153962e-06,+3.137140e-06,+1.442756e-06,+2.461616e-06,-1.514675e-06,-2.032794e-06,-1.900112e-06,+3.598444e-07,+1.077534e-06,-6.946718e-07,
-2.577196e-06,+2.510696e-07,+4.922748e-07,-1.789320e-06,-3.410673e-06,-1.516124e-07,-5.658895e-07,+2.100455e-05,-1.953129e-05,-1.334397e-06,
-2.652713e-07,-2.545341e-06,-1.144762e-06,-1.404342e-05,-1.244544e-05,+2.861585e-06,+3.332886e-06,+3.137140e-06,+1.322267e-05,+3.742704e-06,
+2.897316e-06,+2.190498e-06,+2.205173e-06,+1.143671e-06,+6.378957e-07,+5.430563e-06,+1.286320e-06,+2.526458e-07,+2.177544e-07,+2.828000e-06,
+4.718304e-06,+1.089443e-06,+1.589336e-06,+9.667085e-07,+8.915135e-06,-8.678863e-06,-1.120393e-06,-3.995973e-07,-1.465675e-06,-4.898731e-07,
-5.620170e-06,-5.308212e-06,+1.100077e-06,+1.514102e-06,+1.442756e-06,+3.742704e-06,+3.252501e-06,+1.140401e-06,+9.556437e-07,+1.074267e-06,
+6.025770e-07,+2.229803e-07,+2.229262e-06,+7.350322e-07,-1.697000e-07,-3.558592e-08,+1.125714e-06,+1.675207e-06,+4.637744e-07,+7.091632e-07,
+1.517128e-07,+3.237941e-05,-1.468040e-05,-3.215213e-06,-1.947112e-06,-4.679705e-06,-1.704683e-06,-1.767036e-06,-9.088221e-06,+2.904772e-06,
+3.783409e-06,+2.461616e-06,+2.897316e-06,+1.140401e-06,+4.663585e-06,-6.658251e-07,-6.095897e-07,-2.046872e-06,+6.571694e-07,+7.337406e-07,
+4.691496e-07,-2.068054e-06,+1.187547e-06,+5.255661e-07,-1.547152e-06,-2.647444e-06,+1.106522e-07,+4.997985e-08,-1.370915e-05,-2.999372e-05,
-3.748124e-06,+1.223827e-05,+5.273817e-06,+1.394768e-07,-3.168803e-05,-1.134928e-05,+4.547653e-06,-1.389411e-06,-1.514675e-06,+2.190498e-06,
+9.556437e-07,-6.658251e-07,+5.056387e-05,+5.733135e-06,+4.876198e-06,-7.288427e-07,+8.635461e-06,+6.028157e-06,-1.714284e-06,+1.581347e-08,
+3.218589e-06,+1.109427e-05,+1.044090e-05,+4.586428e-06,+6.845056e-06,-1.869792e-05,-8.443994e-06,-5.266702e-07,+2.273053e-06,+4.392979e-06,
+7.696074e-07,-1.800529e-05,-7.996861e-06,+2.510788e-07,-1.112583e-06,-2.032794e-06,+2.205173e-06,+1.074267e-06,-6.095897e-07,+5.733135e-06,
+2.044822e-05,+5.782109e-06,+3.019946e-07,+6.052292e-06,+4.901258e-06,+3.139930e-06,-7.229161e-08,+2.063215e-06,+9.119608e-06,+5.919809e-06,
+2.508532e-06,+1.952830e-06,-2.818437e-05,-3.243783e-06,+1.385250e-06,+5.392409e-06,+6.707402e-06,+1.714821e-06,-1.292217e-05,-1.762774e-06,
-6.168826e-07,-2.430049e-06,-1.900112e-06,+1.143671e-06,+6.025770e-07,-2.046872e-06,+4.876198e-06,+5.782109e-06,+1.050558e-05,-4.464648e-07,
+3.715678e-06,+2.971496e-06,+4.662022e-06,-2.933217e-07,+1.929764e-06,+8.453756e-06,+6.938162e-06,+2.054123e-06,+2.107770e-06,+4.867242e-06,
+8.457205e-06,+2.594281e-06,-4.957146e-07,-1.777882e-06,-4.269990e-07,-4.453331e-06,-8.772736e-06,-3.769066e-07,+2.748058e-07,+3.598444e-07,
+6.378957e-07,+2.229803e-07,+6.571694e-07,-7.288427e-07,+3.019946e-07,-4.464648e-07,+6.154732e-06,-4.826199e-07,+7.972120e-08,-1.076161e-08,
+3.200496e-07,+5.782623e-08,-2.955595e-07,-2.241002e-06,-1.286211e-07,-1.497812e-08,-1.484061e-06,-3.331677e-05,-2.793294e-06,+6.043096e-07,
+3.340004e-06,+3.351277e-07,-3.319978e-05,-6.480971e-06,+1.891723e-06,+5.139875e-07,+1.077534e-06,+5.430563e-06,+2.229262e-06,+7.337406e-07,
+8.635461e-06,+6.052292e-06,+3.715678e-06,-4.826199e-07,+4.450922e-05,+1.085514e-05,+5.123687e-06,+5.368289e-07,+7.298413e-06,+1.028015e-05,
+4.371002e-06,+3.204427e-06,+2.715323e-06,-9.871850e-07,-2.585784e-05,-3.918135e-06,+1.373694e-06,+3.041857e-06,+5.123849e-07,-1.786611e-05,
-9.015757e-06,+2.750163e-06,+1.300466e-06,-6.946718e-07,+1.286320e-06,+7.350322e-07,+4.691496e-07,+6.028157e-06,+4.901258e-06,+2.971496e-06,
+7.972120e-08,+1.085514e-05,+2.321785e-05,+6.428409e-06,+1.349267e-06,+5.085592e-06,+6.280345e-06,+5.900478e-06,+2.167341e-06,+2.780337e-06,
-2.301505e-05,-9.963355e-06,+3.701831e-06,+7.145170e-06,+1.046557e-05,+2.740856e-06,-1.284912e-05,-9.206922e-06,-8.707644e-07,-2.462923e-06,
-2.577196e-06,+2.526458e-07,-1.697000e-07,-2.068054e-06,-1.714284e-06,+3.139930e-06,+4.662022e-06,-1.076161e-08,+5.123687e-06,+6.428409e-06,
+2.909737e-05,+1.584565e-06,+5.282449e-06,+6.901083e-06,+5.037513e-06,+2.245936e-06,+2.913474e-06,+1.239755e-05,-1.006820e-05,-1.357614e-06,
+4.236828e-08,-1.228042e-06,-4.688568e-07,+4.276601e-07,-5.871784e-06,+1.931014e-06,+1.697307e-06,+2.510696e-07,+2.177544e-07,-3.558592e-08,
+1.187547e-06,+1.581347e-08,-7.229161e-08,-2.933217e-07,+3.200496e-07,+5.368289e-07,+1.349267e-06,+1.584565e-06,+1.202174e-05,+6.348587e-07,
+7.476531e-07,-3.430022e-07,-2.735203e-08,+1.029644e-06,+2.707801e-06,-1.510451e-05,+1.050664e-06,+2.855058e-06,+3.006442e-06,+3.553067e-07,
-1.347161e-05,-5.546947e-06,+1.641571e-06,+1.043461e-06,+4.922748e-07,+2.828000e-06,+1.125714e-06,+5.255661e-07,+3.218589e-06,+2.063215e-06,
+1.929764e-06,+5.782623e-08,+7.298413e-06,+5.085592e-06,+5.282449e-06,+6.348587e-07,+5.694378e-06,+4.881555e-06,+3.197262e-06,+1.468878e-06,
+1.591542e-06,-2.730393e-05,-2.297489e-05,+1.330178e-06,+1.196767e-05,+1.042182e-05,+2.336701e-06,-3.513421e-05,-9.127687e-06,+2.066175e-06,
-7.175271e-07,-1.789320e-06,+4.718304e-06,+1.675207e-06,-1.547152e-06,+1.109427e-05,+9.119608e-06,+8.453756e-06,-2.955595e-07,+1.028015e-05,
+6.280345e-06,+6.901083e-06,+7.476531e-07,+4.881555e-06,+3.021823e-05,+1.387415e-05,+5.421927e-06,+6.312873e-06,-3.439632e-05,-2.822061e-05,
-1.893086e-06,+6.420234e-06,+6.513678e-06,+5.732269e-07,-2.668011e-05,-1.096436e-05,+7.301858e-07,-2.794225e-06,-3.410673e-06,+1.089443e-06,
+4.637744e-07,-2.647444e-06,+1.044090e-05,+5.919809e-06,+6.938162e-06,-2.241002e-06,+4.371002e-06,+5.900478e-06,+5.037513e-06,-3.430022e-07,
+3.197262e-06,+1.387415e-05,+6.129212e-05,+4.943961e-06,+6.180768e-06,-2.850229e-06,-1.058895e-05,-8.355664e-07,+3.212209e-06,+2.545036e-06,
+3.052545e-07,-1.162202e-05,-5.814707e-06,+1.004776e-06,-1.512459e-07,-1.516124e-07,+1.589336e-06,+7.091632e-07,+1.106522e-07,+4.586428e-06,
+2.508532e-06,+2.054123e-06,-1.286211e-07,+3.204427e-06,+2.167341e-06,+2.245936e-06,-2.735203e-08,+1.468878e-06,+5.421927e-06,+4.943961e-06,
+6.807734e-06,+2.987841e-06,-1.212873e-06,-1.523963e-05,-9.429525e-07,+5.828089e-06,+3.286945e-06,+2.869335e-07,-1.245654e-05,-7.572505e-06,
+2.087421e-06,+2.815655e-07,-5.658895e-07,+9.667085e-07,+1.517128e-07,+4.997985e-08,+6.845056e-06,+1.952830e-06,+2.107770e-06,-1.497812e-08,
+2.715323e-06,+2.780337e-06,+2.913474e-06,+1.029644e-06,+1.591542e-06,+6.312873e-06,+6.180768e-06,+2.987841e-06,+9.474784e-06), 27, 27)
colnames(rr) <- rownames(rr) <- c("USD","EUR","GBP","CAD","AUD","NZD","JPY","XAU","NOK","SEK","CZK","PLN","HUF","RON","RUB","TRY","ZAR","ILS","INR","IDR","KRW","TWD","PHP",
"MXN","BRL","CLP","COP")
library("MASS")
set.seed(1234)
data = mvrnorm(260 * 5, rep(0, 27), rr)

Clearly principal component 1 is EUR vs the USD, so EURUSD. This will be easy to replicate in the market with a single transaction, namely, EURUSD:

barplot(eigen(cov(data))$vectors[, 1], names.arg = colnames(data), las = 2)
title("PC1")

enter image description here

PC2 however, which I interpret as "emerging markets + commodity currencies", is more interesting, and would need lots of transactions to replicate accurately:

barplot(eigen(cov(data))$vectors[, 2], names.arg = colnames(data), las = 2)
title("PC2")

enter image description here

So in an attempt to lower the transaction costs, I want a sparse fit for PC2. Enter glmnet:

library("glmnet")
cvmodel <- cv.glmnet(x = data, y = data %*% eigen(cov(data))$vectors[, 2])
barplot(coef(cvmodel, s = "lambda.1se")[-1], names.arg = colnames(data), las = 2, main = "cv.glmnet lambda.1se coefficients")

enter image description here

As you can see I really haven't obtained much sparsity, and if I plot actual PC2 vs my "sparse" PC2, I get an altogether too good fit.

plot(data %*% eigen(cov(data))$vectors[, 2], data %*% coef(cvmodel, s = "lambda.1se")[-1], main = "actual vs cv.glmnet")

enter image description here

So my question is, how do I get a "worse" fit, using much more sparse model, using glmnet, in a rigorous way? I know the fit provides a whole bunch of lambda values but I don't know how to choose a worse one that still suits my purposes. Could I go through each one, create an lm each time, and choose the first with an r-squared for example that's better than 0.8, say? What more efficient strategy can I use here?

Improve this question
$\endgroup$
1
  • 1
    $\begingroup$ See sparse principal components in R packages "PMA" (newer and recommended) and "elasticnet" (older), that will give you a more direct solution. You will find citations of research papers in the documentation of the main functions of these packages. $\endgroup$
    – Richard Hardy
    Commented Aug 9, 2016 at 20:23

2 Answers 2

Reset to default
3
$\begingroup$

The most rigorous way to get a more sparse model is to just increase the size of $\lambda$. That's all. Larger $\lambda$ will shrink the coefficients more towards zero, which causes more of them to be pinned at 0, which is precisely the sparsity property that you want. So pick the value of $\lambda$ that corresponds to how much sparseness you want; sparisty of each solution is reported by the nzero part of the returned cv.glmnet object.

Improve this answer
$\endgroup$
5
  • $\begingroup$ The reason for rotating FX market data is that 75% of variance is explained by 4 principal components (ee <- eigen(cov(data))$values; cumsum(ee)/sum(ee)). There is no point in trading specific currencies if they are heavily affected by the principal components. For good diversification, we want to separate out common risk and specific risk. Enter PCA. $\endgroup$
    – Thomas Browne
    Commented Aug 9, 2016 at 2:21
  • 1
    $\begingroup$ Yes I know that lambda must become greater, but how much greater is the question. Can I relate each lambda back to say, r.squared of the resultant fit, without doing 60 linear models? "lambda.1se" has already given us a hint at how to choose lambda. Is there a "lambda.2se"? "3se"? $\endgroup$
    – Thomas Browne
    Commented Aug 9, 2016 at 2:22
  • 1
    $\begingroup$ glmnet also reports how sparse each value of lambda is. Since you know how much sparsity you want, choose the lambda with the desired amount of sparsity. $\endgroup$
    – Sycorax ♦
    Commented Aug 9, 2016 at 2:24
  • $\begingroup$ Excellent, which variable is that? $\endgroup$
    – Thomas Browne
    Commented Aug 9, 2016 at 2:25
  • 1
    $\begingroup$ It's nzero in the returned cv.glmnet object. This is discussed in the documentation. $\endgroup$
    – Sycorax ♦
    Commented Aug 9, 2016 at 2:27
1
$\begingroup$

Choosing a $\lambda$ value smaller than lambda.min or lambda.1se likely penalizes the large coefficients too heavily. Often, the $\lambda$ value optimal for selection is not optimal for prediction. The so-called 'relaxed lasso' (Meinshausen, 2007) reduces this issue; it 'debiases' the lasso solution by unshrinking coefficients of selected predictors towards their unpenalized (e.g., OLS) values.

how do I get a "worse" fit, using much more sparse model, using glmnet [...] I know the fit provides a whole bunch of lambda values [...]. Could I go through each one, create an lm each time [...] What more efficient strategy can I use here?

By specifying relax = TRUE in the call to cv.glmnet. See also Hastie, Tibshirani and Tibshirani (2020), or the relaxed-lasso vignette (https://cran.r-project.org/package=glmnet/vignettes/relax.pdf).

Can I relate each lambda back to say, r.squared of the resultant fit, without doing 60 linear models?

One could first standardize the response so it has a variance of 1. That way, the cross-validated MSE can be interpreted as a coefficient of determination, or 1-$R^2$.

For your data, one could do:

scaled_y <- scale(data %*% eigen(cov(data))$vectors[, 2])
library("glmnet")
set.seed(42)
cvmodel <- cv.glmnet(x = data, y = scaled_y, relax = TRUE)
plot(cvmodel)

This yields the following:

Plot of mean-squared-error versus log(lambda)

Note that $\gamma$ is the mixing parameter, giving the weight of the penalized fit. The original lasso fit is obtained with $\gamma = 1$, while $\gamma = 0$ returns unbiased OLS coefficients, after selecting variables with the lasso. It appears $\gamma = 0$ is always most beneficial for this data problem.

One can print numerical results to obtain a more precise result:

head(data.frame(cvmodel$relaxed$statlist$`g:0`), 12)

       lambda        cvm        cvsd       cvup       cvlo nzero
s0  0.7189748 0.96510044 0.028974384 0.99407483 0.93612606     0
s1  0.6551030 0.33704936 0.009816319 0.34686568 0.32723304     2
s2  0.5969055 0.33526067 0.010529985 0.34579066 0.32473069     2
s3  0.5438780 0.33526067 0.010529994 0.34579067 0.32473068     2
s4  0.4955614 0.33512809 0.010494051 0.34562214 0.32463404     2
s5  0.4515371 0.30674833 0.013334435 0.32008277 0.29341390     2
s6  0.4114238 0.22634878 0.009249976 0.23559876 0.21709881     4
s7  0.3748740 0.22035701 0.009661487 0.23001850 0.21069552     4
s8  0.3415713 0.12800604 0.007542284 0.13554833 0.12046376     6
s9  0.3112270 0.08427203 0.003880242 0.08815227 0.08039179     7
s10 0.2835785 0.08083322 0.003713613 0.08454684 0.07711961     8
s11 0.2583861 0.07755071 0.004199024 0.08174973 0.07335168     8

To obtain an $R^2$ of about 0.70 (1 - 0.30), one could employ $\lambda$ s5.

Meinshausen, N. (2007). Relaxed lasso. Computational Statistics & Data Analysis, 52(1), 374-393. https://doi.org/10.1016/j.csda.2006.12.019

Hastie, T., Tibshirani, R., & Tibshirani, R. (2020). Best Subset, Forward Stepwise or Lasso? Analysis and Recommendations Based on Extensive Comparisons. Statistical Science, 35(4), 579-592. https://doi.org/10.1214/19-STS733

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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