2 juin 2006 1
Introduction à la régression
cours n°3
Influence d’un prédicteur
ENSM.SE – 1A
Olivier Roustant
2 juin 2006 2
Objectif du cours
Utiliser les résultats théoriques sur
l’estimation des paramètres pour savoir si un
prédicteur d’un modèle linéaire est influent
2 juin 2006 3
Influent ou non influent ?
yi = β0 + β1xi + ei avec e1, …, e4 i.i.d N(0, 0.042)
Imaginer un pèse-
personne avec une
erreur de mesure
de 50 kg !
"ˆ1 (# ) < 0 "ˆ1 (# ') > 0 !!
En fait ici, l’erreur d’estimation sur β1 = |β1| !! (=0.2)
! !
2 juin 2006 4
Influent ou non influent (suite)
Le même ex., mais en planifiant mieux les expériences
Cette fois, l’erreur d’estimation sur la pente = 0.0536
2 juin 2006 5
Exercice
Vous pouvez réaliser n expériences pour estimer un
phénomène linéaire sur [a,b] impliquant 1 prédicteur
Comment répartir les expériences dans le domaine expérimental
[a,b] de façon à ce que l’estimation soit la plus précise possible ?
2 juin 2006 6
Influent ou non influent (suite)
Le même exemple, planification optimale (voir diapo. 26)
L’erreur d’estimation sur la pente = 0.04
2 juin 2006 7
Influent ou non influent ?
Morale de l’exemple
Prendre en compte l’erreur d’estimation d’un
paramètre pour savoir s’il est important ou pas
→ Décision en milieu incertain : test statistique
L’impossibilité de décider peut venir d’une
mauvaise planification des expériences
2 juin 2006 8
Formalisation
Considérons le modèle linéaire
yi = β0 +β1x1,i + … + βpxp,i + ei
avec e1, …, en i.i.d N(0,σ2)
Le prédicteur xi est influent si βi ≠ 0
Test statistique opposant les hypothèses
{βi = 0} et {βi ≠ 0}
2 juin 2006 9
Construction du test statistique
1ère étape : hypothèse H0
On veut contrôler le risque de décider qu’un
prédicteur est influent alors qu’il ne l’est pas.
Quelle est l’hypothèse H0 ?
H 0 = {" i = 0}
Autre raison : pouvoir faire les calculs !
!
2 juin 2006 10
Construction du test (suite)
2ème étape : choix d’une statistique de décision
On part de l’estimateur des moindres carrés (EMC)
de βi
Matriciellement, on vérifie qu’on a :
(Y-Xβ)’(Y-Xβ) minimum ssi X’(Y-Xβ) = 0
ssi (X’X) β = X’Y
D’où
ˆ #1
" = (X' X) X'Y
2 juin 2006 11
Construction du test (suite)
2ème étape (suite) : loi de l’EMC sous H0 ?
Une combinaison linéaire de v.a. de lois normales
indépendantes est encore de loi normale (admis)
Ici, "ˆ = (X' X)#1 X'Y = (X' X)#1 X'(X" + e)
d’où "ˆ = " + (X' X)#1 X'e
! Chaque "ˆ i est donc une combinaison linéaire des e ,
i
! centrée sur β ⇒ loi normale centré sur β
i i
!
2 juin 2006 12
Construction du test (suite)
2ème étape (suite)
Exercice. Pour un vecteur u, n×1, on définit la
matrice de covariance par cov(u)= (cov(ui,uj))1≤i,j≤n
a) Mq Cov(u) = E((u-m)(u-m)’), avec m=E(u) (vect. n×1)
b) Mq Cov(e) = σ2 In
c) Déduire de a) et b) que cov("ˆ ) = # 2 (X' X)$1
puis que var("ˆ i ) := # i 2 = # 2 ((X' X)$1 ) ii
Conclusion : sous
! H , l’EMC est de loi N(0, σ 2)
0 i
!
2 juin 2006 13
Construction du test (suite)
2ème étape (suite)
La loi de l’estimateur dépend des paramètres
Intuitivement, sous H0 :
ˆ
"
"ˆ i ~ N(0,# ((X' X) ) ii ) %
2 $1 i
$1
& N(0,1)
#ˆ ((X' X) ) ii
Résultat exact : remplacer la loi N(0,1) par la loi de
Student tn-p-1.
!
"ˆ i
Conclusion : choix de la statistique T =
#ˆ ((X' X)$1 ) ii
2 juin 2006 14
Construction du test (suite)
2ème étape (résumé)
"ˆ i
Choix de la statistique de décision T =
#ˆ ((X' X)$1 ) ii
Interprétation : estimation du paramètre rapporté à
son écart-type d’estimation
!
Vocabulaire : T est appelé t-ratio (à cause de la loi
de Student, notée t)
– Propriété : T est de loi de Student tn-p-1
→ En pratique, dès que n-p-1≥20, on approche tn-p-1 par
N(0,1)
2 juin 2006 15
Construction du test (suite)
3ème étape : détermination d’un seuil
"ˆ i,obs
Notation : Tobs =
#ˆ obs ((X' X)$1 ) ii
Au niveau 5%, on rejette H0
–!n-p-1≥20 : si Tobs dépasse 1.96 en valeur absolue
– n-p-1<20 : utiliser les tables de la loi de Student
– Mieux (dans tous les cas) : utiliser la p-valeur
2 juin 2006 16
Construction du test (fin)
3ème étape : p-valeur
On appelle p-valeur la probabilité d’obtenir
pire que ce qu’on a :
p - valeur = PH 0 ( T > Tobs )
Permet de ne pas avoir à calculer de seuil :
p-valeur < 5% ⇒ rejet au niveau 5%
!
p-valeur < 1% ⇒ rejet au niveau 1%
…
2 juin 2006 17
Test de signification : pratique
En pratique, les logiciels donnent le tableau
suivant :
erreur
coefficient estimation t " ratio p " valeur
d'estimation
#ˆ i,obs
#i #ˆ i,obs $ˆ i,obs Tobs = PH 0 ( T > Tobs )
$ˆ i,obs
!
2 juin 2006 18
Exemple 1
Données de pollution (cf cours 1)
2 juin 2006 19
Régression avec R
Le fichier de données : NO3 SO4
(format .txt) 0,45 0,78
0,09 0,25
lm : 1,44 2,39
linear model … …
> pollution <- read.table("pollution.txt", header=TRUE, dec=",", sep="\t")
> modele_degre_1 <- lm(log(NO3)~log(SO4), data=pollution)
> summary(modele_degre_1)
> modele_degre_2 <- lm(log(NO3)~log(SO4)+I(log(SO4)^2), data=pollution)
> summary(modele_degre_2)
2 juin 2006 20
Sorties à commenter
Call:
lm(formula = log(NO3) ~ log(SO4),
Comme
datan-p-1>20,
= pollution)
on peut aussi se baser
sur le fait que |t-ratio| > 2
Residuals: ou
Min 1Q Median 3Q Maxque l’erreur d’estimation est
-0.80424 -0.14485 -0.01087 0.16564 0.56666
< la moitié de l’estimation
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.43642 0.03679 -11.86 <2e-16 ***
log(SO4) 0.92168 0.03356 27.47 <2e-16 ***
---
p-valeur
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1<‘ 0.05
’1
⇒ paramètres significatifs au niveau 5%
Residual standard error: 0.2417 on (on
165est mêmeoftrès
degrees large : p=2e-16 !)
freedom
Multiple R-Squared: 0.8205, Adjusted R-squared: 0.8195
F-statistic: 754.4 on 1 and 165 DF, p-value: < 2.2e-16
2 juin 2006 21
Call:
lm(formula = log(NO3) ~ log(SO4)
Comme
+ I(log(SO4)^2),
n-p-1>20, ondata
peut
= pollution)
aussi se baser
sur le fait que |t-ratio| < 2
Residuals: ou
Min 1Q Median 3Q Max que l’erreur d’estimation est
-0.79819 -0.14085 -0.01470 0.16158 0.57136
> la moitié de l’estimation
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.42918 0.03955 -10.852 <2e-16 ***
log(SO4) 0.95337 0.07098 13.432 <2e-16 ***
I(log(SO4)^2) 0.01886 0.03720 0.507 0.613
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2423 on 164 degreesp-valeur > 0.05
of freedom
Multiple R-Squared: 0.8208, paramètreR-squared:
⇒ Adjusted non significatif
0.8186au niveau 5%
F-statistic: 375.7 on 2 and 164 DF, p-value: < 2.2e-16
2 juin 2006 22
Exemple 2
Retour sur les simulations (cf transp. n°3)
yi = β0 + β1xi + ei avec e1, …, e4 i.i.d N(0, 0.042)
2 juin 2006 23
Call:
La t-valeur est > 1.96 en valeur absolue,
lm(formula = ysim ~ experiences)
Pourtant on ne rejette pas H0
Cela est dû au fait qu’on ne peut pas utiliser
Residuals:
l’approximation normale (ici n=4 << 20)
1 2 3 4
La p-valeur est calculée à partir de la loi de Student
0.038612 -0.038612 0.002785 -0.002785
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.1404 0.0987 11.554 0.00741 **
experiences -0.4099 0.1936 -2.118 0.16838
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Moralité : la pente de la droite est négative,
Residual standard error: 0.03871Maison 2l’erreur
degreesd’estimation
of freedom est trop importante
Multiple R-Squared: 0.6916, Et le paramètre
Adjusted est statistiquement
R-squared: 0.5374 non
significatif
F-statistic: 4.485 on 1 and 2 DF, p-value: au niveau 5% …rassurant !
0.1684
2 juin 2006 24
Exemple 2 (suite)
Retour sur les simulations (cf transp. n°6)
yi = β0 + β1xi + ei avec e1, …, e4 i.i.d N(0, 0.042)
2 juin 2006 25
Call:
lm(formula = ysim ~ experiences)
Moralité : la pente de la droite est négative,
Residuals: Cette fois l’erreur d’estimation est assez faible
1 2 3 4 Et le paramètre est statistiquement
0.01300 -0.01300 0.03190 -0.03190
significatif au niveau 5% (mais pas 1%)
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.97956 0.02436 40.22 0.000618 ***
experiences -0.26142 0.03444 -7.59 0.016921 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.03444 on 2 degrees of freedom
Multiple R-Squared: 0.9664, Adjusted R-squared: 0.9497
Remarque : la pente
F-statistic: 57.6 réelle (inconnue)
on 1 and est
2 DF, p-value: 0.01692
- 0.2
2 juin 2006 26
Exercice : planification des
expériences en dimension 1
Vérifier que, en dimension 1
1 1 " x (x%
2
"1 x% (1
X' X = n$ ' puis ( X' X ) = $ '
#x x 2& n x ( x #(x 1 &
2 2
2
# 1
puis var "ˆ1 = cov("ˆ1, "ˆ1 ) = (# 2 (X' X)$1 )11 =
( ) n x2 $ x2
! En déduire que pour minimiser l’erreur d’estimation de la pente, il
faut que la variance empirique des xi soit la plus grande possible
! Prenons n pair, et considérons le domaine expérimental [-1,1]. Montrer que
le maximum est atteint lorsque la moitié des points est placée sur le bord
gauche (en x = -1), et l’autre moitié sur le bord droit (x = 1)