1
$\begingroup$

I was looking to overlay annotation on my phylogeny made,which was solved here question asked where i was looking to compare two different phylogeny which i could do.

Now i want to label the tips into various groups such as if its human or canine etc.One of the possible solution i found is this one not sure how to implement it.

alignment file annotation file

library(dendextend)
library(seqinr)
library(phytools)
library(phangorn)

a<-read.alignment("clean_dup_align_fast.fas", format="fasta") 
a.phydat<-as.phyDat(a)
dist.a.phydat<-dist.dna(as.DNAbin(a.phydat))
upgma.a<-upgma(dist.a.phydat)
parsimony(upgma.a,a.phydat)
pars.a <- optim.parsimony(upgma.a, a.phydat)
pars.a<-acctran(pars.a, a.phydat)
pars.a.rooted<-root(pars.a, outgroup="AAA64460", resolve.root=T)
pars.a.rooted.dd<-as.dendrogram(force.ultrametric(pars.a.rooted))

i want to label according to the start column which is my source in my annotation file Any suggestion or help would be really appreciated

My working solution using this , this works as I did the alignment of sequences using full header name which kind of help to look for the pattern and grep it. The data file

I would still like to see the ggtree way instead of putting lots of label in the plot I would only like to see if I can annotate all the tips which have one organism source as one color rather putting all the labels

pars.a
# "Phylogenetic tree with 516 tips and 514 internal nodes"

# by using the class() function
class(pars.a)
# "phylo"

# or by using the str() structure function
str(tree) 
str(pars.a)
# "List of 4"
# this list includes $edge, $Nnode, $ tip.label and $edge.length
# the tree$tip.label includes family designation
#tree$tip.label  # 516 of these
pars.a$tip.label
# from the Science paper, we have seven kinase families:
# kinase categories... TK, TKL, STE, CK1, AGC, CAMK, CMGC
# with the following colours
# "red", "green", "paleblue", "orange", "yellow", "purple", "pink", "green" 


# by using the grep()function on the tree$tip.label part of the object
# we can find the tip labels that include "TK/" - i.e. tyrosine kinases
#grep("TK/", tree$tip.label)  # gives a list of numbers with "TK/" in tip label
#length(grep("TK/", tree$tip.label))

grep("Human",pars.a$tip.label)
length(grep("Human",ml.a.rooted$tip.label)
)
# thus there are 94 tip labels with that are designated TK (not TKL tyrosine kinase like)

# make a vector for each tip.label called tipcol with black on all of these...
#tipcol <- rep('black', length(tree$tip.label))
tipcol <- rep('black',length(pars.a$tip.label))

# make a vector with our list of kinase categories
#kinaseCats <- c("TK/", "TKL", "STE", "CK1", "AGC", "CAMK", "CMGC", "RGC")
kinaseCats <- c("Human","Rat","Murine","Bovine","Canine","Feline","Porcine","Gallus","Bat","Infectious","Transmissible")

# make a vector of color we want:
colorsList <-c("red", "darkolivegreen3", "blue", "orange", "blueviolet", "purple", "pink", "green","brown","navyblue","purple4")

# replace colours where grep gives "TK" as red, etc in a loop
#for(i in 1:length(kinaseCats)){
 # tipcol[grep(kinaseCats[i], tree$tip.label)] <- colorsList[i]
#}


for(i in 1:length(kinaseCats)){
  tipcol[grep(kinaseCats[i], pars.a$tip.label)] <- colorsList[i]
}



# plot with edge length false to see nodes better
#plot(tree, 
 #    use.edge.length = FALSE,
  #   tip.color=tipcol, 
   #  cex = 0.25)


plot(pars.a, 
     use.edge.length = FALSE,
     tip.color=tipcol, 
     cex = 0.5)
# slow to draw due to text - a bit annoying!
nodelabels(cex=0.4)

figure plot

$\endgroup$

1 Answer 1

1
$\begingroup$

there are nice features in R ggtree package. You may use the geom_cladelabel() function for your purposes. But I haven't find yet a way to automatize clade annotation. To do it, you have to be sure of the clade homogeneity before removing tiplabs. here are some commands to do clade label annotation:

annot<-read.delim("annota.txt", header = T)

ggtr<-as.treedata(pars.a)
g<-ggtree(ggtr)
#find the node number corresponding to the common ancestor of the top clade in you tree
MRCA(ggtr, "AAL40400", "AAF25519")
#here, check that all accessions belong to the same organism/group
subset(annot, label=="AAL40400" | label=="AAF25519")
#in this case, sequences are related to bovine.
g<-g + geom_cladelabel(node=89, label="Bovine", align=F, fontsize=4)
plot(g)

#another example
MRCA(ggtr, "AAA46456", "AAF97738")
subset(annot, label=="AAA46456" | label=="AAF97738")
g<-g+ geom_cladelabel(node=73, label="Rodent", align=F, fontsize=4)
plot(g)

build with ggtree

$\endgroup$
7
  • $\begingroup$ looks bit more complicated due to this "have to be sure of the clade homogeneity" if there is a mix up then its difficult $\endgroup$
    – kcm
    May 28, 2020 at 8:28
  • 1
    $\begingroup$ you may keep tip labels where clade are heterogeneous. Or you should consider smaller clades or consider larger taxonomic groups (eg I summarized mice and rats as rodents in the second example). In any case, i think it deserve a visual/manual inspection to ensure the correctness. $\endgroup$ May 28, 2020 at 8:31
  • 1
    $\begingroup$ yes it is fully customizable. Finding the appropriate number is the more fastidious step I think. There is probably a way to automatize it however. $\endgroup$ May 28, 2020 at 8:44
  • 1
    $\begingroup$ nodelabels(cex=0.4) this is how i can label the node in base R which i was doing it is the same for the object which is made using this "g<-ggtree(ggtr)" $\endgroup$
    – kcm
    May 28, 2020 at 8:46
  • 1
    $\begingroup$ this a good thing to identify interesting node numbers! $\endgroup$ May 28, 2020 at 8:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.