From: paradis Date: Tue, 28 Feb 2012 12:57:46 +0000 (+0000) Subject: new alex() X-Git-Url: https://git.donarmstrong.com/?p=ape.git;a=commitdiff_plain;h=21426f51c5940cb37f3198a7853ef59743729b85 new alex() git-svn-id: https://svn.mpl.ird.fr/ape/dev/ape@180 6e262413-ae40-0410-9e79-b911bd7a66b7 --- diff --git a/DESCRIPTION b/DESCRIPTION index 8d34b18..7d6ef89 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: ape -Version: 3.0-1 -Date: 2012-02-17 +Version: 3.0-2 +Date: 2012-02-27 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 57a418c..2ebf72e 100644 --- a/NEWS +++ b/NEWS @@ -1,3 +1,13 @@ + CHANGES IN APE VERSION 3.0-2 + + +NEW FEATURES + + o The new function alex (alignment explorator) zooms in a DNA + alignment and opens the result in a new window. + + + CHANGES IN APE VERSION 3.0-1 @@ -17,6 +27,9 @@ BUG FIXES o read.GenBank() has been updated to work with EFetch 2.0. + o unroot() on trees in "pruningwise" order did not respect this + attribute + CHANGES IN APE VERSION 3.0 diff --git a/R/alex.R b/R/alex.R new file mode 100644 index 0000000..ab846a7 --- /dev/null +++ b/R/alex.R @@ -0,0 +1,42 @@ +## alex.R (2012-03-27) + +## Alignment Explorer With Multiple Devices + +## Copyright 2012 Emmanuel Paradis + +## This file is part of the R-package `ape'. +## See the file ../COPYING for licensing issues. + +alex <- function(x, ...) +{ + n <- nrow(x) + s <- ncol(x) + devmain <- dev.cur() + on.exit(dev.set(devmain)) + NEW <- TRUE + cat("Click on two opposite corners of the zone you want to zoom-in. +Right-click to exit.\n") + repeat { + xy <- locator(2) + if (is.null(xy)) break + xy <- lapply(xy, function(x) sort(round(x))) + i1 <- xy$y[1L]; i2 <- xy$y[2L] + j1 <- xy$x[1L]; j2 <- xy$x[2L] + if (i1 > n || j1 > s) cat("Try again!\n") else { + if (i1 <= 0) i1 <- 1L + if (j1 <= 0) j1 <- 1L + if (i2 > n) i2 <- n + if (j2 > s) j2 <- s + if (NEW) { + dev.new() + devsub <- dev.cur() + NEW <- FALSE + } else dev.set(devsub) + image(x[i1:i2, j1:j2], xaxt = "n", ...) + atx <- axTicks(1) + axis(1, atx, labels = (j1:j2)[atx]) + title(sub = paste("From", sQuote(deparse(substitute(x))))) + dev.set(devmain) + } + } +} diff --git a/R/root.R b/R/root.R index 8c1d9b2..c51184d 100644 --- a/R/root.R +++ b/R/root.R @@ -21,43 +21,73 @@ unroot <- function(phy) { if (!inherits(phy, "phylo")) stop('object "phy" is not of class "phylo"') - if (dim(phy$edge)[1] < 3) + N <- dim(phy$edge)[1] + if (N < 3) stop("cannot unroot a tree with less than three edges.") + ## delete FIRST the root.edge (in case this is sufficient to ## unroot the tree, i.e. there is a multichotomy at the root) if (!is.null(phy$root.edge)) phy$root.edge <- NULL if (!is.rooted(phy)) return(phy) - ## We remove one of the edges coming from the root, and - ## eventually adding the branch length to the other one - ## also coming from the root. - ## In all cases, the node deleted is the 2nd one (numbered - ## nb.tip+2 in 'edge'), so we simply need to renumber the - ## nodes by adding 1, except the root (this remains the - ## origin of the tree). - nb.tip <- length(phy$tip.label) - ROOT <- nb.tip + 1L - EDGEROOT <- which(phy$edge[, 1] == ROOT) - ## j: the target where to stick the edge - ## i: the edge to delete - if (phy$edge[EDGEROOT[1], 2] == ROOT + 1L) { - j <- EDGEROOT[2] - i <- EDGEROOT[1] + + n <- length(phy$tip.label) + ROOT <- n + 1L + +### EDGEROOT[1]: the edge to delete +### EDGEROOT[2]: the target where to stick the edge to delete + +### If the tree is in pruningwise order, then the last two edges +### are those connected to the root; the node situated in +### phy$edge[N - 2L, 1L] will be the new root... + + ophy <- attr(phy, "order") + if (!is.null(ophy) && ophy == "pruningwise") { + NEWROOT <- phy$edge[N - 2L, 1L] + EDGEROOT <- c(N, N - 1L) } else { - j <- EDGEROOT[1] - i <- EDGEROOT[2] + +### ... otherwise, we remove one of the edges coming from +### the root, and eventually adding the branch length to +### the other one also coming from the root. +### In all cases, the node deleted is the 2nd one (numbered +### nb.tip+2 in 'edge'), so we simply need to renumber the +### nodes by adding 1, except the root (this remains the +### origin of the tree). + + EDGEROOT <- which(phy$edge[, 1L] == ROOT) + NEWROOT <- ROOT + 1L } - ## This should work whether the tree is in pruningwise or - ## cladewise order. - phy$edge <- phy$edge[-i, ] - nodes <- phy$edge > ROOT # renumber all nodes except the root - phy$edge[nodes] <- phy$edge[nodes] - 1L + + ## make sure EDGEROOT is ordered as described above: + if (phy$edge[EDGEROOT[1L], 2L] != NEWROOT) + EDGEROOT <- EDGEROOT[2:1] + + phy$edge <- phy$edge[-EDGEROOT[1L], ] + + s <- phy$edge == NEWROOT # renumber the new root + phy$edge[s] <- ROOT + + s <- phy$edge > NEWROOT # renumber all nodes greater than the new root + phy$edge[s] <- phy$edge[s] - 1L + if (!is.null(phy$edge.length)) { - phy$edge.length[j] <- phy$edge.length[j] + phy$edge.length[i] - phy$edge.length <- phy$edge.length[-i] + phy$edge.length[EDGEROOT[2L]] <- + phy$edge.length[EDGEROOT[2L]] + phy$edge.length[EDGEROOT[1L]] + phy$edge.length <- phy$edge.length[-EDGEROOT[1L]] } + phy$Nnode <- phy$Nnode - 1L - if (!is.null(phy$node.label)) - phy$node.label <- phy$node.label[-2] + + if (!is.null(phy$node.label)) { + if (NEWROOT == n + 2L) + phy$node.label <- phy$node.label[-1] + else { + lbs <- phy$node.label + tmp <- lbs[NEWROOT - n] + lbs <- lbs[-c(1, NEWROOT)] + phy$node.label <- c(tmp, lbs) + } + } phy } diff --git a/R/vcv.phylo.R b/R/vcv.phylo.R index 294b767..fb3f2a0 100644 --- a/R/vcv.phylo.R +++ b/R/vcv.phylo.R @@ -1,4 +1,4 @@ -## vcv.phylo.R (2012-02-09) +## vcv.phylo.R (2012-02-21) ## Phylogenetic Variance-Covariance or Correlation Matrix @@ -45,18 +45,20 @@ vcv.phylo <- function(phy, model = "Brownian", corr = FALSE, ...) } } - diag(vcv) <- xx[1:n] + diag.elts <- 1 + 0:(n - 1)*(n + 1) + vcv[diag.elts] <- xx[1:n] if (corr) { - ## This is inspired from the code of `cov2cor' (2005-09-08): - M <- vcv - Is <- sqrt(1/M[1 + 0:(n - 1)*(n + 1)]) - vcv[] <- Is * M * rep(Is, each = n) - vcv[1 + 0:(n - 1)*(n + 1)] <- 1 + ## This is inspired from the code of cov2cor (2005-09-08): + Is <- sqrt(1 / vcv[diag.elts]) + ## below 'vcv[] <- ...' has been changed to 'vcv <- ...' + ## which seems to be twice faster for n = 1000 and + ## respects the additional attributes (2012-02-21): + vcv <- Is * vcv * rep(Is, each = n) + vcv[diag.elts] <- 1 } dimnames(vcv)[1:2] <- list(phy$tip.label) - vcv } diff --git a/man/alex.Rd b/man/alex.Rd new file mode 100644 index 0000000..92176e0 --- /dev/null +++ b/man/alex.Rd @@ -0,0 +1,44 @@ +\name{alex} +\alias{alex} +\title{Alignment Explorer With Multiple Devices} +\description{ + This function helps to explore DNA alignments by zooming in. The user + clicks twice defining the opposite corners of the portion which is + extracted and drawned on a new window. +} +\usage{ +alex(x, ...) +} +\arguments{ + \item{x}{an object of class \code{"DNAbin"}.} + \item{\dots}{further arguments to pass to \code{image.DNAbin}.} +} +\details{ + This function works with a DNA alignment (freshly) plotted on an + interactive graphical device (i.e., not a file) with \code{image}. + After calling \code{alex}, the user clicks twice defining a rectangle + in the alignment, then this portion of the alignment is extacted and + plotted on a \emph{new} window. The user can click as many times on + the alignment. The process is stopped by a right-click. If the user + clicks twice outside the alignment, a message ``Try again!'' is + printed. + + Each time \code{alex} is called, the alignment is plotted on a new + window without closing or deleting those possibly already plotted. + + In all cases, the device where \code{x} is plotted is the active + window after the operation. It should \emph{not} be closed during the + whole process. +} +\value{NULL} +\author{Emmanuel Paradis} +\seealso{ + \code{\link{image.DNAbin}}, \code{\link{trex}} +} +\examples{ +\dontrun{ +data(woodmouse) +image(woodmouse) +alex(woodmouse) +}} +\keyword{hplot} diff --git a/man/dbd.Rd b/man/dbd.Rd index f86b97b..ad2e620 100644 --- a/man/dbd.Rd +++ b/man/dbd.Rd @@ -37,9 +37,9 @@ dbdTime(x, birth, death, t, conditional = FALSE, } \details{ These three functions compute the probabilities to observe \code{x} - species starting from a single one after time units \code{t} (assumed - to be continuous). The first one is a short-cut for the second with - \code{mu = 0} and with default values for the two other parameters. + species starting from a single one after time \code{t} (assumed to be + continuous). The first function is a short-cut for the second one with + \code{mu = 0} and with default values for the two other arguments. \code{dbdTime} is for time-varying \code{lambda} and \code{mu} specified as \R functions. @@ -52,7 +52,7 @@ dbdTime(x, birth, death, t, conditional = FALSE, The returned value is, logically, zero for values of \code{x} out of range, i.e., negative or zero for \code{dyule} or if \code{conditional = TRUE}. However, it is not checked if the values of \code{x} are - non-integers and the probabilities are computed and returned. + positive non-integers and the probabilities are computed and returned. The details on the form of the arguments \code{birth}, \code{death}, \code{BIRTH}, \code{DEATH}, and \code{fast} can be found in the links diff --git a/man/trex.Rd b/man/trex.Rd index e05ed42..7b06b6b 100644 --- a/man/trex.Rd +++ b/man/trex.Rd @@ -2,7 +2,7 @@ \alias{trex} \title{Tree Explorer With Multiple Devices} \description{ - This functions requires a plotted tree: the user is invited to click + This function requires a plotted tree: the user is invited to click close to a node and the corresponding subtree (or clade) is plotted on a new window. } diff --git a/src/dist_dna.c b/src/dist_dna.c index 289be08..ab630a9 100644 --- a/src/dist_dna.c +++ b/src/dist_dna.c @@ -93,14 +93,6 @@ void distDNA_indel(unsigned char *x, int *n, int *s, double *d) #define X(i, j) i - 1 + *n * (j - 1) -#define DINDEX(i, j) n * (i - 1) - i*(i - 1)/2 + j - i - 1 - -int give_index(int i, int j, int n) -{ - if (i > j) return(DINDEX(j, i)); - else return(DINDEX(i, j)); -} - void distDNA_indelblock(unsigned char *x, int *n, int *s, double *d) { int i1, i2, s1, s2, target, N, start_block, end_block;