改变plot.dendrogram中的叶子颜色,就像使用package ape的plot.phylo一样

2019年8月4日 425次阅读

我试图用与使用包’ape’绘制树的相同’风格’来绘制凝聚聚类(UPGMA与Agnes)的结果.我在下图中包含了一个简单的例子

关键问题是我希望能够根据叶子标签中的图案为树形图的叶子着色.我尝试了两种方法：要么使用hc2Newick,要么使用Joris Meys的代码作为对Change Dendrogram leaves的回答.两者都没有给出令人满意的结果.可能是我还没有完全理解树形图的构建方式.可以在https://www.dropbox.com/s/gke9qnvwptltkky/abundance.agnes.ave找到abundance.agnes.ave对象的ASCII保存(从运行的agnes中存储).

当我使用第一个选项(使用来自bioconductor的ctc包的hc2Newick)时,我在使用此代码时会得到下图：

write(hc2Newick(as.hclust(abundance.agnes.ave)),file="all_samples_euclidean.tre")
eucltree<-read.tree(file="all_samples_euclidean.tre")
eucltree.laz<-ladderize(eucltree,FALSE)
tiplabs<-eucltree$tip.label
numbertiplabs<-length(tiplabs)
colourtips<-rep("green",numbertiplabs)
colourtips[grep("II",tiplabs)]<-"red"
plot(eucltree.laz,tip.color=colourtips,adj=1,cex=0.6,use.edge.length=F)
add.scale.bar()

这显然不太理想,情节的“对齐”并不像我想的那样.我认为这与分支长度计算有关,但我没有最模糊的想法如何解决这个问题.当然,与colLab函数的结果相比,它看起来更像我想要报告的树形图样式.此外,在上面的代码中使用use.edge.length = T确实给出了一个未正确“对齐”的聚类：

第二种方法使用Joris Meys的colLab函数和下面的代码给出了下图

clusDendro<-as.dendrogram(as.hclust(abundance.agnes.ave))
labelColors<-c("red","green")
clusMember<-rep(1,length(rownames(abundance.x)))
clusMember[grep("II",rownames(abundance.x))]<-2
names(clusMember)<-rownames(abundance.x)

colLab <- function(n)
{
  if(is.leaf(n)) {
    a <- attributes(n)
    # clusMember - a vector designating leaf grouping
    # labelColors - a vector of colors for the above grouping
    labCol <- labelColors[clusMember[which(names(clusMember) == a$label)]]
    attr(n, "nodePar") <- c(a$nodePar, lab.col = labCol)
  }
  n
}

clusDendro<-dendrapply(clusDendro, colLab)
plot(clusDendro,horiz=T,axes=F)

这个情节越来越接近我想要的,但是我不知道为什么空心圆出现在叶子上以及如何去掉它们.

任何帮助深表感谢.

亲切的问候,

调频

最佳答案 我不久前编写了这段代码,看起来机制中有一些变化.

我使用的plot.dendrogram函数有一个参数nodePar.自从我上次使用该函数以来,行为已经改变,虽然这通常用于内部节点,但它显然也对外部节点有影响.根据帮助文件,pch的默认值现在是1：2.

因此,您需要在添加到colLab函数中的外部节点的属性中专门指定pch = NA.尝试适应它：

colLab <- function(n)
{
  if(is.leaf(n)) {
    a <- attributes(n)
    # clusMember - a vector designating leaf grouping
    # labelColors - a vector of colors for the above grouping
    labCol <- labelColors[clusMember[which(names(clusMember) == a$label)]]

    attr(n, "nodePar") <- 
        if(is.list(a$nodePar)) c(a$nodePar, lab.col = labCol,pch=NA) else
                               list(lab.col = labCol,pch=NA)
  }
  n
}

在我的机器上,这解决了问题.

或者,您可以查看ape包中函数plot.phylo的参数use.edge.length.你将它设置为FALSE,但从你的解释我相信你希望它设置为默认值,为TRUE.

编辑：为了使函数更通用,将labelColors和clusMember作为参数添加到函数中可能是个好主意.我的快速解决方案不是干净代码的最好例子……

还要忘记我所说的使用边长的内容. ape包将其解释为真实的树形图,并将use.edge.length置于TRUE将边长变换为演化时间.因此,树状图的“怪异”概述.

另请注意,如果treeleaf没有nodePar属性,使用c()函数添加额外参数将导致不良影响：如果添加例如lab.cex = 0.6,c()函数将创建一个向量而不是列表,并在参数列表中有字符值时将lab.cex的值转换为字符.在这种情况下,这将是颜色的名称,这解释了您在评论中谈到的错误.

点赞