From 039cecc7ab5174412f1be3164f991f1ddc637f83 Mon Sep 17 00:00:00 2001 From: paradis Date: Tue, 23 Jun 2009 13:10:19 +0000 Subject: [PATCH] fix a bug in drop.tip + add new option 'rooted' git-svn-id: https://svn.mpl.ird.fr/ape/dev/ape@79 6e262413-ae40-0410-9e79-b911bd7a66b7 --- ChangeLog | 5 +++ DESCRIPTION | 2 +- R/drop.tip.R | 94 +++++++++++++++++++++++++++++++------------------ man/drop.tip.Rd | 12 ++++--- 4 files changed, 73 insertions(+), 40 deletions(-) diff --git a/ChangeLog b/ChangeLog index a1c4552..98602f9 100644 --- a/ChangeLog +++ b/ChangeLog @@ -10,6 +10,9 @@ NEW FEATURES o yule.cov() now fits the null model, and its help page has been corrected with respect to this change. + o drop.tip() has a new option 'rooted' to force (or not) a tree + to be treated as (un)rooted. + BUG FIXES @@ -30,6 +33,8 @@ BUG FIXES o write.tree() failed with objects of class "multiPhylo". + o drop.tip(, subtree = TRUE) sometimes shuffled tip labels. + OTHER CHANGES diff --git a/DESCRIPTION b/DESCRIPTION index 0b18830..4a0986b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: ape Version: 2.3-1 -Date: 2009-06-19 +Date: 2009-06-23 Title: Analyses of Phylogenetics and Evolution Author: Emmanuel Paradis, Ben Bolker, Julien Claude, Hoa Sien Cuong, Richard Desper, Benoit Durand, Julien Dutheil, Olivier Gascuel, Gangolf Jobb, Christoph Heibl, Daniel Lawson, Vincent Lefort, Pierre Legendre, Jim Lemon, Yvonnick Noel, Johan Nylander, Rainer Opgen-Rhein, Korbinian Strimmer, Damien de Vienne Maintainer: Emmanuel Paradis diff --git a/R/drop.tip.R b/R/drop.tip.R index a4a3d03..6e497cb 100644 --- a/R/drop.tip.R +++ b/R/drop.tip.R @@ -1,4 +1,4 @@ -## drop.tip.R (2009-05-10) +## drop.tip.R (2009-06-23) ## Remove Tips in a Phylogenetic Tree @@ -72,12 +72,23 @@ extract.clade <- function(phy, node, root.edge = 0) } drop.tip <- - function(phy, tip, trim.internal = TRUE, subtree = FALSE, root.edge = 0) + function(phy, tip, trim.internal = TRUE, subtree = FALSE, + root.edge = 0, rooted = is.rooted(phy)) { if (!inherits(phy, "phylo")) stop('object "phy" is not of class "phylo"') - phy <- reorder(phy) + Ntip <- length(phy$tip.label) + ## find the tips to drop: + if (is.character(tip)) + tip <- which(phy$tip.label %in% tip) + + if (!rooted && subtree) { + phy <- root(phy, (1:Ntip)[-tip][1]) + root.edge <- 0 + } + + phy <- reorder(phy) NEWROOT <- ROOT <- Ntip + 1 Nnode <- phy$Nnode Nedge <- dim(phy$edge)[1] @@ -93,19 +104,26 @@ drop.tip <- edge1 <- phy$edge[, 1] # local copies edge2 <- phy$edge[, 2] # keep <- !logical(Nedge) + ## find the tips to drop: if (is.character(tip)) tip <- which(phy$tip.label %in% tip) + + if (!rooted && subtree) { + phy <- root(phy, (1:Ntip)[-tip][1]) + root.edge <- 0 + } + ## delete the terminal edges given by `tip': keep[match(tip, edge2)] <- FALSE if (trim.internal) { - internals <- edge2 <= Ntip + ints <- edge2 > Ntip ## delete the internal edges that do not have anymore ## descendants (ie, they are in the 2nd col of `edge' but ## not in the 1st one) repeat { - sel <- !(edge2 %in% edge1[keep]) & internals & keep + sel <- !(edge2 %in% edge1[keep]) & ints & keep if (!sum(sel)) break keep[sel] <- FALSE } @@ -137,48 +155,54 @@ drop.tip <- if (!root.edge) phy$root.edge <- NULL - ## upate the tree; 1) drop the edges and tip labels + ## drop the edges phy$edge <- phy$edge[keep, ] if (wbl) phy$edge.length <- phy$edge.length[keep] - phy$tip.label <- phy$tip.label[-tip] - ## 2) renumber the remaining tips now - TIPS <- phy$edge[, 2] <= Ntip - ## keep the ordering so no need to reorder tip.label: - phy$edge[TIPS, 2] <- rank(phy$edge[TIPS, 2]) - ## 3) update node.label if needed - if (!is.null(phy$node.label)) - phy$node.label <- phy$node.label[sort(unique(phy$edge[, 1])) - Ntip] - Ntip <- length(phy$tip.label) # 4) update Ntip - ## make new tip labels if necessary + ## find the new terminal edges (works whatever 'subtree' and 'trim.internal'): + TERMS <- !(phy$edge[, 2] %in% phy$edge[, 1]) + + ## get the old No. of the nodes and tips that become tips: + oldNo.ofNewTips <- phy$edge[TERMS, 2] + + n <- length(oldNo.ofNewTips) # the new number of tips in the tree + + ## assumes that the ordering of tips is unchanged: + phy$edge[TERMS, 2] <- 1:n + + ## make new tip labels if necessary: if (subtree || !trim.internal) { - new.trms <- !(phy$edge[, 2] %in% phy$edge[, 1]) & phy$edge[, 2] > Ntip - node2tip <- phy$edge[new.trms, 2] - if (subtree) - new.lab <- paste("[", N[node2tip], "_tips]", sep = "") - else { - new.lab <- - if (is.null(phy$node.label)) rep("NA", length(node2tip)) - else phy$node.label[node2tip - Ntip] + ## get the logical indices of the tips that are kept within 'oldNo.ofNewTips': + tips.kept <- oldNo.ofNewTips <= Ntip & !(oldNo.ofNewTips %in% tip) + new.tip.label <- character(n) + new.tip.label[tips.kept] <- phy$tip.label[-tip] + ## get the numbers of the nodes that become tips: + node2tip <- oldNo.ofNewTips[!tips.kept] + new.tip.label[!tips.kept] <- if (subtree) { + paste("[", N[node2tip], "_tips]", sep = "") + } else { + if (is.null(phy$node.label)) rep("NA", length(node2tip)) + else phy$node.label[node2tip - Ntip] } - ## change the #'s in the edge matrix - new.tip <- Ntip + 1:length(node2tip) - phy$edge[new.trms, 2] <- new.tip - phy$tip.label[new.tip] <- new.lab - Ntip <- length(phy$tip.label) if (!is.null(phy$node.label)) phy$node.label <- phy$node.label[-(node2tip - Ntip)] - } - phy$Nnode <- dim(phy$edge)[1] - Ntip + 1L # 3) update phy$Nnode + phy$tip.label <- new.tip.label + } else phy$tip.label <- phy$tip.label[-tip] + + ## update node.label if needed: + if (!is.null(phy$node.label)) + phy$node.label <- phy$node.label[sort(unique(phy$edge[, 1])) - Ntip] + + phy$Nnode <- dim(phy$edge)[1] - n + 1L # update phy$Nnode ## The block below renumbers the nodes so that they conform ## to the "phylo" format -- same as in root() - newNb <- integer(Ntip + phy$Nnode) - newNb[NEWROOT] <- Ntip + 1L - sndcol <- phy$edge[, 2] > Ntip + newNb <- integer(n + phy$Nnode) + newNb[NEWROOT] <- n + 1L + sndcol <- phy$edge[, 2] > n ## executed from right to left, so newNb is modified before phy$edge: phy$edge[sndcol, 2] <- newNb[phy$edge[sndcol, 2]] <- - (Ntip + 2):(Ntip + phy$Nnode) + (n + 2):(n + phy$Nnode) phy$edge[, 1] <- newNb[phy$edge[, 1]] storage.mode(phy$edge) <- "integer" collapse.singles(phy) diff --git a/man/drop.tip.Rd b/man/drop.tip.Rd index f0f2aeb..0b2c32a 100644 --- a/man/drop.tip.Rd +++ b/man/drop.tip.Rd @@ -4,7 +4,7 @@ \title{Remove Tips in a Phylogenetic Tree} \usage{ drop.tip(phy, tip, trim.internal = TRUE, subtree = FALSE, - root.edge = 0) + root.edge = 0, rooted = is.rooted(phy)) extract.clade(phy, node, root.edge = 0) } \arguments{ @@ -18,6 +18,9 @@ extract.clade(phy, node, root.edge = 0) \item{root.edge}{an integer giving the number of internal branches to be used to build the new root edge. This has no effect if \code{trim.internal = FALSE}.} + \item{rooted}{a logical indicated whether the tree must be treated as + rooted or not. This allows to force the tree to be considered as + unrooted (see examples).} \item{node}{a node number or label.} } \description{ @@ -77,9 +80,10 @@ tip <- c( plot(drop.tip(bird.families, tip)) plot(drop.tip(bird.families, tip, trim.internal = FALSE)) data(bird.orders) -plot(drop.tip(bird.orders, 6:23, subtree = TRUE), font = 1) -plot(drop.tip(bird.orders, c(1:5, 20:23), subtree = TRUE), font = 1) - +plot(drop.tip(bird.orders, 6:23, subtree = TRUE)) +plot(drop.tip(bird.orders, c(1:5, 20:23), subtree = TRUE)) +plot(drop.tip(bird.orders, c(1:20, 23), subtree = TRUE)) +plot(drop.tip(bird.orders, c(1:20, 23), subtree = TRUE, rooted = FALSE)) ### Examples of the use of `root.edge' tr <- read.tree(text = "(A:1,(B:1,(C:1,(D:1,E:1):1):1):1):1;") drop.tip(tr, c("A", "B"), root.edge = 0) # = (C:1,(D:1,E:1):1); -- 2.39.5