2

使用 R 中的系统发育树,我想创建一个矩阵,指示树的每个分支(B1 到 B8)是否与每个物种(A 到 E)相关联,其中 1 表示分支相关联。(如下图)

R 函数 which.edge() 可用于识别物种的终端分支。但它不能识别与每个物种相关的所有分支。对于每个物种,我可以使用什么函数来识别树中从根部到尖端的所有分支?

示例树

library(ape)
ex.tree <- read.tree(text="(A:4,((B:1,C:1):2,(D:2,E:2):1):1);") 
plot(ex.tree)
edgelabels() #shows branches 1-8

这是我想创建的矩阵(物种 AE 作为列,分支 B1-B8 作为行),但功能简单,而不是手动创建。

B1 <- c(1,0,0,0,0)
B2 <- c(0,1,1,1,1)
B3 <- c(0,1,1,0,0)
B4 <- c(0,1,0,0,0)
B5 <- c(0,0,1,0,0)
B6 <- c(0,0,0,1,1)
B7 <- c(0,0,0,1,0)
B8 <- c(0,0,0,0,1)
Mat <- rbind(B1,B2,B3,B4,B5,B6,B7,B8)   
colnames(Mat) <- c("A","B","C","D","E")
Mat

例如,分支 B2 属于物种 BE,但不属于物种 A。对于物种 E,存在分支 B2、B6、B8。

哪个 R 函数最好?提前致谢!

4

1 回答 1

2

我不知道有任何内置函数可以做到这一点。我写了一个辅助函数,它可以从存储在tree对象中的边缘数据中计算出来。

branchNodeAdjacency <- function(x) {
    m <- matrix(0, ncol=nt, nrow=nrow(x$edge))
    from <- x$edge[,1]
    to <- x$edge[,2]
    g <- seq_along(x$tip.label)
    while (any(!is.na(g))) {
        i <- match(g, to)
        m[cbind(i, seq_along(i))] <- 1
        g <- from[i]
    }
    rownames(m) <- paste0("B", seq.int(nrow(m)))
    colnames(m) <- x$tip.label
    m
}

branchNodeAdjacency(ex.tree)
#    A B C D E
# B1 1 0 0 0 0
# B2 0 1 1 1 1
# B3 0 1 1 0 0
# B4 0 1 0 0 0
# B5 0 0 1 0 0
# B6 0 0 0 1 1
# B7 0 0 0 1 0
# B8 0 0 0 0 1

这个想法是我们跟踪每个内部节点表示哪些叶节点值。

于 2015-11-09T23:42:40.203 回答