]> git.donarmstrong.com Git - ape.git/commitdiff
several bug fixes while in JKT
authorparadis <paradis@6e262413-ae40-0410-9e79-b911bd7a66b7>
Tue, 7 Jul 2009 06:35:45 +0000 (06:35 +0000)
committerparadis <paradis@6e262413-ae40-0410-9e79-b911bd7a66b7>
Tue, 7 Jul 2009 06:35:45 +0000 (06:35 +0000)
git-svn-id: https://svn.mpl.ird.fr/ape/dev/ape@81 6e262413-ae40-0410-9e79-b911bd7a66b7

ChangeLog
DESCRIPTION
R/all.equal.phylo.R
R/chronopl.R
R/dist.topo.R
R/drop.tip.R
R/root.R
R/vcv.phylo.R
Thanks

index 98602f93fd35b93c2965986ebbc415dd68b2f902..dbe6fe09acbc5774976984d7345941ff816c3597 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,23 @@
+               CHANGES IN APE VERSION 2.3-2
+
+
+BUG FIXES
+
+    o all.equal.phylo() did not compare unrooted trees correctly.
+
+    o dist.topo(... method = "PH85") did not treat unrooted trees
+      correctly (thanks to Tim Wallstrom for the fix).
+
+    o root() sometimes failed to test for the monophyly of the
+      outgroup correctly.
+
+    o extract.clade() sometimes included too many edges.
+
+    o vcv.phylo() did not work correctly when the tree is in
+      "pruningwise" order.
+
+
+
                CHANGES IN APE VERSION 2.3-1
 
 
index 4a0986bf7b27518aeaeb3c0e9cdbf3d4ceea0553..0ff3acf22f2c9e272140b0a23d8d325a3bd9261e 100644 (file)
@@ -1,6 +1,6 @@
 Package: ape
-Version: 2.3-1
-Date: 2009-06-23
+Version: 2.3-2
+Date: 2009-07-06
 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 <Emmanuel.Paradis@ird.fr>
index cd0543793d3de7733d1d6d91303b9e80276197c7..79533d5de6f10561b38aec06b78f3feac3792977 100644 (file)
@@ -1,4 +1,4 @@
-## all.equal.phylo.R (2006-09-12)
+## all.equal.phylo.R (2009-07-05)
 ##
 ##     Global Comparison of two Phylogenies
 
@@ -8,18 +8,18 @@
 ## This file is part of the R-package `ape'.
 ## See the file ../COPYING for licensing issues.
 
-### Recherche de la correspondance entre deux arbres
-### Parcours en profondeur et en parallèle des deux arbres (current et target)
-### current, target: les deux arbres à comparer
-### use.edge.length: faut-il comparer les longueurs de branches ?
-### use.tip.label: faut-il comparer les étiquettes de feuilles ou seulement la
-###    topologie des deux arbres ?
-### index.return: si TRUE, retourner la matrice de correspondance entre noeuds
-###    et feuilles, une matrice à deux colonnes (current et target) avec pour
-###    chaque ligne des paires d'identifiants de noeuds/feuilles, tels qu'ils
-###    apparaissent dans l'attribut 'edge' des objets phylo
-### tolerance, scale: paramètres de comparaison des longueurs de branches
-###    (voir 'all.equal')
+## Recherche de la correspondance entre deux arbres
+## Parcours en profondeur et en parallèle des deux arbres (current et target)
+## current, target: les deux arbres à comparer
+## use.edge.length: faut-il comparer les longueurs de branches ?
+## use.tip.label: faut-il comparer les étiquettes de feuilles ou seulement la
+##     topologie des deux arbres ?
+## index.return: si TRUE, retourner la matrice de correspondance entre noeuds
+##     et feuilles, une matrice à deux colonnes (current et target) avec pour
+##     chaque ligne des paires d'identifiants de noeuds/feuilles, tels qu'ils
+##     apparaissent dans l'attribut 'edge' des objets phylo
+## tolerance, scale: paramètres de comparaison des longueurs de branches
+##     (voir 'all.equal')
 all.equal.phylo <- function(target, current,
                         use.edge.length = TRUE,
                         use.tip.label = TRUE,
@@ -66,6 +66,14 @@ all.equal.phylo <- function(target, current,
         root1 <- Ntip1 + 1
         root2 <- Ntip2 + 1
         if (root1 != root2) return(FALSE)
+        ## Fix by EP so that unrooted trees are correctly compared:
+        if (!is.rooted(target) && !is.rooted(current)) {
+            outg <- target$tip.label[1]
+            if (! outg %in% current$tip.label) return(FALSE)
+            target <- root(target, outg)
+            current <- root(current, outg)
+        }
+        ## End
        result <- same.node(root1, root2)
        if (!isTRUE(index.return)) return(!is.null(result))
        if (is.null(result)) return(result)
index 6a0bd89be6064680f5112f7e889b85156ab140ef..570fa2602fbd9166199c5f68dc6e81fb6b73a4f6 100644 (file)
@@ -1,15 +1,16 @@
-## chronopl.R (2008-11-04)
+## chronopl.R (2009-07-06)
 
 ##   Molecular Dating With Penalized Likelihood
 
-## Copyright 2005-2008 Emmanuel Paradis
+## Copyright 2005-2009 Emmanuel Paradis
 
 ## This file is part of the R-package `ape'.
 ## See the file ../COPYING for licensing issues.
 
-chronopl <- function(phy, lambda, age.min = 1, age.max = NULL,
-                     node = "root", S = 1, tol = 1e-8,
-                     CV = FALSE, eval.max = 500, iter.max = 500, ...)
+chronopl <-
+    function(phy, lambda, age.min = 1, age.max = NULL,
+             node = "root", S = 1, tol = 1e-8,
+             CV = FALSE, eval.max = 500, iter.max = 500, ...)
 {
     n <- length(phy$tip.label)
     ROOT <- n + 1
@@ -84,27 +85,23 @@ chronopl <- function(phy, lambda, age.min = 1, age.max = NULL,
         real.edge.length <- ini.time[e[, 1]] - ini.time[e[, 2]]
         while (any(real.edge.length <= 0)) {
             for (i in EDGES) {
-                if (real.edge.length[i] <= 0) {
-                    if (e[i, 1] %in% node) {
-                        ini.time[e[i, 2]] <-
-                            ini.time[e[, 2]] - 2*real.edge.length[i]
-                        next
-                    }
-                    if (e[i, 2] %in% node) {
-                        ini.time[e[i, 1]] <-
-                            ini.time[e[, 1]] + 2*real.edge.length[i]
-                        next
-                    }
-                    ini.time[e[i, 2]] <-
-                        ini.time[e[, 2]] - real.edge.length[i]
-                    ini.time[e[i, 1]] <-
-                        ini.time[e[, 1]] + real.edge.length[i]
+                if (real.edge.length[i] > 0) next
+                if (e[i, 1] %in% node) {
+                    ini.time[e[i, 2]] <- ini.time[e[1, 2]] - 2 * real.edge.length[i]
+                    next
                 }
+                if (e[i, 2] %in% node) {
+                    ini.time[e[i, 1]] <- ini.time[e[1, 1]] + 2 * real.edge.length[i]
+                    next
+                }
+                browser()
+                ini.time[e[i, 2]] <- ini.time[e[1, 2]] - real.edge.length[i]
+                ini.time[e[i, 1]] <- ini.time[e[1, 1]] + real.edge.length[i]
             }
             real.edge.length <- ini.time[e[, 1]] - ini.time[e[, 2]]
+            print(min(real.edge.length))
         }
     }
-
     ## `unknown.ages' will contain the index of the nodes of unknown age:
     unknown.ages <- n + 1:m
 
index 6741dcc8415f5ee1c33c1055aa58aa1e5fa55a49..1248c4983d4ef8126d90bd11fd2aa84ee50dc19a 100644 (file)
@@ -1,4 +1,4 @@
-## dist.topo.R (2009-05-10)
+## dist.topo.R (2009-07-06)
 
 ##      Topological Distances, Tree Bipartitions,
 ##   Consensus Trees, and Bootstrapping Phylogenies
@@ -15,21 +15,27 @@ dist.topo <- function(x, y, method = "PH85")
     n <- length(x$tip.label)
     bp1 <- .Call("bipartition", x$edge, n, x$Nnode, PACKAGE = "ape")
     bp1 <- lapply(bp1, function(xx) sort(x$tip.label[xx]))
-    bp2 <- .Call("bipartition", y$edge, n, y$Nnode, PACKAGE = "ape")
-    bp2 <- lapply(bp2, function(xx) sort(y$tip.label[xx]))
+    ## fix by Tim Wallstrom:
+    bp2.tmp <- .Call("bipartition", y$edge, n, y$Nnode, PACKAGE = "ape")
+    bp2 <- lapply(bp2.tmp, function(xx) sort(y$tip.label[xx]))
+    bp2.comp <- lapply(bp2.tmp, function(xx) setdiff(1:n, xx))
+    bp2.comp <- lapply(bp2.comp, function(xx) sort(y$tip.label[xx]))
+    ## End
     q1 <- length(bp1)
     q2 <- length(bp2)
     if (method == "PH85") {
         p <- 0
         for (i in 1:q1) {
             for (j in 1:q2) {
-                if (identical(all.equal(bp1[[i]], bp2[[j]]), TRUE)) {
+                if (identical(bp1[[i]], bp2[[j]]) |
+                    identical(bp1[[i]], bp2.comp[[j]])) {
                     p <- p + 1
                     break
                 }
             }
         }
-        dT <- if (q1 == q2) 2*(q1 - p) else 2*(min(q1, q2) - p) + abs(q1 - q2)
+        dT <- q1 + q2 - 2 * p # same than:
+        ##dT <- if (q1 == q2) 2*(q1 - p) else 2*(min(q1, q2) - p) + abs(q1 - q2)
     }
     if (method == "BHV01") {
         dT <- 0
index 6e497cb9c0245dd66e39bbd79d291bd484cab0af..acdbe688d7643f6465b9d04c5a30844ce027ffea 100644 (file)
@@ -1,4 +1,4 @@
-## drop.tip.R (2009-06-23)
+## drop.tip.R (2009-07-06)
 
 ##   Remove Tips in a Phylogenetic Tree
 
@@ -29,7 +29,7 @@ extract.clade <- function(phy, node, root.edge = 0)
     root.node <- which(phy$edge[, 2] == node)
     start <- root.node + 1 # start of the clade looked for
     anc <- phy$edge[root.node, 1] # the ancestor of 'node'
-    next.anc <- which(phy$edge[-(1:start), 1] == anc) # find the next occurence of 'anc'
+    next.anc <- which(phy$edge[-(1:start), 1] <= anc) # find the next occurence of 'anc' or an 'older' node
 
     keep <- if (length(next.anc)) start + 0:(next.anc[1] - 1) else start:Nedge
 
index 832ab0f544ba0ef235870701ceb12798552b6b3c..3abdf650e8ab06091ae20a75e91051b80ccc2b8b 100644 (file)
--- a/R/root.R
+++ b/R/root.R
@@ -1,4 +1,4 @@
-## root.R (2009-05-10)
+## root.R (2009-07-06)
 
 ##   Root of Phylogenetic Trees
 
@@ -11,11 +11,10 @@ is.rooted <- function(phy)
 {
     if (!inherits(phy, "phylo"))
       stop('object "phy" is not of class "phylo"')
-    if (!is.null(phy$root.edge)) return(TRUE)
+    if (!is.null(phy$root.edge)) TRUE
     else
       if (tabulate(phy$edge[, 1])[length(phy$tip.label) + 1] > 2)
-        return(FALSE)
-      else return(TRUE)
+        FALSE else TRUE
 }
 
 unroot <- function(phy)
@@ -66,8 +65,7 @@ root <- function(phy, outgroup, node = NULL, resolve.root = FALSE)
 {
     if (!inherits(phy, "phylo"))
       stop('object "phy" is not of class "phylo"')
-    ord <- attr(phy, "order")
-    if (!is.null(ord) && ord == "pruningwise") phy <- reorder(phy)
+    phy <- reorder(phy)
     n <- length(phy$tip.label)
     ROOT <- n + 1
     if (!is.null(node)) {
@@ -88,7 +86,6 @@ root <- function(phy, outgroup, node = NULL, resolve.root = FALSE)
         ## First check that the outgroup is monophyletic--
         ## unless there's only one tip specified of course
         if (length(outgroup) > 1) {
-            msg <- "the specified outgroup is not monophyletic"
             seq.nod <- .Call("seq_root2tip", phy$edge, n,
                              phy$Nnode, PACKAGE = "ape")
             sn <- seq.nod[outgroup]
@@ -108,9 +105,15 @@ root <- function(phy, outgroup, node = NULL, resolve.root = FALSE)
             ## (below is slightly faster than calling "bipartition")
             desc <- which(unlist(lapply(seq.nod,
                                         function(x) any(x %in% newroot))))
-            if (length(outgroup) != length(desc)) stop(msg)
-            ## both vectors below are already sorted:
-            if (!all(outgroup == desc)) stop(msg)
+            msg <- "the specified outgroup is not monophyletic"
+            ingroup <- (1:n)[-outgroup]
+            ## 'outgroup' and 'desc' are already sorted:
+            if (newroot != ROOT) {
+                if (!identical(outgroup, desc) && !identical(ingroup, desc))
+                    stop(msg)
+            } else { # otherwise check monophyly of the ingroup
+                if (!is.monophyletic(phy, ingroup)) stop(msg)
+            }
         } else newroot <- phy$edge[which(phy$edge[, 2] == outgroup), 1]
     }
     N <- Nedge(phy)
@@ -128,7 +131,7 @@ root <- function(phy, outgroup, node = NULL, resolve.root = FALSE)
                 if (!is.null(phy$edge.length))
                 phy$edge.length <-
                     c(phy$edge.length[a], 0, phy$edge.length[b])
-                phy$Nnode <- phy$Nnode + 1
+                phy$Nnode <- phy$Nnode + 1L
                 ## node renumbering (see comments below)
                 newNb <- integer(n + oldNnode)
                 newNb[newroot] <- n + 1L
index a51c6e4c35f243376237dc82cd9ee9c892474a93..01cbff7827fd6c74547e12ae6f32b12545cb586a 100644 (file)
@@ -1,4 +1,4 @@
-## vcv.phylo.R (2009-05-10)
+## vcv.phylo.R (2009-07-06)
 
 ##   Phylogenetic Variance-Covariance or Correlation Matrix
 
@@ -37,6 +37,7 @@ vcv.phylo <- function(phy, model = "Brownian", cor = FALSE)
         }
     }
 
+    phy <- reorder(phy)
     n <- length(phy$tip.label)
     n.node <- phy$Nnode
     N <- n.node + n
diff --git a/Thanks b/Thanks
index 157d4cc8dcee296a0d44a45c1c84655e5ff0d0e0..abbeca8eb86ecef466cd7cd081a514be52419a46 100644 (file)
--- a/Thanks
+++ b/Thanks
@@ -6,8 +6,8 @@ comments, or bug reports: thanks to all of you!
 
 Significant bug fixes were provided by Cécile Ané, James Bullard,
 Éric Durand, Olivier François, Rich FitzJohn Bret Larget, Nick Matzke,
-Elizabeth Purdom, Dan Rabosky, Klaus Schliep, Li-San Wang, Yan Wong,
-and Peter Wragg. Contact me if I forgot someone.
+Elizabeth Purdom, Dan Rabosky, Klaus Schliep, Tim Wallstrom, Li-San
+Wang, Yan Wong, and Peter Wragg. Contact me if I forgot someone.
 
 Kurt Hornik, of the R Core Team, helped in several occasions to
 fix some problems and bugs.