Reprenons ici le code pour récupérer la base de la taille (et du poids) d’un ensemble de 200 personnes.
Davis=read.table( "http://socserv.socsci.mcmaster.ca/jfox/Books/Applied-Regression-2E/datasets/Davis.txt") Davis[12,c(2,3)]=Davis[12,c(3,2)] X=Davis$height n=length(X)
La première chose que l’on avait vu la semaine passée est le PP-plot,
On va ici ‘tester’ une loi normale (de même moyenne et de même variance que notre échantillon… faisons simple)
PP=data.frame(empirique=(1:n)/n, theorique=pnorm(sort(X),mean(X),sd(X)))
Comme dans tous les posts, je mets de lignes de code basées sur
le package ggplot2 (alors qu’en cours, j’utilise les fonctions de
bases de graphiques). Le code le plus simple pour tracer le PP plot est
library(ggplot2) qplot(empirique,theorique,data=PP)
ou, pour un code qui pourra plus facilement être modifié (pour changer les couleurs, les titres, ou rajouter des courbes)
ggplot(PP,aes(x=empirique,y=theorique)) + geom_abline(intercept=0,slope=1,colour="red") + geom_point()
La seconde chose que l’on avait vu la semaine passée est le QQ-plot,
Le point de départ est la table, avec les coordonnées, abscisse et ordonnnée,
QQ=data.frame(empirique=sort(X), theorique=qnorm((1:n)/n,mean(X),sd(X)))
puis on fait un graphique
ggplot(QQ,aes(x=empirique,y=theorique)) + geom_abline(intercept=0,slope=1,colour="red") + geom_point()
Ici, on utilise le quantile de la loi normale, de même moyenne et même variance que notre échantillon. Mais comme on l’avait vu en cours, on peut en fait mettre une loi normale quelconque (disons centrée et réduite), et on a un test de normalité,
QQ=data.frame(empirique=sort(X), theorique=qnorm((1:n)/n))
que l’on visualise avec
ggplot(QQ,aes(x=theorique,y=empirique)) + geom_abline(intercept = mean(X), slope = sd(X),colour="red") + geom_point()
Comme je l’avais dit en cours, on a ici une méthode simple d’ajustement de loi, pour tester
contre
essentiellement visuelle.
Si on veut aller un peu plus loin, en particulier pour avoir un
intervalle de confiance, on peut faire des simulations. Par exemple,
si la loi que l’on cherche est une loi .
PPs=NULL for(s in 1:100){PPs=rbind(PPs,data.frame("empirique"=(1:n)/n,"theorique"=pnorm(sort(rnorm(n))),sample=s))} ggplot(PPs,aes(x=empirique,y=theorique,group=sample))+ geom_line(colour="light blue")+geom_abline(intercept = 0, slope = 1,colour="red")
pour le PP plot, et pour le QQ plot
QQs=NULL for(s in 1:100){QQs=rbind(QQs,data.frame("empirique"=sort(rnorm(n)),"theorique"=qnorm((1:n)/n),sample=s))} ggplot(QQs,aes(x=empirique,y=theorique,group=sample))+ geom_line(colour="light blue")+ geom_abline(intercept = 0, slope = 1,colour="red")
Dans la second partie du cours, on avait parlé d’estimation à noyau de la densité.Classiquement, tout le monde connaît l’histogramme.
dX=data.frame(X=X) H=ggplot(dX,aes(X)) H+geom_histogram(aes(y=..density..),fill="light blue")
On peut considérer plusieurs formes de noyaux. Un noyau rectangulaire
qui donnera un histogramme glissant. Le code est ici
H+geom_density(kernel="rectangular",fill="light blue",colour=NA)
On peut considérer un noyau triangulaire
H+geom_density(kernel="triangular",fill="light blue",colour=NA)
un noyau d’Epanechnikov
H+geom_density(kernel="epanechnikov",fill="light blue",colour=NA)
ou un noyau Gaussien
H+geom_density(kernel="gaussian",fill="light blue",colour=NA)
Ces deux derniers donnant des visualisation très très proches. On peut ensuite jouer sur la fenêtre. Pour l’instant, on avait
H+geom_density(kernel="gaussian",adjust=1,fill="light blue",colour=NA)
mais on peut prendre une fenêtre plus petite
H+geom_density(kernel="gaussian",adjust=.25,fill="light blue",colour=NA)
Ce qui va augmenter la variance, mais baisser le biais (on va en reparler cette semaine, promis).
OpenEdition suggests that you cite this post as follows:
Arthur Charpentier (September 30, 2014). PP plot, QQ plot et estimation à noyau. Freakonometrics. Retrieved November 5, 2024 from https://doi.org/10.58079/oux0