r-查找定义范围之外的所有范围

Question

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

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在这里有用吗？或者还有其他的好方法可以解决这个问题？

Answer 1

使用Bioconductor软件包GenomicRanges，将原始数据转换为GRanges

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

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

gaps <- gaps(gr)

GRanges知道strand。您没有在GRanges构造函数中指定一个strand，所以strand被指定为*。因此，在+、-和*股上存在"缺口"，您只对*股上的缺口感兴趣

> 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操作在数以百万计的10秒内是有效的。

Answer 2

您可以从IRanges包中使用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