This is a test version of Biostars. For the public version, visit https://www.biostars.org.
coverageBed misses some reads with the overlap of one basepair

My region of interest is chr4:88487220-88488363 I run coverageBed with default parameters

coverageBed reports only 3 reads that fall into this interval.

I checked with samtools the reads that fall into this regions and got 7 reads back:

It seems that bedtools did not take into account reads that start at position 88487171. Reads that start at 88487171, are ending at position 88487220, meaning they have an overlap with the region of interest by only 1 bp. So, they should be detected by coverageBed. Why does it happen so?

HWI-ST227:399:C4Y0RACXX:4:1301:3554:28973    16    chr4    88487171    255    50M    *    0    0    GGGACCGTGGGACTGTCCACTATGTTTGAGCCCAAGGTAGACTGGTGCTG    HGIGJIJIIJIHIJJJIHFHJJJJIIJJJJJJJJJJJHHHHHFFFFFCCC    XA:i:1    MD:Z:22C27    NM:i:1
HWI-ST227:399:C4Y0RACXX:4:1308:16026:23049    16    chr4    88487171    255    50M    *    0    0    GGGACCGTGGGACTGTCCACTATGTTTGAGCCCAAGGTAGACTGGTGCTG    IIHFIIJJIJIHIJJIJIIIJJJJJJJJJJJJIJIIIGHHHHFFFFFCCC    XA:i:1    MD:Z:22C27    NM:i:1
HWI-ST227:399:C4Y0RACXX:4:1310:18272:80730    16    chr4    88487171    255    50M    *    0    0    GGGACCGTGGGACTGTCCACTATGTTTGAGCCCAAGGTAGACTGGTGCTG    HFB@GIHF?B:9HGDBGCFEIIIIHBIIGIFEGEADDFDFBBDBADB8@@    XA:i:1    MD:Z:22C27    NM:i:1
HWI-ST227:399:C4Y0RACXX:4:2111:15178:17225    16    chr4    88487171    255    50M    *    0    0    GGGACCGTGGGACTGTCCACTATGTTTGAGCCCAAGGTAGACTGGTGCTG    JJJJIJJJJJIGIJJJJIHJJJJJJJJJJJJJJJJJJHHHHHFFFFFBCB    XA:i:1    MD:Z:22C27    NM:i:1
HWI-ST227:399:C4Y0RACXX:4:1114:17538:80401    16    chr4    88487196    255    50M    *    0    0    TTGAGCCCAAGGTAGACTGGTGCTGACCAGCAAAGCTTGCCGGGCAGGTT    GHIHFIHHGDGCIGJIJIJHGHGIIIJIIIHJJJIHEHHHHHFFFFFCCC    XA:i:0    MD:Z:50    NM:i:0
HWI-ST227:399:C4Y0RACXX:4:2106:6167:96305    16    chr4    88487196    255    50M    *    0    0    TTGAGCCCAAGGTAGACTGGTGCTGACCAGCAAAGCTTGCCGGGCAGGTT    EBCFADDDGBGHBB?:@HG?@HBD?9F??AHEEHGCDFC@FA????D???    XA:i:0    MD:Z:50    NM:i:0
HWI-ST227:399:C4Y0RACXX:4:2111:4985:84167    16    chr4    88487196    255    50M    *    0    0    TTGAGCCCAAGGTAGACTGGTGCTGACCAGCAAAGCTTGCCGGGCAGGTT    IJIJIHBGIIJJJIJJIGIJIJJIIGIJJIIJIGBGEHFHGHFFFFD@?@    XA:i:0    MD:Z:50    NM:i:0
alignment

1 answer

I think the XA:i: tag is telling you that those reads are not mapped uniquely, so coverageBed probably ignores them :)

Log in to answer this question.