Dans ce notebook l'exemple qui nous servira à comprendre la régression multiple avec R est une étude sur le poids de naissance des bébés.
Il s'agit d'une enquête concernant les facteurs de risque associés au faible poids de naissance de nourrissons (données collectées au centre médical de Baystate dans le Massachusetts pendant l'année 1986). Le faible poids de naissance est un événement qui intéresse les médecins depuis plusieurs années en raison du taux de mortalité infantile et du taux d'anomalies infantiles très élevés chez les nourrissons de faible poids. Le comportement d'une femme pendant la grossesse (régime alimentaire, habitudes tabagiques...) peut altérer de façon importante les chances de mener la grossesse à terme, et par conséquent, de donner naissance à un enfant de poids normal. Le fichier de données contient les informations sur 189 femmes venant consulter dans le centre médical. Les variables sont :
Nous considérons le modèle linéaire suivant: $$\text{BWT} = \beta_0 + \beta_1\text{AGE} + \beta_2\text{LWT}+\beta_3\text{SMOKE}+\varepsilon.$$
Commençons par importer les données sous R :
naiss <- read.table('poids_naiss',header=T,sep=':')
head(naiss)
attach(naiss)
Avant d'estimer les paramètres, nous calculons la matrice de corrélation et nous présentons un diagramme de dispersion de toutes les paires de ces variables. Ceci permet de visualiser la relation entre la variable à expliquer et chacune des variables explicatives et de juger de la corrélation entre les variables explicatives.
cor(naiss)
pairs(naiss)
Comme en régression linéaire simple, nous estimons le modèle par la fonction lm:
mod1 <- lm(BWT ~ AGE+LWT+SMOKE)
mod1
Il existe une écriture plus rapide pour définir un modèle linéaire avec toutes les variables d'un dataframe :
mod1 <- lm(BWT~.,data=naiss)
mod1
La dimension des données étant trop élevée, on ne peut plus représenter l'ensemble des données et la fonction de régression dans un seul graphique. Mais on peut faire des représentations graphiques variable par variable.
Notons qu'on peut couper la fenêtre graphique en plusieurs partie par l'instruction par(mfrow=c(m,n)) afin de tracer m*n graphiques dans une fenêtre (organisée dans une matrice de taille mxn)
par(mfrow=c(1,3))
plot(BWT ~ AGE)
abline(coef(mod1)[1],coef(mod1)[2])
plot(BWT ~ LWT)
abline(coef(mod1)[1],coef(mod1)[3])
plot(BWT ~ SMOKE)
abline(coef(mod1)[1],coef(mod1)[4])
Les tests sur les divers paramètres sont obtenus grâce à la fonction summary.
summary(mod1)
Les résultats fournis par summary se présentent de façon identique à ceux de la régression linéaire simple. On y retrouve les estimations des paramètres de régression dans la colonne Estimate.
Les valeurs réalisées des statistiques des tests de Student associés aux hypothèses $H_0 : \beta_i = 0$ versus $H_l : \beta_i \neq 0$ se trouvent dans la colonne t value, les p-valeurs associées dans la colonne Pr(>|t|). Residual standard error fournit l'estimation de $\sigma$ ainsi que le nombre de degrés de liberté associés $n - p - 1$. On trouve enfin le coefficient de détermination $R^2$ (Multiple R-squared) ainsi que la version ajustée $R^2_a$ (Adjusted R-squared). Enfin, on trouve la réalisation du test de Fisher global (F-statistic) ainsi que sa p-valeur (p-value) associée.
Le test de Fisher global permet de tester l'apport global et conjoint de l'ensemble des variables explicatives présentes dans le modèle pour "expliquer" les variations de $Y$. L'hypothèse nulle ici est $H_0 : \beta_1= \dots =\beta_3 = 0$ (l'ensemble des $p$ variables explicatives n'apporte pas une information utile pour la prédiction de $Y$, sous le modèle linéaire). L'assertion d'intérêt est $H_1:$ au moins l'un des coefficients $\beta_j$ est significativement différent de zéro (au moins une des variables explicatives et associée à $Y$ après ajustement sur les autres variables explicatives).
L'instruction suivante produit un certain nombre de graphiques qui peuvent servir à analyser les résidus :
plot(mod1)
Nous reprenons les données de la concentration d'ozone du dernier TP que nous allons analyser davantage. Nous expliquons la concentration d'ozone maxO3 par dix variables :