R/delta.plot.R

   1 ## delta.plot.R (2010-01-12)
   2
   3 ##   Delta Plots
   4
   5 ## Copyright 2010 Emmanuel Paradis
   6
   7 ## This file is part of the R-package `ape'.
   8 ## See the file ../COPYING for licensing issues.
   9
  10 delta.plot <- function(X, k = 20, plot = TRUE, which = 1:2)
  11 {
  12     if (is.matrix(X)) X <- as.dist(X)
  13     n <- attr(X, "Size")
  14     if (n < 4) stop("need at least 4 observations")
  15
  16     ## add a category for the cases delta = 1
  17     ans <- .C("delta_plot", as.double(X), as.integer(n),
  18               as.integer(k), integer(k + 1), double(n),
  19               NAOK = TRUE, DUP = FALSE, PACKAGE = "ape")
  20     counts <- ans[[4]]
  21     ## add the counts of delta=1 to the last category:
  22     counts[k] <- counts[k] + counts[k + 1]
  23     counts <- counts[-(k + 1)]
  24
  25     delta.bar <- ans[[5]]/choose(n - 1, 3)
  26
  27     if (plot) {
  28         if (length(which) == 2) layout(matrix(1:2, 1, 2))
  29         if (1 %in% which) {
  30             barplot(counts, space = 0, xlab = expression(delta[q]))
  31             a <- axTicks(1)
  32             axis(1, at = a, labels = a/k)
  33         }
  34         if (2 %in% which)
  35             plot(delta.bar, type = "h", ylab = expression(bar(delta)))
  36     }
  37
  38     invisible(list(counts = counts, delta.bar = delta.bar))
  39 }