Cross validation after LASSO in complex survey data

Question

I am trying to do model selection on some candidate predictors using LASSO with a continuous outcome. The goal is to select the optimal model with the best prediction performance, which usually can be done by K-fold cross validation after obtaining a solution path of the tuning parameters from LASSO. The issue here is that the data are from a complex multi-stage survey design (NHANES), with cluster sampling and stratification. The estimation part is not hard since glmnet in R can take sampling weights. But the cross validation part is less clear to me since observations now are not i.i.d anymore, and how can the procedure account for sampling weights representing a finite population?

So my questions are:

1) How to carry out K-fold cross validation with complex survey data to select the optimal tuning parameter? More specifically, how to appropriately partition the sample data into training and validation sets? And how to define the estimate of prediction error?

2) Is there an alternative way to select the optimal tuning parameter?

Answer 1

A recent paper demonstrates how K-fold cross-validation can be applied to complex designs (featuring clustering, stratification, and unequal selection probabilities). The paper is open-access, fairly short, and quite readable, so check it out!

Wieczorek, Guerin, and McMahon 2022. "K-fold cross-validation for complex sample surveys." https://doi.org/10.1002/sta4.454

The following excerpt provides a quick summary of their recommendations:

Nonetheless, Wolter's advice for taking subsamples also applies to forming CV folds: For SRS sampling, partition the observations in the dataset completely at random into K equal-sized folds, as usual.

For cluster sampling, partition the data at the level of the PSUs. All elements from a given PSU should be placed in the same fold, so that the folds are a random partition of PSUs rather than of elementary sampling units. (Note: with multistage sampling, after a first-stage cluster sample there is further subsampling in the selected clusters. Since there is no straightforward way for CV to mimic this subsampling, we simply form folds at the PSU level even in multistage samples.)

For stratified sampling, make each fold a stratified sample of units from each stratum. Create SRS CV folds separately within each stratum, then combine them across strata.

For unequal probability samples, such as probability proportional to size (PPS), make each fold an SRS of the data. However, see below for advice on using sampling weights to estimate $\hat{𝐸𝑟𝑟}_{𝐶𝑉}(𝑓)$ as a weighted mean and to fit each $\hat{𝑓}_{𝑡𝑟𝑎𝑖𝑛_𝑗}$ .

For more complex sample designs that combine several of the design features above, combine these rules as needed.

The authors of that paper created an R package, named 'surveyCV', which implements the K-fold cross-validation method they proposed. That R package is available on CRAN.

https://github.com/ColbyStatSvyRsch/surveyCV

Answer 2

I don't have a detailed answer, just some pointers to work I've been meaning to read:

You could take a look at McConville (2011) on complex-survey LASSO, to be sure your use of LASSO is appropriate for your data. But maybe it's not a big deal if you're doing LASSO only for variable selection, then fitting something else to the remaining variables.

For cross-validation with complex survey data (though not LASSO), McConville also cites Opsomer & Miller (2005) and You (2009). But their methods seems to use leave-one-out CV, not K-fold.

Leave-one-out should be simpler to implement with complex surveys---there's less concern about how to partition the data appropriately. (On the other hand, it can take longer to run than K-fold. And if your goal is model selection, it's known that leave-one-out can be worse than K-fold for large samples.)

Answer 3

EDIT by OP: Not applicable to complex survey data.

The cv.glmet function could help you to perform the cross validation required. The lambda.min value is the value of λ where the CV error is minimal. The lambda.1se represents the value of λ in the search that was simpler than the best model (lambda.min), but which has error within 1 standard error of the best model.

1. Choose a grid of values you could choose from for both alpha and lambda

grid <- expand.grid(.alpha = (1:10) * 0.1, .lambda = (1:10) * 0.1)

2. Setup the control parameters of your model. The below train control does repeatedcv for 10 iterations. Go over the methods available and pick the one that would fit your current scenario.

cv.glmmod <-cv.glmnet(xTrain,y=yTrain,alpha=grid$.alpha, > = T,lambda = grid$.lambda)

The lambda.min value could be accessed from the model itself as shown below.

cv.glmmod$lambda.min