How can I plot a tree (and squirrels) in R?

Question

Here is my tree:

tree = data.frame(branchID = c(1,11,12,111,112,1121,1122), length = c(32, 21, 19, 5, 12, 6, 2))

> tree
  branchID length
1        1     32
2       11     21
3       12     19
4      111      5
5      112     12
6     1121      6
7     1122      2

This tree is in 2D and is made of branches. Each branch has an ID. 1 is the trunk. Then the trunk bifurcate into two branches, 11 on the left and 12 on the right. 11 bifurcates as well in the branches called 111 (going toward the left) and 112 (going toward the right). etc.. Each branch has a certain length.

On this tree there are squirrels:

squirrels = data.frame(branchID = c(1,11,1121,11,111), PositionOnBranch = c(23, 12, 4, 2, 1), name=c("FluffyTail", "Ginger", "NutCracker", "SuperSquirrel", "ChipnDale"))

> squirrels
  branchID PositionOnBranch          name
1        1               23    FluffyTail
2       11               12        Ginger
3     1121                4    NutCracker
4       11                2 SuperSquirrel
5      111                1     ChipnDale

Each squirrel is found on a specific branch. For example the FluffyTail is on the trunk at position 23 (the total length of the trunk being 32). ChipnDale is on the branch 111 at position 1 (the total length of the branch 111 is 5). The position is taken relatively to the lower extremity of the branch.

How can I plot my tree and my squirrels?

Maybe also consider some modification of `rpart` plotting to suit your tree structure? — thelatemail, Jan 27 '15 at 06:09

jbaums · Accepted Answer · 2015-08-26T10:09:47.953

I put a bit more thought/time into this, and have packaged up some horticultural functions in package trees, here.

With trees, you can:

generate a random tree design (a random seed, so to speak) with seed();
sow the seed to give rise to a magnificent tree with germinate();
add randomly-located leaves (or squirrels) with foliate();
add squirrels (for example) to specified locations with squirrels(); and
prune() the tree.

# Install the package and set the RNG state
devtools::install_github('johnbaums/trees')
set.seed(1)

Let's fertilise a seed and grow a tree

# Create a tree seed    
s <- seed(70, 10, min.branch.length=0, max.branch.length=4,
          min.trunk.height=5, max.trunk.height=8)

head(s, 10)

#       branch    length
# 1          0 6.3039785
# 2          L 2.8500587
# 3         LL 1.5999775
# 4        LLL 1.3014086
# 5       LLLL 3.0283486
# 6      LLLLL 0.8107690
# 7     LLLLLR 2.8444849
# 8    LLLLLRL 0.4867677
# 9   LLLLLRLR 0.9819541
# 10 LLLLLRLRR 0.5732175

# Germinate the seed
g <- germinate(s, col='peachpuff4')

enter image description here

And add some leaves

leafygreens <- colorRampPalette(paste0('darkolivegreen', c('', 1:4)))(100)
foliate(g, 5000, 4, pch=24:25, col=NA, cex=1.5, bg=paste0(leafygreens, '30'))

enter image description here

Or some squirrels

plot(g, col='peachpuff4')
squirrels(g, 
          branches=c("LLLLRRRL", "LRLRR", "LRRLRLLL", "LRRRLL", "RLLLLLR", 
                     "RLLRL", "RLLRRLRR", "RRRLLRL", "RRRLLRR", "RRRRLR"),
          pos=c(0.22, 0.77, 0.16, 0.12, 0.71, 0.23, 0.18, 0.61, 0.8, 2.71),
          pch=20, cex=2.5)

enter image description here

Plotting @Remi.b's tree and squirrels

g <- germinate(list(trunk.height=32, 
                   branches=c(1, 2, 11, 12, 121, 122),
                   lengths=c(21, 19, 5, 12, 6, 2)), 
              left='1', right='2', angle=40)

xy <- squirrels(g, c(0, 1, 121, 1, 11), pos=c(23, 12, 4, 2, 1), 
               left='1', right='2', pch=21, bg='white', cex=3, lwd=2)
text(xy$x, xy$y, labels=seq_len(nrow(xy)), font=2)
legend('bottomleft', bty='n',
      legend=paste(seq_len(nrow(xy)), 
                   c('FluffyTail', 'Ginger', 'NutCracker', 'SuperSquirrel', 
                     'ChipnDale'), sep='. '))

enter image description here

EDIT:

Following @baptiste's hot tip about @ScottChamberlain's rphylopic package, it's time to upgrade those dots to squirrels (though they may resemble coffee beans).

library(rphylopic)
s <- seed(50, 10, min.branch.length=0, max.branch.length=5,
          min.trunk.height=5, max.trunk.height=8)
g <- germinate(s, trunk.width=15, col='peachpuff4')
leafygreens <- colorRampPalette(paste0('darkolivegreen', c('', 1:4)))(100)
foliate(g, 2000, 4, pch=24:25, col=NA, cex=1.2, bg=paste0(leafygreens, '50'))
xy <- foliate(g, 2, 2, 4, xy=TRUE, plot=FALSE)

# snazzy drop shadow
add_phylopic_base(
    image_data("5ebe5f2c-2407-4245-a8fe-397466bb06da", size = "64")[[1]], 
    1, xy$x, xy$y, ysize = 2.3, col='black')
add_phylopic_base(
    image_data("5ebe5f2c-2407-4245-a8fe-397466bb06da", size = "64")[[1]], 
    1, xy$x, xy$y, ysize = 2, col='darkorange3')

Ok this is insane! Did you actually created a package for my needs?! These trees look soo cool. Thank you so much. — Remi.b, Feb 12 '15 at 05:18
@Remi.b: I created a package because it was a fun question! `:-)` — jbaums, Feb 12 '15 at 05:25
This might be the greatest answer I have ever seen in the R tag. Brilliant! — thelatemail, May 06 '15 at 02:23
you should try adding squirrels with this package https://github.com/sckott/rphylopic — baptiste, Aug 25 '15 at 20:07
`grid::grid.raster(rphylopic::image_data("5ebe5f2c-2407-4245-a8fe-397466bb06da", size = "thumb")[[1]])` to plot one — baptiste, Aug 26 '15 at 01:12

score 12 · Answer 2 · answered Jan 28 '15 at 17:08

I probably over-thought this, but... squirrels.

get.coords <- function(a, d, x0, y0) {
  a <- ifelse(a <= 90, 90 - a, 450 - a)
  data.frame(x = x0 + d * cos(a / 180 * pi), 
             y = y0+ d * sin(a / 180 * pi))
}


tree$angle <- sapply(gsub(2, '+45', gsub(1, '-45', tree$branchID)), 
                     function(x) eval(parse(text=x)))
tree$tipy <- tree$tipx <- tree$basey <- tree$basex <- NA

for(i in seq_len(nrow(tree))) {
  if(tree$branchID[i] == 0) {
    tree$basex[i] <- tree$basey[i] <- tree$tipx[i] <- 0
    tree$tipy[i] <- tree$length[i]
    next
  } else if(tree$branchID[i] %in% 1:2) {
    parent <- 0
  } else {
    parent <- substr(tree$branchID[i], 1, nchar(tree$branchID[i])-1)
  }
  tree$basex[i] <- tree$tipx[which(tree$branchID==parent)]
  tree$basey[i] <- tree$tipy[which(tree$branchID==parent)]
  tip <- get.coords(tree$angle[i], tree$length[i], tree$basex[i], tree$basey[i])
  tree$tipx[i] <- tip[, 1]
  tree$tipy[i] <- tip[, 2]
}  

squirrels$nesty <- squirrels$nestx <- NA
for (i in seq_len(nrow(squirrels))) {
  b <- tree[tree$branchID == squirrels$branchID[i], ]
  nest <- get.coords(b$angle, squirrels$PositionOnBranch[i], b$basex, b$basey)
  squirrels$nestx[i] <- nest[1]
  squirrels$nesty[i] <- nest[2]
}

And now we plot.

plot.new()
plot.window(xlim=range(tree$basex, tree$tipx), 
            ylim=range(tree$basey, tree$tipy), asp=1)
with(tree, segments(basex, basey, tipx, tipy, lwd=pmax(10/nchar(branchID), 1)))
points(squirrels[, c('nestx', 'nesty')], pch=21, cex=3, bg='white', lwd=2)
text(squirrels[, c('nestx', 'nesty')], labels=seq_len(nrow(squirrels)), font=2)
legend('bottomleft', legend=paste(seq_len(nrow(squirrels)), squirrels$name), bty='n')

squizzerl

And for kicks we will simulate a bigger tree (and put some apples on it like in Farmville):

twigs <- replicate(50, paste(rbinom(5, 1, 0.5) + 1, collapse=''))
branches <- sort(unique(c(sapply(twigs, function(x) sapply(seq_len(nchar(x)), function(y) substr(x, 1, y))))))
tree <- data.frame(branchID=c(0, branches), length=c(30, sample(10, length(branches), TRUE)), 
                   stringsAsFactors=FALSE)


tree$angle <- sapply(gsub(2, '+45', gsub(1, '-45', tree$branchID)), 
                     function(x) eval(parse(text=x)))
tree$tipy <- tree$tipx <- tree$basey <- tree$basex <- NA

for(i in seq_len(nrow(tree))) {
  if(tree$branchID[i] == 0) {
    tree$basex[i] <- tree$basey[i] <- tree$tipx[i] <- 0
    tree$tipy[i] <- tree$length[i]
    next
  } else if(tree$branchID[i] %in% 1:2) {
    parent <- 0
  } else {
    parent <- substr(tree$branchID[i], 1, nchar(tree$branchID[i])-1)
  }
  tree$basex[i] <- tree$tipx[which(tree$branchID==parent)]
  tree$basey[i] <- tree$tipy[which(tree$branchID==parent)]
  tip <- get.coords(tree$angle[i], tree$length[i], tree$basex[i], tree$basey[i])
  tree$tipx[i] <- tip[, 1]
  tree$tipy[i] <- tip[, 2]
}  

plot.new()
plot.window(xlim=range(tree$basex, tree$tipx), 
            ylim=range(tree$basey, tree$tipy), asp=1)
par(mar=c(0, 0, 0, 0))
with(tree, segments(basex, basey, tipx, tipy, lwd=pmax(20/nchar(branchID), 1)))

apple_branches <- sample(branches, 10)
sapply(apple_branches, function(x) {
  b <- tree[tree$branchID == x, ]
  apples <- get.coords(b$angle, runif(sample(2, 1), 0, b$length), b$basex, b$basey)
  points(apples, pch=20, col='tomato2', cex=2)
})

enter image description here

I forgot to mention that I took the liberty of giving the trunk `branchID=0`, and dropping the first `1` from everything else (so the first two boughs splitting of the trunk are `1` and `2`. — jbaums, Jan 28 '15 at 17:19
Ha ha this is awesome! Great job. I actually just wrote a code to produce tree stochastically. But your code seems more performant than mine and the plot is awesome. I will probably forget mine and take yours. Can you please just point me how, in your code, you define the distribution from which the branch lengths are taken and where the nb of branches is defined and how the apples are randomly placed on the tree. Thanks a lot! — Remi.b, Jan 28 '15 at 19:28
@Remi.b I'll go through this answer to tidy it up and answer your questions when I get a moment. Briefly, though, I took a random (uniform) sample of fifty 6th tier branches (`paste(rbinom(5, 1, 0.5) + 1, collapse='')` gives one such branch, and I replicated this 50 times), and reduced this sample to a unique set. The number of branches is not fully controllable with this approach, but is manipulated to some extent by `replicate(50, ...)`. I then added all the parent branches to the set, including the trunk. — jbaums, Jan 28 '15 at 22:24
For apples, I chose a random sample of 10 branches (`apple_branches <- sample(branches, 10)`), and for each of these branches I randomly selected whether to place 1 or 2 apples, with location of these along the branch random between 0 (the start of the branch) and `b$length` (the end of the branch). — jbaums, Jan 28 '15 at 22:26

score 10 · Answer 3 · edited May 23 '17 at 11:46

Well, you could convert your data to define a "tree" as defined by the ape package. Here's a function that can convert your data.frame to the correct format.

library(ape)

to.tree <- function(dd) {
    dd$parent <- dd$branchID %/% 10

    root <- subset(dd, parent==0)
    dd <- subset(dd, parent!=0)

    ids <- unique(c(dd$parent, dd$branchID))
    tip <- !(ids %in% dd$parent)
    lvl <- ids[order(!tip, ids)]
    edg <- sapply(dd[,c("parent","branchID")], 
        function(x) as.numeric(factor(x, levels=lvl)))

    x<-list(
        edge=edg,
        edge.length=dd$length,
        tip.label=head(lvl, sum(tip)),
        node.label=tail(lvl, length(tip)-sum(tip)),
        Nnode = length(tip)-sum(tip),
        root.edge=root$length[1]
    )
    class(x)<-"phylo"
    reorder(x)    
}

Then we can plot it somewhat easily

xx <- to.tree(tree)
plot(xx, show.node.label=TRUE, root.edge=TRUE)

Now, if we want to add the squirrel information, we need to know where each branch is located. I'm going to borrow getphylo_x and getphylo_y from this answer. Then I can run

sx<-Vectorize(getphylo_x, "node")(xx, as.character(squirrels$branchID)) -
    tree$length[match(squirrels$branchID, tree$branchID)] +
    squirrels$PositionOnBranch
sy<-Vectorize(getphylo_y, "node")(xx, as.character(squirrels$branchID))

points(sx,sy)
text(sx,sy, squirrels$name, pos=3)

to add the squirrel information to the plot. The final result is

enter image description here

It's not perfect but it's not a bad start.

Looks really nice! +1. Thank you. However, I fail to reproduce your tree. The function `reorder` in `to_tree` throws an error. I know nothing about plotting phylogenies. Do you know if we can easily modify your graph so that branches are straight (does not bend) and leave the bifurcation with some angle (say, the new branches form an angle of 80 degrees for example)? — Remi.b, Jan 27 '15 at 14:28
I forgot to include the line, but make sure you load the `ape` library. That library does have other options for drawing the tree, but none of the code I wrote to place the squirrels would work with diagonal branches. — MrFlick, Jan 27 '15 at 17:03
Yes, I loaded the `ape` package. However, the `reorder` function needs two arguments, it returns `Error in tapply(X = X, INDEX = x, FUN = FUN, ...) : argument "X" is missing, with no default` when it is run in the `to.tree` function. I tried with `reorder.phylo` but then the plot fails throwing the message `'x' is a list, but does not have components 'x' and 'y'`. Is there something I get it wrong? (I am using `R 3.1.2 GUI 1.65 Mavericks build`). — Remi.b, Jan 27 '15 at 17:24
Oh alright. I should have spotted the mistake by myself. Thanks a lot! — Remi.b, Jan 27 '15 at 17:42

score 5 · Answer 4 · answered Jan 27 '15 at 06:06

The reshaping of this might take a while, but this is broadly possible. E.g., rejigging your data representation so it looks like:

library(igraph)
dat <- read.table(text="1 1n2
1n2 1.1
1n2 1.2
1.1 1.1.1
1.1 1.1.2
1.1.2 1.1.2.1
1.1.2 1.1.2.2",header=FALSE)

g <- graph.data.frame(dat)
tkplot(g)

And manually moving the tree parts around in tkplot, you can get:

enter image description here

Doing this automatically is a whole different story admittedly.

Christoph Glur · Answer 5 · 2016-05-15T06:17:56.570

A version that supports trees with more than two branches. A bit work is required to convert to a data.tree structure, and to add the squirrels to it. But once you're there, the plotting is straight forward.

df <- data.frame(branchID = c(1,11,12,13, 14, 111,112,1121,1122), length = c(32, 21, 12, 8, 19, 5, 12, 6, 2))
squirrels <- data.frame(branchID = c(1,11,1121,11,111), PositionOnBranch = c(23, 12, 4, 2, 1), squirrel=c("FluffyTail", "Ginger", "NutCracker", "SuperSquirrel", "ChipnDale"), stringsAsFactors = FALSE)

library(magrittr)

#derive pathString from branchID, so we can convert it to data.tree structure
df$branchID %>%
  as.character %>%
  sapply(function(x) strsplit(x, split = "")) %>%
  sapply(function(x) paste(x, collapse = "/")) ->
  df$pathString

df$type <- "branch"

library(data.tree)

tree <- FromDataFrameTable(df)

#climb, little squirrels!
for (i in 1:nrow(squirrels)) {
  squirrels[i, 'branchID'] %>%
    as.character %>%
    strsplit(split = "") %>%
    extract2(1) %>%
    extract(-1) -> path
  if (length(path) > 0) branch <- tree$Climb(path)
  else branch <- tree
  #actually, we add the squirrels as branches to our tree
  #What a symbiotic coexistence!
  #advantage: Our SetCoordinates can be re-used as is
  #disadvantage: may be confusing, and it requires us
  #to do some filtering later
  branch$AddChild(squirrels[i, 'squirrel'],
                 length = squirrels[i, 'PositionOnBranch'],
                 type = "squirrel")
}



SetCoordinates <- function(node, branch) {
  if (branch$isRoot) {
    node$x0 <- 0
    node$y0 <- 0
  } else {
    node$x0 <- branch$parent$x1
    node$y0 <- branch$parent$y1
  }

  #let's hope our squirrels didn't flunk in trigonometry ;-)
  angle <- branch$position / (sum(Get(branch$siblings, "type") == "branch") + 2)
  x <- - node$length * cospi(angle)
  y <- sqrt(node$length^2 - x^2)
  node$x1 <- node$x0 + x
  node$y1 <- node$y0 + y
}

#let it grow!
tree$Do(function(node) {
          SetCoordinates(node, node)
          node$lwd <- 10 * (node$root$height - node$level + 1) / node$root$height
        }, filterFun = function(node) node$type == "branch")
tree$Do(function(node) SetCoordinates(node, node$parent), filterFun = function(node) node$type == "squirrel")

Looking at the data:

print(tree, "type", "length", "x0", "y0", "x1", "y1")

This prints like so:

                    levelName     type length        x0       y0         x1       y1
1  1                            branch     32   0.00000  0.00000   0.000000 32.00000
2   ¦--1                        branch     21   0.00000 32.00000 -16.989357 44.34349
3   ¦   ¦--1                    branch      5 -16.98936 44.34349 -19.489357 48.67362
4   ¦   ¦   °--ChipnDale      squirrel      1 -16.98936 44.34349 -17.489357 45.20952
5   ¦   ¦--2                    branch     12 -16.98936 44.34349 -10.989357 54.73580
6   ¦   ¦   ¦--1                branch      6 -10.98936 54.73580 -13.989357 59.93195
7   ¦   ¦   ¦   °--NutCracker squirrel      4 -10.98936 54.73580 -12.989357 58.19990
8   ¦   ¦   °--2                branch      2 -10.98936 54.73580  -9.989357 56.46785
9   ¦   ¦--Ginger             squirrel     12   0.00000 32.00000  -9.708204 39.05342
10  ¦   °--SuperSquirrel      squirrel      2   0.00000 32.00000  -1.618034 33.17557
11  ¦--2                        branch     12   0.00000 32.00000  -3.708204 43.41268
12  ¦--3                        branch      8   0.00000 32.00000   2.472136 39.60845
13  ¦--4                        branch     19   0.00000 32.00000  15.371323 43.16792
14  °--FluffyTail             squirrel     23   0.00000  0.00000   0.000000 23.00000

Once we're here, plotting is also easy:

plot(c(min(tree$Get("x0")), max(tree$Get("x1"))),
     c(min(tree$Get("y0")), max(tree$Get("y1"))),
     type='n', asp=1, axes=FALSE, xlab='', ylab='')

tree$Do(function(node) segments(node$x0, node$y0, node$x1, node$y1, lwd = node$lwd),
        filterFun = function(node) node$type == "branch")

tree$Do(function(node) {
          points(node$x1, node$y1, lwd = 8, col = "saddlebrown")
          text(node$x1, node$y1, labels = node$name, pos = 2, cex = 0.7)
        },
        filterFun = function(node) node$type == "squirrel")

How can I plot a tree (and squirrels) in R?

5 Answers5

Linked