## delta.plot.R (2010-01-12)
## Delta Plots
## Copyright 2010 Emmanuel Paradis
## This file is part of the R-package `ape'.
## See the file ../COPYING for licensing issues.
delta.plot <- function(X, k = 20, plot = TRUE, which = 1:2)
{
if (is.matrix(X)) X <- as.dist(X)
n <- attr(X, "Size")
if (n < 4) stop("need at least 4 observations")
## add a category for the cases delta = 1
ans <- .C("delta_plot", as.double(X), as.integer(n),
as.integer(k), integer(k + 1), double(n),
NAOK = TRUE, DUP = FALSE, PACKAGE = "ape")
counts <- ans[[4]]
## add the counts of delta=1 to the last category:
counts[k] <- counts[k] + counts[k + 1]
counts <- counts[-(k + 1)]
delta.bar <- ans[[5]]/choose(n - 1, 3)
if (plot) {
if (length(which) == 2) layout(matrix(1:2, 1, 2))
if (1 %in% which) {
barplot(counts, space = 0, xlab = expression(delta[q]))
a <- axTicks(1)
axis(1, at = a, labels = a/k)
}
if (2 %in% which)
plot(delta.bar, type = "h", ylab = expression(bar(delta)))
}
invisible(list(counts = counts, delta.bar = delta.bar))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.