Graphique de régression complexe dans R

10

J'ai besoin de dessiner un graphique complexe pour l'analyse visuelle des données. J'ai 2 variables et un grand nombre de cas (> 1000). Par exemple (le nombre est 100 si la dispersion est moins "normale"):

x <- rnorm(100,mean=95,sd=50)
y <- rnorm(100,mean=35,sd=20)
d <- data.frame(x=x,y=y)

1) J'ai besoin de tracer des données brutes avec une taille en points, correspondant à la fréquence relative des coïncidences, ce plot(x,y)n'est donc pas une option - j'ai besoin de tailles en points. Que faut-il faire pour y parvenir?

2) Sur le même tracé, je dois tracer une ellipse d'intervalle de confiance à 95% et une ligne représentant le changement de corrélation (je ne sais pas comment le nommer correctement) - quelque chose comme ceci:

library(corrgram)
corrgram(d, order=TRUE, lower.panel=panel.ellipse, upper.panel=panel.pts)

corrélogramme

mais avec les deux graphiques sur une seule parcelle.

3) Enfin, je dois dessiner un modèle de régression linéaire résultant en plus de tout cela:

r<-lm(y~x, data=d)
abline(r,col=2,lwd=2)

mais avec une plage d'erreur ... quelque chose comme sur QQ-plot:

QQ-plot

mais pour les erreurs d'ajustement, si cela est possible.

La question est donc:

Comment réaliser tout cela en un seul graphique?

Yuriy Petrovskiy
la source

Réponses:

29

L'image ci-dessous ressemble-t-elle à ce que vous souhaitez réaliser?

entrez la description de l'image ici

Voici le code R mis à jour , suite à vos commentaires:

do.it <- function(df, type="confidence", ...) {
  require(ellipse)
  lm0 <- lm(y ~ x, data=df)
  xc <- with(df, xyTable(x, y))
  df.new <- data.frame(x=seq(min(df$x), max(df$x), 0.1))
  pred.ulb <- predict(lm0, df.new, interval=type)
  pred.lo <- predict(loess(y ~ x, data=df), df.new)
  plot(xc$x, xc$y, cex=xc$number*2/3, xlab="x", ylab="y", ...)
  abline(lm0, col="red")
  lines(df.new$x, pred.lo, col="green", lwd=1.5)
  lines(df.new$x, pred.ulb[,"lwr"], lty=2, col="red")
  lines(df.new$x, pred.ulb[,"upr"], lty=2, col="red")    
  lines(ellipse(cor(df$x, df$y), scale=c(sd(df$x),sd(df$y)), 
        centre=c(mean(df$x),mean(df$y))), lwd=1.5, col="green")
  invisible(lm0)
}

set.seed(101)
n <- 1000
x <- rnorm(n, mean=2)
y <- 1.5 + 0.4*x + rnorm(n)
df <- data.frame(x=x, y=y)

# take a bootstrap sample
df <- df[sample(nrow(df), nrow(df), rep=TRUE),]

do.it(df, pch=19, col=rgb(0,0,.7,.5))

Et voici la version ggplotized

entrez la description de l'image ici

produit avec le code suivant:

xc <- with(df, xyTable(x, y))
df2 <- cbind.data.frame(x=xc$x, y=xc$y, n=xc$number)
df.ell <- as.data.frame(with(df, ellipse(cor(x, y), 
                                         scale=c(sd(x),sd(y)), 
                                         centre=c(mean(x),mean(y)))))
library(ggplot2)

ggplot(data=df2, aes(x=x, y=y)) + 
  geom_point(aes(size=n), alpha=.6) + 
  stat_smooth(data=df, method="loess", se=FALSE, color="green") + 
  stat_smooth(data=df, method="lm") +
  geom_path(data=df.ell, colour="green", size=1.2)

Il pourrait être personnalisé un peu plus en ajoutant des indices d'ajustement de modèle, comme la distance de Cook, avec un effet d'ombrage des couleurs.

chl
la source
1
@chl +1, joli graphique et code court.
mpiktas
@mpiktas Merci. Cela m'a amené à réaliser que je n'avais pas travaillé avec le bon échantillon, en fait :-)
chl
df.new <- data.frame(x = seq(min(x), max(x), 0.1))s size is also strange (too small). Also tryed X,Flibrary(car) cr.plots(m0)
(X,y)car::dataEllipseellipse
2
@Tal L'interprétation de l'ellipse est la même que dans le corrgrampackage: elle montre 95% de région de confiance par paire en supposant une distribution normale bivariée centrée sur la moyenne et mise à l'échelle par SD (x) et SD (y). Je ne suis pas un grand fan de cela lorsqu'il est utilisé dans un nuage de points, cependant. Mais voir Murdoch & Chow, A graphical display of large correlation matrices , Am Stat (1996) 50: 178, ou Friendly, Corrgrams: Exploratory displays for correlation matrices , Am Stat (2002) 56: 316.
chl
2

Pour le point 1, utilisez simplement le cexparamètre sur le tracé pour définir la taille du point.

Par exemple

x = rnorm(100)
plot(x, pch=20, cex=abs(x))

Pour avoir plusieurs graphiques dans un même tracé, utilisez par(mfrow=c(numrows, numcols))une disposition régulièrement espacée ou layoutcréez des graphiques plus complexes.

Nico
la source
1
+1 pour l'astuce cex, mais je pense que l'OP veut que tout soit sur la même région de traçage, pas sur des régions distinctes.
chl
Ahh ... maintenant je comprends la question. Eh bien, alors il peut simplement utiliser curveou pointspour superposer les trois graphiques;)
nico