M2 Pro Ingénierie Mathématique Université d'Angers, Université de Nantes
Année 2011-2012
TD de régression linéaire simple
Exercice 1 : EMV et EMC On dispose d'un échantillon de n couples (xi , yi ), i = 1, . . . , n satisfaisant y i = β0 + β1 x i + i ,
où l'on suppose que les i , i = 1 . . . n sont des variables aléatoires i.i.d. de loi normale N (0, σ 2 ), σ 2 inconnue. Le but de cet exercice est de comparer la méthode des moindres carrés et la méthode du maximum de vraisemblance dans ce modèle. 1. Calculer la vraisemblance de l'échantillon. 2. Calculer les estimateurs de β0 , β1 et σ 2 par la méthode du maximum de vraisemblance. 3. Calculer les estimateurs de β0 , β1 et σ 2 à l'aide de la méthode des moindres carrés. 4. Comparer les résultats obtenus. Exercice 2 : Modèle de croissance humaine Un père a deux garçons, et s'inquiète de la croissance de son cadet qu'il trouve petit. Il décide de faire un modèle familial à partir des mesures de taille en fonction de l'age de l'aîné :
âge taille
3 96
4 104.8
5 110.3
6 115.3
7 121.9
8 127.4
9 130.8
10 136
11 139.7
12 144.5
1. Représenter les données sur un graphique et justier l'utilisation d'un modèle de régression linéaire simple. Discuter les hypothèses nécessaires. 2. Estimer les coecients de la régression et tracez sur le graphique la droite de régression estimée. 3. Calculer le R2 et représenter les résidus. La régression semble-t'elle valable ? Pour information, les données proviennent des études auxologique du Docteur Sempé dont une partie a été publiée par Abidi et al (1996). Ces données mesurées sur des milliers d'enfants (de 1 mois à 19 ans) ont permis d'établir un modèle de croissance humaine qui fournit les prédictions du carnet de santé. Il s'écrit de la manière suivante : Y = θ1 1 −
1 + ((X + θ8 )/θ2 )θ3
1 , + ((X + θ8 )/θ4 )θ5 + ((x + θ8 )/θ6 )θ7
où θ1 représente la taille adulte, θ8 le temps de grossesse, et les couples (θ2 , θ3 ), (θ4 , θ5 ) et (θ6 , θ7 ) permettent de modéliser respectivement la phase de croissance initiale (juste après la naissance), la phase de croissance centrale (pré-adolescente) et la phase nale. Exercice 3 : Hauteur des arbres Nous souhaitons exprimer la hauteur Y d'un arbre en fonction de son diamètre X à 1m30 du sol. Pour cela, nous avons mesuré 20 couples diamètre-hauteur et les résultats ci-dessous sont disponibles : 20
x ¯ = 34.9; 20
1 X (xi − x ¯)2 = 28.29; y¯ = 18.34 20 i=1 20
1 X 1 X (yi − y¯)2 = 2.85; (xi − x ¯)(yi − y¯) = 6.26. 20 i=1 20 i=1
1
1. On note Yˆ = βˆ0 + βˆ1 X l'estimation de la droite de régression. Donner l'expression de βˆ0 et βˆ1 en fonction des statistiques élémentaires ci-dessus. Calculer βˆ0 et βˆ1 . 2. Donner une mesure de qualité d'ajustement des données au modèle. Exprimer cette mesure à l'aide des statistiques élémentaires. Calculer et commenter. 3. Testez H0 : ”βj = 0” contre H1 : ”βj 6= 0” pour j = 0, 1. Commentez. Exercice 4 : Natalité en Amérique La tableau suivant contient la liste de 14 pays d'Amérique du Nord et d'Amérique Centrale, dont la population dépassait le million d'habitants en 1985. Pour chaque pays, on mesure le taux de natalité yi (nombre de naissances annuel pour 1000 habitants) ainsi que le taux d'urbanisation xi (pourcentage de la population vivant dans des villes de plus de 100000 habitants). On fait l'hypothèse d'un modèle de regréssion linéaire simple du type yi = β0 + β1 xi + i , c'est-à-dire que le taux de natalité dépend linéairement du taux d'urbanisation.
Observations 1 2 3 4 5 6 7 8 9 10 11 12 13 14
pays Canada Costa-Rica Cuba USA El Salvador Guatemala Haïti Honduras Jamaïque Mexique Nicaragua Trinidade/Tobago Panama Rep. Dominicaine
taux d'urbanisation 55.0 27.3 33.3 56.5 11.5 14.2 13.9 19.0 33.1 43.2 28.5 6.8 37.7 37.1
1. 2. 3. 4. 5.
taux de natalité 16.2 30.5 16.9 16.0 40.2 38.4 41.3 43.9 28.3 33.9 44.2 24.6 28.0 33.1
Représenter graphiquement les données. Estimer les paramètres β0 et β1 du modèle et tracer la droite de régression correspondante. Calculer la somme des résidus. Calculer SCtot , SCreg et SCres puis R2 . Tester l'hypothèse H0 : ”β1 = 0” contre H1 : ”β1 6= 0” et donner un intervalle de conance à 95% pour β1 . 6. Tester l'hypothèse H0 : ”β0 = 0” contre H1 : ”β0 6= 0” et donner un intervalle de conance à 95% pour β0 . 7. Représenter graphiquement un intervalle de conance de 95% autour de la droite de régression à l'aide d'une grille de 10 points.
2
M2 Pro Ingénierie Mathématique Université d'Angers, Université de Nantes
Année 2011-2012
TP de régression linéaire simple
Exercice 1 : concentration en ozone Nous allons traiter les 50 données journalières de la concentration en ozone en fonction de la température. Les données se trouvent dans le chier "ozone.txt". La variable à expliquer est la concentration en ozone, notée "maxO3", et la variable explicative est la température à midi, notée "T12". 1. Commencer par représenter les données à l'aide des commandes suivantes : >ozone<-read.table("ozone.txt",header=T) >plot(maxO3~T12,data=ozone)
Une regression linéaire simple semble-t'elle justiée graphiquement ? 2. Eectuer la régression linéaire à l'aide de la commande >reg<-lm(maxO3~T12,data=ozone)
et consulter les résultats à l'aide de la commande >resume<-summary(reg)
Que représente les coecients de la matrice coecients ? 3. Tracer l'estimation de la droite de régression, ainsi qu'un intervalle de conance à 95% de celle-ci grâce aux commandes suivantes : >plot(maxO3~T12,data=ozone) >T12=seq(min(ozone[,"T12"]),max(ozone[,"T12"]),length=100) >grille<-data.frame(T12) >ICdte<-predict(reg,new=grille,interval="confidence",level=0.95) >matlines(grille$T12,cbind(ICdte),lty=c(1,2,2),col=1)
Ce graphique permet de vérier visuellement l'ajustement des données au modèle de régression proposé. Que remarquez-vous ? Représentez le vecteur des résidus grâce aux commandes : >res<-rstudent(reg) >plot(res,pch=15,ylab=Résidus,ylim=c(-3,3)) >abline(h=c(-2,0,2),lty=c(2,1,2)).
4. On s'intéresse à présent à la qualité de prévision du modèle. Pour cela, on va tracer un intervalle de conance des prévisions de la manière suivante : >plot(maxO3~T12,data=ozone) >T12=seq(min(ozone[,"T12"]),max(ozone[,"T12"]),length=100) >grille<-data.frame(T12) >ICprev<-predict(reg,new=grille,interval="pred",level=0.95) >matlines(grille$T12,cbind(ICprev),lty=c(1,2,2),col=1)
5. On va maintenant calculer les intervalles de conances des coecients β0 et β1 du modèle de régression. Pour cela, on utilise la fonction coef() qui permet d'extraire les estimateurs de β0 et β1 et leurs écarts types empiriques. >seuil<-qt(0.975,df=reg$df.res) >beta0min<-coef(resume)[1,1]-seuil*coef(resume)[1,2] >beta0max<-coef(resume)[1,1]+seuil*coef(resume)[1,2] >beta1min<-coef(resume)[2,1]-seuil*coef(resume)[2,2] >beta1max<-coef(resume)[2,1]+seuil*coef(resume)[2,2]
3
Que remarquez-vous sur l'intervalle de conance de β0 ? Comment l'expliquez-vous ? 6. Pour être plus précis et tenir compte de la dépendance entre β0 et β1 , on peut aussi construire une région de conance pour βˆ. Les commandes suivantes permettent de visualiser la diérence entre le rectangle de conance, simple juxtaposition des deux intervalles de conance et la région de conance. Elles nécessitent l'installation du package ellipse. >library(ellipse) >plot(ellipse(reg,level=0.95),type="l",xlab="beta0",ylab="beta1") >points(coef(reg)[1],coef(reg)[2],pch=3) >lines(c(beta0min,beta0min,beta0max,beta0max,beta0min),c(beta1min,beta1max,beta1max, beta1min,beta1min),lty=2) >plot(ellipse(reg,level=0.95),type="l",xlab="beta0",ylab="beta1") >points(coef(reg)[1],coef(reg)[2],pch=3) >lines(c(beta0min,beta0min,beta0max,beta0max,beta0min),c(beta1min,beta1max,beta1max, beta1min,beta1min),lty=2)
Exercice 2 : Hauteur des eucalyptus On veut expliquer la hauteur des eucalyptus en fonction de leur circonférence à partir d'une régression linéaire simple. On dispose de 1737 couples circonférence-hauteur qui se trouvent dans le chier "eucalyptus.txt". 1. Extraire et représenter les données dans le plan. 2. Eectuer la régression et commenter les résultats obtenus. 3. Tracer l'estimation de la droite de régression et un intervalle de conance à 95% de celle-ci. Que déduisez-vous de la qualité de l'estimation ? 4. Calculer les intervalles de conance des coecients β0 et β1 du modèle de régression et tracer le rectangle de conance associé. Faites de même pour la région de conance du couple β = (β0 , β1 ). Commenter. 5. On veut à présent prédire la taille d'une nouvelle série d'eucalyptus de circonférence 50, 100, 200 puis 500. Donner les estimateurs de la taille de chacun d'entre eux et les intervalles de conances associés. 6. Que se passe-t'il pour les faibles valeurs de circonférences ? Proposer une amélioration possible du modèle pour tenir compte de ce phénomène. Cette amélioration sera traitée dans le prochain TP de régression multiple. Exercice 3 : Modèle quadratique Au vu de la représentation de la concentration d'ozone en fonction de la température à midi de l'Exercice 1, nous souhaitons modéliser l'ozone par la température au carré. 1. Ecrire le modèle et estimer les paramètres. 2. Comparer ce modèle au modèle de régression linéaire.
4