时间戳

首页 > 解决方案 > 如何根据条件查找两个数据帧之间的重叠区域

r - 如何根据条件查找两个数据帧之间的重叠区域

问题描述

我有两个数据框，一个叫做strain_1，另一个叫做strain_2. 每个数据框有 4 列（st_A, ed_A, st_B, ed_B: 表示“开始”和“结束”位置），但行数不同。st_A,ed_A和st_B,分别是block_A和block_Bed_B的“开始”和“结束”位置（参见图 1和下面的示例）。

我正在寻找和 之间的常见重叠块。strain_1strain_2

以图 1 为例：

strain_1 <- data.frame(st_A=c(7,25,35,48,89), ed_A=c(9,28,38,51,91),
                       st_B=c(123,97,140,73, 13), ed_B=c(127,98,145,76,16))

strain_2 <- data.frame(st_A=c(5,20,36,49) ,  ed_A=c(8,25,39,50),    
                       st_B=c(124,95,141,105) ,  ed_B=c(129,100,147,110))

从这个例子中，我们看到三个重叠区域（图 1）：

重叠区域定义为：block_A 和 block_B 的最小值( or )和( or st_A )st_B最大值（见图2：绿色框 = 公共区域）。ed_Aed_B

data frame目标是用这些公共区域（一对块）创建一个新的

## result_desired 
result_desired <- data.frame(st_A=c(5,20,35), ed_A=c(9,28,39),
                             st_B=c(123,95,140), ed_B=c(129,100,147))

根据每个块的大小，有 16 种可能的组合（见图3 ）。

有没有快速的方法来做到这一点？知道我有几千行的数据。

我尝试了一些代码，基于@Gregor 评论，但我无法得到想要的结果：

## require(dplyr)
require(dplyr)

## data 
strain_1 <- data.frame(st_A=c(7,25,35,48,89), ed_A=c(9,28,38,51,91),
                       st_B=c(123,97,140,73, 13), ed_B=c(127,98,145,76,16))

strain_2 <- data.frame(st_A=c(5,20,36,49) ,  ed_A=c(8,25,39,50),    
                       st_B=c(124,95,141,105) ,  ed_B=c(129,100,147,110))

# merge data to get cross join 
cj_data <-merge(strain_1,strain_2, by = NULL)

# Check block1 and block2
cj_filtered <- cj_data %>% mutate(c_block1= case_when(st_A.x <= st_A.y & ed_A.x <= ed_A.y | 
                                                      st_A.x >= st_A.y & ed_A.x >= ed_A.y |
                                                      st_A.x <= st_A.y & ed_A.x >= ed_A.y | 
                                                      st_A.x >= st_A.y & ed_A.x <= ed_A.y ~ "overlap_OK",
                                                      TRUE ~ "NO"),

                                  c_block2= case_when(st_B.x <= st_B.y & ed_B.x <= ed_B.y | 
                                                      st_B.x >= st_B.y & ed_B.x >= ed_B.y |
                                                      st_B.x <= st_B.y & ed_B.x >= ed_B.y | 
                                                      st_B.x >= st_B.y & ed_B.x <= ed_B.y ~ "overlap_OK",
                                                      TRUE ~ "NO"))

## cj_filtered:
st_A.x  ed_A.x  st_B.x  ed_B.x  st_A.y  ed_A.y  st_B.y  ed_B.y  c_block1    c_block2
7       9       123     127     5       8       124     129     overlap_OK  overlap_OK
25      28      97      98      5       8       124     129     overlap_OK  overlap_OK
35      38      140     145     5       8       124     129     overlap_OK  overlap_OK
48      51      73      76      5       8       124     129     overlap_OK  overlap_OK
89      91      13      16      5       8       124     129     overlap_OK  overlap_OK
7       9       123     127     20      25      95      100     overlap_OK  overlap_OK
25      28      97      98      20      25      95      100     overlap_OK  overlap_OK
35      38      140     145     20      25      95      100     overlap_OK  overlap_OK
48      51      73      76      20      25      95      100     overlap_OK  overlap_OK
89      91      13      16      20      25      95      100     overlap_OK  overlap_OK
7       9       123     127     36      39      141     147     overlap_OK  overlap_OK
25      28      97      98      36      39      141     147     overlap_OK  overlap_OK
35      38      140     145     36      39      141     147     overlap_OK  overlap_OK
48      51      73      76      36      39      141     147     overlap_OK  overlap_OK
89      91      13      16      36      39      141     147     overlap_OK  overlap_OK
7       9       123     127     49      50      105     110     overlap_OK  overlap_OK
25      28      97      98      49      50      105     110     overlap_OK  overlap_OK
35      38      140     145     49      50      105     110     overlap_OK  overlap_OK
48      51      73      76      49      50      105     110     overlap_OK  overlap_OK
89      91      13      16      49      50      105     110     overlap_OK  overlap_OK

谢谢你的帮助。

标签： rdplyr

解决方案

以下是使用 2 个选项data.table：

1a) 使用非 equi 连接：

cols <- c(paste0("x.", names(strain_1)), paste0("i.", names(strain_2)))
DT <- rbindlist(list(
    strain_1[strain_2, on=.(st_A>=st_A, st_A<=ed_A), nomatch=0L, mget(cols)],
    strain_1[strain_2, on=.(st_A<=st_A, ed_A>=st_A), nomatch=0L, mget(cols)]
))

1b）使用foverlaps：

setkey(strain_1, st_A, ed_A)
setkey(strain_2, st_A, ed_A)
foverlaps(strain_1, strain_2, nomatch=0L)

然后是另一个步骤 2 以获得所需的输出：

DT[between(x.st_B, i.st_B, i.ed_B) | between(i.st_B, x.st_B, x.ed_B), 
    .(st_A=pmin(x.st_A, i.st_A), 
        ed_A=pmax(x.ed_A, i.ed_A), 
        st_B=pmin(x.st_B, i.st_B), 
        ed_B=pmax(x.ed_B, i.ed_B))]

输出：

   st_A ed_A st_B ed_B
1:    5    9  123  129
2:   20   28   95  100
3:   35   39  140  147

数据：

library(data.table)
strain_1 <- data.frame(st_A=c(7,25,35,48,89), ed_A=c(9,28,38,51,91),
    st_B=c(123,97,140,73, 13), ed_B=c(127,98,145,76,16))

strain_2 <- data.frame(st_A=c(5,20,36,49) ,  ed_A=c(8,25,39,50),    
    st_B=c(124,95,141,105) ,  ed_B=c(129,100,147,110))

result_desired <- data.frame(st_A=c(5,20,35), ed_A=c(9,28,39),
    st_B=c(123,95,140), ed_B=c(129,100,147))

setDT(strain_1)
setDT(strain_2)
setDT(result_desired)

ps：Bioconducter中应该也有IRanges的东西。

推荐阅读