On cherche à expliquer une variable binaire \(Y\) par deux variables quantitatives \(X_1\) et \(X_2\) à l’aide du jeu de données suivant

n <- 2000
set.seed(12345)
X1 <- runif(n)
set.seed(5678)
X2 <- runif(n)
set.seed(9012)
R1 <- X1<=0.25
R2 <- (X1>0.25 & X2>=0.75)
R3 <- (X1>0.25 & X2<0.75)
Y <- rep(0,n)
Y[R1] <- rbinom(sum(R1),1,0.25)
Y[R2] <- rbinom(sum(R2),1,0.25)
Y[R3] <- rbinom(sum(R3),1,0.75)
donnees <- data.frame(X1,X2,Y)
donnees$Y <- as.factor(donnees$Y)
  1. Séparer le jeu de données en un échantillon d’apprentissage de taille 1500 et un échantillon test de taille 500.
indapp <- 1:1500
dapp <- donnees[indapp,]
dtest <- donnees[-indapp,]
  1. On considère la régle de classification des \(k\) plus proches voisins. Pour un entier \(k\) plus petit que \(n\) et un nouvel individu \(x\), cette règle affecte à \(x\) le label majoritaire des \(k\) plus proches voisins de \(x\). Sur R on utilise la fonction knn du package class. On peut par exemple obtenir les prévisions des individus de l’échantillon test de la règle des 3 plus proches voisins avec
library(class)
knn3 <- knn(dapp[,1:2],dtest[,1:2],cl=dapp$Y,k=3)
head(knn3)
## [1] 1 1 1 0 0 1
## Levels: 0 1

Calculer l’erreur de classification de la règle des 3 plus proches voisins sur les données test.

mean(knn3!=dtest$Y)
## [1] 0.308
  1. Expliquer la fonction knn.cv

On prédit le groupe de chaque individu par validation croisée leave-one-out :

\[\widehat y_i=g_{k,i}(x_i),\quad i=1,\dots,n\]

où \(g_{k,i}\) désigne la règle de \(k\) plus proche voisins construites à partir de l’échantillon amputé de la \(i\)ème observation.

  1. Calculer l’erreur de classification de la règle des 3 plus proches voisins par validation croisée leave-one-out.
prev_cv <- knn.cv(donnees[,-3],cl=donnees$Y,k=3)

On peut alors estimer l’erreur de la règle des 10 ppv par \[\frac{1}{n}\sum_{i=1}^n1_{g_{k,i}(x_i)\neq y_i}.\]

mean(prev_cv!=donnees$Y)
## [1] 0.319
  1. On considère le vecteur de plus proches voisins suivant :
K_cand <- seq(1,500,by=20)

Proposer 2 façons de choisir une valeur de \(k\) dans ce vecteur.

err.ho <- rep(0,length(K_cand))
for (i in 1:length(K_cand)){
  knni <- knn(dapp[,1:2],dtest[,1:2],cl=dapp$Y,k=K_cand[i])
  err.ho[i] <- mean(knni!=dtest$Y)
}

Puis on choisir la valeur de \(k\) pour laquelle l’erreur est minimale.

K_cand[which.min(err.ho)]
## [1] 41
err.cv <- rep(0,length(K_cand))
for (i in 1:length(K_cand)){
  knni <- knn.cv(donnees[,-3],cl=donnees$Y,k=K_cand[i])
  err.cv[i] <- mean(knni!=donnees$Y)
} 
K_cand[which.min(err.cv)]
## [1] 41

On souhaite maintenant utiliser le package caret pour estimer des critères d’erreur et sélectionner des paramètres. On garde le même cadre que précédemment où on cherche à sélectionner le paramètre \(k\) de la règle des plus proches voisins. On pourra consulter l’url http://topepo.github.io/caret/index.html

  1. Expliquer les sorties des commandes
library(caret)
ctrl1 <- trainControl(method="LGOCV",number=1,index=list(1:1500))
KK <- data.frame(k=K_cand)
ee1 <- train(Y~.,data=donnees,method="knn",trControl=ctrl1,tuneGrid=KK)
ee1
## k-Nearest Neighbors 
## 
## 2000 samples
##    2 predictor
##    2 classes: '0', '1' 
## 
## No pre-processing
## Resampling: Repeated Train/Test Splits Estimated (1 reps, 75%) 
## Summary of sample sizes: 1500 
## Resampling results across tuning parameters:
## 
##   k    Accuracy  Kappa    
##     1  0.620     0.2298338
##    21  0.772     0.5349141
##    41  0.778     0.5459306
##    61  0.766     0.5213863
##    81  0.770     0.5305816
##   101  0.770     0.5310869
##   121  0.778     0.5469092
##   141  0.776     0.5430735
##   161  0.770     0.5315911
##   181  0.774     0.5387454
##   201  0.774     0.5382478
##   221  0.776     0.5415923
##   241  0.774     0.5362466
##   261  0.768     0.5231596
##   281  0.766     0.5182648
##   301  0.768     0.5221224
##   321  0.754     0.4924570
##   341  0.748     0.4786581
##   361  0.742     0.4653692
##   381  0.736     0.4514354
##   401  0.720     0.4136862
##   421  0.716     0.4046521
##   441  0.706     0.3799143
##   461  0.702     0.3700747
##   481  0.694     0.3509909
## 
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 121.
plot(ee1)

La probabilité d’erreur est estimée par validation hold out (même résultat que err.ho).

  1. Utiliser caret pour sélectionner \(k\) par validation croisée leave-one-out.
ctrl2 <- trainControl(method="LOOCV",number=1)
ee2 <- train(Y~.,data=donnees,method="knn",trControl=ctrl2,tuneGrid=KK)
ee2
## k-Nearest Neighbors 
## 
## 2000 samples
##    2 predictor
##    2 classes: '0', '1' 
## 
## No pre-processing
## Resampling: Leave-One-Out Cross-Validation 
## Summary of sample sizes: 1999, 1999, 1999, 1999, 1999, 1999, ... 
## Resampling results across tuning parameters:
## 
##   k    Accuracy  Kappa    
##     1  0.6230    0.2375953
##    21  0.7305    0.4503096
##    41  0.7360    0.4610133
##    61  0.7340    0.4560060
##    81  0.7325    0.4529966
##   101  0.7305    0.4491411
##   121  0.7300    0.4482950
##   141  0.7320    0.4524981
##   161  0.7330    0.4547727
##   181  0.7305    0.4498428
##   201  0.7305    0.4497260
##   221  0.7350    0.4590865
##   241  0.7340    0.4571606
##   261  0.7355    0.4601644
##   281  0.7340    0.4572758
##   301  0.7295    0.4470971
##   321  0.7255    0.4380849
##   341  0.7255    0.4378456
##   361  0.7260    0.4388098
##   381  0.7210    0.4279598
##   401  0.7205    0.4262618
##   421  0.7150    0.4145343
##   441  0.7105    0.4043357
##   461  0.7065    0.3945508
##   481  0.6990    0.3776774
## 
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 41.
plot(ee2)

  1. Faire de même pour la validation croisée 10 blocs.
ctrl3 <- trainControl(method="cv",number=10)
ee3 <- train(Y~.,data=donnees,method="knn",trControl=ctrl3,tuneGrid=KK)
ee3
## k-Nearest Neighbors 
## 
## 2000 samples
##    2 predictor
##    2 classes: '0', '1' 
## 
## No pre-processing
## Resampling: Cross-Validated (10 fold) 
## Summary of sample sizes: 1800, 1799, 1799, 1801, 1800, 1801, ... 
## Resampling results across tuning parameters:
## 
##   k    Accuracy   Kappa    
##     1  0.6239518  0.2385112
##    21  0.7294508  0.4479972
##    41  0.7344509  0.4576559
##    61  0.7314784  0.4508143
##    81  0.7314659  0.4512200
##   101  0.7314559  0.4519329
##   121  0.7294533  0.4472738
##   141  0.7289608  0.4466366
##   161  0.7329784  0.4548592
##   181  0.7344659  0.4580714
##   201  0.7329634  0.4554590
##   221  0.7329509  0.4549839
##   241  0.7334434  0.4557615
##   261  0.7309408  0.4504602
##   281  0.7279408  0.4433117
##   301  0.7269383  0.4413434
##   321  0.7234507  0.4336784
##   341  0.7219633  0.4300231
##   361  0.7174607  0.4198864
##   381  0.7154706  0.4155530
##   401  0.7134657  0.4103142
##   421  0.7039606  0.3893561
##   441  0.7004680  0.3806363
##   461  0.6914828  0.3607878
##   481  0.6804700  0.3354355
## 
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 181.
plot(ee3)

Les validations croisés peuvent se révéler couteuses en temps de calcul. On utlise souvent des techniques de parallélisation pour améliorer les performances computationnelles. Ces techniques sont relativement facile à mettre en oeuvre avec caret, on peut par exemple utiliser la librairie doParallel :

library(doParallel)
cl <- makePSOCKcluster(1)
registerDoParallel(cl)
system.time(ee3 <- train(Y~.,data=donnees,method="knn",trControl=ctrl3,tuneGrid=KK))
##    user  system elapsed 
##  12.405   0.013  12.426
stopCluster(cl)
cl <- makePSOCKcluster(4)
registerDoParallel(cl)
system.time(ee3 <- train(Y~.,data=donnees,method="knn",trControl=ctrl3,tuneGrid=KK))
##    user  system elapsed 
##   0.570   0.014   5.752
stopCluster(cl)