Chapitre 6 Faire de la régression sur R

Les problèmes de régression et de classification supervisée consistent à expliquer et/ou prédire une sortie \(y\in\mathcal Y\) avec

  • \(\mathcal Y=\mathbb R\) pour la régression
  • \(\mathcal Y\) de cardinal fini pour la classification supervisée,

par des entrées \(x\in\mathbb R^p\). Il s’agit donc de trouver une fonction \[m:\mathbb R^p\to\mathcal Y\] à partir de données \((X_1,Y_1),\dots,(X_n,Y_n)\).

Ces données sont souvent collectées dans un dataframe df de la forme

\(Y\) \(X_1\) \(X_2\) \(\dots\) \(X_p\)
\(y_1\) \(x_{1,1}\) \(x_{1,2}\) \(\dots\) \(x_{1,p}\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(y_n\) \(x_{n,1}\) \(x_{n,2}\) \(\dots\) \(x_{n,p}\)

Le protocole pour construire un algorithme de régression sur R est toujours le même. Il faut spécifier :

  • la méthode (ou l’algorithme)
  • la variable à expliquer
  • les variables explicatives
  • le jeu de données
  • les éventuelles options de la méthode considérée.

Par exemple la commande

method(Y~X1+X3,data=df,...)

ajustera le modèle method pour expliquer \(Y\) par \(X_1\) et \(X_3\) avec les données dans df (les points représentent d’éventuelles options). Voici quelques exemples de méthodes :

fonction R algorithme Package Problème
lm modèle linéaire Reg
glm modèle logistique Class
lda analyse discriminante linéaire MASS Class
svm Support Vector Machine e1071 Class
knn.reg plus proches voisins FNN Reg
knn plus proches voisins class Class
rpart arbres rpart Reg et Class
glmnet ridge et lasso glmnet Reg et Class
gbm boosting gbm Reg et Class
randomForest forêts aléatoires randomForest Reg et Class

Remarque: pour glmnet, on ne peut pas utiliser de formule de la forme Y~.. Il faut spécifier une matrice pour les \(X\) et un vecteur pour \(Y\). La fonction model.matrix peut se révéler très utile pour calculer la matrice des \(X\).

Puisqu’il existe un grand nombre d’algorithmes pour répondre à un même problème de régression, il est important de définir des critères de performance afin de les comparer. Ces critères sont généralement inconnus et doivent être estimés à l’aide de procédure de type apprentissage/validation ou validation croisée. On a souvent besoin d’utiliser la fonction predict pour calculer ces critères. Cette fonction est une fonction générique : on peut utiliser predict pour une régression linéaire, logistique, un arbre, une forêt aléatoire… Pour obtenir l’aide de cette fonction pour

  • la régression linéaire : taper help(predict.lm)
  • la régression logisitque : taper help(predict.glm)
  • les régressions pénalisées : taper help(predict.glmnet)
  • les arbres : taper help(predict.rpart)
  • les forêts aléatoires : taper help(predict.randomForest)
  • …

Dans la suite on suppose que \(\mathcal Y=\mathbb R\) et on considère le modèle de régression \[Y=m(X)+\varepsilon.\] La performance d’un estimateur \(\widehat{m}\) de \(m\) sera mesurée par son erreur quadratique de prédiction : \[E[(Y-\widehat m(X))^2].\]

6.1 Modèle linéaire : fonctions lm et predict

On considère le modèle de régression linéaire \[Y=\beta_0+\beta_1X_1+\dots+\beta_pX_p+\varepsilon\] où \(X_1,\dots,X_p\) sont les variables explicatives, \(Y\) la variable à expliquer et \(\varepsilon\) le terme d’erreur. On fixe \(p=5\) et on considère les données suivantes :

n <- 1000
p <- 5
set.seed(1234)
X.mat <- matrix(rnorm(n*p),ncol=p)
eps <- rnorm(n,mean = 0,sd=0.5)
df <- data.frame(X.mat,eps)
df <- df %>% mutate(Y=X1+X2+X3+X4+X5+eps) %>% select(-eps)

  1. Construire un modèle linaire permettant d’expliquer \(Y\) par \(X_1,\dots,X_5\) (utiliser la fonction lm) et afficher les estimateurs de \(\beta_0,\dots,\beta_5\) (on pourra utiliser les fonctions coef et summary).

  2. On considère le jeu de données test suivant.

    m <- 500
p <- 5
set.seed(12345)
X.mat <- matrix(rnorm(m*p),ncol=5)
eps <- rnorm(m,mean = 0,sd=0.5)
df.test <- data.frame(X.mat,eps)
df.test <- df.test %>% mutate(Y=X1+X2+X3+X4+X5+eps) %>% select(-eps)

    Calculer, pour chaque individu de ce nouveau jeu de données, les prédictions faites par le modèle de la question précédente (utiliser la fonction predict avec l’option newdata).

  3. Créer un nouveau dataframe qui contiennent les valeurs prédites \(\widehat y_i\) à la question précédente sur une colonne et les valeurs observées \(y_i\) du jeu de données df.test sur une autre colonne.

  4. A l’aide du verbe summarize, calculer l’erreur quadratique moyenne (estimée) du modèle linéaire : \[\frac{1}{m}\sum_{i\in test}(\widehat y_i-y_i)^2.\]

6.2 Sélection de variables

On considère les données suivantes

n <- 1000
p <- 105
set.seed(1234)
X.mat <- matrix(rnorm(n*p),ncol=p)
eps <- rnorm(n,mean = 0,sd=0.5)
df <- data.frame(X.mat,eps)
df <- df %>% mutate(Y=X1+X2+X3+X4+X5+eps) %>% select(-eps)

issues du modèle \[Y=\beta_0+\beta_1X_1+\dots+\beta_pX_p+\varepsilon\] avec \(p=105\). On remarquera que seules les variables \(X_1,\dots,X_5\) sont explicatives.

  1. Ajuster un modèle linéaire (fonction lm) sur df et afficher les estimateurs de \(\beta_0,\dots,\beta_{105}\).

  2. On propose d’utiliser une procédure de sélection de variables backward à partir du critère BIC. Effectuer cette procédure à l’aide de la fonction step (on pourra utiliser les options direction=“backward” et k=log(n)). On appellera ce modèle mod.step.

On a sélectionné un modèle avec 8 variables : les 5 explicatives et 3 variables de bruit.

  1. Calculer les erreurs quadratiques de prévision \[\frac{1}{m}\sum_{i\in test}(\widehat y_i-y_i)^2\] des deux modèles (le modèle complet et le modèle sélectionné) en utilisant le jeu de données test suivant.

    m <- 300
p <- 105
set.seed(12345)
X.mat <- matrix(rnorm(m*p),ncol=p)
eps <- rnorm(m,mean = 0,sd=0.5)
df.test <- data.frame(X.mat,eps)
df.test <- df.test %>% mutate(Y=X1+X2+X3+X4+X5+eps) %>% select(-eps)

6.3 Régression logistique et arbre

On considère le jeu de données spam disponible ici

library(kernlab)
data(spam)

Le problème est d’expliquer la variable type (un email est un spam ou non) par les 57 autres variables.

  1. Séparer les données en un échantillon d’apprentissage dapp de taille 3000 et un échantillon test dtest de taille 1601. On pourra utiliser la fonction sample.

  2. Construire un modèle logistique permettant de répondre au problème en utilisant uniquement les données d’apprentissage. On utilisera la fonction glm avec l’option family="binomial".

  3. A l’aide de la fonction step, effectuer une sélection backward (ça peut prendre quelques minutes).

  4. A l’aide de la fonction rpart du package rpart, construire un arbre de régression (toujours sur les données d’apprentissage) pour répondre au problème. On utilisera les paramètres par défaut de la fonction.

  5. Visualiser l’arbre construit à l’aide des fonctions rpart.plot et visTree des packages rpart.plot et visNetwork

  6. Pour les 3 modèles construits (logistique, backward et arbre) calculer les prédictions de la variable type pour les individus de l’échantillon dtest. On pourra regrouper ces prévisions dans un data-frame à 3 colonnes.

  7. Ajouter au data-frame précédent une colonne où on mettra les valeurs observées de la variable à expliquer.

  8. A l’aide de summarize_at calculer les erreurs de classification des 3 modèles.

  9. Représenter les courbes ROC et calculer les AUC. On pourra consulter les pages 346 et 347 dans Cornillon et al. (2018) pour le tracé de courbes ROC sur R.

Références

Cornillon, P. A., A. Guyader, F. Husson, N. Jégou, J. Josse, N. Klutchnikoff, E. Le Pennec, E. Matzner-Løber, L. Rouvière, and B. Thieurmel. 2018. R Pour La Statistique et La Science Des Données. PUR. https://r-stat-sc-donnees.github.io.