Avant de parler davantage du premier devoir, un mot pour expliquer pourquoi dans la première partie il est important d’éviter d’utiliser des données temporelles (le cours présente des outils pour modéliser des données individuelles). Considérons un exemple simple (évoqué voilà plusieurs mois ici même), où on cherche à prédire un tour de poitrine conditionnellement à un tour de taille (pour une personne de sexe féminin, cf. le mannequin ci-contre). L’idée intuitive est qu’il doit exister une relation positive, possiblement linéaire entre les deux mensurations. A partir de données, on devrait pouvoir quantifier cette relation. On dispose pour cela de données, en ligne sur
> mensurations<-read.table(" + http://freakonometrics.free.fr/mensurations.csv", + header=TRUE,sep=";") > head(mensurations[,-1]) POITRINE TAILLE HANCHE HAUTEUR POIDS AGE 1 NA NA NA NA NA NA 2 94 58 89 165 51 22.24230 3 NA NA NA NA NA NA 4 94 58 89 165 51 22.40383 5 91 58 97 165 55 NA 6 NA NA NA NA NA NA
Si on regarde rapidement la corrélation entre la taille et la poitrine, cette dernière est faible,
> Z=mensurations[,c("TAILLE","POITRINE")] > I=apply(Z,1,function(x) sum(is.na(x))) > cor(Z[I==0,]) TAILLE POITRINE TAILLE 1.00000000 0.09700002 POITRINE 0.09700002 1.00000000
… mais significativement non nulle,
> cor.test(Z[,1],Z[,2],,method = "pearson") Pearsons product-moment correlation data: Z[, 1] and Z[, 2] t = 2.4559, df = 635, p-value = 0.01432 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: 0.0194634 0.1733769 sample estimates: cor 0.09700002
Avec une -value de l’ordre de 1.5%, on va rejeter l’hypothèse
de nullité de la corrélation (c’est ici un test basé sur la transformation de Fisher). On peut aussi faire une régression simple,
> reg=lm(POITRINE~TAILLE,data=mensurations) > summary(reg) Call: lm(formula = POITRINE ~ TAILLE, data = mensurations) Residuals: Min 1Q Median 3Q Max -9.7253 -3.4424 -0.2091 1.5576 13.8579 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 82.67798 2.82953 29.220 <2e-16 *** TAILLE 0.11663 0.04749 2.456 0.0143 * --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 4.133 on 635 degrees of freedom (22 observations deleted due to missingness) Multiple R-squared: 0.009409, Adjusted R-squared: 0.007849 F-statistic: 6.031 on 1 and 635 DF, p-value: 0.01432
(on notera que le test sur la corrélation, le test de significativité de la – seule – variable explicative, et le test de significativité globale donnent ici exactement la même -value de l’ordre de 1.5%). Là encore, on voit que le tour de taille a un effet significativement non nul… mais marginal. Graphiquement, la prédiction avec le tour de taille est le trait rouge,
La valeur moyenne est en trait pointillé. On peut voir que la différence est infime,
> base=data.frame(TAILLE=seq(50,70,by=5)) > base$Pred.POITRINE=predict(reg,newdata=base) > base$Moy.POITRINE=mean(POITRINE,na.rm=TRUE) > base TAILLE Pred.POITRINE Moy.POITRINE 1 50 88.50935 89.61538 2 55 89.09248 89.61538 3 60 89.67562 89.61538 4 65 90.25876 89.61538 5 70 90.84189 89.61538
Le soucis est que cette base est en fait une base sur données individuelles. En fait, il s’agit des mensurations des demoiselles qui ont posé comme playmate dans la version U.S. du magasine Playboy. Si on se restreint sur des sous-période, par exemple les années 70, on a des prédictions relativement différentes,
> submensurations=mensurations[abs(DATE-1975)<=5,] > subreg=lm(POITRINE~TAILLE,data=submensurations) > summary(subreg) Call: lm(formula = POITRINE ~ TAILLE, data = submensurations) Residuals: Min 1Q Median 3Q Max -8.655 -1.871 -0.655 1.750 13.155 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 66.1485 7.3839 8.958 6.12e-15 *** TAILLE 0.4053 0.1230 3.296 0.0013 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 3.574 on 117 degrees of freedom (2 observations deleted due to missingness) Multiple R-squared: 0.08496, Adjusted R-squared: 0.07714 F-statistic: 10.86 on 1 and 117 DF, p-value: 0.001299 > base$Pred70s.POITRINE=predict(subreg,newdata=base) > base TAILLE Pred.POITRINE Moy.POITRINE Pred70s.POITRINE 1 50 88.50935 89.61538 86.41274 2 55 89.09248 89.61538 88.43917 3 60 89.67562 89.61538 90.46559 4 65 90.25876 89.61538 92.49201 5 70 90.84189 89.61538 94.51844
que l’on peut visualiser ci-dessous,
L’explication est que la morphologie des femmes (ou en tous les cas de celles qui ont posé nues dans Playboy) a beaucoup changé en 60 ans, comme évoquéprécédemment. On peut le voir en regardant l’évolution de la pente, par date (5ans)
(le trait horizontal est la pente obtenue sur l’ensemble de la période). On peut d’ailleurs intégrer l’année dans notre modèle, en rajoutant le temps comme variable explicative dans la régression,
> reg=lm(POITRINE~TAILLE+DATE,data=mensurations) > summary(reg) Call: lm(formula = POITRINE ~ TAILLE + DATE, data = mensurations) Residuals: Min 1Q Median 3Q Max -9.7304 -2.2184 -0.5405 2.0742 12.1041 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 368.296266 18.008455 20.451 < 2e-16 *** TAILLE 0.323319 0.042140 7.673 6.39e-14 *** DATE -0.150301 0.009393 -16.002 < 2e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 3.491 on 634 degrees of freedom (22 observations deleted due to missingness) Multiple R-squared: 0.2944, Adjusted R-squared: 0.2922 F-statistic: 132.3 on 2 and 634 DF, p-value: < 2.2e-16
Visuellement, on aurait la prédiction suivante,
 |
 |
Si on se place autour de 1970,
> submensurations=mensurations[abs(DATE-1970)<=5,] > subreg=lm(POITRINE~TAILLE,data=submensurations) > predict(subreg,newdata=data.frame(TAILLE=c(55,70))) 1 2 89.14361 95.75476
alors qu’autour de 1990,
> submensurations=mensurations[abs(DATE-1990)<=5,] > subreg=lm(POITRINE~TAILLE,data=submensurations) > predict(subreg,newdata=data.frame(TAILLE=c(55,70))) 1 2 87.69853 92.70228
Notons qu’on fait ici des régressions locales (dans ces deux derniers calculs), ce qui n’est pas ce qui est représenté visuellement… Si on veut une meilleure modélisation locale, il faudrait quitter le monde des modèles linéaires pour aller un peu plus loin. On verra en fin de cours si on a le temps d’en parler
 |
 |
Cela dit, globalement, on avait un effet positif, mais il aurait été tout à fait possible d’avoir un effet (global négatif). C’est ce qui arrive dans le paradoxe de Simpson.