# R/qqplot-das.R In StatDA: Statistical Analysis for Environmental Data

```"qqplot.das" <-
function (x, distribution = "norm", ylab = deparse(substitute(x)),
xlab = paste(distribution, "quantiles"), main = "", las = par("las"),
datax=FALSE, envelope = 0.95, labels = FALSE, col = palette()[2], lwd = 2,
pch = 1, line = c("quartiles", "robust", "none"), cex=1,xaxt="s",
add.plot=FALSE, xlim = NULL, ylim=NULL, ...)
{
# Modified by PF, 06.11.2009
# original function: qq.plot from library(car)
# If datax=TRUE: x and y axes are exchanged
result <- NULL
line <- match.arg(line)
good <- !is.na(x)
ord <- order(x[good])
ord.x <- x[good][ord]
q.function <- eval(parse(text = paste("q", distribution,
sep = "")))
d.function <- eval(parse(text = paste("d", distribution,
sep = "")))
n <- length(ord.x)
P <- ppoints(n)
z <- q.function(P, ...)
if (datax)
plot(ord.x, z, xlab = ylab, ylab = xlab, main = main,
las = las, col = col, pch = pch, cex = cex, xaxt = xaxt,
ylim = ylim,xlim=xlim)
else plot(z, ord.x, xlab = xlab, ylab = ylab, main = main,
las = las, col = col, pch = pch, cex = cex, xaxt = xaxt,
ylim = ylim,xlim=xlim)
}
else {
if (datax)
points(ord.x, z, xlab = ylab, ylab = xlab, main = main,
las = las, col = col, pch = pch, cex = cex, xaxt = xaxt)
else points(z, ord.x, xlab = xlab, ylab = ylab, main = main,
las = las, col = col, pch = pch, cex = cex, xaxt = xaxt)

}

if (line == "quartiles") {
Q.x <- quantile(ord.x, c(0.25, 0.75))
Q.z <- q.function(c(0.25, 0.75), ...)
if (datax) {
b <- (Q.z[2] - Q.z[1])/(Q.x[2] - Q.x[1])
a <- Q.z[1] - b * Q.x[1]
}
else {
b <- (Q.x[2] - Q.x[1])/(Q.z[2] - Q.z[1])
a <- Q.x[1] - b * Q.z[1]
}
abline(a, b, col = col, lwd = lwd)
}
if (line == "robust") {
#if (!require("MASS"))
#    stop("MASS package not available")
if (datax)
coef <- coefficients(rlm(z ~ ord.x))
else coef <- coefficients(rlm(ord.x ~ z))
a <- coef[1]
b <- coef[2]
abline(a, b)
}
if (line != "none" & envelope != FALSE) {
zz <- qnorm(1 - (1 - envelope)/2)
if (datax) {
b <- (Q.x[2] - Q.x[1])/(Q.z[2] - Q.z[1])
a <- Q.x[1] - b * Q.z[1]
SE <- (b/d.function(z, ...)) * sqrt(P * (1 - P)/n)
fit.value <- a + b * z
upper <- fit.value + zz * SE
lower <- fit.value - zz * SE
lines(upper, z, lty = 2, lwd = lwd/2, col = col)
lines(lower, z, lty = 2, lwd = lwd/2, col = col)
}
else {
SE <- (b/d.function(z, ...)) * sqrt(P * (1 - P)/n)
fit.value <- a + b * z
upper <- fit.value + zz * SE
lower <- fit.value - zz * SE
lines(z, upper, lty = 2, lwd = lwd/2, col = col)
lines(z, lower, lty = 2, lwd = lwd/2, col = col)
}
}
if (labels[1] == TRUE & length(labels) == 1)
labels <- seq(along = z)
if (labels[1] != FALSE) {
selected <- identify(z, ord.x, labels[good][ord])
result <- seq(along = x)[good][ord][selected]
}
if (is.null(result))
invisible(result)
else sort(result)
}
```

## Try the StatDA package in your browser

Any scripts or data that you put into this service are public.

StatDA documentation built on March 13, 2020, 2:42 a.m.