Dans cet exercice nous allons illustrer l’utilisation du package \(\texttt{glmnet}\) sur un problème de classification.
Pour cela nous travaillerons sur le jeu de données South Africa Heart Disease, utilisé à fins illstratives dans le livre Elements of Statistical Learning que l’on peut télécharger sur cette page.
Nous nous limiterons ici à un problème de classification binaire, mettant donc en jeu un modèle de régression logistique, mais cette vignette illustre l’utilisation du package de manière bien plus complète.
## Warning: package 'knitr' was built under R version 3.5.2
Le jeu de données est contenu dans le fichier SAheart.data. Notons que la première colonne du fichier contient le nom des lignes et qu’il convient de les interpréter en tant que telle.
tab = read.csv("datasets/SAheart.data", row.names = 1)
extraire la variable réponse, qui est contenue dans le champ chd
transformer le descripteur qualitatif famhist en descripteur(s) quantitatifs
standardiser les descripteurs
# extract outcome #
#-----------------#
y = tab$chd
# convert to a factor
y = factor(y)
# convert famhist to numeric #
#----------------------------#
# (NB: only two levels, so can use a single variable)
tab$famhist = as.numeric(tab$famhist) - 1
# extract matrix #
#----------------#
# discard outcome
X = tab[, -which(colnames(tab) == "chd")]
# convert to matrix
X = as.matrix(X)
# scale columns
X = scale(X)
On construit le modèle avec la fonction \(\texttt{glmnet}\) et on représente le chemin de régularisation obtenu avec la fonction \(\texttt{plot.glmnet}\).
Notons que par défaut la fonction \(\texttt{glmnet}\) considère 100 valeurs du paramètre de régularisation définies automatiquement (se référer à la documentation pour davantage de précisions).
# load package
library(glmnet)
## Loading required package: Matrix
## Loading required package: foreach
## Loaded glmnet 2.0-16
# fit model
fit.lasso = glmnet(x = X, y = y, family = "binomial")
# plot
plot(fit.lasso)
Il suffit pour cela de modifier le paramètre \(\alpha\) qui définit la pénalité elastic-net: \[ \Omega(w) = \alpha ||w||_1 + \frac{1-\alpha}{2} ||w||_2^2\].
Le paramètre \(\alpha\) vaut par défaut zéro, ce qui correspond à un modèle lasso. Se référer à la documentation pour davantage de précision.
# fit model
fit.ridge = glmnet(x = X, y = y, family = "binomial", alpha = 0)
# plot
plot(fit.ridge)