R/bind.tree.R

   1 ## bind.tree.R (2010-02-12)
   2
   3 ##    Bind Trees
   4
   5 ## Copyright 2003-2010 Emmanuel Paradis
   6
   7 ## This file is part of the R-package `ape'.
   8 ## See the file ../COPYING for licensing issues.
   9
  10 bind.tree <- function(x, y, where = "root", position = 0)
  11 {
  12     nb.tip <- length(x$tip.label)
  13     nb.node <- x$Nnode
  14     ROOT <- nb.tip + 1
  15     if (where == 0 || where == "root")
  16       where <- ROOT
  17     if (position < 0) position <- 0
  18     if (where > nb.tip + nb.node) stop("node number out of range for tree 'x'")
  19     nb.edge <- dim(x$edge)[1]
  20     yHasNoRootEdge <- is.null(y$root.edge)
  21     xHasNoRootEdge <- is.null(x$root.edge)
  22
  23     ## check whether both trees have branch lengths:
  24     wbl <- TRUE
  25     noblx <- is.null(x$edge.length)
  26     nobly <- is.null(y$edge.length)
  27     if (noblx && nobly) wbl <- FALSE
  28     if (xor(noblx, nobly)) {
  29         if (nobly) x$edge.length <- NULL
  30         else y$edge.length <- NULL
  31         wbl <- FALSE
  32         warning("one tree has no branch lengths, they will be ignored")
  33     }
  34
  35     ## To avoid problems with tips or nodes with indentical
  36     ## labels we substitute the one where `y' is grafted:
  37     if (where <= nb.tip) {
  38         Tip.Label.where <- x$tip.label[where]
  39         x$tip.label[where] <- "TheTipWhereToGraftY"
  40     }
  41     if (where > ROOT) {
  42         xHasNoNodeLabel <- TRUE
  43         if (is.null(x$node.label)) {
  44             x$node.label <- paste("NODE", 1:nb.node, sep = "")
  45             x$node.label[where - nb.tip] <- "TheNodeWhereToGraftY"
  46         } else {
  47             Node.Label.where <- x$node.label[where - nb.tip]
  48             x$node.label[where - nb.tip] <- "TheNodeWhereToGraftY"
  49             xHasNoNodeLabel <- FALSE
  50         }
  51     }
  52
  53     ## if we bind `y' under a node or tip of `y', we first
  54     ## adjust the edge lengths if needed
  55     if (position && wbl) {
  56         if (where == ROOT) {
  57             if (xHasNoRootEdge) stop("tree 'x' has no root edge")
  58             if (x$root.edge < position)
  59               stop("argument 'position' is larger than the root edge.")
  60             x$root.edge <- x$root.edge - position
  61         } else {
  62             i <- which(x$edge[, 2] == where)
  63             if (x$edge.length[i] < position)
  64               stop("argument 'position' is larger than the specified edge.")
  65             x$edge.length[i] <- x$edge.length[i] - position
  66         }
  67         y$root.edge <- if (yHasNoRootEdge) position else y$root.edge + position
  68     }
  69
  70     if (is.null(y$root.edge) && where > nb.tip) y$root.edge <- 0
  71
  72     X <- write.tree(x)
  73     Y <- write.tree(y)
  74     Y <- substr(Y, 1, nchar(Y) - 1)
  75
  76     if (where <= nb.tip) {
  77         if (position)
  78           X <- gsub("TheTipWhereToGraftY",
  79                     paste("(", "TheTipWhereToGraftY", ",", Y, ")",
  80                           sep = ""), X)
  81         else
  82           X <- gsub("TheTipWhereToGraftY", Y, X)
  83     }
  84     if (where == ROOT) {
  85         rmvx <- if (xHasNoRootEdge) "\\);$" else ";$"
  86         X <- gsub(rmvx, "", X)
  87         Y <- gsub("^\\(", "", Y)
  88         if (!xHasNoRootEdge) X <- paste("(", X, sep = "")
  89         X <- paste(X, ",", Y, ";", sep = "")
  90     }
  91     if (where > ROOT) {
  92         if (position) {
  93             ## find where is the node in `X':
  94             ## below 19 is: nchar("TheNodeWhereToGraftY") - 1
  95             for (i in 1:nchar(X)) {
  96                 if ("TheNodeWhereToGraftY" == substr(X, i, i + 19))
  97                   break
  98                 i <- i + 1
  99             }
 100             ## now go back to find the left matching parentheses
 101             n.paren <- 1
 102             i <- i - 2
 103             while (n.paren > 0) {
 104                 if (substr(X, i, i) == ")") n.paren <- n.paren + 1
 105                 if (substr(X, i, i) == "(") n.paren <- n.paren - 1
 106                 i <- i - 1
 107             }
 108             ## insert the left parenthesis:
 109             ## here 21 is: nchar("TheNodeWhereToGraftY") + 1
 110             X <- paste(substr(X, 1, i - 1), "(",
 111                        substr(X, i, 21), sep = "")
 112             ## and insert `y':
 113             X <- gsub("TheNodeWhereToGraftY",
 114                       paste("TheNodeWhereToGraftY", ",", Y,
 115                             sep = ""), X)
 116         } else {
 117             xx <- paste(")", "TheNodeWhereToGraftY", sep = "")
 118             X <- gsub(xx, paste(",", Y, xx, sep = ""), X)
 119         }
 120     }
 121     phy <- read.tree(text = X)
 122     ## restore the labels:
 123     if (where <= nb.tip)
 124       phy$tip.label[which(phy$tip.label == "TheTipWhereToGraftY")] <-
 125         Tip.Label.where
 126     if (where > ROOT) {
 127         if (xHasNoNodeLabel) phy$node.label <- NULL
 128         else
 129           phy$node.label[which(phy$node.label == "TheNodeWhereToGraftY")] <-
 130             Node.Label.where
 131     }
 132     phy
 133 }