intersect result is wrong compare to bedtools intersect

[StayHungryStayFool]

I compare R function intersect in GenomicRanges package(terms RI), linux bedtools intersect(LI) and your tools, find that RI is the same as LI, but some wrong in your tools. I can not send you an example pic in here.

[StayHungryStayFool]

[StayHungryStayFool]

I want to get the intersection region of A, B, C. Orange track from your tools. can you explain what happened?

[asntech]

When you intersect with intervene the order of your files does matter and you will end up having slightly different results with a different order. Also, make sure that for Intervene we are using pybedtoost with u=True or v=True. I hope this helps,

[canulef3]

Hi! I'm having the same issue and i would like to know how to add the option u=True or v=True in my command line.

[amizeranschi]

Hi @asntech

Regarding your earlier comment that When you intersect with intervene the order of your files does matter and you will end up having slightly different results with a different order.

Why is this so? Set intersection should be a commutative operation, even when evaluating the overlaps between sets of genomic regions.

Based on this older remark (https://github.com/daler/pybedtools/issues/45#issuecomment-2543863), it would make sense to run bedtools with u=False, because reporting the actual overlapping regions correctly is more informative than simply counting the number of overlaps between regions. Similarly, I don't see the logic for setting v=True.

The way I see it, genomic regions are also sets themselves, with the unit element being the base pair. Thus, intersecting (sets of) genomic regions should accurately output the base-pairs common to those regions, even if the resulting number of overlapping elements won't match with the numbers of regions from the original BED files.

[Rohit-Satyam]

Hi!!

I find this tool fairly helpful. I tested intervene on my files and it works fine for me. I actually manually calculated the overlapping fraction of regions between two files and intervene calculated the correct percentage. When finding intersecting regions what important is to know which file will act as your reference. Above the diagonal of the matrix file that intervene produce are the scenarios where the files representing the column names are taken as reference (i.e. -a $column_header_files).

Below diagonals are the scenarios where the row names represent the files taken as a reference (-a). It is fair enough to calculate what fraction of regions of the reference file overlaps with the other bed file. Using -u option will prevent the overcounting of features present in the reference files that overlap with other bed files. For two files at a time intervene performs two intersections. As indicated by (daler/pybedtools#45 (comment)) Not using u=True results in another issue -- the total number of features for overlapping with b is greater than the number of features in a in the first place:

If the amount of overlap is your concern you can use 50% of overlap as threshold and -r option of bedtools intersect which does not count the overlaps discussed by (daler/pybedtools#45 (comment)) by using --bedtools-options f=0.50,r function

[amizeranschi]

@Rohit-Satyam

Thanks for your thoughts. Keep in mind, I'm not interested in counting the overlapping regions as separate entities. Instead, I'm interested in getting the total length (in bp) of the overlapping regions (1 bp minimum overlap) from two or more BED files. From this point of view, the intersection operation should be commutative, i.e. the order of input BED files shouldn't matter (see the example below).

I've also compared Intervene's results with those of two other tools and the results are different (even though the other tools "agree" with each other). Have a look: https://github.com/asntech/intervene/issues/34.

I would like to somehow get Intervene to give me the same results as bedops and bedr as shown in the previous link. According to https://github.com/daler/pybedtools/issues/45#issuecomment-2543863, running bedtools within Intervene with u=False should give me the right results. Using the example from the previous link:

a.bed        -------------        --------------------------
b.bed            -----------               -----     ------
Overlaps:        ---------                 -----     ------

However, I don't see a way of doing this with Intervene right now.

Out of curiosity, have you used bedops before? Its --intersect operation is much simpler (IMO) than that of bedtools. It also, supposedly, runs quite faster, when the input files are sorted. https://bedops.readthedocs.io/en/latest/content/reference/set-operations/bedops.html#intersect-i-intersect.

It seems that Intervene's idea of intersecting BED files is more akin to the --element-of operation from bedops: https://bedops.readthedocs.io/en/latest/content/reference/set-operations/bedops.html#element-of-e-element-of.

[asntech]

Thanks @Rohit-Satyam for your interest and input.

@amizeranschi Intervene is not designed to give you a list of overlapping regions, rather it spits out how many of the individual regions in a.bed overlap with b.bed and vise vera. To get a list of overlapping regions you need to use bedtools, bedops, or bedr.

For the u=True; v=True I just posted here https://github.com/asntech/intervene/issues/35#issuecomment-658593662

You can also provide additional bedtools arguments using --bedtools-options

[amizeranschi]

@asntech I was hoping to use Intervene for my use case, due to its user-friendly way of creating UpSet diagrams. While bedr does what I need when setting feature = "bp", it only outputs Venn diagrams for up to 5 BED files.

[Rohit-Satyam]

Hi @asntech

I have some more doubts regarding the fractions of overlap reported by intervene. I will try to talk this with an example:

I had two-bed files adipose-tissue.csv and blood.csv

and the following is the number of features (coordinates) in them

adipose-tissue.csv 978 blood.csv 3326 I tried doing what intervene would do given two files bedtools intersect -a blood.csv -b adipose-tissue.csv -f 0.50 -u | wc -l

Gives me 354 overlapping features. If I had to calculate the fraction of features of adipose-tissue.csv overlapping with the blood.csv, I would calculate 354/978 which will give me 0.3619. However, Intervene is calculating 354/3326 which I think is quite misleading. Can you help me understand why Intervene does that?

[whiteorchid]

Thanks so much for the great tool of intervene!

I have some doubts and wish can have your guidance!

It seems the --bedtools-options r ("bedtools intersect -r ") is default from the venn plot result, which in the Venn plot may has less overlap, while 90% of bed1 file may has overlap with 30% of bed2, in the -r mode, it may only display 10% overlap in the Venn plot.

So which is better to display in the Venn plot, with or without the -r mode?

Thanks so much for your guidance! @asntech Appreciate!