1

この問題を解決するために findOverlap を使用しようとしましたが、条件なしでオーバーラップ領域のみを見つけたので、データを選択する条件があれば。どうすればいいですか?

以下のような2つのデータフレームがあるとしましょう

データフレーム

Sample, start, stop, event, probe, length, length/probe, region
CNV1234,  2000,  3000,  CN gain,  23,  235, 9, intron
CNV1534,  1200,  1800,  CN loss,  60,  600  10, exon

データフレームb

Sample, start, stop, event, probe, length, length/probe, region
CNV234,  2500,  3500,  CN gain,  23,  235, 9, exon 
CNV34,  1200,  1800,  CN loss,  60,  600  10, intron

2つの質問があります

まず、CNVの長さが50%を超えており、このオーバーラップがイントロン領域にあるこの2つのデータフレーム間のオーバーラップを見つけたい

次に、オーバーラップ領域の長さを知りたいです。

結果を次のようなデータフレームにしたい

Sample, start, stop, event, probe, length, length/probe, region, overlap, length of overlap
CNV1234,  2000,  3000,  CN gain,  23,  235, 9, intron, CNV234, 500
4

1 に答える 1

1

これがあなたのデータです

a <- read.csv(textConnection(
    "Sample, start, stop, event, probe, length, length/probe, region
     CNV1234,  2000,  3000,  CN gain,  23,  235, 9, intron
     CNV1534,  1200,  1800,  CN loss,  60,  600  10, exon"))

b <- read.csv(textConnection(
    "Sample, start, stop, event, probe, length, length/probe, region
     CNV234,  2500,  3500,  CN gain,  23,  235, 9, exon 
     CNV34,  1200,  1800,  CN loss,  60,  600  10, intron"))

GenomicRangesパッケージをロードします (データは実際には複数の染色体から取得されていると想定しており、すべての染色体でこれを実行したいと考えています。「A」は染色体名です)。

library(GenomicRanges)
gr1 <- with(a, GRanges("A", IRanges(start, stop - 1L),
                       Sample=Sample, event=event))
gr2 <- with(b, GRanges("A", IRanges(start, stop - 1L, names=Sample),
                       Region=Sample))

GRanges が範囲 (開始座標と終了座標を含む 1 ベース) を表す方法に注意してください。これらのオブジェクト間のすべてのオーバーラップを見つけます (最小幅の 1/2 より短いオーバーラップなど、一部のオーバーラップを除外するには min.overlaps を使用できます)。

h <- findOverlaps(gr1, gr2)

幅の「50%」が何であるかは明確ではありません -- 幅は? b? -- そこで、すべてのオーバーラップの幅を計算し、幅が「a」の幅の 1/2 より大きいものを保持します。

wd <- width(pintersect(gr1[queryHits(h)], gr2[subjectHits(h)]))
ok <- wd > width(gr1[queryHits(h)]) / 2
h <- h[ok]

最後に、オーバーラップ基準を満たすクエリを選択し、メタデータ列 ( mcols) とオーバーラップする領域のオーバーラップ幅を追加して、結果を組み立てます。

result <- gr1[queryHits(h)]
mcols(result) <- cbind(mcols(result), mcols(gr2[subjectHits(h)]))
result$`width of overlap` <- wd[ok]

結果は強制的にデータ フレームに戻される可能性がありas.data.frame(result)ます。または、ダウンストリーム分析は GRanges インフラストラクチャを使用して自然に行われる可能性があります。

Bioconductorメーリング リストでBioconductorパッケージについて質問することをお勧めします(購読は必要ありません)。これを行うためのより効率的な方法がある可能性が高く、そのメーリング リストのメンバーがこれらの解決策を提供します。

于 2013-08-13T06:12:29.057 に答える