查找定义的范围集之外的所有范围答案

【问题标题】：Find all ranges outside of defined set of ranges查找定义的范围集之外的所有范围
【发布时间】：2013-11-28 19:34:58
【问题描述】：

我想知道定义给定范围集未涵盖的所有范围的最佳方法是什么。例如，如果我有一组已知坐标的基因：

dtGenes <- fread(
  "id,start,end
 1,1000,1300
 2,1200,1500
 3,1600,2600
 4,3000,4000
")

假设我知道染色体的总长度（为简单起见，假设它们都在同一条染色体上）是 10000。所以，最后我希望有以下基因间区域列表：

"startR,endR
    0,1000
 1500,1600
 2600,3000
 4000,10000
"

Bioconductor 的IRange 在这里有用吗？还是有其他好办法解决这个问题？

【问题讨论】：

标签： r bioconductor iranges

【解决方案1】：

使用 Bioconductor 包 GenomicRanges，将您的原始数据转换为 GRanges

library(GenomicRanges)
gr <- with(dtGenes, GRanges("chr1", IRanges(start, end, names=id),
                            seqlengths=c(chr1=10000)))

然后找出你的基因之间的差距

gaps <- gaps(gr)

GRanges 知道链。您没有在GRanges 构造函数中指定链，因此链被分配了*。因此 +、- 和 * 链上存在“间隙”，您只对 * 链上的那些感兴趣

> gaps[strand(gaps) == "*"]
GRanges with 4 ranges and 0 metadata columns:
      seqnames        ranges strand
         <Rle>     <IRanges>  <Rle>
  [1]     chr1 [   1,   999]      *
  [2]     chr1 [1501,  1599]      *
  [3]     chr1 [2601,  2999]      *
  [4]     chr1 [4001, 10000]      *
  ---
  seqlengths:
    chr1
   10000

注意 Bioconductor 约定，染色体从 1 开始，并且范围是封闭的——start 和 end 坐标包含在范围内。在gr 上使用shift 和narrow 使您的范围与Bioconductor 约定一致。 GRanges 操作在数以百万计的范围内非常有效。

【讨论】：

【解决方案2】：

您可以使用IRanges 包中的reduce

reduce 首先将 x 中的范围从左到右排序，然后合并重叠或相邻的。

library(IRanges)
dat <- read.csv(text="id,start,end
 1,1000,1300
 2,1200,1500
 3,1600,2600
 4,3000,4000
")

ir <- IRanges(dat$start,dat$end)
rir <- reduce(ir)
IRanges of length 3
    start  end width
[1]  1000 1500   501
[2]  1600 2600  1001
[3]  3000 4000  1001

【讨论】：

感谢指点！ reduce 并没有完全解决寻找外部区域的任务，但它做了重要的第一步；在其他一些情况下它肯定对我有用，但是对于这个给定的任务，我非常喜欢gaps 的解决方案，它通过单个操作完成所有工作。
@VasilyA 我不使用很多 iranges 包。因此，MartinMorgan 解决方案当然是这里的必经之路！