使用 data.table 识别与第二个表的相关距离

Question

我有四个来自系统发育的分类单元。我已将系统发育划分为“世代”，以识别特定切点（进化枝）的相关分类群。下面的数据表中可以看到这样的示例：

library(data.table)
clades = data.table(generation = rep(c("g0", "g4"), each = 4),
                    taxa = rep(c(LETTERS[1:4]), times = 2),
                    clade = c(1:4, 1, 1, 2, 2)
                    )

在时间 0，每个人彼此无关，在时间 4，A 和 B 一起处于一个进化枝，C 和 D 一起处于一个进化枝。

我还知道每个类群之间的地理距离，以长格式存储在这里：

distance = data.table(
    row = c("A", "A", "B", "A", "B", "C"),
    col = c("B", "C", "C", "D", "D", "D"), 
    value = c(1.41, 5.65, 5.83, 5.83, 5.65, 1.41)
  )

我想有效地计算，在每一代中，在每个进化枝内，一定距离内的类群数量（在本例中为 2），以及总共有多少个比较 - 给出两个类群在一定距离内的概率彼此之间有一定的距离。

distance$cut = distance$value < 2

在此示例中，在第 0 代中，应该有四个进化枝，每个进化枝在截止范围内没有对，也没有对（因为所有进化枝都包含单个分类群）。在第四个时间点，我们应该有两个进化枝，并且每个进化枝内应该有一对在距离内，并且一对进行比较。如下表：

result = data.table(generation = c("g0", "g0", "g0", "g0", "g4", "g4"),
                    clade = c(1, 2, 3, 4, 1, 2),
                    numerator = c(0, 0, 0, 0, 1, 1),
                    denominator = c(0, 0, 0, 0, 1, 1))

我已经设法使用 *apply 函数来做到这一点，但它非常慢，我需要大规模地做到这一点。如果有人可以建议一种快速方法，我将不胜感激。

generations = unique(clades$generation)
output = lapply(generations, function(g) {
  generation_df = clades[clades$generation == g]
  clade_sets <- unique(generation_df$clade)
  
  clade_result <- sapply(clade_sets, function(cs) {
    clade_taxa = generation_df$taxa[generation_df$clade == cs]
    trait = distance$trait[distance$row %in% clade_taxa & distance$col %in% clade_taxa]
    
    if (length(trait) >= 1) {
      numerator = sum(trait)
      denominator = length(trait)
    } else {
      numerator = 0
      denominator = 0
    }
    list(numerator = numerator, denominator = denominator)
  })
  
  # Combine clade result with generation and clade information
  data.frame(
    generation = g,
    clade = clade_sets,
    numerator = unlist(clade_result[1, ]),
    denominator = unlist(clade_result[2, ]),
    stringsAsFactors = FALSE
  )
})
dplyr::bind_rows(output)

Answer 1

我认为这得益于使用

distance

作为快速查找的矩阵：

distance_m <- dcast(row ~ col, value.var = "value", data = distance)
distance_m <- as.matrix(distance_m[,-1]) |>
  `dimnames<-`(list(distance_m[[1]], colnames(distance_m)[-1])) 
distance_m[lower.tri(distance_m)] <- distance_m[upper.tri(distance_m)]
distance_m
#      B    C    D
# A 1.41 5.65 5.83
# B 5.65 5.83 5.65
# C 5.83 5.65 1.41

从这里开始，

out <- clades[, if (.N == 1L) list(num=0, den=0L) else {
  m <- combn(taxa, 2)
  list(num=apply(m, 2, function(z) distance_m[z[1], z[2]]), den=ncol(m))
}, by = .(generation, clade)]
out
#    generation clade   num   den
#        <char> <num> <num> <int>
# 1:         g0     1  0.00     0
# 2:         g0     2  0.00     0
# 3:         g0     3  0.00     0
# 4:         g0     4  0.00     0
# 5:         g4     1  1.41     1
# 6:         g4     2  1.41     1

从这里可以轻松地将

num

减少为“计数低于 2”：

out[den > 0, num := +(num < 2)]
#    generation clade   num   den
#        <char> <num> <num> <int>
# 1:         g0     1     0     0
# 2:         g0     2     0     0
# 3:         g0     3     0     0
# 4:         g0     4     0     0
# 5:         g4     1     1     1
# 6:         g4     2     1     1