From: paradis Date: Sat, 31 Mar 2012 10:09:02 +0000 (+0000) Subject: various bug fixes X-Git-Url: https://git.donarmstrong.com/?p=ape.git;a=commitdiff_plain;h=756ad52b92dc1ac2922cf62ce882469ad4cacc2c various bug fixes git-svn-id: https://svn.mpl.ird.fr/ape/dev/ape@183 6e262413-ae40-0410-9e79-b911bd7a66b7 --- diff --git a/DESCRIPTION b/DESCRIPTION index 903c8b1..a048969 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: ape Version: 3.0-2 -Date: 2012-03-22 +Date: 2012-03-30 Title: Analyses of Phylogenetics and Evolution Author: Emmanuel Paradis, Ben Bolker, Julien Claude, Hoa Sien Cuong, Richard Desper, Benoit Durand, Julien Dutheil, Olivier Gascuel, Christoph Heibl, Daniel Lawson, Vincent Lefort, Pierre Legendre, Jim Lemon, Yvonnick Noel, Johan Nylander, Rainer Opgen-Rhein, Andrei-Alin Popescu, Klaus Schliep, Korbinian Strimmer, Damien de Vienne Maintainer: Emmanuel Paradis @@ -9,16 +9,19 @@ Suggests: gee Imports: gee, nlme, lattice ZipData: no Description: ape provides functions for reading, writing, plotting, - and manipulating phylogenetic trees, analyses of comparative data - in a phylogenetic framework, analyses of diversification and - macroevolution, computing distances from allelic and nucleotide - data, reading nucleotide sequences, and several tools such as - Mantel's test, computation of minimum spanning tree, generalized - skyline plots, estimation of absolute evolutionary rates and - clock-like trees using mean path lengths, non-parametric rate - smoothing and penalized likelihood. Phylogeny estimation can be - done with the NJ, BIONJ, ME, MVR, SDM, and triangle methods, and - several methods handling incomplete distance matrices (NJ*, BIONJ*, - MVR*, and the corresponding triangle method). + and manipulating phylogenetic trees, analyses of comparative data in + a phylogenetic framework, ancestral character analyses, analyses of + diversification and macroevolution, computing distances from allelic + and nucleotide data, reading and writing nucleotide sequences, and + several tools such as Mantel's test, minimum spanning tree, + generalized skyline plots, graphical exploration of phylogenetic + data (alex, trex, kronoviz), estimation of absolute evolutionary + rates and clock-like trees using mean path lengths and penalized + likelihood. Phylogeny estimation can be done with the NJ, BIONJ, ME, + MVR, SDM, and triangle methods, and several methods handling + incomplete distance matrices (NJ*, BIONJ*, MVR*, and the + corresponding triangle method). Some functions call external + applications (PhyML, Clustal, T-Coffee, Muscle) whose results are + returned into R. License: GPL (>= 2) URL: http://ape.mpl.ird.fr/ diff --git a/NEWS b/NEWS index de2da10..acbbbf9 100644 --- a/NEWS +++ b/NEWS @@ -17,14 +17,27 @@ BUG FIXES o mltt.plot(, backward = FALSE) did not set the x-axis correctly. - o A bug was introduced in prop.clades() with ape 3.0. + o A bug was introduced in prop.clades() with ape 3.0. The help page + has been clarified relative to the use of the 'rooted' option. o mantel.test() printed a useless warning message. o plot.phylo(, direction = "downward") ignored 'y.lim'. o is.monophyletic() did not work correctly if 'tips' was not stored - as integers + as integers. + + o prop.part() could make R crash if the first tree had many + multichotomies. + + o njs(), bionjs(), and mvrs() now return an error if 'fs < 1'. + + +OTHER CHANGES + + o The DESCRIPTION file has been updated. + + o The option 'original.data' of write.nexus() has been removed. diff --git a/R/SDM.R b/R/SDM.R index 94a8f31..faa263e 100644 --- a/R/SDM.R +++ b/R/SDM.R @@ -18,7 +18,7 @@ SDM <- function(...) tot <- tot + length(rownames(st[[i]])) } sp <- mat.or.vec(1,k) - for (i in c(k+1:k)) + for (i in c(k + 1:k)) sp[i - k] <- st[[i]] astart <- mat.or.vec(1, tot) # start of aip, astart[p] is start of aip @@ -65,7 +65,7 @@ SDM <- function(...) dijp <- -1 if (l == p) { # calculate alpha_p for (i in 1:n) { - for (j in 1:n) { #check if {i,j}\subset L_l + for (j in 1:n) { # check if {i,j}\subset L_l d <- st[[l]] if (i != j & is.element(labels[i],rownames(d)) & is.element(labels[j],colnames(d))) { dij <- d[rownames(d) == labels[i], colnames(d) == labels[j]] @@ -79,7 +79,7 @@ SDM <- function(...) } } else { for (i in 1:n) { - for (j in 1:n) { #check if {i,j}\subset L_l + for (j in 1:n) { # check if {i,j}\subset L_l d <- st[[l]] if (i != j & is.element(labels[i], rownames(d)) & is.element(labels[j], colnames(d)) & is.element(labels[i], rownames(d_p)) & is.element(labels[j], colnames(d_p))) { dij <- d[rownames(d) == labels[i], colnames(d) == labels[j]] @@ -160,11 +160,11 @@ SDM <- function(...) } } a <- solve(Q, col, 1e-19) - for(i in 1:n) { - for(j in 1:n) { + for (i in 1:n) { + for (j in 1:n) { sum <- 0 sumv <- 0 - for(p in 1:k) { + for (p in 1:k) { d <- st[[p]] if (is.element(labels[i], rownames(d)) & is.element(labels[j], rownames(d))) { ipos <- which(rownames(d) == labels[i]) diff --git a/R/ltt.plot.R b/R/ltt.plot.R index 520294f..98f188d 100644 --- a/R/ltt.plot.R +++ b/R/ltt.plot.R @@ -99,7 +99,7 @@ mltt.plot <- mc <- as.character(match.call())[-(1:2)] nms <- mc[1:n] for (i in 1:n) { - if (class(dts[[i]]) == "phylo") { + if (inherits(dts[[i]], "phylo")) { a <- list(ltt.plot.coords(dts[[i]], backward, tol)) names(a) <- nms[i] } else { # a list of trees diff --git a/R/mvr.R b/R/mvr.R index a2218fc..8c4c8b3 100644 --- a/R/mvr.R +++ b/R/mvr.R @@ -1,4 +1,4 @@ -## mvr.R (2011-10-11) +## mvr.R (2012-03-30) ## Minimum Variance Reduction @@ -29,6 +29,9 @@ mvr <- function(X, V) mvrs <- function(X, V, fs = 15) { + if (fs < 1) + stop("argument 'fs' must be a non-zero positive integer") + if (is.matrix(X)) X <- as.dist(X) if (is.matrix(V)) V <- as.dist(V) diff --git a/R/njs.R b/R/njs.R index 1f7ee8a..4345bfe 100644 --- a/R/njs.R +++ b/R/njs.R @@ -1,4 +1,4 @@ -## njs.R (2011-10-11) +## njs.R (2012-03-30) ## Tree Reconstruction from Incomplete Distances With NJ* or bio-NJ* @@ -9,6 +9,8 @@ .NJS_BIONJS <- function(libs, X, fs = 15) { + if (fs < 1) + stop("argument 'fs' must be a non-zero positive integer") if (is.matrix(X)) X <- as.dist(X) X[is.na(X)] <- -1 X[X < 0] <- -1 diff --git a/R/write.nexus.R b/R/write.nexus.R index 00fcffc..ac383f7 100644 --- a/R/write.nexus.R +++ b/R/write.nexus.R @@ -1,4 +1,4 @@ -## write.nexus.R (2012-03-02) +## write.nexus.R (2012-03-30) ## Write Tree File in Nexus Format @@ -7,7 +7,7 @@ ## This file is part of the R-package `ape'. ## See the file ../COPYING for licensing issues. -write.nexus <- function(..., file = "", translate = TRUE, original.data = NULL) +write.nexus <- function(..., file = "", translate = TRUE) { obj <- list(...) ## We insure that all trees are in a list, even if there is a single one: @@ -22,9 +22,6 @@ write.nexus <- function(..., file = "", translate = TRUE, original.data = NULL) cat(paste("[R-package APE, ", date(), "]\n\n", sep = ""), file = file, append = TRUE) - if (!is.null(original.data)) - warning("the option 'original.data' is deprecated and will be removed soon. Please update your code.") - N <- length(obj[[1]]$tip.label) cat("BEGIN TAXA;\n", file = file, append = TRUE) diff --git a/man/boot.phylo.Rd b/man/boot.phylo.Rd index fcbf1f3..8ac8668 100644 --- a/man/boot.phylo.Rd +++ b/man/boot.phylo.Rd @@ -90,6 +90,10 @@ prop.clades(phy, ..., part = NULL, rooted = FALSE) passed as \code{\dots} fulfills conditions (i) and (ii) above, then it might be faster to first call, e.g., \code{pp <- prop.part(...)}, then use the option \code{part}: \code{prop.clades(phy, part = pp)}. + + You have to be careful that by default \code{prop.clades} considers + the trees as unrooted and this may result in spurious results if the + trees are rooted (see examples). } \value{ \code{prop.part} returns an object of class \code{"prop.part"} which @@ -138,6 +142,13 @@ pp20 <- prop.part(TR) length(pp20) plot(pp10, pch = "x", col = 2) plot(pp20, pch = "x", col = 2) + +set.seed(1) +tr <- rtree(10) # rooted by default +prop.clades(tr, tr) # clearly wrong +prop.clades(tr, tr, rooted = TRUE) +tr <- rtree(10, rooted = FALSE) +prop.clades(tr, tr) # correct } \keyword{manip} \keyword{htest} diff --git a/man/mvr.Rd b/man/mvr.Rd index 68f0bc9..b6984fd 100644 --- a/man/mvr.Rd +++ b/man/mvr.Rd @@ -12,7 +12,8 @@ mvrs(X, V, fs = 15) \arguments{ \item{X}{a distance matrix.} \item{V}{a variance matrix.} - \item{fs}{agglomeration criterion parameter.} + \item{fs}{agglomeration criterion parameter: it is coerced as an + integer and must at least equal to one.} } \details{ The MVR method can be seen as a version of BIONJ which is not diff --git a/man/njs.Rd b/man/njs.Rd index 90a3444..04bafb9 100644 --- a/man/njs.Rd +++ b/man/njs.Rd @@ -12,15 +12,16 @@ bionjs(X, fs = 15) } \arguments{ \item{X}{a distance matrix.} - \item{fs}{argument \emph{s} of the agglomerative criterion.} + \item{fs}{argument \emph{s} of the agglomerative criterion: it is + coerced as an integer and must at least equal to one.} } \details{ Missing values represented by either \code{NA} or any negative number. Basically, the Q* criterion is applied to all the pairs of leaves, and the \emph{s} highest scoring ones are chosen for further analysis by - the agglomeration criteria that better handle handle missing - distances (see references for details). + the agglomeration criteria that better handle missing distances (see + references for details). } \value{ an object of class \code{"phylo"}. diff --git a/man/write.nexus.Rd b/man/write.nexus.Rd index faf2ef6..d4f145f 100644 --- a/man/write.nexus.Rd +++ b/man/write.nexus.Rd @@ -2,7 +2,7 @@ \alias{write.nexus} \title{Write Tree File in Nexus Format} \usage{ -write.nexus(..., file = "", translate = TRUE, original.data = NULL) +write.nexus(..., file = "", translate = TRUE) } \arguments{ \item{\dots}{either (i) a single object of class \code{"phylo"}, (ii) a @@ -14,7 +14,6 @@ write.nexus(..., file = "", translate = TRUE, original.data = NULL) \item{translate}{a logical, if \code{TRUE} (the default) a translation of the tip labels is done which are replaced in the parenthetic representation with tokens.} - \item{original.data}{deprecated; will be removed soon.} } \description{ This function writes trees in a file with the NEXUS format. diff --git a/src/bipartition.c b/src/bipartition.c index 884fc6c..d29d554 100644 --- a/src/bipartition.c +++ b/src/bipartition.c @@ -1,6 +1,6 @@ -/* bipartition.c 2011-10-22 */ +/* bipartition.c 2012-03-26 */ -/* Copyright 2005-2011 Emmanuel Paradis, and 2007 R Development Core Team */ +/* Copyright 2005-2012 Emmanuel Paradis, and 2007 R Development Core Team */ /* This file is part of the R-package `ape'. */ /* See the file ../COPYING for licensing issues. */ @@ -160,8 +160,8 @@ SEXP prop_part(SEXP TREES, SEXP nbtree, SEXP keep_partitions) INTEGER(nbtip)[0] = Ntip; INTEGER(nbnode)[0] = Nnode; - if (KeepPartition) Npart = Ntree*(Nnode - 1) + 1; - else Npart = Nnode; + if (KeepPartition) Npart = Ntree * (Ntip - 2) + 1; + else Npart = Ntip - 1; PROTECT(number = allocVector(INTSXP, Npart)); no = INTEGER(number); /* copy the pointer */