r - 如何找到到最近的非重叠元素的距离？

Question

我有一个如下表所示的表，其中每个簇(第 1 列)包含小区域中不同元素(第 4 列)的注释，并具有起始坐标(第 2 列)和结束坐标(第 3 列)。对于每个条目，我想添加一列，该列对应于到该簇中最近的其他元素的距离。但我想排除簇中一对元素具有相同的开始/结束坐标或重叠区域的情况。如何为此类数据框生成额外的 nearest_distance 列？

cluster-47593-walk-0125 252     306     AR    
cluster-47593-walk-0125 6       23      ZNF148
cluster-47593-walk-0125 357     381     CEBPA 
cluster-47593-walk-0125 263     276     CEBPB 
cluster-47593-walk-0125 246     324     NR3C1 
cluster-47593-walk-0125 139     170     HMGA1 
cluster-47593-walk-0125 139     170     HMGA2 
cluster-47593-walk-0125 207     227     IRF8  
cluster-47593-walk-0125 207     227     IRF1  
cluster-47593-walk-0125 207     245     IRF2  
cluster-47593-walk-0125 207     227     IRF3  
cluster-47593-walk-0125 207     227     IRF4  
cluster-47593-walk-0125 207     227     IRF5  
cluster-47593-walk-0125 207     227     IRF6  
cluster-47593-walk-0125 204     245     IRF7  
cluster-47593-walk-0125 13      36      PATZ1 
cluster-47593-walk-0125 14      143     PAX4  
cluster-47593-walk-0125 4       25      RREB1 
cluster-47593-walk-0125 73      87      SMAD1 
cluster-47593-walk-0125 73      87      SMAD2 
cluster-47593-walk-0125 73      87      SMAD3 
cluster-47593-walk-0125 71      89      SMAD4 
cluster-47593-walk-0125 11      40      SP1   
cluster-47593-walk-0125 11      38      SP2   
cluster-47593-walk-0125 7       38      SP3   
cluster-47593-walk-0125 11      38      SP4   
cluster-47593-walk-0125 13      33      GTF2I 
cluster-47593-walk-0125 281     352     YY1   
cluster-47586-walk-0222 252     306     AR    
cluster-47586-walk-0222 6       23      ZNF148
[...]

Answer 1

首先，一些列名称

names(data) <- c("cluster", "start", "end", "element")
data
                   cluster start end element
1  cluster-47593-walk-0125   252 306      AR
2  cluster-47593-walk-0125     6  23  ZNF148
3  cluster-47593-walk-0125   357 381   CEBPA
4  cluster-47593-walk-0125   263 276   CEBPB

现在创建新列

data$nearest_distance <- apply(data, 1, function(x) 
 {
     cluster <- x[1]
     start <- as.numeric(x[2])
     end <- as.numeric(x[3])
     elem <- x[4]
     posb <- data[data$cluster == cluster & data$element != elem & 
                  ((data$start > end) | (data$end < start)), ]
     startDist <- as.matrix(dist(c(end, posb$start)))[, 1]
     endDist <- as.matrix(dist(c(start, posb$end)))[, 1]
     best.dist <- min(startDist[startDist > 0], endDist[endDist > 0])
     return(best.dist)
  }
)

我真的不喜欢至少这个函数的开头，但我无法想出更好的解决方案..所以我们有

                   cluster start end element nearest_distance
1  cluster-47593-walk-0125   252 306      AR                7
2  cluster-47593-walk-0125     6  23  ZNF148               48
3  cluster-47593-walk-0125   357 381   CEBPA                5
4  cluster-47593-walk-0125   263 276   CEBPB                5
5  cluster-47593-walk-0125   246 324   NR3C1                1 
.....

编辑:修复system.time()测试后发现这是一种非常低效的方法。显然，计算整个 dist() 矩阵是多余的，因此我们可以将这两行更改为

startDist <- abs(end-posb$start)
endDist <- abs(start-posb$end)

另一个小变化是我们可以删除约束 data$element != elem 因为后面有 > 0。在 1000 个簇(每个簇有 30 行)上测试这个函数花了三分钟多的时间。仍然存在子集问题，所以我尝试将数据拆分成一个列表，这允许我们使用矩阵而不是数据框(因为簇的约束消失了) ，这也提高了效率。这次我们有 10 000 个簇，每个簇 30 行

data <- data[rep(1:30, each = 10000), ]
data$cluster <- factor(rep(1:10000, 30))

spl <- split(data[, c(2:3)], data$cluster)
spl <- lapply(spl, data.matrix)

system.time({
x = lapply(spl, function(z) {
     apply(z, 1, function(x) {
       start <- x[1]
       end <- x[2]
       posb <- z[z[,1] > end | z[,2] < start, , drop = FALSE]
       startDist <- abs(end-posb[, 1])
       endDist <- abs(start-posb[, 2])
       best.dist <- min(startDist[startDist > 0], endDist[endDist > 0])
       return(best.dist)
     })
  })
})
data$nearest_distance = unsplit(x, data$cluster)


user  system elapsed 
18.16    0.00   18.35