qqnorm y qqline en ggplot2

Question

Say tiene un modelo lineal LM que quiero una parcela qq de los residuos. Normalmente me gustaría utilizar los gráficos de base R:

qqnorm(residuals(LM), ylab="Residuals")
qqline(residuals(LM))

Me puede encontrar la manera de obtener la parte qqnorm de la trama, pero me parece que no puede gestionar el qqline:

ggplot(LM, aes(sample=.resid)) +
    stat_qq()

sospecho que me falta algo bastante básico, pero parece que tendría que haber una manera fácil de hacer esto.

EDIT: Muchas gracias por la solución a continuación. He modificado el código (muy ligeramente) para extraer la información del modelo lineal, de modo que la trama funciona como la trama de conveniencia en el paquete de gráficos de base R.

ggQQ <- function(LM) # argument: a linear model
{
    y <- quantile(LM$resid[!is.na(LM$resid)], c(0.25, 0.75))
    x <- qnorm(c(0.25, 0.75))
    slope <- diff(y)/diff(x)
    int <- y[1L] - slope * x[1L]
    p <- ggplot(LM, aes(sample=.resid)) +
        stat_qq(alpha = 0.5) +
        geom_abline(slope = slope, intercept = int, color="blue")

    return(p)
}

Answer 1

El siguiente código le dará la trama que desea. No parece el paquete ggplot a contener código para el cálculo de los parámetros de la qqline, así que no sé si es posible lograr una gráfica de este tipo en un (comprensible) de una sola línea.

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)

}

Answer 2

También puede agregar intervalos de confianza / bandas de confianza con esta función (Partes del código copiado de car:::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$label <- ifelse(df$ord.x > df$upper | 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
}

Ejemplo:

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

Answer 3

La Q-Q estándar de diagnóstico para modelos gráficos lineales cuantiles de los estandarizado residuos vs. cuantiles teóricos de N (0,1). @ GgQQ función de Peter gráficamente los residuos. El siguiente fragmento de código que modifica y añade algunos cambios cosméticos para hacer la trama más como la que se obtendría a partir plot(lm(...)).

ggQQ = function(lm) {
  # extract standardized residuals from the fit
  d <- data.frame(std.resid = rstandard(lm))
  # calculate 1Q/4Q line
  y <- quantile(d$std.resid[!is.na(d$std.resid)], c(0.25, 0.75))
  x <- qnorm(c(0.25, 0.75))
  slope <- diff(y)/diff(x)
  int <- y[1L] - slope * x[1L]

  p <- ggplot(data=d, aes(sample=std.resid)) +
    stat_qq(shape=1, size=3) +           # open circles
    labs(title="Normal Q-Q",             # plot title
         x="Theoretical Quantiles",      # x-axis label
         y="Standardized Residuals") +   # y-axis label
    geom_abline(slope = slope, intercept = int, linetype="dashed")  # dashed reference line
  return(p)
}

Ejemplo de uso:

# sample data (y = x + N(0,1), x in [1,100])
df <- data.frame(cbind(x=c(1:100),y=c(1:100+rnorm(100))))
ggQQ(lm(y~x,data=df))

Answer 4

Desde la versión 2.0, ggplot2 tiene una interfaz bien documentada para la extensión; por lo que ahora puede escribir una nueva estadística para el qqline por sí mismo (que yo he hecho, por primera vez, por lo que las mejoras son bienvenidos ):

qq.line <- 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, ...))
}

Esto también generaliza sobre la distribución (exactamente como stat_qq hace), y se puede utilizar como sigue:

> test.data <- data.frame(sample=rnorm(100, 10, 2)) # normal distribution
> test.data.2 <- data.frame(sample=rt(100, df=2))   # t distribution
> 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))

(Por desgracia, ya que el qqline está en una capa separada, no pude encontrar una manera de "reutilización" los parámetros de la distribución, pero que sólo debería ser un problema menor.)

Answer 5

¿Por qué no lo siguiente?

Dado algún vector, por ejemplo,

myresiduals <- rnorm(100) ^ 2

ggplot(data=as.data.frame(qqnorm( myresiduals , plot=F)), mapping=aes(x=x, y=y)) + 
    geom_point() + geom_smooth(method="lm", se=FALSE)

Pero parece extraño que tenemos que utilizar una gráfica tradicional funcionar para apuntalar ggplot2.

¿No podemos obtener el mismo efecto de alguna manera, comenzando con el vector para los que queremos la trama cuantil y luego aplicar la apropiada "stat" y "funciones" en geom ggplot2?

¿El Hadley Wickham controlar estos mensajes? Tal vez nos puede mostrar una mejor manera.

Answer 6

Con la versión más reciente ggplot2 (> = 3.0), se implementó una nueva función stat_qq_line ( https://github.com/tidyverse/ggplot2/blob/master/NEWS.md ) y una línea qq para residuos del modelo se puede añadir con:

library(ggplot2)
model <- lm(mpg ~ wt, data=mtcars)
ggplot(model, aes(sample = rstandard(model))) + geom_qq() + stat_qq_line()

Se necesita

rstandard(model) para obtener el residuo estandarizado. (@Jlhoward crédito y @qwr)

Si recibe un 'Error en stat_qq_line (): no se pudo encontrar la función 'stat_qq_line'', su versión ggplot2 es demasiado viejo y lo puede solucionar mediante la actualización de su paquete ggplot2:. install.packages("ggplot2")

Answer 7

Se podría robar una página de los veteranos que hicieron estas cosas con papel de probabilidad normal. Una mirada cuidadosa a una ggplot () + stat_qq (gráfico) sugiere que una línea de referencia puede ser añadido con geom_abline (), como este

df <- data.frame( y=rpois(100, 4) )

ggplot(df, aes(sample=y)) +
  stat_qq() +
  geom_abline(intercept=mean(df$y), slope = sd(df$y))

Answer 8

v.3.0.0 ggplot2 tiene ahora una estadística qqline. Desde la página de ayuda:

df <- data.frame(y = rt(200, df = 5))
p <- ggplot(df, aes(sample = y))
p + stat_qq() + stat_qq_line()

! Ggplot2 v3.0.0 Ejemplo stats equivalente a qqnorm más abline ] 1