Exemplo

Com a multiplicidade de distribuições que podem ser ajustadas, a capacidade de comparar os ajustes de vários modelos com os dados é importante. Um gráfico QQ-plot pode ser usado para comparar duas distribuições, traçando os quantis “Q” um contra o outro.

Se as distribuições forem iguais, devemos esperar que os pontos estejam aproximadamente na reta \(x = y\). No exemplo a seguir, a resposta \(y1\) segue uma distribuição de Birnbaum-Saunders e \(y2\) é distribuída por Fréchet. Ajustando uma distribuição de Birnbaum-Saunders a ambos, os gráficos QQ-plot resultantes confirmam que o ajuste a \(y1\) é razoável, mas o ajuste a \(y2\) não é, ver Figura 1. 

library(VGAM)
n <- 100; set.seed(1)
bdata <- data.frame(x1 = rep(1, n)) # bdata has n rows
bdata <- transform(bdata, Scale = exp(1), Shape = exp(0.5)) # Recycling
bdata <- transform(bdata, y1 = rbisa(n, scale = Scale, shape = Shape),
y2 = rfrechet(n, scale = Scale, shape = Shape))
fit1 <- vglm(y1 ~ 1, bisa, data = bdata, trace = FALSE)
fit2 <- vglm(y2 ~ 1, bisa, data = bdata, trace = FALSE)
pfit1 <- predict(fit1, untransform = TRUE)
pfit2 <- predict(fit2, untransform = TRUE)
bdata <- transform(bdata, scale1.hat = pfit1[, "scale"], 
                   shape1.hat = pfit1[, "shape"], scale2.hat = pfit2[, "scale"], 
                   shape2.hat = pfit2[, "shape"])
pvector <- ppoints(n) # n equally-spaced points on (0, 1)
par(mfrow=c(2,2), pch=19,mar=c(4,4,2,1))
plot(qbisa(pvector, scale = scale1.hat, shape = shape1.hat) ~ sort(y1), 
     data = bdata, col = "blue", log = "xy", ylab = "Estimated quantiles", 
     main = "(a); q-type function")
abline(a = 0, b = 1, col = "orange", lty = "dashed"); grid()
plot(sort(pbisa(y1, scale = scale1.hat, shape = shape1.hat)) ~ pvector, 
     data = bdata, col = "blue", ylab = "Estimated CDF from the model", 
     main = "(b); p-type function")
abline(a = 0, b = 1, col = "orange", lty = "dashed"); grid()
plot(qbisa(pvector, scale = scale2.hat, shape = shape2.hat) ~ sort(y2), 
     data = bdata, col = "blue", log = "xy", ylab = "Estimated quantiles", 
     main = "(c); q-type function")
abline(a = 0, b = 1, col = "orange", lty = "dashed");grid()
plot(sort(pbisa(y2, scale = scale2.hat, shape = shape2.hat)) ~ pvector, 
     data = bdata, col = "blue", ylab = "Estimated CDF from the model", 
     main = "(d); p-type function")
abline(a = 0, b = 1, col = "orange", lty = "dashed");grid()

Figura 1: Gráficos QQ-plot de 2 conjuntos de dados simulados. Uma distribuição Birnbaum-Saunders é adequada para ambos. (a)–(b) Dados de uma distribuição Birnbaum-Saunders; (c)–(d) dados de uma distribuição de Fréchet. Os gráficos estão em escala logarítmica. Uma linha \(x = y\) aparece em todos os gráficos.


Teoria

Para explicar alguns detalhes por trás disso, a função R ppoints(n) é usada para formar um vetor igualmente espaçado de comprimento \(n\) no intervalo (0,1), por exemplo, as posições mostradas gráficamente são algo semelhante a \[ p_i=\dfrac{i + 12}{1 + n}\cdot \]

Se \(Y_1,\cdots,Y_n\) são uma amostra aleatória de alguma distribuição contínua com função de distribuição \(F_Y(y;\theta)\), então \(U_i = F_Y (y_i ;\theta)\) tem uma distribuição uniforme padrão.

Seja \(Y_{(i)}\) a \(i\)-ésima estatística de ordem, de modo que \(Y_{(1)}\leq Y_{(2)}\leq \cdots \leq Y_{(n)}\). Então \[ F_Y (Y_{(i)} ; \theta) = U_{(i)}\approx p_i\cdot \]

Os gráficos (a) e (c) verificam se \(F^{-1}(p_i;\widehat{\theta})\approx y_{(i)}\), enquanto os gráficos (b) e (d) verificam se \(F(y_{(i)};\widehat{\theta})\approx p_i\). Para o leitor interessado, as ideias aqui são emprestadas de, por exemplo, qqnorm() e qqplot().

Um erro comum entre analistas de dados amadores é a interpretação exagerada dos gráficos QQ-plot. Quando \(n\) é pequeno, há uma variação estatística substancial entre as estatísticas de ordem, portanto, não é razoável esperar que os pontos traçados fiquem muito próximos da linha \(x = y\). É comum quantificar a concordância de uma distribuição ajustada com os dados observados citando um coeficiente de correlação.