R/evonet.R

   1 ## evonet.R (2011-06-09)
   2
   3 ##   Evolutionary Networks
   4
   5 ## Copyright 2011 Emmanuel Paradis
   6
   7 ## This file is part of the R-package `ape'.
   8 ## See the file ../COPYING for licensing issues.
   9
  10 evonet <- function(phy, from, to = NULL)
  11 {
  12     if (!inherits(phy, "phylo"))
  13         stop('object "phy" is not of class "phylo".')
  14     if (!is.rooted(phy))
  15         warning("the tree is unrooted")
  16     x <- phy
  17
  18     if (is.null(to)) {
  19         if (is.data.frame(from))
  20             from <- as.matrix(from)
  21         if (!is.matrix(from))
  22             stop("'from' must be a matrix or a data frame if 'to' is not given")
  23         if (ncol(from) > 2) {
  24             warning("'from' has more than two columns: only the first two will be used.")
  25             ret <- from[, 1:2]
  26         } else if (ncol(from) < 2) {
  27             stop("'from' must have at least two columns")
  28         } else ret <- from
  29     } else {
  30         from <- as.vector(from)
  31         to <- as.vector(to)
  32         if (length(from) != length(to))
  33             stop("'from' and 'to' not of the same length after coercing as vectors")
  34         ret <- cbind(from, to)
  35     }
  36
  37     ## check that values are not out of range:
  38     storage.mode(ret) <- "integer"
  39     if (any(is.na(ret)))
  40         stop("some values are NA's after coercing as integers")
  41     if (any(ret < 0) || any(ret > Ntip(phy) + phy$Nnode))
  42         stop("some values are out of range")
  43
  44     x$reticulation <- ret
  45     class(x) <- c("evonet", "phylo")
  46     x
  47 }
  48
  49 as.phylo.evonet <- function(x, ...)
  50 {
  51     x$reticulation <- NULL
  52     class(x) <- "phylo"
  53     x
  54 }
  55
  56 plot.evonet <- function(x, col = "blue", lty = 1, lwd = 1, alpha = 0.5,
  57                         arrows = 0, arrow.type = "classical", ...)
  58 {
  59     plot.phylo(as.phylo(x), ...)
  60     ape::edges(x$reticulation[, 1], x$reticulation[, 2],
  61                col = rgb(t(col2rgb(col)), alpha = 255 * alpha,
  62                maxColorValue = 255), lty = lty, lwd = lwd,
  63                arrows = arrows, type = arrow.type)
  64 }
  65
  66 as.networx.evonet <- function(x, weight = NA, ...)
  67 {
  68     if (any(x$reticulation <= Ntip(x)))
  69         stop("some tips are involved in reticulations: cannot convert to \"networx\"")
  70     x <- reorder(x, "p")
  71     ned <- Nedge(x)
  72     nrt <- nrow(x$reticulation)
  73     x$edge <- rbind(x$edge, x$reticulation)
  74     colnames(x$edge) <- c("oldNodes", "newNodes")
  75     x$reticulation <- NULL
  76     x$edge.length <- c(x$edge.length, rep(weight, length.out = nrt))
  77     x$split <- c(1:ned, 1:nrt)
  78     class(x) <- c("networx", "phylo")
  79     x
  80 }
  81
  82 as.network.evonet <- function(x, directed = TRUE, ...)
  83 {
  84     class(x) <- NULL
  85     x$edge <- rbind(x$edge, x$reticulation)
  86     as.network.phylo(x, directed = directed, ...)
  87 }
  88
  89 as.igraph.evonet <- function(x, directed = TRUE, use.labels = TRUE, ...)
  90 {
  91     class(x) <- NULL
  92     x$edge <- rbind(x$edge, x$reticulation)
  93     as.igraph.phylo(x, directed = directed, use.labels = use.labels, ...)
  94 }
  95
  96 print.evonet <- function(x, ...)
  97 {
  98     nr <- nrow(x$reticulation)
  99     cat("\n    Evolutionary network with", nr, "reticulation")
 100     if (nr > 1) cat("s")
 101     cat("\n\n               --- Base tree ---")
 102     print.phylo(as.phylo(x))
 103 }