optimization - R 优化 : How can I avoid a for loop in this situation?

Question

我正在尝试在 R 中做一个简单的基因组轨道交叉，并遇到主要的性能问题，可能与我使用 for 循环有关。

在这种情况下，我以 100bp 的间隔预定义了窗口，并且我正在尝试计算 mylist 中的注释覆盖每个窗口的多少。从图形上看，它看起来像这样:

          0    100   200    300    400   500   600  
windows: |-----|-----|-----|-----|-----|-----|

mylist:    |-|   |-----------|

所以我写了一些代码来做到这一点，但它相当慢并且已经成为我代码中的瓶颈:

##window for each 100-bp segment    
windows <- numeric(6)

##second track
mylist = vector("list")
mylist[[1]] = c(1,20)
mylist[[2]] = c(120,320)


##do the intersection
for(i in 1:length(mylist)){
  st <- floor(mylist[[i]][1]/100)+1
  sp <- floor(mylist[[i]][2]/100)+1
  for(j in st:sp){       
    b <- max((j-1)*100, mylist[[i]][1])
    e <- min(j*100, mylist[[i]][2])
    windows[j] <- windows[j] + e - b + 1
  }
}

print(windows)
[1]  20  81 101  21   0   0

自然，这用于比我在此处提供的示例大得多的数据集。通过一些分析，我可以看到瓶颈在 for 循环中，但是我使用 *apply 函数对其进行矢量化的笨拙尝试导致代码运行速度慢了一个数量级。

我想我可以用 C 写一些东西，但如果可能的话，我想避免这种情况。任何人都可以提出另一种可以加快计算速度的方法吗？

Answer 1

“正确”的做法是使用生物导体IRanges包，它使用 IntervalTree 数据结构来表示这些范围。

拥有您自己的两个对象 IRanges对象，然后您将使用 findOverlaps功能取胜。

在这里获取:

http://www.bioconductor.org/packages/release/bioc/html/IRanges.html

顺便说一下，包的内部是用 C 编写的，所以它的速度非常快。

编辑

再想一想，它不像我所建议的那样扣篮(单线)，但是如果您完全使用基因组间隔(或其他类型)，您绝对应该开始使用这个库......你可能需要做一些设置操作和东西。抱歉，虽然没有时间提供确切的答案。

我只是认为向您指出这个库很重要。