r - 计算数据串中的某些化学 react

Question

最有效的计数方法是什么？一串分子式中的甲基化(CH2 的差异)。考虑一串 100 种不同的分子式，我有一些分子式，例如 C6H14O3、C5H12O3 和 C4H10O3。第一个和第二个因 CH2 不同，第二个和第三个因 CH2 不同，因此该字符串中有 2 个甲基化。当字符串非常大时，这会变得更加复杂。所以我想计算有多少个分子式仅因字符串中的一个 CH2 而不同。

考虑我有:

DT<- data.frame(formula=c("C6H12O7S1","C6H10O8S1","C7H4O2N4","C8H12O5S1","C8H16O5S1","C8H12O3N2","C8H14O4S1","C9H7O3N1S1","C9H11O6N1S1","C9H9O6N1S1","C9H12O5S1","C9H18O5S1","C9H14O5","C9H20O5S1","C9H9O4N1S1","C9H9O5N1S1","C9H14O6","C10H11O5N1S1","C10H14O6","C10H16O6S1","C10H17O5N1","C10H20O7S1","C10H14O4","C10H12O7N2","C10H16O6","C10H14O6N2","C10H7O4N1S1","C10H18O6S1","C10H16O5","C10H13O6N1S1","C10H18O7S1","C11H18O6S1","C11H15O6N1","C11H22O7S1","C11H16O6S1","C11H16O6","C11H18O6"))

我想计算这个字符串中有多少个如上所述的 CH2 差异。

有人现在有一种简单的方法可以做到这一点，特别是对于非常大的字符串吗？

非常感谢。

Answer 1

此解决方案使用 thelatemail 的建议作为分割字符串的起点:

重现原始数据集(将 stringsAsFactors 设置为 F)

DT <- data.frame(formula=c("C6H12O7S1","C6H10O8S1","C7H4O2N4","C8H12O5S1","C8H16O5S1",
                          "C8H12O3N2","C8H14O4S1","C9H7O3N1S1","C9H11O6N1S1","C9H9O6N1S1",
                          "C9H12O5S1","C9H18O5S1","C9H14O5","C9H20O5S1","C9H9O4N1S1",
                          "C9H9O5N1S1","C9H14O6","C10H11O5N1S1","C10H14O6","C10H16O6S1",
                          "C10H17O5N1","C10H20O7S1","C10H14O4","C10H12O7N2","C10H16O6",
                          "C10H14O6N2","C10H7O4N1S1","C10H18O6S1","C10H16O5","C10H13O6N1S1",
                          "C10H18O7S1","C11H18O6S1","C11H15O6N1","C11H22O7S1","C11H16O6S1","C11H16O6",
                          "C11H18O6"), stringsAsFactors = F) %>% arrange(formula)

> head(DT)
       formula
1 C10H11O5N1S1
2   C10H12O7N2
3 C10H13O6N1S1
4     C10H14O4
5     C10H14O6
6   C10H14O6N2

展开列，使每个元素用一列来计数(缺少的按 0 计)

DT2 <- strsplit(DT$formula, "(?<=[0-9])(?=[A-Z])|(?<=[A-Z])(?=[0-9])", perl=TRUE) %>%
  lapply(function(x){structure(list(formula = rep(paste0(x, collapse = ""), length(x)/2),
                                    element = x[seq(from = 1, to = length(x), by = 2)],
                                    count = as.integer(x[seq(from = 2, to = length(x), by = 2)])),
                               .Names = c("formula", "element", "count"),
                               row.names = c(NA, -length(x)/2), 
                               class = "data.frame")}) %>%
  data.table::rbindlist() %>% 
  spread(element, count, fill = 0)

>head(DT2)
        formula  C  H N O S
1: C10H11O5N1S1 10 11 1 5 1
2:   C10H12O7N2 10 12 2 7 0
3: C10H13O6N1S1 10 13 1 6 1
4:     C10H14O4 10 14 0 4 0
5:     C10H14O6 10 14 0 6 0
6:   C10H14O6N2 10 14 2 6 0

获取所有可能的成对组合的列表并扩展数据集

pairwise.combos <- combn(nrow(DT2), m = 2)

DT3 <- rbind(DT2[pairwise.combos[1,],],
             DT2[pairwise.combos[2,],])
DT3$pair <- rep(seq.int(dim(pairwise.combos)[2]),2)

> head(DT3)
        formula  C  H N O S pair
1: C10H11O5N1S1 10 11 1 5 1    1
2: C10H11O5N1S1 10 11 1 5 1    2
3: C10H11O5N1S1 10 11 1 5 1    3
4: C10H11O5N1S1 10 11 1 5 1    4
5: C10H11O5N1S1 10 11 1 5 1    5
6: C10H11O5N1S1 10 11 1 5 1    6

根据 OP 的说明编辑最后一步，即一对公式应仅在 CH2 上有所不同，并且其他元素的计数应相同。

检查每对CH2是否不同(此步骤可以修改以检查其他差异)

DT4 <- DT3 %>% group_by(pair) %>%
  arrange(C, H) %>%
  summarise(CH2.diff = (diff(C) == 1) && (diff(H) == 2) && 
              (diff(N) == 0) & (diff(O) == 0) & (diff(S) == 0)) %>%
  ungroup() %>%
  filter(CH2.diff == 1) %>% select(pair)
DT4 <- right_join(DT3, DT4)

# total count of CH2 in pairwise comparisons
> length(unique(DT4$pair))
[1] 9

# check which pairs differ by CH2
> head(DT4)
       formula  C  H N O S pair
1 C10H11O5N1S1 10 11 1 5 1   35
2   C9H9O5N1S1  9  9 1 5 1   35
3 C10H13O6N1S1 10 13 1 6 1   96
4  C9H11O6N1S1  9 11 1 6 1   96
5     C10H14O6 10 14 0 6 0  149
6     C11H16O6 11 16 0 6 0  149