Régression: données individuelles ou temporelles ?

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 https://latex.codecogs.com/gif.latex?p-value de l’ordre de 1.5%, on va rejeter l’hypothèse https://latex.codecogs.com/gif.latex?H_0 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 https://latex.codecogs.com/gif.latex?p-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 (https://latex.codecogs.com/gif.latex?\pm5ans)

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

La recherche, l’enseignement et les élèves

Depuis quelques temps, lorsqu’on me demande ce que je fais dans la vie, je réponds “enseignant-chercheur” et je trouve que le terme me convient bien…

I.− Part. prés. de enseigner*: a) [L’obj. désigne le destinataire] Transmettre un savoir à. Synon.

Prononc. et Orth. : [ʃ ε ʀ ʃ œ:ʀ], fém. [-ø:z]. Ds Ac. 1694-1932. Étymol. et Hist. 1. le plus souvent avec compl. introduit par de a) 1538 cercheur « personne qui cherche, qui s’informe » (Est., s.v. conquisitores); 1636 chercheur (Monet); b) 1840 adj. « qui cherche (qqc.) » (Proudhon, Qu’est-ce que la propriété, p. 133); 2. emploi abs. 1835 adj. chercheuse « avide de connaissances, qui essaie de découvrir » (Balzac, loc. cit.);

Si on regarde les définitions, on voit qu’un chercheur s’informe (et reçoit du savoir) alors que l’enseignant le transmet (et communique son savoir). D’aucuns prétendront que le chercheur aussi transmet, mais quand on voit l’impact des articles académiques, on peut se demander à qui on transmet. Sans chercher à polémiquer, je trouve d’ailleurs surprenant qu’en France on ait pu séparer les deux, en créant des postes de chercheurs qui n’étaient pas tenu de partager leur savoir (mais j’ouvre une parenthèse que je ne tiens pas à ouvrir davantage)…. Bref, je suis très content d’enseigner car c’est cette activité qui me donne une utilité sociale.
Sauf que comme tous les enseignants, je déteste faire des examens. Pire encore, les corriger…
Lors de mon dernier examen, j’avais fait un questionnaire, avec des choix multiples (4 réponses, une seule bonne). Ça demande un temps fou à préparer, mais la correction est beaucoup plus rapide que lire des dizaines de pages (activité que je me laisse pour les heures à venir, car en plus des questionnaires, il y avait deux gros devoirs sur des bases de données). Et en plus, avec les choix multiples la correction est impartiale !
Dans un monde idéal, j’aimerais pouvoir dire “si vous pensez avoir réussi, passez au cours suivant, sinon reprenez le…” Alors j’ai tenté l’expérience. J’ai demandé aux étudiants d’écrire sur la copie le nombre de bonnes réponses qu’ils pensaient avoir. L’exercice était libre, et sans influence sur la note finale. Enfin, presque. Suite à une suggestion, j’ai accepté de donner des points à toute personne qui prédirait exactement son score !

> ACT6420=read.table("ACT6420-H2012.csv",
> base=ACT6420[ACT6420$obs>0,]
> attach(base)

Tout d’abord, comme l’exercice n’était pas obligatoire, certains (plusieurs) n’ont pas fait de prédiction… et comme souvent, la “non-réponse” est (très) informative,

> reponse=(is.na(pred)==FALSE)
> boxplot(obs~reponse,horizontal=TRUE)

Ceux qui n’ont pas répondu ont moins bien réussi que ceux qui ont répondu. Et sans vouloir revenir sur un point vu en cours (corrélation n’est pas causalité), je suis toujours surpris de fait que les personnes qui ne répondent pas aient un comportement aussi clairement différent de ceux qui répondent.

Maintenant regardons ce qui a été répondu par ceux qui ont joué le jeu, et qui ont fait une prédiction,

> plot(pred,obs,xlim=c(0,40),ylim=c(0,40),
+ xlab="Nombre de bonne réponses prédit",
+ ylab="Nombre de bonne réponses observé")
> abline(a=0,b=1,lty=1,col="red")
> abline(lm(obs~pred),lty=2)

En abscisse la prédiction, et en ordonnée, la valeur observée (un peu comme dans le cours de prévision justement). La courbe en rouge est une prédiction parfaite. Et la courbe en pointillé, la droite de régression. Les meilleurs se sous-estiment (en moyenne), alors que les moins bons se surestiment. Comme le note @AlbertineP on voit que les moins bons sont en fait un peu meilleurs qu’ils ne le pensent, alors que les meilleurs se surestiment ! Un effet Dunning-Kruger inversé, pour reprendre la terminologie proposée par @Tournyol. Damned, je ne pourrais donc pas demander aux étudiants de s’auto-évaluer, il va falloir que je continue à corriger mes copies ! (et je vais renouveler l’expérience, pour voir sur plus de copies ce que ça donne !)

ACT6420 compléments de lecture

Histoire de revoir, d’approfondir certaines notions abordées dans la première partie du cours, quelques liens vers de documents, ici ou là… Parmi les ouvrages cités dans le plan de cours, il y a le site de Jed Frees qui propose des codes R et des bases de données, utilisés dans le livre. De même pour le livre de William Greene, un site permet d’avoir accès aux données par exemple. Enfin, Russell Davidson et James MacKinnon ont mis en ligne l’intégralité de leur livre d’estimation et d’inférence en économétrie. Bien que le livre soit plus poussé que ce que l’on devrait faire, il apportera de nombreux éclaircissements, et une autre manière de voir le cours. Je peux aussi renvoyer aux notes de cours de Bruno Crepon et de Nicolas Jacquemet, d’économétrie linéaire appliquée.

  • statistique mathématique
La section 3. des notes de cours de Didier Concordet rappelle les principaux résultats à connaître sur l’inférence statistique, ou sinon, je peux renvoyer vers les notes de Jérôme Dupuis.
  • tests statistiques
Une référence un peu générale, le livre de Gaël Millot, en ligne sur Google Books . Sinon, les sections 4. et 5. des notes de cours de Didier Concordet reprennent l’essentiel des résultats à connaître sur les tests pour suivre le cours. Je renverrais aussi vers les transparents de Wissem Maazoun sur les tests usuels. Enfin, sur la théorie des tests, des compléments peuvent être trouvés dans les notes de cours écrites par Paul Doukhan, ou les notes de Nicolas Chopin, Emmanuelle Cretois et Jean Diebolt. (ces deux documents vont beaucoup plus loin que ce qui est demandé comme prérequis).
  • optimisation et aspects numériques
On reviendra plus longuement sur ces aspects là sur la section de séries temporelles, mais les transparents de Grégory Vial peuvent être un bon rappel (en particulier sur la descente de gradient). Je peux aussi mentionner une revue d’algorithme dans un document de Joël Pothier (qui va beaucoup plus loin que ce qui pourrait servir dans ce cours, avec page 10 la descente de gradient).
  • calcul différentiel et calcul matriciel
Sur les rappels d’algèbre matricielle, je peux renvoyer à un court document synthétique de Xavier Joutard. Et sur le calcul différentiel appliqué à des matrices (que l’on utilisera lorsque l’on abordera la régression multiple) je renvoie à la première section du document de Thomas Minka.
  • espérance conditionnelle et géométrie
Pour des rappels sur l’espérance conditionnelle, et l’interprétation géométrique en terme de projection, je revois vers des transparents de Samy Tindel, qui reprennent tout ce qu’il faut savoir. Sinon Massimiliano Gubinelli a mis en ligne un document de 4 pages, un cours document de Vincent Bansaye, ou la Section 1.4. des notes de cours d’Arnaud Guyader permettent d’apporter des éclaircissements.
  • un peu de programmation (en R)
Quelques liens vers des livres de R, d’introduction à la programmation, dont l’introduction d’Emmanuel Paradis, le brise glace d’Andrew Robinson et Arnaud Schloesing, ou le livre de Vincent Goulet. Bien qu’aucun prérequis en informatique ne soit demandé pour le cours, je m’attends à ce que l’on connaissent les principes généraux de programmation. Par exemple connaitre la différence entre une variable locale et une variable globale dans un code.
> y=cars$dist 
> f=function(z) {return(mean(log(z)))} 
> meanlog=f(y)
D’autres compléments seront mis en ligne, en particulier sur les processus, dès que nous aborderons les séries temporelles.

J’ai probablement raté des documents intéressants, les commentaires sont ouverts si quelqu’un souhaite partager des documents…

Actuariat 1, ACT2121, premier cours

Le premier cours de l’automne d’ACT2121 actuariat I  aura lieu lundi après midi, au local SH3220. Comme il s’agit de préparer l’examen P de la SOA, nous allons passer les séances à faire des exercices. Le plan de cours sera mis en ligne prochainement. L’ouvrage de référence est le livre de Jacques Labelle, exercices de probabilités appliquées. A chaque séance, je proposerai une série de 90 exercices, et on en fera autant que possible (dans les conditions de l’examen, avec préparation de 6 minutes par exercices). Pour la première séance, les exercices sont en ligne,

Je rajouterai au fur et à mesure des compléments de cours, voire des corrections…