4

特定の範囲のセットでカバーされていないすべての範囲を定義する最良の方法は何だろうと思っています。たとえば、既知の座標を持つ一連の遺伝子がある場合:

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
"

バイオコンダクターはIRangeここで役に立ちますか? またはこれを解決する他の良い方法はありますか？

4

2 に答える 2

4

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

染色体が 1 から始まり、範囲が閉じているという生体伝導体の規則に注意してください。startとend座標は範囲に含まれます。shiftおよびnarrowonを使用grして、範囲を生体伝導体の慣習と一致させます。GRanges 操作は、数千万の範囲で効率的です。

于 2013-11-28T20:09:39.473 に答える
1

パッケージreduceからご利用いただけます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
于 2013-11-28T20:10:11.593 に答える