From: paradis Date: Thu, 10 Nov 2011 10:23:02 +0000 (+0000) Subject: a few changes... X-Git-Url: https://git.donarmstrong.com/?p=ape.git;a=commitdiff_plain;h=a02ce8c6e9fd80d3d7a749cc24699366fb8e54b6 a few changes... git-svn-id: https://svn.mpl.ird.fr/ape/dev/ape@172 6e262413-ae40-0410-9e79-b911bd7a66b7 --- diff --git a/DESCRIPTION b/DESCRIPTION index 5882a71..2ac8960 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: ape -Version: 2.8 -Date: 2011-10-20 +Version: 2.8-1 +Date: 2011-11-10 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 diff --git a/NEWS b/NEWS index 17f28e2..39b3905 100644 --- a/NEWS +++ b/NEWS @@ -1,3 +1,19 @@ + CHANGES IN APE VERSION 2.8-1 + + +BUG FIXES + + o mantel.test() could return P-values > 1 with the default + two-tailed test. + + +OTHER CHANGES + + o The code of yule() has been simplified and is now much faster for + big trees. + + + CHANGES IN APE VERSION 2.8 @@ -22,6 +38,11 @@ NEW FEATURES the new function ltt.plot.coords. +BUG FIXES + + o prop.part() crashed if some trees had some multichotomies. + + CHANGES IN APE VERSION 2.7-3 diff --git a/R/mantel.test.R b/R/mantel.test.R index 9704152..a23bc32 100644 --- a/R/mantel.test.R +++ b/R/mantel.test.R @@ -1,4 +1,4 @@ -## mantel.test.R (2011-06-22) +## mantel.test.R (2011-11-10) ## Mantel Test for Similarity of Two Matrices @@ -18,20 +18,19 @@ mant.zstat <- function(m1, m2) sum(lower.triang(m1 * m2)) lower.triang <- function(m) { - d <- dim(m) - if (d[1] != d[2]) print("Warning: non-square matrix") + if (!diff(dim(m))) print("Warning: non-square matrix") m[col(m) <= row(m)] } -mantel.test <- function (m1, m2, nperm = 1000, graph = FALSE, - alternative = "two.sided", ...) +mantel.test <- function(m1, m2, nperm = 1000, graph = FALSE, + alternative = "two.sided", ...) { alternative <- match.arg(alternative, c("two.sided", "less", "greater")) n <- nrow(m1) realz <- mant.zstat(m1, m2) nullstats <- replicate(nperm, mant.zstat(m1, perm.rowscols(m2, n))) pval <- switch(alternative, - "two.sided" = 2 * sum(abs(nullstats) > abs(realz)), + "two.sided" = 2 * min(sum(nullstats > realz), sum(nullstats < realz)), "less" = sum(nullstats < realz), "greater" = sum(nullstats > realz)) pval <- pval / nperm diff --git a/R/yule.R b/R/yule.R index 041bda5..abd719c 100644 --- a/R/yule.R +++ b/R/yule.R @@ -1,11 +1,11 @@ -## yule.R (2009-06-08) +## yule.R (2011-11-03) ## Fits Yule Model to a Phylogenetic Tree ## yule: standard Yule model (constant birth rate) ## yule.cov: Yule model with covariates -## Copyright 2003-2009 Emmanuel Paradis +## Copyright 2003-2011 Emmanuel Paradis ## This file is part of the R-package `ape'. ## See the file ../COPYING for licensing issues. @@ -13,18 +13,17 @@ yule <- function(phy, use.root.edge = FALSE) { if (!is.binary.tree(phy)) - stop("tree must be dichotomous to fit the Yule model.") - bt <- rev(sort(branching.times(phy))) # branching times from past to present - ni <- cumsum(rev(table(bt))) + 1 + stop("tree must be dichotomous to fit the Yule model.") + X <- sum(phy$edge.length) nb.node <- phy$Nnode - if (!is.null(phy$root.edge) && use.root.edge) { - X <- X + phy$root.edge - ni <- c(1, ni) - } else nb.node <- nb.node - 1 + + if (!is.null(phy$root.edge) && use.root.edge) X <- X + phy$root.edge + else nb.node <- nb.node - 1 + lambda <- nb.node/X se <- lambda/sqrt(nb.node) - loglik <- -lambda*X + lfactorial(phy$Nnode) + nb.node*log(lambda) + loglik <- -lambda * X + lfactorial(phy$Nnode) + nb.node * log(lambda) obj <- list(lambda = lambda, se = se, loglik = loglik) class(obj) <- "yule" obj diff --git a/Thanks b/Thanks index 0d9b473..98cf253 100644 --- a/Thanks +++ b/Thanks @@ -6,7 +6,7 @@ comments, or bug reports: thanks to all of you! Significant bug fixes were provided by Cécile Ané, Jeremy Beaulieu, James Bullard, Otto Cordero, Éric Durand, Olivier François, Rich -FitzJohn, , Jos Käfer, Bret Larget, Naim Matasci, Nick Matzke, Michael +FitzJohn, Jos Käfer, Bret Larget, Naim Matasci, Nick Matzke, Michael Phelan, Elizabeth Purdom, Dan Rabosky, Filipe Vieira, Tim Wallstrom, Li-San Wang, Yan Wong, Peter Wragg, Janet Young, and Jinlong Zhang. Contact me if I forgot someone. diff --git a/man/bionj.Rd b/man/bionj.Rd index 261542a..18b06b4 100644 --- a/man/bionj.Rd +++ b/man/bionj.Rd @@ -35,9 +35,10 @@ x <- c(7, 8, 11, 13, 16, 13, 17, 5, 8, 10, 13, 10, 14, 5, 7, 10, 7, 11, 8, 11, 8, 12, 5, 6, 10, 9, 13, 8) M <- matrix(0, 8, 8) -M[row(M) > col(M)] <- x -M[row(M) < col(M)] <- x -rownames(M) <- colnames(M) <- 1:8 +M[lower.tri(M)] <- x +M <- t(M) +M[lower.tri(M)] <- x +dimnames(M) <- list(1:8, 1:8) tr <- bionj(M) plot(tr, "u") ### a less theoretical example diff --git a/man/fastme.Rd b/man/fastme.Rd index d29d305..9174b48 100644 --- a/man/fastme.Rd +++ b/man/fastme.Rd @@ -42,9 +42,10 @@ x <- c(7, 8, 11, 13, 16, 13, 17, 5, 8, 10, 13, 10, 14, 5, 7, 10, 7, 11, 8, 11, 8, 12, 5, 6, 10, 9, 13, 8) M <- matrix(0, 8, 8) -M[row(M) > col(M)] <- x -M[row(M) < col(M)] <- x -rownames(M) <- colnames(M) <- 1:8 +M[lower.tri(M)] <- x +M <- t(M) +M[lower.tri(M)] <- x +dimnames(M) <- list(1:8, 1:8) tr <- fastme.bal(M) plot(tr, "u") ### a less theoretical example diff --git a/man/nj.Rd b/man/nj.Rd index 612d780..88d2258 100644 --- a/man/nj.Rd +++ b/man/nj.Rd @@ -35,9 +35,10 @@ x <- c(7, 8, 11, 13, 16, 13, 17, 5, 8, 10, 13, 10, 14, 5, 7, 10, 7, 11, 8, 11, 8, 12, 5, 6, 10, 9, 13, 8) M <- matrix(0, 8, 8) -M[row(M) > col(M)] <- x -M[row(M) < col(M)] <- x -rownames(M) <- colnames(M) <- 1:8 +M[lower.tri(M)] <- x +M <- t(M) +M[lower.tri(M)] <- x +dimnames(M) <- list(1:8, 1:8) tr <- nj(M) plot(tr, "u") ### a less theoretical example diff --git a/src/bipartition.c b/src/bipartition.c index a2b34d3..884fc6c 100644 --- a/src/bipartition.c +++ b/src/bipartition.c @@ -1,4 +1,4 @@ -/* bipartition.c 2011-06-23 */ +/* bipartition.c 2011-10-22 */ /* Copyright 2005-2011 Emmanuel Paradis, and 2007 R Development Core Team */ @@ -189,6 +189,11 @@ SEXP prop_part(SEXP TREES, SEXP nbtree, SEXP keep_partitions) /* We start on the 2nd tree: */ for (k = 1; k < Ntree; k++) { + +/* in case there are trees with multichotomies: */ + nbnode = getListElement(VECTOR_ELT(TREES, k), "Nnode"); + Nnode = INTEGER(nbnode)[0]; + PROTECT(bp = bipartition(getListElement(VECTOR_ELT(TREES, k), "edge"), nbtip, nbnode)); for (i = 1; i < Nnode; i++) {