Machine learning

L. Rouvière
[Link]@[Link]
Octobre 2021

Table des matières

I Apprentissage : contexte et formalisation 4

1 Motivations 4

2 Quelques exemples 6

3 Cadre statistique pour l’apprentissage supervisé 8

4 Exemples de fonction de perte 10

5 Estimation du risque 15

6 Le sur-apprentissage 17

7 Le package tidymodels 19

8 Annexe : le package caret 22

9 Bibliographie 25

II Algorithmes linéaires 27
1 Estimation par moindres carrés 28

2 Sélection de variables 30

3 Régularisation 33
3.1 Régression ridge . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.2 Régression Lasso . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
3.3 Variantes de ridge/lasso . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.4 Discrimination binaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

4 Support vector machine 42

4.1 SVM - cas séparable . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.2 SVM : cas non séparable . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.3 SVM non linéaire : astuce du noyau . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.4 Scores et probabilités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
4.5 Compléments : SVM multi-classes et SVR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
4.5.1 SVM multiclasses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
4.5.2 Support vector regression (SVR) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59

5 Bibliographie 62

III Algorithmes non linéaires 64

1
1 Arbres 64
1.1 Arbres binaires . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
1.2 Choix des coupures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
1.2.1 Cas de la régression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
1.2.2 Cas de la classification supervisée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
1.3 Elagage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
1.4 Importance des variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75

2 Réseaux de neurones 76
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
2.2 Le perceptron simple . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
2.3 Perceptron multicouches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
2.4 Estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
2.5 Choix des paramètres et surapprentissage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84

3 Bibliographie 86

IV Agrégation 87
1 Bagging et forêts aléatoires 87
1.1 Bagging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
1.2 Forêts aléatoires . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89
1.2.1 Algorithme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
1.2.2 Choix des paramètres . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91
1.2.3 Erreur OOB et importance des variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

2 Boosting 95
2.1 Algorithme de gradient boosting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96
2.2 Choix des paramètres . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97
2.3 Compléments/conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101

3 Bibliographie 102

Présentation
— Objectifs : comprendre les aspects théoriques et pratiques des algorithmes machine learning de référence.
— Pré-requis : théorie des probabilités, modélisation statistique, régression (linéaire et logistique). R, niveau
avancé.
— Enseignant : Laurent Rouvière [Link]@[Link]
— Recherche : statistique non paramétrique, apprentissage statistique
— Enseignements : statistique et probabilités (Université, école d’ingénieur et de commerce, formation
continue).
— Consulting : energie, finance, marketing, sport.

Programme
— Matériel :
— slides : https: // lrouviere. github. io/ machine_ learning/
— Tutoriel long : https: // lrouviere. github. io/ TUTO_ ML/
— Tutoriel court : https: // lrouviere. github. io/ machine_ learning/ tuto_ court_ ml_ sans_
correc. html

— 4 parties :
1. Machine Learning : cadre, objectif, risque...
2. Algorithmes linéaires : MCO, régularisation (ridge, lasso), SVM
3. Algorithmes non linéaires : arbres et réseaux de neurones
4. Agrégation : forêts aléatoires et boosting

2
Objectifs/questions
— Buzzword : machine learning, big data, data mining, intelligence artificielle...

— Machine learning versus statistique (traditionnelle)

— Risque =⇒ calcul ou estimation : ré-échantillonnage, validation croisée...
— Algorithmes versus estimateurs...
— Classification des algorithmes. Tous équivalents ? Cadre propice...
— ...

3
Première partie
Apprentissage : contexte et formalisation
1 Motivations
Apprentissage statistique ?
Plusieurs "définitions"
1. "... explores way of estimating functional dependency from a given collection of data" [Vapnik, 2000].
2. "...vast set of tools for modelling and understanding complex data" [James et al., 2015]

Wikipedia
L’apprentissage automatique (en anglais : machine learning), apprentissage artificiel ou apprentissage statistique
est un champ d’étude de l’intelligence artificielle qui se fonde sur des approches mathématiques et statistiques pour
donner aux ordinateurs la capacité d’apprendre à partir de donnée...
=⇒ Interface : Mathématiques-statistiques/informatique.
Constat
— Le développement des moyens informatiques fait que l’on est confronté à des données de plus en plus com-
plexes.
— Les méthodes traditionnelles se révèlent souvent peu efficaces face à ce type de données.
— Nécessité de proposer des algorithmes/modèles statistiques qui apprennent directement à partir des données.

Un peu d’histoire - voir [Besse, 2018]

Période Mémoire Ordre de grandeur

1940-70 Octet n = 30, p ≤ 10
1970 kO n = 500, p ≤ 10
1980 MO Machine Learning
1990 GO Data-Mining
2000 TO p > n, apprentissage statistique
2010 PO n explose, cloud, cluster...
2013 serveurs Big data
2017 ?? Intelligence artificielle...

Conclusion
Capacités informatiques =⇒ Data Mining =⇒ Apprentissage statistique =⇒ Big Data =⇒ Intelligence artificielle...

Approche statistique
Objectif =⇒ expliquer
— notion de modèle ;

— retrouver des lois de probabilités ;

— décisions prises à l’aide de tests statistiques, intervalles de confiance.

Exemples
— Tests indépendance/adéquation...
— Modèle linéaire : estimation, sélection de variables, analyse des résidus...
— Régression logistique...
— Séries temporelles...

4
Approche machine learning
Objectif =⇒ prédire
— notion d’algorithmes de prévision ;
— critères d’erreur de prévision ;
— calibration de paramètres (tuning).

Exemples
— Algorithmes linéaires (moindres carrés, régularisation, "SVM") ;
— Arbres, réseaux de neurones ;
— Agrégation : boosting, bagging (forêts aléatoires) ;
— Deep learning (apprentissage profond).

Statistique vs apprentissage
— Les objectifs diffèrent :
— recherche de complexité minimale en statistique =⇒ le modèle doit être interprétable !
— complexité moins importante en machine learning =⇒ on veut "juste bien prédire".
— Approches néanmoins complémentaires :
— bien expliquer =⇒ bien prédire ;
— "récentes" évolutions d’aide à l’interprétation des algorithmes ML =⇒ scores d’importance des va-
riables...
— un bon algorithme doit posséder des bonnes propriétés statistiques (convergence, biais, variance...).

Conclusion
Ne pas dissocier les deux approches.

Problématiques associées à l’apprentissage

— Apprentissage supervisé : prédire une sortie y ∈ Y à partir d’entrées x ∈ X ;

— Apprentissage non supervisé : établir une typologie des observations ;

— Règles d’association : identifier des liens entre différents produits ;

— Systèmes de recommendation : identifier les produits susceptibles d’intéresser des consommateurs.

Nombreuses applications
finance, économie, marketing, biologie, médecine...

Théorie de l’apprentissage statistique

Approche mathématique
— Ouvrage fondateur : [Vapnik, 2000]
— voir aussi [Bousquet et al., 2003].

5
The Elements of Statistical Learning [Hastie et al., 2009, James et al., 2015]

— Disponibles (avec jeux de données, codes...) aux url :

[Link] [Link]

Wikistat
— Page de cours et tutoriels très bien faits sur la statistique classique et moderne.
— On pourra notamment regarder les vignettes sur la partie apprentissage :
— [Wikistat, 2020a]
— [Wikistat, 2020b]
— ...
— Plusieurs parties de ce cours sont inspirées de ces vignettes.

2 Quelques exemples
Reconnaissance de l’écriture
Apprentissage statistique
Comprendre et apprendre un comportement à partir d’exemples.

Qu’est-ce qui est écrit ? 0, 1, 2... ?

Reconnaissance de la parole

6
 0.5 0.5
YES NO

−0.5 −0.5
0 1 0 1
1 1
BOAT GOAT

−1 −1
0 1 0 1
0.2 1
SH AO

−0.2 −1
0 1 0 1

Apprentissage sur les réseaux

Prévision de pics d’ozone

— On a mesuré pendant 366 jours la concentration maximale en ozone (V4) ;
— On dispose également d’autres variables météorologiques (température, nébulosité, vent...).
> head(Ozone)
## V1 V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 V12 V13
## 1 1 1 4 3 5480 8 20 NA NA 5000 -15 30.56 200
## 2 1 2 5 3 5660 6 NA 38 NA NA -14 NA 300
## 3 1 3 6 3 5710 4 28 40 NA 2693 -25 47.66 250
## 4 1 4 7 5 5700 3 37 45 NA 590 -24 55.04 100
## 5 1 5 1 5 5760 3 51 54 45.32 1450 25 57.02 60
## 6 1 6 2 6 5720 4 69 35 49.64 1568 15 53.78 60

Question
Peut-on prédire la concentration maximale en ozone du lendemain à partir des prévisions météorologiques ?

Détection de spam
— Sur 4 601 mails, on a pu identifier 1813 spams.

7
 — On a également mesuré sur chacun de ces mails la présence ou absence de 57 mots.
> spam %>% select(c(1:8,58)) %>% head()
## make address all num3d our over remove internet type
## 1 0.00 0.64 0.64 0 0.32 0.00 0.00 0.00 spam
## 2 0.21 0.28 0.50 0 0.14 0.28 0.21 0.07 spam
## 3 0.06 0.00 0.71 0 1.23 0.19 0.19 0.12 spam
## 4 0.00 0.00 0.00 0 0.63 0.00 0.31 0.63 spam
## 5 0.00 0.00 0.00 0 0.63 0.00 0.31 0.63 spam
## 6 0.00 0.00 0.00 0 1.85 0.00 0.00 1.85 spam

Question
Peut-on construire à partir de ces données une méthode de détection automatique de spam ?

3 Cadre statistique pour l’apprentissage supervisé

Régression vs classification
— Données de type entrée-sortie : dn = (x1 , y1 ), . . . , (xn , yn ) où xi ∈ X représente l’entrée et yi ∈ Y la sortie.

Objectifs
1. Expliquer le(s) méchanisme(s) liant les entrée xi aux sorties yi ;
2. Prédire « au mieux » la sortie y associée à une nouvelle entrée x ∈ X .

Vocabulaire
— Lorsque la variable à expliquer est quantitative (Y ⊆ R), on parle de régression.
— Lorsqu’elle est qualitative (Card(Y) fini), on parle de classification (supervisée).

Exemples
— La plupart des problèmes présentés précédemment peuvent être appréhendés dans un contexte d’apprentissage
supervisé : on cherche à expliquer une sortie y par des entrées x :

yi xi
Chiffre image Classification
Mot courbe Classification
Spam présence/absence de mots Classification
C. en O3 données météo. Régression

Remarque
La nature des variables associées aux entrées xi est variée (quanti, quali, fonctionnelle...).

Un début de formalisation mathématique

— Etant données des observations dn = {(x1 , y1 ), . . . , (xn , yn )} on cherche à expliquer/prédire les sorties yi ∈ Y
à partir des entrées xi ∈ X .
— Il s’agit donc de trouver une fonction de prévision f : X → Y telle que

f (xi ) ≈ yi , i = 1, . . . , n.

— Nécessité de se donner un critère qui permette de mesurer la qualité des fonctions de prévision f .
— Le plus souvent, on utilise une fonction de perte ` : Y × Y → R+ telle que

`(y, y 0 ) = 0 si y = y 0

`(y, y 0 ) > 0 si y 6= y 0 .

8
Vision statistique
— On suppose que les données dn = {(x1 , y1 ), . . . , (xn , yn )} sont des réalisations d’un n-échantillon Dn =
{(X1 , Y1 ), . . . , (Xn , Yn )} de loi inconnue.
— Les Xi sont des variables aléatoires à valeurs dans X , les Yi dans Y.
— Le plus souvent on supposera que les couples (Xi , Yi ), i = 1, . . . , n sont i.i.d de loi P.

Performance d’une fonction de prévision

— Etant donné une fonction de perte ` : Y × Y → R+ , la performance d’une fonction de prévision f : X → Y
est mesurée par
R(f ) = E[`(Y, f (X))]
où (X, Y ) est indépendant des (Xi , Yi ) et de même loi P .
— R(f ) est appelé risque ou erreur de généralisation de f .

Fonction de prévision optimale

— R(f ) =⇒ "Erreur moyenne" de f par rapport à la loi des données.
— Idée : trouver f qui a la plus petite erreur.

Aspect théorique
— Pour une fonction de perte ` : Y × Y → R+ donnée, le problème théorique consiste à trouver

f ? ∈ argmin R(f ) ⇐⇒ R(f ? ) ≤ R(f ) ∀f

f

— Une telle fonction f ? (si elle existe) est appelée fonction de prévision optimale pour la perte `.

Aspect pratique
— La fonction de prévision optimale f ? dépend le plus souvent de la loi P des (X, Y ) qui est en pratique
inconnue.
— Le job du statisticien est de trouver un estimateur fn = fn (., Dn ) tel que R(fn ) ≈ R(f ? ).

Définition
Un algorithme de prévision est représenté par une suite (fn )n d’applications (mesurables) telles que pour n ≥ 1,
fn : X × (X × Y)n → Y.

Propriétés statistiques d’un algorithme

— 1 un algorithme : 1 estimateur fn (.) = fn (., Dn ) de f ? .

Propriétés statistiques
— Biais : E[fn (x)] − f ? (x) =⇒ prévisions "en moyenne" ;
— Variance : V[fn (x)] =⇒ stabilité des prévisions ;
— Consistance : limn→∞ R(fn ) = R(f ? ) =⇒ précision quand n augmente ;
— ...

9
4 Exemples de fonction de perte
Choix de la fonction de perte
— Le cadre mathématique développé précédemment sous-entend qu’une fonction est performante (voire opti-
male) vis-à-vis d’un critère (représenté par la fonction de perte `)).

— Un algorithme de prévision performant pour un critère ne sera pas forcément performant pour un autre.

Conséquence pratique
Avant de s’attacher à construire un algorithme de prévision, il est capital de savoir mesurer la performance d’un
algorithme de prévision.
— Une fonction de perte ` : Y × Ye → R+ dépend de l’espace des observations Y et de celui des prévisions Y.
e

— On distingue 3 catégories de fonction de perte en fonction de ces espaces :

1. Prévisions numériques : problème de régression où on cherche à prédire la valeur de Y : ` : R × R → R+ ;

2. Prévision de groupes : problème de classification où on veut prédire un label : ` : {1, . . . , K} ×

{1, . . . , K} → R+ ;

3. Prévision de probabilités : problème de classification où on veut prédire les probabilités P(Y = k|X = x) :
` : {1, . . . , K} × RK → R+ .

Régression
— Y = R, une prévision = un réel =⇒ m : X → R ;
— Une perte = une distance entre deux nombres, par exemple la perte quadratique :

` : R × R → R+
(y, y 0 ) 7→ (y − y 0 )2

— Le risque (risque quadratique) est alors donné par

R(m) = E[(Y − m(X)]2

— et la fonction optimale (inconnue), appelée fonction de régression, par

m? (x) = E[Y |X = x].

Classification
— Y = {1, . . . , K}, une prévision = un groupe =⇒ g : X → {1, . . . , K} ;
— Une perte = 1 coût pour une mauvaise prévision, par exemple la perte indicatrice

` : {1, . . . , K} × {1, . . . , K} → R+
(y, y 0 ) 7→ 1y6=y0

— Le risque (erreur de classification) est alors donné par

R(g) = E[1g(X)6=Y ] = P(g(X) 6= Y ).

— et la fonction optimale (inconnue), appelée règle de Bayes, par

g ? (x) = argmax P(Y = k|X = x).

k

10
Classification binaire
— Y = {−1, 1}, une prévision = un groupe =⇒ g : X → {−1, 1}.
— Ce cadre permet une analyse plus fine des différents types d’erreur.
— En effet, seules 4 situations peuvent se produire
g(x) = −1 g(x) = 1
y = −1 VN FP
y=1 FN VP
— On peut les quantifier en terme de probabilités.
Pour aller plus vite

Erreurs binaires
— Spécificité =⇒ bien prédire les négatifs :

sp(g) = P(g(X) = −1|Y = −1),

— Sensibilité =⇒ bien prédire les positifs :

se(g) = P(g(X) = 1|Y = 1),

— Taux de faux négatifs =⇒ prédire négatif à tort :

fn(g) = P(g(X) = −1|Y = 1),

— Taux de faux positifs =⇒ prédire positif à tort :

fp(g) = P(g(X) = 1|Y = −1).

Critères binaires
De nombreux critères s’obtiennent en combinant ces probabilités :

EC(g) = P(g(X) 6= Y ) = fp(g)P(Y = −1) + fn(g)P(Y = 1).

Quelques critères binaires

— Balanced Accuracy :
1 1 1
P(g(X) = −1|Y = −1) + P(g(X) = 1|Y = 1) = (se(g) + sp(g)).
2 2 2
— F1 -score :
Precision × Recall
2 ,
Precision + Recall
avec
Precision P(Y = 1|g(X) = 1) et Recall = P(g(X) = 1|Y = 1).
— Kappa de Cohen...

Remarque
Mieux adapté que l’erreur de classification au cas de données déséquilibrées.

Classification (pour des probabilités)

— Y = {1, . . . , K}, une prévision = K − 1 probabilités pk (x) = P(Y = k|X = x), k = 1, . . . , K − 1.

— Les fonctions de perte sont généralement définies comme généralisation de pertes spécifiques au problème de
classification binaire.

— Classification binaire avec Y = {−1, 1} et S : X → R (S(x) = P(Y = 1|X = x) ou une transformation

bijective de cette probabilité) =⇒ fonction de score.

11
 P(Y = 1|X = x) petit P(Y = 1|X = x) grand
S(x)
s

Fonction de score
— Objectif d’un score : ordonner
— avant (d’éventuellement) classer en fixant un seuil s ∈ R :

1 si S(x) > s
gs (x) =
−1 sinon.

— Pour un seuil s donné, on a les erreurs (FP et FN)

α(s) = P(S(X) > s|Y = −1) = 1 − sp(s)
et
β(s) = P(S(X) ≤ s|Y = 1) = 1 − se(s).

Courbe ROC
— Idée : s’affranchir du choix de s en visualisant les erreurs α(s) et β(s) sur un graphe 2D pour toutes les
valeurs de s.

Définition
La courbe ROC de S est la courbe paramétrée par les valeurs de seuil s dont les abscisses et ordonnées sont définies
par 
x(s) = P(S(X) > s|Y = −1) = α(s)
y(s) = P(S(X) > s|Y = 1) = 1 − β(s).

Visualisation
— Abscisses : les faux positifs ou la spécificité ;
— Ordonnées : les faux négatifs ou la sensibilité.

Analyse de la courbe ROC

— Une proba est entre 0 et 1 =⇒ ROC vit dans le carré [0, 1]2 .
— x(−∞) = y(−∞) = 1 et x(+∞) = y(+∞) = 1 =⇒ ROC part du point (1, 1) pour arriver en (0, 0).
— ROC parfaite : il existe s? tel que α(s? ) = β(s? ) = 0 =⇒ ROC est définie par l’union des segments
[(1, 1); (0, 1)] et [(0, 1); (0, 0)].

— Mauvaise ROC : S(X) et Y sont indépendantes =⇒ x(s) = y(s) pour tout s ∈ R et ROC correspond à la
première bissectrice.

Visualisation
1.00

s=s∗
s=−∞
0.75

Score
1−fn

0.50 Aléatoire

Parfait

0.25

0.00 s=+∞
0.00 0.25 0.50 0.75 1.00
fp

12
Interprétation
On évalue la performance d’un score par sa capacité à se rapprocher le plus vite possible de la droite y = 1.

Exemple 1 Exemple 2
1.00

0.75
Score
1−fn 0.50 S1

S2

0.25

0.00
0.00 0.25 0.50 0.75 1.000.00 0.25 0.50 0.75 1.00
fp

Comparaison
— Exemple 1 : S2 meilleur que S1.
— Exemple 2 : il y a débat...
— Idée : utiliser l’aire sous la courbe.

AUC
Définition
On appelle AUC l’aire sous la courbe ROC de S.

Propriété
— 0.5 ≤ AUC(S) ≤ 1.
— Plus l’AUC est grand, meilleur est le score.

Interprétation de l’AUC
Propriété
Soit (X1 , Y1 ) et (X2 , Y2 ) indépendants et de même loi que (X, Y ), on a

AUC(S) =P(S(X1 ) > S(X2 )|Y1 = 1, Y2 = −1)

1
+ P(S(X1 ) = S(X2 )|Y1 = 1, Y2 = −1).
2
En particulier si S(X) est continue alors

AUC(S) = P(S(X1 ) ≥ S(X2 )|Y1 = 1, Y2 = −1).

Interprétation
— L’AUC correspond à la probabilité que le score ordonne correctement deux observations prélevées aléatoire-
ment dans les groupes -1 et 1.
— AUC(S) = 0.9 =⇒ dans 90% des cas, le score d’un individu positif sera plus grand que le score d’un individu
négatif.

13
Perte AUC et score optimal
— Remarquons que
AUC(S) = E[1S(X1 )>S(X2 ) |Y1 = 1, Y2 = −1].
— L’AUC peut donc s’écrire comme l’espérance d’une fonction de perte particulière

`((y1 , y2 ), (S(x1 ), S(x2 ))) = 1S(x1 )>S(x2 ) avec y1 = 1 et y2 = −1.

Proposition
Le score optimal par rapport à l’AUC est

S ? (x) = P(Y = 1|X = x).

En effet pour tout score S : X → R on a

AUC(S ? ) ≥ AUC(S).

Résumé

Perte `(y, f (x)) Risque R(f ) Champion f ?

Régression (y − m(x))2 E[Y − m(X)]2 E[Y |X = x]

Classif. binaire 1y6=g(x) P(Y 6= g(X)) Bayes

Scoring 1S(x1 )>S(x2 ) AUC(S) P(Y = 1|X = x)

Le package yardstick
— Nous verrons dans la section suivante que ces critères se calculent (ou plutôt s’estiment) en confrontant les
valeurs observées yi aux valeurs prédites d’un algorithme. Par exemple
> head(tbl)
## # A tibble: 6 x 3
## obs proba class
## <fct> <dbl> <fct>
## 1 0 0.117 0
## 2 0 0.288 0
## 3 1 0.994 1
## 4 0 0.528 1
## 5 0 0.577 1
## 6 1 0.997 1

— Le package yardstick contient un ensemble de fonctions qui permettent de calculer les critères :
https: // yardstick. tidymodels. org/ articles/ metric-types. html

Exemples
— Erreur de classification (ou plutôt accuracy) avec accuracy :
> library(yardstick)
> tbl %>% accuracy(truth=obs,estimate=class)
## # A tibble: 1 x 3
## .metric .estimator .estimate
## <chr> <chr> <dbl>
## 1 accuracy binary 0.834

— AUC avec roc_auc

14
 > tbl %>% roc_auc(truth=obs,estimate=proba,event_level="second")
## # A tibble: 1 x 3
## .metric .estimator .estimate
## <chr> <chr> <dbl>
## 1 roc_auc binary 0.926

— On peut aussi définir plusieurs critères :

> multi_metric <- metric_set(accuracy,bal_accuracy,f_meas,kap)
> tbl %>% multi_metric(truth=obs,estimate=class,event_level="second")
## # A tibble: 4 x 3
## .metric .estimator .estimate
## <chr> <chr> <dbl>
## 1 accuracy binary 0.834
## 2 bal_accuracy binary 0.834
## 3 f_meas binary 0.832
## 4 kap binary 0.668

— et tracer des courbes ROC avec roc_curve et autoplot

> tbl %>% roc_curve(truth=obs,estimate=proba,event_level="second") %>%
+ autoplot()

1.00

0.75
sensitivity

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
1 − specificity

5 Estimation du risque
Rappels
— n observations (X1 , Y1 ), . . . , (Xn , Yn ) i.i.d à valeurs dans X × Y.

Objectif
Etant donnée une fonction de perte ` : Y × Y → R+ , on cherche un algorithme de prévision fn (x) = fn (x, Dn ) qui
soit "proche" de l’oracle f ? défini par
f ? ∈ argmin R(f )
f

où R(f ) = E[`(Y, f (X))].

Question
Etant donné un algorithme fn , que vaut son risque R(fn ) ?

Risque empirique
— La loi de (X, Y ) étant inconnue en pratique, il est impossible de calculer R(fn ) = E[`(Y, fn (X))].

— Première approche : R(fn ) étant une espérance, on peut l’estimer (LGN) par sa version empirique
n
1X
Rn (fn ) = `(Yi , fn (Xi )).
n i=1

15
Problème
— L’échantillon Dn a déjà été utilisé pour construire l’algorithme de prévision fn =⇒ La LGN ne peut donc
s’appliquer !
— Conséquence : Rn (fn ) conduit souvent à une sous-estimation de R(fn ).

Une solution
Utiliser des méthodes de type validation croisée ou bootstrap.

Apprentissage - Validation ou Validation hold out

— Elle consiste à séparer l’échantillon Dn en :
1. un échantillon d’apprentissage Dn,app pour construire fn ;
2. un échantillon de validation Dn,test utilisé pour estimer le risque de fn .

Algorithme
Entrée : {A, T } une partition de {1, . . . , n} en deux parties.
1. Ajuster l’algorithme de prévision en utilisant uniquement les données d’apprentissage Dapp = {(xi , yi ) : i ∈
A}. On désigne par fapp (., Dapp ) l’algorithme obtenu.
2. Calculer les valeurs prédites fapp (xi , Dapp ) par l’algorithme pour chaque observation de l’échantillon test
Dtest = {(xi , yi ) : i ∈ T }
Retourner :
1 X
`(yi , fapp (xi , Dapp )).
|T |
i∈T

X1 Xp Y

Dapp fapp (., Dapp )

Dtest

Commentaires
Nécessite d’avoir un nombre suffisant d’observations dans
1. Dapp pour bien ajuster l’algorithme de prévision ;
2. Dtest pour bien estimer l’erreur de l’algorithme.

Validation croisée K-blocs

— Principe : répéter la hold out sur différentes partitions.

Algorithme - CV
Entrée : {B1 , . . . , BK } une partition de {1, . . . , n} en K blocs.
Pour k = 1, . . . , K :
1. Ajuster l’algorithme de prévision en utilisant l’ensemble des données privé du k e bloc, c’est-à-dire Bk =
{(xi , yi ) : i ∈ {1, . . . , n}\Bk }. On désigne par fk (.) = fk (., Bk ) l’algorithme obtenu.

16
 2. Calculer la valeur prédite par l’algorithme pour chaque observation du bloc k : fk (xi ), i ∈ Bk et en déduire
le risque sur le bloc k :
b k) = 1
X
R(f `(yi , fk (xi )).
|Bk |
i∈Bk

1
PK
Retourner : K k=1 R(f
b k ).

B1 B1 B1
B2 B2 B2
B3 B3 B3

Commentaires
— Le choix de K doit être fait par l’utilisateur (souvent K = 10).
— Avantage : plus adapté que la technique apprentissage/validation =⇒ plus stable et précis.
— Inconvénient : plus couteux en temps de calcul.

Leave one out

— Lorsque K = n, on parle de validation croisée leave one out ;
— Le risque est alors estimé par
n
b n (fn ) = 1
X
R `(Yi , fni (Xi ))
n i=1

où fni désigne l’algorithme de prévision construit sur Dn amputé de la i-ème observation.

=⇒ recommandé uniquement lorsque n est petit.

Autres approches
— Estimation par pénalisation : critère ajustement/complexité, Cp de Mallows, AIC-BIC...

— Validation croisée Monte-Carlo : répéter plusieurs fois la validation hold out ;

— Bootstrap : notamment Out Of Bag ;

— voir [Wikistat, 2020b].

6 Le sur-apprentissage
— La plupart des modèles statistiques renvoient des estimateurs qui dépendent de paramètres λ à calibrer.

Exemples
— nombres de variables dans un modèle linéaire ou logistique.
— paramètre de pénalités pour les régressions pénalisées.
— profondeur des arbres.
— nombre de plus proches voisins.

17
 — nombre d’itérations en boosting.
— ...

Remarque importante
Le choix de ces paramètres est le plus souvent crucial pour la performance de l’estimateur sélectionné.
— Le paramètre λ à sélectionner représente la complexité du modèle :

Complexité =⇒ compromis biais/variance

— λ petit =⇒ modèle peu flexible =⇒ mauvaise adéquation sur les données =⇒ biais %, variance &.
— λ grand =⇒ modèle trop flexible =⇒ sur-ajustement =⇒ biais &, variance %.

Overfitting
Sur-ajuster signifie que le modèle va (trop) bien ajuster les données d’apprentissage, il aura du mal à s’adapter à
de nouveaux individus.
Erreur

Ajustement

Prévision

Sur−apprentissage

Complexité

Overfitting en régression

0
Y

−1

−4 0 4
X

18
 1

10ppv
0

Y
199ppv

1ppv

−1

−4 0 4
X

Overfitting en classification supervisée

1.00

0.9
0.75

0.6
X2

X2

0.50

0.3 0.25

0.0 0.00
0.0 0.3 0.6 0.9 0.00 0.25 0.50 0.75 1.00
X1 X1

1ppv 400ppv
1.00

0.75
X2

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00
X1

Application shiny
[Link]

7 Le package tidymodels
Présentation du package
— Successeur de caret pour conduire des projets machine learning sur R.

— Meta package qui inclut

19
 — rsample : pour ré-échantilloner
— yardstick : pour les fonctions de perte
— recipe : pour les recettes de préparation... des données
— tune : pour calibrer les algorithme
— ...

— Tutoriel : [Link]

Calibrer des paramètres

— Tous les algorithmes dépendent de paramètres θ que l’utilisateur doit sélectionner.
— Le procédé est toujours le même et peut se résumer dans l’algorithme suivant.

Choix de paramètres par minimisation du risque (grid search)

Entrées :
— Une grille [Link] de valeurs pour θ ;
— Un risque de prévision R ;
— un algorithme d’estimation du risque.
Pour chaque θ dans [Link] :
• Estimer R(fn,θ ) par l’algorithme choisi =⇒ R(f
b n,θ )
Retourner : θb une valeur de θ qui minimise R(f
b n,θ ).

— Ce procédé est automatisé dans tidymodels.

— Il faut spécifier les différents paramètres :

— la méthode (logistique, ppv, arbre, randomForest...)
— Une grille pour les paramètres (nombre de ppv...)
— Le critère de performance (erreur de classification, AUC, risque quadratique...)
— La méthode d’estimation du critère (apprentissage validation, validation croisée, bootstrap...)

— Nous l’illustrons à travers le choix du nombre de voisins de l’algorithme des k-ppv.

Les données
— Une variable binaire à expliquer par 2 variables continues
> head(don.2D.500)
## # A tibble: 6 x 3
## X1 X2 Y
## <dbl> <dbl> <fct>
## 1 0.721 0.209 0
## 2 0.876 0.766 1
## 3 0.761 0.842 1
## 4 0.886 0.934 0
## 5 0.456 0.676 0
## 6 0.166 0.859 1

Le workflow
— On commence par renseigner l’algorithme et la manière dont on va choisir les paramètres.
> library(tidymodels)
> tune_spec <-
+ nearest_neighbor(neighbors=tune(),weight_func="rectangular") %>%
+ set_mode("classification") %>%
+ set_engine("kknn")

— On créé ensuite la workflow :

> ppv_wf <- workflow() %>%
+ add_model(tune_spec) %>%
+ add_formula(Y ~ .)

20
Ré-échantillonnage et grille de paramètres
— On spécifie ensuite la méthode de ré-échantillonnage, ici une validation croisée 10 blocs
> [Link](12345)
> re_ech_cv <- vfold_cv(don.2D.500,v=10)
> re_ech_cv %>% head()
## # A tibble: 6 x 2
## splits id
## <list> <chr>
## 1 <split [450/50]> Fold01
## 2 <split [450/50]> Fold02
## 3 <split [450/50]> Fold03
## 4 <split [450/50]> Fold04
## 5 <split [450/50]> Fold05
## 6 <split [450/50]> Fold06

— Puis vient la grille de paramètres

> grille_k <- tibble(neighbors=1:100)

=⇒ consulter https: // www. tidymodels. org/ find/ parsnip/ pour trouver les identifiants des algorithmes
et de leurs paramètres.

Estimation du risque
— Fonction tune_grid
> tune_grid(...,resamples=...,grid=...,metrics=...)

— Calcul du risque pour chaque valeur de la grille :

> [Link] <- ppv_wf %>%
+ tune_grid(
+ resamples = re_ech_cv,
+ grid = grille_k,
+ metrics=metric_set(accuracy))

— On lit les résultats avec collect_metrics :

> [Link] %>% collect_metrics() %>% select(1:5) %>% head()
## # A tibble: 6 x 5
## neighbors .metric .estimator mean n
## <int> <chr> <chr> <dbl> <int>
## 1 1 accuracy binary 0.618 10
## 2 2 accuracy binary 0.618 10
## 3 3 accuracy binary 0.672 10
## 4 4 accuracy binary 0.672 10
## 5 5 accuracy binary 0.69 10
## 6 6 accuracy binary 0.69 10

Visualisation des erreurs

> tbl <- [Link] %>% collect_metrics()
> ggplot(tbl)+aes(x=neighbors,y=mean)+geom_line()+ylab("Accuracy")

0.725

0.700
Accuracy

0.675

0.650

0.625

0 25 50 75 100
neighbors

21
Sélection du meilleur paramètre
— On visualise les meilleures valeurs de paramètres :
> [Link] %>% show_best() %>% select(1:6)
## # A tibble: 5 x 6
## neighbors .metric .estimator mean n std_err
## <int> <chr> <chr> <dbl> <int> <dbl>
## 1 13 accuracy binary 0.72 10 0.0255
## 2 14 accuracy binary 0.72 10 0.0255
## 3 35 accuracy binary 0.72 10 0.0207
## 4 36 accuracy binary 0.72 10 0.0207
## 5 39 accuracy binary 0.718 10 0.0199

— et on choisit celle qui maximise l’accuracy :

> best_k <- [Link] %>% select_best()
> best_k
## # A tibble: 1 x 2
## neighbors .config
## <int> <chr>
## 1 13 Preprocessor1_Model013

Algorithme final et prévision

— L’algorithme final s’obtient en entrainant la méthode sur toutes les données pour la valeur de paramètre
sélectionné :
> final_ppv <-
+ ppv_wf %>%
+ finalize_workflow(best_k) %>%
+ fit(data = don.2D.500)

— On peut maintenant prédire de nouveaux individus :

> newx <- tibble(X1=0.3,X2=0.8)
> predict(final_ppv,new_data=newx)
## # A tibble: 1 x 1
## .pred_class
## <fct>
## 1 0

Conclusion
— Les choix de l’utilisateur sont des paramètres de la procédure.

— =⇒ facilement personnalisable.

— Aisé de changer le critère, la méthode de ré-échantillonnage...

8 Annexe : le package caret

Le package caret
— Il permet d’évaluer la performance de plus de 230 méthodes : [Link]

— Il suffit d’indiquer :
— la méthode (logistique, ppv, arbre, randomForest...)
— Une grille pour les paramètres (nombre de ppv...)
— Le critère de performance (erreur de classification, AUC, risque quadratique...)
— La méthode d’estimation du critère (apprentissage validation, validation croisée, bootstrap...)

Apprentissage-validation

22
> library(caret)
> K_cand <- [Link](k=seq(1,500,by=20))
> library(caret)
> ctrl1 <- trainControl(method="LGOCV",number=1,index=list(1:1500))
> e1 <- train(Y~.,data=donnees,method="knn",trControl=ctrl1,tuneGrid=K_cand)
> e1
## 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.2382571
## 21 0.718 0.4342076
## 41 0.722 0.4418388

## 61 0.718 0.4344073
## 81 0.720 0.4383195
## 101 0.714 0.4263847
## 121 0.716 0.4304965
## 141 0.718 0.4348063
## 161 0.718 0.4348063
## 181 0.718 0.4348063
## 201 0.720 0.4387158
## 221 0.718 0.4350056
## 241 0.718 0.4350056
## 261 0.722 0.4428232
## 281 0.714 0.4267894
## 301 0.714 0.4269915
## 321 0.710 0.4183621
## 341 0.696 0.3893130
## 361 0.696 0.3893130
## 381 0.688 0.3727988
## 401 0.684 0.3645329
## 421 0.686 0.3686666
## 441 0.686 0.3679956
## 461 0.684 0.3638574
## 481 0.680 0.3558050
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 261.

Validation croisée

> library(doMC)
> registerDoMC(cores = 3)
> ctrl2 <- trainControl(method="cv",number=10)
> e2 <- train(Y~.,data=dapp,method="knn",trControl=ctrl2,tuneGrid=K_cand)
> e2
## k-Nearest Neighbors
##
## 1500 samples
## 2 predictor
## 2 classes: ’0’, ’1’
##
## No pre-processing
## Resampling: Cross-Validated (10 fold)
## Summary of sample sizes: 1350, 1350, 1350, 1350, 1350, 1350, ...
## Resampling results across tuning parameters:
##
## k Accuracy Kappa
## 1 0.6240000 0.2446251
## 21 0.7393333 0.4745290
## 41 0.7306667 0.4570024
## 61 0.7340000 0.4636743

23
## 81 0.7333333 0.4632875
## 101 0.7313333 0.4593480
## 121 0.7326667 0.4624249
## 141 0.7333333 0.4640787
## 161 0.7366667 0.4708178
## 181 0.7313333 0.4602309
## 201 0.7326667 0.4626618
## 221 0.7293333 0.4559741
## 241 0.7306667 0.4585960
## 261 0.7353333 0.4676751
## 281 0.7286667 0.4537842
## 301 0.7253333 0.4463516
## 321 0.7173333 0.4294524
## 341 0.7113333 0.4168003
## 361 0.7080000 0.4099303
## 381 0.7140000 0.4213569
## 401 0.7073333 0.4073761
## 421 0.7100000 0.4126434
## 441 0.7066667 0.4054984
## 461 0.6966667 0.3844183
## 481 0.6860000 0.3612515
##
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 21.

Validation croisée répétée

> ctrl3 <- trainControl(method="repeatedcv",repeats=5,number=10)

> e3 <- train(Y~.,data=dapp,method="knn",trControl=ctrl3,tuneGrid=K_cand)
> e3
## k-Nearest Neighbors
##
## 1500 samples
## 2 predictor
## 2 classes: ’0’, ’1’
##
## No pre-processing
## Resampling: Cross-Validated (10 fold, repeated 5 times)
## Summary of sample sizes: 1350, 1350, 1350, 1350, 1350, 1350, ...
## Resampling results across tuning parameters:
##
## k Accuracy Kappa
## 1 0.6232000 0.2438066
## 21 0.7354667 0.4665640
## 41 0.7314667 0.4585144
## 61 0.7317333 0.4592608
## 81 0.7302667 0.4568784
## 101 0.7310667 0.4589567

## 121 0.7320000 0.4609326

## 141 0.7322667 0.4616077
## 161 0.7336000 0.4643374
## 181 0.7340000 0.4649895
## 201 0.7332000 0.4632905
## 221 0.7325333 0.4620114
## 241 0.7316000 0.4600484
## 261 0.7305333 0.4578098
## 281 0.7286667 0.4536040
## 301 0.7238667 0.4434101
## 321 0.7189333 0.4330787
## 341 0.7136000 0.4215865
## 361 0.7122667 0.4183400
## 381 0.7098667 0.4131761
## 401 0.7090667 0.4112403
## 421 0.7058667 0.4043164
## 441 0.7001333 0.3920207
## 461 0.6952000 0.3811374
## 481 0.6872000 0.3636126
##
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was k = 21.

24
Critère AUC
> donnees1 <- donnees
> names(donnees1)[3] <- c("Class")
> levels(donnees1$Class) <- c("G0","G1")
> ctrl11 <- trainControl(method="LGOCV",number=1,index=list(1:1500),
+ classProbs=TRUE,summary=twoClassSummary)
> e4 <- train(Class~.,data=donnees1,method="knn",trControl=ctrl11,
+ metric="ROC",tuneGrid=K_cand)
> e4
## k-Nearest Neighbors
##
## 2000 samples
## 2 predictor
## 2 classes: ’G0’, ’G1’
##
## No pre-processing
## Resampling: Repeated Train/Test Splits Estimated (1 reps, 75%)
## Summary of sample sizes: 1500
## Resampling results across tuning parameters:
##

## k ROC Sens Spec

## 1 0.6190866 0.5983264 0.6398467
## 21 0.7171484 0.6903766 0.7432950
## 41 0.7229757 0.6861925 0.7547893
## 61 0.7200500 0.6945607 0.7394636
## 81 0.7255567 0.6945607 0.7432950
## 101 0.7319450 0.6903766 0.7356322
## 121 0.7382452 0.6945607 0.7356322
## 141 0.7353757 0.7029289 0.7318008
## 161 0.7308549 0.7029289 0.7318008
## 181 0.7351272 0.7029289 0.7318008
## 201 0.7340050 0.7029289 0.7356322
## 221 0.7324099 0.7071130 0.7279693
## 241 0.7349028 0.7071130 0.7279693
## 261 0.7365780 0.7071130 0.7356322
## 281 0.7349749 0.6987448 0.7279693
## 301 0.7356963 0.7029289 0.7241379
## 321 0.7341493 0.6861925 0.7318008
## 341 0.7343898 0.6527197 0.7356322
## 361 0.7306385 0.6527197 0.7356322
## 381 0.7301816 0.6359833 0.7394636
## 401 0.7270957 0.6276151 0.7356322
## 421 0.7255487 0.6317992 0.7356322

## 441 0.7258933 0.6192469 0.7471264

## 461 0.7220619 0.6150628 0.7471264
## 481 0.7236330 0.6108787 0.7432950
##
## ROC was used to select the optimal model using the largest value.
## The final value used for the model was k = 121.

9 Bibliographie
Références

Biblio1

[Besse, 2018] Besse, P. (2018). Science des données - Apprentissage Statistique. INSA - Toulouse. [Link]
[Link]/~besse/pub/Appren_stat.pdf.
[Bousquet et al., 2003] Bousquet, O., Boucheron, S., and Lugosi, G. (2003). Introduction to Statistical Learning
Theory, chapter Advanced Lectures on Machine Learning. Springer.
[Clémençon et al., 2008] Clémençon, S., Lugosi, G., and Vayatis, N. (2008). Ranking and empirical minimization
of u-statistics. The Annals of Statistics, 36(2) :844–874.
[Hastie et al., 2009] Hastie, T., Tibshirani, R., and Friedman, J. (2009). The Elements of Statistical Learning :
Data Mining, Inference, and Prediction. Springer, second edition.

25
[James et al., 2015] James, G., Witten, D., Hastie, T., and Tibshirani, R. (2015). The Elements of Statistical
Learning : Data Mining, Inference, and Prediction. Springer.
[Vapnik, 2000] Vapnik, V. (2000). The Nature of Statistical Learning Theory. Springer, second edition.
[Wikistat, 2020a] Wikistat (2020a). Apprentissage machine — introduction. [Link]
[Link].
[Wikistat, 2020b] Wikistat (2020b). Qualité de prévision et risque. [Link]
pdf.

26
Deuxième partie
Algorithmes linéaires
— Rappel : une fonction de prévision f : Rd → R.

Fonction de prévision linéaire

Une fonction de prévision est dite linéaire si elle se met sous la forme
f (x) = fβ (x) = β0 + β1 x1 + . . . + βd xd .

Remarque
— Possibilité d’inclure des effets non linéaires :
fβ (x) = β0 + β11 x1 + β12 x21 + β21 x2 + β22 x22 + β12 x1 x2 + β31 x3 + β32 exp(x3 ) . . .

— Variables qualitatives codées en indicatrices :

fβ (x) = β0 + β1 1x1 =A + β2 1x1 =B + β3 1x1 =C + . . .
muni d’une contrainte identifiante, par exemple β1 = 0.

Régression
— Y à valeurs dans R.
— On utilise souvent le terme modèle linéaire :
yi = β0 + β1 xi1 + . . . + βd xid + εi
où les εi sont i.i.d tels que E[εi ] = 0 et V[εi ] = σ 2 .
— Fonction de prévision :
mβ (x) = E[Y |X = x] = β0 + β1 x1 + . . . + βd xd .

Classification binaire
— Y à valeurs dans {0, 1}.
— La classification s’effectue à partir de la probabilité
p(x) = P(Y = 1|X = x).

— Frontière entre les deux classes :

 
p(x)
{x : p(x) = 1 − p(x)} = x : log =0 .
1 − p(x)
— La frontière est linéaire si
p(x)
log = β0 + β1 x 1 + . . . + βd x d .
1 − p(x)
— =⇒ Modèle logistique.

Questions
1. Comment calculer (ou plutôt estimer) les βj ?
— MCO-vraisemblance
— Approches régularisées =⇒ ridge-lasso...
— Machines à support vecteur (SVM).

2. Comment choisir la combinaison linéaire ?

— Sélection de variables
— Régression sur composantes =⇒ PCR-PLS...
— Transformation de variables =⇒ résidus partiels, modèle additifs...

27
1 Estimation par moindres carrés
Minimiser les erreurs
— Les données : (x1 , y1 ), . . . , (xn , yn ) à valeurs dans Rd × R.
— Le modèle
yi = β0 + β1 xi1 + . . . + βd xid + εi
— εi représente l’écart (ou l’erreur) entre la prévision du modèle β et la valeur observée.

Idée
Choisir β de manière à minimiser ces erreurs.

Estimateurs des moindres carrés

Définition
On appelle critère des moindres carrés ordinaires ou somme des carrés résiduelles la fonction de β :
n
X
SCR(β) = (yi − (β0 + β1 xi1 + · · · + βd xid ))2 = kY − Xβk2
i=1

avec    
y1 1 x11 ... x1d
Y =  ... 
 .. .. ..  .
et X = .
 
. . 
yn 1 xn1 ... xnd

Propriété
Si X est de plein rang alors l’estimateur des MCO βb = (Xt X)−1 Xt Y minimise SCR(β).

Exemple
— Données Hitters, 263 individus, 20 variables
> Hitters %>% select(c(1:5,19)) %>% head()
## AtBat Hits HmRun Runs RBI Salary
## 1 315 81 7 24 38 475.0
## 2 479 130 18 66 72 480.0
## 3 496 141 20 65 78 500.0
## 4 321 87 10 39 42 91.5
## 5 594 169 4 74 51 750.0
## 6 185 37 1 23 8 70.0

— Problème : Expliquer/prédire le salaire (Salary) par les autres variables.

— Calcul des estimateurs MCO avec lm :
> mod <- lm(Salary~.,data=Hitters)
> coef(mod)[1:5]
## (Intercept) AtBat Hits HmRun Runs
## 163.103588 -1.979873 7.500768 4.330883 -2.376210

— Prévision du salaire de nouveaux individus

> xnew %>% select(1:5)
## AtBat Hits HmRun Runs RBI
## 1 585 139 31 93 94

avec predict :
> predict(mod,newdata=xnew)
## 1
## 1129.376

28
Modèle gaussien
— En supposant de plus que les erreurs εi suivent une loi Gaussienne, on obtient la loi des estimateurs

βbj − βj
∼ Tn−(d+1) .
σ
bβbj

— On en déduit des procédures de test :

> broom::tidy(mod) %>% head()
## # A tibble: 6 x 5
## term estimate [Link] statistic [Link]
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 163. 90.8 1.80 0.0736
## 2 AtBat -1.98 0.634 -3.12 0.00201
## 3 Hits 7.50 2.38 3.15 0.00181
## 4 HmRun 4.33 6.20 0.698 0.486
## 5 Runs -2.38 2.98 -0.797 0.426
## 6 RBI -1.04 2.60 -0.402 0.688

— Ainsi que des intervalles de confiance pour les paramètres :

> confint(mod) %>% head()
## 2.5 % 97.5 %
## (Intercept) -15.709647 341.9168228
## AtBat -3.228667 -0.7310792
## Hits 2.817562 12.1839734
## HmRun -7.884569 16.5463352
## Runs -8.247625 3.4952055
## RBI -6.168102 4.0781779

— ou pour les prévisions :

> predict(mod,newdata=xnew,interval="confidence")
## fit lwr upr
## 1 1129.376 889.2244 1369.528

Cas du modèle logistique

— Toutes ces notions se généralisent (assez) rapidement au modèle logistique

pβ (x)
log = β0 + β1 x 1 + . . . + βd x d .
1 − pβ (x)

— Le critère des MCO est remplacé par la log-vraisemblance (à maximiser) :

n
X
L(y1 , . . . , yn ; β) = [yi xti β − log(1 + exp(xti β))].
i=1

— Pas de solution explicite mais de (bons) algorithmes qui convergent vers le max.

Exemple
— On considère les données SAheart :
> head(SAheart)
## sbp tobacco ldl adiposity famhist typea obesity alcohol age chd
## 1 160 12.00 5.73 23.11 Present 49 25.30 97.20 52 1
## 2 144 0.01 4.41 28.61 Absent 55 28.87 2.06 63 1
## 3 118 0.08 3.48 32.28 Present 52 29.14 3.81 46 0
## 4 170 7.50 6.41 38.03 Present 51 31.99 24.26 58 1
## 5 134 13.60 3.50 27.78 Present 60 25.99 57.34 49 1
## 6 132 6.20 6.47 36.21 Present 62 30.77 14.14 45 0

— Problème : expliquer/prédire la variable binaire chd par les autres variables.

— On obtient les estimateurs avec glm

29
 > logit <- glm(chd~.,data=SAheart,family="binomial")
> broom::tidy(logit)
## # A tibble: 10 x 5
## term estimate [Link] statistic [Link]
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) -6.15 1.31 -4.70 0.00000258
## 2 sbp 0.00650 0.00573 1.14 0.256
## 3 tobacco 0.0794 0.0266 2.98 0.00285
## 4 ldl 0.174 0.0597 2.92 0.00355
## 5 adiposity 0.0186 0.0293 0.635 0.526
## 6 famhistPresent 0.925 0.228 4.06 0.0000490
## 7 typea 0.0396 0.0123 3.21 0.00131
## 8 obesity -0.0629 0.0442 -1.42 0.155
## 9 alcohol 0.000122 0.00448 0.0271 0.978
## 10 age 0.0452 0.0121 3.73 0.000193

— Les prévisions de la probabilité de l’évènement {chd=1} pour de nouveaux individus

> xnew
## sbp tobacco ldl adiposity famhist typea obesity alcohol age
## 1 146 0 6.62 25.69 Absent 60 28.07 8.23 63

— s’obtiennent avec predict :

> predict(logit,newdata=xnew,type="response")
## 1
## 0.4719671

Conclusion
Remarque
La qualité de ces modèles (et donc des prévisions) reposent sur deux postulats :
1. le modèle est bon : Y s’explique bien par une combinaison linéaire des X ;
2. les estimateurs sont bons : ils possèdent de bonnes propriétés statistiques.
— La qualité du modèle est toujours difficile à vérifier =⇒ ajouter d’autres effets dans la combinaison linéaire
(quadratique, interactions...).
— On en sait plus sur la performance des estimateurs :
1. Trop de variables =⇒ % de la variance (sur-ajustement).
2. Colinéarités =⇒ % de la variance (sur-ajustement).

2 Sélection de variables
— Une approche naturelle pour répondre aux 2 problèmes évoqués précédemment est de sélectionner des va-
riables explicatives parmi {X1 , . . . , Xd }.

Idée
Supprimer les variables
— qui n’expliquent pas Y .
— dont l’effet est déjà expliqué par d’autres variables
=⇒ ce n’est pas parce qu’une variable n’est pas sélectionnée qu’elle n’est pas liée à Y !

Best subset selection

— d variables explicatives =⇒ 2d modèles concurrents.
— Idée : construire les 2d modèles et les comparer.

Algorithme BSS
Entrée : un critère de choix de modèle (AIC, BIC. . .).
Pour j = 0, . . . , d :
 
d
1. Construire les modèles linéaires à j variables ;
j
2. Choisir parmi ces modèles celui qui a la plus petite SCR. On note Mj le modèle sélectionné.
Retourner : le meilleur modèle parmi M0 , M1 , . . . , Md au sens du critère de choix de modèle.

30
Exemples de critères (voir [Cornillon et al., 2019])
— AIC : Akaike Information Criterion
−2Ln (β̂) + 2d.
— BIC : Bayesian Information Criterion
−2Ln (β̂) + log(n)d.
— R2 ajusté :
n−1 SSR kŶ − Ȳ1k2
Ra2 = 1 − (1 − R2 ) où R2 = = .
n−d+1 SST kY − Ȳ1k2
— Cp de Mallow : !
n
1 X
Cp = (Yi − Ŷi )2 + 2dσ̂ 2 .
n i=1

Ajustement/complexité
— Ces critères sont constitués de deux parties :
1. une qui mesure la qualité d’ajustement du modèle ;
2. une autre qui mesure sa complexité.

Exemple AIC
— −2Ln (β̂) mesure l’ajustement ;
— 2p mesure la complexité.

=⇒ l’idée est de choisir un modèle de complexité minimale qui ajuste bien les données.

Le coin R
— On peut utiliser les packages leaps et bestglm.
— On propose de présenter bestglm qui fait appel à leaps pour la régression et fonctionne également pour le
modèle logistique.
> Hitters1 <- Hitters[,c(1:18,20,19)]
> [Link] <- bestglm(Hitters1)
> [Link]$Subsets %>% select(c(1:5,22)) %>% head()
## (Intercept) AtBat Hits HmRun Runs BIC
## 0 TRUE FALSE FALSE FALSE FALSE 3213.768
## 1 TRUE FALSE FALSE FALSE FALSE 3117.350
## 2 TRUE FALSE TRUE FALSE FALSE 3079.270
## 3 TRUE FALSE TRUE FALSE FALSE 3072.569
## 4 TRUE FALSE TRUE FALSE FALSE 3066.387
## 5 TRUE TRUE TRUE FALSE FALSE 3064.125

— On obtient le modèle sélectionné avec :

> [Link]$BestModel %>% broom::tidy()
## # A tibble: 7 x 5
## term estimate [Link] statistic [Link]
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 91.5 65.0 1.41 1.60e- 1
## 2 AtBat -1.87 0.527 -3.54 4.70e- 4
## 3 Hits 7.60 1.66 4.57 7.46e- 6
## 4 Walks 3.70 1.21 3.06 2.49e- 3
## 5 CRBI 0.643 0.0644 9.98 5.05e-20
## 6 DivisionW -123. 39.8 -3.09 2.24e- 3
## 7 PutOuts 0.264 0.0748 3.53 4.84e- 4

Remarque
— L’approche exhaustive peut se révéler coûteuse en temps de calcul lorsque d > 50.
— On utilise généralement des méthodes pas à pas dans ce cas.

31
Pas à pas ascendant
Algorithme forward
Entrée : un critère de choix de modèle (AIC, BIC. . .)
1. Construire M0 le modèle linéaire qui contient uniquement la constante ;
2. Pour j = 0, . . . , d − 1 :
(a) Construire les d − j modèles linéaires en ajoutant une variable, parmi les variables non utilisées, à Mj ;
(b) Choisir, parmi ces d − j modèles, celui qui minimise la SCR → Mj+1 .
Retourner : le meilleur modèle parmi M0 , M1 , . . . , Md au sens du critère de choix de modèle.

Le coin R
Utiliser method=forward dans bestglm.

Pas à pas descendant

Algorithme backward
Entrée : un critère de choix de modèle (AIC, BIC. . .)
1. Construire Md le modèle linéaire complet (avec toutes les variables explicatives) ;
2. Pour j = d, . . . , 1 :
(a) Construire les j modèles linéaires en supprimant une variable, parmi les variables non utilisées, à Mj ;
(b) Choisir, parmi ces j modèles, celui qui minimise la SCR → Mj−1 .
Retourner : le meilleur modèle parmi M0 , M1 , . . . , Md au sens du critère de choix de modèle.

Le coin R
Utiliser method=backward dans bestglm.

32
3 Régularisation

Erreur

Ajustement

Prévision

Sur−apprentissage

Complexité

Complexité linéaire
Le nombre de variables est une mesure de la complexité des algorithmes linéaires.

Illustration numérique
— On génère des données (xi , yi ), i = 1, . . . , 500 selon le modèle

yi = 1xi1 + 0xi2 + . . . + 0xiq + εi

où x1 , . . . , xq , ε sont i.i.d. de loi N (0, 1).

— Seule X1 est explicative, les q − 1 autres variables peuvent être vues comme du bruit.
— On calcule l’estimateur de MCO de β1 sur 1000 répétitions. On trace les boxplot de ces estimateurs pour q = 0
et q = 400.

1.2
beta1

1.0

0.8

1 400
nb_var

Conclusion
Plus de variance (donc moins de précision) lorsque le nombre de variables inutiles augmente.
— Lorsque le nombre de variables d est grand, les estimateurs des moindres carrés du modèle linéaire

Y = β1 X1 + . . . + βd Xd + ε

possèdent généralement une grande variance.

Idée des méthodes pénalisés

— Contraindre la valeur des estimateurs des moindres carrés de manière à réduire la variance (quitte à augmenter
un peu le biais).

33
 — Comment ? En imposant une contrainte sur la valeur des estimateurs des moindres carrés :
 2
n
X d
X
β̂ pen = argmin yi − xij βj 
β i=1 j=1

sous la contrainte kβk? ≤ t.

Questions
— Quelle norme utiliser pour la contrainte ?

— Existence/unicité des estimateurs ? Solutions explicites du problème d’optimisation ?

— Comment choisir t ?
— t petit =⇒ estimateurs contraints (proche de 0) ;
— t grand =⇒ estimateurs des moindres carrés (non pénalisés).

3.1 Régression ridge

— La régression ridge consiste à minimiser le critère des moindres carrés pénalisé par la norme 2 des coefficients.

Définition
1. Les estimateurs ridge β̂ R s’obtiennent en minimisant
 2
n
X d
X d
X
yi − β0 − xij βj  sous la contrainte βj2 ≤ t (1)
i=1 j=1 j=1

2. ou de façon équivalente
  2 

X n Xd d
X 

β̂ R = argmin yi − β0 − xij βj  + λ βj2 . (2)
β 
 i=1 j=1 j=1



Quelques remarques
— Les définitions (1) et (2) sont équivalentes dans le sens où pour tout t il existe un unique λ tels que les
solutions aux deux problèmes d’optimisation coïncident.

— La constante β0 n’entre généralement pas dans la pénalité.

— L’estimateur dépend bien entendu du paramètre t (ou λ) : β̂ R = β̂ R (t) = β̂ R (λ).

— Les variables explicatives sont le plus souvent réduites pour éviter les problèmes d’échelle dans la pénalité.

Exemple avec les données Hitters

— Il existe plusieurs fonctions et packages qui permettent de faire de la régression pénalisée sur R. Nous présen-
tons ici glmnet.
— glmnet n’accepte pas d’objet formule. Il faut spécifier la matrice des X et le vecteur des Y :
> Hitters.X <- [Link](Salary~.,data=Hitters)[,-1]

Ridge avec glmnet

34
> library(glmnet)
> [Link] <- glmnet(Hitters.X,Hitters$Salary,alpha=0)
> par(mfrow=c(1,2))
> plot([Link],lwd=2)
> plot([Link],lwd=2,xvar="lambda")

19 19 19 19 19 19 19 19 19 19

50

50
0

0
Coefficients

Coefficients
−50

−50
−100

−100
0 50 100 150 200 4 6 8 10 12

L1 Norm Log Lambda

Propriétés des estimateurs ridge

Propriétés
1. Lorsque les variables explicatives sont centrée-réduites, l’estimateur Ridge solution de (2) s’écrit

β̂ R = β̂ R (λ) = (Xt X + λI)−1 Xt Y.

2. On déduit
biais(β̂ R ) = −λ(Xt X + λI)−1 β
et
V(β̂ R ) = σ 2 (Xt X + λI)−1 Xt X(Xt X + λI)−1 .

Commentaires
— Si λ = 0, on retrouve le biais et la variance de l’estimateur des MCO.
— λ % =⇒ biais % et variance & et réciproquement lorsque λ &.

Choix de λ
— Il est crucial : si λ ≈ 0 alors β̂ R ≈ β̂ M CO , si λ "grand" alors β̂ R ≈ 0.

— Le choix de λ se fait le plus souvent de façon "classique" :

1. Estimation d’un critère de choix de modèle pour toutes les valeurs de λ ;

2. Choix du λ qui minimise le critère estimé.

— Exemple : la fonction [Link] choisit la valeur de λ qui minimise l’erreur quadratique moyenne

E[(Y − mβ̂ R (λ) (X))2 ]

estimée par validation croisée.

> [Link](321)
> [Link] <- [Link](Hitters.X,Hitters$Salary,alpha=0,
+ lambda=exp(seq(0,10,length=100)))
> bestlam <- [Link]$[Link]
> bestlam
## [1] 233.8186
> plot([Link])

35
 19 19 19 19 19 19 19 19 19 19 19 19 19 19 19

220000
180000
Mean−Squared Error

140000
100000
0 2 4 6 8 10

Log(λ)

3.2 Régression Lasso

— La régression lasso consiste à minimiser le critère des moindres carrés pénalisé par la norme 1 des coefficients.

Définition [Tibshirani, 1996]

1. Les estimateurs lasso β̂ L s’obtiennent en minimisant
 2
X n Xd d
X
Yi − β0 − Xij βj  sous la contrainte |βj | ≤ t (3)
i=1 j=1 j=1

2. ou de façon équivalente
  2 

X n Xd Xd 

β̂ L = argmin Yi − β0 − Xij βj  + λ |βj | . (4)
β 
 i=1 j=1 j=1



Comparaison Ridge-Lasso
— Dans le cas où la matrice X est orthonormée, on a une écriture explicite pour les estimateurs ridge et lasso.

Propriété
Si la matrice de design X est orthonormée, alors
β̂j
β̂jR = et β̂jL = signe(β̂j )(|β̂j | − λ)+
1+λ
où β̂j est l’estimateur MCO de βj .

Commentaires
— Ridge "diminue" l’estimateur MCO de façon proportionnelle ;
— Lasso translate et tronque l’estimateur MCO (lorsque ce dernier est petit).

2 Est
Lasso
0
MCO
Ridge
−2

−4
−4 −2 0 2 4

36
Conclusion
Le lasso va avoir tendance à "mettre" des coefficients à 0 et donc à faire de la sélection de variables.

β2 β2

β̂ β̂

β1 β1

Remarque
Ces approches reviennent (d’une certaine façon) à projeter l’estimateur des MCO sur les boules unités associées à
1. la norme 2 pour la régression ridge ;
2. la norme 1 pour le lasso.

Quelques remarques
— Comme pour la régression ridge :

— on préfère souvent réduire la matrice de design avant d’effectuer la régression lasso ;

— Le choix de λ est crucial (il est le plus souvent sélectionné en minimisant un critère empirique).

— λ % =⇒ biais % et variance & et réciproquement lorsque λ &.

— MAIS, contrairement à ridge : λ % =⇒ le nombre de coefficients nuls augmente ([Bühlmann and van de Geer, 2011]).

Le coin R

> [Link] <- glmnet(Hitters.X,Hitters$Salary,alpha=1)

> par(mfrow=c(1,2))
> plot([Link],lwd=2)
> plot([Link],lwd=2,xvar="lambda")

0 6 6 11 17 19 17 12 6
50

50
0

0
Coefficients

Coefficients
−50

−50
−100

−100

0 50 100 150 200 −2 0 2 4

L1 Norm Log Lambda

37
Sélection de λ
> [Link](321)
> [Link] <- [Link](Hitters.X,Hitters$Salary,alpha=1)
> bestlam <- [Link]$[Link]
> bestlam
## [1] 17.19108
> plot([Link])

19 18 18 18 17 17 13 13 12 9 7 6 6 6 6 5 4 2

220000
Mean−Squared Error

180000
140000
100000

−2 0 2 4

Log(λ)

Résolution numérique
— Il existe plusieurs façons de résoudre le problème numérique d’optimisation lasso (ou ridge).

— Un des plus utilisé est l’algorithme de descente de coordonnées [Hastie et al., 2015].

— On considère le problème lasso

  2 

X n Xd Xd 

β̂ L = argmin Yi − β0 − Xij βj  + λ |βj |
β 
 i=1 j=1 j=1



avec les variables explicatives centrées-réduites (pour simplifier).

Descente de coordonnées
1. Initialisation : β̂0 = ȳ, β̂j = ..., j = 1, . . . , d.

2. Répéter jusqu’à convergence : Pour j = 1, . . . , d :

(j) P
(a) Calculer les résidus partiels ri = yi − k6=j xik β̂k ;
(j)
(b) Faire la régression simple des yi contre ri =⇒ β̃j ;
(c) Mettre à jour β̂j = signe(β̃j )(|β̃j | − λ)+

3. Retourner : β̂j , j = 1, . . . , d.

3.3 Variantes de ridge/lasso

Différentes pénalités
— Les approches ridge et lasso diffèrent uniquement au niveau de la pénalité ajoutée au critère des moindres
carrés.

— Norme 2 pour ridge et norme 1 pour le lasso.

— Il existe tout un tas d’autres stratégies de pénalisations.

— Nous en présentons quelques unes dans cette partie.

— On pourra consulter [Hastie et al., 2015] pour plus de détails.

38
Elastic net
— [Zou and Hastie, 2005] ont proposé de combiner les approches ridge et lasso en proposant une pénalité (appelée
elastic net) de la forme
Xd
λ ((1 − α)βj2 + α|βj |)
j=1

où α ∈ [0, 1].

— Le paramètre α définit le compromis ridge/lasso :

— α = 1 =⇒ Lasso ;
— α = 0 =⇒ Ridge ;
— Ce paramètre correspond (évidemment) à l’argument alpha de la fonction glmnet.
— Avantage : on a plus de flexibilité car la pénalité elastic net propose une gamme de modèles beaucoup plus
large que lasso et ridge ;
— Inconvénient : en plus du λ il faut aussi sélectionner le α !

Group Lasso
— Dans certaines applications, les variables explicatives appartiennent à des groupes de variables prédéfinis.
— Nécessité de "shrinker" ou sélectionner les variables par groupe.

Exemple : variables qualitatives

— 2 variables explicatives qualitatives X1 et X2 et une variable explicative continue X3 .
— Le modèle s’écrit

Y =β0 + β1 1X1 =A + β2 1X1 =B + β3 1X1 =C

+ β4 1X2 =D + β5 1X2 =E + β6 1X2 =F + β7 1X2 =G + β8 X3 + ε

muni des contraintes β1 = β4 = 0.

— 3 groupes : X1 = (1X1 =B , 1X1 =C ), X2 = (1X2 =E , 1X2 =F , 1X2 =G ) et X3 = X3 .

Définition
En présence de d variables réparties en L groupes X1 , . . . , XL de cardinal d1 , . . . , dL . On note β` , ` = 1, . . . , L le
vecteur des coefficients associé au groupe X` . Les estimateurs group-lasso s’obtiennent en minimisant le critère
n L
!2 L p
X X X
y i − β0 − Xi` β` +λ d` kβ` k2
i=1 `=1 `=1

Remarque
Puisque kβ` k2 = 0 ssi β`1 = . . . = β`d` = 0, cette procédure encourage la mise à zéro des coefficients d’un même
groupe.

Le coin R
— La fonction gglasso du package gglasso permet de faire du groupe lasso sur R.
> summary(donnees)
## X1 X2 X3 Y
## Length:200 Length:200 Min. :0.009496 Min. :-3.23315
## Class :character Class :character 1st Qu.:0.237935 1st Qu.:-0.50404
## Mode :character Mode :character Median :0.485563 Median : 0.16759
## Mean :0.483286 Mean : 0.09792
## 3rd Qu.:0.734949 3rd Qu.: 0.66918
## Max. :0.998741 Max. : 3.04377
> D <- [Link](Y~.,data=donnees)[,-1]
> model <- glmnet(D,Y,alpha=1)
> groupe <- c(1,1,2,2,2,3)
> library(gglasso)
> model1 <- gglasso(D,Y,group=groupe)
> plot(model1)

39
 6 6 5 4

0.5
0.5
4

0.4
0.4
3

0.3
0.3
5

Coefficients

Coefficients

0.2
0.2

0.1
0.1

0.0
0.0
2

−0.2

−0.2
6

−10 −8 −6 −4 −11 −9 −7 −5

Log Lambda Log Lambda

Remarque
Les coefficients s’annulent par groupe lorsque λ augmente (graphe de droite).

Sparse group lasso

— La norme 2 de la pénalité group-lasso implique que, généralement, tous les coefficients d’un groupe sont tous
nuls ou tous non nuls.
— Dans certains cas, il peut être intéressant de mettre de la sparsité dans les groupes aussi. Comment ?
— En ajoutant la norme 1 dans la pénalité.
Pénalité sparse group lasso
L
X
λ [(1 − α)kβ` k2 + αkβ` k1 ] .
`=1

— Sur R : package SGL.

Fused lasso
— Utile pour prendre en compte la spatialité des données.
— Idée : deux coefficients successifs doivent être proches.
Pénalité fused lasso
d
X d
X
λ1 |βj | + λ2 |βj+1 − βj |
j=1 j=2

qui peut se re-paramétrer en

d
X
λ |βj+1 − βj |.
j=2

— Sur R : package genlasso.

3.4 Discrimination binaire

Discrimination binaire
— Les méthodes ridge et lasso ont été présentées dans un cadre de régression linéaire.

— Ces techniques d’adaptent directement à la régression logistique Y = {−1, 1}.

— Les pénalités sont identiques.

— Seul changement : le critère moindre carré est remplacé par la déviance =⇒ ce qui revient à minimiser l’opposé
de la vraisemblance plus la pénalité.

40
Lasso et Ridge pour la logistique
Définition
On note ỹi = (yi + 1)/2.
— On appelle estimateur ridge en régression logistique l’estimateur
 
 X n d
X 
β̂ R = argmin − (ỹi xti β − log(1 + exp(xti β))) + λ βj2 .
β 
i=1

j=1

— On appelle estimateur lasso en régression logistique l’estimateur

 
 X n d
X 
β̂ L = argmin − (ỹi xti β − log(1 + exp(xti β))) + λ |βj | .
β 
i=1

j=1

Le coin R
— Pour faire du ridge ou lasso en logistique, il suffit d’ajouter l’argument family=binomial dans glmnet.

— Tout reste identique pour le reste (tracé du chemin des coefficients, choix du λ...).
— Exemple : données SAheart
> head(SAheart)
## sbp tobacco ldl adiposity famhist typea obesity alcohol age chd
## 1 160 12.00 5.73 23.11 Present 49 25.30 97.20 52 1
## 2 144 0.01 4.41 28.61 Absent 55 28.87 2.06 63 1
## 3 118 0.08 3.48 32.28 Present 52 29.14 3.81 46 0
## 4 170 7.50 6.41 38.03 Present 51 31.99 24.26 58 1
## 5 134 13.60 3.50 27.78 Present 60 25.99 57.34 49 1
## 6 132 6.20 6.47 36.21 Present 62 30.77 14.14 45 0

— On obtient les chemins de régularisation ridge et lasso avec les commandes suivantes :
> SAheart.X <- [Link](chd~.,data=SAheart)
> [Link] <- glmnet(SAheart.X,SAheart$chd,family="binomial",alpha=0)
> [Link] <- glmnet(SAheart.X,SAheart$chd,family="binomial",alpha=1)
> plot([Link],xvar="lambda")

9 9 9 9 9 8 8 7 6 5 1
0.8

0.8
0.6

0.6
Coefficients

Coefficients
0.4

0.4
0.2

0.2
0.0
0.0

−4 −2 0 2 4 −7 −6 −5 −4 −3 −2

Log Lambda Log Lambda

41
4 Support vector machine
Cadre et notation
— Discrimination binaire : Y à valeurs dans {−1, 1} et X = (X1 , . . . , Xd ) dans Rd .
— Les extensions multi-classes et régression seront présentées à la fin de cette partie.

Objectif
— Estimer la fonction de score S(x) = P(Y = 1|X = x) ;
— En déduire une règle de classification g : Rd → {−1, 1}.

Règles linéaires
— Elles consistent à séparer l’espace des X par un hyperplan.
— On classe ensuite 1 d’un coté de l’hyperplan, -1 de l’autre coté.

Mathématiquement
— On cherche une combinaison linéaire des variables w1 x1 + . . . + wd xd .
— Règle associée : 
1 si w1 x1 + . . . + wd xd ≥ 0
g(x) =
−1 sinon.

Exemple 1 : régression logistique

— Modèle :
p(x)
logit = β0 + β1 x 1 + . . . + βd x d
1 − p(x)
où p(x) = P(Y = 1|X = x).
— Règle de classification : 
1 si p(x) ≥ 0.5
g(x) =
−1 sinon.
— équivalent à 
1 si β0 + β1 x1 + . . . + βd xd ≥ 0
g(x) =
−1 sinon.

Exemple 2 : LDA
— Modèle : L(X|Y = k) = N (µk , Σ), k = 0, 1.
— Règle de classification : 
1 si p(x) ≥ 0.5
g(x) =
−1 sinon.
— équivalent à
si c + xt Σ−1 (µ1 − µ0 ) ≥ 0

1
g(x) =
−1 sinon.

Illustration avec p = 2
— On considère les données suivantes :

42
 1.00

0.75

X2
0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1
— On compare les prévisions logistique et lda.

lda logit
1.00

0.75
X2

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00
X1

Remarque
On retrouve bien la linéarité (et la proximité) de ces deux méthodes.
— Ces approches linéaires s’obtiennent à partir d’un modèle statistique

— sur la loi de Y sachant X pour la logistique ;

— sur la loi de X sachant Y pour la discriminante linéaire.

— L’approche SVM repose sur le calcul direct du "meilleur" hyperplan séparateur qui sera déterminé à partir
d’algorithmes d’optimisation.

4.1 SVM - cas séparable

Bibliographie
En plus des documents cités précédemment, cette partie s’appuie sur les diapos de cours de
— Magalie Fromont, Apprentissage statistique, Université Rennes 2 ([Fromont, 2015]).

— Jean-Philippe Vert, Support vector machines and applications in computational biology, disponible à l’url
http: // cbio. ensmp. fr/ ~jvert/ svn/ kernelcourse/ slides/ kernel2h/ kernel2h. pdf

Remarque
Les aspects techniques ne seront pas présentés ici, on pourra en trouver dans la partie 2.4 du tutoriel.

43
Présentation
— L’approche SVM [Vapnik, 2000] peut être vue comme une généralisation de "recherche d’hyperplan optimal".

Cas simple
Les données (x1 , y1 ), . . . , (xn , yn ) sont dites linéairement séparables si il existe (w, b) ∈ Rd × R tel que pour tout i :
— yi = 1 si hw, xi i + b = wt xi + b > 0 ;
— yi = −1 si hw, xi i + b = wt xi + b < 0.

hw, xi + b = 0

Vocabulaire
— L’équation hw, xi + b définit un hyperplan séparateur de vecteur normal w.
— La fonction signe(hw, xi + b) est une règle de discrimination potentielle.

Problème
Il existe une infinité d’hyperplans séparateurs donc une infinité de règles de discrimination potentielles.

Solution
[Vapnik, 2000] propose de choisir l’hyperplan ayant la marge maximale.

44
 hw, xi + b = 0

1
M = kwk

marge

1
M = kwk

Le problème d’optimisation
— On veut trouver l’hyperplan de marge maximale qui sépare les groupes.

Hyperplan séparateur optimal

Solution du problème d’optimisation sous contrainte :
— Version 1 :
max M
w,b,kwk=1

sous les contraintes yi (wt xi + b) ≥ M, i = 1, . . . , n.

— Version 2 :
1
min kwk2
w,b 2

sous les contraintes yi (wt xi + b) ≥ 1, i = 1, . . . , n.

Solutions
— On obtient
n
X
w? = αi? yi xi .
i=1

où les αi? sont des constantes positives qui s’obtiennent en résolvant le dual du problème précédent.
— De plus, b? s’obtient en résolvant
αi? [yi (xti w + b) − 1] = 0
pour un αi? non nul.

Remarque
w? s’écrit omme une combinaison linéaire des xi .

Vecteurs supports
Propriété (conditions KKT)
αi? [yi (xti w? + b) − 1] = 0, i = 1, . . . , n.

Conséquence (importante)
— Si αi? 6= 0 alors yi (xti w? + b) = 1 et xi est sur la marge.
— w? se calcule uniquement à partir de ces points là.
— Ces points sont appelés les vecteurs supports de la SVM.

45
Représentation
hw? , xi + b? = 0

M = kw1? k

α?
i > 0

marge

M = kw1? k

Le coin R
— La fonction svm du package e1071 permet d’ajuster des SVM.
1.00

0.75
X1

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75
X2
> library(e1071)
> [Link] <- svm(Y~.,data=df,kernel="linear",cost=10000000000)

La fonction svm
— Les vecteurs supports :
> [Link]$index
## [1] 6 14 12

— [Link]$coefs = αi? yi pour chaque vecteur support

> [Link]$coefs
## [,1]
## [1,] 1.898982
## [2,] 1.905497
## [3,] -3.804479

— On peut en déduire l’hyperplan séparateur

> w <- apply([Link]$coefs*df[[Link]$index,2:3],2,sum)
> b <- -[Link]$rho
> w
## X1 X2
## -0.5470382 0.5427583
> b
## [1] -0.4035113

On peut ainsi visualiser

— les vecteurs supports ;
— l’hyperplan séparateur ;

46
 — la marge.
1.00

0.75

X1
0.50 0

0.25

0.00
0.00 0.25 0.50 0.75
X2

— La fonction plot donne aussi une représentation de l’hyperplan séparateur.

> plot([Link],data=df,fill=TRUE,grid=100)

SVM classification plot

o o o
o oo
0.8 o

1
o
0.6
o ox
X1

x o
0.4 o
o
o

0
o
0.2
o
x o
0.2 0.4 0.6 0.8

X2

4.2 SVM : cas non séparable

Problème
Dans la vraie vie, les données ne sont (quasiment) jamais linéairement séparables...

47
 hw, xi + b = 0

1
M = kwk

marge
Idée
M = 1 Autoriser certains points
kwk
1. à être bien classés mais à l’intérieur de la
hw, xi + b = 0
marge ;
2. et/ou à être mal classés.

1
M = kwk

marge

1
M = kwk

Slack variables
Rappel : cas séparable
1
min kwk2
w,b 2

sous les contraintes yi (wt xi + b) ≥ 1, i = 1, . . . , n.

— Les contraintes yi (wt xi + b) ≥ 1 signifient que tous les points se trouvent en dehors de la frontière définie par
la marge ;
— Cas non séparable : le problème ci-dessus n’admet pas de solution !

Variables ressorts
On introduit des variables ressorts (slack variables) positives ξ1 , . . . , ξn telles que yi (wt xi + b) ≥ 1 − ξi . 2 cas sont
à distinguer :
1. ξi ∈ [0, 1] =⇒ bien classé mais dans la région définie par la marge ;
2. ξi > 1 =⇒ mal classé.

— Bien entendu, on souhaite avoir le maximum de variables ressorts ξi nulles ;

— Lorsque ξi > 0, on souhaite que ξi soit le plus petit possible.

Cas non séparable : problème d’optimisation (primal)

— Il s’agit de minimiser en (w, b, ξ)
n
1 X
kwk2 +C ξi
2 i=1

yi (wt xi + b) ≥ 1−ξi

sous les contraintes
ξi ≥ 0, i = 1, . . . , n.

48
 — C > 0 est un paramètre à calibrer (paramètre de coût).
— Le cas séparable correspond à C → ∞.

— Les solutions de ce nouveau problème d’optimisation s’obtiennent de la même façon que dans le cas séparable
(Lagrangien, problème dual...).
— L’hyperplan optimal est défini par
Xn
?
w = αi? yi xi
i=1

et b est solution de yi (hw , xi i + b ) = 1 pour tout i tel que 0 < αi? < C.
? ? ?

Vecteurs supports
— Les xi tels que αi? > 0 sont les vecteurs supports ;
— On en distingue 2 types :
1. ceux sur la frontière définie par la marge : ξi? = 0 ;
2. ceux en dehors : ξi? > 0 et αi? = C.
— Les vecteurs non supports vérifient αi? = 0 et ξi? = 0.

hw? , xi + b? = 0

ξ1?

ξ3? 1
M =
ξ2? kw? k

ξ5?

ξ4? 0 < α?i ≤ C, ξi? = 0

marge

1
M = kw? k

Le coin R
— On utilise la même fonction que dans le cas séparable (svm du package e1071) ;
— L’argument cost correspond à la constante de régularisation C.
1.00

0.75

Y
X1

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75
X2

> mod.svm1 <- svm(Y~.,data=df1,kernel="linear",cost=1000)

> mod.svm1$index
## [1] 6 13 14 10 12 15

Visualisation de l’hyperplan séparateur

49
> plot(mod.svm1,data=df1,fill=TRUE,grid=100)

SVM classification plot

o o o
o oo
0.8 o

1
x
0.6
o ox

X1
x x
0.4 o
o
o

0
o
0.2
x
x o
0.2 0.4 0.6 0.8

X2

Choix de C
Ce paramètre régule le compromis biais/variance de la svm :
— C & : la marge est privilégiée et les ξi % =⇒beaucoup d’observations dans la marge ou mal classées (et donc
beaucoup de vecteurs supports).
— C %=⇒ ξi & donc moins d’observations mal classées =⇒ meilleur ajustement mais petite marge =⇒ risque
de surajustement.

Conclusion
Il est donc très important de bien choisir ce paramètre.
— Le choix est souvent effectué de façon "classique" :

1. On se donne un critère de performance (taux de mal classés par exemple) ;

2. On estime la valeur du critère pour différentes valeurs de C ;

3. On choisit la valeur de C pour laquelle le critère estimé est minimum.

— La fonction [Link] permet de choisir C en estimant le taux de mal classés par validation croisée. On peut
aussi (bien entendu) utiliser la fonction train du package caret.

Un exemple

1.00

0.75

Y
V1

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.00
V2

> mod.svm1 <- svm(Y~.,data=df3,kernel="linear",cost=0.000001)

> mod.svm2 <- svm(Y~.,data=df3,kernel="linear",cost=0.1)
> mod.svm3 <- svm(Y~.,data=df3,kernel="linear",cost=5)

50
 SVM classification plot SVM classification plot SVM classification plot
o o oooooooooo ox x xx xxxxxxxxxxxxxxxxxxxxx xxxxx xxxx x o o oooooooooo oo o x ooooooxoooo x
ooo xxxx xxxx x
o o oooooooooo oo o x ooooooxoooo x
ooo oxxx xxxx x
ooooo
oooo x
xx xx xxxxxxxxxx xxxxx xxxxxxxxxxx xxxx xxxxxx ooooo
oooo oo oxo o
o oo oooooooooxo oo ooooo
oooo oo oxo o
o oo oooooooooxo oo
xxo oo xxxx xxxxxx xoo oo ooooo oooooo oo xox oo oo oooooooxxx xxxx xxxxx xoo oo ooooo oooooo oo xox oo oo ooooooooxx xxxx xxxxx
ooo ooooxxxxxxx x xx x x xxxxxx xxxxxxxx xxxxx
x x
ooo ooooooooo o xo x xx xx ooo ooooooooo o xo ox xx xx
o x xx xxx x xx x xxx o o ooo oooo o o oo ooo xxxxxxxxxx x xx xx o o ooo oooo o o oo ooo xxxxxxxxxx x xx xx
0.8 o xx x x x xxx
oo xxxxxxxxxxxxxxxxx xx xx xxx x xxxxxxxxxxx xxxxxxx x xxxxx 0.8 o o o oooo
oo oxoooooooo x x o oooooo oxxxxxx xxxxx xxxx x
o
0.8 o o o oooo
oo oxoooooooo x x o oooooo ooxoxxxx xxxxx xxxx o
o
o x x xxxx xx xxxxxxxxxxxx x xx x xxxxx oooooxoo o ooo o ooxxxxxxxxxxx xx xx xxxxoooo oooooxoo o ooo o ooxooxxoxxxxx xx xx xxoooooo

1
o o
x x x x xx x x xxxxxx xxx x o o o oo o oo o ooo xxxxxxxxxxxxxxxx o oo o o o o oo o oo o ooo xxxxxxxxxxxxxxxx o oo o
xx xxx xx xx xxxx xxx xxxxx xxxxxxx x xx xxx oo xoo oo oooooxo ox xxxxx xxxxxxx x o o oo oo xoo oo oooooxo ox oxxxx xxxxxxx x o o oo
x xxxxxxxx xxxxxx xxxxxxxxxxxxxxxxxxxxx xxxx xxxx xxxx o oxoooooo ooooooooo oooo o xx xx x x x x o oxoooooo ooooooooo oooo o ooxxx x x xxx ooo oo oo
0.6 xxx x xxxxxxxxxxxxx xxxxxx x xxxxxxxxxx xxx xxxxxx xxx x
x 0.6 ooo oooo o oooxoooxoxxxxxxxxxx xxxxx xxxx xxo oo oooo 0.6 ooo oooo o oooxoooooooxx xxxx xxx o xo
oxooooo o oo xxxx x xxxx xx xooooooo o oo oo oxooooo o oo xxxx x xxxx xxxxxoo ooooo o o ooo
oo
x xxxxx xxxx xxxxx xxxx xxxx xxxxxxxxx xxxxxxxxxx o ooo
o o ox x
xx x o o x o
oo o o ooo
o o ox x
oo x o
oo x oo o o
xxxxxxxxxxxx xx xx xx x xxxxxxxxxxxxx xxx xxx xxx x xxx x ooo oxoooo x xxxxxxxxxxxxxxxx xxx oooo ooo oooxo o
oooooxoooo ooo oxoooo x xxxxxxxxxxxxxxxx xoo oooo ooo oooxo o
oooooxoooo
V1

V1

V1
o o
xxxx xxxxxxxx xx x x x xxxx xx x xx ooo oooooooo xx x x x oxxo oo ooo o ooo oooooooo ox x x x ooxo oo ooo o
0.4 xxxx xxx xxxxxxxx xxxxxx xxxxx xxxxx x x xxxxxx xxx xxx 0.4 oooo ooo oxoxxxxx xxxxxx xxxxx xxxxx o o ooo o xx oo 0.4 oooo ooo oooooxxx xxxxxx xxxxx xxoxx o o ooo o xx oo
xx x xx xx x xxxxx x xx xxxxxx xx o ooo xxxxxx x xx xxx o o ooo o o o ooo ox o ooo xxxxxx x xx oxx o o ooo o o o oo o ox
xx x xxxxx x x o o o o o
o
o o o o o o o
oo o
xxxxx xxxxxxx xxx xxxxxx xxxxxxxxxxxx xx xxxxxxxxxxxxxxxx xxx x x o
oooo o xx x
oox x x xxxxxxoo ooo x ooooooooooooo ooo o o o
oooo o xx x
oox x x xxxoxooo ooo x ooooooooooooo ooo o o
x x xx x x xxxxxxxxx x xxxxx xxxxxxxx x oxx xxxxx xxx xxx xxxo o o
oo xooo o o o o
oo oxooo o oxooxxxx xxx xxx oxoo o o
ooooxooo o oxo o oooo o
oo ooxo
x x ooox x xx xxx x o o x ooo o oooo x xx xxx x o

0
xxx x xx xx x
xx xx x x xx xxx xxxxxxxxx xxxxxxxxx xx x ooo ooo ooo oo o o
0.2 x x xxxx xxx xxxx xxxx x xxxx x x x 0.2 o xxxxxxxx xxxxxxxx xooo o oooooooo ooooooooo oooo 0.2 o xxxxxxxx xxxxxxxx xooo o oooooooo ooooooooo oooo
xxx xx x x x x x x x x xxx oo oo oo oxxooooo o o xxx oo oo oo oxxooooo o o
xxxxxxxxx xxxxxx xxxxxxxxx xxxxxxx xxxxxxxxxxx xxxx xxx xx xxxxxxxxx xxxxxx xooooo o
ooo o ooo xo oooooooo o o xxxxxxxxx xxxxxo xooooo o
ooo o ooo xo oooooooo o
oooo oo
o xoo o oo o oooo oo
o xoo o oo o o o
xxxx xxxxx x x xxx x x x x xx xx x
xx ooxooo o o o oooooo oo
o o o o xx o o
xx ooxooo o o o oooooo o o o o
xxx xxx xxxx xxxxxx xxxxxxxxxxx xxxxxxxx xxxx xxxxx xxxxxxx xxx xxx xxxx xo ooooo ooo oo xo oo
ooxo ooo oo oooo oo
ooooo xxx xxx xxxx ooo ooo ooo oo xo oo
ooxo ooo oo oooo oo
ooooo

0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8

V2 V2 V2

> mod.svm1$nSV
## [1] 480 480
> mod.svm2$nSV
## [1] 190 190
> mod.svm3$nSV
## [1] 166 165

Choix de C avec tune

> [Link](1234)
> [Link] <- tune(svm,Y~.,data=df3,kernel="linear",
+ ranges=list(cost=c(0.001,0.01,1,10,100,1000)))
> summary([Link])
## Parameter tuning of 'svm':
## - sampling method: 10-fold cross validation
## - best parameters:
## cost
## 1
##
## - best performance: 0.075
##
## - Detailed performance results:
## cost error dispersion
## 1 1e-03 0.127 0.07087548
## 2 1e-02 0.080 0.03944053
## 3 1e+00 0.075 0.03439961
## 4 1e+01 0.075 0.03439961
## 5 1e+02 0.075 0.03439961
## 6 1e+03 0.075 0.03439961

> bestmod <- [Link]$[Link]

> summary(bestmod)
##
## Call:
## [Link](method = svm, train.x = Y ~ ., data = df3, ranges = list(cost = c(0.001,
## 0.01, 1, 10, 100, 1000)), kernel = "linear")
##
##
## Parameters:
## SVM-Type: C-classification
## SVM-Kernel: linear
## cost: 1
##
## Number of Support Vectors: 336
##
## ( 168 168 )
##
##
## Number of Classes: 2
##
## Levels:
## 0 1

Approche tune_grid de tidymodels

1. Initialisation du workflow :

51
 > library(tidymodels)
> tune_spec <-
+ svm_poly(cost=tune(),degree=1,scale_factor=1) %>%
+ set_mode("classification") %>%
+ set_engine("kernlab")
> svm_wf <- workflow() %>%
+ add_model(tune_spec) %>%
+ add_formula(Y ~ .)

2. Ré-échantillonnage et grille de paramètres :

> [Link](12345)
> re_ech_cv <- vfold_cv(df3,v=10)
> grille_C <- tibble(cost=c(0.001,0.01,1,10,100,1000))

3. Calcul des erreurs :

> [Link](123)
> [Link] <- svm_wf %>%
+ tune_grid(
+ resamples = re_ech_cv,
+ grid = grille_C,
+ metrics=metric_set(accuracy))

4. Visualisation des résultats :

> [Link] %>% collect_metrics() %>% dplyr::select(-7)
## # A tibble: 6 x 6
## cost .metric .estimator mean n std_err
## <dbl> <chr> <chr> <dbl> <int> <dbl>
## 1 0.001 accuracy binary 0.864 10 0.0163
## 2 0.01 accuracy binary 0.915 10 0.00778
## 3 1 accuracy binary 0.927 10 0.00700
## 4 10 accuracy binary 0.927 10 0.00775
## 5 100 accuracy binary 0.929 10 0.00809
## 6 1000 accuracy binary 0.929 10 0.00809

5. Sélection du meilleur paramètre

> best_C <- [Link] %>% select_best()
> best_C
## # A tibble: 1 x 2
## cost .config
## <dbl> <chr>
## 1 100 Preprocessor1_Model5

6. Ajustement de l’algorithme final :

> final_svm <-
+ svm_wf %>%
+ finalize_workflow(best_C) %>%
+ fit(data = df3)

4.3 SVM non linéaire : astuce du noyau

— Les solutions linéaires ne sont pas toujours intéressantes.

brutes transf
1.0 1.00

0.5 0.75
Y
X2

0.0 0.50 0

−0.5 0.25

−1.0 0.00
−1.0 −0.5 0.0 0.5 0.00 0.25 0.50 0.75
X1

52
Idée
Trouver une transformation des données telle que les données transformées soient linéairement séparables.

Noyau
Définition 4.1. Soit Φ : X → H une application qui va de l’espace des observations X dans un Hilbert H. Le
noyau K entre x et x0 associé à Φ est le produit scalaire entre Φ(x) et Φ(x0 ) :

K :X ×X →R
(x, x0 ) 7→ hΦ(x), Φ(x0 )iH .

Exemple
Si X = H = R2 et ϕ(x1 , x2 ) = (x21 , x22 ) alors

K(x, x0 ) = (x1 x01 )2 + (x2 x02 )2 .

L’astuce noyau
— L’astuce consiste donc à envoyer les observations xi dans un espace de Hilbert H appelé espace de représen-
tation ou feature space...

— en espérant que les données (Φ(x1 ), y1 ), . . . , (Φ(xn ), yn ) soient (presque) linéairement séparables de manière
à appliquer une svm sur ces données transformées.

Remarque
1. Beaucoup d’algorithmes linéaires (en particulier les SVM) peuvent être appliqués sur Φ(x) sans calculer
explicitement Φ ! Il suffit de pouvoir calculer le noyau K(x, x0 ) ;
2. On n’a pas besoin de connaitre l’espace H ni l’application Φ, il suffit de se donner un noyau K !

SVM dans l’espace original

— Le problème dual consiste à maximiser
n n n
X 1 XX
LD (α) = αi − αi αk yi yk hxi , xk i
i=1
2 i=1
k=1

0P≤ αi ≤ C, i = 1, . . . , n
sous les contraintes n
i=1 αi yi = 0.

— La règle de décision s’obtient en calculant le signe de

n
X
f (x) = αi? yi hxi , xi + b? .
i=1

SVM dans le feature space

— Le problème dual consiste à maximiser
n n n
X 1 XX
LD (α) = αi − αi αk yi yk hΦ(xi ), Φ(xk )i
i=1
2 i=1
k=1

0P≤ αi ≤ C, i = 1, . . . , n
sous les contraintes n
i=1 αi yi = 0.

— La règle de décision s’obtient en calculant le signe de

n
X
f (x) = αi? yi hΦ(xi ), Φ(x)i + b? .
i=1

53
SVM dans le feature space avec un noyau
— Le problème dual consiste à maximiser
n n n
X 1 XX
LD (α) = αi − αi αk yi yk K(xi , xk )
i=1
2 i=1
k=1

0P≤ αi ≤ C, i = 1, . . . , n
sous les contraintes n
i=1 αi yi = 0.
— La règle de décision s’obtient en calculant le signe de
n
X
f (x) = αi? yi K(xi , x) + b? .
i=1

Conclusion
— Pour calculer la svm, on n’a pas besoin de connaitre H ou Φ, il suffit de connaitre K !
Questions
Qu’est-ce qu’un noyau ? Comment construire un noyau ?
Théorème 4.1 ([Aronszajn, 1950]). Une fonction K : X × X → R est un noyau si et seulement si elle est
(symétrique) définie positive, c’est-à-dire ssi
1. K(x, x0 ) = K(x0 , x) ∀(x, x0 ) ∈ X 2 ;
2. ∀(x1 , . . . , xN ) ∈ X N et ∀(a1 , . . . , aN ) ∈ RN
N X
X N
ai aj K(xi , xj ) ≥ 0.
i=1 j=1

Exemple
brutes transf
1.0 1.00

0.5 0.75
Y
X2

0.0 0.50 0

−0.5 0.25

−1.0 0.00
−1.0 −0.5 0.0 0.5 0.00 0.25 0.50 0.75
X1

Si
Φ: R2 → R3
√
(x1 , x2 ) 7→ (x21 , 2x1 x2 , x22 )
alors K(x, x0 ) = (xt x0 )2 (noyau polynomial de degré 2).

Exemples de noyau
1. Linéaire (vanilladot) : K(x, x0 ) = hx, x0 i.
2. Polynomial (polydot) : K(x, x0 ) = (scalehx, x0 i + offset)degree .
3. Gaussien (Gaussian radial basis function ou RBF - rbfdot)
K(x, x0 ) = exp(−sigmakx − x0 k2 .
4. Laplace (sur R) : K(x, x0 ) = exp(−sigmakx − x0 k).
5. ...
Remarque
Les paramètres correspondent aux noyaux proposés par la fonction ksvm de kernlab (voir [Karatzoglou et al., 2004]).

54
Commentaires
1. En l’absence d’information a priori le noyau radial est préconisé.
2. Procédure d’optimisation pour sigma proposé dans ksvm.

Remarque
N’importe quelle fonction définie positive fait l’affaire... Possibilité de construire des noyaux (et donc de faire des
svm) sur des objets plus complexes (courbes, images, séquences de lettres...).

Le coin R - exemple 1
— Argument kernel dans la fonction svm.
> svm(Y~.,data=donnees,cost=1,kernel="linear")
> svm(Y~.,data=donnees,cost=1,kernel="polynomial",degree=2)
> svm(Y~.,data=donnees,cost=1,kernel="radial",gamma=1)
SVM classification plot SVM classification plot SVM classification plot
o o o o o o
x x o o x o o o o
o x o x o x
o o o
x x o o o x
0.5 x 0.5 x 0.5 x
x x x
1

1
x o o
x x o x x x x x o
x o o x o o
x x x o o x x o x
0.0 x x x 0.0 x o o 0.0 o o x
X1

X1

X1
x x o o o o
x x x x x x o x x o o x
x xx o o xx o x xx o
x x o o o o
−0.5 −0.5 −0.5
0

0
x x o o x o x x o
o o o
xx oo ox
x o o o o o
o o o
−0.5 0.0 0.5 −0.5 0.0 0.5 −0.5 0.0 0.5

X2 X2 X2

Le coin R - exemple 2

Y
X1

0 0

−1

−1 0 1 2
X2

> svm(Y~.,data=donnees,kernel="linear",cost=1)
> svm(Y~.,data=donnees,kernel="polynomial",degree=2,cost=1)

SVM classification plot SVM classification plot

oo o o oo o o
o o
1.5 o ooooo 1.5 o ooooo
o o o o
o oo ooooooo oo o o o
o oo o oo ooooooo oo o o o
o oo
o o oooo ooo o o o oooo ooo o
1.0 o oooo o o oooooooooo 1.0 o oooo o o oooooooooo
o oo o oo o oo o o oo o oo o oo o
1

oo o ooo oo o oo o ooo oo o
ooo oo o oo o o ooo oo o oo o o
0.5 o o 0.5 o o
oo oo
o o
X1

X1

0.0 0.0
o o
o o
−0.5
o o o o oo
o oo
−0.5
o x o o oo
o oo
x xxx x xx oo o o ox ooo x oo oo o o
oxox oxo ooooo ooooo oxox ooo ooooo ooooo
−1.0 x x x xxxxxo xoo o oooo oo −1.0 o o x oxxxoo xoo o oooo oo
0

x xo xo xx o o oo o o oo ooo xo o o oo o
xx x o o oooo oo o ox x o o oooo oo o
xx oo oo xo oo oo
−1.5 o ox o −1.5 o xx o
o o
o o o o

−1.5 −0.5 0.5 1.5 −1.5 −0.5 0.5 1.5

X2 X2

Le package kernlab
— Il propose un choix plus large de noyau.
> library(kernlab)
> [Link] <- ksvm(Y~.,data=donnees,kernel="polydot",
+ kpar=list(degree=2),C=0.001)
> plot([Link])

55
 SVM classification plot

1.5
1.0
1.0

0.5 0.5

X1
0.0
0.0
−0.5

−1.0 −0.5

−1.5
−1.0
−1.5 −1.0 −0.5 0.0 0.5 1.0 1.5

X2

4.4 Scores et probabilités

— Jusqu’à présent nous avons utiliser la SVM uniquement pour classer :
Pn
— 1 si on est d’un coté de l’hyperplan =⇒ i=1 αi? yi K(xi , x) + b? ≥ 0 ;
Pn
— -1 si on est de l’autre coté =⇒ i=1 αi? yi K(xi , x) + b? < 0.
— Rappel : dans le cas linéaire la fonction
n
X
f (x) = αi? yi K(xi , x) + b?
i=1

mesure la distance entre x et l’hyperplan séparateur.

— Conclusion : cette fonction peut être utilisée comme un score, puisque sa valeur (absolue) traduit une confiance
que l’on a dans la prévision.

Probabilités
— La valeur de f (x) est difficilement interprétable en tant que telle.
— Il peut être intéressant de la "ramener” entre 0 et 1 pour l’interpréter comme une estimation de P(Y = 1|X =
x).

Une solution
— Considérer un modèle logistique :
1
P(Y = 1|X = x) =
1 + exp(af (x) + b)
— et d’estimer a et b par maximum de vraisemblance sur les données (f (xi ), yi ), i = 1, . . . , n.

Le coin R
> [Link] <- svm(Y~.,data=df,kernel="linear",cost=10000000000,probability=TRUE)
> plot([Link],data=df,fill=TRUE,grid=100)

SVM classification plot

o
o

0.8
o
1

o o
0.6
x o
o
o
X1

0.4 o
o
o
o
0

0.2 x
o o o
o

x o
0.2 0.4 0.6 0.8

X2

56
 — Nouvelle observation :
> newX <- [Link](X1=0.2,X2=0.6)

— Calcul du score et de la proba :

> predict([Link],newdata=newX,[Link] = TRUE,probability=TRUE)
## 1
## 0
## attr(,"[Link]")
## 0/1
## 1 36.44796
## attr(,"probabilities")
## 0 1
## 1 0.9770206 0.02297939
## Levels: 0 1

— On peut retrouver cette proba avec :

> a <- [Link]$probA
> b <- [Link]$probB
> 1/(1+exp(a*36.44796+b))
## [1] 0.9770206

4.5 Compléments : SVM multi-classes et SVR

Cible multi-classes ou quantitative
— On a abordé ici uniquement le problème de la classification binaire : yi ∈ {−1, 1}.

— Les SVM se généralisent aux cas multi-classes : yi ∈ {1, . . . , M }

— et à la régression : yi ∈ R.

4.5.1 SVM multiclasses

— On suppose ici que yi ∈ {1, . . . , M }
— Il existe plusieurs approches pour généraliser les SVM à ce contexte, notamment :
— One against one
Idée
Faire une SVM binaire sur toutes les paires (j, k) ∈ {1, . . . , M }2 avec j 6= k et choisir le groupe qui gagne le
plus souvent.
— One against all
Idée
Faire une SVM binaire de chaque groupe contre les autres et choisir le groupe qui a la "plus belle victoire”.

One against one

Algorithme
1. Pour chaque paire (j, k), faire la SVM binaire avec uniquement les individus des groupes k et j ;
2. On obtient ainsi M (M − 1)/2 règles "linéaires" fj,k (x).
3. On calcule pour j = 1, . . . , M X
V (j) = signe (fj,k (x))
k6=j

qui représente le nombre de fois où on a voté j contre les autres groupes.

4. On classe un nouvel individu x dans le groupe qui a remporté le plus de suffrage :

f (x) = argmax V (j).

j

57
One against all
Algorithme
1. Faire une SVM binaire avec tous les individus de chaque groupe contre les autres.
2. On obtient ainsi M règles "linéaires" fj (x) (groupe j contre les autres).
3. On classe un nouvel individu dans la classe qui a le score le plus élevé :
f (x) = argmax fj (x).
j

Comparaison
— M SVM binaire avec l’approche one against all contre M (M − 1)/2 avec le one against one mais...
— moins d’individus dans les one against one.
— Risque de déséquilibre plus fort avec le one against all (mais généralement plus rapide).
— Comme dans le cas binaire, il faut sélectionner le paramètre de complexité, le noyau, les paramètres du
noyau...

Le coin R
— L’approche one against one est plus souvent utilisée.
— C’est le cas par défaut avec svm de e1071 et ksvm de kernlab.

Exemple
1.0

0.5
Y
1
X1

0.0
2

3
−0.5

−1.0
−1.0 −0.5 0.0 0.5 1.0
X2

SVM linéaire multi classes

> multi1 <- svm(Y~.,data=df,cost=10,kernel="linear")
> plot(multi1,data=df,grid=100)

SVM classification plot

o o o o o o o oo
o o
x x oo o o o
o o oo o o o o
o o o o
oo o o o o
3

xo
o o o x xx o o oo o
0.5 o o o
o o o
ooo o x o o
o o oo o o o
o
o o
o
o o o x o o o
o oo o x o
o o o
o o
X1

0.0 x
2

oo o o x oo o
o o o
o o oo o x x xo x x
o o oo o o
o o x o o x
x o o o
−0.5 o
oo x oo o o oo x x o oo o
o o o x
o o o oo o
o ox x o o o oo
1

oo oo
x
o oo x x
xx o o o o o oo o o
o ooo o o o o o
o

−0.5 0.0 0.5

X2

58
SVM non linéaire multi classes
— Il "suffit" d’utiliser un noyau.
> multi2 <- svm(Y~.,data=df,cost=0.1,kernel="sigmoid")
> plot(multi2,data=df,grid=100)

SVM classification plot

o x x x x x x xx
x x
x x xx x x o
o o ox x x x x
x x xo
oo o o x x

3
xx
o o x x xx x x oo o
0.5 o o o
x o o
ooo o x x o
o o xx x x o
o
o o
oo o x x x x x x x
x x x x
x x x
x x
X1
0.0 x

2
xx x x x xx x
x x x
x x xx x x x
xx x x x x xx x x
x x x x x x
x x x x
−0.5 x
xx x xx x x xx x x x xx x
x
x x x x xx xx x
x xx x x x x xx

1
xx oo
x
o oo x x
xx x x x x x xx x x
x ooo o x x x x
x

−0.5 0.0 0.5

X2

4.5.2 Support vector regression (SVR)

— On suppose ici que les yi sont dans R.

— On ne va plus chercher l’hyperplan qui sépare au mieux les groupes mais

— l’hyperplan (w, b) qui "approche au mieux” les valeurs yi

|hw, xi i + b − yi | petits.

Comparaison avec les MCO

— Approche MCO (rappel) : on cherche (w, b) qui minimise
n
X
(yi − hw, xi i − b)2
i=1

— Approche SVR : on veut

1. tous les points à distance de moins de ε de (w, b) ;
2. (w, b) de marge maximale (kwk minimale).

Optimisation SVR
On va chercher à minimiser la norme de w en se fixant comme contrainte que les yi ne soient pas "trop loin” de
l’hyperplan :
1
min kwk2
w,b 2

sous les contraintes |yi − hw, xi i − b| ≤ ε, i = 1, . . . , n,

où ε > 0 est un paramètre à calibrer par l’utilisateur.

Un exemple en dimension 1

59
 2.0

1.5

Y 1.0

0.00 0.25 0.50 0.75 1.00

X

Un exemple en dimension 1

wx + b + ε
2.0

wx + b
1.5
Y

wx + b − ε
1.0

0.0 0.5 1.0

X

Influence de ε
— Il contrôle le niveau de tolérance que l’on se donne.

2.0

epsilon
1.5
0.5
Y

0.65

1.0

0.00 0.25 0.50 0.75 1.00

X

Alléger la contrainte...
— La contrainte nécessite souvent de prendre des grandes valeurs de ε...

60
 3.0

2.5

2.0

Y
1.5

1.0

0.00 0.25 0.50 0.75 1.00

X

— Clairement pas satisfaisant de prendre ε trop grand.

Idée
— Comme pour la SVM binaire, autoriser des observations à se situer en dehors de la marge !
— Comment ? En introduisant des slack variables !

SVR cas général

Le problème d’optimisation
On cherche (w, b, ξ, ξ ? ) qui minimise
n
1 X
kwk2 + C (ξi + ξi? )
2 i=1

 yi − hw, xi i − b ≤ ε + ξi , i = 1, . . . , n,
sous les contraintes hw, xi i + b − yi ≤ ε + ξi? , i = 1, . . . , n
ξi ≥ 0, ξi? ≥ 0, i = 1, . . . , n


Slack variables en régression

2.0
ξi

1.5
Y

1.0

ξi?

0.00 0.25 0.50 0.75 1.00

X

61
Rien ne change après...
— Les solutions s’obtiennent en résolvant le problème dual =⇒ αi , αi? .

— Les données (les X) sont généralement centrées-réduites pour éviter les problèmes d’échelle.

— Les observations vérifiant αi? − αi 6= 0 sont les vecteurs supports.

— L’hyperplan optimal se déduit des vecteurs supports :

n
X
w? = (αi? − αi )xi .
i=1

— L’astuce du noyau reste d’actualité pour prendre en compte de la non linéarité.

— Il faut sélectionner C, le noyau (et ses paramètres) ainsi que ε...

Le coin R
— Là aussi, pas grand chose ne change.
> svm(Y~.,data=df,kernel="linear",epsilon=0.5,cost=100)
##
## Call:
## svm(formula = Y ~ ., data = df, kernel = "linear", epsilon = 0.5,
## cost = 100)
##
##
## Parameters:
## SVM-Type: eps-regression
## SVM-Kernel: linear
## cost: 100
## gamma: 1
## epsilon: 0.5
##
##
## Number of Support Vectors: 11

Conclusion
— Algorithme machine learning pouvant être utilisé en régression et en classification supervisée.

— Méthode linéaire mais prise en compte possible de la non linéarité grâce à l’astuce du noyau.

— Calibration difficile : beaucoup de paramètres

1. paramètre de cout C
2. noyau
3. paramètres du noyau
4. seuil de tolérance ε pour la régression
— et souvent peu d’information a priori sur la valeur de ces paramètres...

5 Bibliographie
Références

Biblio2

[Aronszajn, 1950] Aronszajn, N. (1950). Theory of reproducing kernels. Transactions of the American Mathematical
Society, 68 :337–404.
[Bühlmann and van de Geer, 2011] Bühlmann, P. and van de Geer, S. (2011). Statistics for high-dimensional data.
Springer.
[Cornillon et al., 2019] Cornillon, P., Hengartner, N., Matzner-Løber, E., and Rouvière, L. (2019). Régression avec
R. EDP Sciences.

62
[Fromont, 2015] Fromont, M. (2015). Apprentissage statistique. Université Rennes 2, diapos de cours.
[Hastie et al., 2009] Hastie, T., Tibshirani, R., and Friedman, J. (2009). The Elements of Statistical Learning :
Data Mining, Inference, and Prediction. Springer, second edition.
[Hastie et al., 2015] Hastie, T., Tibshirani, R., and Wainwright, M. (2015). Statistical Learning with Sparsity : The
Lasso and Generalizations. CRC Press. [Link]
[Link].
[Karatzoglou et al., 2004] Karatzoglou, A., Smola, A., Hornik, K., and Zeileis, A. (2004). kernlab – an s4 package
for kernel methods in r. Journal of Statistical Software, 11(9).
[Tibshirani, 1996] Tibshirani, R. (1996). Regression shrinkage and selection via the lasso. Journal of the Royal
Statistical Society, Series B, 58 :267–288.
[Zou and Hastie, 2005] Zou, H. and Hastie, T. (2005). Regularization and variable selection via the elastic net.
Journal of the Royal Statistical Society, Series B, 67 :301–320.

63
Troisième partie
Algorithmes non linéaires
— Algorithmes linéaires :
f (x) = fβ (x) = β0 + β1 x1 + . . . + βd xd .
— Problème : tous les problèmes ne sont pas linéaires.
— Possible d’ajouter de la non linéarité dans les algorithmes linéaires : effets quadratiques, interaction...
— Difficile pour l’utilisateur de trouver quels effets ajouter ! Surtout lorsque d est grand.

Dans cette partie

Présentation de quelques algorithmes non linéaires :
— Méthodes par arbres.
— Réseaux de neurones.

1 Arbres
Présentation
— Les arbres sont des algorithmes de prédiction qui fonctionnent en régression et en discrimination.

— Il existe différentes variantes permettant de construire des prédicteurs par arbres.

— Nous nous focalisons dans cette partie sur la méthode CART [Breiman et al., 1984] qui est la plus utilisée.

1.1 Arbres binaires

Notations
— On cherche à expliquer une variable Y par d variables explicatives X1 , . . . , Xd .

— Y peut admettre un nombre quelconque de modalités et les variables X1 , . . . , Xd peuvent être qualitatives et/ou
quantitatives.

— Néanmoins, pour simplifier on se place dans un premier temps en discrimination binaire : Y admet 2 modalités
(-1 ou 1). On suppose de plus que l’on a simplement 2 variables explicatives quantitatives.

Représentation des données

— On dispose de n observations (x1 , y1 ), . . . , (xn , yn ) où xi ∈ R2 et yi ∈ {0, 1}.
1.00

0.75

Y
X2

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

Approche par arbres

Trouver une partition des observations qui sépare "au mieux" les points rouges des points bleus.

64
Arbres binaires
— La méthode CART propose de construire une partition basée sur des divisions successives parallèles aux axes.
— 2 exemples de partition :

Arbre 1 Arbre 2
1.00

0.75
X2

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00
X1

— A chaque étape, la méthode cherche une nouvelle division : une variable et un seuil de coupure.

1.00

0.75
X2

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

1.00

0.75
X2

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

65
 1.00

0.75

X2
0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

1.00

0.75
X2

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

1.00

0.75
X2

0.50

0.25

0.00
0.00 0.25 0.50 0.75 1.00
X1

Représentation de l’arbre
1.00

1
0.57
0.75 100%
yes X2 < 0.39 no

0 1
0.32 0.75
X2

0.50 42% 58%

X1 >= 0.58 X2 >= 0.79

1 0
0.25 0.61 0.17
21% 15%
X1 < 0.23 X1 >= 0.31

0.00
0 0 1 0 1 1
0.00 0.25 0.50 0.75 1.00 0.03 0.08 1.00 0.00 0.57 0.95
X1 21% 9% 12% 11% 5% 43%

66
Remarque
Visuel de droite plus pertinent :
— Plus d’information.
— Généralisation à plus de deux dimensions.

Vocabulaire
— Chaque coupure divise une partie de Rd en deux parties appelées nœuds.
— Le premier nœud, qui contient toutes les observations, est le nœud racine.
— Une coupure divise en noeud en deux nœuds fils :

Xj ≥ s Xj < s

N1 (j, s) N2 (j, s)

— Les nœuds qui ne sont pas découpés (en bas de l’arbre) sont les nœuds terminaux ou feuilles de l’arbre.

Arbre et algorithme de prévision

— L’arbre construit, les prévisions se déduisent à partir de moyennes faites dans les feuilles.
— On note N (x) la feuille de l’arbre qui contient x ∈ Rd , les prévisions s’obtiennent selon :
1. Régression =⇒ moyenne des yi de la feuille
1 X
mn (x) = yi
|N (x)|
i:xi ∈N (x)

2. Classification (classe) =⇒ vote à la majorité :

X
gn (x) = argmax 1yi =k
k
i:xi ∈N (x)

3. Classification (proba) =⇒ proportion d’obs. du groupe k :

1 X
Sk,n (x) = 1yi =k .
|N (x)|
i:xi ∈N (x)

Questions
1. Comment découper un nœud ?
=⇒ si on dispose d’un algorithme pour découper un nœud, il suffira de le répéter.
2. Comment choisir la profondeur de l’arbre ?
— Profondeur maximale ? (on découpe jusqu’à ne plus pouvoir) sur-ajustement ?
— Critère d’arrêt ?
— Élagage ? (on construit un arbre profond et on enlève des branches "inutiles"...).

1.2 Choix des coupures

— Une coupure = un couple (j, s) ∈ {1, . . . , d} × R.
— Idée : définir un critère mesure la performance d’une coupure et choisir celle qui optimise le critère.
— Coupure performante =⇒ les deux nœuds fils sont homogènes vis-à-vis de Y .

Fonction d’impureté
— Objectif : mesurer l’homogénéité d’un nœud.
— Intérêt : choisir la coupure qui maximise la pureté des nœuds fils.

67
Critère de découpe
— L’impureté I d’un nœud doit être :
1. faible lorsque un nœud est homogène : les valeurs de Y dans le nœud sont proches.
2. élevée lorsque un nœud est hétérogène : les valeurs de Y dans le nœud sont dispersées.

L’idée
Une fois I définie, on choisira le couple (j, s) qui maximise le gain d’impureté :

∆(j, s) = p(N )I(N ) − (p(N1 (j, s))I(N1 (j, s)) + p(N2 (j, s))I(N2 (j, s)))

où p(N ) représente la proportion d’observations dans le nœud N .

1.2.1 Cas de la régression

— Une mesure naturelle de l’impureté d’un nœud N en régression est la variance du nœud :
1 X
I(N ) = (yi − ȳN )2 ,
|N |
i:xi ∈N

où ȳN désigne la moyenne des Yi dans N .

12345 12345
6(2.92) 6(2.92)

234 12 45
1 5(4)
6(2.19) 3(0.67) 6(0.67)

=⇒ coupure de droite plus performante.

Exemple

Coupure 1 Coupure 2

2
Y

0.00 0.25 0.50 0.75 1.000.00 0.25 0.50 0.75 1.00

X

I(N ) I(N1 ) I(N2 ) ∆

Gauche 1.05 0.01 0.94 0.34
Droite 1.05 0.04 0.01 1.02

Pour aller plus vite

68
1.2.2 Cas de la classification supervisée
— Les Yi , i = 1, . . . , n sont à valeurs dans {1, . . . , K}.

— On cherche une fonction I telle que I(N ) soit

— petite si un label majoritaire se distingue clairement dans N ;

— grande sinon.

Impureté
L’impureté d’un nœud N en classification se mesure selon
K
X
I(N ) = f (pj (N ))
j=1

où
— pj (N ) représente la proportion d’observations de la classe j dans le nœud N .
— f est une fonction (concave) [0, 1] → R+ telle que f (0) = f (1) = 0.

Exemples de fonctions f
— Si N est pur, on veut I(N ) = 0 =⇒ c’est pourquoi f doit vérifier f (0) = f (1) = 0.

— Les 2 mesures d’impureté les plus classiques sont :

1. Gini : f (p) = p(1 − p) ;
2. Information : f (p) = −p log(p).

Cas binaire
Dans ce cas on a
1. I(N ) = 2p(1 − p) pour Gini
2. I(N ) = −p log p − (1 − p) log(1 − p) pour Information
où p désigne la proportion de 1 (ou 0) dans N .

Impureté dans le cas binaire

0.6

Indice
valeurs

0.4
gini

information
0.2

0.0
0.00 0.25 0.50 0.75 1.00
p

69
Exemple 1

N•N N•N
•N•(0.5) •N•(0.5)

NNN
N•N(0.44) •N•(0.44) •••(0)
(0)

=⇒ coupure de droite plus performante.

Exemple 2

Coupure 1 Coupure 2
1.00

0.75
Y
X2

0.50 0

0.25

0.00
0.00 0.25 0.50 0.75 1.000.00 0.25 0.50 0.75 1.00
X1

I(N ) I(N1 ) I(N2 ) ∆

Gauche 0.50 0.34 0.35 0.16
Droite 0.50 0.43 0.50 0.02

1.3 Elagage
Pourquoi élaguer ?
— Les coupures permettent de séparer les données selon Y =⇒ plus on coupe mieux on ajuste !
— Risque de sur-ajustement si on coupe trop !
Erreur

Ajustement

Prévision

Sur−apprentissage

Complexité

Complexité d’un arbre

Représentée par son nombre de coupures ou sa profondeur.

70
Comment faire ?
— Tester tous les arbres ? =⇒ possible uniquement sur de petits échantillons !
— Critère d’arrêt : ne plus découper si une certaine condition est vérifiée.=⇒ possible mais... une coupure peut
ne pas être pertinente alors que des coupures plus basses le seront !

Élaguer
1. Considérer un arbre (trop) profond =⇒ qui sur-ajuste ;
2. Supprimer les branches peu utiles.

Élagage CART
— Tester tous les sous-arbres d’un arbre très profond se révèlent souvent trop couteux en temps de calcul.
— [Breiman et al., 1984] propose une stratégie d’élagage qui permet de se ramener à une suite d’arbres emboités

Tmax = T0 ⊃ T1 ⊃ . . . ⊃ TK .

de taille raisonnable (plus petite que n).

— Il est ensuite possible de choisir un arbre dans cette suite par des méthodes traditionnelles :
1. choix d’un risque ;
2. optimisation de ce risque (par validation croisée par exemple).
Pour aller plus vite

Construction de la suite de sous arbres

— Soit T un arbre à |T | nœuds terminaux N1 , . . . , N|T | .
— Soit R(N ) un risque (d’ajustement) dans le nœud N :
— Régression :
1 X
Rm (T ) = (yi − ȳNm )2
Nm
i:xi ∈Nm

— Classification :
1 X
Rm (T ) = 1yi 6=yNm
Nm
i:xi ∈Nm

Définition
Soit α ≥ 0, le critère coût/complexité est défini par :
|T |
X
Cα (T ) = Nm Rm (T ) + α|T |.
m=1

Idée
— Cα (T ) est un critère qui prend en compte l’adéquation d’un arbre et sa complexité.
— L’idée est de chercher un arbre Tα qui minimise Cα (T ) pour une valeur de α bien choisie.

Remarque
— α = 0 =⇒ Tα = T0 = Tmax .
— α = +∞ =⇒ Tα = T+∞ = Troot arbre sans coupure.

Question (a priori difficile)

Comment calculer Tα qui minimise Cα (T ) ?

71
Deux lemmes
Lemme 1
Si T1 et T2 sont deux sous-arbres de Tmax avec Rα (T1 ) = Rα (T2 ). Alors T1 ⊂ T2 ou T2 ⊂ T1
=⇒ garantit une unique solution de taille minimale.
Lemme 2
Si α > α0 alors Tα = Tα0 ou Tα ⊂ Tα0 .
=⇒ garantit une stabilité des solutions lorsque α parcourt R+ =⇒ elles vont être emboîtées les unes dans les autres.

Théorème [Breiman et al., 1984]

Il existe une suite finie α0 = 0 < α1 < · · · < αM avec M ≤ |Tmax | et une suite associée d’arbres emboîtés (Tαm )m

Tmax = Tα0 ⊃ Tα1 ⊃ · · · ⊃ TαM = Troot

telle que ∀α ∈ [αm , αm+1 [

Tm ∈ argmin Cα (T ).
T ⊆Tmax

Tα0 Tα1 ... TαM

[ [ [ [
α0 = 0 α1 α2 ... αM

Commentaires
— Nombre de minimiseurs de Cα (T ) est "petit".
— Ils s’obtiennent en élaguant : en supprimant des branches.

Exemple
— On visualise la suite de sous-arbres avec la fonction printcp ou dans l’objet rpart :
> library(rpart)
> [Link](123)
> arbre <- rpart(Y~.,data=[Link],cp=0.0001,minsplit=2)
> arbre$cptable
## CP nsplit rel error xerror xstd
## 1 0.353846154 0 1.00000000 1.0000000 0.09336996
## 2 0.230769231 1 0.64615385 0.7076923 0.08688336
## 3 0.138461538 2 0.41538462 0.5076923 0.07805324
## 4 0.061538462 4 0.13846154 0.2153846 0.05481185
## 5 0.015384615 5 0.07692308 0.1846154 0.05111769
## 6 0.007692308 6 0.06153846 0.2461538 0.05816388
## 7 0.000100000 14 0.00000000 0.2153846 0.05481185

Sorties printcp
— Suite de 7 arbres emboités.
— CP : complexity parameter, il mesure la complexité de l’arbre : CP &=⇒ complexité %.
— nsplit : nombre de coupures de l’arbre.
— [Link] : erreur (normalisée) calculée sur les données d’apprentissage =⇒ erreur d’ajustement.
— xerror : erreur (normalisée) calculée par validation croisée 10 blocs =⇒ erreur de prévision (voir diapos
suivantes).
— xstd : écart-type associé à l’erreur de validation croisée.

Visualisation
— On peut les visualiser en combinant prune (extraction) et [Link] (tracé) :

72
 > arbre1 <- prune(arbre,cp=0.01)
> arbre2 <- prune(arbre,cp=0.1)
> library([Link])
> [Link](arbre1);[Link](arbre2)

1 1
0.57
100%
0.57
100%
yes X2 < 0.39 no
yes X2 < 0.39 no
0 1
0.32 0.75
42% 58%
0 1
0.32 0.75
X1 >= 0.58 X2 >= 0.79
42% 58%
1 0
0.61 0.17 X1 >= 0.58 X2 >= 0.79
21% 15%
1
X1 < 0.23 X1 >= 0.18
0.61
0
0.08 21%
9% X1 < 0.23
X1 >= 0.016

0 0 1 0 1
0 0 1 1 0 1 1
0.03 0.00 1.00 1.00 0.00 1.00 0.95 0.03 0.08 1.00 0.17 0.95
21% 8% 1% 12% 13% 3% 43% 21% 9% 12% 15% 43%

Choix de l’arbre final

— Choisir un arbre dans la suite revient à choisir une valeur de α.
— Ce choix s’effectue généralement de façon classique :
1. Choix d’un risque.
2. Estimation du risque par ré-échantillonnage (CV par exemple) pour tous les αm .
3. Sélection du αm qui minimise le risque estimé.

Remarque
La fonction rpart effectue par défaut une validation croisée 10 blocs en prenant :
— le risque quadratique en régression.
— l’erreur de classification en classification.

Validation croisée rpart

√ √
1. Calculer β0 = 0, β1 = α1 α2 , ... βM −1 = αM −1 αM , βM = +∞.
2. Pour k = 1, . . . , K
(a) Construire l’arbre maximal sur l’ensemble des données privé du k e bloc, c’est-à-dire B −k = {(xi , yi ) :
i ∈ {1, . . . , n}\Bk }.
(b) Appliquer l’algorithme d’élagage à cet arbre maximal, puis extraire les arbres qui correspondent aux
valeurs βm , m = 0, . . . , M =⇒ Tβm (., B −k ).
(c) Calculer les valeurs prédites par chaque arbre sur le bloc k : Tβm (xi , B −k ), i ∈ Bk .
3. En déduire les erreurs pour chaque βm :
K X
b m) = 1
X
R(β `(yi , Tβm (xi , B −k )).
n
k=1 i∈Bk

Retourner : une valeur αm telle que R(β

b m ) est minimum.

— Les erreurs de validation croisée se trouvent dans la colonne xerror de l’élément cptable.
— On peut les visualiser avec plotcp :
> plotcp(arbre)

73
 size of tree

1 2 3 5 6 7 15

1.2
1.0
X−val Relative Error

0.8
0.6
0.4
0.2
0.0
Inf 0.29 0.18 0.092 0.031 0.011 0.00088

cp

— Il reste à choisir l’arbre qui minimise l’erreur de prévision :

> cp_opt <- as_tibble(arbre$cptable) %>% arrange(xerror) %>%
+ slice(1) %>% select(CP) %>% [Link]()
> cp_opt
## [1] 0.01538462

— et à le visualiser :
> arbre_final <- prune(arbre,cp=cp_opt)
> [Link](arbre_final)

1
0.57
100%
yes X2 < 0.39 no

0 1
0.32 0.75
42% 58%
X1 >= 0.58 X2 >= 0.79
1 0
0.61 0.17
21% 15%
X1 < 0.23 X1 >= 0.18

0 0 1 0 1 1
0.03 0.08 1.00 0.00 1.00 0.95
21% 9% 12% 13% 3% 43%

— 2 variables explicatives =⇒ on peut visualiser l’arbre final

— en coloriant le carré [0, 1]2 en fonction des valeurs prédites.

1.00

0.75

Y
X2

0.50 0

0.25

0.00

0.00 0.25 0.50 0.75 1.00

X1

Prévision
— Nouvel individu :

74
 > xnew <- tibble(X1=0.4,X2=0.5)

— Prévision de la classe :
> predict(arbre_final,newdata=xnew,type="class")
## 1
## 1
## Levels: 0 1

— Prévision des probabilités :

> predict(arbre_final,newdata=xnew,type="prob")
## 0 1
## 1 0.046875 0.953125

1.4 Importance des variables

— La visualisation de l’arbre peut donner une idée sur l’importance des variables dans l’algorithme.
— Pas suffisant ! Il se peut en effet que des variables possèdent une grande importance sans pour autant apparaitre
explicitement dans l’arbre !
— Difficile de quantifier l’importance juste en regardant l’arbre !
— Il se peut en effet que des variables possèdent une grande importance sans pour autant apparaitre en
haut de l’arbre !

Mesure d’importance d’un arbre

Basée sur le gain d’impureté des nœuds internes.

Xj1 , ∆j1

— Nœuds internes =⇒ Nt , t = 1, . . . , J − 1 ;
— Variables de coupure =⇒ Xjt ; Xj2 , ∆j2
— Gain d’impureté =⇒ i2jt .

Xj3 , ∆j3

Mesure d’impureté de la variable `

|T |−1
X
I` (T ) = ∆t 1jt =` .
t=1

Exemple
0
0.35
100%
yes age < 51 no

0 1
0.22 0.56
63% 37%
age < 31 famhist = Absent
0 0 1
0.31 0.40 0.70
39% 18% 19%
typea < 69 tobacco < 7.6 ldl < 5
1
0.54
8%
adiposity >= 28
0
0.35
4%
tobacco < 4.2

0 0 1 0 1 0 1 1 1
0.07 0.27 0.83 0.28 0.71 0.10 0.60 0.74 0.82
23% 37% 3% 13% 5% 2% 2% 4% 11%

— Visualisation des importance à l’aide de vip :

75
 > library(vip)
> vip(arbre)

age

tobacco

adiposity

typea

ldl

sbp

famhist

obesity

alcohol

0 10 20 30
Importance

Bilan
1. Avantages :
— Méthode « simple » relativement facile à mettre en œuvre.
— Fonctionne en régression et en classification.
— Résultats interprétables (à condition que l’arbre ne soit pas trop profond).
2. Inconvénients :
— Performances prédictives limitées.
— méthode connue pour être instable, sensible à de légères perturbations de l’échantillon. =⇒ Cet inconvé-
nient sera un avantage pour des agrégations bootstrap =⇒ forêts aléatoires.

2 Réseaux de neurones
2.1 Introduction
Bibliographie
— Wikistat : Neural natworks and introduction to deep learning
— Eric Rakotomalala : Deep learning : Tensorflow et Keras sous R
— Rstudio : R interface to Keras

Historique
— Modélisation du neurone formel [McCulloch and Pitts, 1943].
— Concept mis en réseau avec une couche d’entrée et une sortie [Rosenblatt, 1958].
— Origine du perceptron
— Approche connexioniste (atteint ses limites technologiques et théoriques au début des années 70)
— Relance de l’approche connexioniste au début des années 80 avec l’essor technologique et quelques avancées
théoriques
— Estimation du gradient par rétro-propagation de l’erreur [Rumelhart et al., 1986].
— Développement considérable (au début des années 90)
— Remis en veilleuse au milieu des années 90 au profit d’autres algorithmes d’apprentissage : boosting, support
vector machine...
— Regain d’intérêt dans les années 2010, énorme battage médiatique sous l’appellation d’apprentissage pro-
fond/deep learning.
— Résultats spectaculaires obtenus par ces réseaux en reconnaissance d’images, traitement du langage naturel...

76
Différentes architectures
Il existe différents types de réseaux neuronaux :
— perceptron multicouches : les plus anciens et les plus simples ;

— réseaux de convolution : particulièrement efficaces pour le traitement d’images ;

— réseaux récurrents : adaptés à des données séquentielles (données textuelles, séries temporelles).

Dans cette partie

nous nous intéresserons uniquement au perceptron multicouches.

Neurone : vision biologique

Définition : neurone biologique

Un neurone biologique est une cellule qui se caractérise par
— des synapses : les points de connexion avec les autres neurones ;
— dentrites : entrées du neurones ;
— les axones ou sorties du neurone vers d’autres neurones ;
— le noyau qui active les sorties.

Définition : neurone formel

Un neurone formel est un modèle qui se caractérise par
— des entrées x1 , . . . , xp ;
— des poids w0 , w1 , . . . , wp ;
— une fonction d’activation σ : R → R ;
— une sortie :
ŷ = σ(w0 + w1 x1 + . . . + xp xp ).

2.2 Le perceptron simple

— Le problème : expliquer une sortie y ∈ R par des entrées x = (x1 , . . . , xp ).

Définition
Le perceptron simple est une fonction f des entrées x
— pondérées par un vecteur w = (w1 , . . . , wp ),
— complétées par un neurone de biais w0 ,
— et une fonction d’activation σ : R → R
ŷ = f (x) = σ(w0 + w1 x1 + . . . + xp xp ).

77
Fonction d’activation
Plusieurs fonctions d’activation peuvent être utilisées :
— Identité : σ(x) = x ;
— sigmoïde ou logistique : σ(x) = 1/(1 + exp(−x)) ;
— seuil : σ(x) = 1x≥0 ;
— ReLU (Rectified Linear Unit) : σ(x) = max(x, 0) ;
p
— Radiale : σ(x) = 1/2π exp(−x2 /2).
Remarque
Les poids wj sont estimés à partir des données (voir plus loin).

Représentation graphique

Entrées Poids Neurone Act. Sortie

w0

x1 w1

x2 w2

x3 w3 Σ σ y

..
.
..
.
xp wp

Le coin R
— Plusieurs packages R permettent d’ajuster des réseaux de neurones : nnet, deepnet. . .
— Nous présentons ici le package keras, initialement programmé en Python et qui a été "traduit" récemment en
R.
> library(keras)
> install_keras()

Exemple
— On veut expliquer une variable Y binaire par 4 variables d’entrées X1 , . . . , X4 .
— On dispose d’un échantillon d’apprentissage de taille 300 :
> head(dapp)
## X1 X2 X3 X4 Y
## 1 0.5855288 -1.4203239 1.67751179 -0.1746226 1
## 2 0.7094660 -2.4669386 0.07947405 -0.6706167 1
## 3 -0.1093033 0.4847158 -0.85642750 0.5074258 0
## 4 -0.4534972 -0.9379723 -0.77877729 1.2474343 0
## 5 0.6058875 3.3307333 -0.38093608 -1.2482755 1
## 6 -1.8179560 -0.1629455 -1.89735834 -1.9347187 1

Définition du modèle
— Elle s’effectue à l’aide des fonctions keras_model_sequential et layer_dense.
> model <- keras_model_sequential()
> model %>% layer_dense(units=1,input_shape=c(4),
+ activation="sigmoid")

— units : nombre de neurones souhaités ;

— activation : choix de la fonction d’activation.

78
Summary
— Un summary du modèle permet de visualiser le nombre de paramètres à estimer.

> summary(model)
## ___________________________________________________________________________
## Layer (type) Output Shape Param #
## ===========================================================================
## dense (Dense) (None, 1) 5
## ===========================================================================
## Total params: 5
## Trainable params: 5
## Non-trainable params: 0
## ___________________________________________________________________________

Estimation des paramètres

— On indique dans la fonction compile la fonction de perte pour l’estimation des paramètres du modèle et le
critère de performance
> model %>% compile(
+ loss="binary_crossentropy",
+ optimizer="adam",
+ metrics="accuracy"
+ )

Estimation
— On utilise la fonction fit pour entrainer le modèle
> Xtrain <- [Link](dapp[,1:4])
> Ytrain <- dapp$Y
> model %>% fit(x=Xtrain,y=Ytrain,epochs=300,batch_size=5)

— Et on obtient les poids avec get_weights :

> W <- get_weights(model)
> W
## [[1]]
## [,1]
## [1,] 0.250867128
## [2,] 0.092339918
## [3,] -0.162947521
## [4,] -0.005261241
##
## [[2]]
## [1] 0.1739036

Visualisation du réseau

Entrées Poids Neurone Act. Sortie

0.174

x1 0.251

x2 0.092

x3 −0.163 Σ σ y

x4 −0.005

Estimation
1
Pb(Y = 1|X = x) =
1 + exp(−(0.174 + 0.251x1 + . . . − 0.005x4 ))

79
Prévision
— On calcule la prévision de la probabilité P (Y = 1|X = x) pour le premier individu de l’échantillon test :
> w <- W[[1]]
> w0 <- W[[2]]
> Xtest <- [Link](dtest[,1:4])
> sc1 <- w0+sum(w*Xtest[1,])
> 1/(1+exp(-sc1))
## [1] 0.6209704

— que l’on retrouve avec predict_proba :

> prev <- model %>% predict_proba(Xtest)
> prev[1]
## [1] 0.6209704

2.3 Perceptron multicouches

Constat
— Règle de classification : le preceptron simple affecte un individu dans le groupe 1 si

P(Y = 1|X = x) ≥ 0.5 ⇐⇒ w0 + w1 x1 + . . . + wp xp ≥ 0.

— Il s’agit donc d’une règle linéaire.

=⇒ peu efficace pour représenter des phénomènes "complexes".

Idée
Conserver cette structure de réseau en considérant plusieurs couches de plusieurs neurones.

Perceptron simple

Entrées C. Sortie Sortie

x1

x2

x3
σ y
..
.
..
.
xp

Une couche cachée

Deux couches cachées

Commentaires
— Les neurones de la première couche (cachée) calculent des combinaisons linéaires des entrées.
— Ces combinaisons linéaires sont ensuite activées par une fonction d’activation, produisant une sortie par
neurone.
— Chaque neurone de la deuxième couche (cachée) est une combinaison linaire des sorties de la couche précé-
dente...
— activées par une fonction d’activation, produisant une sortie par neurone...

80
Entrées C. cachée C. Sortie Sortie

x1
σ1,1
x2

x3 σ1,2
σ y
.. ..
. .
..
.
σ1,k1
xp

Entrées C. cachée 1 C. cachée 2 C. Sortie Sortie

x1
σ1,1 σ2,1
x2

x3 σ1,2 σ2,2
σ y
.. .. ..
. . .
..
.
σ1,k1 σ2,k2
xp

Remarque
Le nombre de neurones dans la couche finale est définie par la dimension de la sortie y :
— Régression ou classification binaire =⇒ 1 neurone.
— Classification multiclasse (K) =⇒ K (ou K − 1) neurones.

Le coin R
— L’ajout de couches cachées dans keras est relativement simple.
— Il suffit de définir ces couches au moment de la spécification du modèle.
— Par exemple, pour deux couches cachées avec 10 et 5 neurones, on utilisera :
> model <- keras_model_sequential()
> model %>% layer_dense(units=10,input_shape=c(4),activation="sigmoid") %>%
+ layer_dense(units=5,activation="sigmoid") %>%
+ layer_dense(units=1,activation="sigmoid")

> summary(model)
## ___________________________________________________________________________
## Layer (type) Output Shape Param #
## ===========================================================================
## dense_1 (Dense) (None, 10) 50
## ___________________________________________________________________________
## dense_2 (Dense) (None, 5) 55
## ___________________________________________________________________________
## dense_3 (Dense) (None, 1) 6
## ===========================================================================
## Total params: 111
## Trainable params: 111
## Non-trainable params: 0
## ___________________________________________________________________________

81
2.4 Estimation
— L’utilisateur doit choisir le nombre de couches, le nombre de neurones par couche, les fonctions d’activation
de chaque neurone.
— Une fois ces paramètres choisis, il faut calculer (estimer) tous les vecteurs de poids dans tous les neurones.

L’approche
— On désigne par θ l’ensemble des paramètres à estimer =⇒ f (x, θ) la règle associée au réseau.
— Minimisation de risque empirique : minimiser
n
1X
Rn (θ) = `(yi , f (xi , θ))
n i=1

où ` est une fonction de perte (classique).

Fonctions de perte
— Erreur quadratique (régression) :
`(y, f (x)) = (y − f (x))2 .
— Cross-entropy ou log-vraisemblance négative (classification binaire 0/1) :

`(y, p(x)) = −(y log(p(x)) + (1 − y) log(1 − p(x)))

où p(x) = P(Y = 1|X = x).

— Cross-entropy ou log-vraisemblance négative (classification multi-classes) :
K
X
`(y, p(x)) = − 1y=k log(pk (x))
k=1

où pk (x) = P(Y = k|X = x).

Descente de gradient
— La solution s’obtient à l’aide de methodes de type descente de gradient :

θnew = θold − ε∇θ Rn (θold ).

— Le réseau étant structuré en couches, la mise à jour des paramètres n’est pas directe.

Algorithme de rétropropagation (voir ici)

1. Etape forward : calculer tous les poids associés à θold et stocker toutes les valeurs intermédiaires.
2. Etape backward :
(a) Calculer le gradient dans la couche de sortie.
(b) En déduire les gradients des couches cachées.

Batch et epoch
— L’algorithme de rétropropagation n’est généralement pas appliqué sur l’ensemble des données, mais sur des
sous-ensemble de cardinaux m appelés batch.
— Cette approche est classique sur les gros volumes de données et permet de prendre en compte des données
séquentielles.
— Pour prendre en compte toutes les données sur une étape de la descente de gradient, on va donc appliquer
n/m fois l’algorithme de rétropropagation.
— Une itération sur l’ensemble des données est appelée epoch.

82
Algorithme de rétropropagation stochastique
Algorithme
Entrées : ε (learning rate), m (taille des batchs), nb (nombre d’epochs).
1. Pour ` = 1 à nb
2. Partitionner aléatoire les données en n/m batch de taille m =⇒ B1 , . . . , Bn/m .
(a) Pour j = 1 à n/m
i. Calculer les gradients sur le batch j avec l’algorithme de rétropropagation : ∇θ .
ii. Mettre à jour les paramètres
θnew = θold − ε∇θold .
Sorties : θnew et f (x, θnew ).

Choix des paramètres

— ε (pas de la descente de gradient), généralement petit. Existence de versions améliorées de l’algorithme
précédent moins sensible à ce paramètre (RMSProp, Adam...).
— m (taille des batch) : généralement petit (pas trop en fonction du temps de calcul). L’utilisateur peut (doit)
faire plusieurs essais.
— nb (nombre d’epoch), proche du nombre d’itérations en boosting =⇒ risque de surapprentissage si trop grand.

En pratique
Il est courant de visualiser l’évolution de la fonction de perte et/ou d’un critère de performance en fonction du
nombre d’epoch.

Un exemple
— On considère un réseau à 2 couches cachées comportant 50 nœuds (2851 paramètres).
> model1 <- keras_model_sequential()
> model1 %>% layer_dense(units=50,input_shape=c(4),
+ activation="sigmoid") %>%
+ layer_dense(units = 50,activation = "sigmoid") %>%
+ layer_dense(units = 1,activation = "sigmoid")

— On utilise
— crossentropy comme perte.
— Adam comme algorithme d’optimisation.
— accuracy (taux de bien classés) comme mesure de performance.
> model1 %>% compile(
+ loss="binary_crossentropy",
+ optimizer="adam",
+ metrics="accuracy"
+ )

— On estime les paramètres avec m = 5 et nb = 1000 et utilise 20% des données dans l’échantillon de validation.
> history <- model1 %>% fit(
+ x=Xtrain,
+ y=Ytrain,
+ epochs=1000,
+ batch_size=5,
+ validation_split=0.2
+ )

Erreur et perte
> plot(history)

— On compare ce nouveau réseau avec le perceptron simple construit précédemment.

83
 ●
●●
●
● ●●●● ●●●
●● ● ●
●● ●
●●
0.75 ●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●●
●
●
●
●●
●
●
●
●●
● ●
●
● ●●● ● ●● ●
●
●
●
●●
●
●
●
●●●
●●●
●●
●●
●
●
●
●●
●
● ●
●
●●
●●●
●
●●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
● ●
●● ●●●
●●●
●
●●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●
●● ●
●
●●
●●●●
●
●
●
●
●
●●
●
●
● ●
●●
●
●
●●●
●
●●
●
●
●
● ●●
●
●●
●●
●
●●
● ●●
●●
●●●
●
●
●
●
● ●●
●
●●
● ●
●
●
●
● ●
●
●●
●
●
●
●
● ●●
●
●
●
●
● ●
●
●●
●
●
●
●
●
●
●
●
●
●
●●●●
●●
●
● ●
●
●
●● ●
●
●
●●
●
●
●● ●
●●●
●
●
●●
●
● ●●●
0.50 ●
●●
●
●
●
●
●
●●

loss
●
●
●●
●
●
● ●●●
●●●
●
●●
●
●
● ●●
●
●
●●
●
●
●
●
●● ● ●
●●
●●
●
●● ●
●
●
●
●
●
● ●
●●
● ●
●
●
●
●
●
● ●●●
●
●
●
●
●●
● ●●
●
●
●
●
●●
● ●
●
●
●●
●
●
●
●
●●
● ●●
●●
●
●
●
● ●
●●
●●
●
●
●● ●●
●
●●
●
●
●
● ●●
●
●●
●
●
●
●
●● ●●
●
●
●
●●
●
●
●
● ●●
●
●●
●
●
●
● ●●
●●
●
● ●
●●
●
●
●
●
●●●●●
●
●
●
●
● ●
●●
●
●●
0.25 ●
●
●
● ●
●●
●
●
●●● ●
●
●
●
●
●●●
●
●
●
●●
●
●●
●
●
●
●
●
●●
●
● ● ●
●
●
●●
●
●●
●●
●●
●
●
●
●
●●
●
●
●●
● ●
●
●
●
●
●●
●
●
●●●
●
●● ●
●
●
●
●
●
●●
●●
●
●
●
●
●
●
●
●
●●
●
● ● ●
● ● ●
● ●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
● ●● ●
● ● ● ●● ● ● ●
●●●
●
●
●●●
●
●●●
● ●●●
●●
●
●
●
●●
●
●●
●●●
●
●
●
●●
●
●
●
●●
●●●
●●●
●● ●●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
● ●
●
● ● ● ● ● ● ●● ● ● ● ● ● ●●●
●
●●●●
●●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●●
●
●
●
●
●
●●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●●●
●
●●
●
●
●
●●
●
●
●
●●●●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
● ●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●●
●●
●●
●
●
●●
●
●●●●
●●● ●●●
●●● ●● ●● ●● ●● ● ●●●●
●●
●
●
●●
●
● ●
●
●
●
● ●
●
●●
●
●●
●
●●
●
●
●
●●
●
●●
●
●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●●
●
●●●●
●● ●● ● ●
● ●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●●
●
●●
●
●●
●●
●
●
●●
●
●●
●
●
●
●●
●●
● ●●
●
●●● ●● ●●●
●
● ● ●●
●
●●● ●
●
●
●
●●
●
●●
●
●●
●●
●
●
●●
●
●
●
●
●●
●●
●
●
●
●●
●
●●
●
●●
●
●●●
● ● ● ●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
● ●
●
●●
●●
●
● ●
●●●
●
●
●●
●●
●
●
●
●●
●
●
●●
●
●●
●
●●
●
●
●●
●●
●
●
●●
●●
●●
●●
●
●
●●
●
●
●●●
●
●●
●
●
●
●●
●
●●
●
●
●●
●●
●
●
●
●
●●
●
●
●●
●
●●
●
●●
●
●
●●
●
●
●
●●●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●●
●
● ● ●
●●
●●
●
●
●●
●
●●
●
●●
●
●●
●●
●
●●●
●
●
●
●● ●
●●
● ●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●●
●
●●
●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●●
●
●
●
●●
●
●
●●●
●
●●●
●
●
●●
● ●
● ●
●●
● ●
● ●
● ●
● ●
● ●
●
● ●
●
●
●
● ●● ●
data
0.00 ● ●●
●●
●
●
●●●
●
●●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●

● training
1.0 ● ●
●●●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●
●● ●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
● ● validation
●
●●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●●●
●●
●
●●
●
●●●
●
●
●●●
●●●
●
●●
●
●
●
●● ●●●● ●
●
●●
●
●
● ●
●
●
●●●
●
●
●●
●
●●
●●
●
●●
●
●
●●
●
●●
●
●●
●
●●
●
●
●●
●
●
●●
●
●
●●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●●
●
●
●●
●
●
●●
●
●
●●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●●
●
●
●
●●
●
●
●
●●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●●
●
●
●
●
●
●●
●●
●
●●
●
●●
●
●
●●
●
●●
●
●
●
●●
●
●
● ●●
●
●●●●
●
●●
●●
●
●●●
●
●
●●
●●
●●●●
●●
●
●●
●●
●
●●
●
●●
●
●
●●
●●●
●
●●
●●
●●
●●
●
●●● ● ●
● ● ● ●● ●● ● ●
●
●
●
●●
●
●
●
●●
●
●●
●
●●
●
●●
●●
●
●●
● ● ●
●
●
●
●●
●
●
●
●●●
●
●●●
●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●
●●
●
●●
●●
●
●●
●
●
●●
●
●
●
●●●
●
●●
●
●●
●
●●
●
●●●
●
●●
●● ●
●●●
●●
●●
●●
●
●
●●●
●●●
●
●●
●
●●
●
●
●
●●
●
●●●
● ●●
●
●●
●●
●●
●
●
●●
● ●●●
●
●●
●●
●
●
●
●●● ●
●●
●
●
●● ●●●●
●●●●
●
●
●●
●
●●
0.8 ●
●
●
●
●
●
●
●●
●
●
●
●
●
●●
●
●
●
●
●
●
●
●●
●
●●
●
●
●
●
●
●●
●
●
●
●●
●
●●
●
●
●
●●● ●
●
●●●●
● ●
●
●●● ●●●
●
●●●
●
●●
●
●
●●
●
●
●●●●●
●●
● ●●
●
●●●●●
●
●
●
●
●
●
●
●
●
●●● ●● ●
● ●
● ●
●
●
●
●
●
●
● ●●
●●● ●
●●
● ●●● ●
●
● ●●●
●
●●
●●
● ●●
●
●●
acc

●
●● ●
●●
●
●
●●
●
●●● ● ●
●
● ●
●● ●
●
0.6 ●
●
●
●
●
● ● ●
●
●
●
● ●
●●
●
●●
●
●
● ● ● ●
●●
●● ●●
●
●
●●
●
●●
●
●
●●
●●
●
●●
●
● ●●●
●
●●
●
●
●
● ●●
●
●
●●
●
●●
●
●
●
●
●
●●
●
●
●
● ●
●
●●
●●
●●
●
●
●●
●
●
●●● ●●
●●
●
●
●● ●
●
●● ●
●●● ● ●
● ●● ●
● ●●
● ●
● ●●
● ●
● ●●

● ●
●●
●
●
●
●
●
●
●
●●
●
● ●●
●
●
0.4 ●
●
●●
●
●
●●
●●
●●
●● ●
●
●
●
●
●
●
●
●●●
●●
●●
●
●
●
●
●
● ●
●●
●●
●●●●
●●
●●
●●
●●
●●
●●
●●
●●●
●
●●
●
●
●
●●
●●●
●
●●●●

0 200 400 600 800 1000

epoch

> Xtest <- [Link](dtest[,1:4])

> Ytest <- dtest$Y
> model %>% evaluate(Xtest,Ytest)
## $loss
## [1] 0.7259337
##
## $acc
## [1] 0.39
> model1 %>% evaluate(Xtest,Ytest)
## $loss
## [1] 0.3290039
##
## $acc
## [1] 0.935

Nombre de couches et de neurones

— A choisir par l’utilisateur.
— Il est généralement mieux d’en avoir trop que pas assez =⇒ plus "facile" de capter des non linéarités complexes
avec beaucoup de couches et de neurones.
— On fait généralement plusieurs essais que l’on compare (avec caret par exemple).
— Voir par exemple l’appli suivante :
[Link]

2.5 Choix des paramètres et surapprentissage

Surapprentissage
— Plusieurs paramètres peuvent causer du surapprentissage, notamment les nombres de couches cachées, de
neurones et d’epoch.

Plusieurs solutions
1. Régularisation de type ridge/lasso :
n
1X
Rn (θ) = `(yi , f (xi , θ)) + λΩ(θ).
n i=1

=⇒ ajouter kernel_regularizer = regularizer_l2(l = 0.001)) dans la fonction layer_dense par exemple.

2. Early stopping : on stoppe l’algorithme lorsque l’ajout d’epoch n’améliore pas suffisamment un critère donné.
3. Dropout : suppression (aléatoire) de certains neurones dans les couches =⇒ souvent la solution privilégiée.

84
C1 C2 C3 S C1 C2 C3 S

Dropout
— A chaque étape de la phase d’entrainement, on supprime un nombre de neurones (selon une Bernoulli de
paramètre p).

Le coin R
— Il suffit d’ajouter layer_dropout après les couches cachées.
> model3 <- keras_model_sequential()
> model3 %>% layer_dense(units=50,input_shape=c(4),activation="sigmoid") %>%
+ layer_dropout(0.5) %>%
+ layer_dense(units = 50,activation = "sigmoid") %>%
+ layer_dropout(0.5) %>%
+ layer_dense(units = 1,activation = "sigmoid")

Sélection avec caret

— On peut sélectionner la plupart des paramètres avec caret.
— On propose par exemple, pour un réseau avec une couche cachée, de choisir
1. le nombre de neurones dans la couche cachée parmi 10, 50, 100
2. la fonction d’activation : sigmoïde ou relu.
— On définit d’abord les paramètres du modèle
> library(caret)
> dapp1 <- dapp
> dapp1$Y <- [Link](dapp1$Y)
> param_grid <- [Link](size=c(10,50,100),
+ lambda=0,batch_size=5,lr=0.001,
+ rho=0.9,decay=0,
+ activation=c("relu","sigmoid"))

— On calcule ensuite les taux de bien classés par validation croisée 5 blocs pour chaque combinaison de para-
mètres.
> caret_mlp <- train(Y~.,data=dapp1,method="mlpKerasDecay",
+ tuneGrid=param_grid,epoch=500,verbose=0,
+ trControl=trainControl(method="cv",number=5))

> caret_mlp
## Multilayer Perceptron Network with Weight Decay
## 300 samples
## 4 predictor
## 2 classes: ’0’, ’1’
## No pre-processing
## Resampling: Cross-Validated (5 fold)
## Summary of sample sizes: 240, 240, 240, 240, 240
## Resampling results across tuning parameters:
## size activation Accuracy Kappa
## 10 relu 0.9200000 0.8394122
## 10 sigmoid 0.8966667 0.7913512
## 50 relu 0.9266667 0.8515286
## 50 sigmoid 0.9066667 0.8127427

85
## 100 relu 0.9366667 0.8722974
## 100 sigmoid 0.9300000 0.8595025
## Tuning parameter ’lambda’ was held constant at a value of 0
## Tuning parameter ’rho’ was held constant at a value of 0.9
## Tuning parameter ’decay’ was held constant at a value of 0
## Accuracy was used to select the optimal model using the largest value.
## The final values used for the model were size = 100, lambda =
## 0, batch_size = 5, lr = 0.001, rho = 0.9, decay = 0 and activation = relu.

Conclusion
— Avantages :

— Méthode connue pour être efficace pour (quasiment) tous les problèmes.

— Plus particulièrement sur des architectures particulières : images, données textuelles.

— Inconvénients :

— Gain plus discutable sur des problèmes standards.

— (Beaucoup) plus difficile à calibrer que les autres algorithmes ML.

— Niveau d’expertise important.

3 Bibliographie
Références

Biblio3

[Breiman et al., 1984] Breiman, L., Friedman, J., Olshen, R., and Stone, C. (1984). Classification and regression
trees. Wadsworth & Brooks.
[McCulloch and Pitts, 1943] McCulloch, W. and Pitts, W. (1943). A logical calculus of ideas immanent in nervous
activity. Bulletin of Mathematical Biophysics, 5 :115–133.
[Rosenblatt, 1958] Rosenblatt, F. (1958). The perceptron : a probabilistic model for information storage and
organization in the brain. Psychological Review, 65 :386–408.
[Rumelhart et al., 1986] Rumelhart, D. E., Hinton, G. E., and R. J. Williams, R. J. (1986). Learning representations
by back-propagating errors. Nature, pages 533–536.

86
Quatrième partie
Agrégation
— Idée : construire un grand nombre d’algorithmes "simples" et les agréger pour obtenir une seule prévision.
Par exemple

Dn,B TB (x, Dn,B )

..
. 1
PB
fn (x) = B k=1 Tk (x, Dn,k )
Dn,2 T2 (x, Dn,2 )

Dn,1 T1 (x, Dn,1 )

Questions
1. Comment choisir les échantillons Dn,b ?
2. Comment choisir les algorithmes ?
3. ...

1 Bagging et forêts aléatoires

Cadre
— Idem que précédemment, on cherche à expliquer une variable Y par d variables explicatives X1 , . . . , Xd .

— Pour simplifier on se place en régression : Y est à valeurs dans R mais tout ce qui va être fait s’étant
directement à la classification binaire ou multiclasses.

— Notations :

— (X, Y ) un couple aléatoire à valeurs dans Rd × R.

— Dn = (X1 , Y1 ), . . . , (Xn , Yn ) un n-échantillon i.i.d. de même loi que (X, Y ).

— Un algorithme de la forme :
B
1 X
fn (x) = Tb (x)
B
b=1

— Hypothèse : les T1 , . . . , Tb sont identiquement distribuées.

Propriété
1 − ρ(x)
E[fn (x)] = E[T1 (x)] et V[fn (x)] = ρ(x)V[T1 (x)] + V[T1 (x)]
B
où ρ(x) = corr(T1 (x), T2 (x)).

Conséquence
— Biais non modifié.
— Variance & si B % et ρ(x) &.

— Ajuster le même algorithme sur les mêmes données n’est d’aucun intérêt.
— Ajuster le même algorithme sur des sous-échantillons disjoints est d’un intérêt limité.
— Utiliser un grand nombre d’algorithmes différents est compliqué...

Idée
Ajuster le même algorithme sur des échantillons bootstraps.

87
 1 2 3 4 5 6 7 8 9 10

3 4 6 10 3 9 10 7 7 1 T1
2 8 6 2 10 10 2 9 5 6 T2
2 9 4 4 7 7 2 3 6 7 T3
6 1 3 3 9 3 8 10 10 1 T4
3 7 10 3 2 8 6 9 10 2 T5
.. ..
. .
7 10 3 4 9 10 10 8 6 1 TB

1.1 Bagging
— Le bagging désigne un ensemble de méthodes introduit par Léo Breiman [Breiman, 1996].

— Bagging : vient de la contraction de Bootstrap Aggregating.

— Idée : plutôt que de constuire un seul estimateur, en construire un grand nombre (sur des échantillons boots-
trap) et les agréger.

Idée : échantillons bootstrap

— Echantillon initial :
— Echantillons bootstrap : tirage de taille n avec remise
— A la fin, on agrège :
B
1 X
fn (x) = Tb (x)
B
b=1

Algorithme bagging
Entrées :
— B un entier positif ;
— T un algorithme de prévision.
Pour b entre 1 et B :
1. Faire un tirage aléatoire avec remise de taille n dans {1, . . . , n}. On note θb l’ensemble des indices sélectionnés
?
et Dn,b = {(xi , yi ), i ∈ θb } l’échantillon bootstrap associé.
?
2. Entraîner l’algorithme T sur Dn,b =⇒ T (., θb , Dn ).
B
Retourner : fn (x) = B1 b=1 T (x, θb , Dn ).
P

Un algorithme pas forcément aléatoire

— L’aléa bootstrap implique que l’algorithme "change" lorsqu’on l’exécute plusieurs fois mais...
B
1 X
lim T (x, θb , Dn ) = Eθ [T (x, θ, Dn )] = f¯n (x, Dn )
B→+∞ B
b=1

Conséquence
— L’algorithme se stabilise (converge) lorsque B %.
— Recommandation : choisir B le plus grand possible.

88
Choix de T

1 − ρ(x)
E[fn (x)] = E[T1 (x)] et V[fn (x)] = ρ(x)V[T1 (x)] + V[T1 (x)].
B
Conclusion
— Bagger ne modifie pas le biais.
— B grand =⇒ V[fn (x)] ≈ ρ(x)V[T1 (x)] =⇒ la variance diminue d’autant plus que la corrélation entre les
prédicteurs diminue.
— Il est donc nécessaire d’agréger des estimateurs sensibles à de légères perturbations de l’échantillon.
— Les arbres sont connus pour posséder de telles propriétés.

1.2 Forêts aléatoires

Rappels sur les arbres
1
0.57
100%
yes X2 < 0.39 no

0 1
0.32 0.75
42% 58%
X1 >= 0.58 X2 >= 0.79

1 0
0.61 0.17
21% 15%
X1 < 0.23 X1 >= 0.31

0 0 1 0 1 1
0.03 0.08 1.00 0.00 0.57 0.95
21% 9% 12% 11% 5% 43%

Complexité
Profondeur
— petite : biais %, variance &
— grande : biais &, variance % (sur-apprentissage).

Définition
— Comme son nom l’indique, une forêt aléatoire est définie à partir d’un ensemble d’arbres.

Définition
Soit Tk (x), k = 1, . . . , B des prédicteurs par arbre (Tk : Rd → R). Le prédicteur des forêts aléatoires est obtenu par
agrégation de cette collection d’arbres :
B
1 X
fn (x) = Tk (x).
B
k=1

Forêts aléatoires
— Forêts aléatoires = collection d’abres.

— Les forêts aléatoires les plus utilisées sont (de loin) celles proposées par Léo Breiman (au début des années
2000).

— Elles consistent à agréger des arbres construits sur des échantillons bootstrap.

— On pourra trouver de la doc à l’url

http: // www. stat. berkeley. edu/ ~breiman/ RandomForests/
et consulter la thèse de Robin Genuer [Genuer, 2010].

89
1.2.1 Algorithme
Coupures "aléatoires"

Xj ≥ s Xj < s

N1 (j, s) N2 (j, s)

Arbres pour forêt

— Breiman propose de sélectionner la "meilleure" variable dans un ensemble composé uniquement de mtry
variables choisies aléatoirement parmi les d variables initiales.
— Objectif : diminuer la corrélation entre les arbres que l’on agrège.

Algorithme forêts aléatoires

Entrées :
— B un entier positif ;
— mtry un entier entre 1 et d ;
— [Link] un entier plus petit que n.
Pour b entre 1 et B :
?
1. Faire un tirage aléatoire avec remise de taille n dans {1, . . . , n}. On note Ib l’ensemble des indices sélectionnés et Dn,b = {(xi , yi ), i ∈ Ib }
l’échantillon bootstrap associé.
?
2. Construire un arbre CART à partir de Dn,b en découpant chaque nœud de la façon suivante :
(a) Choisir mtry variables au hasard parmi les d variables explicatives ;
(b) Sélectionner la meilleure coupure Xj ≤ s en ne considérant que les mtry variables sélectionnées ;
(c) Ne pas découper un nœud s’il contient moins de [Link] observations.
3. On note T (., θb , Dn ) l’arbre obtenu.
1
PB
Retourner : fn (x) = B b=1 T (x, θb , Dn ).

Type de prévision
La prévision dépend de la nature de Y et de ce que l’on souhaite estimer
— Régression : T (x, θb , Dn ) ∈ R et
B
1 X
mn (x) = T (x, θb , Dn ).
B
b=1

— Classification (classe) : T (x, θb , Dn ) ∈ {1, . . . , K} et

B
X
gn (x) ∈ argmax 1T (x,θb ,Dn )=k , k = 1, . . . , K.
k∈{1,...,K} b=1

— Classification (proba) : Tk (x, θb , Dn ) ∈ [0, 1] et

B
1 X
Sn,k (x) = Tk (x, θb , Dn ), k = 1, . . . , K.
B
b=1

Le coin R
— Notamment 2 packages avec à peu près la même syntaxe.
— randomforest : le plus ancien et probablement encore le plus utilisé.
— ranger [Wright and Ziegler, 2017] : plus efficace au niveau temps de calcul (codé en C++).

90
 > library(ranger)
> [Link](12345)
> foret <- ranger(type~.,data=spam)
> foret
## ranger(type ~ ., data = spam)
## Type: Classification
## Number of trees: 500
## Sample size: 4601
## Number of independent variables: 57
## Mtry: 7
## Target node size: 1
## Variable importance mode: none
## Splitrule: gini
## OOB prediction error: 4.59 %

1.2.2 Choix des paramètres

— B réglé =⇒ le plus grand possible. En pratique on pourra s’assurer que le courbe d’erreur en fonction du
nombre d’arbres est stabilisée.
— Pour les autres paramètres on étudie à nouveau :
1 − ρ(x)
E[fn (x)] = E[T1 (x)] et V[fn (x)] = ρ(x)V[T1 (x)] + V[T1 (x)].
B
Conséquence
— Le biais n’étant pas amélioré par "l’agrégation bagging", il est recommandé d’agréger des estimateurs qui
possèdent un biais faible (contrairement au boosting).
— Arbres "profonds", peu d’observations dans les nœuds terminaux.
— Par défaut dans randomForest, [Link] = 5 en régression et 1 en classification.

Choix de mtry
— Il est en relation avec la corrélation entre les arbres ρ(x).
— Ce paramètre a une influence sur le compromis biais/variance de la forêt.
— mtry &
1. tendance à se rapprocher d’un choix "aléatoire" des variables de découpe des arbres =⇒ les arbres sont
de plus en plus différents =⇒ ρ(x) & =⇒ la variance de la forêt diminue.
2. mais... le biais des arbres %=⇒ le biais de la forêt %.
— Inversement lorsque mtry % (risque de sur-ajustement).
Conclusion
— Il est recommandé de comparer les performances de la forêt pour plusieurs valeurs de mtry.
√
— Par défaut mtry = d/3 en régression et d en classification.
— Visualisation d’erreur en fonction de [Link] et mtry

0.09
[Link]
1

0.08 5
Erreur

15
0.07
50

100
0.06
500

0.05
1 3 10 30
mtry

91
Commentaires
[Link] petit et mtry à calibrer.

En pratique
— On peut bien entendu calibrer ces paramètres avec les approches traditionnelles mais...
— les valeurs par défaut sont souvent performantes !
— On pourra quand même faire quelques essais, notamment pour mtry.

Un exemple avec tidymodels

1. Initialisation du workflow :
> tune_spec <- rand_forest(mtry = tune(),min_n= tune()) %>%
+ set_engine("ranger") %>%
+ set_mode("classification")
> rf_wf <- workflow() %>% add_model(tune_spec) %>% add_formula(type ~ .)

2. Ré-échantillonnage et grille de paramètres :

> blocs <- vfold_cv(spam, v = 10,repeats = 5)
> rf_grid <- [Link](mtry=c(seq(1,55,by=5),57),
+ min_n=c(1,5,15,50,100,500))

3. Calcul des erreurs :

> rf_res <- rf_wf %>% tune_grid(resamples = blocs,grid = rf_grid)

4. Visualisation des résultats (AUC et accuracy) :

> rf_res %>% show_best("roc_auc") %>% select(-8)
## # A tibble: 5 x 7
## mtry min_n .metric .estimator mean n std_err
## <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl>
## 1 4 1 roc_auc binary 0.988 50 0.000614
## 2 5 1 roc_auc binary 0.988 50 0.000623
## 3 6 1 roc_auc binary 0.988 50 0.000617
## 4 5 5 roc_auc binary 0.988 50 0.000621
## 5 7 1 roc_auc binary 0.988 50 0.000645

> rf_res %>% show_best("accuracy") %>% select(-8)

## # A tibble: 5 x 7
## mtry min_n .metric .estimator mean n std_err
## <dbl> <dbl> <chr> <chr> <dbl> <int> <dbl>
## 1 4 1 accuracy binary 0.954 50 0.00159
## 2 6 1 accuracy binary 0.954 50 0.00141
## 3 7 1 accuracy binary 0.954 50 0.00149
## 4 5 1 accuracy binary 0.954 50 0.00153
## 5 8 1 accuracy binary 0.953 50 0.00146

Remarque
On retrouve bien [Link] petit et mtry proche de la valeur par défaut (7).
5. Ajustement de l’algorithme final :
> foret_finale <- rf_wf %>%
+ finalize_workflow(list(mtry=7,min_n=1)) %>%
+ fit(data=spam)

1.2.3 Erreur OOB et importance des variables

— Comme pour tous les algorithmes de prévision on peut évaluer la performance des forêts aléatoires en estimant
un risque par ré-échantillonnage.
— Les tirages bootstraps permettent de définir une alternative, souvent moins couteuse en temps de calcul, au
ré-échantillonnage.
— Idée/astuce : utiliser les observations non sélectionnées dans les échantillons bootstraps pour estimer le risque.

92
 3 4 6 10 3 9 10 7 7 1 T1
2 8 6 2 10 10 2 9 5 6 T2
2 9 4 4 7 7 2 3 6 7 T3
6 1 3 3 9 3 8 10 10 1 T4
3 7 10 3 2 8 6 9 10 2 T5
7 10 3 4 9 10 10 8 6 1 T6

OOB illustration
— Les échantillons 2, 3 et 5 ne contiennent pas la première observation, donc
1
ŷ1 = (T2 (x1 ) + T3 (x1 ) + T5 (x1 )).
3
— On fait de même pour toutes les observations =⇒ ŷ2 , . . . , ŷn .
— On calcule l’erreur selon
n n
1X 1X
(ŷi − yi )2 ou 1ŷ 6=y .
n i=1 n i=1 i i

OOB définition
— Pour i = 1, . . . , n on note
OOB(i) = {b ≤ B : i ∈
/ Ib }
l’ensemble des tirages bootstrap qui ne contiennent pas i et
1 X
fn,OOB(i) (xi ) = T (xi , θb , Dn )
|OOB(i)|
b∈OOB(i)

la prévision de la forêt en ne considérant que les arbres pour lesquels i n’est pas dans le tirage bootstrap.
— L’erreur OOB s’obtient en confrontant ces prévisions au valeurs observées, par exemple
n n
1X 1X
(yi − fn,OOB(i) (xi ))2 ou 1f (x )6=y .
n i=1 n i=1 n,OOB(i) i i

=⇒ erreur renvoyée par défaut dans ranger et randomforest.

Importance des variables

Deux mesures sont généralement utilisées.
— Score d’impureté : simplement la moyenne des importances de Xj dans chaque arbre de la forêt :
B
1 X
Ijimp = Ij (Tb ),
B
b=1

voir chapitre sur les arbres pour la définition de Ij (Tb ).

— Importance par permutation : comparer les erreurs de chaque arbre sur l’échantillon
1. OOB de l’arbre ;
2. OOB en permutant les valeurs de la variables j.
=⇒ Idée : Si Xj est importante ces erreurs doivent êtres très différentes.

Importance par permutation

— On présente ce score en régression mais rien ne change pour la classification.
— On note
1 X
Err(OOBb ) = (yi − T (xi , θb , Dn ))2 ,
|OOBb |
i∈OOBb
avec
OOBb = {i ≤ n : i ∈
/ Ib }.
=⇒ Erreur de l’arbre b calculée sur les données OOB.

93
 — On recalcule cette erreur mais sur OOBb où on permute les valeurs de la j e colonne.
   
x11 ... x1j ... x1d x11 ... x3j ... x1d
x21 ... x2j ... x2d  x21 ... x5j ... x2d 
   
x51 ... x3j ... x3d 
 =⇒ x51 ... x1j ... x3d 

 
x41 ... x4j ... x4d  x41 ... x2j ... x4d 
x51 ... x5j ... x5d x51 ... x4j ... x5d

— On note x̃ji les individus de l’échantillon OOBb permuté et on calcule

1 X
Err(OOBjb ) = (yi − T (x̃ji , θb , Dn ))2 .
|OOBb |
i∈OOBb

Importance par permutation

B
1 X
Ijperm = (Err(OOBjb ) − Err(OOBb )).
B
b=1

Le coin R
— On peut calculer et visualiser facilement ces importances avec ranger :
> [Link](1234)
> [Link] <- ranger(type~.,data=spam,importance="impurity")
> [Link] <- ranger(type~.,data=spam,importance="permutation")
> vip([Link]);vip([Link])

charExclamation capitalLong

charDollar hp

remove charExclamation

free remove

capitalAve capitalAve

your capitalTotal

capitalLong charDollar

hp free

capitalTotal your

money george

0 50 100 150 200 250 0.00 0.01 0.02 0.03 0.04

Importance Importance

Conclusion
Beaucoup d’avantages
— Bonnes performances prédictives =⇒ souvent parmi les algorithmes de tête dans les compétitions
[Fernández-Delgado et al., 2014].
— Facile à calibrer.

Assez peu d’inconvénients

Coté boîte noire (mais guère plus que les autres méthodes...)

94
2 Boosting
— Le terme Boosting s’applique à des méthodes générales permettant de produire des décisions précises à partir
de règles faibles (weaklearner).

— Historiquement, le premier algorithme boosting est adaboost [Freund and Schapire, 1996].

— Beaucoup de travaux ont par la suite été développés pour comprendre et généraliser ces algorithmes (voir
[Hastie et al., 2009]) :
— modèle additif
— descente de gradient =⇒ gradient boosting machine, extreme gradient bossting (Xgboost).
— ...

— Dans cette partie =⇒ descente de gradient.

Retour aux sources...

— Machine learning =⇒ objectifs prédictifs =⇒ minimisation de risque.
— Risque d’une fonction de prévision f : Rd → R :

R(f ) = E[`(Y, f (X))].

— R(f ) inconnu =⇒ version empirique

n
1X
Rn (f ) = `(yi , f (xi )).
n i=1

Idée
Minimiser Rn (f ) sur une classe d’algorithmes F.

Choix de F
— Il est bien entendu crucial.
— F riche/complexité élevée =⇒ Rn (f ) & =⇒ f (xi ) ≈ yi , i = 1, . . . , n =⇒ sur-ajustement.
— et réciproquement pour des classes F simple/complexité faible.

Combinaisons d’arbres
— [Friedman, 2001, Friedman, 2002] propose de se restreindre à des combinaisons d’arbres :
( B )
X
F= λb T (x, θb ), λb ∈ R, θb ∈ Θ
b=1

où θb désigne les paramètres de l’arbre (impureté, profondeur)...

— Rappel : un arbre peut s’écrire

L
X
T (x, θb ) = γb` 1x∈Nb`
`=1

où Nb` désigne les feuilles et γb` les prévisions dans les feuilles.
— Les paramètres B, θb définissent la complexité de F.
— Il faudra les calibrer à un moment mais nous les considérons fixés pour l’instant.

Un premier problème
Chercher f ∈ F qui minimise Rn (f ).
— Résolution numérique trop difficile.
— Nécessité de trouver un algorithme qui approche la solution.

95
2.1 Algorithme de gradient boosting
Descentes de gradient
— Définissent des suites qui convergent vers des extrema locaux de fonctions Rp → R.
— Le risque Rn (f ) ne dépend que des valeurs de f aux points xi .
— En notant f = (f (x1 ), . . . , f (xn )) ∈ Rn , on a
n
e n (f ) = 1
X
Rn (f ) = R `(yi , f (xi ))
n i=1

e n : Rn → R.
avec R

Nouveau problème
Minimiser R
e n . =⇒ en gardant en tête que minimiser de Rn (f ) n’est pas équivalent à minimiser R
e n (f ).

— Descente de gradient =⇒ suite (fb )b de vecteurs de Rn qui convergent vers des extrema (locaux) de R
en.
— Suite récursive :
fb = fb−1 − ρb ∇R
fn (fb−1 ),

où ∇Rfn (fb−1 ) désigne le vecteur gradient de R fn évalué en fb−1 . =⇒ vecteur de Rn donc la ie coordonnée
vaut
∂R
fn (f ) ∂`(yi , f (xi ))
(fb−1 ) = (fb−1 (xi )).
∂f (xi ) ∂f (xi )

Exemple
Si `(y, f (x)) = 1/2(y − f (x)2 alors

∂`(yi , f (xi ))
− (fb−1 (xi )) = yi − fb−1 (xi ),
∂f (xi )

=⇒ résidu de fb−1 (xi ).

— Si tout se passe bien... la suite (fb )b doit converger vers un minimum de R
fn .

Deux problèmes
1. Cette suite définit des prévisions uniquement aux points xi =⇒ impossible de prédire en tout x.
2. Les éléments de la suite ne s’écrivent pas comme des combinaisons d’arbres.

Une solution
[Friedman, 2001] propose d’ajuster un arbre sur les valeurs du gradient à chaque étape de la descente.

Algorithme de gradient boosting

1
Pn
1. Initialisation : f0 (.) = argminc n i=1 `(yi , c)
2. Pour b = 1 à B :
(a) Calculer l’opposé du gradient − ∂f ∂(xi ) `(yi , f (xi )) et l’évaluer aux points fb−1 (xi ) :

∂
ui = − `(yi , f (xi )) , i = 1, . . . , n.
∂f (xi ) f (xi )=fb−1 (xi )

(b) Ajuster un arbre de régression à J feuilles sur (xi , ui ), . . . , (xn , un ).

(c) Calculer les valeurs prédites dans chaque feuille
n
X
γjb = argmin `(yi , fb−1 (xi ) + γ).
γ
i:xi ∈Njb

PJ
(d) Mise à jour : fb (x) = fb−1 (x) + j=1 γjb 1x∈Njb .
Retourner : l’algorithme fn (x) = fB (x).

96
Paramètres
Nous donnons les correspondances entre les paramètres et les options de la fonction gbm :
— ` la fonction de perte =⇒ distribution
— B nombre d’itérations =⇒ [Link]
— J le nombre de feuilles des arbres =⇒ [Link] (=J − 1)
— λ le paramètre de rétrécissement =⇒ shrinkage.

Stochastic gradient boosting

[Friedman, 2002] montre qu’ajuster les arbres sur des sous-échantillons (tirage sans remise) améliore souvent les
performances de l’algorithme. =⇒ [Link] : taille des sous-échantillons.

Exemple
— Données sinus

0
Y

−1

−4 0 4
X .
— On entraîne l’algorithme :
> [Link](1234)
> library(gbm)
> boost.5000 <- gbm(Y~.,data=data_sinus,
+ distribution="gaussian",shrinkage=0.1,[Link] = 5000)

— On visualise les prévisions en fonction du nombre d’itérations :

1 50

−1

boost
Y

300 5000
sinus
1

−1

−4 0 4 −4 0 4
X

2.2 Choix des paramètres

Fonction de perte
— Pas vraiment un paramètre...
— Elle doit

97
 1. mesurer un coût (comme d’habitude). =⇒ elle caractérise la fonction de prévision à estimer =⇒ fn est
en effet un estimateur de
f ? ∈ argmin E[`(Y, f (X))].
f :Rd →R

2. être convexe et dérivable par rapport à son second argument (spécificité gradient).

L2 -boosting en régression
— Correspond à la perte quadratique
1
`(y, f (x)) = (y − f (x))2 .
2
— fonction de prévision optimale : f ? (x) = E[Y |X = x].

Remarque
— Avec cette perte, les ui sont donnés par

∂`(yi , f (xi ))
ui = − (fb−1 (xi )) = yi − fb−1 (xi ),
∂f (xi )

— fb s’obtient donc en corrigeant fb−1 avec une régression sur ses résidus.

Version simplifiée du L2 -boosting

La boucle de l’algorithme de gradient boosting peut se réécrire :
1. Calculer les résidus ui = yi − fb−1 (xi ), i = 1, . . . , n ;
2. Ajuster un arbre de régression pour expliquer les résidus ui par les xi ;
3. Corriger fb−1 en lui ajoutant l’arbre construit.

Interprétation
— On "corrige" fb−1 en cherchant à expliquer "l’information restante" qui est contenue dans les résidus.
— Meilleur ajustement lorsque b % =⇒ biais & (mais variance %).

Logitboost
— Classification binaire avec Y dans {−1, 1} et Ỹ = (Y + 1)/2 dans {0, 1}.
— Log-vraisemblance binomiale de la prévision p(x) ∈ [0, 1] par rapport à l’observation ỹ :

L(ỹ, p(x)) = ỹ log p(x) + (1 − ỹ) log(1 − p(x)).

— Soit f : Rd → R telle que

1 p(x) 1
f (x) = log ⇐⇒ p(x) = .
2 1 − p(x) 1 + exp(−2f (x))

=⇒ re-paramétrisation.
— Chercher p(x) qui maximise L(ỹ, p(x)) revient à chercher f (x) qui minimise son opposé :
 
y+1 y+1
−L(y, f (x)) = − log p(x) − 1 − log(1 − p(x))
2 2
y+1
= log(1 + exp(−2f (x)))+
2  
y+1
1− log(1 + exp(2f (x)))
2
= log(1 + exp(−2yf (x))).

98
Remarque
f (x) 7→ log(1 + exp(−2yf (x))) est convexe et dérivable.

Logitboost
Algorithme de gradient boosting avec la fonction de perte

`(y, f (x)) = log(1 + exp(−2yf (x))).

— Fonction optimale
1 P(Y = 1|X = x)
f ? (x) = log .
2 1 − P(Y = 1|X = x)
— fn estimant f ? , on estime P(Y = 1|X = x) avec
1
.
1 + exp(−2fn (x))

Adaboost
— Remarque : f (x) 7→ exp(−yf (x)) est aussi convexe et dérivable.

Adaboost
Algorithme de gradient boosting avec la fonction de perte

`(y, f (x)) = exp(−yf (x)).

Remarque
— Même nom que l’algorithme initial de [Freund and Schapire, 1996] car quasi-similaire [Hastie et al., 2009].
— Même f ? que logitboost.

Adaboost - version 1
Algorithme [Freund and Schapire, 1996]
Entrées : une règle faible, M nombre d’itérations.
1. Initialiser les poids wi = 1/n, i = 1, . . . , n
2. Pour m = 1 à M :
a) Ajuster la règle faible sur l’échantillon dn pondéré par les poids w1 , . . . , wn , on note gm (x) l’estimateur
issu de cet ajustement Pn
b) Calculer le taux d’erreur : wi 1y 6=g (x )
em = i=1Pn i m i .
i=1 wi

c) Calculer : αm = log((1 − em )/em )

d) Réajuster les poids : wi = wi exp(αm 1yi 6=gm (xi ) ), i = 1, . . . , n
PM
Sorties : l’algorithme de prévision m=1 αm gm (x).

Récapitulatif
— Les principales fonctions de perte pour la régression et classification sont résumées dans le tableau :

Y Perte Prév. optimale

L2 -boosting R (y − f (x))2 E[Y |X = x]
1 P(Y =1|X=x)
Logitboost {−1, 1} log(1 + exp(−2yf (x))) 2 log 1−P(Y =1|X=x)
1 P(Y =1|X=x)
Adaboost {−1, 1} exp(−yf (x)) 2 log 1−P(Y =1|X=x)

— Dans gbm on utilise distribution=

— gaussian pour le L2 -boosting.
— bernoulli pour logitboost.
— adaboost pour adaboost.

99
Profondeur des arbres
— [Link] qui correspond au nombre de coupures =⇒ nombre de feuilles J − 1.
— On parle d’interaction car ce paramètre est associé au degrés d’interactions que l’algorithme peut identifier :
X X X
f ? (x) = fj (xj ) + fj,k (xj , xk ) + fj,k,` (xj , xk , x` ) + . . .
1≤j≤d 1≤j,k≤d 1≤j,k,`≤d

=⇒ [Link]=
— 1 =⇒ premier terme
— 2 =⇒ second terme (interactions d’ordre 2)
— ...
— Boosting : réduction de biais.
— Nécessité d’utiliser des arbres biaisés =⇒ peu de coupures.

Recommandation
Choisir [Link] entre 2 et 5.

Nombre d’itérations
— Le nombre d’arbres [Link] mesure la complexité de l’algorithme.
— Plus on itere, mieux on ajuste =⇒ si on itère trop, on sur-ajuste.
— Nécessité de calibrer correctement ce paramètre.

Comment ?
Avec des méthodes classiques d’estimation du risque.

Sélection de [Link] dans gbm

— gbm propose d’estimer le risque associé au paramètre distribution par ré-échantillonnage :
— [Link] pour du Out Of Bag.
— [Link] pour de la validation hold out.
— [Link] pour de la validation croisée.
— La valeur sélectionnée s’obtient avec [Link].

Exemple

> [Link](321)
> boost.5000 <- gbm(Y~.,data=data_sinus,[Link] = 0.75,
+ distribution="gaussian",shrinkage=0.1,[Link] = 5000)
> [Link](boost.5000)
## [1] 1040
0.5
0.4
Squared error loss

0.3
0.2
0.1

0 1000 2000 3000 4000 5000

Iteration

=⇒ Risque quadratique estimé par hold out avec 75% d’observations dans l’échantillon d’apprentissage.

100
Rétrécissement
— shrinkage dans gbm.
— Correspond au pas de la descente de gradient : shrinkage % =⇒ minimisation plus rapide.

Conséquence
shrinkage est lié à [Link] :
— shrinkage % =⇒ [Link] &.
— shrinkage & =⇒ [Link] %.

Illustration
Erreur test train

0.001 0.01

0.4

0.2

0.0
Err

0.1 1

0.4

0.2

0.0
0 100 200 300 400 500 0 100 200 300 400 500
[Link]

Remarque
Le nombre d’itération optimal diminue lorsque shrinkage augmente.

Recommandation
— Pas nécessaire de trop optimiser shrinkage.

— Tester 3 ou 4 valeurs (0.01, 0.1,0.5...) et regarder les courbes de risque.

— S’assurer que le nombre d’itérations optimal se trouve sur un "plateau" pour des raisons de stabilité.

2.3 Compléments/conclusion
Importance des variables
— Similaire aux forêts aléatoires.
— Score d’impureté :
B
1 X
Ijimp = Ij (Tb ).
B
b=1

— Visualisation avec vip.

Comparaison Boosting/Forêts aléatoires

— Deux algorithmes qui agrègent des arbres :
B
X
fn (x) = αb Tb (x).
b=1

— Indépendance pour les forêts =⇒ Tb se construit indépendamment de Tb−1 .

101
 — Récursivité pour le boosting =⇒ Tb se construit à partir de Tb−1 .

Interprétation statistique
— Boosting : réduction de biais =⇒ arbres peu profonds.
— Random Forest : réduction de variance =⇒ arbres très profonds.
=⇒ les arbres sont ajustés de façon différente pour ces deux algorithmes. =⇒ dans les deux cas, il faut des arbres
"mauvais".

3 Bibliographie
Références

Biblio4

[Breiman, 1996] Breiman, L. (1996). Bagging predictors. Machine Learning, 26(2) :123–140.
[Fernández-Delgado et al., 2014] Fernández-Delgado, M., Cernadas, E., Barro, S., and Amorim, D. (2014). Do we
need hundreds of classifiers to solve real world classification problems ? Journal of Machine Learning Research,
15 :3133–3181.
[Freund and Schapire, 1996] Freund, Y. and Schapire, R. (1996). Experiments with a new boosting algorithm. In
Proceedings of the Thirteenth International Conference on Machine Learning.
[Friedman, 2001] Friedman, J. H. (2001). Greedy function approximation : A gradient boosting machine. Annals
of Statistics, 29 :1189–1232.
[Friedman, 2002] Friedman, J. H. (2002). Stochastic gradient boosting. Computational Statistics & Data Analysis,
28 :367–378.
[Genuer, 2010] Genuer, R. (2010). Forêts aléatoires : aspects théoriques, sélection de variables et applications. PhD
thesis, Université Paris XI.
[Hastie et al., 2009] Hastie, T., Tibshirani, R., and Friedman, J. (2009). The Elements of Statistical Learning :
Data Mining, Inference, and Prediction. Springer, second edition.
[Wright and Ziegler, 2017] Wright, M. and Ziegler, A. (2017). ranger : A fast implementation of random forests
for high dimensional data in c++ and r. Journal of Statistical Software, 17(1).

Discussion/comparaison des algorithmes

Linéaire SVM Réseau Arbre Forêt Boosting

Performance
Calibration
Coût calc.
Interprétation

Commentaires
— Résultats pour données tabulaires.
— Différent pour données structurées (image, texte..) =⇒ performance % réseaux pré-entrainés =⇒ apprentis-
sage profond/deep learning.

102