RDl*_*ady 6 plot r phylogeny ape-phylo
我为一个蛋白质家族构建了一个系统发育树,可以将其分成不同的组,按照受体类型或反应类型对每个蛋白质进行分类.树中的节点被标记为受体类型.
在系统发育树中,我可以看到属于相同组或受体类型的蛋白质在同一分支中聚集在一起.所以我想折叠这些具有共同标签的分支,按给定的关键字列表对它们进行分组.
命令将是这样的:
./collapse_tree_by_label -f phylogenetic_tree.newick -l list_of_labels_to_collapse.txt -o collapsed_tree.eps(或pdf)
我的list_of_labels_to_collapse.txt将是这样的:A B C D
我的纽约树将是这样的:(A_1:0.05,A_2:0.03,A_3:0.2,A_4:0.1):0.9,(((B_1:0.05,B_2:0.02,B_3:0.04):0.6,(C_1:0.6) ,C_2:0.08):0.7):0.5,(D_1:0.3,D_2:0.4,D_3:0.5,D_4:0.7,D_5:0.4):1.2)
没有折叠的输出图像是这样的:http: //i.stack.imgur.com/pHkoQ.png
输出图像折叠应该是这样的(collapsed_tree.eps):http: //i.stack.imgur.com/TLXd0.png
三角形的宽度应表示分支长度,三角形的高度必须表示分支中的节点数.
我一直在玩R.中的"猿"包.我能够绘制一个系统发育树,但我仍然无法弄清楚如何通过标签中的关键字来折叠分支:
require("ape")
Run Code Online (Sandbox Code Playgroud)
这将加载树:
cat("((A_1:0.05,A_2:0.03,A_3:0.2,A_4:0.1):0.9,(((B_1:0.05,B_2:0.02,B_3:0.04):0.6,(C_1:0.6,C_2:0.08):0.7):0.5,(D_1:0.3,D_2:0.4,D_3:0.5,D_4:0.7,D_5:0.4):1.2):0.5);", file = "ex.tre", sep = "\n")
tree.test <- read.tree("ex.tre")
Run Code Online (Sandbox Code Playgroud)
这应该是崩溃的代码
这将绘制树:
plot(tree.test)
Run Code Online (Sandbox Code Playgroud)
您存储在 R 中的树已经将提示存储为多分法。只需用代表多分法的三角形绘制树即可。
据我所知,没有任何功能ape可以做到这一点,但是如果你稍微搞乱绘图功能,你就可以完成它
# Step 1: make edges for descendent nodes invisible in plot:
groups <- c("A", "B", "C", "D")
group_edges <- numeric(0)
for(group in groups){
group_edges <- c(group_edges,getMRCA(tree.test,tree.test$tip.label[grepl(group, tree.test$tip.label)]))
}
edge.width <- rep(1, nrow(tree.test$edge))
edge.width[tree.test$edge[,1] %in% group_edges ] <- 0
# Step 2: plot the tree with the hidden edges
plot(tree.test, show.tip.label = F, edge.width = edge.width)
# Step 3: add triangles
add_polytomy_triangle <- function(phy, group){
root <- length(phy$tip.label)+1
group_node_labels <- phy$tip.label[grepl(group, phy$tip.label)]
group_nodes <- which(phy$tip.label %in% group_node_labels)
group_mrca <- getMRCA(phy,group_nodes)
tip_coord1 <- c(dist.nodes(phy)[root, group_nodes[1]], group_nodes[1])
tip_coord2 <- c(dist.nodes(phy)[root, group_nodes[1]], group_nodes[length(group_nodes)])
node_coord <- c(dist.nodes(phy)[root, group_mrca], mean(c(tip_coord1[2], tip_coord2[2])))
xcoords <- c(tip_coord1[1], tip_coord2[1], node_coord[1])
ycoords <- c(tip_coord1[2], tip_coord2[2], node_coord[2])
polygon(xcoords, ycoords)
}
Run Code Online (Sandbox Code Playgroud)
然后你只需要循环遍历组来添加三角形
for(group in groups){
add_polytomy_triangle(tree.test, group)
}
Run Code Online (Sandbox Code Playgroud)
| 归档时间: |
|
| 查看次数: |
3319 次 |
| 最近记录: |