No BWA alignment output files!

[thh32]

Hi,

I have been using Bowtie for the last few months for alignment but now have a file for indexing which is 6Gb and Bowtie cannot deal with so I have moved onto using BWA, the indexing went fine and I got my files.

When reading the documentation is seemed as though "bwa mem" was the best option for the actual alignment however all the output is being sent to the screen and now an actual file (determined as iy is going to the output file from the cluster we have and no new files have been created).

Basically was this the correct choice and did I miss something which has caused the output of aligned and unaligned files to not be created?

here is the code I used:

bwa index -p 50Gb_index -a bwtsw contigs.fa


bwa mem -t 20 50Gb_index ../../NZ_bowtie_nonaligned_hungate_thread.done.gz


THANKS!

[Brian Bushnell]

Use the redirect operator ">" to pipe to a file:

bwa mem -t 20 50Gb_index ../../NZ_bowtie_nonaligned_hungate_thread.done.gz > mapped.sam

[thh32]

Originally Posted by Brian Bushnell

Use the redirect operator ">" to pipe to a file:

bwa mem -t 20 50Gb_index ../../NZ_bowtie_nonaligned_hungate_thread.done.gz > mapped.sam

This would work if I wanted a SAM file output but I am aiming for a similar output to that provided with bowtie where you get a file of the reads that align to the indexed file and a seperate file for reads that dont align to the indexed file.

Is there a way to do this with BWA?

[Brian Bushnell]

I don't see an option for it, though you can split the output with samtools into mapped and unmapped reads.

Alternatively, BBMap supports unlimited reference size and also separate files for mapped and unmapped reads, using the "outm=" and "outu=" flags. The output can be in sam or other formats, depending on the file extension.

[dariober]

Originally Posted by thh32

This would work if I wanted a SAM file output but I am aiming for a similar output to that provided with bowtie where you get a file of the reads that align to the indexed file and a seperate file for reads that dont align to the indexed file.

Is there a way to do this with BWA?

Untested, but it should work in bash:

bwa mem myIndex myReads.fq.gz \
| tee >(samtools view -Sh -f4 - > unmapped.sam) \
| tee >(samtools view -Sh -F4 - > mapped.sam) > /dev/null

This sends unmapped reads (-f4) to unmapped.sam and mapped reads (-F4) to mapped.sam. You could customize the samtools commands in the subshells to decide what goes in each file (e.g. by using -q options). In fact, while you are at it, I would convert to bam and sort the mapped reads.

Dario