This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Identify deletion in multiple bam files compared to the same bacterial reference

Hi, I have done an illumina WGS experiment in many gut samples and mapped those samples back to a reference bacterial genome.

I have identified a small deletion (2bp) in one bam file compared to the reference bacterial genome. I have the coordinates of the region. Is there a quick way to chech if the same deletion is present in all other samples (bam files) compared to the reference??

Thanks

bedtools bam
• 1,427 views • link
written by David David

You can have the hits on this region :

samtools view -b your_file.bam "chr_name:start-end" > hits_file.bam

Then check out by eyes or process the CIGAR

• 1 views • link
written by Bastien Hervé Bastien Hervé

1 answer

Thanks for the answer, i have tried as follows

samtools view -b your_file.bam "sequence1:345000-345510" > hits_file.bam

Then to count how many reads have at least the deletion i have done:

samtools view  17.mapped_sorted.bam  "sequence1:345000-345510" | grep  -cE '[0-9]{1,2}D' 
202

To count all i did:

samtools view  17.mapped_sorted.bam  "sequence1:345000-345510" | wc -l 
410

The idea is to know how many deletions( i got in this region , so here is 202/410*100= 51%

Would that make sense or is there any easy way to do that ?

• 1 views • link
updated by Ram Ram • written by David David

You posted this as an answer. Did you mean to do that or did you mean to add it as a comment to Bastien's comment? Your're pulling information from the CIGAR string with the grep statement, right?

• 1 views • link
written by Kevin Blighe Kevin Blighe

Log in to answer this question.

More posts like this