Specifying colours in phylogenetic ancestral state reconstruction

Question

I am trying to specify the colours used in my ancestral state reconstruction analysis so that I can then use the same colours in other graphs the same data. However, when I used nodelabels from ape the colours I specify do not line up with what is displayed. Does anyone know why this is? And how to work around it? See example of the problem below:

rm(list = ls())

library(phytools)
library(ape)
library(RColorBrewer)

set.seed(1237)

tree = rtree(50)

plot(tree)

variable = rTraitDisc(tree, k = 5)
names(variable) = tree$tip.label

cols = data.frame(type = levels(variable),
                    color = I(brewer.pal(nlevels(variable), name = 'Set1')))
cols_vector = cols$color
names(cols_vector) = cols$type

fit = ace(x = variable, phy = tree, type = 'discrete')

nodelabels(node=1:tree$Nnode+Ntip(tree),
           pie=fit$lik.anc, piecol=cols_vector)
tiplabels(pie=to.matrix(variable,sort(unique(variable))),piecol=cols_vector)
tiplabels(text = variable)

Which gives the image:

Most of the colours line up ok, but E is displayed as purple but has the code #FF7F00 (which is orange). Which is exemplified by this code:

t = table(variable)
barplot(t, col = cols_vector[match(names(t), names(cols_vector))])

Any advice would be appreciated.

Answer 1

The problem is in command to.matrix(variable,sort(unique(variable))). Given your seed, no tips have value D and thus if you construct a matrix with columns representing unique values, only the first four colours from your palette will be used.

head(to.matrix(variable,sort(unique(variable))))
    A B C E
t47 1 0 0 0
t14 1 0 0 0
t2  1 0 0 0
t19 1 0 0 0
t29 1 0 0 0
t5  1 0 0 0

Replacing the seq = argument with the vector of all your possible trait values should solve the problem with non-matching colours.

tiplabels(pie = to.matrix(x = variable, seq = LETTERS[1:5]), piecol = cols_vector)