-
Notifications
You must be signed in to change notification settings - Fork 178
reproduce phyloseq example using ggtree
Guangchuang Yu edited this page Aug 31, 2015
·
5 revisions
This example is demonstrated in phyloseq
vignette.
require(phyloseq)
data(GlobalPatterns)
GP <- prune_taxa(taxa_sums(GlobalPatterns) > 0, GlobalPatterns)
GP.chl <- subset_taxa(GP, Phylum=="Chlamydiae")
plot_tree(GP.chl, color="SampleType", shape="Family", label.tips="Genus", size="Abundance")
require(ggplot2)
require(ggtree)
tree <- GP.chl@phy_tree
p <- ggtree(tree, ladderize = F)
p <- p + geom_text(subset=.(!isTip), aes(label=label), hjust=-.2, size=4)
dd <- psmelt(GP.chl)
dd <- dd[dd$Abundance > 0, ]
data <- merge(p$data, dd, by.x="label", by.y="OTU")
spacing <- 0.02
idx <- with(data, sapply(table(node)[unique(node)], function(i) 1:i)) %>% unlist
hjust <- spacing * idx * max(data$x)
data$xdodge <- data$x + hjust
p <- p + geom_point(data=data, aes(x=xdodge, color=SampleType, shape=Family, size=Abundance), na.rm=T) +
theme(legend.position="right") + scale_size_continuous(trans=log_trans(5))
d2 <- unique(data[, c("x", "y", "Genus")])
p + geom_text(data=d2, aes(label=Genus), hjust=-.3, na.rm=T, size=4)
Although Abundance
points were adjusted to prevent overlap, it's hard to detect pattern in this figure.
We believe using gheatmap
is more easy to interpret the data.
.
We can also use dots to encode more information with aligning them vertically by SampleType
. This is a better way and more easy to interpret.