SEQanswers

Go Back   SEQanswers > Applications Forums > Genomic Resequencing



Similar Threads
Thread Thread Starter Forum Replies Last Post
SAM/BAM format to wiggle format pinki999 Bioinformatics 19 08-12-2015 12:35 AM
Can GSNAP generate SAM/BAM files? efoss Bioinformatics 4 10-16-2011 08:11 PM
Generate consensus sequence from BAM alignment Adam Witney General 2 09-14-2011 02:49 AM
the tophat generate the bam file instead of sam files? dingkai0564 Bioinformatics 1 11-10-2010 07:33 PM
which softwate can generate fastq format file? biocc Illumina/Solexa 1 08-21-2008 06:02 AM

Reply
 
Thread Tools
Old 06-18-2010, 02:56 AM   #1
ForeignMan
Member
 
Location: Germany

Join Date: Jun 2010
Posts: 20
Default Generate subset from BAM format

Greetings everyone,

I have a quite large BAM-File after an alignment with bwa. My question is pretty general:
Is there a tool to extract some reads (by read name or row index) from a BAM-file?

I already did this once on the fastq files, but it would be really helpful to work directly on the BAM-file so I don't have to do all the preceding steps (alignment, conversion from SAM to BAM, sorting, removing of duplicate reads etc.) again.
I could not find any tool yet, that is able to read and write BAM-file. Does anyone have experience with Rsamtools, a Bioconductor package to read BAM-files into R? I couldn't find any way to export the data from R back into a BAM-file. If Rsamtools could write BAM-files it would be perfect. I'm looking for something similar.
I would be very great if someone knows such a tool or can give any advice on how to get a subset of my aligned reads without having to repeat the whole alignment procedure and post-processing.

Thanks in advance.
ForeignMan is offline   Reply With Quote
Old 06-18-2010, 03:07 AM   #2
ffinkernagel
Senior Member
 
Location: Marburg, Germany

Join Date: Oct 2009
Posts: 110
Default

If you want scripting, PySam works well (enough) from python and reads/writes both sam and bam files, and you can filter on whatever you want.

Otherwise 'samtools' with samtools view my_file | grep readname should do for a quick job.

Rsamtools reads bam, but on large files (60M aligned reads) I had it die with 'no negative indices allowed' or such.
ffinkernagel is offline   Reply With Quote
Old 06-20-2010, 10:01 AM   #3
ForeignMan
Member
 
Location: Germany

Join Date: Jun 2010
Posts: 20
Default

Thank you very much for your answer!

Didn't know PySam yet, but I might give it a try, although I got no experience with python.
Your tip concerning samtools view with grep will work only for sam output, right? That would be a compromise, but working just with BAM-files would be the best, because I wouldn't need to convert SAM to BAM again. I've got Illumina paired-end sequencing data, and calling bwa sampe also takes quite some time. My problem is, that want to select a subset several times (maybe 10 times or more). Thus, I try to save as much calculations as possible. And in general, I'm interested in tools to directly manipulate BAM-files.
But your tips sound good, and I guess I will work on SAM-files this time if I can't find another way.

Concerning Rsamtools: I remember I had the same error once, too. But, if I remember it right, it occurred only when reading specific (or all) columns. I was able to read a >1GB BAM-file when I set the parameter to only read qname, flag, or position for example. But that's a serious problem and wouldn't help me here, when I want to read my whole BAM-file (its size is >2GB). Thanks for reminding me of that.
ForeignMan is offline   Reply With Quote
Old 06-20-2010, 10:35 PM   #4
ffinkernagel
Senior Member
 
Location: Marburg, Germany

Join Date: Oct 2009
Posts: 110
Default

samtools: They read (and write) BAM as well.

RSamtools: Filtering the columns was not enough - my BAM is ~2.4 GB.
ffinkernagel is offline   Reply With Quote
Old 02-28-2011, 06:20 AM   #5
Joker!sAce
Member
 
Location: Denmark

Join Date: Feb 2011
Posts: 21
Default

Hi,

I'm working on something similar, almost identical. You could try this script to get all the aligned reads. Pysam can read/ write to SAM files.
Here are two (Documentation 1, Documentation 2 ) links that can help you.

Quote:
import pysam
samfile = pysam.Samfile( "NA06984.454.MOSAIK.SRP000033.2009_11.bam", "rb" )

for alignedread in samfile.fetch():
print alignedread
This is a very simple script to get all the aligned reads.
"NA06984.454.MOSAIK.SRP000033.2009_11.bam" was the name of my BAM file.

Also, Could I have a look at what your data looks like, one or two records would do. I have a free week, I could write a small script that would do the trick.

Best,
Joker!sAce

Last edited by Joker!sAce; 02-28-2011 at 06:43 AM.
Joker!sAce is offline   Reply With Quote
Reply

Tags
bam, tools and techniques

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 08:57 AM.


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