R

qqnorm et qqline dans ggplot2

Qqnorm Qqline Ggplot2



Solution:

Le code suivant vous donnera le tracé que vous voulez. Le package ggplot ne semble pas contenir de code pour calculer les paramètres de la qqline, donc je ne sais pas s'il est possible de réaliser un tel tracé dans un one-liner (compréhensible).

qqplot.data<- function (vec) # argument: vector of numbers { # following four lines from base R's qqline() y <- quantile(vec[!is.na(vec)], c(0.25, 0.75)) x <- qnorm(c(0.25, 0.75)) slope <- diff(y)/diff(x) int <- y[1L] - slope * x[1L] d <- data.frame(resids = vec) ggplot(d, aes(sample = resids)) + stat_qq() + geom_abline(slope = slope, intercept = int) }   

Vous pouvez également ajouter des intervalles de confiance/bandes de confiance avec cette fonction (Parties du code copié à partir devoiture:::qqPlot)



gg_qq<- function(x, distribution = 'norm', ..., line.estimate = NULL, conf = 0.95, labels = names(x)){ q.function <- eval(parse(text = paste0('q', distribution))) d.function <- eval(parse(text = paste0('d', distribution))) x <- na.omit(x) ord <- order(x) n <- length(x) P <- ppoints(length(x)) df <- data.frame(ord.x = x[ord], z = q.function(P, ...)) if(is.null(line.estimate)){ Q.x <- quantile(df$ord.x, c(0.25, 0.75)) Q.z <- q.function(c(0.25, 0.75), ...) b <- diff(Q.x)/diff(Q.z) coef <- c(Q.x[1] - b * Q.z[1], b) } else { coef <- coef(line.estimate(ord.x ~ z)) } zz <- qnorm(1 - (1 - conf)/2) SE <- (coef[2]/d.function(df$z)) * sqrt(P * (1 - P)/n) fit.value <- coef[1] + coef[2] * df$z df$upper <- fit.value + zz * SE df$lower <- fit.value - zz * SE if(!is.null(labels)) df$ord.x < df$lower, labels[ord],'') p <- ggplot(df, aes(x=z, y=ord.x)) + geom_point() + geom_abline(intercept = coef[1], slope = coef[2]) + geom_ribbon(aes(ymin = lower, ymax = upper), alpha=0.2) if(!is.null(labels)) p <- p + geom_text( aes(label = label)) print(p) coef }  

Exemple:

Animaux2<- data(Animals2, package = 'robustbase') mod.lm <- lm(log(Animals2$brain) ~ log(Animals2$body)) x <- rstudent(mod.lm) gg_qq(x)  

entrez la description de l


Depuis la version 3.0, unstat_qq_line équivalent à ce qui suit a trouvé sa place dans la version officiellecode ggplot2.


Ancienne version :

Depuis la version 2.0, ggplot2 dispose d'une interface d'extension bien documentée ; nous pouvons donc maintenant facilement écrire une nouvelle statistique pour la qqline par elle-même (ce que j'ai fait pour la première fois, donc les améliorations sont les bienvenues):

qq.ligne<- function(data, qf, na.rm) { # from stackoverflow.com/a/4357932/1346276 q.sample <- quantile(data, c(0.25, 0.75), na.rm = na.rm) q.theory <- qf(c(0.25, 0.75)) slope <- diff(q.sample) / diff(q.theory) intercept <- q.sample[1] - slope * q.theory[1] list(slope = slope, intercept = intercept) } StatQQLine <- ggproto('StatQQLine', Stat, # http://docs.ggplot2.org/current/vignettes/extending-ggplot2.html # https://github.com/hadley/ggplot2/blob/master/R/stat-qq.r required_aes = c('sample'), compute_group = function(data, scales, distribution = stats::qnorm, dparams = list(), na.rm = FALSE) { qf <- function(p) do.call(distribution, c(list(p = p), dparams)) n <- length(data$sample) theoretical <- qf(stats::ppoints(n)) qq <- qq.line(data$sample, qf = qf, na.rm = na.rm) line <- qq$intercept + theoretical * qq$slope data.frame(x = theoretical, y = line) } ) stat_qqline <- function(mapping = NULL, data = NULL, geom = 'line', position = 'identity', ..., distribution = stats::qnorm, dparams = list(), na.rm = FALSE, show.legend = NA, inherit.aes = TRUE) { layer(stat = StatQQLine, data = data, mapping = mapping, geom = geom, position = position, show.legend = show.legend, inherit.aes = inherit.aes, params = list(distribution = distribution, dparams = dparams, na.rm = na.rm, ...)) }  

Cela se généralise également sur la distribution (exactement commestat_qq le fait) et peut être utilisé comme suit :

> test.data test.data.2 ggplot(test.data, aes(sample=sample)) + stat_qq() + stat_qqline() > ggplot(test.data.2, aes(sample=sample)) + stat_qq(distribution =qt, dparams=list(df=2)) + + stat_qqline(distribution=qt, dparams=list(df=2))

(Malheureusement, étant donné que qqline se trouve sur un calque séparé, je n'ai pas trouvé de moyen de ' réutiliser ' les paramètres de distribution, mais cela ne devrait être qu'un problème mineur.)