Existe-t-il un moyen d'utiliser la validation croisée pour effectuer une sélection de variable / fonctionnalité dans R?

10

J'ai un ensemble de données avec environ 70 variables que j'aimerais réduire. Ce que je cherche à faire, c'est d'utiliser CV pour trouver les variables les plus utiles de la manière suivante.

1) Sélectionnez au hasard disons 20 variables.

2) Utilisez stepwise/ LASSO/ lars/ etc pour choisir les variables les plus importantes.

3) Répétez ~ 50x et voyez quelles variables sont sélectionnées (pas éliminées) le plus fréquemment.

C'est dans le sens de ce que randomForestferait, mais le rfVarSelpackage ne semble fonctionner que pour les facteurs / classification et j'ai besoin de prédire une variable dépendante continue.

J'utilise R donc toutes les suggestions seraient idéalement mises en œuvre là-bas.

screechOwl
la source
Toutes les fonctionnalités sont-elles importantes? Combien d'échantillons avez-vous? Si je comprends bien le problème, vous pouvez essayer de faire une variante du boosting - sélectionner à plusieurs reprises un sous-ensemble d'échantillons et y ajuster toutes les variables et voir celles qui apparaissent le plus souvent.
Ofelia
1
Je pense que votre procédure est peu susceptible d'améliorer le LASSO, dont les implémentations dans R (par exemple glmnet et pénalisé) utilisent par défaut la validation croisée pour trouver le paramètre de régularisation "optimal". Vous pourriez envisager de répéter plusieurs fois la recherche LASSO de ce paramètre pour faire face à la variance potentiellement importante de la validation croisée (CV répété). Bien sûr, aucun algorithme ne peut battre vos connaissances préalables spécifiques au sujet.
miura

Réponses:

9

Je crois que ce que vous décrivez est déjà implémenté dans le caretpackage. Regardez la rfefonction ou la vignette ici: http://cran.r-project.org/web/packages/caret/vignettes/caretSelection.pdf

Maintenant, cela dit, pourquoi avez-vous besoin de réduire le nombre de fonctionnalités? De 70 à 20, ce n'est pas vraiment un ordre de grandeur. Je pense que vous auriez besoin de plus de 70 fonctionnalités avant d'avoir une ferme conviction que certaines fonctionnalités n'ont vraiment pas d'importance. Mais là encore, c'est là qu'intervient un prieur subjectif, je suppose.

Parkes de karité
la source
5

Il n'y a aucune raison pour que la fréquence de sélection des variables fournisse des informations que vous n'avez pas déjà obtenues de l'importance apparente des variables dans le modèle initial. Il s'agit essentiellement d'une reprise de la signification statistique initiale. vous ajoutez également un nouveau niveau d'arbitraire lorsque vous essayez de décider d'un seuil de fréquence de sélection. Le rééchantillonnage de la sélection des variables est gravement endommagé par la colinéarité en plus des autres problèmes.

Frank Harrell
la source
2

J'ai révisé ma réponse plus tôt aujourd'hui. J'ai maintenant généré quelques exemples de données sur lesquelles exécuter le code. D'autres ont suggéré à juste titre que vous envisagiez d'utiliser le package caret, avec lequel je suis d'accord. Dans certains cas, cependant, vous pouvez trouver nécessaire d'écrire votre propre code. Ci-dessous, j'ai tenté de montrer comment utiliser la fonction sample () dans R pour assigner au hasard des observations aux plis de validation croisée. J'utilise également des boucles pour effectuer une présélection variable (en utilisant une régression linéaire univariée avec une valeur de coupure indulgente de 0,1) et une construction de modèle (en utilisant une régression pas à pas) sur les dix ensembles d'apprentissage. Vous pouvez ensuite écrire votre propre code pour appliquer les modèles résultants aux plis de validation. J'espère que cela t'aides!

################################################################################
## Load the MASS library, which contains the "stepAIC" function for performing
## stepwise regression, to be used later in this script
library(MASS)
################################################################################


################################################################################
## Generate example data, with 100 observations (rows), 70 variables (columns 1
## to 70), and a continuous dependent variable (column 71)
Data <- NULL
Data <- as.data.frame(Data)

for (i in 1:71) {
for (j in 1:100) {
Data[j,i]  <- rnorm(1) }}

names(Data)[71] <- "Dependent"
################################################################################


################################################################################
## Create ten folds for cross-validation. Each observation in your data will
## randomly be assigned to one of ten folds.
Data$Fold <- sample(c(rep(1:10,10)))

## Each fold will have the same number of observations assigned to it. You can
## double check this by typing the following:
table(Data$Fold)

## Note: If you were to have 105 observations instead of 100, you could instead
## write: Data$Fold <- sample(c(rep(1:10,10),rep(1:5,1)))
################################################################################


################################################################################
## I like to use a "for loop" for cross-validation. Here, prior to beginning my
## "for loop", I will define the variables I plan to use in it. You have to do
## this first or R will give you an error code.
fit <- NULL
stepw <- NULL
training <- NULL
testing <- NULL
Preselection <- NULL
Selected <- NULL
variables <- NULL
################################################################################


################################################################################
## Now we can begin the ten-fold cross validation. First, we open the "for loop"
for (CV in 1:10) {

## Now we define your training and testing folds. I like to store these data in
## a list, so at the end of the script, if I want to, I can go back and look at
## the observations in each individual fold
training[[CV]] <- Data[which(Data$Fold != CV),]
testing[[CV]]  <- Data[which(Data$Fold == CV),]

## We can preselect variables by analyzing each variable separately using
## univariate linear regression and then ranking them by p value. First we will
## define the container object to which we plan to output these data.
Preselection[[CV]] <- as.data.frame(Preselection[CV])

## Now we will run a separate linear regression for each of our 70 variables.
## We will store the variable name and the coefficient p value in our object
## called "Preselection".
for (i in 1:70) {
Preselection[[CV]][i,1]  <- i
Preselection[[CV]][i,2]  <- summary(lm(Dependent ~ training[[CV]][,i] , data = training[[CV]]))$coefficients[2,4]
}

## Now we will remove "i" and also we will name the columns of our new object.
rm(i)
names(Preselection[[CV]]) <- c("Variable", "pValue")

## Now we will make note of those variables whose p values were less than 0.1.
Selected[[CV]] <- Preselection[[CV]][which(Preselection[[CV]]$pValue <= 0.1),] ; row.names(Selected[[CV]]) <- NULL

## Fit a model using the pre-selected variables to the training fold
## First we must save the variable names as a character string
temp <- NULL
for (k in 1:(as.numeric(length(Selected[[CV]]$Variable)))) {
temp[k] <- paste("training[[CV]]$V",Selected[[CV]]$Variable[k]," + ",sep="")}
variables[[CV]] <- paste(temp, collapse = "")
variables[[CV]] <- substr(variables[[CV]],1,(nchar(variables[[CV]])-3))

## Now we can use this string as the independent variables list in our model
y <- training[[CV]][,"Dependent"]
form <- as.formula(paste("y ~", variables[[CV]]))

## We can build a model using all of the pre-selected variables
fit[[CV]] <- lm(form, training[[CV]])

## Then we can build new models using stepwise removal of these variables using
## the MASS package
stepw[[CV]] <- stepAIC(fit[[CV]], direction="both")

## End for loop
}

## Now you have your ten training and validation sets saved as training[[CV]]
## and testing[[CV]]. You also have results from your univariate pre-selection
## analyses saved as Preselection[[CV]]. Those variables that had p values less
## than 0.1 are saved in Selected[[CV]]. Models built using these variables are
## saved in fit[[CV]]. Reduced versions of these models (by stepwise selection)
## are saved in stepw[[CV]].

## Now you might consider using the predict.lm function from the stats package
## to apply your ten models to their corresponding validation folds. You then
## could look at the performance of the ten models and average their performance
## statistics together to get an overall idea of how well your data predict the
## outcome.
################################################################################

Avant d'effectuer une validation croisée, il est important que vous preniez connaissance de son utilisation appropriée. Ces deux références offrent d'excellentes discussions sur la validation croisée:

  1. Simon RM, Subramanian J, Li MC, Menezes S.Utilisation de la validation croisée pour évaluer l'exactitude prédictive des classificateurs du risque de survie sur la base de données de grande dimension. Bref Bioinform. Mai 2011; 12 (3): 203-14. En ligne du 15 février 2011. Http://bib.oxfordjournals.org/content/12/3/203.long
  2. Richard Simon, Michael D. Radmacher, Kevin Dobbin et Lisa M. McShane. Pièges dans l'utilisation des données de puces à ADN pour la classification diagnostique et pronostique. JNCI J Natl Cancer Inst (2003) 95 (1): 14-18. http://jnci.oxfordjournals.org/content/95/1/14.long

Ces articles s'adressent aux biostatisticiens, mais seraient utiles à tous.

En outre, gardez toujours à l'esprit que l'utilisation de la régression pas à pas est dangereuse (bien que l'utilisation de la validation croisée devrait aider à atténuer le sur-ajustement). Une bonne discussion sur la régression pas à pas est disponible ici: http://www.stata.com/support/faqs/stat/stepwise.html .

Faites-moi savoir si vous avez des questions supplémentaires!

Alexandre
la source
0

Je viens de trouver quelque chose de sympa ici: http://cran.r-project.org/web/packages/Causata/vignettes/Causata-vignette.pdf

Essayez ceci peut-être lorsque vous utilisez le package glmnet

# extract nonzero coefficients
coefs.all <- as.matrix(coef(cv.glmnet.obj, s="lambda.min"))
idx <- as.vector(abs(coefs.all) > 0)
coefs.nonzero <- as.matrix(coefs.all[idx])
rownames(coefs.nonzero) <- rownames(coefs.all)[idx]
Simon Nehls
la source