Variabilidade nos resultados do cv.glmnet

Estou usando cv.glmnetpara encontrar preditores. A configuração que eu uso é a seguinte:

lassoResults<-cv.glmnet(x=countDiffs,y=responseDiffs,alpha=1,nfolds=cvfold)
bestlambda<-lassoResults$lambda.min

results<-predict(lassoResults,s=bestlambda,type="coefficients")

choicePred<-rownames(results)[which(results !=0)]

Garantir que os resultados sejam reprodutíveis set.seed(1). Os resultados são altamente variáveis. Corri exatamente o mesmo código 100 para ver como os resultados eram variáveis. Nas corridas 98/100, sempre havia um preditor em particular selecionado (às vezes por conta própria); outros preditores foram selecionados (o coeficiente era diferente de zero) geralmente 50/100 vezes.

Por isso, diz-me que cada vez que a validação cruzada estiver em execução, provavelmente irá selecionar um melhor lambda diferente, porque a randomização inicial das dobras é importante. Outros viram esse problema ( resultados do CV.glmnet ), mas não há uma solução sugerida.

Eu estou pensando que talvez aquele que aparece 98/100 provavelmente esteja bastante correlacionado com todos os outros? Os resultados não estabilizar, se eu LOOCV apenas correr ( ), mas estou curioso por que eles são tão variáveis quando \ text {nfold} <n .$\text{fold-size} = n$$\text{nfold} < n$

O ponto aqui é que nas cv.glmnetdobras K ("partes") são escolhidas aleatoriamente.

Na validação cruzada de dobras K, o conjunto de dados é dividido em partes partes são usadas para prever a parte K-ésima (isso é feito vezes, usando uma parte diferente a cada vez). Isso é feito para todas as lambdas e é a que fornece o menor erro de validação cruzada.K - 1 K $K$$K-1$$K$$K$lambda.min

É por isso que, quando você usa os resultados não mudam: cada grupo é composto de um, portanto, não há muita escolha para os grupos$nfolds = n$$K$

No cv.glmnet()manual de referência:

Observe também que os resultados do cv.glmnet são aleatórios, pois as dobras são selecionadas aleatoriamente. Os usuários podem reduzir essa aleatoriedade executando o cv.glmnet várias vezes e calculando a média das curvas de erro.

### cycle for doing 100 cross validations
### and take the average of the mean error curves
### initialize vector for final data.frame with Mean Standard Errors
MSEs <- NULL
for (i in 1:100){
                 cv <- cv.glmnet(y, x, alpha=alpha, nfolds=k)  
                 MSEs <- cbind(MSEs, cv$cvm)
             }
  rownames(MSEs) <- cv$lambda
  lambda.min <- as.numeric(names(which.min(rowMeans(MSEs))))

MSEs é o quadro de dados que contém todos os erros para todas as lambdas (para as 100 execuções), lambda.miné sua lambda com erro médio mínimo.

Ultimamente, enfrentei o mesmo problema. Tentei repetir o CV várias vezes, como 100, 200, 1000 no meu conjunto de dados, tentando encontrar o melhor e (estou usando uma rede elástica). Mas, mesmo se eu criar 3 testes cv, cada um com 1000 iterações que calculam a média dos mínimos MSEs para cada , recebo três casais diferentes ( , ).α α λ $\lambda$$\alpha$$\alpha$$\lambda$$\alpha$

Não tocarei no problema aqui, mas decidi que minha melhor solução não é calcular a média dos min MSEs, mas extrair os coeficientes para cada iteração best e depois tratá-los como uma distribuição de valores (uma variável aleatória).$\alpha$$\lambda$

Então, para cada preditor, recebo:

• coeficiente médio
• desvio padrão
• Resumo de 5 números (mediana, quartis, mínimo e máximo)
• A porcentagem de vezes é diferente de zero (ou seja, tem influência)

Dessa forma, recebo uma descrição bastante sólida do efeito do preditor. Depois de ter distribuições para os coeficientes, você pode executar qualquer coisa estatística que considere valer para obter IC, valores de p, etc ... mas ainda não investiguei isso.

Este método pode ser usado com mais ou menos qualquer método de seleção que eu possa pensar.

Vou adicionar outra solução, que lida com o bug no @ Alice devido à falta de lambdas, mas não requer pacotes extras como o @Max Ghenis. Agradecemos a todas as outras respostas - todo mundo faz pontos úteis!

lambdas = NULL
for (i in 1:n)
{
    fit <- cv.glmnet(xs,ys)
    errors = data.frame(fit$lambda,fit$cvm)
    lambdas <- rbind(lambdas,errors)
}
# take mean cvm for each lambda
lambdas <- aggregate(lambdas[, 2], list(lambdas$fit.lambda), mean)

# select the best one
bestindex = which(lambdas[2]==min(lambdas[2]))
bestlambda = lambdas[bestindex,1]

# and now run glmnet once more with it
fit <- glmnet(xy,ys,lambda=bestlambda)

A resposta de Alice funciona bem na maioria dos casos, mas às vezes ocorre erros devido às cv.glmnet$lambdavezes retornando resultados de diferentes tamanhos, por exemplo:

Erro nos nomes de usuários <- (tmp, valor = c (0.135739830284452, 0.12368107787663,: comprimento de 'dimnames' [1] diferente da extensão da matriz).

OptimLambdaabaixo deve funcionar no caso geral e também é mais rápido, aproveitando o mclapplyprocessamento paralelo e evitando loops.

Lambdas <- function(...) {
  cv <- cv.glmnet(...)
  return(data.table(cvm=cv$cvm, lambda=cv$lambda))
}

OptimLambda <- function(k, ...) {
  # Returns optimal lambda for glmnet.
  #
  # Args:
  #   k: # times to loop through cv.glmnet.
  #   ...: Other args passed to cv.glmnet.
  #
  # Returns:
  #   Lambda associated with minimum average CV error over runs.
  #
  # Example:
  #   OptimLambda(k=100, y=y, x=x, alpha=alpha, nfolds=k)
  #
  require(parallel)
  require(data.table)
  MSEs <- data.table(rbind.fill(mclapply(seq(k), function(dummy) Lambdas(...))))
  return(MSEs[, list(mean.cvm=mean(cvm)), lambda][order(mean.cvm)][1]$lambda)
}

Você pode controlar a aleatoriedade se definir explicitamente foldid. Aqui está um exemplo de CV 5 vezes

library(caret)
set.seed(284)
flds <- createFolds(responseDiffs, k = cvfold, list = TRUE, returnTrain = FALSE)
foldids = rep(1,length(responseDiffs))
foldids[flds$Fold2] = 2
foldids[flds$Fold3] = 3
foldids[flds$Fold4] = 4
foldids[flds$Fold5] = 5

Agora execute o cv.glmnet com esses foldids.

lassoResults<-cv.glmnet(x=countDiffs,y=responseDiffs,alpha=1,foldid = foldids)

Você obterá os mesmos resultados todas as vezes.