Using samtools to analyze bwa output

[niteesh.prasad]

I have a query file with 70,000 lines of sequences and when I do the

bwasw

each sam output file is about 35 mb and obviously contains all the sequence lines. I need to perform the alignment on about 1000 files, so you can imagine that I don't want to have 35 gigs worth of output to sort.

I am trying to get rid of the non-aligned reads from the SAM output file ans I searched through forums and tried the following code with no avail.

./samtools view -bS -f 4 testing.sam > output.bam

All my sam files give the following error:

[samopen] SAM header is present: 1 sequences.
[sam_read1] reference '0' is recognized as '*'.
Parse warning at line 3: mapped sequence without CIGAR
Parse error at line 3: missing colon in auxiliary data
Aborted (core dumped)

Here's what a part of my sam files looks like

@SQ	SN:gi|89106884|ref|AC_000091.1|	LN:4646332

	4	*	0	0	*	*	0	0	<CGTAAAAAGGA> 	*
	4	*	0	0	*	*	0	0	<TAGCCGTAGGG>   *

I'd appreciate help with this!

[nilshomer]

There's no read name in the SAM records you list.

[niteesh.prasad]

But this is the output I get from bwa. how should i fix this?

[nilshomer]

Post the BWA help mailing list. It looks like a bug to me. Make sure you include a test case.

[niteesh.prasad]

Before I do that, I want to make sure I am not wrong. Here is the code I used for the bwa:

./bwa bwasw INDEXES/AC_000091.fna QUERY/635plus104_500bp_reads.fasta > results.sam

The input AC_000091.fna was the indexes created with the following command:

./bwa index AC_000091.fna "AC_000091.fna";

and the query file was a FASTA sequence file.

So, in theory this should have worked?