1 ## ltt.plot.R (2011-09-24)
3 ## Lineages Through Time Plot
5 ## Copyright 2002-2011 Emmanuel Paradis
7 ## This file is part of the R-package `ape'.
8 ## See the file ../COPYING for licensing issues.
10 ltt.plot.coords <- function(phy, backward = TRUE, tol = 1e-6)
12 if (is.ultrametric(phy, tol)) {
13 if (is.binary.tree(phy)) {
14 N <- numeric(phy$Nnode + 1)
17 node.order <- tabulate(phy$edge[, 1])
18 N <- node.order[-(1:length(phy$tip.label))] - 1
20 bt <- branching.times(phy)
22 o <- order(bt, decreasing = TRUE)
24 if (!is.binary.tree(phy)) N <- c(1, N[o])
26 if (!is.binary.tree(phy)) phy <- multi2di(phy)
30 event <- time.event <- numeric(n + m)
35 for (i in 1:nrow(phy$edge))
36 time.event[phy$edge[i, 2]] <- time.event[phy$edge[i, 1]] + phy$edge.length[i]
38 present <- max(time.event)
40 event[ROOT:(n + m)] <- 1
42 ## delete the events that are too close to present:
43 past.event <- present - time.event > tol
44 event <- event[past.event]
45 time.event <- time.event[past.event]
48 o <- order(time.event)
49 time.event <- time.event[o]
52 time <- c(time.event - present, 0)
56 if (!is.null(phy$root.edge)) {
57 time <- c(time[1] - phy$root.edge, time)
60 if (!backward) time <- time - time[1]
64 ltt.plot <- function(phy, xlab = "Time", ylab = "N",
65 backward = TRUE, tol = 1e-6, ...)
67 if (!inherits(phy, "phylo"))
68 stop("object \"phy\" is not of class \"phylo\"")
70 xy <- ltt.plot.coords(phy, backward, tol)
72 plot.default(xy, xlab = xlab, ylab = ylab, xaxs = "r",
73 yaxs = "r", type = "S", ...)
76 ltt.lines <- function(phy, backward = TRUE, tol = 1e-6, ...)
78 xy <- ltt.plot.coords(phy, backward, tol)
79 lines(xy, type = "S", ...)
83 function(phy, ..., dcol = TRUE, dlty = FALSE, legend = TRUE,
84 xlab = "Time", ylab = "N", log = "", backward = TRUE, tol = 1e-6)
86 if (inherits(phy, "phylo")) { # if a tree of class "phylo"
87 TREES <- list(ltt.plot.coords(phy, backward, tol))
88 names(TREES) <- deparse(substitute(phy))
89 } else { # a list of trees
90 TREES <- lapply(phy, ltt.plot.coords, backward = backward, tol = tol)
91 names(TREES) <- names(phy)
92 if (is.null(names(TREES)))
94 paste(deparse(substitute(phy)), "-", 1:length(TREES))
99 mc <- as.character(match.call())[-(1:2)]
102 if (class(dts[[i]]) == "phylo") {
103 a <- list(ltt.plot.coords(dts[[i]], backward, tol))
105 } else { # a list of trees
106 a <- lapply(dts[[i]], ltt.plot.coords, backward = backward, tol = tol)
107 names(a) <- names(dts[[i]])
108 if (is.null(names(a)))
109 names(a) <- paste(deparse(substitute(phy)), "-", 1:length(a))
115 xl <- c(min(unlist(lapply(TREES, function(x) min(x[, 1])))), 0)
116 yl <- c(1, max(unlist(lapply(TREES, function(x) max(x[, 2])))))
118 plot.default(NA, type = "n", xlim = xl, ylim = yl, xaxs = "r",
119 yaxs = "r", xlab = xlab, ylab = ylab, log = log)
121 lty <- if (!dlty) rep(1, n) else 1:n
122 col <- if (!dcol) rep(1, n) else topo.colors(n)
125 lines(TREES[[i]], col = col[i], lty = lty[i], type = "S")
128 legend(xl[1], yl[2], legend = names(TREES),
129 lty = lty, col = col, bty = "n")
132 ltt.coplot <- function(phy, backward = TRUE, ...)
134 layout(matrix(1:2, 2))
135 par(mar = c(0, 3, 0.5, 0.5))
136 o <- plot(phy, root.edge = TRUE, ...)
137 par(mar = c(3, 3, 0, 0.5))
138 ltt.plot(phy, xlim = o$x.lim, backward = FALSE, xaxt = "n")
139 if (backward) axisPhylo() else axis(1)