Chapitre II

Régression linéaire multiple

Licence 3 MIASHS - Université de Bordeaux

Marie Chavent

Chapitre 2 Régression linéaire multiple 1/40

Un exemple

On cherche à modéliser la relation entre poids des bébés à naissance et l’âge, le

poids et le statut tabagique de la mère durant la grossesse. On pose :
- y = poids de naissance en grammes (bwt),
- x1 = âge de la mère (age),
- x2 = poids de la mère en kilos (weight),
- x3 = statut tabagique de la mère pendant la grossesse (smoke) codée
1=oui et 0=non.
On suppose que cette relation est linéaire de la forme :

y = β0 + β1 x1 + β2 x2 + β3 x3

On veut estimer cette relation avec un modèle de régression multiple.

On utilise un échantillon de n = 1174 naissances pour lesquelles le poids
du bébé, l’âge, le poids et le statut tabagique de la mère, ont été mesurés.

Chapitre 2 Régression linéaire multiple 2/40

 load("poids.RData")
print(data[1:5,c("bwt","age","weight","smoke")],digits=4)

## bwt age weight smoke

## 1 3402 27 45.36 0
## 2 3203 33 61.23 0
## 3 3629 28 52.16 1
## 4 3062 23 56.70 1
## 5 3856 25 42.18 0

pairs(data) #diagrammes de dispersion

150 250 350 15 25 35 45 40 80

3500

bwt

modele <- lm(bwt~ age+weight+smoke,data=data)

1500
350

modele$coefficients
250

gestation
150

0.4 0.8

parity
0.0
15 25 35 45

age ## (Intercept) age weight smoke

## 3050.56238 -0.91802 7.90266 -254.25425
170

height
140
80

weight
40

0.8

smoke
0.0 0.4

1500 3500 0.0 0.4 0.8 140 170 0.0 0.4 0.8

Chapitre 2 Régression linéaire multiple 3/40

1. Le modèle

On cherche à modéliser la relation entre plus de 2 variables quantitatives.

Un modèle de régression linéaire multiple est de la forme suivante :
p
X
y = β0 + βj xj + ε (1)
j=1

où :
- y est la variable à expliquer (à valeurs dans R) ;
- x1 , . . . , xp sont les variables explicatives (à valeurs dans R) ;
- ε est le terme d’erreur aléatoire du modèle ;
- β0 , β1 , . . . , βp sont les paramètres à estimer.

Commentaires :
- La désignation “multiple” fait référence au fait qu’il y a plusieurs variables
explicatives xj pour expliquer y .
- La désignation “linéaire” correspond au fait que le modèle (1) est linéaire.

Chapitre 2 Régression linéaire multiple 4/40

Pour n observations, on peut écrire le modèle de régression linéaire multiple
sous la forme :
p
X
yi = β0 + βj xij + εi pour i = 1, . . . , n. (2)
j=1

Dans ce chapitre, on suppose que :

- εi est une variable aléatoire, non observée,
- xij est observé et non aléatoire,
- yi est observé et aléatoire.

On fait les trois hypothèses additionnelles suivantes :

(A1) E[εi ] = 0, ∀i = 1, . . . , n,
ou de manière équivalente :
E[yi ] = β0 + pj=1 βj xij , ∀i = 1, . . . , n.
P

Commentaire sur l’hypothèse (A1) : elle indique que les erreurs sont centrées

Chapitre 2 Régression linéaire multiple 5/40

(A2) V(εi ) = σ 2 , ∀i = 1, . . . , n,
ou de manière équivalente :
V(yi ) = σ 2 , ∀i = 1, . . . , n.
Commentaires sur l’hypothèse (A2) :
- On parle d’hypothèse d’homoscédasticité (' homogénéité des variances).
- Cette variance σ 2 est un paramètre du modèle qu’il faudra estimer.

(A3) Cov(εi , εi 0 ) = 0, ∀i 6= i 0
ou de manière équivalente :
Cov(yi , yi 0 ) = 0, ∀i 6= i 0 .

Commentaire sur l’hypothèse (A3) :

- Sous cette hypothèse, les termes d’erreur εi sont non corrélés.

Chapitre 2 Régression linéaire multiple 6/40

 On peut écrire matriciellement le modèle (2) de la manière suivante :

Y = Xβ +  (3)

où
       
y1  1 x11 ... x1p  β0   ε1 
       
       
y  1 x21 ... x2p  β  ε 
 2    1  2
 . ,
Y =  X =
.
,
 . ,
β=  et  . .
= 
. . .. ..  . .
. . . .  . .

       
       
yn 1 xn,1 ... xnp βp εn

- Y désigne le vecteur à expliquer de taille n,

- X la matrice explicative de taille n × (p + 1),
-  le vecteur d’erreurs de taille n.
Exercice : Touver X et Y pour les données sur les appartements.

Chapitre 2 Régression linéaire multiple 7/40

Les hypothèses peuvent alors s’écrire sous forme matricielle :

(A1’) E() = 0n
ou de manière équivalente :
E(Y ) = X β ∈ Rn .

(A2’) V() = σ 2 In
ou de manière équivalente :
V(Y ) = σ 2 In .

Dans la suite de ce chapitre, on suppose que

n > (p + 1) et rang(X ) = p + 1

On a donc plus d’observations que de variables et il n’existe pas de liaison

linéaire entre les variables explicatives xj c’est à dire pas de multicolinéarité.

Chapitre 2 Régression linéaire multiple 8/40

 Remarque.

Il est important de bien faire la différence entre

- l’expression E(yi ) = β0 + pj=1 βj xij (qui désigne l’espérance d’une

P
variable aléatoire scalaire), et l’expression E(Y ) = X β (qui désigne
l’espérance d’une variable aléatoire vectorielle) : on obtient dans un cas
un scalaire, dans l’autre cas un vecteur de Rn .

- l’expression V(yi ) = σ 2 (qui désigne la variance d’une variable aléatoire

scalaire), et l’expression V(Y ) = σ 2 In (qui désigne la covariance d’une
variable aléatoire vectorielle) : on obtient dans un cas un scalaire (σ 2 ),
dans l’autre cas une matrice carrée (σ 2 In ) de dimension n × n.

Chapitre 2 Régression linéaire multiple 9/40

2. Estimation des paramètres β0 , β1 , . . . , βp , et σ 2

A partir de l’echantillon (aléatoire) de n observations

{(xi1 , . . . , xip , yi ), i = 1, . . . , n},

on veut estimer les paramètres

β0 , β1 , . . . , βp et σ 2 .

- Pour estimer β = (β0 , β1 , . . . , βp ), on peut utiliser la méthode des

moindres carrés qui ne nécessite pas d’hypothèse supplémentaire sur la
distribution de εi , contrairement à la méthode du maximum de
vraisemblance qui est fondée sur la normalité de εi .
- La méthode des moindres carrés ne fournit pas un estimateur de σ 2 .

Chapitre 2 Régression linéaire multiple 10/40

Estimation de β par les moindres carrés

On cherche βb ∈ Rp+1 qui minimise la somme des erreurs quadratiques

ε2i = (yi − β0 − β1 xi1 − . . . − βp xip )2

On doit donc résoudre le problème d’optimisation suivant :

n p
X X
βb = arg min [yi − (β0 + βj xij )]2 . (4)
β∈Rp+1
i=1 j=1

Vocabulaire :
- ŷi = βb0 + p
P
j=1 βbj xij est appelé la valeur prédite.
- ε̂i = yi − ŷi est appelé le résidu.

En notant xiT = (1, xi1 , . . . , xip ), la valeur prédite ŷi s’écrit

ŷi = xiT β.
b

Chapitre 2 Régression linéaire multiple 11/40

Résolution du problème d’optimisation

Le problème d’optimisation est :

min F (β),
β∈Rp+1

avec
n p
X X
F (β) = [yi − (β0 + βj xij )]2
i=1 j=1

= (Y − X β)T (Y − X β)
= Y T Y − 2Y T X β + β T X T X β

Le minimum est atteint pour

∂F (β)
= 0.
∂β

Rappels. Soient a et x deux vecteurs de dimension K , et soit A une matrice de

dimension K × K . On a :

∂aT x ∂x T a ∂x T Ax
= =a et = 2Ax si A est symétrique,.
∂x ∂x ∂x

Chapitre 2 Régression linéaire multiple 12/40

Solution du problème d’optimisation
On en déduit après quelques manipulations :

βb = (X T X )−1 X T Y , (5)
sous réserve que X T X soit inversible.

Commentaires
Pn
- Le minimum de F est égal à i=1 ε̂2i . Ce minimum est appelé la somme
des carrés des résidus (SCR).
Pp
- La valeur prédite ybi estime E[yi ] = β0 + j=1 βj xij et non pas yi . Une
meilleure notation serait E[y
di ].

- Aucune des hypothèses n’a été utilisée ici pour obtenir β.

b

Chapitre 2 Régression linéaire multiple 13/40

Propriétés de βb

Sous les hypothèses (A1’) et (A2’), on peut montrer que

- E[β]
b = β,
b = σ 2 (X T X )−1
- V(β)

Commentaires
- L’estimateur βb est sans biais.
- Il est aussi de variance minimale parmi tous les estimateurs linéaires par
rapport à Y ) sans biais (propriété dite de Gauss-Markov).

Chapitre 2 Régression linéaire multiple 14/40

Estimation de σ 2
Le paramètre σ 2 est défini par

σ 2 = V(εi ) = V(yi )= E (yi − E[yi ])2 .

 

En prenant ŷi = xiT βb comme estimateur de E[yi ], il apparaît naturel d’estimer

σ 2 par

n Pn
2 1 X 2 εi )2
i=1 (b SCR
s = (yi − ŷi ) = = .
n − (p + 1) i=1 n−p−1 n−p−1

Commentaires

- s 2 est un estimateur sans biais de σ 2

- La perte de p + 1 degrés de liberté dans l’expression de s 2 est le “coût” de
l’estimation de β0 , β1 , . . . , βp nécessaire pour obtenir les ŷi .

Chapitre 2 Régression linéaire multiple 15/40

Sorties R des données poids de naissance

summary(modele)

##
## Call:
## lm(formula = bwt ~ age + weight + smoke, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1961 -308 11 309 1487
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3050.562 108.861 28.02 < 2e-16 ***
## age -0.918 2.535 -0.36 0.72
## weight 7.903 1.568 5.04 5.4e-07 ***
## smoke -254.254 29.939 -8.49 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 499 on 1170 degrees of freedom
## Multiple R-squared: 0.081,Adjusted R-squared: 0.0786
## F-statistic: 34.4 on 3 and 1170 DF, p-value: <2e-16

Chapitre 2 Régression linéaire multiple 16/40

3. Tests d’hypothèses et intervalles pour les paramètres βj

On veux maintenant tester la nullité des coefficients βj du modèle de

régression.

Pour faire ces tests, il est nécessaire de faire une hypothèse supplémentaire :

(A3)’  ∼ N (0n , σ 2 In )
ou de manière équivalente
Y ∼ N (X β, σ 2 In ).

Commentaire. L’unique "nouveauté" ici est la normalité.

Chapitre 2 Régression linéaire multiple 17/40

Test de signification du modèle

Typiquement, on commence par tester :

H0 : “β1 = . . . = βp = 0” contre H1 : “∃ j ∈ {1, . . . , p}, βj 6= 0”.

On utilise la statistique suivante :

Pn 2
i=1 (ŷi − ȳn ) /p SCE /p
Fn = Pn 2
=
i=1 (yi − ŷi ) /(n − p − 1) SCR/(n − p − 1)

qui est distribuée sous H0 selon une loi de Fisher à p et n − p − 1 degrés de

libertés. On rejette H0 avec un risque 0 ≤ α ≤ 1 si

Fn ≥ f1−α (p, n − p − 1)

où f1−α (p, n − p − 1) est le fractile d’ordre 1 − α de la loi F (p, n − p − 1).

Chapitre 2 Régression linéaire multiple 18/40

Table d’analyse de la variance (ANOVA) :

Source de variation Somme des carrés ddl carré moyen F

− ȳn )2 − ȳn )2
Pn 1 Pn SCE /p
régression (expliquée) SCE = i=1 (b
yi p p i=1 (b
yi SCR/(n−p−1)

− ybi )2 1 − ybi )2
Pn Pn
Résiduelle SCR = i=1 (yi n-(p+1) n−p−1 i=1 (yi

− ȳn )2 1 − ȳn )2
Pn Pn
Totale SCT = i=1 (yi n-1 n−1 i=1 (yi

Remarques :
- On retrouve la statistique dite de Fisher Fn qui permet de tester
l’ajustement du modèle.
- On retrouve la propriété fondamentale SCT = SCE + SCT qui permet de
mesurer l’ajustement du modèle par le coefficient de détermination
SCE SCR
R2 = =1− .
SCT SCT

- Le coefficient R 2 donne la proportion de variabilité de y qui est expliquée

par le modèle. Plus le R 2 est proche de 1, meilleure est l’adéquation du
modèle aux données.

Chapitre 2 Régression linéaire multiple 19/40

Test de significativité d’un paramètre βj

On désire maintenant tester :

H0 : “βj = 0” contre H1 : “βj 6= 0”

Nouvelles propriétés pour les estimateurs βbj et s 2

Sous les hypothèses (A1’)-(A3’), on a :
(a) βbj ∼ N βj , σ 2 cjj

où cjj est le terme (j + 1, j + 1) de la matrice (X T X )−1

(n − p − 1)s 2
(b) 2
∼ χ2 (n − p − 1)
σ
(c) βbj et s 2 sont indépendants

Un rappel de probabilité
U
Si U ∼ N (0, 1), V ∼ χ2 (ν) et U est indépendant de V , alors q ∼ T (ν).
V
ν

Chapitre 2 Régression linéaire multiple 20/40

On déduit alors des propriétés (a)-(c) que

βbj − βj
p
σ 2 cjj βbj − βj
s = √ ∼ T (n − p − 1).
(n−p−1)s 2 s cjj
σ2
n−p−1

On utilisera donc la statistique suivante :

βbj − βj
Tn = √ ,
s cjj

qui est distribuée selon une loi de Student à n − p − 1 degrés de libertés.

Chapitre 2 Régression linéaire multiple 21/40

Test de H0 contre H1

Sous l’hypothèse H0 : “βj = 0”, on a

βbj
Tn = √ ∼ T (n − p − 1). (6)
s cjj

Pour une hypothèse alternative H1 : “βj 6= 0” bilatérale, on rejette H0 avec un

risque 0 ≤ α ≤ 1 si
|t| ≥ t1−α/2 (n − p − 1)
où t est la réalisation de Tn et t1−α/2 (n − p − 1) est le fractile d’ordre 1 − α/2
de la loi T (n − p − 1).

Chapitre 2 Régression linéaire multiple 22/40

Remarques.
Pour réaliser ce test, on peut également :
- regarder la p-valeur aussi appelée niveau de signification du test : si
p-valeur ≤ α, on rejette H0 . Dans le cas d’un test bilatéral
(H1 : “β1 6= 0”), on a :

p-valeur = P(|Tn | > |t| / H0 ). (7)

On rejette H0 si p-valeur ≤ α.

- construire l’intervalle de confiance de βj :

√
[βbj ± t1−α/2 (n − p − 1)s cjj ].

On rejette H0 si 0 n’appartient pas à cet intervalle.

Chapitre 2 Régression linéaire multiple 23/40

Rejeter H0 signifie :

- que le coefficient βj est significativement non nul,

- que βj s’interprète comme le taux d’accroissement moyen de y en
fonction d’une variation de xj lorsque tous les autres régresseurs
x1 , . . . , xj−1 , xj+1 , . . . , xp restent fixés.

Exemple des données poids de naissance.

summary(modele)$coefficients

## Estimate Std. Error t value Pr(>|t|)

## (Intercept) 3050.56238 108.8612 28.0225 1.2569e-132
## age -0.91802 2.5353 -0.3621 7.1734e-01
## weight 7.90266 1.5675 5.0414 5.3511e-07
## smoke -254.25425 29.9388 -8.4925 6.0621e-17

confint(modele)

## 2.5 % 97.5 %
## (Intercept) 2836.9775 3264.1473
## age -5.8922 4.0562
## weight 4.8271 10.9782
## smoke -312.9940 -195.5145

Chapitre 2 Régression linéaire multiple 24/40

Contribution jointe d’un ensemble de régresseurs
On peut maintenant tester la nullité de q ≤ p paramètres :

H0 : “β1 = . . . = βq = 0” contre H1 : “∃ j ∈ {1, . . . , q}, βj 6= 0”.

Cela revient à comparer deux modèles :

- le modèle complet à p regresseurs (modèle 1) pour lequel on évalue la
somme des carrés des résidus SCR1 ,
- le modèle réduit à p − q regresseurs (modèle 0) pour lequel on évalue la
somme des carrés des résidus SCR0 .

On peut montrer que sous H0 :

(SCR0 − SCR1 )/q
∼ F (q, n − p − 1).
SCR1 /(n − p − 1)

Chapitre 2 Régression linéaire multiple 25/40

La zone de rejet associée à cette statistique de test est donc :

R = ]f1−α (q, n − p − 1), +∞[.

Rejeter H0 signifie qu’au moins un des q coefficients est non nul.

Exemple des données poids de naissance.

modele0 <- lm(bwt~ smoke,data=data)

modele1 <- lm(bwt~ age+weight+smoke,data=data)
anova(modele0,modele1)

## Analysis of Variance Table

##
## Model 1: bwt ~ smoke
## Model 2: bwt ~ age + weight + smoke
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1172 297411671
## 2 1170 291055628 2 6356043 12.8 3.2e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Chapitre 2 Régression linéaire multiple 26/40

3. Prévision d’une valeur ultérieure

On désire prévoir à l’aide du modèle la valeur de la variable y pour des

observations futures (x1,0 , . . . , xp,0 ) des p variables explicatives.
Posons
x0 = (1, x1,0 , . . . , xp,0 )T ∈ Rp+1
D’après le modèle on a :
y0 = x0 T β + ε0 ,
et la prédiction est :
[ Tb
yb0 = E[y0 ] = x0 β.

L’erreur de prédiction est définie par yb0 − y0 et on peut montrer que sous les
hypothèses du modèle (incluant l’hypothèse de normalité), on a :
  
2 T T −1
yb0 − y0 ∼ N 0, σ 1 + x0 (X X ) x0 . (8)

Chapitre 2 Régression linéaire multiple 27/40

On en déduit que :
y0 − ŷ0
p ∼ N (0, 1).
σ 1 + x0T (X T X )−1 x0
On peut montrer que :
y0 − ŷ0
p ∼ T (n − p − 1).
s 1 + x0T (X T X )−1 x0

On utilise ce résultat pour construire un intervalle de prédiction pour y0 , c’est à

dire l’intervalle [A, B] tel que

P(A ≤ y0 ≤ B) = 1 − α.

Ici, y0 est une variable aléatoire et non pas un paramètre. L’intervalle de

prédiction est donc un intervalle dans lequel une future observation y0 va
tomber avec une certaine probabilité (différent d’un intervalle de confiance).

Chapitre 2 Régression linéaire multiple 28/40

On en déduit l’intervalle de prédiction pour y0 au niveau de confiance 1 − α
suivant :  q 
ŷ0 ± t1−α/2 (n − p − 1)s 1 + x0T (X T X )−1 x0

On peut aussi construire un intervalle de confiance de la valeur moyenne

E[y0 ] = x0T β,

qui est cette fois un paramètre. On va donc chercher l’intervalle aléatoire [A, B]
tel que
P(A ≤ E[y0 ] ≤ B) = 1 − α.

Chapitre 2 Régression linéaire multiple 29/40

Pour construire cet intervalle, on montre que :

 
ŷ0 ∼ N x00 β, σ 2 x0T (X T X )−1 x0 ,
ŷ0 − x0T β
p ∼ T (n − p − 1).
s x0T (X T X )−1 x0

On en déduit l’intervalle de confiance de E[y0 ] = x0T β suivant :

 q 
yb0 ∓ t1−α/2 (n − p − 1) s x0T (X T X )−1 x0 .

Exemple des données poids de naissance.

#prevision de l'age du bebe d'une femme de 30 ans, 50 kg et fumeuse
predict(modele,data.frame(age=30,weight=50,smoke=1),interval="prediction")

## fit lwr upr

## 1 3163.9 2183.8 4144

predict(modele,data.frame(age=30,weight=50,smoke=1),interval="confidence")

## fit lwr upr

## 1 3163.9 3109.1 3218.7

Chapitre 2 Régression linéaire multiple 30/40

3. Analyse des résidus

Exemple des données poids de naissance.

#On calcule les residus studentises

residus=rstudent(modele)

#On calcule les valeurs predites

val.predites <- predict(modele)

#Graphique predictions-residus
plot(val.predites ,residus, xlab="val.predites", ylab="Residus")
abline(h=c(-2,0,2), lty=c(2,1,2),col=c(1,2,1))

3
2
1
0
Residus

−1
−2
−3
−4

3200 3400 3600 3800

val.predites

Chapitre 2 Régression linéaire multiple 31/40

#normalite des residus

shapiro.test(residus)

##
## Shapiro-Wilk normality test
##
## data: residus
## W = 0.993, p-value = 0.000022

hist(residus)

Histogram of residus
250
200
150
Frequency

100
50
0

−4 −3 −2 −1 0 1 2 3

residus

Chapitre 2 Régression linéaire multiple 32/40

3. Sélection de variables

Il s’agit maintenant de sélectionner parmi les p variables explicatives, les q ≤ p

variables qui donnent le “meilleur" modèle pour prédire y .

Il faut donc :
- un critère de qualité d’un modèle afin de comparer deux modèles n’ayant
pas nécessairement le même nombre de variables explicatives.
- une procédure qui permet de choisir parmi tous les modèles, le meilleur
au sens de ce critère. On parle de procédure de choix de modèle.

Un problème de complexité :
p
X
- Le nombre de modèles à considérer est Cpq = 2p − 1. Ce nombre croît
q=1
exponentiellement avec p. Par exemple, si p = 30, on devrait considérer
230 = 109 modèles...
- En pratique, on utilise donc des heuristiques dont les plus simples sont les
procédures pas à pas ascendante ou descendante.

Chapitre 2 Régression linéaire multiple 33/40

Les critères R 2 et R 2 ajusté

- Le coefficient R 2 = 1 − SCR
SCT

- mesure l’ajustement du modèle aux données,

- augmente lorsque le nombre de variables incluses dans le
modèle augmente,
- permet de comparer des modèles ayant le même nombre de
variables.
2 SCR/(n−p−1)
- Le coefficient Rajuste =1− SCT /(n−1)
2
2 V(ε) σ
- estime le Rpopulation = 1 − V(Y ) = 1 − σY2 ,
- n’augmente pas forcément lorsque le nombre de variables
introduites dans le modèle augmente,
- permet de comparer des modèles ayant un nombre de variables
différent.

Chapitre 2 Régression linéaire multiple 34/40

Les critères AIC et BIC .

Ce sont deux critères de vraisemblance pénalisées définis par :

- AIC = −2ln(L) + 2k : Akaike Information Criterion
- BIC = −2ln(L) + kln(n) : Bayesian Information Criterion
où L est la vraisemblance maximisée et k est le nombre de paramètres libres du
modèle.

En régression multiple :
- il y a q + 2 paramètres β0 , β1 , . . . , βq , σ et une equation donc k = q + 1
paramètres libres.
- la vraisemblance est définie comme la densité conjointe des yi et son
expression est
 
2 1 1 T
L(β, σ ) = exp − 2 (Y − X β) (Y − X β)
(2πσ 2 )n/2 2σ

Les estimateurs du maximum de vraisemblance sont β̃ = (X T X )−1 X T Y

et σ̃ 2 = SCR
n
. La vraisemblance maximisée est L = L(β̃, σ̃ 2 ) et on obtient :

−2ln(L) = n(ln(2πσ̃ 2 ) + 1)

Chapitre 2 Régression linéaire multiple 35/40

Ecriture simplifiée en régression multiple.

AIC = n ln(SCR) + 2k + cste

BIC = n ln(SCR) + k ln(n) + cste

Ces critères doivent être minimisés dans une procédure de choix de modèle.

Chapitre 2 Régression linéaire multiple 36/40

Procédure pas à pas ascendante (forward stepwise).
- On part du modèle nul sans variable.
- On effectue p régressions linéaires simples et on sélectionne le modèle qui
minimise le critère AIC .
- On effectue p − 1 régressions linéaires avec 2 variables explicatives et on
sélectionne le modèle qui minimise le critère AIC .
- On recommence jusqu’à ce que le critère AIC ne diminue plus.

Procédure pas à pas descendante (backward stepwise).

On part cette fois du modèle complet à p variables explicatives et on supprime
pas à pas les variables. Le test d’arrêt et le critère sont les mêmes que pour la
procédure ascendante.

Chapitre 2 Régression linéaire multiple 37/40

Exemple des données poids de naissance.

full <- lm(bwt ~ gestation + age + weight + smoke, data=data)

null <- lm(bwt ~ 1, data=data)
back <- step(full, direction="backward")

## Start: AIC=14379
## bwt ~ gestation + age + weight + smoke
##
## Df Sum of Sq RSS AIC
## - age 1 60465 242763748 14377
## <none> 242703282 14379
## - weight 1 5352463 248055745 14402
## - smoke 1 14379595 257082877 14444
## - gestation 1 48352346 291055628 14590
##
## Step: AIC=14377
## bwt ~ gestation + weight + smoke
##
## Df Sum of Sq RSS AIC
## <none> 242763748 14377
## - weight 1 5637978 248401726 14402
## - smoke 1 14556053 257319800 14443
## - gestation 1 48324498 291088245 14588

formula(back)

## bwt ~ gestation + weight + smoke

Chapitre 2 Régression linéaire multiple 38/40

forw <- step(null, scope=list(lower=null,upper=full), direction="forward", trace = 1)

## Start: AIC=14683
## bwt ~ 1
##
## Df Sum of Sq RSS AIC
## + gestation 1 52601391 264100599 14472
## + smoke 1 19290318 297411671 14611
## + weight 1 7699680 309002310 14656
## <none> 316701990 14683
## + age 1 230584 316471406 14684
##
## Step: AIC=14472
## bwt ~ gestation
##
## Df Sum of Sq RSS AIC
## + smoke 1 15698873 248401726 14402
## + weight 1 6780798 257319800 14443
## + age 1 754996 263345603 14471
## <none> 264100599 14472
##
## Step: AIC=14402
## bwt ~ gestation + smoke
##
## Df Sum of Sq RSS AIC
## + weight 1 5637978 242763748 14377
## <none> 248401726 14402
## + age 1 345981 248055745 14402
##
## Step: AIC=14377
## bwt ~ gestation + smoke + weight
##
## Df Sum of Sq RSS AIC
## <none> 242763748 14377
## + age 1 60465 242703282 14379

Chapitre 2 Régression linéaire multiple 39/40

formula(forw)

## bwt ~ gestation + smoke + weight

summary(lm(forw,data=data)) # idem forw

##
## Call:
## lm(formula = forw, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1471.9 -305.0 -7.9 276.2 1455.9
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -499.702 246.393 -2.03 0.043 *
## gestation 12.703 0.832 15.26 < 2e-16 ***
## smoke -229.004 27.341 -8.38 < 2e-16 ***
## weight 7.386 1.417 5.21 2.2e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 456 on 1170 degrees of freedom
## Multiple R-squared: 0.233,Adjusted R-squared: 0.231
## F-statistic: 119 on 3 and 1170 DF, p-value: <2e-16

Chapitre 2 Régression linéaire multiple 40/40