SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
How to use samtools idxstats <bamfile> serhatg Bioinformatics 2 07-02-2014 10:28 PM
read counts from samtools idxstats witty Bioinformatics 3 01-30-2014 10:35 AM
Why does samtools idxstats not show the stats kgulukota Bioinformatics 5 04-30-2012 12:18 AM
samtools idxstats Irina Pulyakhina Bioinformatics 4 10-19-2011 05:05 AM
idxstats option for samtools GenoMax Bioinformatics 0 08-10-2011 06:42 AM

Reply
 
Thread Tools
Old 08-19-2014, 10:05 AM   #1
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default Issue with new samtools 1.0 index and idxstats

When using samtools 1.19, I would index my .bam, and use idxstatst to make a file of read counts per sequence in my reference, and the unmapped reads would be at the bottom of the list.

The line would look like

* 0 0 1027820

To show 1027820 reads did not align

When using samtools 1.0 to index, now, the idxstats file always claims there is exactly 1 read that mapped to nothing, which is clearly wrong.

Can anyone else replicate this? This can't be the desired behavior, can it?
swbarnes2 is offline   Reply With Quote
Old 08-19-2014, 01:43 PM   #2
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

I just tested on an old BAM file (reindexed and ran idxstats) and got 0, as appropriate. Perhaps this only happens if there are unmapped reads in the file? I should note that you should ensure that both htslib and samtools are on the 1.0 tag release. When I just recloned them (just to ensure that any local changes weren't in there) I couldn't even compile samtools since the most current commit breaks compilation with the version of htslib I had from earlier in the day.
dpryan is offline   Reply With Quote
Old 08-19-2014, 02:15 PM   #3
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

I went into the htslib-1.0 folder of my samtools 1.0 folder, and "make" seems to finish without errors, but samtools 1.0 still does the same thing: idxstats run from the v 1.0 indexing shows only 1 unmapped read. Rebuilding the index with 0.1.19 and rerunning idxstats shows the correct number.

Last edited by swbarnes2; 08-19-2014 at 02:32 PM.
swbarnes2 is offline   Reply With Quote
Old 08-20-2014, 05:33 AM   #4
dariober
Senior Member
 
Location: Cambridge, UK

Join Date: May 2010
Posts: 311
Default

Quote:
Originally Posted by swbarnes2 View Post
Can anyone else replicate this? This can't be the desired behavior, can it?
I just upgraded to samtools-1.0.

Indeed the same happens to me:

Code:
samtools --version
samtools 1.0
Using htslib 1.0
Copyright (C) 2014 Genome Research Ltd.

samtools index fk042_F5_14_CHEM1.bam
samtools idxstats fk042_F5_14_CHEM1.bam | tail -n1
*	0	0	1

## Now with samtools 0.18.0:
~/applications/samtools/samtools-0.1.18/samtools index fk042_F5_14_CHEM1.bam
samtools idxstats fk042_F5_14_CHEM1.bam | tail -n1
*	0	0	8882
dariober is offline   Reply With Quote
Old 08-20-2014, 05:43 AM   #5
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

@dariober: In fact this is the issue I alluded to over on biostars :P

It turns out that this only happens if there are unmapped read in the BAM file (otherwise, it correctly reports 0).
dpryan is offline   Reply With Quote
Old 08-20-2014, 05:50 AM   #6
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

I've submitted a bug report.
dpryan is offline   Reply With Quote
Old 08-20-2014, 09:44 AM   #7
swbarnes2
Senior Member
 
Location: San Diego

Join Date: May 2008
Posts: 912
Default

Quote:
Originally Posted by dpryan View Post
I've submitted a bug report.
Thank you. I didn't want to do that until I had some confirmation that other people were seeing the problem.
swbarnes2 is offline   Reply With Quote
Old 08-26-2014, 04:20 AM   #8
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

I created a pull request to fix this.
dpryan is offline   Reply With Quote
Old 08-26-2014, 04:38 AM   #9
dariober
Senior Member
 
Location: Cambridge, UK

Join Date: May 2010
Posts: 311
Default

Quote:
Originally Posted by dpryan View Post
I created a pull request to fix this.
Thanks a lot for taking care of this issue. Here's my comment to your post on github:

Quote:
I'm quite keen on querying the index to get the number of unmapped reads, so I favour the second solution. (With the disclaimer the I'm not the one going to fix it though!)
dariober is offline   Reply With Quote
Reply

Tags
samtools-1.0

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 03:29 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO