Efficiently finding a relative position inside a GRanges object
1
0
Entering edit mode
@aliaksei-holik-4992
Last seen 8.1 years ago
Spain/Barcelona/Centre for Genomic Regu…

Dear BioConductors,

I have two GRanges objects, one of which contains several sequences (let's say motifs) and the other is a single bp (let's say it's a variant):

library(GenomicRanges)
variant.position <- 63552595
variant.GR <- GRanges("chr17", IRanges(variant.position, width=1), "+")
motifs.GR <- GRanges("chr17", IRanges(c(63552587, 63552589), width = 10), c("+"), motif=c("MYC", "TCF4"))

I'd like to get the relative position of the variant inside the motifs that it matches. I know I can use findOverlaps and then subtract the start position of the motif from the variant position:

hits <- findOverlaps(variant.GR, motifs.GR)
relative.position <- variant.position - start(motifs.GR[subjectHits(hits)]) + 1
names(relative.position) <- motifs.GR[subjectHits(hits)]$motif

But I'm aware that this approach is using three steps, i.e. subsetting the GRanges object, extracting the start position and subtracting it from the variant position. Besides, I will have to reverse the calculation for the motifs on the negative strand. Given that I'm likely to run it on millions of variants and hundreds of thousands of motifs, I'm concerned about the computational overheads and wondering if there's a more efficient way to get that relative position. Ideally, the function would take the above objects and return the index of the motif containing the variant and the position of the variant inside the motif, ideally, in a strand-specific manner.

Many thanks!

 

genomicranges • 1.1k views
ADD COMMENT
2
Entering edit mode
@michael-lawrence-3846
Last seen 2.3 years ago
United States

Probably could just use mapToTranscripts(). The motifs are not transcripts, obviously, but it will give you the answer you want.

ADD COMMENT
0
Entering edit mode

It seems to be just what I needed. Thanks a lot!

ADD REPLY

Login before adding your answer.

Traffic: 760 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6