0% ont trouvé ce document utile (0 vote)
135 vues84 pages

Cours de Machine Learning: Théorie et Pratique

Transféré par

SEMEVO fleur
Copyright
© © All Rights Reserved
Nous prenons très au sérieux les droits relatifs au contenu. Si vous pensez qu’il s’agit de votre contenu, signalez une atteinte au droit d’auteur ici.
Formats disponibles
Téléchargez aux formats PDF, TXT ou lisez en ligne sur Scribd
0% ont trouvé ce document utile (0 vote)
135 vues84 pages

Cours de Machine Learning: Théorie et Pratique

Transféré par

SEMEVO fleur
Copyright
© © All Rights Reserved
Nous prenons très au sérieux les droits relatifs au contenu. Si vous pensez qu’il s’agit de votre contenu, signalez une atteinte au droit d’auteur ici.
Formats disponibles
Téléchargez aux formats PDF, TXT ou lisez en ligne sur Scribd

Machine learning

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

Novembre 2019

Présentation

— Objectifs : comprendre les aspects théoriques et pratiques des quelques algorithmes machine learning.

— 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.

Programme

— Matériel : slides + Notebook R.

— 5 parties :
1. Contexte mathématique pour l’apprentissage (rappels ?) 2h00.
2. Support vector machine : 4h.
3. Agrégation : forêts aléatoires et boosting+arbres (rappels ?) : 4h.
4. Réseau de neurones et introduction au deep learning : 2h.

Table des matières

I Apprentissage : contexte et formalisation 4

1 Motivations 4

2 Cadre mathématique pour l’apprentissage supervisé 7

3 Exemples de fonction de perte 8

4 Estimation du risque 10

5 Le sur-apprentissage 12

6 Le package caret 14

1
7 Annexe : compléments sur les scores 17

8 Bibliographie 19

II Support vector machine 20

1 SVM - cas séparable 21

2 SVM : cas non séparable 25

3 SVM non linéaire : astuce du noyau 31

III Arbres 35

1 Arbres binaires 35

2 Choix des découpes 38


2.1 Cas de la régression . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
2.2 Cas de la classification supervisée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

3 Elagage 40

4 Annexe 1 : impureté, cas multiclasses 45

5 Annexe 2 : algorithme élagage 46

6 Annexe 3 : arbres Chaid 49


6.1 Regroupement des modalités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
6.2 Division d’un nœud . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
6.3 Choix des paramètres . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52

7 Bibliographie 56

IV Agrégation 57

1 Bagging et forêts aléatoires 57


1.1 Bagging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
1.2 Forêts aléatoires . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60

2 Boosting 65
2.1 Algorithmes de gradient boosting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
2.2 Choix des paramètres . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67

3 Bibliographie 71

2
V Réseaux de neurones 73

1 Introduction 73

2 Le perceptron simple 74

3 Perceptron multicouches 77

4 Estimation 79

5 Choix des paramètres et surapprentissage 82

6 Bibliographie 83

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]
3. Comprendre et apprendre un comportement à partir d’exemples.

Constat
— Le développement des moyens informatiques fait que l’on est confronté à des données de plus en plus complexes.
— 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 and Laurent, ]

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 ?? Big data
2017 ?? Intelligence artificielle...

Conclusion
Moyens informatiques =⇒ Data Mining =⇒ Apprentissage statistique =⇒ Intelligence artificielle...

Reconnaissance de l’écriture

Apprentissage statistique
Comprendre et apprendre un comportement à partir d’exemples.

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

4
Reconnaissance de la parole

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 ?

5
Détection de spam

— Sur 4 601 mails, on a pu identifier 1813 spams.


— 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 ?

Problématiques associées à l’apprentissage

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

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

— Règles d’association : mesurer le lien 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].

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

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

[Link] [Link]

6
2 Cadre mathématique pour l’apprentissage supervisé
Régression vs Discrimination

— 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 discrimination ou 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 Discri.
Mot courbe Discri.
Spam présence/absence de mots Discri.
C. en O3 données météo. Régression

Remarque 2.1. — 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 sortie 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 .

7
Approche 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 (mesurable) 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

Aspect théorique
— Pour une fonction de perte ` : Y × Y → R+ donnée, le problème théorique consiste à trouver
f ? ∈ argmin R(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 2.1. — 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.
— On dit que la suite (fn )n est universellement consistante si ∀P
lim R(fn ) = R(f ? ).
n→∞

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.

3 Exemples de fonction de perte


Régression

— Dans un contexte de régression (Y = R), la perte quadratique est la plus souvent utilisée. Elle est définie par :
` : R × R → R+
(y, y 0 ) 7→ (y − y 0 )2
— Le risque pour une fonction de prévision ou régresseur m : X → R est alors donné par
R(m) = E((Y − m(X))2 .

8
Classification binaire

— Dans un contexte de discrimination binaire (Y = {−1, 1}), la perte indicatrice est la plus souvent utilisée.
Elle est définie par :
` : {−1, 1} × {−1, 1} → R+
(y, y 0 ) 7→ 1y6=y0
— Le risque pour une fonction de prévision ou règle de prévision g : X → {−1, 1} est alors donné par
R(g) = E(1g(X)6=Y ) = P(g(X) 6= Y ).

Fonction de score

— On reste dans un cadre de classification binaire (Y = {−1, 1}).


— Mais... plutôt que de chercher une règle de prévision g : X → {−1, 1}, on cherche une fonction S : X → R
telle que
P(Y = 1|X = x) faible P(Y = 1|X = x) élevée

S(x)

— Une telle fonction est appelée fonction de score : plutôt que de prédire directement le groupe d’un nouvel
individu x ∈ X , on lui donne une note S(x)
— élevée si il a des "chances" d’être dans le groupe 1 ;
— faible si il a des "chances" d’être dans le groupe -1 ;

Courbe ROC et AUC

— On utilise souvent la courbe ROC pour visualiser la performance d’un score :



x(s) = α(s) = 1 − sp(s) = P(S(X) > s|Y = −1)
y(s) = 1 − β(s) = se(s) = P(S(X) ≥ s|Y = 1)
— On déduit de ce critère un risque pour les scores en considérant l’aire sous la courbe ROC (AUC) :
R(S) = AUC(S).

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

0.75
true_positive_fraction

Scores
Perfect

0.50 random
S1
S2

0.25

0.00

0.00 0.25 0.50 0.75 1.00


false_positive_fraction

> library(pROC)
> df1 %>% group_by(Scores) %>% summarize(auc(D,M))
## # A tibble: 4 x 2
## Scores ‘auc(D, M)‘
## <chr> <dbl>
## 1 Perfect 1
## 2 random 0.5
## 3 S1 0.896
## 4 S2 0.699

9
Perte `(y, f (x)) Risque E[`(Y, f (X))] Champion f ?

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

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

Scoring AU C(S) P(Y = 1|X = x)

Résumé

4 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

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.

10
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ées. Dn : données, {A, V} : partition de {1, . . . , n}.
1. Construire l’algorithme de prédiction sur Dn,app = {(Xi , Yi ) : i ∈ A}, on le note fn,app ;
b n (fn ) = 1 P
2. Calculer R |V| i∈V `(Yi , fn,app (Xi )).

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

Validation croisée K-blocs

— Principe : répéter l’algorithme apprentissage/validation sur différentes partitions.

Algorithme - CV
Entrées. Dn : données, K un entier qui divise n ;
1. Construire une partition {I1 , . . . , IK } de {1, . . . , n} ;
2. Pour k = 1, . . . , K
(a) Iapp = {1, . . . , n}\Ik et Itest = Ik ;
(b) Construire l’algorithme de prédiction sur Dn,app = {(Xi , Yi ) : i ∈ Iapp }, on le note fn,k ;
(c) En déduire fn (Xi ) = fn,k (Xi ) pour i ∈ Itest ;
3. Retourner
n
b n (fn ) = 1
X
R `(Yi , fn (Xi )).
n i=1

Commentaires

— Plus adapté que la technique apprentissage/validation lorsqu’on a peu d’observations.


— Le choix de K doit être fait par l’utilisateur (souvent K = 10).

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.

11
5 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.
— 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 le plus souvent 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 sur les données d’apprentissage, il aura du mal à s’adapter
à de nouveaux individus.

Train error
Test error

OVERFITTING

Complexity (λ)

Overfitting en regression

12



● ● ●


● ● ●
● ●
1.0 ● ●
●●


●●
● ● ● ● ●● ● ●
● ●
● ● ●
● ● ●
● ● ●

● ● ● ●
● ● ●
● ●
●● ●
● ●
● ● ●

● ● ●
● ●


0.5 ● ● ●
● ● ●

● ●● ●
● ● ●
● ● ● ●

● ●
● ●
● ● ●
● ●


●● ● ●

0.0 ● ● ● ●
● ●

Y
● ●
● ●
● ● ●
● ● ● ●

● ● ● ● ●
●● ● ●


● ●
● ● ●

−0.5 ● ● ●
● ●
● ● ● ●


●● ● ●

● ●
●●
● ● ●●● ●
●●
●●● ●
●●
● ● ●●
● ●●
●● ●
●● ● ●


−1.0 ●
● ●
● ●

● ●●
● ●
● ●


−1.5

−4 0 4
X

1.0

0.5

Estimate
0.0 Over
y

Good
Under

−0.5

−1.0

−1.5

−4 0 4
x

Overfitting en classification supervisée

1.00 ● ● ● ●●● ● ●● ● ● ● ● ●
● ● ● ● ● ●●● ● ● ● ● ●

● ●● ● ● ● ●
● ● ● ●●
●●
● ●
● ● ●● ● ● ● ● ●●
● ● ● ● ●● ● ● ● ● ● ●
●● ● ● ● ● ●● ● ● ● ● ● ●
● ● ●● ● ●●●● ●

●● ● ● ●● ● ● ●● ●●● ● ● ● ● ●●
● ● ● ●●
● ●
● ● ● ● ● ● ●
●●
● ●● ● ● ● ●
● ●● ● ●●
● ● ●●● ● ●● ● ● ●● ● ● ● ●●
● ● ●
●●
● ●●●
● ● ● ●● ● ● ● ● ● ● ● ● ● ● ●●
● ● ●● ● ● ● ● ● ● ● ● ● ● ●● ●
● ● ●● ● ● ●● ● ● ●● ● ●● ●
● ●● ● ● ● ●
● ●● ●
● ● ● ● ● ●
● ● ● ● ● ●
● ● ● ● ● ● ● ●● ● ●



● ● ● ● ● ●
● ●
● ● ● ●
● ● ●●●
● ● ● ● ● ● ● ● ● ●● ● ● ●●● ●●
● ● ● ●
●● ●●
● ● ● ● ● ● ●
● ● ●
●● ● ● ● ● ● ● ●● ●
● ● ●

● ● ●●

● ● ● ● ● ● ● ● ● ● ● ●
● ●● ●● ● ● ● ● ●●
●● ● ●● ●● ● ● ● ●●
●● ● ● ● ● ●● ● ● ● ● ●

● ● ●
● ● ● ●● ● ● ● ●
● ●● ● ● ● ● ● ●● ● ● ● ●●

● ● ●● ●● ●● ● ● ● ●
● ● ● ● ● ●● ● ● ● ● ● ● ●●●
●● ●

● ● ●● ● ●
● ● ● ●
● ● ● ●● ● ●

●●● ●
● ●● ● ● ● ● ●● ●
●● ● ● ●●● ●
●● ● ●● ● ●
0.75 ● ● ● ● ●●
● ● ● ●●
● ●
● ● ●

●● ●● ● ●
● ●● ● ● ●● ●
● ● ● ●
●● ● ●●●
● ●●●● ● ● ● ●
● ● ● ●● ●

●● ● ● ● ● ●● ● ● ● ● ● ●
● ●● ● ●●●

● ●●● ● ● ● ● ● ●● ● ●●
● ● ● ● ●● ●
● ● ● ● ●
● ● ● ●● ● ● ● ● ● ● ● ●
● ● ● ●
● ● ●● ●● ● ● ● ● ●● ●
● ● ● ● ●●


●● ● ● ● ● ● ● ● ● ● ●
● ●● ● ●
● ● ●●● ● ● ●● ● ● ● ●
● ●●● ● ● ● ● ●●
● ●
●● ● ● ●●
●● ● ● ● ●

● ●● ●● ● ● ● ● ● ● ● ● ● ●
●●
●● ● ● ● ●

●● ● ● ● ● ● ●● ● ●● ●

● ●● ● ● ● ●
● ● ●● ● ●



● ● ● ● ● ● ● ●● ●
● ● ● ● ●●
● ● ● ● ● ●

● ●● ●
● ● ● ● ●
● ● ● ● ●● ● ●● ●


● ● ● ● ● ●
● ● ● ●● ● ● ● ●
● ● ●●
●● ● ● ●
●●
● ● ● ●

● ● ● ● ● ● ● ●
● ● ● ● ● ● ● ● ● ●● ● ● ● ● ● ● ●

● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ●

●● ●● ●● ●

● ●
● ●●

●●

●● ●
● ● ●●

●●●
● ●● ●
● ●
● ●
● ●● ● Y
● ● ● ● ●● ● ●
● ● ● ●● ● ● ● ● ● ● ● ●
● ● ● ● ● ● ● ● ● ● ● ● ● ●● ●● ●
● ● ● ● ● ●
X2

● ● ●
0.50 ●● ● ●● ● ● ● ● ● ● ●
● ● ● ●● ● ●
● ●
● ● ● 0
● ● ●●
● ● ●
● ● ●
● ●● ● ● ● ● ● ●● ● ●
● ● ● ● ● ●● ● ●
●●
● ●
● ● ●● ● ● ● ●●
● ● ● ● ●●●
● ●
● ● ● ●●● ● ● ● ● ● ●● ● ● ●●
● ●
●●
● ●
●● ● ● ●
● ● ●● ● ● ● ● ● 1
● ● ●● ●● ● ● ● ● ● ● ● ●● ● ●● ● ● ● ● ●
●● ● ●
● ●● ● ● ● ● ●●●● ● ●

● ●
● ●
●● ●
● ●
● ●

● ●

● ● ● ●● ● ●● ●● ●●● ● ●●
● ●●●
● ● ● ●●
● ● ● ● ●
● ● ● ● ● ●● ● ● ● ● ●


●● ● ●
● ●
● ● ●
●● ● ● ● ● ● ●● ●● ●●
●● ● ● ● ●● ●●●● ●
● ●
● ● ● ● ● ● ● ● ● ● ● ●
● ● ● ● ● ● ● ●● ●● ● ● ●
● ●● ●● ● ●●● ● ●
●● ●

●●●● ● ● ● ● ● ● ●● ● ●● ●
● ● ● ●
● ●● ● ● ●● ●● ● ● ●
● ●

● ● ● ●●● ●● ● ● ● ●
● ●●
● ● ●● ●
●● ● ●● ● ●● ● ● ●●● ● ●● ● ●
● ● ● ●
●● ● ●●● ● ● ●

● ●● ●

● ●

● ● ● ●● ● ● ● ●● ● ●●● ●
●●●● ●
● ● ● ●● ●
● ● ●
● ● ● ●● ● ● ● ● ●● ● ●●
● ● ● ●

● ● ●

● ●● ●
● ● ● ● ● ● ● ● ●

● ● ●● ● ● ● ● ● ● ● ● ●● ●● ●● ● ●
●● ● ● ● ● ●
● ● ● ● ●● ● ● ● ●
● ● ●
● ● ● ●● ●● ●●● ● ● ●
● ● ●● ●
0.25
● ● ● ● ●● ● ● ● ● ● ●
● ●●
● ●● ● ● ● ● ●● ●
● ● ● ● ● ●● ●
●● ● ●
● ●●
●●

● ● ● ● ●

●● ●● ● ●● ●
● ● ● ● ● ●● ● ●●
● ●

● ●
●● ●● ● ● ● ● ● ● ● ● ●
● ● ● ● ● ● ● ●

● ● ●
● ●● ●
●● ● ● ● ● ● ● ● ●● ● ● ● ● ● ●
● ●● ● ● ● ●


● ● ● ● ● ●
● ● ● ● ●● ●●
● ● ● ● ● ●●
● ● ● ●● ● ●● ● ●
● ● ● ●

●● ● ● ●
●● ●
● ● ●● ● ● ●● ● ● ● ● ● ●
● ● ● ● ●● ● ● ● ●●
● ●● ● ● ●● ● ● ● ●● ● ● ●
●● ● ● ● ●
● ● ● ● ● ● ●●
● ● ●● ●
● ● ● ● ● ● ●● ● ● ● ● ● ● ●
● ● ● ●● ● ● ●● ● ●● ●● ● ● ●
●● ● ● ●● ● ●● ● ● ● ● ●● ● ●
● ●● ● ● ●● ● ●
● ● ● ●● ●● ● ●●
●● ●●● ● ● ● ● ●


●●● ●
●● ● ● ● ●
● ●● ● ● ● ● ● ● ● ● ● ● ●
● ● ●● ● ● ● ● ● ● ●
● ●● ● ● ● ● ●
● ● ● ●● ● ● ● ●●
●● ● ● ● ● ●●● ●● ● ● ●
● ●● ●●● ●● ●● ● ●
● ● ● ● ● ● ● ● ● ● ●●

● ● ●● ● ●● ●● ●● ●● ●● ●
● ● ● ● ● ● ● ● ● ● ● ● ● ●●
●● ● ● ● ● ● ● ● ●● ● ●
● ● ● ●● ●●
● ● ● ●●●
● ●● ● ● ● ●
● ● ● ● ●● ● ● ●● ●●● ● ● ● ●
● ●● ●
● ●●

● ● ● ●●●

● ●● ● ● ●
● ●
●● ● ● ● ● ●●
● ●● ●● ● ●
● ● ● ● ● ● ● ● ● ● ●
0.00 ● ● ●

0.00 0.25 0.50 0.75 1.00


X1

13
1.00 1.00

0.75 0.75

label label

X2

X2
0.50 0 0.50 0
1 1

0.25 0.25

0.00 0.00

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

1.00 1.00

0.75 0.75

label
label
X2

X2
0.50 0 0.50
1
1

0.25 0.25

0.00 0.00

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

6 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
> 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

14
## 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.

> plot(e1)

● ●
0.72 ● ●
● ● ● ● ● ● ●

● ● ●

0.70
Accuracy (Repeated Train/Test Splits)

● ●

● ●
● ●

0.68 ●

0.66

0.64

0.62 ●

0 100 200 300 400 500

#Neighbors

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

## 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

15
## 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.

Critère AUC

16
> 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.

7 Annexe : compléments sur les scores


Scores parfait et aléatoire

S ●● ●●
● ● ● ● ●●● ●

● ● ● ● ● ●
●● ● ● ●
● ●
● ●●● ● ●● ●● ● ●● ●

● ●●●●●

Y
random ●● ●
● ●● ●●● ●
● ●
●●
● ● ●
● ●
● ●
● ●● ● ● ●● ●
● ● ●●
●●●● ●●
●●
● ●● ● ● ● ● ● ● 0
● 1

Perfect ●●●●●
●● ● ● ● ●●● ●

● ● ● ● ● ●● ● ● ●
● ●
● ●●● ●
● ●● ●● ● ●● ●

● ● ●●

0.00 0.25 0.50 0.75 1.00


Scores

Définition 7.1. — Score parfait : il est tel qu’il existe un seuil s? tel que
P(Y = 1|S(X) ≥ s? ) = 1 et P(Y = −1|S(X) < s? ) = 1.
— Score aléatoire : il est tel que S(X) et Y sont indépendantes.

17
Lien score/règle de prévision

— Etant donné un score S, on peut déduire une règle de prévision en fixant un seuil s (la réciproque n’est pas
vraie) : 
1 si S(x) ≥ s
gs (x) =
−1 sinon.
— Cette règle définit la table de confusion
gs (X) = −1 gs (X) = 1
Y = −1 OK E1
Y =1 E2 OK
— Pour chaque seuil s, on distingue deux types d’erreur
α(s) = P(gs (X) = 1|Y = −1) = P(S(X) ≥ s|Y = −1)
et
β(s) = P(gs (X) = −1|Y = 1) = P(S(X) < s|Y = 1).

On définit également
— Spécificité : sp(s) = P(S(X) < s|Y = −1) = 1 − α(s)
— Sensibilité : se(s) = P(S(X) ≥ s|Y = 1) = 1 − β(s)

Performance d’un score


Elle se mesure généralement en visualisant les erreurs α(s) et β(s) et/ou la spécificité et la sensibilité pour tous les
seuils s.

Courbe ROC

— Idée : représenter sur un graphe 2d les deux types d’erreur pour tous les seuils s.
Définition 7.2. C’est une courbe paramétrée par le seuil :

x(s) = α(s) = 1 − sp(s) = P(S(X) > s|Y = −1)
y(s) = 1 − β(s) = se(s) = P(S(X) ≥ s|Y = 1)

Remarque
— La courbe ROC d’un score parfait passe par le point (0,1).
— La courbe ROC d’un score aléatoire correspond à la première bissectrice.

1.00

0.75
true_positive_fraction

Scores
Perfect

0.50 random
S1
S2

0.25

0.00

0.00 0.25 0.50 0.75 1.00


false_positive_fraction

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

18
AUC
Définition 7.3. — L’aire sous la courbe ROC d’un score S, notée AU C(S) est souvent utilisée pour mesurer
sa performance.
— Pour un score parfait on a AU C(S) = 1, pour un score aléatoire AU C(S) = 1/2.
Proposition 7.1. — Etant données deux observations (X1 , Y1 ) et (X2 , Y2 ) indépendantes et de même loi que
(X, Y ), on a
AU C(S) = P(S(X1 ) ≥ S(X2 )|(Y1 , Y2 ) = (1, −1)).

AUC

> library(pROC)
> df1 %>% group_by(Scores) %>% summarize(auc(D,M))
## # A tibble: 4 x 2
## Scores ‘auc(D, M)‘
## <chr> <dbl>
## 1 Perfect 1
## 2 random 0.5
## 3 S1 0.896
## 4 S2 0.699

Score optimal

— Le critère AU C(S) peut être interprété comme une fonction de perte pour un score S ;
— Se pose donc la question d’existence d’un score optimal S ? vis-à-vis de ce critère.
Théorème 7.1 ([Clémençon et al., 2008]). Soit S ? (x) = P(Y = 1|X = x), on a alors pour toutes fonctions de
score S
AU C(S ? ) ≥ AU C(S).

Conséquence
Le problème pratique consistera à trouver un "bon" estimateur Sn (x) = Sn (x, Dn ) de

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

8 Bibliographie
Références

Biblio1

[Besse and Laurent, ] Besse, P. and Laurent, B. Apprentissage Statistique modeélisation, preévision, data mining.
INSA - Toulouse. [Link]
[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.
[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.

19
Deuxième partie

Support vector machine


Cadre et notation

— Discrimination binaire : Y à valeurs dans {−1, 1} et X = (X1 , . . . , Xp ) dans Rp .

Objectif
— Estimer la fonction de score S(x) = P(Y = 1|X = x) ;
— En déduire une règle de classification g : Rp → {−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 + . . . + wp Xp .
— Règle associée : 
1 si w1 X1 + . . . + wp Xp ≥ 0
g(x) =
−1 sinon.

Exemple 1 : régression logistique

— Modèle :
p(x)
logit = β0 + β1 x 1 + . . . + βp x p
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 + . . . + βp xp ≥ 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 + x0 Σ−1 (µ1 − µ0 ) ≥ 0

1
g(x) =
−1 sinon.

Illustration avec p = 2

— 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.

20
1.00 ● 1.00 ●
● ●

● ●

0.75 ● 0.75 ●

● ●

● ● ● ● ● ●
● ●
● ● ● ●
● ●
● label ● label
X2

X2
● ● ● 0 ● ● ● 0
0.50 ● 0.50 ●
● ● ● ●
● 1 ● 1
● ● ● ●
● ● ● ●

● ●
● ● ● ●
● ●
● ● ● ●

0.25 ● 0.25 ●
● ●
● ●
● ●

● ● ● ● ● ● ● ●
● ● ● ● ● ●
● ●
● ●
● ● ● ●

● ● ● ●
● ●
0.00 0.00

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

Figure 1 – Règles logistique (gauche) et lda (droite).

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
[Link]

Remarque
Les aspects techniques ne seront pas présentés ici, on pourra les trouver à l’url [Link]
p5awah6ij8ykiil/slides_apprentissage.pdf?dl=0

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) ∈ Rp × 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.

21
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.

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.

22
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? = 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.

Représentation

hw? , xi + b? = 0

M = kw1? k

α?
i > 0

marge

M = kw1? k

23
Le coin R

— La fonction svm du package e1071 permet d’ajuster des SVM.

1.00

● ●



0.75




Y

X1
0.50 ● 0

● ●
● 1

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 ;
— la marge.

1.00

● ●



0.75




Y
X1

0.50 ● 0

● ●
● 1

0.25 ●

● ●

0.00
0.00 0.25 0.50 0.75
X2

24
— La fonction plot donne aussi une représentation de l’hyperplan séparateur.
> plot([Link],data=df,fill=TRUE,grid=500)

SVM classification plot


o
o o
o
oo

0.8 o

1
o

0.6
x
o
o
X1

x o

0.4 o

o
o

0
o
0.2

o
x o
0.2 0.4 0.6 0.8

X2

2 SVM : cas non séparable

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

25
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)

26
— 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.
— 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
1


0.25 ●

● ●
0.00
0.00 0.25 0.50 0.75
X2

27
> 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

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

SVM classification plot


o
o o
o
oo

0.8 o

1
x

0.6
x
o
o
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

28
1.00 ●
● ● ● ● ●

● ●
● ●
● ● ●●
●● ● ● ● ●● ● ● ●
●● ●●● ● ● ● ● ● ● ●● ● ● ● ●

● ● ● ● ●
● ●
● ● ● ●
●● ● ● ●● ●
● ●● ● ● ● ●
●● ● ● ●
● ● ●● ●
● ●
●● ● ●
● ●● ● ●
● ● ● ●●

● ●● ● ● ● ●● ●●
● ● ● ● ● ● ● ●●
● ● ●
● ●
● ● ● ●
● ● ●
● ● ● ●● ● ●● ●
● ● ●
● ● ●● ● ● ●
● ● ● ● ● ● ●● ●
● ● ● ● ●
● ● ● ● ●
● ● ● ● ●● ● ●

● ● ● ● ● ● ●●●
● ● ● ● ● ● ●

● ● ● ●● ● ●●

● ● ● ●
● ● ● ●● ● ● ● ● ●
● ●

● ● ● ●● ●
● ● ● ●
● ● ● ●● ● ●
●● ● ● ● ●● ● ●
● ●● ● ● ● ● ● ● ●

● ●
0.75 ●

● ●●●● ●
● ● ●

● ● ● ●● ● ● ●

● ● ● ● ● ●●
● ●

● ● ● ●● ●
● ● ● ●
● ● ● ● ●
● ● ●

● ●● ● ●
● ● ● ● ● ●
● ●
● ● ● ●● ● ●● ●
●● ● ● ● ● ● ●
● ●
● ●● ●
●● ● ●● ●● ●
● ● ● ● ●
● ● ● ● ● ● ● ● ● ●
● ● ● ● ● ● ●
● ● ● ● ●● ● ● ●● ● ● ●
●● ● ● ●
● ● ●● ● ● ●


● ● ● ●● ● ● ● ● ●● ● ●●
● ● ● ● ●
● ● ● ● ● ●
● ●
● ● ● ● ●● ●
●●● ● ● ● ● ● ●● ● ● ● ●
● ●● ●
● ● ● ●
●● ●
● ●
● ● ●
● ● ●●
● ●● ● ● ●●

● ● ●

Y
● ● ● ●
● ● ● ● ● ● ● ●
● ●● ● ● ● ● ●●

V1
● ● ● ●● ● ● ● ● ●
● ●● ● ● 0
0.50 ●● ● ● ●
● ● ● ●



● ● ● ● ●● ●
●● ● ● ●●
●●● ● ● ●
●● ● ● ● ● ● ●

● ● 1
● ● ●
● ● ●● ● ● ● ●
● ●
●● ● ● ●
● ● ● ● ●
● ●● ● ●
● ● ● ● ● ●● ●
● ●
● ●
● ● ●● ● ●
● ● ● ● ● ●●
● ● ● ● ● ● ● ● ● ●
● ● ●
● ● ● ● ● ● ● ● ● ● ●
● ● ● ● ●
● ● ● ●
●●
● ●
● ● ● ● ● ●● ●
● ● ●● ● ● ●
● ● ● ● ●
● ●● ● ●●
● ● ● ● ● ●●
● ●
● ●●● ● ● ● ● ● ● ●● ●●●
● ●
● ● ●
● ● ● ● ●
● ●
●● ● ● ● ● ● ● ● ●
● ● ● ● ● ● ● ●
● ●
● ●● ●
● ●● ● ●● ●
● ●



● ● ●● ● ●
● ●● ●
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)

SVM classification plot SVM classification plot SVM classification plot


o o o oooooo oooo x x x x x o o o oooooo oooo o o o o x o o o oooooo oooo o o o o x
x x xxx xxxxx x xx xx x xx x o oooooox ooooox ooooo xx x xx x o oooooox ooooox ooooo ox x xx
o
x ooo oox x x xxx xx x xx x xxxx x xxxxx x x
xxxx x xx x x
o
o ooo ooo o o o ooooo ox o oo
o oooo x
xxxx x xx x x
o
o ooo ooo o o o ooooo ox o oo
o oooo x
xxxx x xx x x
x x x
x ooo o xxx xx xx xx xxx x xx x xxxx xx xxx x xxx x ooo o ooo oo x oo ooooo x oo o oooo xx xxx x xxx x ooo o ooo oo x oo ooooo x oo o oooo ox xxx x xxx
oo x oo o oo oo o oo
o x x x xx x x xxx
xx x o o xo o x xxx
xx x o o xo o o xxx
xx x
o o ooo o x xx xxx x x x xxx x xx x o o ooo o ooo ooo o o o ooo x xx x o o ooo o ooo ooo o o o ooo x xx x
o xx x x x x x xx xxx x xxx x x x o oo o o o o o xx xxx x xxx x x x o oo o o o o o xx xxx x xxx x x x
o x x x x xxxx xx x x x xx x
x x x xxx x xxx x xxx x x x o o o o o oooo ox o x o oo o
o o o xxx x xxx x xxx x x x o o o o o oooo ox o x o oo o
o o o oxo x xxx x xxx x x o
0.8 o o x x xx xxx x x xx x xx 0.8 o o x o oo ooo o o oo x xx 0.8 o o x o oo ooo o o oo x xx
x xx x x xx x xx x xxx xxx x x o oo o ooo o xo o xxx xxx x o o oo o ooo o xo o xox xxx x o
xxxx xx x xx xx xx xx xx xx x xx xx oooo ooo oo xx xx xx xx xx x xxoo oooo ooo oo xo ox xx xx xx x oooo
1

1
o x x xxx xx x x xx o o x xxx xx x x oo o o x xxx xx x x oo
x x x x xx x xxx xx x o o o o oo o ooo xx o o o o o oo o ooo xx o
x x x xx x xxx x x x xxx x xx x xx x x o o o ox x xxx x x x xxx x xx x oo o o o o o ox o xxx x x x xxx x xx x oo o o
xxxx xxx x xx x xx xx x x x x ooxo xoo o oo x xx xx x x x o ooxo xoo o oo x xx xx x x x o
x xx xx o oo ox o oo ox
x xx xxx x x xx x x xxx xxxx xxxx x xx x x x xx x x x o oo oooo o oo o o ooox oooo xxxx x xx x x x xx x o o o oo oooo o oo o o ooox oooo ooxx x xx x x x xx o o o
x
x x x x x x x xx x x o
o o o x x x x xx o o o
o o o x x x x oo o o
0.6 x x x x xxxx x xx x x x x xx x x xx 0.6 o oo o oooo oo o ooo oox xx x x x x xxx x xx x oo oo 0.6 o oo o oooo oo o ooo ooooo x x x x xxx xox o oo oo
xx x x xx x x x x xxx x x x x xxx x x x xxx xxx oo x o o x
x x x x x oo o o oo x o o x
x x x x x oo o o
xx x x x xxxx x x oo o x x x oooo o o oo o x x o oooo o o
x x x xxx xxx x x x x x x xx o o o oox oxo x x x o x o oo o o o oox oxo x x x o x o oo
x xxx xx x x x x xx xx x x x xx xxx x x xxxxx x o ooo ox o x x x xx xx x x x oo ooo o ooooxo o o ooo ox o o o x xx xx x x o oo ooo o ooooxo o
V1

V1

V1

xx xxx xxx x x x x x xxx x x xxx xx x xx x x oo oxo ooo o o o x x xxx x x xxx xx x oo o o oo oxo ooo o o o x x xxx x x xxx xx o oo o o
xx x xx x x x x xx xx x xx oo o oo o x x x xx oo o oo oo o oo o o x x ox oo o oo
x x xxxxxxx xx x x x o o ooooooo oo x o o o o ooooooo oo x o o
x xx x xxx xx xx xx x x x x x
xxxx xxx xxx o oo o xxx xx xx xx x x x x o
oooo oxx ooo o oo o oxx xx xx xx x x x x o
oooo oxx ooo
0.4 xxx x xx x x 0.4 ooo x xx o o 0.4 ooo x xo o o
x x xx xx x x x x xxxx xx xx x x x o o oo xx x x x x xxxx xx xo o o x o o oo oo x x x x xxxx xx xo o o x
x x x x xxx x x x xx x x o o o x xxx o o o ooo o o o o x xxx o o o ooo o
x x x xx xx x xxx x x x x xx xxx x x xx x x o o ooo xx x xxx x x x o oo ooo o o oo o o o o ooo xx x xxx x x o o oo ooo o o oo o o
xxxx xx x xx x xx x x x xx x x xx x x x x oooo ox x xx x xx x o o oo o o oo o o o o oooo ox x xx x xx x o o oo o o oo o o o o
x x x x xxxx xxx x x x x xxxxx x x x x x xxx x x x o ox x x xxxx xxx x x x x xooxo o o o o o ooo o o o o ox o o xxxx xxx x x xooooooo o o o o o ooo o o o
x xx xx x x xxx x xxxx x x x oo oo o o ooo o xooo o o o oo oo o o ooo o xooo o o
xxx xx x x xx xx xx x o oo xo xx x o oo xo
0

ooo
0

ooo

0
x x x o x o
xxxxxx x x x xx xx x x x x x x
x xx x x
xx
x oooooo x x x xx xx x x o o o o
o oo o o o oooooo o x x xx xx x x o o o o
o oo o o o
0.2 x xxx xxx xx x x xxx x x xxx xx x xx x xx xx xxxx 0.2 o xxx xxx xx x x xxx x o ooo ooo oo o oo oxoooo o ooo 0.2 o xxx xxx xx x x xxx x o ooo ooo oo o oo oxoooo o ooo
x x x xx xx x x x x x x x x x xx oo o o o o x x x xx oo o o o o
x x xx x x x x o o oo o o o o o o oo o o o o
x xx x xx xx x xxxxx x xx x xx x xx x xx x xx xx x ooooo o oo x oo o oo x xx x xx xx x ooooo o oo x oo o oo
x xx x xxx xx xxxx xx x xx xxx xx x x xx xxx x xx x xxx xx xooo oo o oo xooo oo x o oo ooo x xx x xxx oo xooo oo o oo xooo oo x o oo ooo
x xxxx x xx xx xx xx x x x xx x x xx x xxxx x oo oo oo oo o o o oo o o oo x xxxx o oo oo oo oo o o o oo o o oo
x x x x x x x x x x x x x o x o o o o o o oo x o o x o o o o o o oo
xxx x xxx xxxx x x x xxx xx x
x
xx x x xxxx xxx x xxx oooo o o o ooo xo o oo o o oooo xxx x xxx oooo o o o ooo xo o oo o o oooo
x x x x x xx x x x x x o o o oo o o o x x o o o oo o o o
xx xx xx x x xo oo oo o o xo oo oo o o
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

Un autre exemple

— n = 1000 observations.

1.00 ● ●●
●● ● ● ●
●●
●● ● ● ● ● ●


● ●●● ● ● ●● ●
●● ● ● ●
● ● ● ●
●● ● ● ●

● ● ● ●
● ● ● ●
● ●
●●● ● ● ● ●●
● ● ● ● ● ● ● ● ●
● ● ● ● ●● ● ● ● ● ● ● ● ●●
● ●
● ● ●
● ● ● ● ●●
● ●
●● ● ● ●●
● ● ● ●●
● ● ●
● ● ● ●
● ●

● ●
● ●
● ● ● ● ●
● ● ●
● ●
● ●
● ● ●

0.75 ●
● ● ●

● ●
● ●

● ● ●
● ● ●
● ● ● ● ●
● ●
● ● ● ●
● ●
● ● ●

● ● ● ●
●● ● ●
● ●
● ● ●
● ● ● ●
● ● ●
● ●
● ● ● ●

● ●
● ● ●
● ●
● Y
● ● ●

● ●
X2

0.50 ● ● ● ● ● ● 0
● ● ● ●
● ●
● ●
● ● ● 1
● ● ● ● ● ● ●
● ●● ●
● ●
● ●
● ● ● ● ● ●
● ● ●
● ● ●
● ● ●● ● ●
● ●

● ●● ● ●
● ●
● ●
● ● ● ● ●
● ● ● ●

● ●
● ● ●
● ● ● ● ●●
● ●
● ●
0.25 ● ● ● ●
● ●

● ● ● ●
● ● ●
● ●
● ● ●
● ●
● ●● ● ●
● ● ●
● ● ●
● ● ● ● ●
● ●
●●
● ● ● ●
● ● ● ●
● ●
● ● ● ● ●

●● ● ●●
● ● ●
● ● ●
● ● ● ● ●
● ● ●
● ● ● ●
● ● ● ●
● ● ●●
● ●
0.00

0.00 0.25 0.50 0.75 1.00


X1

29
> model1 <- svm(Y~.,data=donnees,cost=0.001,kernel="radial",gamma=5)
> model2 <- svm(Y~.,data=donnees,cost=1,kernel="radial",gamma=5)
> model3 <- svm(Y~.,data=donnees,cost=100000,kernel="radial",gamma=5)

1.00 1.00 1.00

0.75 0.75 0.75

0.50 0.50 0.50

0.25 0.25 0.25

0.00 0.00 0.00

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

Choix de C avec tune

> [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.071
## - Detailed performance results:
## cost error dispersion
## 1 1e-03 0.142 0.03675746
## 2 1e-02 0.084 0.03373096
## 3 1e+00 0.071 0.02766867
## 4 1e+01 0.072 0.02820559
## 5 1e+02 0.072 0.02820559
## 6 1e+03 0.071 0.02766867

> 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
## gamma: 0.5
##
## Number of Support Vectors: 336
##
## ( 168 168 )
##
## Number of Classes: 2
##
## Levels:
## 0 1

Choix de C avec caret

> library(caret)
> gr <- [Link](C=c(0.001,0.01,1,10,100,1000))
> ctrl <- trainControl(method="repeatedcv",number=10,repeats=5)
> train(Y~.,data=df3,method="svmLinear",trControl=ctrl,tuneGrid=gr)

30
## Support Vector Machines with Linear Kernel
## No pre-processing
## Resampling: Cross-Validated (10 fold, repeated 5 times)
## Summary of sample sizes: 900, 900, 900, 900, 900, 900, ...
## Resampling results across tuning parameters:
## C Accuracy Kappa
## 1e-03 0.8700 0.7377051
## 1e-02 0.9188 0.8369121
## 1e+00 0.9304 0.8604317
## 1e+01 0.9292 0.8580356
## 1e+02 0.9294 0.8584333
## 1e+03 0.9294 0.8584333
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was C = 1.

3 SVM non linéaire : astuce du noyau


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

0.75

X1^2
0.50

1.0

0.25


● ●●
0.5

● ●
● 0.00 ●● ● ●● ●●

0.00 0.25 0.50 0.75 1.00


● ●
X2^2
● ● ●
X1

0.0 ●● ● ●


1.00

● ● ●

−0.5

0.75

−1.0

−1.0 −0.5 0.0 0.5 1.0


X2
X1^2

0.50

0.25


● ●●


0.00 ●● ●● ●● ●●

0.00 0.25 0.50 0.75 1.00


X2^2

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

Noyau
Définition 3.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

31
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 3.1. 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

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

32
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 3.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
1.00
1.0

0.75
0.5

● ●

X1^2
● ● 0.50
X1

0.0 ●● ● ●

● ● ●

0.25
−0.5


● ●●


−1.0 0.00 ●● ●● ●● ●●

−1.0 −0.5 0.0 0.5 1.0 0.00 0.25 0.50 0.75 1.00
X2 X2^2

— 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 (sur Rd ) : K(x, x0 ) = xt x0 .


2. Polynomial (sur Rd ) : K(x, x0 ) = (xt x0 + 1)d .
3. Gaussien (Gaussian radial basis function ou RBF) (sur Rd )
kx − x0 k
 
K(x, x0 ) = exp − .
2σ 2

4. Laplace (sur R) : K(x, x0 ) = exp(−γ|x − x0 |).


5. Noyau min (sur R+ ) : K(x, x0 ) = min(x, x0 ).

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.

33
> 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 o o o
x o o
x o x o x o
o o o

0.5 o 0.5 o 0.5 o

1
x x o x x o
x x x x x o
x x o x o x

x x o x x x x o o
x x x o x o o o x
x o o
X1

X1

X1
0.0 xx x x 0.0 xx o o 0.0 xo o o
x x x
x o o

x x x x o x x x x o x x x x o
x x o x x x
o x o
x o o
0

0
−0.5 x x −0.5 x o −0.5 o o
x x o x o x
x x x
o x o

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

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

SVM classification plot SVM classification plot


x o

2 o x x o oo

1.5 xx xx x
x 1.5 oo o
o
oo
x x o o
xx oo

● o x xxx xx o o ooo oo
x xxxx xx x o oooo oo o
●● ●
ooxx x x x x oooo o o o o


●● x xxx xx x x x x x x x xx o ooo oo o o o o o o o oo
1.0 xx 1.0 oo
1

1
●● x x x x xx x xxx x x x xx
o o o o oo o
o oo o o o oo

●●
● ● ● ●
x x x xx x x x
x xx x o o o oo o o x
o oo o
● ● ● ● xx xx x oo xx o
1 ●●
● ● ●● xx x x oo x o
● ● ●●
● ● ● ● x o
● 0.5 x o 0.5 x o
● ●● ● ● x x
●● ● x o

● xx xx
X1

X1


0.0 0.0
●●
X1

x x x x
0
x x x x
−0.5 x x x x −0.5 o o o o
● x x x x x o
x x x xx x o o o oo o
x x x o oo o
x xx
0

0
● xx x x x oo o x o
x x x x x o o o o o
● ● x x x o o o
● −1.0
xx xx xx xx xx x x −1.0
oo oo oo oo ooo o
● ● ● x xxx xx x x x xx x x o ooo oo o o o xx o o
x x xx xx o o oo oo


● o x x x x x o oo o
● ● x x x xx x
xxx o o o oo o
ooo
x x xxxx o o oooo
−1 ●
● ● x xx xx o oo oo
● ● ●
● ● ●●●● ● ● ● o x x o o o
● ●
●●●●
● −1.5 x x −1.5 o o
● ● ● ● ● ●
●● ●
●●

xx oo


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

−1 0 1 X2 X2
X2

Le package kernlab

— Il propose un choix plus large de noyaux.

> library(kernlab)
> [Link] <- ksvm(Y~.,data=donnees,kernel="polydot",
+ kpar=list(degree=2),C=0.001)
> plot([Link])

SVM classification plot

1.0


1.5 ● ●
● ●


● ● ●
● ●

1.0 ●● ● ●
● ● 0.5
● ●
●●
●● ●
●● ●
● ●
● ● ●
●● ● ●

● ●
●● ●
0.5 ●
● ●


0.0
X1

0.0

−0.5 ●
● ●●
● ● ● −0.5
● ● ● ●
● ●● ● ●
● ● ●

● ●●●
● ●●
−1.0 ● ● ●●


● ● ●
● ● ●
● ●● ●
● ●
● ●
● −1.0
−1.5

−1.5 −1.0 −0.5 0.0 0.5 1.0 1.5

X2

34
Troisième partie

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.
La méthode CHAID est proposée en annexe.

1 Arbres binaires
Notations

— On cherche à expliquer une variable Y par p variables explicatives X1 , . . . , Xp .

— Y peut admettre un nombre quelconque de modalités et les variables X1 , . . . , Xp 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 obervations (X1 , Y1 ), . . . , (Xn , Yn ) où Xi ∈ R2 et Yi ∈ {−1, 1}.

Approche par arbres


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

Définitions

Arbre binaire
Un arbre binaire de décision CART est
— un algorithme de moyennage local par partition (moyenne ou vote à la majorité sur les éléments de la
partition),
— dont la partition est construite par divisions successives au moyen d’hyperplans orthogonaux aux axes de Rp ,
dépendant des données (Xi , Yi ).

35
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 :

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

s1

36
s2

s1

s2

s1 s3

s2

s4

s1 s3

Représentation de l’arbre

s2 X1 < s 1 X1 ≥ s1

X 2 ≥ s2 X2 < s 2

X1 < s 3 X 1 ≥ s3
s4

X2 < s 4 X 2 ≥ s4

s1 s3

37
s2 X1 < s 1 X1 ≥ s1

X 2 ≥ s2 X2 < s 2

X1 < s 3 X 1 ≥ s3
s4

X2 < s 4 X 2 ≥ s4

s1 s3

Règle de classification
On effectue un vote à la majorité dans les nœuds terminaux de l’arbre.

Définitions

Définition
— Les éléments de la partition d’un arbre sont appelés les nœuds terminaux ou les feuilles de l’arbre.
— L’ensemble Rp constitue le nœud racine.
— Chaque division définit deux nœuds, les nœuds fils à gauche et à droite.

2 Choix des découpes

Questions
1. Comment choisir les découpes ?
2. Faut-il stopper les découpes ? Si oui, quand ?

— A chaque étape, on cherche un couple (j, s) qui split un noeud N en deux nœuds fils :

N1 (j, s) = {X ∈ N |Xj ≤ s} et N2 (j, s) = {X ∈ N |Xj > s}.

— La sélection du couple (j, s) s’effectue en optimisant un critère qui mesure l’(im)pureté ou l’hétérogénité des
deux nœuds fils.

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é :

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

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

38
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 .

Découpe en régression
A chaque étape, on choisit le couple (j, s) qui minimise
X X
(Yi − Ȳ1 )2 + (Yi − Ȳ2 )2
Xi ∈N1 (j,s) Xi ∈N2 (j,s)

1
P
où Ȳk = |Nk (j,s)| Xi ∈Nk (j,s) Yi , k = 1, 2.

Exemple
2.5

2.5
2.0

2.0
1.5

1.5
1.0

1.0
0.5

0.5
0.0

0.0

0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0
2.5

2.5
2.0

2.0
1.5

1.5
1.0

1.0
0.5

0.5
0.0

0.0

0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0
2.5

2.5
2.0

2.0

0,25 76,95 2,71


1.5

1.5
1.0

1.0

0.81
0.5

0.5
0.0

0.0

0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0

Sélection
On choisira le seuil de droite.

39
2.2 Cas de la classification supervisée
— On se place ici dans le cas binaire, Y dans {0, 1} (voir Annexe pour le cas multiclasse).
— Un nœud est pur si
— il contient beaucoup de 0 et peu de 1 (ou l’inverse) ;
— la proportion de 1 est proche de 1 (ou de 0).

Impureté de Gini
I(N ) = 2p(N )(1 − p(N ))
où p(N ) représente la proportion de 1 dans N .

Exemple

I(N ) = 0.4872
1.00 1.00
● ●
● ●
● ●
● ● ● ●

● ●

0.75 0.75

● ●
● ●
● ●

● ●
Y Y
0.50 ● ● 0 0.50 ● ● 0
1 1
● ●
● ●
● ● ● ●
● ●
● ● ● ●
● ●
● ●
0.25 0.25
● ● ● ●
● ●
● ●
● ●
● ●
● ●
● ●
0.00 0.00

0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00

1.00 1.00
● ●
● ●
● ●
● ● ● ●

● ●

0.75 0.75

● ●
● ●
● ●

● ●
Y Y
0.50 ● ● 0 0.50 ● ● 0
1 1
● ●
● ●
● ● ● ●
● ●
● ● ● ●
● ●
● ●
0.25 0.25
● ● ● ●
● ●
● ●
● ●
● ●
● ●
● ●
0.00 0.00

0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00

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


Gauche 0.287 0.137 0.281
Droite 0.488 0.437 0.031

Conclusion
On choisira la découpe de gauche.

3 Elagage
Questions

— Comment construire un "bon" arbre ?


— Construire l’arbre maximal ? (on découpe les nœuds jusqu’à ce qu’on ne puisse plus).
— Faut-il se donner un critère d’arrêt ?
— Faut-il construire un arbre grand et choisir un sous-arbre de ce dernier ?

40
1.00 ● ● ●
● ●
● ●
● ● ●
● ● ● ● ●
● ●
● ● ● ●


● ●● ● ● ● ●
● ● ●
● ● ● ● ●● ● ●
● ●
● ●
● ● ●●
● ● ● ●
● ● ● ● ●
● ●
● ● ● ● ●●
● ● ● ● ●
● ● ●● ●




0.75 ● ●
● ●


● ●

● ●


Y
0.50 ● 0
● 1

● ●

● ● ●
● ●
● ● ●

● ● ●
● ●

● ● ●
●●● ●

● ●
●● ● ●
0.25 ● ● ● ●
● ● ●

● ● ●

● ●
● ● ●
● ● ● ● ●
● ●
● ●
● ●● ●
● ● ●● ● ●
● ●
● ●
● ● ●
● ●
● ● ● ●
● ●● ● ●
● ● ●
● ● ●● ●
● ● ● ●●
● ●

● ● ● ● ●

0.00

0.00 0.25 0.50 0.75 1.00

yes X1 >= 0.64 no

X2 >= 0.78
X2 < 0.4

X2 >= 0.79 X2 < 0.017


X1 >= 0.75 X1 >= 0.22

X1 < 0.74 X1 >= 0.99 X1 < 0.095 X2 < 0.23


0 X1 < 0.98 X1 < 0.64 X1 >= 0.26

1 1 X2 >= 0.62 1 1 1 X1 < 0.22


X2 < 0.3 0 0 X1 >= 0.43 X1 >= 0.079 0 X1 >= 0.6

X2 >= 0.32 X2 < 0.42 X1 < 0.42 X1 < 0.038 X2 >= 0.18 X1 >= 0.49
X1 >= 0.66 X1 < 0.77 0 0 0 X1 >= 0.2

X1 < 0.66 1 X1 >= 0.95 1 1 1 X1 < 0.12 X2 < 0.43 X2 >= 0.59
0 0 0 X2 >= 0.42 X2 >= 0.8 X1 >= 0.028 X2 < 0.2 X1 < 0.2 X1 < 0.49

1 X1 < 0.8 1 X1 >= 0.35 1 X1 < 0.038 1 X1 >= 0.21 1 X2 < 0.31 1
0 X1 < 0.96 0 X2 < 0.92 0 0 X2 < 0.1 0 X2 >= 0.37 0 X2 < 0.62

1 1 X2 >= 0.92 1 X1 >= 0.54 1 1 X2 < 0.27 X1 >= 0.55 1


0 X1 >= 0.78 0 0 0 X2 >= 0.035 0 X1 < 0.16 X1 >= 0.57 0

1 1 X2 >= 0.23 1 1 1 1 1
0 X1 < 0.36 0 X1 < 0.058 0 X1 < 0.1 0 X1 < 0.55

X1 >= 0.37 1 X1 >= 0.091 1 X2 >= 0.56


0 0 0 0 0

1 1 1
0 0 X2 < 0.61

1
0

Un exemple en discrimination

Arbre optimal ?
Intuitivement, on a envie de faire à peu près 5 classes.

Arbre « maximal »

> library(rpart)
> library([Link])
> arbre1 <- rpart(Y~.,data=donnees[1:350,],cp=0.0001,minsplit=2)
> prp(arbre1)

Un arbre plus petit

> arbre2 <- rpart(Y~.,data=donnees[1:350,])


> prp(arbre2)

Comparaison des deux arbres

— On compare les performances des deux arbres en estimant leur probabilité de mauvais classement sur un
échantillon test :

41
yes X1 >= 0.64 no

X2 < 0.4 X2 >= 0.78

0 X2 >= 0.79 X1 >= 0.22 1

0 1 0 1

> prev <- [Link](arbre1=predict(arbre1,newdata=donnees[351:500,],type="class"),


+ arbre2=predict(arbre2,newdata=donnees[351:500,],type="class"),
+ obs=donnees[351:500,]$Y)
> prev %>% summarize_at(1:2,funs(mean(.!=obs))) %>% round(3)
## arbre1 arbre2
## 1 0.133 0.127

Conclusion
La performance n’augmente pas forcément avec la profondeur.

Sur-ajustement pour les arbres

Train error
Test error

OVERFITTING

Complexity (λ)

Remarque
La complexité d’un arbre est mesurée par sa taille ou profondeur.

Biais et variance
La profondeur régule le compromis biais/variance :
1. Peu de découpes (arbres peu profonds) =⇒ arbres stables =⇒ peu de variance... mais... beaucoup de biais.

42
2. Beaucoup de découpes (arbres profonds) =⇒ arbres instables =⇒ peu de biais... mais... beaucoup de variance
(surapprentissage).

Principe d’élagage [Breiman et al., 1984]


Plutôt que de choisir « quand couper » on raisonne en 3 temps :
1. On construit un arbre maximal (très profond) Tmax ;
2. On sélectionne une suite d’arbres emboités :

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

3. On sélectionne un arbre dans cette sous-suite.

— La construction de la suite de sous-arbres emboitées est détaillée en Annexe.


— Sur R, on obtient cette sous-suite à l’aide de la fonction printcp :
> arbre <- rpart(Y~.,data=donnees,cp=0.0001,minsplit=2)
> printcp(arbre)
##
## Classification tree:
## rpart(formula = Y ~ ., data = donnees, cp = 1e-04, minsplit = 2)
##
## Variables actually used in tree construction:
## [1] X1 X2
##
## Root node error: 204/500 = 0.408
##
## n= 500
##
## CP nsplit rel error xerror xstd
## 1 0.2941176 0 1.000000 1.00000 0.053870
## 2 0.1225490 1 0.705882 0.72059 0.049938
## 3 0.0931373 3 0.460784 0.51471 0.044646
## 4 0.0637255 4 0.367647 0.42647 0.041555
## 5 0.0122549 5 0.303922 0.35294 0.038483
## 6 0.0098039 7 0.279412 0.35294 0.038483
## 7 0.0049020 9 0.259804 0.35784 0.038704
## 8 0.0040107 25 0.181373 0.39216 0.040184
## 9 0.0036765 41 0.112745 0.39706 0.040386
## 10 0.0032680 49 0.083333 0.40196 0.040586
## 11 0.0024510 52 0.073529 0.41667 0.041174
## 12 0.0001000 82 0.000000 0.45098 0.042473

Sorties printcp

— Suite de 12 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.
— xstd : écart-type associé à l’erreur de validation croisée.

Choix de l’arbre final.


On choisira l’arbre qui a la plus petite erreur de prévision (calculée par validation croisée).

Tracé de l’arbre final


> cp_opt <- arbre$cptable %>% [Link]() %>%
+ filter(xerror==min(xerror)) %>% dplyr::select(CP) %>%
+ slice(1) %>% [Link]()
> cp_opt
## [1] 0.0122549
> arbre_final <- prune(arbre,cp = cp_opt)
> [Link](arbre_final)

43
1
0.59
100%
yes X1 >= 0.57 no

0 1
0.36 0.76
42% 58%
X2 < 0.37 X2 >= 0.8

1 0
0.56 0.34
25% 12%
X2 >= 0.77 X1 >= 0.22

0 0 1 0 1 1
0.07 0.09 0.80 0.10 0.81 0.88
17% 9% 16% 8% 4% 45%

Règle de classification et score par arbre

— L’arbre final T renvoie une partition de Rp en |T | nœuds terminaux N1 , . . . , N|T | .

— Règle de classification :
 P P
1 si i:Xi ∈N (x) 1Yi =1 ≥ i:Xi ∈N (x) 1Yi =0
ĝ(x) =
0 sinon,

où N (x) désigne le nœud terminal qui contient x.

— Score :
1 X
Ŝ(x) = P̂(Y = 1|X = x) = 1Yi =1 .
n
i:Xi ∈N (x)

Fonction predict

— La fonction predict ([Link]) permet d’estimer la classe ou le score :

> x_new <- [Link](X1=0.5,X2=0.85)


> predict(arbre_final,newdata=x_new)
## 0 1
## 1 0.9 0.1
> predict(arbre_final,newdata=x_new,type="class")
## 1
## 0
## Levels: 0 1

Bilan

— Méthode « simple » relativement facile à mettre en œuvre.

— Fonctionne en régression et en discrimination.

— Résultats interprétables (à condition que l’arbre ne soit pas trop profond).

— Un inconvénient : 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.

44
4 Annexe 1 : impureté, cas multiclasses
— 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


— 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 -1) dans N .

Impureté dans le cas binaire

0.6

0.4
indice
impurete

gini
information

0.2

0.0

0.00 0.25 0.50 0.75 1.00


p

45
Découpe en classification supervisée

— On rappelle que pour un nœud N donné et un couple (j, s), on note


N1 (j, s) = {X ∈ N |Xj ≤ s} et N2 (j, s) = {X ∈ N |Xj > s}.

Choix de (j, s)
Pour une mesure d’impureté I donnée, on choisira le couple (j, s) qui maximise le gain d’impureté :
∆(I) = P(N )I(N ) − (P(N1 )I(N1 (j, s)) + P(N2 )I(N2 (j, s)).

5 Annexe 2 : algorithme élagage


Construction de la suite de sous arbres

— Soit T un arbre à |T | nœuds terminaux N1 , . . . , N|T | .


— Soit R(N ) le risque (l’erreur) dans le nœud N :
— Régression :
1 X
R(N ) = (Yi − ȲN )2 .
|N |
i:Xi ∈N

— Classification binaire :
1 X
R(N ) = 1Yi 6=YN .
|N |
i:Xi ∈N

Définition
Soit α > 0, on pose
|T |
X
Cα (T ) = Nm R(Nm ) + α|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+∞ =arbre sans coupure.
— α est appelé paramètre de complexité et Cα (T ) le cout de l’arbre T .
Théorème 5.1 ([Breiman et al., 1984]). Il existe une sous-suite finie α0 = 0 < α1 < . . . < αM avec M < |Tmax |
et une suite associée d’arbres emboités
Tmax = T0 ⊃ T1 ⊃ . . . ⊃ TM
telles que ∀α ∈ [αm , αm+1 [
Tm = argmin Cα (T ).
T

T0 T1 TM

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

Conséquences
— On se ramène à une sous-suite finie d’arbres (emboités).
— Il reste à choisir un arbre (ou une valeur de α).

46
Exemple
> printcp(arbre)
Classification tree:
rpart(formula = Y ~ ., data = donnees, cp = 1e-04, minsplit = 2)
Variables actually used in tree construction:
[1] X1 X2
Root node error: 204/500 = 0.408
n= 500
CP nsplit rel error xerror xstd
1 0.2941176 0 1.000000 1.00000 0.053870
2 0.1225490 1 0.705882 0.71569 0.049838
3 0.0931373 3 0.460784 0.49020 0.043844
4 0.0637255 4 0.367647 0.43627 0.041928
5 0.0122549 5 0.303922 0.34314 0.038034
6 0.0098039 7 0.279412 0.34314 0.038034
7 0.0049020 9 0.259804 0.36275 0.038923
8 0.0040107 25 0.181373 0.34804 0.038260
9 0.0036765 41 0.112745 0.39216 0.040184
10 0.0032680 49 0.083333 0.40196 0.040586
11 0.0024510 52 0.073529 0.41176 0.040980
12 0.0001000 82 0.000000 0.43137 0.041742

> arbre1 <- prune(arbre,cp=0.005)


> prp(arbre)
> prp(arbre1)

yes X1 >= 0.57 no

X2 >= 0.8

X2 < 0.37

X2 >= 0.77 X2 < 0.017

X1 < 1 X1 >= 0.22

1 X2 < 0.4 X1 < 0.095 X2 < 0.27

X1 >= 0.75 X1 < 0.98 X2 < 0.91 X1 >= 0.26

yes X1 >= 0.57 no

X1 < 0.74 1 X1 >= 0.99 X2 >= 0.92 1 1 X2 >= 0.61

0 X2 < 0.9 X1 >= 0.68 0 X2 >= 0.85 0 X1 < 0.18

X2 < 0.37 X2 >= 0.8


1 X2 >= 0.9 1 X2 >= 0.62 X1 >= 0.3 1 X2 >= 0.18 X1 < 0.2

X1 >= 0.6 0 0 0 X2 >= 0.95 X2 < 0.94 X1 >= 0.015 X2 < 0.62

X1 < 0.59 1 X1 < 0.59 X2 < 0.95 1 X1 >= 0.079 1 1 X1 >= 0.54 X1 >= 0.55
0 X2 >= 0.77 X1 >= 0.22 X2 < 0.017
X2 < 0.3 X1 >= 0.66 X1 >= 0.66 0 X1 < 0.35 0 X2 >= 0.027 X2 < 0.19 0 X1 >= 0.2

X2 >= 0.32 1 X1 < 0.64 1 X2 < 0.42 1 1 1 1 X1 >= 0.54 X1 < 0.3 X2 < 0.43 1

X1 >= 0.66 0 0 X1 < 0.8 X1 >= 0.58 0 0 0 X2 < 0.07 0 X1 < 0.56 0 X1 < 0.55
0 X2 < 0.4 0 1 0 1

X1 < 0.66 1 1 X1 >= 0.95 1 1 X1 >= 0.16 X2 >= 0.23 1 X2 >= 0.78 1 1

0 0 0 X2 >= 0.74 0 X2 >= 0.42 0 0 0 X1 >= 0.28 X2 >= 0.41 0

X1 >= 0.68 X1 >= 0.99


1 X2 < 0.71 1 1 X2 >= 0.19 1 X1 >= 0.19 1 X2 >= 0.37

0 0 X1 < 0.96 0 0 X2 < 0.23 0 X2 < 0.79 0

0 1 0 1
1 1 X2 < 0.1 X1 >= 0.41 1 1 1

0 0 X2 < 0.22 0 X1 < 0.22 0 X1 < 0.15

X2 >= 0.26 1 1 1 1

0 X1 < 0.13 X1 < 0.45 0 0

1 1 1

0 X2 >= 0.074 0

Choix d’un arbre


Il reste à sélectionner un arbre dans la suite

Tmax = T0 ⊃ T1 ⊃ . . . ⊃ TM

Sélection d’un arbre

Choix d’un risque


La sélection de l’arbre final s’effectue en choisissant l’élément de la suite qui minimise le risque moyen
E[R(Y, Tm (X)]. Par exemple,
1. l’erreur quadratique E[(Y − Tm (X))2 ] en régression ;
2. la probabilité d’erreur P(Y 6= Tm (X)) en discrimination binaire.
Ce risque (inconnu) est estimé par validation croisée.

Choix de l’arbre final


L’approche consiste à
1. estimer le risque pour chaque αm .
2. choisir le αm qui minimise le risque estimé =⇒ Tαm .

47
Elagage/pruning - Algorithme

Algorithme
1. Calculer la suite α0 = 0 < α1 < . . . < αM et poser
√ √
β1 = 0, β2 = α1 α2 , β3 = α2 α3 , ..., βM +1 = ∞.

2. Séparer les données en K blocs G1 , . . . , Gk de taille k/n. Pour i = 1, . . . , k :


(a) Construire les arbres Tβ1 , . . . , TβM +1 sur l’ensemble des observations privé du ième bloc.
(b) En déduire pour tout j ∈ Gi et tout m ≤ M + 1, Ŷj (βm ) = Tβm (Xj ).
Pn
3. Calculer R(m) = n1 i=1 R(Yi , Ŷi (βm )) pour m = 1, . . . , M + 1.
4. Choisir αm? tel que βm? +1 = argminm≤M +1 R(m).

— Les estimations R(m) se trouvent dans la colonne xerror de la fonction printcp :


CP nsplit rel error xerror xstd
1 0.2941176 0 1.000000 1.00000 0.053870
2 0.1225490 1 0.705882 0.71569 0.049838
3 0.0931373 3 0.460784 0.49020 0.043844
4 0.0637255 4 0.367647 0.43627 0.041928
5 0.0122549 5 0.303922 0.34314 0.038034
6 0.0098039 7 0.279412 0.34314 0.038034
7 0.0049020 9 0.259804 0.36275 0.038923

— On peut représenter les erreurs en fonction des αm à l’aide de plotcp


> plotcp(arbre3)

size of tree

1 2 4 5 6 8 10 26 42 50 53 83
1.0


X−val Relative Error

0.8


0.6

● ●
0.4

● ●


● ● ●
0.2

Inf 0.19 0.11 0.028 0.0069 0.0038 0.0028

cp

On choisira l’arbre à 5 coupures.

Tracé de l’arbre final


> alpha_opt <- arbre$cptable[[Link](arbre$cptable[,"xerror"]),"CP"]
> arbre_final <- prune(arbre,cp=alpha_opt)
> prp(arbre_final)

yes X1 >= 0.57 no

X2 < 0.37 X2 >= 0.8

0 X2 >= 0.77 X1 >= 0.22 1

0 1 0 1

48
6 Annexe 3 : arbres Chaid
— CHAID : Chi2 Automatic Interaction Detection [Kass, 1980].

— 2 étapes χ2 dans le procédé de division d’un nœud :


— regrouper les modalités peu discriminantes de chaque variable explicative Xj ;

— choisir la variable à utiliser pour scinder le nœud.

χ2 d’indépendance : rappel

— Soient X et Y deux variables aléatoires à valeurs dans E et F . On souhaite tester au niveau α les hypothèses
H0 : "X et Y sont indépendantes" contre H1 : "X et Y ne sont pas indépendantes".
— On se donne (E1 , . . . , EI ) et (F1 , . . . , FJ ) deux partitions de E et F .
— On dispose de n mesures du couple (X, Y ) et on désigne par Nij l’effectif observé dans la classe Ei × Fj .

F1 ... Fj ... FJ Total


E1 N11 ... N1j ... N1J N1•
.. ..
. .
Ei Ni1 ... Nij ... NiJ Ni•
.. ..
. .
EI NI1 ... NIj ... NIJ NI•
Total N•1 ... N•j ... N•J n

Le test

Propriété
Sous H0 la statistique
I X
J NI• N•J
2
X
n − Nij
Xn = NI• N•J
i=1 j=1 n

converge en loi vers la loi χ2(I−1)(J−1) .

Conséquence
— Au niveau α, on rejettera l’hypothèse H0 si Xobs est supérieure au quantile d’ordre 1−α de la loi du χ2(I−1)(J−1) .
— Une forte valeur de Xobs (ou une faible valeur de la probabilité critique) signifiera un lien fort entre les deux
variables.

Chaid : le principe

— On suppose dans un premier temps que toutes les variables explicatives Xj , j = 1, . . . , p sont qualitatives à
Mj modalités.

Division d’un nœud


1. Regroupement des modalités peu discriminantes de chaque variable Xj ;
2. Choix de la variable Xj la plus discriminante
3. Le nœud est alors divisé en un nombre de nœuds fils égal au nombre de modalités créées à l’étape 1.

49
6.1 Regroupement des modalités
1. On se place dans un nœud N et on considère une variable Xj à Mj modalités ;
2. Les observations dans le nœud définissent la table de contingence suivante
M1 . . . Mj
1
..
.
K
3. ∀(Mi , M` ) ∈ {M1 , . . . , Mj } , on calcule la statistique du χ2 croisant Y et les modalités (Mi , M` ) =⇒
2

χ2 (Mi , M` ) et p(Mi , M` ) la probabilité critique associée.


Remarque
— 2 modalités discriminantes =⇒ dépendance forte dans le test avec Y =⇒ "Fort rejet" de H0 =⇒ χ2
élevé ou pc faible ;
— Regrouper les modalités peu discriminantes revient donc à regrouper celles qui ont un χ2 faible ou une
pc grande.

4. On choisit la paire de modalités qui minimise le χ2 :

(M̃i , M̃` ) = argmin χ2 (Mi , M` ) = argmax p(Mi , M` ).


(Mi ,M` )∈{M1 ,...,Mj }2 (Mi ,M` )∈{M1 ,...,Mj }2

5. Si p(M̃i , M̃` ) > α2 (α2 ∈]0, 1[ fixé par l’utilisateur) alors on regroupe les modalités M̃i et M̃` et on retourne
à l’étape 2 avec le tableau à Mj − 1 modalités
M1 . . . Mj − 1
1
..
.
K
Sinon, on stoppe les regroupements.

Exemple

— On considère la variable marstat :


> aa <- table(USvoteS$vote3,USvoteS$marstat)
> aa

married widowed divorced never married


Gore 246 57 82 111
Bush 315 44 48 60

— On calcule les probabilités critiques pour les 6 croisements :


> res <- matrix(0,nrow=4,ncol=4)
> rownames(res) <- levels(USvoteS$marstat)
> colnames(res) <- levels(USvoteS$marstat)
> for (i in 1:3)
+ for (j in (i+1):4)
+ res[i,j] <- [Link](aa[,c(i,j)])$[Link]
+
+
> res
married widowed divorced never married
married 0 0.0194 7.64e-05 1.41e-06
widowed 0 0.0000 3.06e-01 1.65e-01
divorced 0 0.0000 0.00e+00 7.42e-01
never married 0 0.0000 0.00e+00 0.00e+00

Exemple de regroupement
Les modalités divorced et never married sont regroupées (si α2 < 0.742).
— En effet
> ctrl <- chaid_control(minsplit = 20,alpha2=0.74)
> a1 <- chaid(vote3~marstat,data=USvoteS,control = ctrl)
> plot(a1)

50
1
marstat

married widowed
divorced, never married

Node 2 (n = 561) Node 3 (n = 101) Node 4 (n = 301)


1 1 1
Gore

Gore

Gore
0.8 0.8 0.8

0.6 0.6 0.6

0.4 0.4 0.4

0.2 0.2 0.2


Bush

Bush

Bush

0 0 0

> ctrl <- chaid_control(minsplit = 20,alpha2=0.75)


> a2 <- chaid(vote3~marstat,data=USvoteS,control = ctrl)
> plot(a2)

1
marstat

married widowed divorcednever married

Node 2 (n = 561) Node 3 (n = 101) Node 4 (n = 130) Node 5 (n = 171)


1 1 1 1
Gore

Gore

Gore

Gore

0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6

0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2


Bush

Bush

Bush

Bush

0 0 0 0

Variables continues et ordinales

— Variables ordinale : le traitement est identique. Seules les modalités contiguës peuvent être regroupées.

— Variables continues : traitées comme des variables ordinales. Penser à utiliser [Link] sur R.

6.2 Division d’un nœud

Un autre χ2 pour choisir la variable

— La phase regroupement effectuée, il faut choisir une variable parmi les p variables regroupées pour diviser le
nœud.

51
— Idée : faire un χ2 pour chaque variable :
(X1 , M1 ) . . . (X1 , M1j ) (X2 , M1 ) ... (X2 , M2j ) ...
1
..
.
K
=⇒ p probabilités critiques p(X1 ), . . . , p(Xp ) et
— Xj discriminante =⇒ rejet de H0 =⇒ p(Xj ) petite.
— On choisit la variable j qui possède la plus petite probabilité critique.

Correction de Bonferroni

— Tendance à favoriser les variables ayant subi le plus de regroupements (erreur de type 1).

— Pour rééquilibrer, les probabilités critiques sont multipliées par le coefficient de Bonferroni :

p0 (Xj ) = bj p(Xj )

où bj correspond au nombre de manières les regrouper les Mj modalités initiales de Xj en M̃j modalités
finales.

— Variable qualitative et ordinale :


M̃j −1
− i)Mj
 
X
i (M̃j Mj − 1
bj = (−1) bj = .
i!(M̃j − i)! M̃j − 1
i=0

— On choisira la variable j ? qui minimise p0 (Xj )...

— à condition que p0 (Xj ) soit plus petit qu’un certain seuil α4 fixé par l’utilisateur.

— Le nœud sera scindé en autant de groupes que Xj possède de modalités (après la phase de regroupement).

Critère d’arrêt

Un nœud ne sera pas divisé si :


— p0 (Xj ) > α4 pour tout j = 1, . . . , p.
— le nœud est pur ou quasiment pur.
— le nœud contient trop peu d’observations...

Remarque
Sur R, on pourra regarder la fonction [Link] :

chaid_control(alpha2 = 0.05, alpha3 = -1, alpha4 = 0.05,


minsplit = 20, minbucket = 7, minprob = 0.01,
stump = FALSE, maxheight = -1)

6.3 Choix des paramètres


— En plus des paramètres associés au critère d’arrêt, deux paramètres sont à calibrer pour construire l’arbre :
les niveaux α2 et α4 .

— Il en existe un troisième (α3 ) qui concerne le remise en cause des regroupements des modalités.

Choix de α4
Degrés d’exigence pour couper un nœud :
— petit : très exigeant =⇒ arbres peu profonds (beaucoup de biais et peu de variance) ;

52
— grand : peu exigeant =⇒ arbres profonds (beaucoup de variance et peu de biais).

Choix de α2
Degrés d’exigence pour regrouper des modalités :
— petit : peu exigeant =⇒ beaucoup de regroupements (on se rapproche des arbres binaires) ;
— grand : très exigeant =⇒ peu de regroupements.

Illustration α4
> ctrl <- chaid_control(minsplit = 20,alpha4=0.0005)
> a1 <- chaid(vote3~.,data=USvoteS,control=ctrl)
> plot(a1)

1
marstat

married widowed, divorced, never married

Node 2 (n = 560) Node 3 (n = 401)


1 1
Gore

Gore
0.8 0.8

0.6 0.6

0.4 0.4

0.2 0.2
Bush

Bush

0 0

> ctrl <- chaid_control(minsplit = 20,alpha4=0.25)


> a2 <- chaid(vote3~.,data=USvoteS,control=ctrl)
> plot(a2)

1
marstat

married widowed, divorced, never married


2 15
educr gender

<HS, HS,
College,
>HS Post Coll male
female
4 16 25
empstat empstat ager

yes, retired no yes, retired


18−24,
no 25−34, 35−44,
55−64,45−54
65+
5 12 17 27
empstat gender educr educr

yes, no
retired < HS
HS, >HS, College, Post Coll <HS,
>HS,
HSCollege, Post Coll
7 19 28 33
educr empstat educr marstat

<HS, HS, >HS,


Post Coll
College
male
female yes,
retired
no HS, >HS,<
College,
HSmarried,
Postwidowed,
Coll
never married
divorced
8 20 30 34
gender marstat empstat ager

male
female divorced,
married, widowed
never married 18−24,
yes,
25−34,
retired
no 35−44,
65+ 45−54, 55−64

Node
Node
3 (n Node
=
6 311)
(n Node
=9189)
(nNode
=
1023)
(n
Node
11
= Node
9)
(n13
= 6)
(n
Node
14= (n
3)
Node
18
= Node
19)
(n21
= Node
8)
(n
22=(n
Node
3)23
= 104)
Node
(n24
= 26)
(n
Node
26= (n
18)
Node
29
= 127)
(n
Node
31
= 17)
(n
Node
32
= 16)
(n
Node
35
= 26)
(n
Node
36
= 10)
(n37
= 32)
(n = 14)
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

Bush Gore

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

Illustration α2
> ctrl <- chaid_control(minsplit = 20,alpha2=0.005)
> a3 <- chaid(vote3~.,data=USvoteS,control=ctrl)
> plot(a3)

53
1
marstat

married
widowed, divorced, never married

2 5
educr gender

male female

7
<HS, HS,
College,
>HS Post Coll ager

18−24, 25−34, 35−44,


55−64,
45−54
65+

Node 3 (n = 311) Node 4 (n = 249) Node 6 (n = 159) Node 8 (n = 127) Node 9 (n = 115)
1 1 1 1 1

Gore

Gore

Gore

Gore

Gore
0.8 0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6 0.6

0.4 0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2 0.2


Bush

Bush

Bush

Bush

Bush
0 0 0 0 0

> ctrl <- chaid_control(minsplit = 20,alpha2=0.5)


> a4 <- chaid(vote3~.,data=USvoteS,control=ctrl)
> plot(a4)

1
marstat

marrieddivorced,
widowed never married
2 18 21
educr ager gender

< HSPost
HS, >HS College
Coll male female
4 14 22 25
ager empstat empstat educr

18−24,
35−44
45−54
25−34
55−64
65+ 18−24, 25−34,
55−64,
35−44,
65+ 45−54 > HS
<HS,
College
Post Coll
6 10 29
educr gender no,
yesretired no,
yesretired empstat

>HS, <HS,
College,
HS Post Coll
male
female yes,
retired
no

Node
Node
3 (n
Node
5
= (n
26)
Node
7
= (n
33)
Node
8
= (n
29)
Node
9
= (n
55)
Node
11
= 71)
Node
(n
12=Node
(n
19)
13= Node
(n
26)
15=(n
Node
52)
16
= 113)
(n
Node
17
= Node
(n
53)
19
=Node
83)
(n
20=Node
(n
23
9)=(n
Node
92)
24
= 107)
Node
(n
26= Node
(n
29)
27= Node
(n
37)
28= (n
58)
Node
30
= (n
45)
31
= 20)
(n = 4)
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
Gore
0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8 0.8
0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6 0.6
0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
Bush
0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

En pratique...

— L’influence de ces deux paramètres est bien entendu conjointe.

— Il n’est pas facile de les calibrer simultanément.

— Approche classique : évaluer les performances (erreur de classification AUC...) pour plusieurs valeurs de
(α2 , α4 ) sur un échantillon test ou par validation croisée.

Exemple

— On veut expliquer avec un arbre CHAID la variable chd par les autres variables du jeu de données SAheart.
> donnees <- SAheart
> donnees$chd <- [Link](donnees$chd)
> for (i in c(1:4,6:9)){donnees[,i] <- [Link](donnees[,i])}

— On va séparer l’échantillon en 2 et estimer l’erreur de classification sur une grille de valeur de α2 et α4 :


> alpha2 <- seq(0.01,0.35,by=0.05)
> alpha4 <- seq(0.01,0.35,by=0.05)
> [Link] <- [Link](alpha2,alpha4)
> names([Link]) <- c("alpha2","alpha4")
> [Link]$perf <- 0
> [Link](1234)
> perm <- sample(nrow(SAheart))
> dapp <- donnees[perm[1:300],]
> dtest <- donnees[-perm[1:300],]

— On estime l’erreur de classification sur les données test :


> for (i in 1:nrow([Link])){
> ctrl <- chaid_control(alpha2=[Link][i,1],alpha4=[Link][i,2])
> a <- chaid(chd~.,data=dapp,control=ctrl)
> prev <- predict(a,newdata = dtest)
> [Link]$perf[i] <- mean(prev!=dtest$chd)
}

— On récupère les valeurs de α2 et α4 qui minimisent l’erreur estimée :

54
> alpha_opt <- [Link][[Link]([Link]$perf),]
> alpha_opt
alpha2 alpha4 perf
1 0.01 0.01 0.2716049

— On peut tracer l’arbre sélectionné :


> ctrl <- chaid_control(alpha2=alpha_opt[1],alpha4=alpha_opt[2])
> arbre_final <- chaid(chd~.,data=donnees,control=ctrl)
> plot(arbre_final)

1
age

25, 26, 27, 28, 29, 30, 31,15,


32,16,
33,17,
34,18,
35,19,
51,
36,20,
52,
37,21,
53,
38,23,
54,
39,24
55,
40,56,
41,57,
42,58,
43,59,
44,60,
45,61,
46,62,
47,63,
48,64
49, 50

3 6
typea famhist

13, 20, 25, 26, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44,
69,45,
70,46,
71,47,
72,48,
73,49,
74,50,
75,51,
77,52,
7853, 54, Absent
55, 56, 57,
Present
58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68

Node 2 (n = 67) Node 4 (n = 209) Node 5 (n = 14) Node 7 (n = 82) Node 8 (n = 90)
1 1 1 1 1
0

0
0.8 0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6 0.6

0.4 0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2 0.2


1

1
0 0 0 0 0

Avec Caret

— On peut faire la même chose avec caret (en plus efficace) :


> grille <- [Link][,1:2]
> grille$alpha3 <- -1
> library(doMC)
> registerDoMC(cores = 3)
> bb <- train(donnees[,-10],donnees$chd,method="chaid",
tuneGrid=grille,trControl=ctrl1,metric="Accuracy")
> bb
CHi-squared Automated Interaction Detection

462 samples
9 predictor
2 classes: ’0’, ’1’

No pre-processing
Resampling: Repeated Train/Test Splits Estimated (1 reps, 75%)
Summary of sample sizes: 300
Resampling results across tuning parameters:

alpha2 alpha4 Accuracy Kappa


0.01 0.01 0.7283951 0.3847747
0.01 0.06 0.7283951 0.3847747
0.01 0.11 0.7283951 0.3847747
0.01 0.16 0.7283951 0.3847747
0.01 0.21 0.7283951 0.3847747
0.01 0.26 0.6851852 0.2528486
0.01 0.31 0.6851852 0.2528486
0.06 0.01 0.6851852 0.2528486
0.06 0.06 0.6851852 0.2528486
0.06 0.11 0.6851852 0.2528486
0.06 0.16 0.6851852 0.2528486
0.06 0.21 0.6851852 0.2528486
0.06 0.26 0.6728395 0.3284843
0.06 0.31 0.6728395 0.2302313
0.11 0.01 0.6419753 0.2394366
0.11 0.06 0.6419753 0.2839506
0.11 0.11 0.6419753 0.2839506
0.11 0.16 0.6419753 0.2839506
0.11 0.21 0.6419753 0.2839506
0.11 0.26 0.6296296 0.2646391
0.11 0.31 0.6419753 0.2839506
0.16 0.01 0.6419753 0.2394366

55
0.16 0.06 0.6419753 0.2839506
0.16 0.11 0.6419753 0.2839506
0.16 0.16 0.6419753 0.2839506
0.16 0.21 0.6419753 0.2839506
0.16 0.26 0.6296296 0.2646391
0.16 0.31 0.6419753 0.2839506
0.21 0.01 0.6419753 0.2394366
0.21 0.06 0.6419753 0.2394366
0.21 0.11 0.6419753 0.2394366
0.21 0.16 0.6419753 0.2394366
0.21 0.21 0.6419753 0.2394366
0.21 0.26 0.6419753 0.2394366
0.21 0.31 0.6419753 0.2394366
0.26 0.01 0.6419753 0.2394366
0.26 0.06 0.6419753 0.2394366
0.26 0.11 0.6419753 0.2394366
0.26 0.16 0.6419753 0.2394366
0.26 0.21 0.6419753 0.2394366
0.26 0.26 0.6419753 0.2394366
0.26 0.31 0.6419753 0.2394366
0.31 0.01 0.6419753 0.2394366
0.31 0.06 0.6419753 0.2394366
0.31 0.11 0.6419753 0.2394366
0.31 0.16 0.6419753 0.2394366
0.31 0.21 0.6419753 0.2394366
0.31 0.26 0.6419753 0.2394366
0.31 0.31 0.6419753 0.2394366

Tuning parameter ’alpha3’ was held constant at a value of -1


Accuracy was used to select the optimal model using the largest value.
The final values used for the model were alpha2 = 0.01, alpha3 = -1
and alpha4 = 0.21.

7 Bibliographie
Références

Biblio2

[Breiman et al., 1984] Breiman, L., Friedman, J., Olshen, R., and Stone, C. (1984). Classification and regression
trees. Wadsworth & Brooks.
[Cornillon and Matzner-Løber, 2011] Cornillon, P. and Matzner-Løber, E. (2011). Régression avec R. Springer.
[Devroye et al., 1996] Devroye, L., Györfi, L., and Lugosi, G. (1996). A Probabilistic Theory of Pattern Recognition.
Springer.
[Devroye and Krzyżak, 1989] Devroye, L. and Krzyżak, A. (1989). An equivalence theorem for l1 convergence of
the kernel regression estimate. Journal of statistical Planning Inference, 23 :71–82.
[Fahrmeir and Kaufmann, 1985] Fahrmeir, L. and Kaufmann, H. (1985). Consistency and asymptotic normality
of the maximum likelihood estimator in generalized linear models. The Annals of Statistics, 13 :342–368.
[Grob, 2003] Grob, J. (2003). Linear regression. Springer.
[Györfi et al., 2002] Györfi, L., Kohler, M., Krzyzak, A., and Harro, W. (2002). A Distribution-Free Theory of
Nonparametric Regression. Springer.
[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.
[Kass, 1980] Kass, G. (1980). An exploratory technique for investigating large quantities of categorical data.
Applied Statistics, 29(2) :119–127.
[Stone, 1977] Stone, C. J. (1977). Consistent nonparametric regression. Annals of Statistics, 5 :595–645.

56
Quatrième partie

Agrégation
Les approches que nous allons étudier sont basées sur l’agrégation :
1. construire un grand nombre de classifieurs "simples" g1 , . . . , gB
2. que l’on agrège
B
1 X
ĝ(x) = gk (x).
B
k=1

Questions
1. Intérêt d’agréger ?
2. Comment construire les gk pour que ĝ soit performant ?

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 ).

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
bootstrap) et les agréger.

Pourquoi agréger ?

— On se place dans le modèle de régression.


Y = m(X) + ε.
— On note
B
1 X
m̂B (x) = mk (x)
B
k=1

un estimateur de m obtenu en agrégeant B estimateurs m1 , . . . , mB .


— Rappels : m̂B (x) = m̂B (x; (X1 , Y1 ), . . . , (Xn , Yn )) et mk (x) = mk (x; (X1 , Y1 ), . . . , (Xn , Yn )) sont des variables
aléatoires.

— On peut mesurer l’intérêt d’agréger en comparant les performances de m̂B (x) à celles des mk (x), k = 1, . . . , B
(en comparant, par exemple, le biais et la variance de ces estimateurs).

57
1 2 3 4 5 6 7 8 9 10

3 4 6 10 3 9 10 7 7 1 m1
2 8 6 2 10 10 2 9 5 6 m2
2 9 4 4 7 7 2 3 6 7 m3
6 1 3 3 9 3 8 10 10 1 m4
3 7 10 3 2 8 6 9 10 2 m5
.. ..
. .
7 10 3 4 9 10 10 8 6 1 mB

Biais et variance

— Hypothèse : les variables aléatoires m1 , . . . , mB sont i.i.d.


— Biais :
E[m̂B (x)] = E[mk (x)].

Conclusion
Agréger ne modifie pas le biais.
— Variance :
1
V[m̂B (x)] = V[mk (x)].
B
Conclusion
Agréger tue la variance.

— Les conclusions précédentes sont vraies sous l’hypothèse que les variables aléatoires m1 , . . . , mB sont i.i.d.

— Les estimateurs m1 , . . . , mB étant construits sur le même échantillon, l’hypothèse d’indépendance n’est clai-
rement pas raisonnable !

Idée
Atténuer la dépendance entre les estimateurs mk , k = 1, . . . , B en introduisant de nouvelles sources d’aléa.

Idée : échantillons bootstrap

— Echantillon initial :
— Echantillons bootstrap :
— A la fin, on agrège :
B
1 X
m̂B (x) = mk (x).
B
k=1

Bagging

— Les mk ne vont pas être construits sur l’échantillon Dn = (X1 , Y1 ), . . . , (Xn , Yn ), mais sur des échantillons
bootstrap de Dn .

Bagging
Entrées :
— x ∈ Rd l’observation à prévoir, Dn l’échantillon
— un régresseur (arbre CART, 1 plus proche voisin...)
— B le nombre d’estimateurs que l’on agrège.
Pour k = 1, . . . , B :
1. Tirer un échantillon bootstrap dans Dn
2. Ajuster le régresseur sur cet échantillon bootstrap : mk (x)
PB
Sortie : L’estimateur m̂B (x) = B1 k=1 mk (x).

58
Tirage de l’échantillon bootstrap

— Les tirages bootstrap sont représentés par B variables aléatoires θk , k = 1, . . . , B.


— Les tirages bootstrap sont généralement effectués selon la même loi et de façon indépendante : θ1 , . . . , θB sont
i.i.d. de même loi que θ.
— 2 techniques sont généralement utilisées :

1. tirage de n observations avec remise ;

2. tirage de ` < n observation sans remise.

Conséquence
Les estimateurs agrégés contiennent 2 sources d’aléa (échantillon et tirage bootstrap) :

mk (x) = m(x, θk , Dn ).

Choix du nombre d’itérations

— Deux paramètres sont à choisir : le nombre d’itérations B et le régresseur.

— On a d’après la loi des grands nombres


B B
1 X 1 X
lim m̂B (x) = lim mk (x) = lim m(x, θk , Dn )
B→∞ B→∞ B B→∞ B
k=1 k=1
=Eθ [m(x, θ, Dn )] = m̄(x, Dn ) p.s|Dn .

— Lorsque B est grand, m̂B se "stabilise" vers l’estimateur bagging m̄(x, Dn ).

Conséquence importante
Le nombre d’itérations B n’est pas un paramètre à calibrer, il est conseillé de le prendre le plus grand possible en
fonction du temps de calcul.

Choix du régresseur

Propriété : biais et variance


On a E[m̂B (x)] = E[mk (x, θk , Dn )] et

1 − ρ(x)
V[m̂B (x)] = ρ(x)V[m(x, θk , Dn )] + V[m(x, θk , Dn )]
B
où ρ(x) = corr(m(x, θk , Dn ), m(x, θk0 , Dn ))) pour k 6= k 0 .

Conclusion
— Bagger ne modifie pas le biais.
— B grand =⇒ V[m̂B (x)] ≈ ρ(x)V[m̂k (x, θk (Dn )] =⇒ 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.

59
1.2 Forêts aléatoires

Rappels sur les arbres

s2 X1 < s 1 X1 ≥ s1

X 2 ≥ s2 X2 < s 2

X1 < s 3 X 1 ≥ s3
s4

X2 < s 4 X 2 ≥ s4

s1 s3

Paramètre à calibrer
Profondeur
— petite : biais %, variance &
— grande : biais &, variance %

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
T̂B (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


[Link]
et consulter la thèse de Robin Genuer [Genuer, 2010].

60
s2 X1 < s 1 X1 ≥ s1

X 2 ≥ s2 X2 < s 2

X1 < s 3 X 1 ≥ s3
s4

X2 < s 4 X 2 ≥ s4

s1 s3

Arbres pour forêt


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

Algorithme : randomforest
Entrées :
— x ∈ Rd l’observation à prévoir, Dn l’échantillon ;
— B nombre d’arbres ; nmax nombre max d’observations par nœud
— m ∈ {1, . . . , d} le nombre de variables candidates pour découper un nœud.
Pour k = 1, . . . , B :
1. Tirer un échantillon bootstrap dans Dn
2. Construire un arbre CART sur cet échantillon bootstrap, chaque coupure est sélectionnée en minimisant la
fonction de coût de CART sur un ensemble de m variables choisies au hasard parmi les d. On note T (., θk , Dn )
l’arbre construit.
PB
Sortie : l’estimateur TB (x) = B1 k=1 T (x, θk , Dn ).

Commentaires

— Si on est en discrimination (Y qualitative), l’étape d’agrégation consiste à faire voter les arbres à la majorité.

— Il y a deux sources d’aléa présentes dans θk : le tirage bootstrap et les m variables sélectionnées à chaque
étape de la construction de l’arbre.

— Méthode simple à mettre en oeuvre et déjà implémentée sur la plupart des logiciels statistiques (sur R, il suffit
de lancer la fonction randomForest du package randomForest).

— Estimateur connu pour fournir des estimations précises sur des données complexes (beaucoup de variables,
données manquantes...).

— Estimateur peu sensible au choix de ses paramètres (B, nmax , m...)

Choix des paramètres

— B : réglé... le plus grand possible.

61
Intérêt du bagging (rappel)
Diminuer la variance des estimateurs qu’on agrège :
1 − ρ(x)
V[T̂B (x)] = ρ(x)V[T (x, θk , Dn )] + V[T (x, θk , Dn )]
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, nmax = 5 en régression et 1 en classification.

Choix de m

— 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.
— m&
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 m %.

Conclusion
— Il est recommandé de comparer les performances de la forêt pour plusieurs valeurs de m.

— Par défaut m = d/3 en régression et d en classification.

Application sur les données spam


> library(randomForest)
> foret1 <- randomForest(type~.,data=spam)
> foret1
##
## Call:
## randomForest(formula = type ~ ., data = spam)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 7
##
## OOB estimate of error rate: 4.52%
## Confusion matrix:
## nonspam spam [Link]
## nonspam 2708 80 0.02869440
## spam 128 1685 0.07060121

Mesure de performance

— Comme pour les autres classifieurs et régresseurs il convient de définir des critères qui permettent de mesurer
la performance des forêts aléatoires.
— Exemples :
— Erreur de prédiction : E[(Y − T̂B (X))2 ] en régression ;
— Probabilité d’erreur : P(Y 6= T̂B (X)) en classification.
— Comme pour les autres méthodes, ces critères peuvent être évalués par apprentissage/validation ou validation
croisée.
— La phase bootstrap des algorithme bagging permet de définir une nouvelle méthode d’estimation de ces critères :
méthode OOB (Out Of Bag).

62
3 4 6 10 3 9 10 7 7 1 m1
2 8 6 2 10 10 2 9 5 6 m2
2 9 4 4 7 7 2 3 6 7 m3
6 1 3 3 9 3 8 10 10 1 m4
3 7 10 3 2 8 6 9 10 2 m5
7 10 3 4 9 10 10 8 6 1 m6

Erreur Ouf Of Bag

— Pour chaque observation (Xi , Yi ) de Dn , on désigne par IB l’ensemble des arbres de la forêt qui ne contiennent
pas cette observation dans leur échantillon bootstrap.
— La prévision de Y au point Xi se fait selon
1 X
Ŷi = T (Xi , θk , Dn ).
|IB |
k∈IB

Estimateurs Our Of Bag


1
Pn
— L’erreur de prédiction est estimée par (Ŷi − Yi )2 .
n
Pi=1
n
— La probabilité d’erreur est estimée par n1 i=1 1Ŷi 6=Yi .

Exemple

— Les échantillons 2, 3 et 5 ne contiennent pas la première observation, donc


1
Ŷ1 = (m2 (X1 ) + m3 (X1 ) + m5 (X1 )).
3

— On fait de même pour toutes les observations =⇒ Ŷ2 , . . . , Ŷn .


— On estime l’erreur selon
n
1X
(Ŷi − Yi )2 .
n i=1

Exemple

— On construit la forêt avec m = 1 :


> foret2 <- randomForest(type~.,data=spam,mtry=1)
> foret2
##
## Call:
## randomForest(formula = type ~ ., data = spam, mtry = 1)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 1
##
## OOB estimate of error rate: 8.06%
## Confusion matrix:
## nonspam spam [Link]
## nonspam 2725 63 0.02259684
## spam 308 1505 0.16988417

Remarque
L’erreur OOB est de 8.06%, elle est de 4.52% lorsque m = 7.

63
Importance des variables

— Un des reproches souvent fait aux forêts est l’aspect boîte noire et manque d’interprétabilité par rapport aux
modèles paramétriques tels que le modèle logistique.

— Il existe un indicateur qui permet de mesurer l’importance des variables présentes dans le modèle.

— Comme l’erreur OOB, ce critère est basé sur le fait que toutes les observations ne sont pas utilisées pour
construire les arbres de la forêt.

— Soit OOBk l’échantillon Out Of Bag associé au k eme arbre : il contient les observations qui ne sont pas dans
le k eme échantillon bootstrap.
— Soit EOOBk l’erreur de prédiction de l’arbre k mesurée sur cet échantillon :
1 X
EOOBk = (T (Xi , θk , Dn ) − Yi )2 .
|OOBk |
i∈OOBk

— Soit OOBkj l’échantillon OOBk dans lequel on a perturbé aléatoirement les valeurs de la variable j et EOOB j
k
l’erreur de prédiction de l’arbre k mesurée sur cet échantillon :
j 1 X
EOOB k
= j
(T (Xij , θk , Dn ) − Yi )2 ,
|OOBk | j
i∈OOBk

Définition
L’importance de la j eme variable est définie par
B
1 X j
Imp(Xj ) = (EOOBk − EOOBk ).
B
k=1

Exemple

— L’importance s’obtient facilement avec le package randomForest


> foret <- randomForest(type~.,data=spam,importance=TRUE)
> Imp <- importance(foret,type=1) %>% [Link]() %>%
+ mutate(variable=names(spam)[-58]) %>% arrange(desc(MeanDecreaseAccuracy))
> head(Imp)
## MeanDecreaseAccuracy variable
## 1 52.88382 charExclamation
## 2 48.69346 remove
## 3 42.21059 capitalAve
## 4 38.86023 charDollar
## 5 38.85976 hp
## 6 37.77500 capitalLong

> ggplot(Imp[1:15,]) + aes(x=reorder(variable,MeanDecreaseAccuracy),


+ y=MeanDecreaseAccuracy)+
+ geom_bar(stat="identity")+coord_flip()+xlab("")+theme_classic()

charExclamation

capitalAve

remove

hp

charDollar

edu

free

capitalLong

capitalTotal

george

your

our

num1999

re

you

0 10 20 30 40 50
MeanDecreaseAccuracy

64
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].

— Il a ensuite été montré que cet algorithme peut-être vu comme un cas particulier d’algorithmes de descente
de gradient =⇒ gradient boosting.

— C’est cette famille d’algorithmes que nous présentons dans cette partie.

2.1 Algorithmes de gradient boosting


— (X, Y ) couple aléatoire à valeurs dans Rd × Y. Etant donnée G une famille de règles, on se pose la question
de trouver la meilleure règle dans G.
— Choisir la règle qui minimise une fonction de perte, par exemple

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

Problème : la fonction de perte n’est pas calculable.


— Idée : choisir la règle qui minimise la version empirique de la fonction de perte :
n
1X
Rn (g) = `(Yi , g(Xi )).
n i=1

— On considère G l’ensemble des arbres binaires et on veut trouver la meilleure combinaison linéaire d’abres
binaires.

Un premier problème
PM
Trouver g(x) = m=1 αm hm (x) ∈ G qui minimise
n
1X
Rn (g) = `(Yi , g(Xi )).
n i=1

sans solution...
— Pas de solution explicite.
— Nécessité de trouver un algorithme pour approcher la solution.

Descente de gradient

— Idée : utiliser un algorithme de descente de gradient de type Newton-Raphson.


— Algorithme récursif :
— itération m : gm
— itération m + 1 : on ajoute à gm un nouvel arbre binaire hm+1 tel que le risque Rn (gm + λhm+1 ) diminue
le plus fortement (λ ∈ R+ petit).
— Approche classique : utiliser l’opposé du gradient fm+1 (x) = −∇Rn (gm )(x)

gm+1 (x) = gm (x) + λfm+1 (x).

Une restriction
Chaque élément de la suite doit être une combinaison d’arbres et fm+1 n’est pas (forcément) un arbre.

65
— Pour trouver l’arbre le plus proche du gradient fm+1 , on cherche un arbre h qui minimise
n
1X
(fm+1 (Xi ) − h(Xi ))2 .
n i=1

— Si on désigne par

Ui = fm+1 (Xi ) = −∇Rn (gm )(Xi ) = − `(yi , g(xi )) ,
∂g(xi ) g(xi )=gm−1 (xi )

la solution hm+1 s’obtient en ajustant un arbre sur l’échantillon (X1 , U1 ) . . . , (Xn , Un ).

Mise à jour de l’algorithme


gm+1 (x) = gm (x) + λhm+1 (x).

Gradient Boost Algorithm [Friedman, 2001] : FGD


Entrées :

— dn = (x1 , y1 ), . . . , (xn , yn ) l’échantillon, λ un paramètre de régularisation tel que 0 < λ ≤ 1.


— M le nombre d’itérations.
— paramètres de l’arbre (nombre de coupures...)
1
Pn
1. Initialisation : g0 (.) = argminc n i=1 `(yi , c)
2. Pour m = 1 à M :

a) Calculer l’opposé du gradient − ∂g(x i)
`(yi , g(xi )) et l’évaluer aux points gm−1 (xi ) :


Ui = − `(yi , g(xi )) , i = 1, . . . , n.
∂g(xi ) g(xi )=gm−1 (xi )

b) Ajuster un arbre sur l’échantillon (x1 , U1 ), . . . , (xn , Un ), on note hm l’arbre ainsi défini.
c) Mise à jour : gm (x) = gm−1 (x) + λhm (x).
3. Sortie : la suite de règles (gm (x))m .

Commentaires

— Il est facile de voir que chaque règle gm s’écrit


m
X
gm (x) = g0 + λ hk (x)
k=1

où les hk sont des arbres binaires. =⇒ ces règles sont donc bien des combinaisons d’arbres.

— Plusieurs paramètres sont à choisir ou à calibrer par l’utilisateur :


— fonction de perte ` ;
— nombre d’itérations M ;
— paramètre de régularisation λ ;
— paramètres de l’arbre (nombre de coupures notamment).

66
2.2 Choix des paramètres

Choix de la fonction de perte

— La fonction de perte `(y, g(x)) doit remplir 2 conditions :


1. mesurer une erreur : prendre des valeurs élevées lorsque y est loin de g(x) et faibles dans le cas inverse.
2. être régulière : notamment dérivable et convexe en son second argument.

Exemple
— Classification binaire avec Y dans {−1, 1} :
1. `(y, g(x)) = exp(−yg(x)) =⇒ adaboost ;
2. `(y, g(x)) = log(1 + exp(−2yg(x))) =⇒ logitboost ;
— Régression avec Y dans R : `(y, g(x)) = 0.5(y − g(x))2 =⇒ L2 -boosting.

Remarque

— Pour le L2 -boosting, les Ui de l’étape 2.a de l’algorithme FGD s’écrivent



Ui = − `(yi , g(xi )) = yi − gm−1 (xi ).
∂g(xi ) g(xi )=gm−1 (xi )

— Ces quantités correspondent aux résidus du régresseur à l’étape m − 1.

Interprétation
— L’estimateur à l’étape m est construit en faisant une régression sur les résidus correpondants à l’estimateur
à l’étape m − 1.
— On "corrige" gm−1 en cherchant à expliquer "l’information restante" qui est contenue dans les résidus.

— L’algorithme L2 -boosting (simplifié) peut alors s’écrire.

L2 -boosting - autre version


1. Initialisation g0 .
2. Pour m = 1 à M :
a) Calculer les résidus Ui = yi − gm−1 (xi ), i = 1, . . . , n.
b) Ajuster la règle faible sur (x1 , U1 ), . . . , (xn , Un ) =⇒ hm .
c) Mise à jour : gm (x) = gm−1 (x) + λhm (x).

Remarque importante [Bühlmann and Yu, 2003, Bühlmann and Hothorn, 2007, Cornillon et al., 2014]
Il a été montré (sous certaines hypothèses) qu’à chaque itération :
— Le biais diminue : B(gm ) ≤ B(gm−1 ).
— La variance augmente V (gm ) ≥ V (gm−1 ).
— D’où l’importance d’utiliser des règles faibles : beaucoup de biais et peu de variance (des arbres avec peu de
nœuds terminaux par exemple).

67
Illustration

On considère l’échantillon suivant



● ● ●
● ●

● ● ●
● ●
1.0 ● ●
●●●

● ●
● ● ● ● ● ● ● ●
● ● ●
● ●
● ● ●
● ● ●●
● ● ● ● ●
● ●
● ●
● ● ● ●

● ● ● ●
● ● ●
● ●
● ●
0.5 ● ● ●
● ● ●

● ●● ●
● ● ●
● ● ● ●

● ●

● ● ●

● ●
● ●
●● ● ●

0.0 ●
● ● ● ●
● ●
● ●

● ● ●
● ● ● ●

● ●

●●
● ● ●●

● ●

● ● ● ●

−0.5 ● ● ●
● ●
● ● ● ●●

●● ● ●


●● ● ●● ● ●●
● ● ●
● ●● ●
●●
● ● ●●
● ● ●
●● ●
● ● ● ●
● ●
−1.0 ●
● ●
● ●
● ●●
● ●
● ●


−1.5

−4 0 4

On ajuste un premier arbre simple :



● ● ●
● ●

● ● ●
● ●
1.0 ● ●
●●●

● ●
● ● ● ● ● ● ● ●
● ● ●
● ●
● ● ●
● ● ●●
● ● ● ● ●
● ●
● ●
● ● ● ●

● ● ● ●
● ● ●
● ●
● ●
0.5 ● ● ●
● ● ●

● ●● ●
● ● ●
● ● ● ●

● ●

● ● ●

● ●
● ●
●● ● ●

0.0 ●
● ● ● ●
● ●
● ●

● ● ●
● ● ● ●

● ●

●●
● ● ●●

● ●

● ● ● ●

−0.5 ● ● ●
● ●
● ● ● ●●

●● ● ●


●● ● ●● ● ●●
● ● ●
● ●● ●
●●
● ● ●●
● ● ●
●● ●
● ● ● ●
● ●
−1.0 ●
● ●
● ●
● ●●
● ●
● ●


−1.5

−4 0 4

On calcule les résidus :



1.0




●● ● ●

● ●
● ●

● ● ●
● ●

0.5 ● ●



● ●● ● ●● ●
● ● ●
● ● ●
● ●

● ●● ●
● ●
● ● ● ●
● ● ● ● ●
● ● ●
● ● ● ●

● ●
●● ● ● ●
● ● ● ● ●
● ● ●

● ●
● ● ●
● ● ● ● ●
● ●
0.0 ● ●
● ● ●


● ● ● ● ●
● ●
● ● ●
● ● ● ●
● ●
● ● ● ● ● ●
● ● ●
● ●
● ●● ● ● ●
● ●
● ●
● ●
●● ●● ●
● ● ●
● ● ●
● ● ● ●
● ●● ● ●
● ● ● ●
● ● ● ● ●
● ● ●
● ● ●●
● ●●

−0.5 ● ●

● ●

● ● ●

● ●





● ●

−4 0 4

On ajuste un nouvel arbre sur les résidus :

68


1.0




●● ● ●

● ●
● ●

● ● ●
● ●

0.5 ● ●



● ●● ● ●● ●
● ● ●
● ● ●
● ●

● ●● ●
● ●
● ● ● ●
● ● ● ● ●
● ● ●
● ● ● ●

● ●
●● ● ● ●
● ● ● ● ●
● ● ●

● ●
● ● ●
● ● ● ● ●
● ●
0.0 ● ●
● ● ●


● ● ● ● ●
● ●
● ● ●
● ● ● ●
● ●
● ● ● ● ● ●
● ● ●
● ●
● ●● ● ● ●
● ●
● ●
● ●
●● ●● ●
● ● ●
● ● ●
● ● ● ●
● ●● ● ●
● ● ● ●
● ● ● ● ●
● ● ●
● ● ●●
● ●●

−0.5 ● ●

● ●

● ● ●

● ●





● ●

−4 0 4

On obtient ainsi 2 arbres :

1.0

0.5

0.0

−0.5

−1.0

−4 0 4

Que l’on ajoute pour déduire un nouvel estimateur (2 itérations)...

1.0

0.5

0.0

−0.5

−1.0

−4 0 4

— Pour 1,000 et 500,000 itérations, on obtient :

69
M=1000 M=500000

● ●

● ●
● ●
● ● ● ● ● ●
● ● ● ●

● ● ● ● ● ●
● ● ● ● ● ●
1.0 ● ●●


●● ● ●●


●●
● ●● ● ●●● ● ● ●● ● ●●● ●
● ● ● ●
● ●● ● ●●
●● ● ●● ●
● ● ● ● ● ●
● ●
● ● ● ● ● ● ● ●
● ● ● ● ● ●
● ● ● ●
● ● ● ●
● ●● ● ●●
● ● ● ● ● ●
● ●
● ● ● ● ● ●
● ● ● ●
● ● ● ●
0.5 ● ● ● ● ● ●
● ●● ● ●●
● ● ● ●
● ● ● ● ● ●
● ●● ● ●●
● ● ●● ● ● ●●
● ●
● ● ● ●
● ● ● ●
●● ● ●● ●
● ● ● ●
● ● ● ●
● ● ● ●
●● ● ●● ●
0.0 ●
●● ● ● ●
●● ● ●
● ● ● ●
● ● ● ●
● ●
● ● ● ● ● ●
● ●● ● ● ●● ●
● ●
● ● ● ● ● ● ● ● ● ●

● ●● ●
● ●●
● ●
● ●
● ● ● ●
● ● ● ● ● ● ● ●

−0.5 ● ● ● ● ● ●
● ● ● ●
● ● ● ● ● ● ● ●
● ●
● ●
●● ● ● ●● ● ●
● ●
● ● ● ●
●●●
● ● ●● ●●●
● ● ●●
● ● ●● ● ● ●●
●● ● ●● ●
●● ●●
●● ●

● ●● ●● ●

● ●●

● ● ●
● ●
●●● ● ●●● ●
● ● ● ●
−1.0 ● ●
● ● ● ●
● ● ● ●

● ● ● ●
● ●
● ● ● ●
● ● ● ●

● ●

● ●
−1.5

−4 0 4 −4 0 4

Remarque importante
L’algorithme sur-ajuste si le nombre d’itérations est (trop) grand.

Choix de m et λ

— Le choix du nombre d’itérations est crucial pour les estimateurs boosting.


— Si m est trop grand on sur-ajuste (estimateurs avec peu de biais mais beaucoup de variance) et réciproquement
si m est trop petit.
— Le paramètre de régularisation λ représente le pas de la descente de gradient.
— Ce paramètre est lié à m : un λ grand nécessitera peu d’itérations et réciproquement.

En pratique
— On considère 2 ou 3 valeurs (petites) pour λ (0.1,0.01) ;
— Pour chaque λ, on choisit le meilleur m en utilisant des techniques de type validation croisée.
— L’algorithme étant construit pour une fonction de perte ` donnée, il est d’usage d’utiliser la même fonction
de perte pour sélectionner m.
— On va donc chercher le nombre d’itérations qui minimise la perte :
m̂ = argmin E[`(Y, gm (X))].
m≤M

— L’espérance ci-dessus étant inconnue en pratique, elle est approchée par des algorithmes de validation croisée.

Le coin R

— La fonction gbm du package gbm [Ridgeway, 2006] permet de faire du gradient boosting. Elle admet notam-
ment comme paramètres :
1. fonction de perte (distribution)
2. nombre d’itérations maximal M ([Link])
3. nombre de nœuds terminaux des arbres plus 1 ([Link])
4. paramètre de régularisation λ (shrinkage)
5. paramètres pour la validation croisée ([Link] pour la validation hold out ou [Link] pour la
validation croisée).

— La fonction [Link] permet de sélectionner le nombre d’itérations.

70
Exemple
> library(gbm)
> [Link](1234)
> spam1 <- spam
> spam1$type <- [Link](spam1$type)-1
> ada <- gbm(type~.,data=spam1,distribution="adaboost",[Link]=5,
+ [Link]=1000,shrinkage=0.05)
> mopt <- [Link](ada)
> mopt
## [1] 440

0.9
0.8
AdaBoost exponential bound

0.7
0.6
0.5
0.4
0.3

0 200 400 600 800 1000

Iteration

Conclusion

— Les algorithmes randomforest et boosting agrègent des arbres :


m
X
ĝm (x) = αk hk (x).
k=1

— Pour être efficace les arbres hk doivent être des règles faibles (weaklearner), donc des arbres peu performants :
— randomforest : arbres très profonds avec beaucoup de variance et peu de biais ;
— boosting : arbres peu profonds avec peu de variance et beaucoup de biais.

Résumé
— Agrégation RF : réduction de variance ;
— Agrégation boosting : réduction de biais.

3 Bibliographie
Références

Biblio4

[Aronszajn, 1950] Aronszajn, N. (1950). Theory of reproducing kernels. Transactions of the American Mathematical
Society, 68 :337–404.
[Breiman, 1996] Breiman, L. (1996). Bagging predictors. Machine Learning, 26(2) :123–140.
[Bühlmann and Hothorn, 2007] Bühlmann, P. and Hothorn, T. (2007). Boosting algorithms : regularization, pre-
diction and model fitting. Statistical Science, 22 :477–505.
[Bühlmann and Yu, 2003] Bühlmann, P. and Yu, B. (2003). Boosting with the l2 loss : regression and classification.
Journal of American Statistical Association, 98 :324–339.
[Cornillon et al., 2014] Cornillon, P., Hengartner, N., and Matzner-Lø ber, E. (2014). Recursive bias estimation
for multivariate regression smoothers. ESAIM : Probability and Statistics, 18(483-502).
[Devroye and Krzyżak, 1989] Devroye, L. and Krzyżak, A. (1989). An equivalence theorem for l1 convergence of
the kernel regression estimate. Journal of statistical Planning Inference, 23 :71–82.

71
[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.
[Freund and Schapire, 1997] Freund, Y. and Schapire, R. (1997). A decision-theoritic generalization of online
learning and an application to boosting. Journal of Computer and System Sciences, 55 :119–139.
[Freund and Schapire, 1999] Freund, Y. and Schapire, R. (1999). A short introduction to boosting. Journal of
Japanese Society for Artificial Intelligence, 14(5) :771–780.
[Friedman, 2001] Friedman, J. H. (2001). Greedy function approximation : A gradient boosting machine. Annals
of Statistics, 29 :1189–1232.
[Fromont, 2015] Fromont, M. (2015). Apprentissage statistique. Université Rennes 2, diapos de cours.
[Genuer, 2010] Genuer, R. (2010). Forêts aléatoires : aspects théoriques, sélection de variables et applications. PhD
thesis, Université Paris XI.
[Györfi et al., 2002] Györfi, L., Kohler, M., Krzyzak, A., and Harro, W. (2002). A Distribution-Free Theory of
Nonparametric Regression. Springer.
[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.
[Ridgeway, 2006] Ridgeway, G. (2006). Generalized boosted models : A guide to the gbm package.
[Stone, 1977] Stone, C. J. (1977). Consistent nonparametric regression. Annals of Statistics, 5 :595–645.
[Vapnik, 2000] Vapnik, V. (2000). The Nature of Statistical Learning Theory. Springer, second edition.
[Vert, 2014] Vert, J. (2014). Support vector machines and applications in computational biology. disponible à l’url
[Link]

72
Cinquième partie

Réseaux de neurones
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...

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

73
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 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 ).

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).

74
Représentation graphique

Entrées Poids Neuronne 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.

Summary

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

75
> 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 Neuronne 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 ))

76
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

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

77
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

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...

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")

78
> 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
## ___________________________________________________________________________

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.

79
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.

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.

80

●●

● ●●●● ●●●
●● ● ●
●● ●
●●
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

> 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.

> 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]

81
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.

Dropout

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

C1 C2 C3 S C1 C2 C3 S

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

82
> 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
## 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.

6 Bibliographie
Références

83
Biblio5

[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.

84

Vous aimerez peut-être aussi