behaviour of regioneR permTest function
1
0
Entering edit mode
shoko • 0
@shoko-21942
Last seen 5.2 years ago

Hello,

I would be grateful if you could help me understand the functions better.

I created bed file A and B and genome, all of which have 3 peaks in exactly the same region. They look like this. Chr1 1 10
Chr2 100 150
Chr3 2 3
A and B are the ChIP peaks between which I would like to compare the similarity, and genome is the known TF map-able region.

The number of overlaps are 3 as expected, but when I run AvsB <- permTest(A=A, B=B, ntimes=10, randomize.function=randomizeRegions, evaluate.function=numOverlaps, count.once=TRUE, genome=genome) I get 0s and 1s and 2s for AvsB[["numOverlaps"]][["permuted"]]. I thought since reads in A could only map onto "genome" region and that is where peaks from B are, it would always be 3.

Do you have any ideas why this could be? Thank you so much!

regioneR • 819 views
ADD COMMENT
0
Entering edit mode
bernatgel ▴ 150
@bernatgel-7226
Last seen 13 hours ago
Spain

Hi @shoko

I think the numbers you are seeing are the expected result (mostly). Randomize regions will randomly place the regions into the available genome, possibly overlapping (see ?randomizeRegions). This means that the 3 regions might end up in chr2, and in fact this is what happens most of the time if you test it with

A <- toGRanges(c("chr1:1-10", "chr2:100-150", "chr3:2-3"))
randomizeRegions(A=A, genome=A)

then, since you have set "count.once=TRUE", only one overlap per region in B will be reported, so most of the time it will be a 1.

This will explain the 1's, 2's and 3's, but it does not explain the 0's. In this case, this is due to a limitation of randomizeRegions, which can only work with genomes whose chromosomes start at 1. Therefore, in this case some of the random regions might appear in the chr2:1-99 region, which in theory is out of the genome. Take into account that this situation is only found in strange cases like yours (custom toy mini-genome) and some advanced usages but will not affect most use cases.

In any case, we've added a note in randomizeRegions documentation stating that and we'll try to fix it in next release of regioneR.

Hope this helps

ADD COMMENT

Login before adding your answer.

Traffic: 698 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